X-Git-Url: https://pintos-os.org/cgi-bin/gitweb.cgi?a=blobdiff_plain;f=src%2Flanguage%2Fstats%2Fcrosstabs.q;h=b71524a11a41175d132ea56978d6a1e4d7ff9400;hb=c6fe58a22249f4f486b42f35fd8bd537c91e8e6e;hp=e57d9343c724ec96b2cb5badf2ec5eaa1fed48c1;hpb=43b1296aafe7582e7dbe6c2b6a8b478d7d9b0fcf;p=pspp-builds.git diff --git a/src/language/stats/crosstabs.q b/src/language/stats/crosstabs.q index e57d9343..b71524a1 100644 --- a/src/language/stats/crosstabs.q +++ b/src/language/stats/crosstabs.q @@ -1,5 +1,5 @@ /* PSPP - a program for statistical analysis. - Copyright (C) 1997-9, 2000, 2006 Free Software Foundation, Inc. + Copyright (C) 1997-9, 2000, 2006, 2009 Free Software Foundation, Inc. This program is free software: you can redistribute it and/or modify it under the terms of the GNU General Public License as published by @@ -46,13 +46,10 @@ #include #include #include -#include #include #include #include #include -#include -#include #include #include #include @@ -61,6 +58,8 @@ #include #include "minmax.h" +#include "xalloc.h" +#include "xmalloca.h" #include "gettext.h" #define _(msgid) gettext (msgid) @@ -178,9 +177,10 @@ static struct pool *pl_col; /* For column data. */ static int internal_cmd_crosstabs (struct lexer *lexer, struct dataset *ds); static void precalc (struct casereader *, const struct dataset *); -static void calc_general (struct ccase *, const struct dataset *); -static void calc_integer (struct ccase *, const struct dataset *); -static void postcalc (void); +static void calc_general (const struct ccase *, const struct dataset *); +static void calc_integer (const struct ccase *, const struct dataset *); +static void postcalc (const struct dataset *); + static void submit (struct tab_table *); static void format_short (char *s, const struct fmt_spec *fp, @@ -191,11 +191,16 @@ int cmd_crosstabs (struct lexer *lexer, struct dataset *ds) { int result = internal_cmd_crosstabs (lexer, ds); + int i; free (variables); pool_destroy (pl_tc); pool_destroy (pl_col); + for (i = 0; i < nxtab; i++) + free (xtab[i]); + free (xtab); + return result; } @@ -301,20 +306,20 @@ internal_cmd_crosstabs (struct lexer *lexer, struct dataset *ds) grouper = casegrouper_create_splits (input, dataset_dict (ds)); while (casegrouper_get_next_group (grouper, &group)) { - struct ccase c; + struct ccase *c; precalc (group, ds); - for (; casereader_read (group, &c); case_destroy (&c)) + for (; (c = casereader_read (group)) != NULL; case_unref (c)) { if (mode == GENERAL) - calc_general (&c, ds); + calc_general (c, ds); else - calc_integer (&c, ds); + calc_integer (c, ds); } casereader_destroy (group); - postcalc (); + postcalc (ds); } ok = casegrouper_destroy (grouper); ok = proc_commit (ds) && ok; @@ -356,7 +361,7 @@ crs_custom_tables (struct lexer *lexer, struct dataset *ds, struct cmd_crosstabs goto done; if (xalloc_oversized (nx, by_nvar[n_by])) { - msg (SE, _("Too many crosstabulation variables or dimensions.")); + msg (SE, _("Too many cross-tabulation variables or dimensions.")); goto done; } nx *= by_nvar[n_by]; @@ -514,12 +519,14 @@ static unsigned hash_table_entry (const void *, const void *); static void precalc (struct casereader *input, const struct dataset *ds) { - struct ccase c; + struct ccase *c; - if (!casereader_peek (input, 0, &c)) - return; - output_split_file_values (ds, &c); - case_destroy (&c); + c = casereader_peek (input, 0); + if (c != NULL) + { + output_split_file_values (ds, c); + case_unref (c); + } if (mode == GENERAL) { @@ -547,7 +554,7 @@ precalc (struct casereader *input, const struct dataset *ds) sorted_tab = xnrealloc (sorted_tab, n_sorted_tab + count, sizeof *sorted_tab); - v = local_alloc (sizeof *v * x->nvar); + v = xmalloca (sizeof *v * x->nvar); for (j = 2; j < x->nvar; j++) v[j] = get_var_range (x->vars[j])->min; for (j = 0; j < count; j++) @@ -581,7 +588,7 @@ precalc (struct casereader *input, const struct dataset *ds) break; } } - local_free (v); + freea (v); } sorted_tab = xnrealloc (sorted_tab, @@ -593,7 +600,7 @@ precalc (struct casereader *input, const struct dataset *ds) /* Form crosstabulations for general mode. */ static void -calc_general (struct ccase *c, const struct dataset *ds) +calc_general (const struct ccase *c, const struct dataset *ds) { /* Missing values to exclude. */ enum mv_class exclude = (cmd.miss == CRS_TABLE ? MV_ANY @@ -611,7 +618,7 @@ calc_general (struct ccase *c, const struct dataset *ds) struct crosstab *x = xtab[t]; const size_t entry_size = (sizeof (struct table_entry) + sizeof (union value) * (x->nvar - 1)); - struct table_entry *te = local_alloc (entry_size); + struct table_entry *te = xmalloca (entry_size); /* Construct table entry for the current record and table. */ te->table = t; @@ -632,12 +639,14 @@ calc_general (struct ccase *c, const struct dataset *ds) te->values[j].f = case_num (c, x->vars[j]); else { - memcpy (te->values[j].s, case_str (c, x->vars[j]), - var_get_width (x->vars[j])); + size_t n = var_get_width (x->vars[j]); + if (n > MAX_SHORT_STRING) + n = MAX_SHORT_STRING; + memcpy (te->values[j].s, case_str (c, x->vars[j]), n); /* Necessary in order to simplify comparisons. */ memset (&te->values[j].s[var_get_width (x->vars[j])], 0, - sizeof (union value) - var_get_width (x->vars[j])); + sizeof (union value) - n); } } } @@ -660,12 +669,12 @@ calc_general (struct ccase *c, const struct dataset *ds) } next_crosstab: - local_free (te); + freea (te); } } static void -calc_integer (struct ccase *c, const struct dataset *ds) +calc_integer (const struct ccase *c, const struct dataset *ds) { bool bad_warn = true; @@ -768,7 +777,7 @@ hash_table_entry (const void *a_, const void *aux UNUSED) hash = a->table; for (i = 0; i < xtab[a->table]->nvar; i++) - hash ^= hsh_hash_bytes (&a->values[i], sizeof a->values[i]); + hash = hash_bytes (&a->values[i], sizeof a->values[i], hash); return hash; } @@ -781,12 +790,13 @@ static void enum_var_values (struct table_entry **entries, int entry_cnt, int var_idx, union value **values, int *value_cnt); static void output_pivot_table (struct table_entry **, struct table_entry **, + const struct dictionary *, double **, double **, double **, int *, int *, int *); -static void make_summary_table (void); +static void make_summary_table (const struct dictionary *); static void -postcalc (void) +postcalc (const struct dataset *ds) { if (mode == GENERAL) { @@ -794,7 +804,7 @@ postcalc (void) sorted_tab = (struct table_entry **) hsh_sort (gen_tab); } - make_summary_table (); + make_summary_table (dataset_dict (ds)); /* Identify all the individual crosstabulation tables, and deal with them. */ @@ -811,7 +821,8 @@ postcalc (void) if (pe == NULL) break; - output_pivot_table (pb, pe, &mat, &row_tot, &col_tot, + output_pivot_table (pb, pe, dataset_dict (ds), + &mat, &row_tot, &col_tot, &maxrows, &maxcols, &maxcells); pb = pe; @@ -822,13 +833,25 @@ postcalc (void) } hsh_destroy (gen_tab); + if (mode == INTEGER) + { + int i; + for (i = 0; i < n_sorted_tab; i++) + { + free (sorted_tab[i]->u.data); + free (sorted_tab[i]); + } + free (sorted_tab); + } } -static void insert_summary (struct tab_table *, int tab_index, double valid); +static void insert_summary (struct tab_table *, int tab_index, + const struct dictionary *, + double valid); /* Output a table summarizing the cases processed. */ static void -make_summary_table (void) +make_summary_table (const struct dictionary *dict) { struct tab_table *summary; @@ -867,7 +890,7 @@ make_summary_table (void) break; while (cur_tab < (*pb)->table) - insert_summary (summary, cur_tab++, 0.); + insert_summary (summary, cur_tab++, dict, 0.); if (mode == GENERAL) for (valid = 0.; pb < pe; pb++) @@ -888,13 +911,13 @@ make_summary_table (void) valid += *data++; } } - insert_summary (summary, cur_tab++, valid); + insert_summary (summary, cur_tab++, dict, valid); pb = pe; } while (cur_tab < nxtab) - insert_summary (summary, cur_tab++, 0.); + insert_summary (summary, cur_tab++, dict, 0.); submit (summary); } @@ -902,15 +925,20 @@ make_summary_table (void) /* Inserts a line into T describing the crosstabulation at index TAB_INDEX, which has VALID valid observations. */ static void -insert_summary (struct tab_table *t, int tab_index, double valid) +insert_summary (struct tab_table *t, int tab_index, + const struct dictionary *dict, + double valid) { struct crosstab *x = xtab[tab_index]; + const struct variable *wv = dict_get_weight (dict); + const struct fmt_spec *wfmt = wv ? var_get_print_format (wv) : & F_8_0; + tab_hline (t, TAL_1, 0, 6, 0); /* Crosstabulation name. */ { - char *buf = local_alloc (128 * x->nvar); + char *buf = xmalloca (128 * x->nvar); char *cp = buf; int i; @@ -923,7 +951,7 @@ insert_summary (struct tab_table *t, int tab_index, double valid) } tab_text (t, 0, 0, TAB_LEFT, buf); - local_free (buf); + freea (buf); } /* Counts and percentages. */ @@ -938,7 +966,7 @@ insert_summary (struct tab_table *t, int tab_index, double valid) for (i = 0; i < 3; i++) { - tab_float (t, i * 2 + 1, 0, TAB_RIGHT, n[i], 8, 0); + tab_double (t, i * 2 + 1, 0, TAB_RIGHT, n[i], wfmt); tab_text (t, i * 2 + 2, 0, TAB_RIGHT | TAT_PRINTF, "%.1f%%", n[i] / n[2] * 100.); } @@ -987,11 +1015,11 @@ static double W; /* Grand total. */ static void display_dimensions (struct tab_table *, int first_difference, struct table_entry *); static void display_crosstabulation (void); -static void display_chisq (void); -static void display_symmetric (void); -static void display_risk (void); +static void display_chisq (const struct dictionary *); +static void display_symmetric (const struct dictionary *); +static void display_risk (const struct dictionary *); static void display_directional (void); -static void crosstabs_dim (struct tab_table *, struct outp_driver *); +static void crosstabs_dim (struct tab_table *, struct outp_driver *, void *); static void table_value_missing (struct tab_table *table, int c, int r, unsigned char opt, const union value *v, const struct variable *var); @@ -1002,6 +1030,7 @@ static void delete_missing (void); hold *MAXROWS entries. */ static void output_pivot_table (struct table_entry **pb, struct table_entry **pe, + const struct dictionary *dict, double **matp, double **row_totp, double **col_totp, int *maxrows, int *maxcols, int *maxcells) { @@ -1050,7 +1079,7 @@ output_pivot_table (struct table_entry **pb, struct table_entry **pe, /* Title. */ { - char *title = local_alloc (x->nvar * 64 + 128); + char *title = xmalloca (x->nvar * 64 + 128); char *cp = title; int i; @@ -1116,7 +1145,7 @@ output_pivot_table (struct table_entry **pb, struct table_entry **pe, strcpy (cp, "]."); tab_title (table, "%s", title); - local_free (title); + freea (title); } tab_offset (table, 0, 2); @@ -1408,11 +1437,11 @@ output_pivot_table (struct table_entry **pb, struct table_entry **pe, if (cmd.miss == CRS_REPORT) delete_missing (); if (chisq) - display_chisq (); + display_chisq (dict); if (sym) - display_symmetric (); + display_symmetric (dict); if (risk) - display_risk (); + display_risk (dict); if (direct) display_directional (); @@ -1496,14 +1525,14 @@ submit (struct tab_table *t) tab_box (t, -1, -1, -1, TAL_GAP, 0, tab_t (t), tab_l (t) - 1, tab_nr (t) - 1); tab_vline (t, TAL_2, tab_l (t), 0, tab_nr (t) - 1); - tab_dim (t, crosstabs_dim); + tab_dim (t, crosstabs_dim, NULL); tab_submit (t); } /* Sets the widths of all the columns and heights of all the rows in table T for driver D. */ static void -crosstabs_dim (struct tab_table *t, struct outp_driver *d) +crosstabs_dim (struct tab_table *t, struct outp_driver *d, void *aux UNUSED) { int i; @@ -1657,7 +1686,7 @@ enum_var_values (struct table_entry **entries, int entry_cnt, int var_idx, if (mode == GENERAL) { - int width = var_get_width (v); + int width = MIN (var_get_width (v), MAX_SHORT_STRING); int i; *values = xnmalloc (entry_cnt, sizeof **values); @@ -1844,7 +1873,6 @@ display_crosstabulation (void) tab_offset (table, -1, tab_row (table) - num_cells * n_rows); for (r = 0; r < n_rows; r++) { - char suffix = 0; bool mark_missing = false; if (cmd.miss == CRS_REPORT @@ -1853,6 +1881,7 @@ display_crosstabulation (void) for (i = 0; i < num_cells; i++) { + char suffix = 0; double v; switch (cells[i]) @@ -1861,7 +1890,7 @@ display_crosstabulation (void) v = row_tot[r]; break; case CRS_CL_ROW: - v = 100.; + v = 100.0; suffix = '%'; break; case CRS_CL_COLUMN: @@ -1899,7 +1928,6 @@ display_crosstabulation (void) { double ct = c < n_cols ? col_tot[c] : W; bool mark_missing = false; - char suffix = 0; int i; if (cmd.miss == CRS_REPORT && c < n_cols @@ -1908,13 +1936,13 @@ display_crosstabulation (void) for (i = 0; i < num_cells; i++) { + char suffix = 0; double v; switch (cells[i]) { case CRS_CL_COUNT: v = ct; - suffix = '%'; break; case CRS_CL_ROW: v = ct / W * 100.; @@ -1937,7 +1965,7 @@ display_crosstabulation (void) NOT_REACHED (); } - format_cell_entry (table, c, i, v, suffix, mark_missing); + format_cell_entry (table, c, i, v, suffix, mark_missing); } last_row = i; } @@ -1953,8 +1981,11 @@ static void calc_chisq (double[N_CHISQ], int[N_CHISQ], double *, double *); /* Display chi-square statistics. */ static void -display_chisq (void) +display_chisq (const struct dictionary *dict) { + const struct variable *wv = dict_get_weight (dict); + const struct fmt_spec *wfmt = wv ? var_get_print_format (wv) : & F_8_0; + static const char *chisq_stats[N_CHISQ] = { N_("Pearson Chi-Square"), @@ -1984,22 +2015,22 @@ display_chisq (void) tab_text (chisq, 0, 0, TAB_LEFT, gettext (chisq_stats[i])); if (i != 2) { - tab_float (chisq, 1, 0, TAB_RIGHT, chisq_v[i], 8, 3); - tab_float (chisq, 2, 0, TAB_RIGHT, df[i], 8, 0); - tab_float (chisq, 3, 0, TAB_RIGHT, - gsl_cdf_chisq_Q (chisq_v[i], df[i]), 8, 3); + tab_double (chisq, 1, 0, TAB_RIGHT, chisq_v[i], NULL); + tab_double (chisq, 2, 0, TAB_RIGHT, df[i], wfmt); + tab_double (chisq, 3, 0, TAB_RIGHT, + gsl_cdf_chisq_Q (chisq_v[i], df[i]), NULL); } else { chisq_fisher = 1; - tab_float (chisq, 4, 0, TAB_RIGHT, fisher2, 8, 3); - tab_float (chisq, 5, 0, TAB_RIGHT, fisher1, 8, 3); + tab_double (chisq, 4, 0, TAB_RIGHT, fisher2, NULL); + tab_double (chisq, 5, 0, TAB_RIGHT, fisher1, NULL); } tab_next_row (chisq); } tab_text (chisq, 0, 0, TAB_LEFT, _("N of Valid Cases")); - tab_float (chisq, 1, 0, TAB_RIGHT, W, 8, 0); + tab_double (chisq, 1, 0, TAB_RIGHT, W, wfmt); tab_next_row (chisq); tab_offset (chisq, 0, -1); @@ -2010,8 +2041,11 @@ static int calc_symmetric (double[N_SYMMETRIC], double[N_SYMMETRIC], /* Display symmetric measures. */ static void -display_symmetric (void) +display_symmetric (const struct dictionary *dict) { + const struct variable *wv = dict_get_weight (dict); + const struct fmt_spec *wfmt = wv ? var_get_print_format (wv) : & F_8_0; + static const char *categories[] = { N_("Nominal by Nominal"), @@ -2059,17 +2093,17 @@ display_symmetric (void) } tab_text (sym, 1, 0, TAB_LEFT, gettext (stats[i])); - tab_float (sym, 2, 0, TAB_RIGHT, sym_v[i], 8, 3); + tab_double (sym, 2, 0, TAB_RIGHT, sym_v[i], NULL); if (sym_ase[i] != SYSMIS) - tab_float (sym, 3, 0, TAB_RIGHT, sym_ase[i], 8, 3); + tab_double (sym, 3, 0, TAB_RIGHT, sym_ase[i], NULL); if (sym_t[i] != SYSMIS) - tab_float (sym, 4, 0, TAB_RIGHT, sym_t[i], 8, 3); - /*tab_float (sym, 5, 0, TAB_RIGHT, normal_sig (sym_v[i]), 8, 3);*/ + tab_double (sym, 4, 0, TAB_RIGHT, sym_t[i], NULL); + /*tab_double (sym, 5, 0, TAB_RIGHT, normal_sig (sym_v[i]), NULL);*/ tab_next_row (sym); } tab_text (sym, 0, 0, TAB_LEFT, _("N of Valid Cases")); - tab_float (sym, 2, 0, TAB_RIGHT, W, 8, 0); + tab_double (sym, 2, 0, TAB_RIGHT, W, wfmt); tab_next_row (sym); tab_offset (sym, 0, -1); @@ -2079,8 +2113,11 @@ static int calc_risk (double[], double[], double[], union value *); /* Display risk estimate. */ static void -display_risk (void) +display_risk (const struct dictionary *dict) { + const struct variable *wv = dict_get_weight (dict); + const struct fmt_spec *wfmt = wv ? var_get_print_format (wv) : & F_8_0; + char buf[256]; double risk_v[3], lower[3], upper[3]; union value c[2]; @@ -2121,14 +2158,14 @@ display_risk (void) } tab_text (risk, 0, 0, TAB_LEFT, buf); - tab_float (risk, 1, 0, TAB_RIGHT, risk_v[i], 8, 3); - tab_float (risk, 2, 0, TAB_RIGHT, lower[i], 8, 3); - tab_float (risk, 3, 0, TAB_RIGHT, upper[i], 8, 3); + tab_double (risk, 1, 0, TAB_RIGHT, risk_v[i], NULL); + tab_double (risk, 2, 0, TAB_RIGHT, lower[i], NULL); + tab_double (risk, 3, 0, TAB_RIGHT, upper[i], NULL); tab_next_row (risk); } tab_text (risk, 0, 0, TAB_LEFT, _("N of Valid Cases")); - tab_float (risk, 1, 0, TAB_RIGHT, W, 8, 0); + tab_double (risk, 1, 0, TAB_RIGHT, W, wfmt); tab_next_row (risk); tab_offset (risk, 0, -1); @@ -2241,12 +2278,12 @@ display_directional (void) } } - tab_float (direct, 3, 0, TAB_RIGHT, direct_v[i], 8, 3); + tab_double (direct, 3, 0, TAB_RIGHT, direct_v[i], NULL); if (direct_ase[i] != SYSMIS) - tab_float (direct, 4, 0, TAB_RIGHT, direct_ase[i], 8, 3); + tab_double (direct, 4, 0, TAB_RIGHT, direct_ase[i], NULL); if (direct_t[i] != SYSMIS) - tab_float (direct, 5, 0, TAB_RIGHT, direct_t[i], 8, 3); - /*tab_float (direct, 6, 0, TAB_RIGHT, normal_sig (direct_v[i]), 8, 3);*/ + tab_double (direct, 5, 0, TAB_RIGHT, direct_t[i], NULL); + /*tab_double (direct, 6, 0, TAB_RIGHT, normal_sig (direct_v[i]), NULL);*/ tab_next_row (direct); } @@ -2439,7 +2476,7 @@ calc_r (double *X, double *Y, double *r, double *ase_0, double *ase_1) for (sum_Xr = sum_X2r = 0., i = 0; i < n_rows; i++) { sum_Xr += X[i] * row_tot[i]; - sum_X2r += X[i] * X[i] * row_tot[i]; + sum_X2r += pow2 (X[i]) * row_tot[i]; } Xbar = sum_Xr / W; @@ -2451,11 +2488,11 @@ calc_r (double *X, double *Y, double *r, double *ase_0, double *ase_1) Ybar = sum_Yc / W; S = sum_XYf - sum_Xr * sum_Yc / W; - SX = sum_X2r - sum_Xr * sum_Xr / W; - SY = sum_Y2c - sum_Yc * sum_Yc / W; + SX = sum_X2r - pow2 (sum_Xr) / W; + SY = sum_Y2c - pow2 (sum_Yc) / W; T = sqrt (SX * SY); *r = S / T; - *ase_0 = sqrt ((sum_X2Y2f - (sum_XYf * sum_XYf) / W) / (sum_X2r * sum_Y2c)); + *ase_0 = sqrt ((sum_X2Y2f - pow2 (sum_XYf) / W) / (sum_X2r * sum_Y2c)); { double s, c, y, t; @@ -2545,9 +2582,9 @@ calc_symmetric (double v[N_SYMMETRIC], double ase[N_SYMMETRIC], Dr = Dc = W * W; for (r = 0; r < n_rows; r++) - Dr -= row_tot[r] * row_tot[r]; + Dr -= pow2 (row_tot[r]); for (c = 0; c < n_cols; c++) - Dc -= col_tot[c] * col_tot[c]; + Dc -= pow2 (col_tot[c]); } { @@ -2720,8 +2757,8 @@ calc_symmetric (double v[N_SYMMETRIC], double ase[N_SYMMETRIC], /* Spearman correlation, Pearson's r. */ if (cmd.a_statistics[CRS_ST_CORR]) { - double *R = local_alloc (sizeof *R * n_rows); - double *C = local_alloc (sizeof *C * n_cols); + double *R = xmalloca (sizeof *R * n_rows); + double *C = xmalloca (sizeof *C * n_cols); { double y, t, c = 0., s = 0.; @@ -2760,8 +2797,8 @@ calc_symmetric (double v[N_SYMMETRIC], double ase[N_SYMMETRIC], calc_r (R, C, &v[6], &t[6], &ase[6]); t[6] = v[6] / t[6]; - local_free (R); - local_free (C); + freea (R); + freea (C); calc_r ((double *) rows, (double *) cols, &v[7], &t[7], &ase[7]); t[7] = v[7] / t[7]; @@ -3056,10 +3093,10 @@ calc_directional (double v[N_DIRECTIONAL], double ase[N_DIRECTIONAL], } for (sum_ri2 = 0., i = 0; i < n_rows; i++) - sum_ri2 += row_tot[i] * row_tot[i]; + sum_ri2 += pow2 (row_tot[i]); for (sum_cj2 = 0., j = 0; j < n_cols; j++) - sum_cj2 += col_tot[j] * col_tot[j]; + sum_cj2 += pow2 (col_tot[j]); v[3] = (W * sum_fij2_ci - sum_ri2) / (W * W - sum_ri2); v[4] = (W * sum_fij2_ri - sum_cj2) / (W * W - sum_cj2); @@ -3149,9 +3186,9 @@ calc_directional (double v[N_DIRECTIONAL], double ase[N_DIRECTIONAL], for (sum_Xr = sum_X2r = 0., i = 0; i < n_rows; i++) { sum_Xr += rows[i].f * row_tot[i]; - sum_X2r += rows[i].f * rows[i].f * row_tot[i]; + sum_X2r += pow2 (rows[i].f) * row_tot[i]; } - SX = sum_X2r - sum_Xr * sum_Xr / W; + SX = sum_X2r - pow2 (sum_Xr) / W; for (SXW = 0., j = 0; j < n_cols; j++) { @@ -3159,7 +3196,7 @@ calc_directional (double v[N_DIRECTIONAL], double ase[N_DIRECTIONAL], for (cum = 0., i = 0; i < n_rows; i++) { - SXW += rows[i].f * rows[i].f * mat[j + i * n_cols]; + SXW += pow2 (rows[i].f) * mat[j + i * n_cols]; cum += rows[i].f * mat[j + i * n_cols]; } @@ -3176,7 +3213,7 @@ calc_directional (double v[N_DIRECTIONAL], double ase[N_DIRECTIONAL], for (sum_Yc = sum_Y2c = 0., i = 0; i < n_cols; i++) { sum_Yc += cols[i].f * col_tot[i]; - sum_Y2c += cols[i].f * cols[i].f * col_tot[i]; + sum_Y2c += pow2 (cols[i].f) * col_tot[i]; } SY = sum_Y2c - sum_Yc * sum_Yc / W; @@ -3186,7 +3223,7 @@ calc_directional (double v[N_DIRECTIONAL], double ase[N_DIRECTIONAL], for (cum = 0., j = 0; j < n_cols; j++) { - SYW += cols[j].f * cols[j].f * mat[j + i * n_cols]; + SYW += pow2 (cols[j].f) * mat[j + i * n_cols]; cum += cols[j].f * mat[j + i * n_cols]; }