Fixed return type of design_matrix_get_case_count.
[pspp-builds.git] / src / language / stats / aggregate.c
1 /* PSPP - a program for statistical analysis.
2    Copyright (C) 1997-9, 2000, 2006, 2008 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-ordering.h>
23 #include <data/case.h>
24 #include <data/casegrouper.h>
25 #include <data/casereader.h>
26 #include <data/casewriter.h>
27 #include <data/dictionary.h>
28 #include <data/file-handle-def.h>
29 #include <data/format.h>
30 #include <data/procedure.h>
31 #include <data/settings.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 case_ordering *sort;         /* Sort criteria (break variable). */
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   case_nullify (&agr.break_case);
191
192   agr.dict = dict_create ();
193   agr.src_dict = dict;
194   dict_set_label (agr.dict, dict_get_label (dict));
195   dict_set_documents (agr.dict, dict_get_documents (dict));
196
197   /* OUTFILE subcommand must be first. */
198   if (!lex_force_match_id (lexer, "OUTFILE"))
199     goto error;
200   lex_match (lexer, '=');
201   if (!lex_match (lexer, '*'))
202     {
203       out_file = fh_parse (lexer, FH_REF_FILE | FH_REF_SCRATCH);
204       if (out_file == NULL)
205         goto error;
206     }
207
208   /* Read most of the subcommands. */
209   for (;;)
210     {
211       lex_match (lexer, '/');
212
213       if (lex_match_id (lexer, "MISSING"))
214         {
215           lex_match (lexer, '=');
216           if (!lex_match_id (lexer, "COLUMNWISE"))
217             {
218               lex_error (lexer, _("while expecting COLUMNWISE"));
219               goto error;
220             }
221           agr.missing = COLUMNWISE;
222         }
223       else if (lex_match_id (lexer, "DOCUMENT"))
224         copy_documents = true;
225       else if (lex_match_id (lexer, "PRESORTED"))
226         presorted = true;
227       else if (lex_match_id (lexer, "BREAK"))
228         {
229           int i;
230
231           lex_match (lexer, '=');
232           agr.sort = parse_case_ordering (lexer, dict,
233
234                                           &saw_direction);
235           if (agr.sort == NULL)
236             goto error;
237           case_ordering_get_vars (agr.sort,
238                                   &agr.break_vars, &agr.break_var_cnt);
239
240           for (i = 0; i < agr.break_var_cnt; i++)
241             dict_clone_var_assert (agr.dict, agr.break_vars[i],
242                                    var_get_name (agr.break_vars[i]));
243
244           /* BREAK must follow the options. */
245           break;
246         }
247       else
248         {
249           lex_error (lexer, _("expecting BREAK"));
250           goto error;
251         }
252     }
253   if (presorted && saw_direction)
254     msg (SW, _("When PRESORTED is specified, specifying sorting directions "
255                "with (A) or (D) has no effect.  Output data will be sorted "
256                "the same way as the input data."));
257
258   /* Read in the aggregate functions. */
259   lex_match (lexer, '/');
260   if (!parse_aggregate_functions (lexer, dict, &agr))
261     goto error;
262
263   /* Delete documents. */
264   if (!copy_documents)
265     dict_clear_documents (agr.dict);
266
267   /* Cancel SPLIT FILE. */
268   dict_set_split_vars (agr.dict, NULL, 0);
269
270   /* Initialize. */
271   agr.case_cnt = 0;
272
273   if (out_file == NULL)
274     {
275       /* The active file will be replaced by the aggregated data,
276          so TEMPORARY is moot. */
277       proc_cancel_temporary_transformations (ds);
278       proc_discard_output (ds);
279       output = autopaging_writer_create (dict_get_next_value_idx (agr.dict));
280     }
281   else
282     {
283       output = any_writer_open (out_file, agr.dict);
284       if (output == NULL)
285         goto error;
286     }
287
288   input = proc_open (ds);
289   if (agr.sort != NULL && !presorted)
290     {
291       input = sort_execute (input, agr.sort);
292       agr.sort = NULL;
293     }
294
295   for (grouper = casegrouper_create_vars (input, agr.break_vars,
296                                           agr.break_var_cnt);
297        casegrouper_get_next_group (grouper, &group);
298        casereader_destroy (group))
299     {
300       struct ccase c;
301
302       if (!casereader_peek (group, 0, &c))
303         {
304           casereader_destroy (group);
305           continue;
306         }
307       initialize_aggregate_info (&agr, &c);
308       case_destroy (&c);
309
310       for (; casereader_read (group, &c); case_destroy (&c))
311         accumulate_aggregate_info (&agr, &c);
312       dump_aggregate_info (&agr, output);
313     }
314   if (!casegrouper_destroy (grouper))
315     goto error;
316
317   if (!proc_commit (ds))
318     {
319       input = NULL;
320       goto error;
321     }
322   input = NULL;
323
324   if (out_file == NULL)
325     {
326       struct casereader *next_input = casewriter_make_reader (output);
327       if (next_input == NULL)
328         goto error;
329
330       proc_set_active_file (ds, next_input, agr.dict);
331       agr.dict = NULL;
332     }
333   else
334     {
335       ok = casewriter_destroy (output);
336       output = NULL;
337       if (!ok)
338         goto error;
339     }
340
341   agr_destroy (&agr);
342   fh_unref (out_file);
343   return CMD_SUCCESS;
344
345 error:
346   if (input != NULL)
347     proc_commit (ds);
348   casewriter_destroy (output);
349   agr_destroy (&agr);
350   fh_unref (out_file);
351   return CMD_CASCADING_FAILURE;
352 }
353
354 /* Parse all the aggregate functions. */
355 static bool
356 parse_aggregate_functions (struct lexer *lexer, const struct dictionary *dict,
357                            struct agr_proc *agr)
358 {
359   struct agr_var *tail; /* Tail of linked list starting at agr->vars. */
360
361   /* Parse everything. */
362   tail = NULL;
363   for (;;)
364     {
365       char **dest;
366       char **dest_label;
367       size_t n_dest;
368       struct string function_name;
369
370       enum mv_class exclude;
371       const struct agr_func *function;
372       int func_index;
373
374       union agr_argument arg[2];
375
376       const struct variable **src;
377       size_t n_src;
378
379       size_t i;
380
381       dest = NULL;
382       dest_label = NULL;
383       n_dest = 0;
384       src = NULL;
385       function = NULL;
386       n_src = 0;
387       arg[0].c = NULL;
388       arg[1].c = NULL;
389       ds_init_empty (&function_name);
390
391       /* Parse the list of target variables. */
392       while (!lex_match (lexer, '='))
393         {
394           size_t n_dest_prev = n_dest;
395
396           if (!parse_DATA_LIST_vars (lexer, &dest, &n_dest,
397                                      PV_APPEND | PV_SINGLE | PV_NO_SCRATCH))
398             goto error;
399
400           /* Assign empty labels. */
401           {
402             int j;
403
404             dest_label = xnrealloc (dest_label, n_dest, sizeof *dest_label);
405             for (j = n_dest_prev; j < n_dest; j++)
406               dest_label[j] = NULL;
407           }
408
409
410
411           if (lex_token (lexer) == T_STRING)
412             {
413               struct string label;
414               ds_init_string (&label, lex_tokstr (lexer));
415
416               ds_truncate (&label, 255);
417               dest_label[n_dest - 1] = ds_xstrdup (&label);
418               lex_get (lexer);
419               ds_destroy (&label);
420             }
421         }
422
423       /* Get the name of the aggregation function. */
424       if (lex_token (lexer) != T_ID)
425         {
426           lex_error (lexer, _("expecting aggregation function"));
427           goto error;
428         }
429
430       exclude = MV_ANY;
431
432       ds_assign_string (&function_name, lex_tokstr (lexer));
433
434       ds_chomp (&function_name, '.');
435
436       if (lex_tokid(lexer)[strlen (lex_tokid (lexer)) - 1] == '.')
437         exclude = MV_SYSTEM;
438
439       for (function = agr_func_tab; function->name; function++)
440         if (!strcasecmp (function->name, ds_cstr (&function_name)))
441           break;
442       if (NULL == function->name)
443         {
444           msg (SE, _("Unknown aggregation function %s."),
445                ds_cstr (&function_name));
446           goto error;
447         }
448       ds_destroy (&function_name);
449       func_index = function - agr_func_tab;
450       lex_get (lexer);
451
452       /* Check for leading lparen. */
453       if (!lex_match (lexer, '('))
454         {
455           if (func_index == N)
456             func_index = N_NO_VARS;
457           else if (func_index == NU)
458             func_index = NU_NO_VARS;
459           else
460             {
461               lex_error (lexer, _("expecting `('"));
462               goto error;
463             }
464         }
465       else
466         {
467           /* Parse list of source variables. */
468           {
469             int pv_opts = PV_NO_SCRATCH;
470
471             if (func_index == SUM || func_index == MEAN || func_index == SD)
472               pv_opts |= PV_NUMERIC;
473             else if (function->n_args)
474               pv_opts |= PV_SAME_TYPE;
475
476             if (!parse_variables_const (lexer, dict, &src, &n_src, pv_opts))
477               goto error;
478           }
479
480           /* Parse function arguments, for those functions that
481              require arguments. */
482           if (function->n_args != 0)
483             for (i = 0; i < function->n_args; i++)
484               {
485                 int type;
486
487                 lex_match (lexer, ',');
488                 if (lex_token (lexer) == T_STRING)
489                   {
490                     arg[i].c = ds_xstrdup (lex_tokstr (lexer));
491                     type = VAL_STRING;
492                   }
493                 else if (lex_is_number (lexer))
494                   {
495                     arg[i].f = lex_tokval (lexer);
496                     type = VAL_NUMERIC;
497                   }
498                 else
499                   {
500                     msg (SE, _("Missing argument %zu to %s."),
501                          i + 1, function->name);
502                     goto error;
503                   }
504
505                 lex_get (lexer);
506
507                 if (type != var_get_type (src[0]))
508                   {
509                     msg (SE, _("Arguments to %s must be of same type as "
510                                "source variables."),
511                          function->name);
512                     goto error;
513                   }
514               }
515
516           /* Trailing rparen. */
517           if (!lex_match (lexer, ')'))
518             {
519               lex_error (lexer, _("expecting `)'"));
520               goto error;
521             }
522
523           /* Now check that the number of source variables match
524              the number of target variables.  If we check earlier
525              than this, the user can get very misleading error
526              message, i.e. `AGGREGATE x=SUM(y t).' will get this
527              error message when a proper message would be more
528              like `unknown variable t'. */
529           if (n_src != n_dest)
530             {
531               msg (SE, _("Number of source variables (%zu) does not match "
532                          "number of target variables (%zu)."),
533                     n_src, n_dest);
534               goto error;
535             }
536
537           if ((func_index == PIN || func_index == POUT
538               || func_index == FIN || func_index == FOUT)
539               && (var_is_numeric (src[0])
540                   ? arg[0].f > arg[1].f
541                   : str_compare_rpad (arg[0].c, arg[1].c) > 0))
542             {
543               union agr_argument t = arg[0];
544               arg[0] = arg[1];
545               arg[1] = t;
546
547               msg (SW, _("The value arguments passed to the %s function "
548                          "are out-of-order.  They will be treated as if "
549                          "they had been specified in the correct order."),
550                    function->name);
551             }
552         }
553
554       /* Finally add these to the linked list of aggregation
555          variables. */
556       for (i = 0; i < n_dest; i++)
557         {
558           struct agr_var *v = xzalloc (sizeof *v);
559
560           /* Add variable to chain. */
561           if (agr->agr_vars != NULL)
562             tail->next = v;
563           else
564             agr->agr_vars = v;
565           tail = v;
566           tail->next = NULL;
567           v->moments = NULL;
568
569           /* Create the target variable in the aggregate
570              dictionary. */
571           {
572             struct variable *destvar;
573
574             v->function = func_index;
575
576             if (src)
577               {
578                 v->src = src[i];
579
580                 if (var_is_alpha (src[i]))
581                   {
582                     v->function |= FSTRING;
583                     v->string = xmalloc (var_get_width (src[i]));
584                   }
585
586                 if (function->alpha_type == VAL_STRING)
587                   destvar = dict_clone_var (agr->dict, v->src, dest[i]);
588                 else
589                   {
590                     assert (var_is_numeric (v->src)
591                             || function->alpha_type == VAL_NUMERIC);
592                     destvar = dict_create_var (agr->dict, dest[i], 0);
593                     if (destvar != NULL)
594                       {
595                         struct fmt_spec f;
596                         if ((func_index == N || func_index == NMISS)
597                             && dict_get_weight (dict) != NULL)
598                           f = fmt_for_output (FMT_F, 8, 2);
599                         else
600                           f = function->format;
601                         var_set_both_formats (destvar, &f);
602                       }
603                   }
604               } else {
605                 struct fmt_spec f;
606                 v->src = NULL;
607                 destvar = dict_create_var (agr->dict, dest[i], 0);
608                 if (func_index == N_NO_VARS && dict_get_weight (dict) != NULL)
609                   f = fmt_for_output (FMT_F, 8, 2);
610                 else
611                   f = function->format;
612                 var_set_both_formats (destvar, &f);
613               }
614
615             if (!destvar)
616               {
617                 msg (SE, _("Variable name %s is not unique within the "
618                            "aggregate file dictionary, which contains "
619                            "the aggregate variables and the break "
620                            "variables."),
621                      dest[i]);
622                 goto error;
623               }
624
625             free (dest[i]);
626             if (dest_label[i])
627               var_set_label (destvar, dest_label[i]);
628
629             v->dest = destvar;
630           }
631
632           v->exclude = exclude;
633
634           if (v->src != NULL)
635             {
636               int j;
637
638               if (var_is_numeric (v->src))
639                 for (j = 0; j < function->n_args; j++)
640                   v->arg[j].f = arg[j].f;
641               else
642                 for (j = 0; j < function->n_args; j++)
643                   v->arg[j].c = xstrdup (arg[j].c);
644             }
645         }
646
647       if (src != NULL && var_is_alpha (src[0]))
648         for (i = 0; i < function->n_args; i++)
649           {
650             free (arg[i].c);
651             arg[i].c = NULL;
652           }
653
654       free (src);
655       free (dest);
656       free (dest_label);
657
658       if (!lex_match (lexer, '/'))
659         {
660           if (lex_token (lexer) == '.')
661             return true;
662
663           lex_error (lexer, "expecting end of command");
664           return false;
665         }
666       continue;
667
668     error:
669       ds_destroy (&function_name);
670       for (i = 0; i < n_dest; i++)
671         {
672           free (dest[i]);
673           free (dest_label[i]);
674         }
675       free (dest);
676       free (dest_label);
677       free (arg[0].c);
678       free (arg[1].c);
679       if (src && n_src && var_is_alpha (src[0]))
680         for (i = 0; i < function->n_args; i++)
681           {
682             free (arg[i].c);
683             arg[i].c = NULL;
684           }
685       free (src);
686
687       return false;
688     }
689 }
690
691 /* Destroys AGR. */
692 static void
693 agr_destroy (struct agr_proc *agr)
694 {
695   struct agr_var *iter, *next;
696
697   case_ordering_destroy (agr->sort);
698   free (agr->break_vars);
699   case_destroy (&agr->break_case);
700   for (iter = agr->agr_vars; iter; iter = next)
701     {
702       next = iter->next;
703
704       if (iter->function & FSTRING)
705         {
706           size_t n_args;
707           size_t i;
708
709           n_args = agr_func_tab[iter->function & FUNC].n_args;
710           for (i = 0; i < n_args; i++)
711             free (iter->arg[i].c);
712           free (iter->string);
713         }
714       else if (iter->function == SD)
715         moments1_destroy (iter->moments);
716
717       var_destroy (iter->subject);
718       var_destroy (iter->weight);
719
720       free (iter);
721     }
722   if (agr->dict != NULL)
723     dict_destroy (agr->dict);
724 }
725 \f
726 /* Execution. */
727
728 /* Accumulates aggregation data from the case INPUT. */
729 static void
730 accumulate_aggregate_info (struct agr_proc *agr, const struct ccase *input)
731 {
732   struct agr_var *iter;
733   double weight;
734   bool bad_warn = true;
735
736   weight = dict_get_case_weight (agr->src_dict, input, &bad_warn);
737
738   for (iter = agr->agr_vars; iter; iter = iter->next)
739     if (iter->src)
740       {
741         const union value *v = case_data (input, iter->src);
742         int src_width = var_get_width (iter->src);
743
744         if (var_is_value_missing (iter->src, v, iter->exclude))
745           {
746             switch (iter->function)
747               {
748               case NMISS:
749               case NMISS | FSTRING:
750                 iter->dbl[0] += weight;
751                 break;
752               case NUMISS:
753               case NUMISS | FSTRING:
754                 iter->int1++;
755                 break;
756               }
757             iter->saw_missing = true;
758             continue;
759           }
760
761         /* This is horrible.  There are too many possibilities. */
762         switch (iter->function)
763           {
764           case SUM:
765             iter->dbl[0] += v->f * weight;
766             iter->int1 = 1;
767             break;
768           case MEAN:
769             iter->dbl[0] += v->f * weight;
770             iter->dbl[1] += weight;
771             break;
772           case MEDIAN:
773             {
774               double wv ;
775               struct ccase cout;
776               case_create (&cout, 2);
777
778               case_data_rw (&cout, iter->subject)->f =
779                 case_data (input, iter->src)->f;
780
781               wv = dict_get_case_weight (agr->src_dict, input, NULL);
782
783               case_data_rw (&cout, iter->weight)->f = wv;
784
785               iter->cc += wv;
786
787               casewriter_write (iter->writer, &cout);
788               case_destroy (&cout);
789             }
790             break;
791           case SD:
792             moments1_add (iter->moments, v->f, weight);
793             break;
794           case MAX:
795             iter->dbl[0] = MAX (iter->dbl[0], v->f);
796             iter->int1 = 1;
797             break;
798           case MAX | FSTRING:
799             if (memcmp (iter->string, v->s, src_width) < 0)
800               memcpy (iter->string, v->s, src_width);
801             iter->int1 = 1;
802             break;
803           case MIN:
804             iter->dbl[0] = MIN (iter->dbl[0], v->f);
805             iter->int1 = 1;
806             break;
807           case MIN | FSTRING:
808             if (memcmp (iter->string, v->s, src_width) > 0)
809               memcpy (iter->string, v->s, src_width);
810             iter->int1 = 1;
811             break;
812           case FGT:
813           case PGT:
814             if (v->f > iter->arg[0].f)
815               iter->dbl[0] += weight;
816             iter->dbl[1] += weight;
817             break;
818           case FGT | FSTRING:
819           case PGT | FSTRING:
820             if (memcmp (iter->arg[0].c, v->s, src_width) < 0)
821               iter->dbl[0] += weight;
822             iter->dbl[1] += weight;
823             break;
824           case FLT:
825           case PLT:
826             if (v->f < iter->arg[0].f)
827               iter->dbl[0] += weight;
828             iter->dbl[1] += weight;
829             break;
830           case FLT | FSTRING:
831           case PLT | FSTRING:
832             if (memcmp (iter->arg[0].c, v->s, src_width) > 0)
833               iter->dbl[0] += weight;
834             iter->dbl[1] += weight;
835             break;
836           case FIN:
837           case PIN:
838             if (iter->arg[0].f <= v->f && v->f <= iter->arg[1].f)
839               iter->dbl[0] += weight;
840             iter->dbl[1] += weight;
841             break;
842           case FIN | FSTRING:
843           case PIN | FSTRING:
844             if (memcmp (iter->arg[0].c, v->s, src_width) <= 0
845                 && memcmp (iter->arg[1].c, v->s, src_width) >= 0)
846               iter->dbl[0] += weight;
847             iter->dbl[1] += weight;
848             break;
849           case FOUT:
850           case POUT:
851             if (iter->arg[0].f > v->f || v->f > iter->arg[1].f)
852               iter->dbl[0] += weight;
853             iter->dbl[1] += weight;
854             break;
855           case FOUT | FSTRING:
856           case POUT | FSTRING:
857             if (memcmp (iter->arg[0].c, v->s, src_width) > 0
858                 || memcmp (iter->arg[1].c, v->s, src_width) < 0)
859               iter->dbl[0] += weight;
860             iter->dbl[1] += weight;
861             break;
862           case N:
863           case N | FSTRING:
864             iter->dbl[0] += weight;
865             break;
866           case NU:
867           case NU | FSTRING:
868             iter->int1++;
869             break;
870           case FIRST:
871             if (iter->int1 == 0)
872               {
873                 iter->dbl[0] = v->f;
874                 iter->int1 = 1;
875               }
876             break;
877           case FIRST | FSTRING:
878             if (iter->int1 == 0)
879               {
880                 memcpy (iter->string, v->s, src_width);
881                 iter->int1 = 1;
882               }
883             break;
884           case LAST:
885             iter->dbl[0] = v->f;
886             iter->int1 = 1;
887             break;
888           case LAST | FSTRING:
889             memcpy (iter->string, v->s, src_width);
890             iter->int1 = 1;
891             break;
892           case NMISS:
893           case NMISS | FSTRING:
894           case NUMISS:
895           case NUMISS | FSTRING:
896             /* Our value is not missing or it would have been
897                caught earlier.  Nothing to do. */
898             break;
899           default:
900             NOT_REACHED ();
901           }
902     } else {
903       switch (iter->function)
904         {
905         case N_NO_VARS:
906           iter->dbl[0] += weight;
907           break;
908         case NU_NO_VARS:
909           iter->int1++;
910           break;
911         default:
912           NOT_REACHED ();
913         }
914     }
915 }
916
917 /* Writes an aggregated record to OUTPUT. */
918 static void
919 dump_aggregate_info (struct agr_proc *agr, struct casewriter *output)
920 {
921   struct ccase c;
922
923   case_create (&c, dict_get_next_value_idx (agr->dict));
924
925   {
926     int value_idx = 0;
927     int i;
928
929     for (i = 0; i < agr->break_var_cnt; i++)
930       {
931         const struct variable *v = agr->break_vars[i];
932         size_t value_cnt = var_get_value_cnt (v);
933         memcpy (case_data_rw_idx (&c, value_idx),
934                 case_data (&agr->break_case, v),
935                 sizeof (union value) * value_cnt);
936         value_idx += value_cnt;
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
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             if (var_is_alpha (i->dest))
953               memset (v->s, ' ', var_get_width (i->dest));
954             else
955               v->f = SYSMIS;
956
957             casewriter_destroy (i->writer);
958
959             continue;
960           }
961
962         switch (i->function)
963           {
964           case SUM:
965             v->f = i->int1 ? i->dbl[0] : SYSMIS;
966             break;
967           case MEAN:
968             v->f = i->dbl[1] != 0.0 ? i->dbl[0] / i->dbl[1] : SYSMIS;
969             break;
970           case MEDIAN:
971             {
972               struct casereader *sorted_reader;
973               struct order_stats *median = percentile_create (0.5, i->cc);
974
975               sorted_reader = casewriter_make_reader (i->writer);
976
977               order_stats_accumulate (&median, 1,
978                                       sorted_reader,
979                                       i->weight,
980                                       i->subject,
981                                       i->exclude);
982
983               v->f = percentile_calculate ((struct percentile *) median,
984                                            PC_HAVERAGE);
985
986               statistic_destroy ((struct statistic *) median);
987             }
988             break;
989           case SD:
990             {
991               double variance;
992
993               /* FIXME: we should use two passes. */
994               moments1_calculate (i->moments, NULL, NULL, &variance,
995                                  NULL, NULL);
996               if (variance != SYSMIS)
997                 v->f = sqrt (variance);
998               else
999                 v->f = SYSMIS;
1000             }
1001             break;
1002           case MAX:
1003           case MIN:
1004             v->f = i->int1 ? i->dbl[0] : SYSMIS;
1005             break;
1006           case MAX | FSTRING:
1007           case MIN | FSTRING:
1008             if (i->int1)
1009               memcpy (v->s, i->string, var_get_width (i->dest));
1010             else
1011               memset (v->s, ' ', var_get_width (i->dest));
1012             break;
1013           case FGT:
1014           case FGT | FSTRING:
1015           case FLT:
1016           case FLT | FSTRING:
1017           case FIN:
1018           case FIN | FSTRING:
1019           case FOUT:
1020           case FOUT | FSTRING:
1021             v->f = i->dbl[1] ? i->dbl[0] / i->dbl[1] : SYSMIS;
1022             break;
1023           case PGT:
1024           case PGT | FSTRING:
1025           case PLT:
1026           case PLT | FSTRING:
1027           case PIN:
1028           case PIN | FSTRING:
1029           case POUT:
1030           case POUT | FSTRING:
1031             v->f = i->dbl[1] ? i->dbl[0] / i->dbl[1] * 100.0 : SYSMIS;
1032             break;
1033           case N:
1034           case N | FSTRING:
1035             v->f = i->dbl[0];
1036             break;
1037           case NU:
1038           case NU | FSTRING:
1039             v->f = i->int1;
1040             break;
1041           case FIRST:
1042           case LAST:
1043             v->f = i->int1 ? i->dbl[0] : SYSMIS;
1044             break;
1045           case FIRST | FSTRING:
1046           case LAST | FSTRING:
1047             if (i->int1)
1048               memcpy (v->s, i->string, var_get_width (i->dest));
1049             else
1050               memset (v->s, ' ', var_get_width (i->dest));
1051             break;
1052           case N_NO_VARS:
1053             v->f = i->dbl[0];
1054             break;
1055           case NU_NO_VARS:
1056             v->f = i->int1;
1057             break;
1058           case NMISS:
1059           case NMISS | FSTRING:
1060             v->f = i->dbl[0];
1061             break;
1062           case NUMISS:
1063           case NUMISS | FSTRING:
1064             v->f = i->int1;
1065             break;
1066           default:
1067             NOT_REACHED ();
1068           }
1069       }
1070   }
1071
1072   casewriter_write (output, &c);
1073 }
1074
1075 /* Resets the state for all the aggregate functions. */
1076 static void
1077 initialize_aggregate_info (struct agr_proc *agr, const struct ccase *input)
1078 {
1079   struct agr_var *iter;
1080
1081   case_destroy (&agr->break_case);
1082   case_clone (&agr->break_case, input);
1083
1084   for (iter = agr->agr_vars; iter; iter = iter->next)
1085     {
1086       iter->saw_missing = false;
1087       iter->dbl[0] = iter->dbl[1] = iter->dbl[2] = 0.0;
1088       iter->int1 = iter->int2 = 0;
1089       switch (iter->function)
1090         {
1091         case MIN:
1092           iter->dbl[0] = DBL_MAX;
1093           break;
1094         case MIN | FSTRING:
1095           memset (iter->string, 255, var_get_width (iter->src));
1096           break;
1097         case MAX:
1098           iter->dbl[0] = -DBL_MAX;
1099           break;
1100         case MAX | FSTRING:
1101           memset (iter->string, 0, var_get_width (iter->src));
1102           break;
1103         case MEDIAN:
1104           {
1105             struct case_ordering *ordering = case_ordering_create ();
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             case_ordering_add_var (ordering, iter->subject, SRT_ASCEND);
1114
1115             iter->writer = sort_create_writer (ordering, 2);
1116             iter->cc = 0;
1117           }
1118           break;
1119         case SD:
1120           if (iter->moments == NULL)
1121             iter->moments = moments1_create (MOMENT_VARIANCE);
1122           else
1123             moments1_clear (iter->moments);
1124           break;
1125         default:
1126           break;
1127         }
1128     }
1129 }