Examine vs. Boxplots: Avoid labels overlapping one another
[pspp] / src / language / stats / examine.c
1 /*
2   PSPP - a program for statistical analysis.
3   Copyright (C) 2012, 2013  Free Software Foundation, Inc.
4   
5   This program is free software: you can redistribute it and/or modify
6   it under the terms of the GNU General Public License as published by
7   the Free Software Foundation, either version 3 of the License, or
8   (at your option) any later version.
9
10   This program is distributed in the hope that it will be useful,
11   but WITHOUT ANY WARRANTY; without even the implied warranty of
12   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
13   GNU General Public License for more details.
14   
15   You should have received a copy of the GNU General Public License
16   along with this program.  If not, see <http://www.gnu.org/licenses/>.
17 */
18
19 #include <config.h>
20
21 #include <math.h>
22 #include <gsl/gsl_cdf.h>
23
24 #include "libpspp/assertion.h"
25 #include "libpspp/message.h"
26 #include "libpspp/pool.h"
27
28
29 #include "data/dataset.h"
30 #include "data/dictionary.h"
31 #include "data/casegrouper.h"
32 #include "data/casereader.h"
33 #include "data/casewriter.h"
34 #include "data/caseproto.h"
35 #include "data/subcase.h"
36
37
38 #include "data/format.h"
39
40 #include "math/interaction.h"
41 #include "math/box-whisker.h"
42 #include "math/categoricals.h"
43 #include "math/chart-geometry.h"
44 #include "math/histogram.h"
45 #include "math/moments.h"
46 #include "math/np.h"
47 #include "math/sort.h"
48 #include "math/order-stats.h"
49 #include "math/percentiles.h"
50 #include "math/tukey-hinges.h"
51 #include "math/trimmed-mean.h"
52
53 #include "output/charts/boxplot.h"
54 #include "output/charts/np-plot.h"
55 #include "output/charts/spreadlevel-plot.h"
56 #include "output/charts/plot-hist.h"
57
58 #include "language/command.h"
59 #include "language/lexer/lexer.h"
60 #include "language/lexer/value-parser.h"
61 #include "language/lexer/variable-parser.h"
62
63 #include "output/tab.h"
64
65 #include "gettext.h"
66 #define _(msgid) gettext (msgid)
67 #define N_(msgid) msgid
68
69 static void 
70 append_value_name (const struct variable *var, const union value *val, struct string *str)
71 {
72   var_append_value_name (var, val, str);
73   if ( var_is_value_missing (var, val, MV_ANY))
74     ds_put_cstr (str, _(" (missing)"));
75 }
76
77 enum bp_mode
78   {
79     BP_GROUPS,
80     BP_VARIABLES
81   };
82
83
84 /* Indices for the ex_proto member (below) */
85 enum
86   {
87     EX_VAL,  /* value */
88     EX_ID,   /* identity */
89     EX_WT    /* weight */
90   };
91
92
93 struct examine
94 {
95   struct pool *pool;
96
97   /* A caseproto used to contain the data subsets under examination,
98      see (enum above)   */
99   struct caseproto *ex_proto;
100
101   size_t n_dep_vars;
102   const struct variable **dep_vars;
103
104   size_t n_iacts;
105   struct interaction **iacts;
106
107   enum mv_class dep_excl;
108   enum mv_class fctr_excl;
109
110   const struct dictionary *dict;
111
112   struct categoricals *cats;
113
114   /* how many extremities to display */
115   int disp_extremes;
116   int calc_extremes;
117   bool descriptives;
118
119   double conf;
120
121   bool missing_pw;
122
123   /* The case index of the ID value (or -1) if not applicable */
124   size_t id_idx;
125   int id_width;
126
127   enum pc_alg pc_alg;
128   double *ptiles;
129   size_t n_percentiles;
130   
131   bool npplot;
132   bool histogram;
133   bool boxplot;
134   bool spreadlevel;
135   int sl_power;
136
137   enum bp_mode boxplot_mode;
138
139   const struct variable *id_var;
140
141   const struct variable *wv;
142 };
143
144 struct extremity
145 {
146   /* The value of this extremity */
147   double val;
148
149   /* Either the casenumber or the value of the variable specified
150      by the /ID subcommand which corresponds to this extremity */
151   union value identity;
152 };
153
154 struct exploratory_stats
155 {
156   double missing;
157   double non_missing;
158
159   struct moments *mom;
160
161   /* Most operations need a sorted reader/writer */
162   struct casewriter *sorted_writer;
163   struct casereader *sorted_reader;
164
165   struct extremity *minima;
166   struct extremity *maxima;
167
168   /* 
169      Minimum should alway equal mimima[0].val.
170      Likewise, maximum should alway equal maxima[0].val.
171      This redundancy exists as an optimisation effort.
172      Some statistics (eg histogram) require early calculation
173      of the min and max
174   */
175   double minimum;
176   double maximum;
177
178   struct trimmed_mean *trimmed_mean;
179   struct percentile *quartiles[3];
180   struct percentile **percentiles;
181
182   struct tukey_hinges *hinges;
183
184   /* The data for the NP Plots */
185   struct np *np;
186
187   struct histogram *histogram;
188
189   /* The data for the box plots */
190   struct box_whisker *box_whisker;
191
192   /* Total weight */
193   double cc;
194
195   /* The minimum weight */
196   double cmin;
197 };
198
199
200 /* Returns an array of (iact->n_vars) pointers to union value initialised to NULL.
201    The caller must free this array when no longer required. */
202 static const union value **
203 previous_value_alloc (const struct interaction *iact)
204 {
205   int ivar_idx;
206
207   const union value **prev_val = xcalloc (iact->n_vars, sizeof (*prev_val));
208
209   for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
210     prev_val[ivar_idx] = NULL;
211
212   return prev_val;
213 }
214
215 /* Set the contents of PREV_VAL to the values of C indexed by the variables of IACT */
216 static int
217 previous_value_record (const struct interaction *iact, const struct ccase *c, const union value **prev_val)
218 {
219   int ivar_idx;
220   int diff_idx = -1;
221
222   for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
223     {
224       const struct variable *ivar = iact->vars[ivar_idx];
225       const int width = var_get_width (ivar);
226       const union value *val = case_data (c, ivar);
227                   
228       if (prev_val[ivar_idx])
229         if (! value_equal (prev_val[ivar_idx], val, width))
230           {
231             diff_idx = ivar_idx;
232             break;
233           }
234     }
235
236   for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
237     {
238       const struct variable *ivar = iact->vars[ivar_idx];
239       const union value *val = case_data (c, ivar);
240       
241       prev_val[ivar_idx] = val;
242     }
243   return diff_idx;
244 }
245
246
247 static void
248 show_boxplot_grouped (const struct examine *cmd, int iact_idx)
249 {
250   int v;
251
252   const struct interaction *iact = cmd->iacts[iact_idx];
253   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
254
255   for (v = 0; v < cmd->n_dep_vars; ++v)
256     {
257       double y_min = DBL_MAX;
258       double y_max = -DBL_MAX;
259       int grp;
260       struct boxplot *boxplot;
261       struct string title;
262       ds_init_empty (&title);
263
264       if (iact->n_vars > 0)
265         {
266           struct string istr;
267           ds_init_empty (&istr);
268           interaction_to_string (iact, &istr);
269           ds_put_format (&title, _("Boxplot of %s vs. %s"),
270                          var_to_string (cmd->dep_vars[v]),
271                          ds_cstr (&istr));
272           ds_destroy (&istr);
273         }
274       else
275         ds_put_format (&title, _("Boxplot of %s"), var_to_string (cmd->dep_vars[v]));
276       
277       for (grp = 0; grp < n_cats; ++grp)
278         {
279           const struct exploratory_stats *es =
280             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
281
282           if ( y_min > es[v].minimum)
283             y_min = es[v].minimum;
284
285           if ( y_max < es[v].maximum)
286             y_max = es[v].maximum;
287         }
288       
289       boxplot = boxplot_create (y_min, y_max, ds_cstr (&title));
290
291       ds_destroy (&title);
292
293       for (grp = 0; grp < n_cats; ++grp)
294         {
295           int ivar_idx;
296           struct string label;
297
298           const struct ccase *c =
299             categoricals_get_case_by_category_real (cmd->cats,  iact_idx, grp);
300
301           const struct exploratory_stats *es =
302             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
303
304           ds_init_empty (&label);
305           for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
306             {
307               struct string l;
308               const struct variable *ivar = iact->vars[ivar_idx];
309               const union value *val = case_data (c, ivar);
310               ds_init_empty (&l);
311
312               append_value_name (ivar, val, &l);
313               ds_ltrim (&l, ss_cstr (" "));
314
315               ds_put_substring (&label, l.ss);
316               if (ivar_idx < iact->n_vars - 1)
317                 ds_put_cstr (&label, "; ");
318
319               ds_destroy (&l);
320             }
321
322           boxplot_add_box (boxplot, es[v].box_whisker, ds_cstr (&label));
323
324           ds_destroy (&label);
325         }
326       
327       boxplot_submit (boxplot);
328     }
329 }
330
331 static void
332 show_boxplot_variabled (const struct examine *cmd, int iact_idx)
333 {
334   int grp;
335   const struct interaction *iact = cmd->iacts[iact_idx];
336   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
337
338   for (grp = 0; grp < n_cats; ++grp)
339     {
340       struct boxplot *boxplot;
341       int v;
342       double y_min = DBL_MAX;
343       double y_max = -DBL_MAX;
344
345       const struct ccase *c =
346         categoricals_get_case_by_category_real (cmd->cats,  iact_idx, grp);
347
348       struct string title;
349       ds_init_empty (&title);
350
351       for (v = 0; v < cmd->n_dep_vars; ++v)
352         {
353           const struct exploratory_stats *es =
354             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
355
356           if ( y_min > es[v].minimum)
357             y_min = es[v].minimum;
358
359           if ( y_max < es[v].maximum)
360             y_max = es[v].maximum;
361         }
362
363       if ( iact->n_vars == 0)
364         ds_put_format (&title, _("Boxplot"));
365       else
366         {
367           int ivar_idx;
368           struct string label;
369           ds_init_empty (&label);
370           for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
371             {
372               const struct variable *ivar = iact->vars[ivar_idx];
373               const union value *val = case_data (c, ivar);
374               
375               ds_put_cstr (&label, var_to_string (ivar));
376               ds_put_cstr (&label, " = ");
377               append_value_name (ivar, val, &label);
378               ds_put_cstr (&label, "; ");
379             }
380
381           ds_put_format (&title, _("Boxplot of %s"),
382                          ds_cstr (&label));
383
384           ds_destroy (&label);
385         }
386
387       boxplot = boxplot_create (y_min, y_max, ds_cstr (&title));
388
389       ds_destroy (&title);
390
391       for (v = 0; v < cmd->n_dep_vars; ++v)
392         {
393           const struct exploratory_stats *es =
394             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
395
396           boxplot_add_box (boxplot, es[v].box_whisker, 
397                            var_to_string (cmd->dep_vars[v]));
398         }
399
400       boxplot_submit (boxplot);
401     }
402 }
403
404
405 static void
406 show_npplot (const struct examine *cmd, int iact_idx)
407 {
408   const struct interaction *iact = cmd->iacts[iact_idx];
409   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
410
411   int v;
412
413   for (v = 0; v < cmd->n_dep_vars; ++v)
414     {
415       int grp;
416       for (grp = 0; grp < n_cats; ++grp)
417         {
418           struct chart_item *npp, *dnpp;
419           struct casereader *reader;
420           struct np *np;
421
422           int ivar_idx;
423           const struct ccase *c =
424             categoricals_get_case_by_category_real (cmd->cats,
425                                                     iact_idx, grp);
426
427           const struct exploratory_stats *es =
428             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
429
430           struct string label;
431           ds_init_cstr (&label, 
432                         var_to_string (cmd->dep_vars[v]));
433
434           if ( iact->n_vars > 0)
435             {
436               ds_put_cstr (&label, " (");
437               for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
438                 {
439                   const struct variable *ivar = iact->vars[ivar_idx];
440                   const union value *val = case_data (c, ivar);
441                   
442                   ds_put_cstr (&label, var_to_string (ivar));
443                   ds_put_cstr (&label, " = ");
444                   append_value_name (ivar, val, &label);
445                   ds_put_cstr (&label, "; ");
446                   
447                 }
448               ds_put_cstr (&label, ")");
449             }
450           
451           np = es[v].np;
452           reader = casewriter_make_reader (np->writer);
453           np->writer = NULL;
454
455           npp = np_plot_create (np, reader, ds_cstr (&label));
456           dnpp = dnp_plot_create (np, reader, ds_cstr (&label));
457
458           if (npp == NULL || dnpp == NULL)
459             {
460               msg (MW, _("Not creating NP plot because data set is empty."));
461               chart_item_unref (npp);
462               chart_item_unref (dnpp);
463             }
464           else
465             {
466               chart_item_submit (npp);
467               chart_item_submit (dnpp);
468             }
469           casereader_destroy (reader);
470
471           ds_destroy (&label);
472         }
473     }
474 }
475
476 static void
477 show_spreadlevel (const struct examine *cmd, int iact_idx)
478 {
479   const struct interaction *iact = cmd->iacts[iact_idx];
480   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
481
482   int v;
483
484   /* Spreadlevel when there are no levels is not useful */
485   if (iact->n_vars == 0)
486     return;
487
488   for (v = 0; v < cmd->n_dep_vars; ++v)
489     {
490       int grp;
491       struct chart_item *sl;
492
493       struct string label;
494       ds_init_cstr (&label, 
495                     var_to_string (cmd->dep_vars[v]));
496
497       if (iact->n_vars > 0)
498         {
499           ds_put_cstr (&label, " (");
500           interaction_to_string (iact, &label);
501           ds_put_cstr (&label, ")");
502         }
503       
504       sl = spreadlevel_plot_create (ds_cstr (&label), cmd->sl_power);
505
506       for (grp = 0; grp < n_cats; ++grp)
507         {
508           const struct exploratory_stats *es =
509             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
510
511           double median = percentile_calculate (es[v].quartiles[1], cmd->pc_alg);
512
513           double iqr = percentile_calculate (es[v].quartiles[2], cmd->pc_alg) -
514             percentile_calculate (es[v].quartiles[0], cmd->pc_alg);
515
516           spreadlevel_plot_add (sl, iqr, median);
517         }
518
519       if (sl == NULL)
520         msg (MW, _("Not creating spreadlevel chart for %s"), ds_cstr (&label));
521       else 
522         chart_item_submit (sl);
523
524       ds_destroy (&label);
525     }
526 }
527
528
529 static void
530 show_histogram (const struct examine *cmd, int iact_idx)
531 {
532   const struct interaction *iact = cmd->iacts[iact_idx];
533   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
534
535   int v;
536
537   for (v = 0; v < cmd->n_dep_vars; ++v)
538     {
539       int grp;
540       for (grp = 0; grp < n_cats; ++grp)
541         {
542           double n, mean, var;
543           int ivar_idx;
544           const struct ccase *c =
545             categoricals_get_case_by_category_real (cmd->cats,
546                                                     iact_idx, grp);
547
548           const struct exploratory_stats *es =
549             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, grp);
550
551           struct string label;
552
553           if (es[v].histogram == NULL)
554             continue;
555
556           ds_init_cstr (&label, 
557                         var_to_string (cmd->dep_vars[v]));
558
559           if ( iact->n_vars > 0)
560             {
561               ds_put_cstr (&label, " (");
562               for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
563                 {
564                   const struct variable *ivar = iact->vars[ivar_idx];
565                   const union value *val = case_data (c, ivar);
566                   
567                   ds_put_cstr (&label, var_to_string (ivar));
568                   ds_put_cstr (&label, " = ");
569                   append_value_name (ivar, val, &label);
570                   ds_put_cstr (&label, "; ");
571                   
572                 }
573               ds_put_cstr (&label, ")");
574             }
575
576
577           moments_calculate (es[v].mom, &n, &mean, &var, NULL, NULL);
578
579           chart_item_submit
580             ( histogram_chart_create (es[v].histogram->gsl_hist,
581                                       ds_cstr (&label), n, mean,
582                                       sqrt (var), false));
583
584           
585           ds_destroy (&label);
586         }
587     }
588 }
589
590 static void
591 percentiles_report (const struct examine *cmd, int iact_idx)
592 {
593   const struct interaction *iact = cmd->iacts[iact_idx];
594   int i, v;
595   const int heading_columns = 1 + iact->n_vars + 1;
596   const int heading_rows = 2;
597   struct tab_table *t;
598
599   const size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
600
601   const int rows_per_cat = 2;
602   const int rows_per_var = n_cats * rows_per_cat;
603
604   const int nr = heading_rows + cmd->n_dep_vars * rows_per_var;
605   const int nc = heading_columns + cmd->n_percentiles;
606
607   t = tab_create (nc, nr);
608   tab_title (t, _("Percentiles"));
609
610   tab_headers (t, heading_columns, 0, heading_rows, 0);
611
612   /* Internal Vertical lines */
613   tab_box (t, -1, -1, -1, TAL_1,
614            heading_columns, 0, nc - 1, nr - 1);
615
616   /* External Frame */
617   tab_box (t, TAL_2, TAL_2, -1, -1,
618            0, 0, nc - 1, nr - 1);
619
620   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
621   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
622
623   tab_joint_text (t, heading_columns, 0,
624                   nc - 1, 0,
625                   TAT_TITLE | TAB_CENTER,
626                   _("Percentiles")
627                   );
628
629   tab_hline (t, TAL_1, heading_columns, nc - 1, 1);
630
631
632   for (i = 0; i < cmd->n_percentiles; ++i)
633     {
634       tab_text_format (t, heading_columns + i, 1,
635                        TAT_TITLE | TAB_CENTER,
636                        _("%g"), cmd->ptiles[i]);
637     }
638
639   for (i = 0; i < iact->n_vars; ++i)
640     {
641       tab_text (t,
642                 1 + i, 1,
643                 TAT_TITLE,
644                 var_to_string (iact->vars[i])
645                 );
646     }
647
648
649
650   if (n_cats > 0)
651     {
652       tab_vline (t, TAL_1, heading_columns - 1, heading_rows, nr - 1);
653
654       for (v = 0; v < cmd->n_dep_vars; ++v)
655         {
656           const union value **prev_vals = previous_value_alloc (iact);
657
658           int ivar_idx;
659           if ( v > 0 )
660             tab_hline (t, TAL_1, 0, nc - 1, heading_rows + v * rows_per_var);
661         
662           tab_text (t,
663                     0, heading_rows + v * rows_per_var,
664                     TAT_TITLE | TAB_LEFT,
665                     var_to_string (cmd->dep_vars[v])
666                     );
667
668           for (i = 0; i < n_cats; ++i)
669             {
670               const struct ccase *c =
671                 categoricals_get_case_by_category_real (cmd->cats,
672                                                         iact_idx, i);
673
674               const struct exploratory_stats *ess =
675                 categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, i);
676
677               const struct exploratory_stats *es = ess + v;
678
679               int diff_idx = previous_value_record (iact, c, prev_vals);
680
681               double hinges[3];
682               int p;
683
684               for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
685                 {
686                   const struct variable *ivar = iact->vars[ivar_idx];
687                   const union value *val = case_data (c, ivar);
688
689                   if (( diff_idx != -1 && diff_idx <= ivar_idx)
690                       || i == 0)
691                     {              
692                       struct string str;
693                       ds_init_empty (&str);
694                       append_value_name (ivar, val, &str);
695               
696                       tab_text (t,
697                                 1 + ivar_idx,
698                                 heading_rows + v * rows_per_var + i * rows_per_cat,
699                                 TAT_TITLE | TAB_LEFT,
700                                 ds_cstr (&str)
701                                 );
702                   
703                       ds_destroy (&str);
704                     }
705                 }
706
707               if ( diff_idx != -1 && diff_idx < iact->n_vars)
708                 {
709                   tab_hline (t, TAL_1, 1 + diff_idx, nc - 1,
710                              heading_rows + v * rows_per_var + i * rows_per_cat
711                              );
712                 }
713
714               tab_text (t, heading_columns - 1, 
715                         heading_rows + v * rows_per_var + i * rows_per_cat,
716                         TAT_TITLE | TAB_LEFT,
717                         gettext (ptile_alg_desc [cmd->pc_alg]));
718
719               tukey_hinges_calculate (es->hinges, hinges);
720
721               for (p = 0; p < cmd->n_percentiles; ++p)
722                 {
723                   tab_double (t, heading_columns + p, 
724                               heading_rows + v * rows_per_var + i * rows_per_cat,
725                               0,
726                               percentile_calculate (es->percentiles[p], cmd->pc_alg),
727                               0);
728               
729                   if (cmd->ptiles[p] == 25.0)
730                     {
731                       tab_double (t, heading_columns + p, 
732                                   heading_rows + v * rows_per_var + i * rows_per_cat + 1,
733                                   0,
734                                   hinges[0],
735                                   0);
736                     }
737                   else if (cmd->ptiles[p] == 50.0)
738                     {
739                       tab_double (t, heading_columns + p, 
740                                   heading_rows + v * rows_per_var + i * rows_per_cat + 1,
741                                   0,
742                                   hinges[1],
743                                   0);
744                     }
745                   else if (cmd->ptiles[p] == 75.0)
746                     {
747                       tab_double (t, heading_columns + p, 
748                                   heading_rows + v * rows_per_var + i * rows_per_cat + 1,
749                                   0,
750                                   hinges[2],
751                                   0);
752                     }
753                 }
754
755
756               tab_text (t, heading_columns - 1, 
757                         heading_rows + v * rows_per_var + i * rows_per_cat + 1,
758                         TAT_TITLE | TAB_LEFT,
759                         _("Tukey's Hinges"));
760           
761             }
762
763           free (prev_vals);
764         }
765     }
766   tab_submit (t);
767 }
768
769 static void
770 descriptives_report (const struct examine *cmd, int iact_idx)
771 {
772   const struct interaction *iact = cmd->iacts[iact_idx];
773   int i, v;
774   const int heading_columns = 1 + iact->n_vars + 2;
775   const int heading_rows = 1;
776   struct tab_table *t;
777
778   size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
779
780   const int rows_per_cat = 13;
781   const int rows_per_var = n_cats * rows_per_cat;
782
783   const int nr = heading_rows + cmd->n_dep_vars * rows_per_var;
784   const int nc = 2 + heading_columns;
785
786   t = tab_create (nc, nr);
787   tab_title (t, _("Descriptives"));
788
789   tab_headers (t, heading_columns, 0, heading_rows, 0);
790
791   /* Internal Vertical lines */
792   tab_box (t, -1, -1, -1, TAL_1,
793            heading_columns, 0, nc - 1, nr - 1);
794
795   /* External Frame */
796   tab_box (t, TAL_2, TAL_2, -1, -1,
797            0, 0, nc - 1, nr - 1);
798
799   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
800   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
801
802
803   tab_text (t, heading_columns, 0, TAB_CENTER | TAT_TITLE,
804             _("Statistic"));
805
806   tab_text (t, heading_columns + 1, 0, TAB_CENTER | TAT_TITLE,
807             _("Std. Error"));
808
809   for (i = 0; i < iact->n_vars; ++i)
810     {
811       tab_text (t,
812                 1 + i, 0,
813                 TAT_TITLE,
814                 var_to_string (iact->vars[i])
815                 );
816     }
817
818   for (v = 0; v < cmd->n_dep_vars; ++v)
819     {
820       const union value **prev_val = previous_value_alloc (iact);
821
822       int ivar_idx;
823       if ( v > 0 )
824         tab_hline (t, TAL_1, 0, nc - 1, heading_rows + v * rows_per_var);
825         
826       tab_text (t,
827                 0, heading_rows + v * rows_per_var,
828                 TAT_TITLE | TAB_LEFT,
829                 var_to_string (cmd->dep_vars[v])
830                 );
831
832       for (i = 0; i < n_cats; ++i)
833         {
834           const struct ccase *c =
835             categoricals_get_case_by_category_real (cmd->cats,
836                                                     iact_idx, i);
837
838           const struct exploratory_stats *ess =
839             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, i);
840
841           const struct exploratory_stats *es = ess + v;
842
843           const int diff_idx = previous_value_record (iact, c, prev_val);
844
845           double m0, m1, m2, m3, m4;
846           double tval;
847
848           moments_calculate (es->mom, &m0, &m1, &m2, &m3, &m4);
849
850           tval = gsl_cdf_tdist_Qinv ((1.0 - cmd->conf) / 2.0, m0 - 1.0);
851
852           for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
853             {
854               const struct variable *ivar = iact->vars[ivar_idx];
855               const union value *val = case_data (c, ivar);
856
857               if (( diff_idx != -1 && diff_idx <= ivar_idx)
858                   || i == 0)
859                 {              
860                   struct string str;
861                   ds_init_empty (&str);
862                   append_value_name (ivar, val, &str);
863               
864                   tab_text (t,
865                             1 + ivar_idx,
866                             heading_rows + v * rows_per_var + i * rows_per_cat,
867                             TAT_TITLE | TAB_LEFT,
868                             ds_cstr (&str)
869                             );
870                   
871                   ds_destroy (&str);
872                 }
873             }
874
875           if ( diff_idx != -1 && diff_idx < iact->n_vars)
876             {
877               tab_hline (t, TAL_1, 1 + diff_idx, nc - 1,
878                          heading_rows + v * rows_per_var + i * rows_per_cat
879                          );
880             }
881
882           tab_text (t,
883                     1 + iact->n_vars,
884                     heading_rows + v * rows_per_var + i * rows_per_cat,
885                     TAB_LEFT,
886                     _("Mean")
887                     );
888
889           tab_double (t,
890                       1 + iact->n_vars + 2,
891                       heading_rows + v * rows_per_var + i * rows_per_cat,
892                       0, m1, 0);
893
894           tab_double (t,
895                       1 + iact->n_vars + 3,
896                       heading_rows + v * rows_per_var + i * rows_per_cat,
897                       0, calc_semean (m2, m0), 0);
898
899           tab_text_format (t,
900                            1 + iact->n_vars,
901                            heading_rows + v * rows_per_var + i * rows_per_cat + 1,
902                            TAB_LEFT,
903                            _("%g%% Confidence Interval for Mean"),
904                            cmd->conf * 100.0
905                            );
906           
907           tab_text (t,
908                     1 + iact->n_vars + 1,
909                     heading_rows + v * rows_per_var + i * rows_per_cat + 1,
910                     TAB_LEFT,
911                     _("Lower Bound")
912                     );
913
914           tab_double (t,
915                       1 + iact->n_vars + 2,
916                       heading_rows + v * rows_per_var + i * rows_per_cat + 1,
917                       0, m1 - tval * calc_semean (m2, m0), 0);
918
919
920           tab_text (t,
921                     1 + iact->n_vars + 1,
922                     heading_rows + v * rows_per_var + i * rows_per_cat + 2,
923                     TAB_LEFT,
924                     _("Upper Bound")
925                     );
926
927           tab_double (t,
928                       1 + iact->n_vars + 2,
929                       heading_rows + v * rows_per_var + i * rows_per_cat + 2,
930                       0, m1 + tval * calc_semean (m2, m0), 0);
931
932
933           tab_text (t,
934                     1 + iact->n_vars,
935                     heading_rows + v * rows_per_var + i * rows_per_cat + 3,
936                     TAB_LEFT,
937                     _("5% Trimmed Mean")
938                     );
939
940           tab_double (t,
941                       1 + iact->n_vars + 2,
942                       heading_rows + v * rows_per_var + i * rows_per_cat + 3,
943                       0,
944                       trimmed_mean_calculate (es->trimmed_mean),
945                       0);
946
947           tab_text (t,
948                     1 + iact->n_vars,
949                     heading_rows + v * rows_per_var + i * rows_per_cat + 4,
950                     TAB_LEFT,
951                     _("Median")
952                     );
953           
954           tab_double (t,
955                       1 + iact->n_vars + 2,
956                       heading_rows + v * rows_per_var + i * rows_per_cat + 4,
957                       0,
958                       percentile_calculate (es->quartiles[1], cmd->pc_alg),
959                       0);
960
961
962           tab_text (t,
963                     1 + iact->n_vars,
964                     heading_rows + v * rows_per_var + i * rows_per_cat + 5,
965                     TAB_LEFT,
966                     _("Variance")
967                     );
968
969           tab_double (t,
970                       1 + iact->n_vars + 2,
971                       heading_rows + v * rows_per_var + i * rows_per_cat + 5,
972                       0, m2, 0);
973
974           tab_text (t,
975                     1 + iact->n_vars,
976                     heading_rows + v * rows_per_var + i * rows_per_cat + 6,
977                     TAB_LEFT,
978                     _("Std. Deviation")
979                     );
980
981           tab_double (t,
982                       1 + iact->n_vars + 2,
983                       heading_rows + v * rows_per_var + i * rows_per_cat + 6,
984                       0, sqrt (m2), 0);
985
986           tab_text (t,
987                     1 + iact->n_vars,
988                     heading_rows + v * rows_per_var + i * rows_per_cat + 7,
989                     TAB_LEFT,
990                     _("Minimum")
991                     );
992
993           tab_double (t,
994                       1 + iact->n_vars + 2,
995                       heading_rows + v * rows_per_var + i * rows_per_cat + 7,
996                       0, 
997                       es->minima[0].val,
998                       0);
999
1000           tab_text (t,
1001                     1 + iact->n_vars,
1002                     heading_rows + v * rows_per_var + i * rows_per_cat + 8,
1003                     TAB_LEFT,
1004                     _("Maximum")
1005                     );
1006
1007           tab_double (t,
1008                       1 + iact->n_vars + 2,
1009                       heading_rows + v * rows_per_var + i * rows_per_cat + 8,
1010                       0, 
1011                       es->maxima[0].val,
1012                       0);
1013
1014           tab_text (t,
1015                     1 + iact->n_vars,
1016                     heading_rows + v * rows_per_var + i * rows_per_cat + 9,
1017                     TAB_LEFT,
1018                     _("Range")
1019                     );
1020
1021           tab_double (t,
1022                       1 + iact->n_vars + 2,
1023                       heading_rows + v * rows_per_var + i * rows_per_cat + 9,
1024                       0, 
1025                       es->maxima[0].val - es->minima[0].val,
1026                       0);
1027
1028           tab_text (t,
1029                     1 + iact->n_vars,
1030                     heading_rows + v * rows_per_var + i * rows_per_cat + 10,
1031                     TAB_LEFT,
1032                     _("Interquartile Range")
1033                     );
1034
1035
1036           tab_double (t,
1037                       1 + iact->n_vars + 2,
1038                       heading_rows + v * rows_per_var + i * rows_per_cat + 10,
1039                       0,
1040                       percentile_calculate (es->quartiles[2], cmd->pc_alg) - 
1041                       percentile_calculate (es->quartiles[0], cmd->pc_alg),
1042                       0);
1043
1044
1045
1046
1047           tab_text (t,
1048                     1 + iact->n_vars,
1049                     heading_rows + v * rows_per_var + i * rows_per_cat + 11,
1050                     TAB_LEFT,
1051                     _("Skewness")
1052                     );
1053
1054           tab_double (t,
1055                       1 + iact->n_vars + 2,
1056                       heading_rows + v * rows_per_var + i * rows_per_cat + 11,
1057                       0, m3, 0);
1058
1059           tab_double (t,
1060                       1 + iact->n_vars + 3,
1061                       heading_rows + v * rows_per_var + i * rows_per_cat + 11,
1062                       0, calc_seskew (m0), 0);
1063
1064           tab_text (t,
1065                     1 + iact->n_vars,
1066                     heading_rows + v * rows_per_var + i * rows_per_cat + 12,
1067                     TAB_LEFT,
1068                     _("Kurtosis")
1069                     );
1070
1071           tab_double (t,
1072                       1 + iact->n_vars + 2,
1073                       heading_rows + v * rows_per_var + i * rows_per_cat + 12,
1074                       0, m4, 0);
1075
1076           tab_double (t,
1077                       1 + iact->n_vars + 3,
1078                       heading_rows + v * rows_per_var + i * rows_per_cat + 12,
1079                       0, calc_sekurt (m0), 0);
1080         }
1081
1082       free (prev_val);
1083     }
1084   tab_submit (t);
1085 }
1086
1087
1088 static void
1089 extremes_report (const struct examine *cmd, int iact_idx)
1090 {
1091   const struct interaction *iact = cmd->iacts[iact_idx];
1092   int i, v;
1093   const int heading_columns = 1 + iact->n_vars + 2;
1094   const int heading_rows = 1;
1095   struct tab_table *t;
1096
1097   size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
1098
1099   const int rows_per_cat = 2 * cmd->disp_extremes;
1100   const int rows_per_var = n_cats * rows_per_cat;
1101
1102   const int nr = heading_rows + cmd->n_dep_vars * rows_per_var;
1103   const int nc = 2 + heading_columns;
1104
1105   t = tab_create (nc, nr);
1106   tab_title (t, _("Extreme Values"));
1107
1108   tab_headers (t, heading_columns, 0, heading_rows, 0);
1109
1110   /* Internal Vertical lines */
1111   tab_box (t, -1, -1, -1, TAL_1,
1112            heading_columns, 0, nc - 1, nr - 1);
1113
1114   /* External Frame */
1115   tab_box (t, TAL_2, TAL_2, -1, -1,
1116            0, 0, nc - 1, nr - 1);
1117
1118   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
1119   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
1120
1121
1122   if ( cmd->id_var ) 
1123     tab_text (t, heading_columns, 0, TAB_CENTER | TAT_TITLE,
1124               var_to_string (cmd->id_var));
1125   else
1126     tab_text (t, heading_columns, 0, TAB_CENTER | TAT_TITLE,
1127               _("Case Number"));
1128
1129   tab_text (t, heading_columns + 1, 0, TAB_CENTER | TAT_TITLE,
1130             _("Value"));
1131
1132   for (i = 0; i < iact->n_vars; ++i)
1133     {
1134       tab_text (t,
1135                 1 + i, 0,
1136                 TAT_TITLE,
1137                 var_to_string (iact->vars[i])
1138                 );
1139     }
1140
1141   for (v = 0; v < cmd->n_dep_vars; ++v)
1142     {
1143       const union value **prev_val = previous_value_alloc (iact);
1144
1145       int ivar_idx;
1146       if ( v > 0 )
1147         tab_hline (t, TAL_1, 0, nc - 1, heading_rows + v * rows_per_var);
1148         
1149       tab_text (t,
1150                 0, heading_rows + v * rows_per_var,
1151                 TAT_TITLE,
1152                 var_to_string (cmd->dep_vars[v])
1153                 );
1154
1155       for (i = 0; i < n_cats; ++i)
1156         {
1157           int e;
1158           const struct ccase *c =
1159             categoricals_get_case_by_category_real (cmd->cats, iact_idx, i);
1160
1161           const struct exploratory_stats *ess =
1162             categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, i);
1163
1164           const struct exploratory_stats *es = ess + v;
1165
1166           int diff_idx = previous_value_record (iact, c, prev_val);
1167
1168           for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
1169             {
1170               const struct variable *ivar = iact->vars[ivar_idx];
1171               const union value *val = case_data (c, ivar);
1172
1173               if (( diff_idx != -1 && diff_idx <= ivar_idx)
1174                   || i == 0)
1175                 {              
1176                   struct string str;
1177                   ds_init_empty (&str);
1178                   append_value_name (ivar, val, &str);
1179               
1180                   tab_text (t,
1181                             1 + ivar_idx,
1182                             heading_rows + v * rows_per_var + i * rows_per_cat,
1183                             TAT_TITLE | TAB_LEFT,
1184                             ds_cstr (&str)
1185                             );
1186                   
1187                   ds_destroy (&str);
1188                 }
1189             }
1190
1191           if ( diff_idx != -1 && diff_idx < iact->n_vars)
1192             {
1193               tab_hline (t, TAL_1, 1 + diff_idx, nc - 1,
1194                          heading_rows + v * rows_per_var + i * rows_per_cat
1195                          );
1196             }
1197           
1198           tab_text (t,
1199                     heading_columns - 2,
1200                     heading_rows + v * rows_per_var + i * rows_per_cat,
1201                     TAB_RIGHT,
1202                     _("Highest"));
1203
1204
1205           tab_hline (t, TAL_1, heading_columns - 2, nc - 1,
1206                      heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes
1207                      );
1208
1209           tab_text (t,
1210                     heading_columns - 2,
1211                     heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes,
1212                     TAB_RIGHT,
1213                     _("Lowest"));
1214
1215           for (e = 0 ; e < cmd->disp_extremes; ++e)
1216             {
1217               tab_double (t,
1218                           heading_columns - 1,
1219                           heading_rows + v * rows_per_var + i * rows_per_cat + e,
1220                           TAB_RIGHT,
1221                           e + 1,
1222                           &F_8_0);
1223
1224               /* The casenumber */
1225               if (cmd->id_var)
1226                 tab_value (t,
1227                            heading_columns,
1228                            heading_rows + v * rows_per_var + i * rows_per_cat + e,
1229                            TAB_RIGHT,
1230                            &es->maxima[e].identity,
1231                            cmd->id_var,
1232                            NULL);
1233               else 
1234                 tab_double (t,
1235                           heading_columns,
1236                             heading_rows + v * rows_per_var + i * rows_per_cat + e,
1237                             TAB_RIGHT,
1238                             es->maxima[e].identity.f,
1239                             &F_8_0);
1240
1241               tab_double (t,
1242                          heading_columns + 1,
1243                          heading_rows + v * rows_per_var + i * rows_per_cat + e,
1244                          0,
1245                          es->maxima[e].val,
1246                          var_get_print_format (cmd->dep_vars[v]));
1247                          
1248
1249               tab_double (t,
1250                           heading_columns - 1,
1251                           heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes + e,
1252                           TAB_RIGHT,
1253                           e + 1,
1254                           &F_8_0);
1255
1256               /* The casenumber */
1257               if (cmd->id_var)
1258                 tab_value (t,
1259                            heading_columns,
1260                            heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes + e,
1261                            TAB_RIGHT,
1262                            &es->minima[e].identity,
1263                            cmd->id_var,
1264                            NULL);
1265               else
1266                 tab_double (t,
1267                             heading_columns,
1268                             heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes + e,
1269                             TAB_RIGHT,
1270                             es->minima[e].identity.f,
1271                             &F_8_0);
1272
1273               tab_double (t,
1274                           heading_columns + 1,
1275                           heading_rows + v * rows_per_var + i * rows_per_cat + cmd->disp_extremes + e,
1276                           0,
1277                           es->minima[e].val,
1278                           var_get_print_format (cmd->dep_vars[v]));
1279             }
1280         }
1281       free (prev_val);
1282     }
1283
1284   tab_submit (t);
1285 }
1286
1287
1288 static void
1289 summary_report (const struct examine *cmd, int iact_idx)
1290 {
1291   const struct interaction *iact = cmd->iacts[iact_idx];
1292   int i, v;
1293   const int heading_columns = 1 + iact->n_vars;
1294   const int heading_rows = 3;
1295   struct tab_table *t;
1296
1297   const struct fmt_spec *wfmt = cmd->wv ? var_get_print_format (cmd->wv) : &F_8_0;
1298
1299   size_t n_cats =  categoricals_n_count (cmd->cats, iact_idx);
1300
1301   const int nr = heading_rows + n_cats * cmd->n_dep_vars;
1302   const int nc = 6 + heading_columns;
1303
1304   t = tab_create (nc, nr);
1305   tab_title (t, _("Case Processing Summary"));
1306
1307   tab_headers (t, heading_columns, 0, heading_rows, 0);
1308
1309   /* Internal Vertical lines */
1310   tab_box (t, -1, -1, -1, TAL_1,
1311            heading_columns, 0, nc - 1, nr - 1);
1312
1313   /* External Frame */
1314   tab_box (t, TAL_2, TAL_2, -1, -1,
1315            0, 0, nc - 1, nr - 1);
1316
1317   tab_hline (t, TAL_2, 0, nc - 1, heading_rows);
1318   tab_vline (t, TAL_2, heading_columns, 0, nr - 1);
1319
1320   tab_joint_text (t, heading_columns, 0,
1321                   nc - 1, 0, TAB_CENTER | TAT_TITLE, _("Cases"));
1322   tab_joint_text (t,
1323                   heading_columns, 1,
1324                   heading_columns + 1, 1,
1325                   TAB_CENTER | TAT_TITLE, _("Valid"));
1326
1327   tab_joint_text (t,
1328                   heading_columns + 2, 1, 
1329                   heading_columns + 3, 1,
1330                   TAB_CENTER | TAT_TITLE, _("Missing"));
1331
1332   tab_joint_text (t,
1333                   heading_columns + 4, 1,
1334                   heading_columns + 5, 1,
1335                   TAB_CENTER | TAT_TITLE, _("Total"));
1336
1337   for (i = 0; i < 3; ++i)
1338     {
1339       tab_text (t, heading_columns + i * 2, 2, TAB_CENTER | TAT_TITLE,
1340                 _("N"));
1341       tab_text (t, heading_columns + i * 2 + 1, 2, TAB_CENTER | TAT_TITLE,
1342                 _("Percent"));
1343     }
1344
1345   for (i = 0; i < iact->n_vars; ++i)
1346     {
1347       tab_text (t,
1348                 1 + i, 2,
1349                 TAT_TITLE,
1350                 var_to_string (iact->vars[i])
1351                 );
1352     }
1353
1354   if (n_cats > 0)
1355     for (v = 0; v < cmd->n_dep_vars; ++v)
1356       {
1357         int ivar_idx;
1358         const union value **prev_values = previous_value_alloc (iact);
1359
1360         if ( v > 0 )
1361           tab_hline (t, TAL_1, 0, nc - 1, heading_rows + v * n_cats);
1362
1363         tab_text (t,
1364                   0, heading_rows + n_cats * v,
1365                   TAT_TITLE,
1366                   var_to_string (cmd->dep_vars[v])
1367                   );
1368
1369
1370         for (i = 0; i < n_cats; ++i)
1371           {
1372             double total;
1373             const struct exploratory_stats *es;
1374
1375             const struct ccase *c =
1376               categoricals_get_case_by_category_real (cmd->cats,
1377                                                       iact_idx, i);
1378             if (c)
1379               {
1380                 int diff_idx = previous_value_record (iact, c, prev_values);
1381
1382                 if ( diff_idx != -1 && diff_idx < iact->n_vars - 1)
1383                   tab_hline (t, TAL_1, 1 + diff_idx, nc - 1,
1384                              heading_rows + n_cats * v + i );
1385
1386                 for (ivar_idx = 0; ivar_idx < iact->n_vars; ++ivar_idx)
1387                   {
1388                     const struct variable *ivar = iact->vars[ivar_idx];
1389                     const union value *val = case_data (c, ivar);
1390
1391                     if (( diff_idx != -1 && diff_idx <= ivar_idx)
1392                         || i == 0)
1393                       {              
1394                         struct string str;
1395                         ds_init_empty (&str);
1396                         append_value_name (ivar, val, &str);
1397               
1398                         tab_text (t,
1399                                   1 + ivar_idx, heading_rows + n_cats * v + i,
1400                                   TAT_TITLE | TAB_LEFT,
1401                                   ds_cstr (&str)
1402                                   );
1403                   
1404                         ds_destroy (&str);
1405                       }
1406                   }
1407               }
1408
1409
1410             es = categoricals_get_user_data_by_category_real (cmd->cats, iact_idx, i);
1411   
1412           
1413             total = es[v].missing + es[v].non_missing;
1414             tab_double (t, 
1415                         heading_columns + 0,
1416                         heading_rows + n_cats * v + i,
1417                         0,
1418                         es[v].non_missing,
1419                         wfmt);
1420
1421
1422             tab_text_format (t, 
1423                              heading_columns + 1,
1424                              heading_rows + n_cats * v + i,
1425                              0,
1426                              "%g%%",
1427                              100.0 * es[v].non_missing / total
1428                              );
1429
1430
1431             tab_double (t, 
1432                         heading_columns + 2,
1433                         heading_rows + n_cats * v + i,
1434                         0,
1435                         es[v].missing,
1436                         wfmt);
1437
1438             tab_text_format (t, 
1439                              heading_columns + 3,
1440                              heading_rows + n_cats * v + i,
1441                              0,
1442                              "%g%%",
1443                              100.0 * es[v].missing / total
1444                              );
1445             tab_double (t, 
1446                         heading_columns + 4,
1447                         heading_rows + n_cats * v + i,
1448                         0,
1449                         total,
1450                         wfmt);
1451
1452             /* This can only be 100% can't it? */
1453             tab_text_format (t, 
1454                              heading_columns + 5,
1455                              heading_rows + n_cats * v + i,
1456                              0,
1457                              "%g%%",
1458                              100.0 * (es[v].missing + es[v].non_missing)/ total
1459                              );
1460           }
1461         free (prev_values);
1462       }
1463
1464   tab_hline (t, TAL_1, heading_columns, nc - 1, 1);
1465   tab_hline (t, TAL_1, heading_columns, nc - 1, 2);
1466
1467   tab_submit (t);
1468 }
1469
1470 /* Attempt to parse an interaction from LEXER */
1471 static struct interaction *
1472 parse_interaction (struct lexer *lexer, struct examine *ex)
1473 {
1474   const struct variable *v = NULL;
1475   struct interaction *iact = NULL;
1476   
1477   if ( lex_match_variable (lexer, ex->dict, &v))
1478     {
1479       iact = interaction_create (v);
1480
1481       while (lex_match (lexer, T_BY))
1482         {
1483           if (!lex_match_variable (lexer, ex->dict, &v))
1484             {
1485               interaction_destroy (iact);
1486               return NULL;
1487             }
1488           interaction_add_variable (iact, v);
1489         }
1490       lex_match (lexer, T_COMMA);
1491     }
1492   
1493   return iact;
1494 }
1495
1496
1497 static void *
1498 create_n (const void *aux1, void *aux2 UNUSED)
1499 {
1500   int v;
1501   
1502   const struct examine *examine = aux1;
1503   struct exploratory_stats *es = pool_calloc (examine->pool, examine->n_dep_vars, sizeof (*es));
1504   struct subcase ordering;
1505   subcase_init (&ordering, 0, 0, SC_ASCEND);
1506
1507   for (v = 0; v < examine->n_dep_vars; v++)
1508     {
1509       es[v].sorted_writer = sort_create_writer (&ordering, examine->ex_proto);
1510       es[v].sorted_reader = NULL;
1511
1512       es[v].mom = moments_create (MOMENT_KURTOSIS);
1513       es[v].cmin = DBL_MAX;
1514
1515       es[v].maximum = -DBL_MAX;
1516       es[v].minimum =  DBL_MAX;
1517     }
1518
1519   subcase_destroy (&ordering);
1520   return es;
1521 }
1522
1523 static void
1524 update_n (const void *aux1, void *aux2 UNUSED, void *user_data,
1525           const struct ccase *c, double weight)
1526 {
1527   int v;
1528   const struct examine *examine = aux1;
1529   struct exploratory_stats *es = user_data;
1530
1531   for (v = 0; v < examine->n_dep_vars; v++)
1532     {
1533       struct ccase *outcase ;
1534       const struct variable *var = examine->dep_vars[v];
1535       const double x = case_data (c, var)->f;
1536       
1537       if (var_is_value_missing (var, case_data (c, var), examine->dep_excl))
1538         {
1539           es[v].missing += weight;
1540           continue;
1541         }
1542
1543       outcase = case_create (examine->ex_proto);
1544
1545       if (x > es[v].maximum)
1546         es[v].maximum = x;
1547
1548       if (x < es[v].minimum)
1549         es[v].minimum =  x;
1550
1551       es[v].non_missing += weight;
1552
1553       moments_pass_one (es[v].mom, x, weight);
1554
1555       /* Save the value and the ID to the writer */
1556       assert (examine->id_idx != -1);
1557       case_data_rw_idx (outcase, EX_VAL)->f = x;
1558       value_copy (case_data_rw_idx (outcase, EX_ID),
1559                   case_data_idx (c, examine->id_idx), examine->id_width);
1560
1561       case_data_rw_idx (outcase, EX_WT)->f = weight;
1562       
1563       es[v].cc += weight;
1564
1565       if (es[v].cmin > weight)
1566         es[v].cmin = weight;
1567
1568       casewriter_write (es[v].sorted_writer, outcase);
1569     }
1570 }
1571
1572 static void
1573 calculate_n (const void *aux1, void *aux2 UNUSED, void *user_data)
1574 {
1575   int v;
1576   const struct examine *examine = aux1;
1577   struct exploratory_stats *es = user_data;
1578
1579   for (v = 0; v < examine->n_dep_vars; v++)
1580     {
1581       int i;
1582       casenumber imin = 0;
1583       casenumber imax;
1584       struct casereader *reader;
1585       struct ccase *c;
1586
1587       if (examine->histogram)
1588         {
1589           /* Sturges Rule */
1590           double bin_width = fabs (es[v].minimum - es[v].maximum)
1591             / (1 + log2 (es[v].cc))
1592             ;
1593
1594           es[v].histogram =
1595             histogram_create (bin_width, es[v].minimum, es[v].maximum);
1596         }
1597
1598       es[v].sorted_reader = casewriter_make_reader (es[v].sorted_writer);
1599       es[v].sorted_writer = NULL;
1600
1601       imax = casereader_get_case_cnt (es[v].sorted_reader);
1602
1603       es[v].maxima = pool_calloc (examine->pool, examine->calc_extremes, sizeof (*es[v].maxima));
1604       es[v].minima = pool_calloc (examine->pool, examine->calc_extremes, sizeof (*es[v].minima));
1605       for (i = 0; i < examine->calc_extremes; ++i)
1606         {
1607           value_init_pool (examine->pool, &es[v].maxima[i].identity, examine->id_width) ;
1608           value_init_pool (examine->pool, &es[v].minima[i].identity, examine->id_width) ;
1609         }
1610       
1611       for (reader = casereader_clone (es[v].sorted_reader);
1612            (c = casereader_read (reader)) != NULL; case_unref (c))
1613         {
1614           const double val = case_data_idx (c, EX_VAL)->f;
1615           const double wt = case_data_idx (c, EX_WT)->f;
1616
1617           moments_pass_two (es[v].mom, val, wt);
1618
1619           if (es[v].histogram)
1620             histogram_add (es[v].histogram, val, wt);
1621
1622           if (imin < examine->calc_extremes)
1623             {
1624               int x;
1625               for (x = imin; x < examine->calc_extremes; ++x)
1626                 {
1627                   struct extremity *min = &es[v].minima[x];
1628                   min->val = val;
1629                   value_copy (&min->identity, case_data_idx (c, EX_ID), examine->id_width);
1630                 }
1631               imin ++;
1632             }
1633
1634           imax --;
1635           if (imax < examine->calc_extremes)
1636             {
1637               int x;
1638
1639               for (x = imax; x < imax + 1; ++x)
1640                 {
1641                   struct extremity *max;
1642
1643                   if (x >= examine->calc_extremes) 
1644                     break;
1645
1646                   max = &es[v].maxima[x];
1647                   max->val = val;
1648                   value_copy (&max->identity, case_data_idx (c, EX_ID), examine->id_width);
1649                 }
1650             }
1651         }
1652       casereader_destroy (reader);
1653
1654       if (examine->calc_extremes > 0)
1655         {
1656           assert (es[v].minima[0].val == es[v].minimum);
1657           assert (es[v].maxima[0].val == es[v].maximum);
1658         }
1659
1660       {
1661         const int n_os = 5 + examine->n_percentiles;
1662         struct order_stats **os ;
1663         es[v].percentiles = pool_calloc (examine->pool, examine->n_percentiles, sizeof (*es[v].percentiles));
1664
1665         es[v].trimmed_mean = trimmed_mean_create (es[v].cc, 0.05);
1666
1667         os = xcalloc (n_os, sizeof *os);
1668         os[0] = &es[v].trimmed_mean->parent;
1669
1670         es[v].quartiles[0] = percentile_create (0.25, es[v].cc);
1671         es[v].quartiles[1] = percentile_create (0.5,  es[v].cc);
1672         es[v].quartiles[2] = percentile_create (0.75, es[v].cc);
1673
1674         os[1] = &es[v].quartiles[0]->parent;
1675         os[2] = &es[v].quartiles[1]->parent;
1676         os[3] = &es[v].quartiles[2]->parent;
1677
1678         es[v].hinges = tukey_hinges_create (es[v].cc, es[v].cmin);
1679         os[4] = &es[v].hinges->parent;
1680
1681         for (i = 0; i < examine->n_percentiles; ++i)
1682           {
1683             es[v].percentiles[i] = percentile_create (examine->ptiles[i] / 100.00, es[v].cc);
1684             os[5 + i] = &es[v].percentiles[i]->parent;
1685           }
1686
1687         order_stats_accumulate_idx (os, n_os,
1688                                     casereader_clone (es[v].sorted_reader),
1689                                     EX_WT, EX_VAL);
1690
1691         free (os);
1692       }
1693
1694       if (examine->boxplot)
1695         {
1696           struct order_stats *os;
1697
1698           es[v].box_whisker = box_whisker_create (es[v].hinges, 
1699                                                   EX_ID, examine->id_var);
1700
1701           os = &es[v].box_whisker->parent;
1702           order_stats_accumulate_idx (&os, 1,
1703                                       casereader_clone (es[v].sorted_reader),
1704                                       EX_WT, EX_VAL);
1705         }
1706
1707       if (examine->npplot)
1708         {
1709           double n, mean, var;
1710           struct order_stats *os;
1711
1712           moments_calculate (es[v].mom, &n, &mean, &var, NULL, NULL);
1713           
1714           es[v].np = np_create (n, mean, var);
1715
1716           os = &es[v].np->parent;
1717
1718           order_stats_accumulate_idx (&os, 1,
1719                                       casereader_clone (es[v].sorted_reader),
1720                                       EX_WT, EX_VAL);
1721         }
1722
1723     }
1724 }
1725
1726 static void
1727 cleanup_exploratory_stats (struct examine *cmd)
1728
1729   int i;
1730   for (i = 0; i < cmd->n_iacts; ++i)
1731     {
1732       int v;
1733       const size_t n_cats =  categoricals_n_count (cmd->cats, i);
1734
1735       for (v = 0; v < cmd->n_dep_vars; ++v)
1736         {
1737           int grp;
1738           for (grp = 0; grp < n_cats; ++grp)
1739             {
1740               int q;
1741               const struct exploratory_stats *es =
1742                 categoricals_get_user_data_by_category_real (cmd->cats, i, grp);
1743
1744               struct order_stats *os = &es[v].hinges->parent;
1745               struct statistic  *stat = &os->parent;
1746               stat->destroy (stat);
1747
1748               for (q = 0; q < 3 ; q++)
1749                 {
1750                   os = &es[v].quartiles[q]->parent;
1751                   stat = &os->parent;
1752                   stat->destroy (stat);
1753                 }
1754
1755               for (q = 0; q < cmd->n_percentiles ; q++)
1756                 {
1757                   os = &es[v].percentiles[q]->parent;
1758                   stat = &os->parent;
1759                   stat->destroy (stat);
1760                 }
1761
1762               os = &es[v].trimmed_mean->parent;
1763               stat = &os->parent;
1764               stat->destroy (stat);
1765
1766               os = &es[v].np->parent;
1767               if (os)
1768                 {
1769                   stat = &os->parent;
1770                   stat->destroy (stat);
1771                 }
1772
1773               statistic_destroy (&es[v].histogram->parent);
1774               moments_destroy (es[v].mom);
1775
1776               casereader_destroy (es[v].sorted_reader);
1777             }
1778         }
1779     }
1780 }
1781
1782
1783 static void
1784 run_examine (struct examine *cmd, struct casereader *input)
1785 {
1786   int i;
1787   struct ccase *c;
1788   struct casereader *reader;
1789
1790   struct payload payload;
1791   payload.create = create_n;
1792   payload.update = update_n;
1793   payload.calculate = calculate_n;
1794   payload.destroy = NULL;
1795   
1796   cmd->wv = dict_get_weight (cmd->dict);
1797
1798   cmd->cats
1799     = categoricals_create (cmd->iacts, cmd->n_iacts,  
1800                            cmd->wv, cmd->dep_excl, cmd->fctr_excl);
1801
1802   categoricals_set_payload (cmd->cats, &payload, cmd, NULL);
1803
1804   if (cmd->id_idx == -1)
1805     {
1806       struct ccase *c = casereader_peek (input,  0);
1807
1808       assert (cmd->id_var == NULL);
1809
1810       cmd->id_idx = case_get_value_cnt (c);
1811       input = casereader_create_arithmetic_sequence (input, 1.0, 1.0);
1812
1813       case_unref (c);
1814     }
1815
1816   /* Remove cases on a listwise basis if requested */
1817   if ( cmd->missing_pw == false)
1818     input = casereader_create_filter_missing (input,
1819                                               cmd->dep_vars,
1820                                               cmd->n_dep_vars,
1821                                               cmd->dep_excl,
1822                                               NULL,
1823                                               NULL);
1824
1825   for (reader = input;
1826        (c = casereader_read (reader)) != NULL; case_unref (c))
1827     {
1828       categoricals_update (cmd->cats, c);
1829     }
1830   casereader_destroy (reader);
1831   categoricals_done (cmd->cats);
1832
1833   for (i = 0; i < cmd->n_iacts; ++i)
1834     {
1835       summary_report (cmd, i);
1836
1837       if (cmd->disp_extremes > 0)
1838         extremes_report (cmd, i);
1839
1840       if (cmd->n_percentiles > 0)
1841         percentiles_report (cmd, i);
1842
1843       if (cmd->boxplot)
1844         {
1845           switch (cmd->boxplot_mode)
1846             {
1847             case BP_GROUPS:
1848               show_boxplot_grouped (cmd, i);
1849               break;
1850             case BP_VARIABLES:
1851               show_boxplot_variabled (cmd, i);
1852               break;
1853             default:
1854               NOT_REACHED ();
1855               break;
1856             }
1857         }
1858
1859       if (cmd->histogram)
1860         show_histogram (cmd, i);
1861
1862       if (cmd->npplot)
1863         show_npplot (cmd, i);
1864
1865       if (cmd->spreadlevel)
1866         show_spreadlevel (cmd, i);
1867
1868       if (cmd->descriptives)
1869         descriptives_report (cmd, i);
1870     }
1871
1872   cleanup_exploratory_stats (cmd);
1873   categoricals_destroy (cmd->cats);
1874 }
1875
1876
1877 int
1878 cmd_examine (struct lexer *lexer, struct dataset *ds)
1879 {
1880   int i;
1881   bool nototals_seen = false;
1882   bool totals_seen = false;
1883
1884   struct interaction **iacts_mem = NULL;
1885   struct examine examine;
1886   bool percentiles_seen = false;
1887
1888   examine.missing_pw = false;
1889   examine.disp_extremes = 0;
1890   examine.calc_extremes = 0;
1891   examine.descriptives = false;
1892   examine.conf = 0.95;
1893   examine.pc_alg = PC_HAVERAGE;
1894   examine.ptiles = NULL;
1895   examine.n_percentiles = 0;
1896   examine.id_idx = -1;
1897   examine.id_width = 0;
1898   examine.id_var = NULL;
1899   examine.boxplot_mode = BP_GROUPS;
1900   
1901   examine.ex_proto = caseproto_create ();
1902
1903   examine.pool = pool_create ();
1904
1905   /* Allocate space for the first interaction.
1906      This is interaction is an empty one (for the totals).
1907      If no totals are requested, we will simply ignore this
1908      interaction.
1909   */
1910   examine.n_iacts = 1;
1911   examine.iacts = iacts_mem = pool_zalloc (examine.pool, sizeof (struct interaction *));
1912   examine.iacts[0] = interaction_create (NULL);
1913
1914   examine.dep_excl = MV_ANY;
1915   examine.fctr_excl = MV_ANY;
1916   examine.histogram = false;
1917   examine.npplot = false;
1918   examine.boxplot = false;
1919   examine.spreadlevel = false;
1920   examine.sl_power = 0;
1921   
1922   examine.dict = dataset_dict (ds);
1923
1924   /* Accept an optional, completely pointless "/VARIABLES=" */
1925   lex_match (lexer, T_SLASH);
1926   if (lex_match_id  (lexer, "VARIABLES"))
1927     {
1928       if (! lex_force_match (lexer, T_EQUALS) )
1929         goto error;
1930     }
1931
1932   if (!parse_variables_const (lexer, examine.dict,
1933                               &examine.dep_vars, &examine.n_dep_vars,
1934                               PV_NO_DUPLICATE | PV_NUMERIC))
1935     goto error;
1936
1937   if (lex_match (lexer, T_BY))
1938     {
1939       struct interaction *iact = NULL;
1940       do
1941         {
1942           iact = parse_interaction (lexer, &examine);
1943           if (iact)
1944             {
1945               examine.n_iacts++;
1946               iacts_mem = 
1947                 pool_nrealloc (examine.pool, iacts_mem,
1948                                examine.n_iacts,
1949                                sizeof (*iacts_mem));
1950               
1951               iacts_mem[examine.n_iacts - 1] = iact;
1952             }
1953         }
1954       while (iact);
1955     }
1956
1957
1958   while (lex_token (lexer) != T_ENDCMD)
1959     {
1960       lex_match (lexer, T_SLASH);
1961
1962       if (lex_match_id (lexer, "STATISTICS"))
1963         {
1964           lex_match (lexer, T_EQUALS);
1965
1966           while (lex_token (lexer) != T_ENDCMD
1967                  && lex_token (lexer) != T_SLASH)
1968             {
1969               if (lex_match_id (lexer, "DESCRIPTIVES"))
1970                 {
1971                   examine.descriptives = true;
1972                 }
1973               else if (lex_match_id (lexer, "EXTREME"))
1974                 {
1975                   int extr = 5;
1976                   if (lex_match (lexer, T_LPAREN))
1977                     {
1978                       extr = lex_integer (lexer);
1979
1980                       if (extr < 0)
1981                         {
1982                           msg (MW, _("%s may not be negative. Using default value (%g)."), "EXTREME", 5.0);
1983                           extr = 5;
1984                         }
1985
1986                       lex_get (lexer);
1987                       if (! lex_force_match (lexer, T_RPAREN))
1988                         goto error;
1989                     }
1990                   examine.disp_extremes  = extr;
1991                 }
1992               else if (lex_match_id (lexer, "NONE"))
1993                 {
1994                 }
1995               else if (lex_match (lexer, T_ALL))
1996                 {
1997                   if (examine.disp_extremes == 0)
1998                     examine.disp_extremes = 5;
1999                 }
2000               else
2001                 {
2002                   lex_error (lexer, NULL);
2003                   goto error;
2004                 }
2005             }
2006         }
2007       else if (lex_match_id (lexer, "PERCENTILES"))
2008         {
2009           percentiles_seen = true;
2010           if (lex_match (lexer, T_LPAREN))
2011             {
2012               while (lex_is_number (lexer))
2013                 {
2014                   double p = lex_number (lexer);
2015                   
2016                   if ( p <= 0 || p >= 100.0)
2017                     {
2018                       lex_error (lexer,
2019                                  _("Percentiles must lie in the range (0, 100)"));
2020                       goto error;
2021                     }
2022
2023                   examine.n_percentiles++;
2024                   examine.ptiles =
2025                     xrealloc (examine.ptiles,
2026                               sizeof (*examine.ptiles) *
2027                               examine.n_percentiles);
2028
2029                   examine.ptiles[examine.n_percentiles - 1] = p;
2030
2031                   lex_get (lexer);
2032                   lex_match (lexer, T_COMMA);
2033                 }
2034               if (!lex_force_match (lexer, T_RPAREN))
2035                 goto error;
2036             }
2037
2038           lex_match (lexer, T_EQUALS);
2039
2040           while (lex_token (lexer) != T_ENDCMD
2041                  && lex_token (lexer) != T_SLASH)
2042             {
2043               if (lex_match_id (lexer, "HAVERAGE"))
2044                 {
2045                   examine.pc_alg = PC_HAVERAGE;
2046                 }
2047               else if (lex_match_id (lexer, "WAVERAGE"))
2048                 {
2049                   examine.pc_alg = PC_WAVERAGE;
2050                 }
2051               else if (lex_match_id (lexer, "ROUND"))
2052                 {
2053                   examine.pc_alg = PC_ROUND;
2054                 }
2055               else if (lex_match_id (lexer, "EMPIRICAL"))
2056                 {
2057                   examine.pc_alg = PC_EMPIRICAL;
2058                 }
2059               else if (lex_match_id (lexer, "AEMPIRICAL"))
2060                 {
2061                   examine.pc_alg = PC_AEMPIRICAL;
2062                 }
2063               else if (lex_match_id (lexer, "NONE"))
2064                 {
2065                   examine.pc_alg = PC_NONE;
2066                 }
2067               else
2068                 {
2069                   lex_error (lexer, NULL);
2070                   goto error;
2071                 }
2072             }
2073         }
2074       else if (lex_match_id (lexer, "TOTAL"))
2075         {
2076           totals_seen = true;
2077         }
2078       else if (lex_match_id (lexer, "NOTOTAL"))
2079         {
2080           nototals_seen = true;
2081         }
2082       else if (lex_match_id (lexer, "MISSING"))
2083         {
2084           lex_match (lexer, T_EQUALS);
2085
2086           while (lex_token (lexer) != T_ENDCMD
2087                  && lex_token (lexer) != T_SLASH)
2088             {
2089               if (lex_match_id (lexer, "LISTWISE"))
2090                 {
2091                   examine.missing_pw = false;
2092                 }
2093               else if (lex_match_id (lexer, "PAIRWISE"))
2094                 {
2095                   examine.missing_pw = true;
2096                 }
2097               else if (lex_match_id (lexer, "EXCLUDE"))
2098                 {
2099                   examine.dep_excl = MV_ANY;
2100                 }
2101               else if (lex_match_id (lexer, "INCLUDE"))
2102                 {
2103                   examine.dep_excl = MV_SYSTEM;
2104                 }
2105               else if (lex_match_id (lexer, "REPORT"))
2106                 {
2107                   examine.fctr_excl = MV_NEVER;
2108                 }
2109               else if (lex_match_id (lexer, "NOREPORT"))
2110                 {
2111                   examine.fctr_excl = MV_ANY;
2112                 }
2113               else
2114                 {
2115                   lex_error (lexer, NULL);
2116                   goto error;
2117                 }
2118             }
2119         }
2120       else if (lex_match_id (lexer, "COMPARE"))
2121         {
2122           lex_match (lexer, T_EQUALS);
2123           if (lex_match_id (lexer, "VARIABLES"))
2124             {
2125               examine.boxplot_mode = BP_VARIABLES;
2126             }
2127           else if (lex_match_id (lexer, "GROUPS"))
2128             {
2129               examine.boxplot_mode = BP_GROUPS;
2130             }
2131           else
2132             {
2133               lex_error (lexer, NULL);
2134               goto error;
2135             }
2136         }
2137       else if (lex_match_id (lexer, "PLOT"))
2138         {
2139           lex_match (lexer, T_EQUALS);
2140
2141           while (lex_token (lexer) != T_ENDCMD
2142                  && lex_token (lexer) != T_SLASH)
2143             {
2144               if (lex_match_id (lexer, "BOXPLOT"))
2145                 {
2146                   examine.boxplot = true;
2147                 }
2148               else if (lex_match_id (lexer, "NPPLOT"))
2149                 {
2150                   examine.npplot = true;
2151                 }
2152               else if (lex_match_id (lexer, "HISTOGRAM"))
2153                 {
2154                   examine.histogram = true;
2155                 }
2156               else if (lex_match_id (lexer, "SPREADLEVEL"))
2157                 {
2158                   examine.spreadlevel = true;
2159                   examine.sl_power = 0;
2160                   if (lex_match (lexer, T_LPAREN))
2161                     {
2162                       examine.sl_power = lex_integer (lexer);
2163
2164                       lex_get (lexer);
2165                       if (! lex_force_match (lexer, T_RPAREN))
2166                         goto error;
2167                     }
2168                 }
2169               else if (lex_match_id (lexer, "NONE"))
2170                 {
2171                   examine.histogram = false;
2172                   examine.npplot = false;
2173                   examine.boxplot = false;
2174                 }
2175               else if (lex_match (lexer, T_ALL))
2176                 {
2177                   examine.histogram = true;
2178                   examine.npplot = true;
2179                   examine.boxplot = true;
2180                 }
2181               else 
2182                 {
2183                   lex_error (lexer, NULL);
2184                   goto error;
2185                 }
2186               lex_match (lexer, T_COMMA);
2187             }          
2188         }
2189       else if (lex_match_id (lexer, "CINTERVAL"))
2190         {
2191           if ( !lex_force_num (lexer))
2192             goto error;
2193         
2194           examine.conf = lex_number (lexer);
2195           lex_get (lexer);
2196         }
2197       else if (lex_match_id (lexer, "ID"))
2198         {
2199           lex_match (lexer, T_EQUALS);
2200
2201           examine.id_var = parse_variable_const (lexer, examine.dict);
2202         }
2203       else
2204         {
2205           lex_error (lexer, NULL);
2206           goto error;
2207         }
2208     }
2209
2210
2211   if ( totals_seen && nototals_seen)
2212     {
2213       msg (SE, _("%s and %s are mutually exclusive"),"TOTAL","NOTOTAL");
2214       goto error;
2215     }
2216
2217   /* If totals have been requested or if there are no factors
2218      in this analysis, then the totals need to be included. */
2219   if ( !nototals_seen || examine.n_iacts == 1)
2220     {
2221       examine.iacts = &iacts_mem[0];
2222     }
2223   else
2224     {
2225       examine.n_iacts--;
2226       examine.iacts = &iacts_mem[1];
2227       interaction_destroy (iacts_mem[0]);
2228     }
2229
2230
2231   if ( examine.id_var )
2232     {
2233       examine.id_idx = var_get_case_index (examine.id_var);
2234       examine.id_width = var_get_width (examine.id_var);
2235     }
2236
2237   examine.ex_proto = caseproto_add_width (examine.ex_proto, 0); /* value */
2238   examine.ex_proto = caseproto_add_width (examine.ex_proto, examine.id_width);   /* id */
2239   examine.ex_proto = caseproto_add_width (examine.ex_proto, 0); /* weight */
2240
2241
2242   if (examine.disp_extremes > 0)
2243     {
2244       examine.calc_extremes = examine.disp_extremes;
2245     }
2246
2247   if (examine.descriptives && examine.calc_extremes == 0)
2248     {
2249       /* Descriptives always displays the max and min */
2250       examine.calc_extremes = 1;
2251     }
2252
2253   if (percentiles_seen && examine.n_percentiles == 0)
2254     {
2255       examine.n_percentiles = 7;
2256       examine.ptiles = xcalloc (examine.n_percentiles,
2257                                 sizeof (*examine.ptiles));
2258
2259       examine.ptiles[0] = 5;
2260       examine.ptiles[1] = 10;
2261       examine.ptiles[2] = 25;
2262       examine.ptiles[3] = 50;
2263       examine.ptiles[4] = 75;
2264       examine.ptiles[5] = 90;
2265       examine.ptiles[6] = 95;
2266     }
2267
2268   assert (examine.calc_extremes >= examine.disp_extremes);
2269   {
2270     struct casegrouper *grouper;
2271     struct casereader *group;
2272     bool ok;
2273     
2274     grouper = casegrouper_create_splits (proc_open (ds), examine.dict);
2275     while (casegrouper_get_next_group (grouper, &group))
2276       run_examine (&examine, group);
2277     ok = casegrouper_destroy (grouper);
2278     ok = proc_commit (ds) && ok;
2279   }
2280
2281   caseproto_unref (examine.ex_proto);
2282
2283   free (examine.ptiles);
2284   free (examine.dep_vars);
2285   pool_destroy (examine.pool);
2286
2287   return CMD_SUCCESS;
2288
2289  error:
2290   caseproto_unref (examine.ex_proto);
2291   examine.iacts = iacts_mem;
2292   for (i = 0; i < examine.n_iacts; ++i)
2293     interaction_destroy (examine.iacts[i]);
2294   free (examine.dep_vars);
2295   free (examine.ptiles);
2296   pool_destroy (examine.pool);
2297
2298   return CMD_FAILURE;
2299 }