src/language/stats/aggregate.c

   1 /* PSPP - a program for statistical analysis.
   2    Copyright (C) 1997-9, 2000, 2006, 2008, 2009 Free Software Foundation, Inc.
   3
   4    This program is free software: you can redistribute it and/or modify
   5    it under the terms of the GNU General Public License as published by
   6    the Free Software Foundation, either version 3 of the License, or
   7    (at your option) any later version.
   8
   9    This program is distributed in the hope that it will be useful,
  10    but WITHOUT ANY WARRANTY; without even the implied warranty of
  11    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  12    GNU General Public License for more details.
  13
  14    You should have received a copy of the GNU General Public License
  15    along with this program.  If not, see <http://www.gnu.org/licenses/>. */
  16
  17 #include <config.h>
  18
  19 #include <stdlib.h>
  20
  21 #include <data/any-writer.h>
  22 #include <data/case.h>
  23 #include <data/casegrouper.h>
  24 #include <data/casereader.h>
  25 #include <data/casewriter.h>
  26 #include <data/dictionary.h>
  27 #include <data/file-handle-def.h>
  28 #include <data/format.h>
  29 #include <data/procedure.h>
  30 #include <data/settings.h>
  31 #include <data/subcase.h>
  32 #include <data/sys-file-writer.h>
  33 #include <data/variable.h>
  34 #include <language/command.h>
  35 #include <language/data-io/file-handle.h>
  36 #include <language/lexer/lexer.h>
  37 #include <language/lexer/variable-parser.h>
  38 #include <language/stats/sort-criteria.h>
  39 #include <libpspp/assertion.h>
  40 #include <libpspp/message.h>
  41 #include <libpspp/misc.h>
  42 #include <libpspp/pool.h>
  43 #include <libpspp/str.h>
  44 #include <math/moments.h>
  45 #include <math/sort.h>
  46 #include <math/statistic.h>
  47 #include <math/percentiles.h>
  48
  49 #include "minmax.h"
  50 #include "xalloc.h"
  51
  52 #include "gettext.h"
  53 #define _(msgid) gettext (msgid)
  54
  55 /* Argument for AGGREGATE function. */
  56 union agr_argument
  57   {
  58     double f;                           /* Numeric. */
  59     char *c;                            /* Short or long string. */
  60   };
  61
  62 /* Specifies how to make an aggregate variable. */
  63 struct agr_var
  64   {
  65     struct agr_var *next;               /* Next in list. */
  66
  67     /* Collected during parsing. */
  68     const struct variable *src; /* Source variable. */
  69     struct variable *dest;      /* Target variable. */
  70     int function;               /* Function. */
  71     enum mv_class exclude;      /* Classes of missing values to exclude. */
  72     union agr_argument arg[2];  /* Arguments. */
  73
  74     /* Accumulated during AGGREGATE execution. */
  75     double dbl[3];
  76     int int1, int2;
  77     char *string;
  78     bool saw_missing;
  79     struct moments1 *moments;
  80     double cc;
  81
  82     struct variable *subject;
  83     struct variable *weight;
  84     struct casewriter *writer;
  85   };
  86
  87 /* Aggregation functions. */
  88 enum
  89   {
  90     NONE, SUM, MEAN, MEDIAN, SD, MAX, MIN, PGT, PLT, PIN, POUT, FGT, FLT, FIN,
  91     FOUT, N, NU, NMISS, NUMISS, FIRST, LAST,
  92     N_AGR_FUNCS, N_NO_VARS, NU_NO_VARS,
  93     FUNC = 0x1f, /* Function mask. */
  94     FSTRING = 1<<5, /* String function bit. */
  95   };
  96
  97 /* Attributes of an aggregation function. */
  98 struct agr_func
  99   {
 100     const char *name;           /* Aggregation function name. */
 101     size_t n_args;              /* Number of arguments. */
 102     enum val_type alpha_type;   /* When given ALPHA arguments, output type. */
 103     struct fmt_spec format;     /* Format spec if alpha_type != ALPHA. */
 104   };
 105
 106 /* Attributes of aggregation functions. */
 107 static const struct agr_func agr_func_tab[] =
 108   {
 109     {"<NONE>",  0, -1,          {0, 0, 0}},
 110     {"SUM",     0, -1,          {FMT_F, 8, 2}},
 111     {"MEAN",    0, -1,          {FMT_F, 8, 2}},
 112     {"MEDIAN",  0, -1,          {FMT_F, 8, 2}},
 113     {"SD",      0, -1,          {FMT_F, 8, 2}},
 114     {"MAX",     0, VAL_STRING,  {-1, -1, -1}},
 115     {"MIN",     0, VAL_STRING,  {-1, -1, -1}},
 116     {"PGT",     1, VAL_NUMERIC, {FMT_F, 5, 1}},
 117     {"PLT",     1, VAL_NUMERIC, {FMT_F, 5, 1}},
 118     {"PIN",     2, VAL_NUMERIC, {FMT_F, 5, 1}},
 119     {"POUT",    2, VAL_NUMERIC, {FMT_F, 5, 1}},
 120     {"FGT",     1, VAL_NUMERIC, {FMT_F, 5, 3}},
 121     {"FLT",     1, VAL_NUMERIC, {FMT_F, 5, 3}},
 122     {"FIN",     2, VAL_NUMERIC, {FMT_F, 5, 3}},
 123     {"FOUT",    2, VAL_NUMERIC, {FMT_F, 5, 3}},
 124     {"N",       0, VAL_NUMERIC, {FMT_F, 7, 0}},
 125     {"NU",      0, VAL_NUMERIC, {FMT_F, 7, 0}},
 126     {"NMISS",   0, VAL_NUMERIC, {FMT_F, 7, 0}},
 127     {"NUMISS",  0, VAL_NUMERIC, {FMT_F, 7, 0}},
 128     {"FIRST",   0, VAL_STRING,  {-1, -1, -1}},
 129     {"LAST",    0, VAL_STRING,  {-1, -1, -1}},
 130     {NULL,      0, -1,          {-1, -1, -1}},
 131     {"N",       0, VAL_NUMERIC, {FMT_F, 7, 0}},
 132     {"NU",      0, VAL_NUMERIC, {FMT_F, 7, 0}},
 133   };
 134
 135 /* Missing value types. */
 136 enum missing_treatment
 137   {
 138     ITEMWISE,           /* Missing values item by item. */
 139     COLUMNWISE          /* Missing values column by column. */
 140   };
 141
 142 /* An entire AGGREGATE procedure. */
 143 struct agr_proc
 144   {
 145     /* Break variables. */
 146     struct subcase sort;                /* Sort criteria (break variables). */
 147     const struct variable **break_vars;       /* Break variables. */
 148     size_t break_var_cnt;               /* Number of break variables. */
 149     struct ccase *break_case;           /* Last values of break variables. */
 150
 151     enum missing_treatment missing;     /* How to treat missing values. */
 152     struct agr_var *agr_vars;           /* First aggregate variable. */
 153     struct dictionary *dict;            /* Aggregate dictionary. */
 154     const struct dictionary *src_dict;  /* Dict of the source */
 155     int case_cnt;                       /* Counts aggregated cases. */
 156   };
 157
 158 static void initialize_aggregate_info (struct agr_proc *,
 159                                        const struct ccase *);
 160
 161 static void accumulate_aggregate_info (struct agr_proc *,
 162                                        const struct ccase *);
 163 /* Prototypes. */
 164 static bool parse_aggregate_functions (struct lexer *, const struct dictionary *,
 165                                        struct agr_proc *);
 166 static void agr_destroy (struct agr_proc *);
 167 static void dump_aggregate_info (struct agr_proc *agr,
 168                                  struct casewriter *output);
 169 \f
 170 /* Parsing. */
 171
 172 /* Parses and executes the AGGREGATE procedure. */
 173 int
 174 cmd_aggregate (struct lexer *lexer, struct dataset *ds)
 175 {
 176   struct dictionary *dict = dataset_dict (ds);
 177   struct agr_proc agr;
 178   struct file_handle *out_file = NULL;
 179   struct casereader *input = NULL, *group;
 180   struct casegrouper *grouper;
 181   struct casewriter *output = NULL;
 182
 183   bool copy_documents = false;
 184   bool presorted = false;
 185   bool saw_direction;
 186   bool ok;
 187
 188   memset(&agr, 0 , sizeof (agr));
 189   agr.missing = ITEMWISE;
 190   agr.break_case = NULL;
 191
 192   agr.dict = dict_create ();
 193   agr.src_dict = dict;
 194   subcase_init_empty (&agr.sort);
 195   dict_set_label (agr.dict, dict_get_label (dict));
 196   dict_set_documents (agr.dict, dict_get_documents (dict));
 197
 198   /* OUTFILE subcommand must be first. */
 199   if (!lex_force_match_id (lexer, "OUTFILE"))
 200     goto error;
 201   lex_match (lexer, '=');
 202   if (!lex_match (lexer, '*'))
 203     {
 204       out_file = fh_parse (lexer, FH_REF_FILE | FH_REF_SCRATCH);
 205       if (out_file == NULL)
 206         goto error;
 207     }
 208
 209   /* Read most of the subcommands. */
 210   for (;;)
 211     {
 212       lex_match (lexer, '/');
 213
 214       if (lex_match_id (lexer, "MISSING"))
 215         {
 216           lex_match (lexer, '=');
 217           if (!lex_match_id (lexer, "COLUMNWISE"))
 218             {
 219               lex_error (lexer, _("while expecting COLUMNWISE"));
 220               goto error;
 221             }
 222           agr.missing = COLUMNWISE;
 223         }
 224       else if (lex_match_id (lexer, "DOCUMENT"))
 225         copy_documents = true;
 226       else if (lex_match_id (lexer, "PRESORTED"))
 227         presorted = true;
 228       else if (lex_match_id (lexer, "BREAK"))
 229         {
 230           int i;
 231
 232           lex_match (lexer, '=');
 233           if (!parse_sort_criteria (lexer, dict, &agr.sort, &agr.break_vars,
 234                                     &saw_direction))
 235             goto error;
 236           agr.break_var_cnt = subcase_get_n_fields (&agr.sort);
 237
 238           for (i = 0; i < agr.break_var_cnt; i++)
 239             dict_clone_var_assert (agr.dict, agr.break_vars[i],
 240                                    var_get_name (agr.break_vars[i]));
 241
 242           /* BREAK must follow the options. */
 243           break;
 244         }
 245       else
 246         {
 247           lex_error (lexer, _("expecting BREAK"));
 248           goto error;
 249         }
 250     }
 251   if (presorted && saw_direction)
 252     msg (SW, _("When PRESORTED is specified, specifying sorting directions "
 253                "with (A) or (D) has no effect.  Output data will be sorted "
 254                "the same way as the input data."));
 255
 256   /* Read in the aggregate functions. */
 257   lex_match (lexer, '/');
 258   if (!parse_aggregate_functions (lexer, dict, &agr))
 259     goto error;
 260
 261   /* Delete documents. */
 262   if (!copy_documents)
 263     dict_clear_documents (agr.dict);
 264
 265   /* Cancel SPLIT FILE. */
 266   dict_set_split_vars (agr.dict, NULL, 0);
 267
 268   /* Initialize. */
 269   agr.case_cnt = 0;
 270
 271   if (out_file == NULL)
 272     {
 273       /* The active file will be replaced by the aggregated data,
 274          so TEMPORARY is moot. */
 275       proc_cancel_temporary_transformations (ds);
 276       proc_discard_output (ds);
 277       output = autopaging_writer_create (dict_get_proto (agr.dict));
 278     }
 279   else
 280     {
 281       output = any_writer_open (out_file, agr.dict);
 282       if (output == NULL)
 283         goto error;
 284     }
 285
 286   input = proc_open (ds);
 287   if (!subcase_is_empty (&agr.sort) && !presorted)
 288     {
 289       input = sort_execute (input, &agr.sort);
 290       subcase_clear (&agr.sort);
 291     }
 292
 293   for (grouper = casegrouper_create_vars (input, agr.break_vars,
 294                                           agr.break_var_cnt);
 295        casegrouper_get_next_group (grouper, &group);
 296        casereader_destroy (group))
 297     {
 298       struct ccase *c = casereader_peek (group, 0);
 299       if (c == NULL)
 300         {
 301           casereader_destroy (group);
 302           continue;
 303         }
 304       initialize_aggregate_info (&agr, c);
 305       case_unref (c);
 306
 307       for (; (c = casereader_read (group)) != NULL; case_unref (c))
 308         accumulate_aggregate_info (&agr, c);
 309       dump_aggregate_info (&agr, output);
 310     }
 311   if (!casegrouper_destroy (grouper))
 312     goto error;
 313
 314   if (!proc_commit (ds))
 315     {
 316       input = NULL;
 317       goto error;
 318     }
 319   input = NULL;
 320
 321   if (out_file == NULL)
 322     {
 323       struct casereader *next_input = casewriter_make_reader (output);
 324       if (next_input == NULL)
 325         goto error;
 326
 327       proc_set_active_file (ds, next_input, agr.dict);
 328       agr.dict = NULL;
 329     }
 330   else
 331     {
 332       ok = casewriter_destroy (output);
 333       output = NULL;
 334       if (!ok)
 335         goto error;
 336     }
 337
 338   agr_destroy (&agr);
 339   fh_unref (out_file);
 340   return CMD_SUCCESS;
 341
 342 error:
 343   if (input != NULL)
 344     proc_commit (ds);
 345   casewriter_destroy (output);
 346   agr_destroy (&agr);
 347   fh_unref (out_file);
 348   return CMD_CASCADING_FAILURE;
 349 }
 350
 351 /* Parse all the aggregate functions. */
 352 static bool
 353 parse_aggregate_functions (struct lexer *lexer, const struct dictionary *dict,
 354                            struct agr_proc *agr)
 355 {
 356   struct agr_var *tail; /* Tail of linked list starting at agr->vars. */
 357
 358   /* Parse everything. */
 359   tail = NULL;
 360   for (;;)
 361     {
 362       char **dest;
 363       char **dest_label;
 364       size_t n_dest;
 365       struct string function_name;
 366
 367       enum mv_class exclude;
 368       const struct agr_func *function;
 369       int func_index;
 370
 371       union agr_argument arg[2];
 372
 373       const struct variable **src;
 374       size_t n_src;
 375
 376       size_t i;
 377
 378       dest = NULL;
 379       dest_label = NULL;
 380       n_dest = 0;
 381       src = NULL;
 382       function = NULL;
 383       n_src = 0;
 384       arg[0].c = NULL;
 385       arg[1].c = NULL;
 386       ds_init_empty (&function_name);
 387
 388       /* Parse the list of target variables. */
 389       while (!lex_match (lexer, '='))
 390         {
 391           size_t n_dest_prev = n_dest;
 392
 393           if (!parse_DATA_LIST_vars (lexer, &dest, &n_dest,
 394                                      PV_APPEND | PV_SINGLE | PV_NO_SCRATCH))
 395             goto error;
 396
 397           /* Assign empty labels. */
 398           {
 399             int j;
 400
 401             dest_label = xnrealloc (dest_label, n_dest, sizeof *dest_label);
 402             for (j = n_dest_prev; j < n_dest; j++)
 403               dest_label[j] = NULL;
 404           }
 405
 406
 407
 408           if (lex_token (lexer) == T_STRING)
 409             {
 410               struct string label;
 411               ds_init_string (&label, lex_tokstr (lexer));
 412
 413               ds_truncate (&label, 255);
 414               dest_label[n_dest - 1] = ds_xstrdup (&label);
 415               lex_get (lexer);
 416               ds_destroy (&label);
 417             }
 418         }
 419
 420       /* Get the name of the aggregation function. */
 421       if (lex_token (lexer) != T_ID)
 422         {
 423           lex_error (lexer, _("expecting aggregation function"));
 424           goto error;
 425         }
 426
 427       exclude = MV_ANY;
 428
 429       ds_assign_string (&function_name, lex_tokstr (lexer));
 430
 431       ds_chomp (&function_name, '.');
 432
 433       if (lex_tokid(lexer)[strlen (lex_tokid (lexer)) - 1] == '.')
 434         exclude = MV_SYSTEM;
 435
 436       for (function = agr_func_tab; function->name; function++)
 437         if (!strcasecmp (function->name, ds_cstr (&function_name)))
 438           break;
 439       if (NULL == function->name)
 440         {
 441           msg (SE, _("Unknown aggregation function %s."),
 442                ds_cstr (&function_name));
 443           goto error;
 444         }
 445       ds_destroy (&function_name);
 446       func_index = function - agr_func_tab;
 447       lex_get (lexer);
 448
 449       /* Check for leading lparen. */
 450       if (!lex_match (lexer, '('))
 451         {
 452           if (func_index == N)
 453             func_index = N_NO_VARS;
 454           else if (func_index == NU)
 455             func_index = NU_NO_VARS;
 456           else
 457             {
 458               lex_error (lexer, _("expecting `('"));
 459               goto error;
 460             }
 461         }
 462       else
 463         {
 464           /* Parse list of source variables. */
 465           {
 466             int pv_opts = PV_NO_SCRATCH;
 467
 468             if (func_index == SUM || func_index == MEAN || func_index == SD)
 469               pv_opts |= PV_NUMERIC;
 470             else if (function->n_args)
 471               pv_opts |= PV_SAME_TYPE;
 472
 473             if (!parse_variables_const (lexer, dict, &src, &n_src, pv_opts))
 474               goto error;
 475           }
 476
 477           /* Parse function arguments, for those functions that
 478              require arguments. */
 479           if (function->n_args != 0)
 480             for (i = 0; i < function->n_args; i++)
 481               {
 482                 int type;
 483
 484                 lex_match (lexer, ',');
 485                 if (lex_token (lexer) == T_STRING)
 486                   {
 487                     arg[i].c = ds_xstrdup (lex_tokstr (lexer));
 488                     type = VAL_STRING;
 489                   }
 490                 else if (lex_is_number (lexer))
 491                   {
 492                     arg[i].f = lex_tokval (lexer);
 493                     type = VAL_NUMERIC;
 494                   }
 495                 else
 496                   {
 497                     msg (SE, _("Missing argument %zu to %s."),
 498                          i + 1, function->name);
 499                     goto error;
 500                   }
 501
 502                 lex_get (lexer);
 503
 504                 if (type != var_get_type (src[0]))
 505                   {
 506                     msg (SE, _("Arguments to %s must be of same type as "
 507                                "source variables."),
 508                          function->name);
 509                     goto error;
 510                   }
 511               }
 512
 513           /* Trailing rparen. */
 514           if (!lex_match (lexer, ')'))
 515             {
 516               lex_error (lexer, _("expecting `)'"));
 517               goto error;
 518             }
 519
 520           /* Now check that the number of source variables match
 521              the number of target variables.  If we check earlier
 522              than this, the user can get very misleading error
 523              message, i.e. `AGGREGATE x=SUM(y t).' will get this
 524              error message when a proper message would be more
 525              like `unknown variable t'. */
 526           if (n_src != n_dest)
 527             {
 528               msg (SE, _("Number of source variables (%zu) does not match "
 529                          "number of target variables (%zu)."),
 530                     n_src, n_dest);
 531               goto error;
 532             }
 533
 534           if ((func_index == PIN || func_index == POUT
 535               || func_index == FIN || func_index == FOUT)
 536               && (var_is_numeric (src[0])
 537                   ? arg[0].f > arg[1].f
 538                   : str_compare_rpad (arg[0].c, arg[1].c) > 0))
 539             {
 540               union agr_argument t = arg[0];
 541               arg[0] = arg[1];
 542               arg[1] = t;
 543
 544               msg (SW, _("The value arguments passed to the %s function "
 545                          "are out-of-order.  They will be treated as if "
 546                          "they had been specified in the correct order."),
 547                    function->name);
 548             }
 549         }
 550
 551       /* Finally add these to the linked list of aggregation
 552          variables. */
 553       for (i = 0; i < n_dest; i++)
 554         {
 555           struct agr_var *v = xzalloc (sizeof *v);
 556
 557           /* Add variable to chain. */
 558           if (agr->agr_vars != NULL)
 559             tail->next = v;
 560           else
 561             agr->agr_vars = v;
 562           tail = v;
 563           tail->next = NULL;
 564           v->moments = NULL;
 565
 566           /* Create the target variable in the aggregate
 567              dictionary. */
 568           {
 569             struct variable *destvar;
 570
 571             v->function = func_index;
 572
 573             if (src)
 574               {
 575                 v->src = src[i];
 576
 577                 if (var_is_alpha (src[i]))
 578                   {
 579                     v->function |= FSTRING;
 580                     v->string = xmalloc (var_get_width (src[i]));
 581                   }
 582
 583                 if (function->alpha_type == VAL_STRING)
 584                   destvar = dict_clone_var (agr->dict, v->src, dest[i]);
 585                 else
 586                   {
 587                     assert (var_is_numeric (v->src)
 588                             || function->alpha_type == VAL_NUMERIC);
 589                     destvar = dict_create_var (agr->dict, dest[i], 0);
 590                     if (destvar != NULL)
 591                       {
 592                         struct fmt_spec f;
 593                         if ((func_index == N || func_index == NMISS)
 594                             && dict_get_weight (dict) != NULL)
 595                           f = fmt_for_output (FMT_F, 8, 2);
 596                         else
 597                           f = function->format;
 598                         var_set_both_formats (destvar, &f);
 599                       }
 600                   }
 601               } else {
 602                 struct fmt_spec f;
 603                 v->src = NULL;
 604                 destvar = dict_create_var (agr->dict, dest[i], 0);
 605                 if (func_index == N_NO_VARS && dict_get_weight (dict) != NULL)
 606                   f = fmt_for_output (FMT_F, 8, 2);
 607                 else
 608                   f = function->format;
 609                 var_set_both_formats (destvar, &f);
 610               }
 611
 612             if (!destvar)
 613               {
 614                 msg (SE, _("Variable name %s is not unique within the "
 615                            "aggregate file dictionary, which contains "
 616                            "the aggregate variables and the break "
 617                            "variables."),
 618                      dest[i]);
 619                 goto error;
 620               }
 621
 622             free (dest[i]);
 623             if (dest_label[i])
 624               var_set_label (destvar, dest_label[i]);
 625
 626             v->dest = destvar;
 627           }
 628
 629           v->exclude = exclude;
 630
 631           if (v->src != NULL)
 632             {
 633               int j;
 634
 635               if (var_is_numeric (v->src))
 636                 for (j = 0; j < function->n_args; j++)
 637                   v->arg[j].f = arg[j].f;
 638               else
 639                 for (j = 0; j < function->n_args; j++)
 640                   v->arg[j].c = xstrdup (arg[j].c);
 641             }
 642         }
 643
 644       if (src != NULL && var_is_alpha (src[0]))
 645         for (i = 0; i < function->n_args; i++)
 646           {
 647             free (arg[i].c);
 648             arg[i].c = NULL;
 649           }
 650
 651       free (src);
 652       free (dest);
 653       free (dest_label);
 654
 655       if (!lex_match (lexer, '/'))
 656         {
 657           if (lex_token (lexer) == '.')
 658             return true;
 659
 660           lex_error (lexer, "expecting end of command");
 661           return false;
 662         }
 663       continue;
 664
 665     error:
 666       ds_destroy (&function_name);
 667       for (i = 0; i < n_dest; i++)
 668         {
 669           free (dest[i]);
 670           free (dest_label[i]);
 671         }
 672       free (dest);
 673       free (dest_label);
 674       free (arg[0].c);
 675       free (arg[1].c);
 676       if (src && n_src && var_is_alpha (src[0]))
 677         for (i = 0; i < function->n_args; i++)
 678           {
 679             free (arg[i].c);
 680             arg[i].c = NULL;
 681           }
 682       free (src);
 683
 684       return false;
 685     }
 686 }
 687
 688 /* Destroys AGR. */
 689 static void
 690 agr_destroy (struct agr_proc *agr)
 691 {
 692   struct agr_var *iter, *next;
 693
 694   subcase_destroy (&agr->sort);
 695   free (agr->break_vars);
 696   case_unref (agr->break_case);
 697   for (iter = agr->agr_vars; iter; iter = next)
 698     {
 699       next = iter->next;
 700
 701       if (iter->function & FSTRING)
 702         {
 703           size_t n_args;
 704           size_t i;
 705
 706           n_args = agr_func_tab[iter->function & FUNC].n_args;
 707           for (i = 0; i < n_args; i++)
 708             free (iter->arg[i].c);
 709           free (iter->string);
 710         }
 711       else if (iter->function == SD)
 712         moments1_destroy (iter->moments);
 713
 714       var_destroy (iter->subject);
 715       var_destroy (iter->weight);
 716
 717       free (iter);
 718     }
 719   if (agr->dict != NULL)
 720     dict_destroy (agr->dict);
 721 }
 722 \f
 723 /* Execution. */
 724
 725 /* Accumulates aggregation data from the case INPUT. */
 726 static void
 727 accumulate_aggregate_info (struct agr_proc *agr, const struct ccase *input)
 728 {
 729   struct agr_var *iter;
 730   double weight;
 731   bool bad_warn = true;
 732
 733   weight = dict_get_case_weight (agr->src_dict, input, &bad_warn);
 734
 735   for (iter = agr->agr_vars; iter; iter = iter->next)
 736     if (iter->src)
 737       {
 738         const union value *v = case_data (input, iter->src);
 739         int src_width = var_get_width (iter->src);
 740
 741         if (var_is_value_missing (iter->src, v, iter->exclude))
 742           {
 743             switch (iter->function)
 744               {
 745               case NMISS:
 746               case NMISS | FSTRING:
 747                 iter->dbl[0] += weight;
 748                 break;
 749               case NUMISS:
 750               case NUMISS | FSTRING:
 751                 iter->int1++;
 752                 break;
 753               }
 754             iter->saw_missing = true;
 755             continue;
 756           }
 757
 758         /* This is horrible.  There are too many possibilities. */
 759         switch (iter->function)
 760           {
 761           case SUM:
 762             iter->dbl[0] += v->f * weight;
 763             iter->int1 = 1;
 764             break;
 765           case MEAN:
 766             iter->dbl[0] += v->f * weight;
 767             iter->dbl[1] += weight;
 768             break;
 769           case MEDIAN:
 770             {
 771               double wv ;
 772               struct ccase *cout;
 773
 774               cout = case_create (casewriter_get_proto (iter->writer));
 775
 776               case_data_rw (cout, iter->subject)->f
 777                 = case_data (input, iter->src)->f;
 778
 779               wv = dict_get_case_weight (agr->src_dict, input, NULL);
 780
 781               case_data_rw (cout, iter->weight)->f = wv;
 782
 783               iter->cc += wv;
 784
 785               casewriter_write (iter->writer, cout);
 786             }
 787             break;
 788           case SD:
 789             moments1_add (iter->moments, v->f, weight);
 790             break;
 791           case MAX:
 792             iter->dbl[0] = MAX (iter->dbl[0], v->f);
 793             iter->int1 = 1;
 794             break;
 795           case MAX | FSTRING:
 796             if (memcmp (iter->string, value_str (v, src_width), src_width) < 0)
 797               memcpy (iter->string, value_str (v, src_width), src_width);
 798             iter->int1 = 1;
 799             break;
 800           case MIN:
 801             iter->dbl[0] = MIN (iter->dbl[0], v->f);
 802             iter->int1 = 1;
 803             break;
 804           case MIN | FSTRING:
 805             if (memcmp (iter->string, value_str (v, src_width), src_width) > 0)
 806               memcpy (iter->string, value_str (v, src_width), src_width);
 807             iter->int1 = 1;
 808             break;
 809           case FGT:
 810           case PGT:
 811             if (v->f > iter->arg[0].f)
 812               iter->dbl[0] += weight;
 813             iter->dbl[1] += weight;
 814             break;
 815           case FGT | FSTRING:
 816           case PGT | FSTRING:
 817             if (memcmp (iter->arg[0].c,
 818                         value_str (v, src_width), src_width) < 0)
 819               iter->dbl[0] += weight;
 820             iter->dbl[1] += weight;
 821             break;
 822           case FLT:
 823           case PLT:
 824             if (v->f < iter->arg[0].f)
 825               iter->dbl[0] += weight;
 826             iter->dbl[1] += weight;
 827             break;
 828           case FLT | FSTRING:
 829           case PLT | FSTRING:
 830             if (memcmp (iter->arg[0].c,
 831                         value_str (v, src_width), src_width) > 0)
 832               iter->dbl[0] += weight;
 833             iter->dbl[1] += weight;
 834             break;
 835           case FIN:
 836           case PIN:
 837             if (iter->arg[0].f <= v->f && v->f <= iter->arg[1].f)
 838               iter->dbl[0] += weight;
 839             iter->dbl[1] += weight;
 840             break;
 841           case FIN | FSTRING:
 842           case PIN | FSTRING:
 843             if (memcmp (iter->arg[0].c,
 844                         value_str (v, src_width), src_width) <= 0
 845                 && memcmp (iter->arg[1].c,
 846                            value_str (v, src_width), src_width) >= 0)
 847               iter->dbl[0] += weight;
 848             iter->dbl[1] += weight;
 849             break;
 850           case FOUT:
 851           case POUT:
 852             if (iter->arg[0].f > v->f || v->f > iter->arg[1].f)
 853               iter->dbl[0] += weight;
 854             iter->dbl[1] += weight;
 855             break;
 856           case FOUT | FSTRING:
 857           case POUT | FSTRING:
 858             if (memcmp (iter->arg[0].c,
 859                         value_str (v, src_width), src_width) > 0
 860                 || memcmp (iter->arg[1].c,
 861                            value_str (v, src_width), src_width) < 0)
 862               iter->dbl[0] += weight;
 863             iter->dbl[1] += weight;
 864             break;
 865           case N:
 866           case N | FSTRING:
 867             iter->dbl[0] += weight;
 868             break;
 869           case NU:
 870           case NU | FSTRING:
 871             iter->int1++;
 872             break;
 873           case FIRST:
 874             if (iter->int1 == 0)
 875               {
 876                 iter->dbl[0] = v->f;
 877                 iter->int1 = 1;
 878               }
 879             break;
 880           case FIRST | FSTRING:
 881             if (iter->int1 == 0)
 882               {
 883                 memcpy (iter->string, value_str (v, src_width), src_width);
 884                 iter->int1 = 1;
 885               }
 886             break;
 887           case LAST:
 888             iter->dbl[0] = v->f;
 889             iter->int1 = 1;
 890             break;
 891           case LAST | FSTRING:
 892             memcpy (iter->string, value_str (v, src_width), src_width);
 893             iter->int1 = 1;
 894             break;
 895           case NMISS:
 896           case NMISS | FSTRING:
 897           case NUMISS:
 898           case NUMISS | FSTRING:
 899             /* Our value is not missing or it would have been
 900                caught earlier.  Nothing to do. */
 901             break;
 902           default:
 903             NOT_REACHED ();
 904           }
 905     } else {
 906       switch (iter->function)
 907         {
 908         case N_NO_VARS:
 909           iter->dbl[0] += weight;
 910           break;
 911         case NU_NO_VARS:
 912           iter->int1++;
 913           break;
 914         default:
 915           NOT_REACHED ();
 916         }
 917     }
 918 }
 919
 920 /* Writes an aggregated record to OUTPUT. */
 921 static void
 922 dump_aggregate_info (struct agr_proc *agr, struct casewriter *output)
 923 {
 924   struct ccase *c = case_create (dict_get_proto (agr->dict));
 925
 926   {
 927     int value_idx = 0;
 928     int i;
 929
 930     for (i = 0; i < agr->break_var_cnt; i++)
 931       {
 932         const struct variable *v = agr->break_vars[i];
 933         value_copy (case_data_rw_idx (c, value_idx),
 934                     case_data (agr->break_case, v),
 935                     var_get_width (v));
 936         value_idx++;
 937       }
 938   }
 939
 940   {
 941     struct agr_var *i;
 942
 943     for (i = agr->agr_vars; i; i = i->next)
 944       {
 945         union value *v = case_data_rw (c, i->dest);
 946         int width = var_get_width (i->dest);
 947
 948         if (agr->missing == COLUMNWISE && i->saw_missing
 949             && (i->function & FUNC) != N && (i->function & FUNC) != NU
 950             && (i->function & FUNC) != NMISS && (i->function & FUNC) != NUMISS)
 951           {
 952             value_set_missing (v, width);
 953             casewriter_destroy (i->writer);
 954             continue;
 955           }
 956
 957         switch (i->function)
 958           {
 959           case SUM:
 960             v->f = i->int1 ? i->dbl[0] : SYSMIS;
 961             break;
 962           case MEAN:
 963             v->f = i->dbl[1] != 0.0 ? i->dbl[0] / i->dbl[1] : SYSMIS;
 964             break;
 965           case MEDIAN:
 966             {
 967               struct casereader *sorted_reader;
 968               struct order_stats *median = percentile_create (0.5, i->cc);
 969
 970               sorted_reader = casewriter_make_reader (i->writer);
 971
 972               order_stats_accumulate (&median, 1,
 973                                       sorted_reader,
 974                                       i->weight,
 975                                       i->subject,
 976                                       i->exclude);
 977
 978               v->f = percentile_calculate ((struct percentile *) median,
 979                                            PC_HAVERAGE);
 980
 981               statistic_destroy ((struct statistic *) median);
 982             }
 983             break;
 984           case SD:
 985             {
 986               double variance;
 987
 988               /* FIXME: we should use two passes. */
 989               moments1_calculate (i->moments, NULL, NULL, &variance,
 990                                  NULL, NULL);
 991               if (variance != SYSMIS)
 992                 v->f = sqrt (variance);
 993               else
 994                 v->f = SYSMIS;
 995             }
 996             break;
 997           case MAX:
 998           case MIN:
 999             v->f = i->int1 ? i->dbl[0] : SYSMIS;
1000             break;
1001           case MAX | FSTRING:
1002           case MIN | FSTRING:
1003             if (i->int1)
1004               memcpy (value_str_rw (v, width), i->string, width);
1005             else
1006               value_set_missing (v, width);
1007             break;
1008           case FGT:
1009           case FGT | FSTRING:
1010           case FLT:
1011           case FLT | FSTRING:
1012           case FIN:
1013           case FIN | FSTRING:
1014           case FOUT:
1015           case FOUT | FSTRING:
1016             v->f = i->dbl[1] ? i->dbl[0] / i->dbl[1] : SYSMIS;
1017             break;
1018           case PGT:
1019           case PGT | FSTRING:
1020           case PLT:
1021           case PLT | FSTRING:
1022           case PIN:
1023           case PIN | FSTRING:
1024           case POUT:
1025           case POUT | FSTRING:
1026             v->f = i->dbl[1] ? i->dbl[0] / i->dbl[1] * 100.0 : SYSMIS;
1027             break;
1028           case N:
1029           case N | FSTRING:
1030             v->f = i->dbl[0];
1031             break;
1032           case NU:
1033           case NU | FSTRING:
1034             v->f = i->int1;
1035             break;
1036           case FIRST:
1037           case LAST:
1038             v->f = i->int1 ? i->dbl[0] : SYSMIS;
1039             break;
1040           case FIRST | FSTRING:
1041           case LAST | FSTRING:
1042             if (i->int1)
1043               memcpy (value_str_rw (v, width), i->string, width);
1044             else
1045               value_set_missing (v, width);
1046             break;
1047           case N_NO_VARS:
1048             v->f = i->dbl[0];
1049             break;
1050           case NU_NO_VARS:
1051             v->f = i->int1;
1052             break;
1053           case NMISS:
1054           case NMISS | FSTRING:
1055             v->f = i->dbl[0];
1056             break;
1057           case NUMISS:
1058           case NUMISS | FSTRING:
1059             v->f = i->int1;
1060             break;
1061           default:
1062             NOT_REACHED ();
1063           }
1064       }
1065   }
1066
1067   casewriter_write (output, c);
1068 }
1069
1070 /* Resets the state for all the aggregate functions. */
1071 static void
1072 initialize_aggregate_info (struct agr_proc *agr, const struct ccase *input)
1073 {
1074   struct agr_var *iter;
1075
1076   case_unref (agr->break_case);
1077   agr->break_case = case_ref (input);
1078
1079   for (iter = agr->agr_vars; iter; iter = iter->next)
1080     {
1081       iter->saw_missing = false;
1082       iter->dbl[0] = iter->dbl[1] = iter->dbl[2] = 0.0;
1083       iter->int1 = iter->int2 = 0;
1084       switch (iter->function)
1085         {
1086         case MIN:
1087           iter->dbl[0] = DBL_MAX;
1088           break;
1089         case MIN | FSTRING:
1090           memset (iter->string, 255, var_get_width (iter->src));
1091           break;
1092         case MAX:
1093           iter->dbl[0] = -DBL_MAX;
1094           break;
1095         case MAX | FSTRING:
1096           memset (iter->string, 0, var_get_width (iter->src));
1097           break;
1098         case MEDIAN:
1099           {
1100             struct caseproto *proto;
1101             struct subcase ordering;
1102
1103             proto = caseproto_create ();
1104             proto = caseproto_add_width (proto, 0);
1105             proto = caseproto_add_width (proto, 0);
1106
1107             if ( ! iter->subject)
1108               iter->subject = var_create_internal (0);
1109
1110             if ( ! iter->weight)
1111               iter->weight = var_create_internal (1);
1112
1113             subcase_init_var (&ordering, iter->subject, SC_ASCEND);
1114             iter->writer = sort_create_writer (&ordering, proto);
1115             subcase_destroy (&ordering);
1116             caseproto_unref (proto);
1117
1118             iter->cc = 0;
1119           }
1120           break;
1121         case SD:
1122           if (iter->moments == NULL)
1123             iter->moments = moments1_create (MOMENT_VARIANCE);
1124           else
1125             moments1_clear (iter->moments);
1126           break;
1127         default:
1128           break;
1129         }
1130     }
1131 }