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