session: Fix two memory leaks.
[pspp] / src / language / stats / rank.c
1 /* PSPP - a program for statistical analysis.
2    Copyright (C) 2005, 2006, 2007, 2009, 2010, 2011, 2012, 2013 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 <math.h>
20 #include <gsl/gsl_cdf.h>
21
22 #include "data/case.h"
23 #include "data/casegrouper.h"
24 #include "data/casereader.h"
25 #include "data/dataset.h"
26 #include "data/dictionary.h"
27 #include "data/format.h"
28 #include "data/variable.h"
29 #include "data/subcase.h"
30 #include "data/casewriter.h"
31 #include "data/short-names.h"
32 #include "language/command.h"
33 #include "language/lexer/lexer.h"
34 #include "language/lexer/variable-parser.h"
35 #include "language/stats/sort-criteria.h"
36 #include "math/sort.h"
37 #include "libpspp/assertion.h"
38 #include "libpspp/i18n.h"
39 #include "libpspp/message.h"
40 #include "libpspp/misc.h"
41 #include "libpspp/pool.h"
42 #include "libpspp/string-set.h"
43 #include "libpspp/taint.h"
44
45 #include "output/tab.h"
46
47 #include "gettext.h"
48 #define _(msgid) gettext (msgid)
49 #define N_(msgid) (msgid)
50
51 struct rank;
52
53 typedef double (*rank_function_t) (const struct rank*, double c, double cc, double cc_1,
54                                    int i, double w);
55
56 static double rank_proportion (const struct rank *, double c, double cc, double cc_1,
57                                int i, double w);
58
59 static double rank_normal (const struct rank *, double c, double cc, double cc_1,
60                            int i, double w);
61
62 static double rank_percent (const struct rank *, double c, double cc, double cc_1,
63                             int i, double w);
64
65 static double rank_rfraction (const struct rank *, double c, double cc, double cc_1,
66                               int i, double w);
67
68 static double rank_rank (const struct rank *, double c, double cc, double cc_1,
69                          int i, double w);
70
71 static double rank_n (const struct rank *, double c, double cc, double cc_1,
72                       int i, double w);
73
74 static double rank_savage (const struct rank *, double c, double cc, double cc_1,
75                            int i, double w);
76
77 static double rank_ntiles (const struct rank *, double c, double cc, double cc_1,
78                            int i, double w);
79
80
81 enum rank_func
82   {
83     RANK,
84     NORMAL,
85     PERCENT,
86     RFRACTION,
87     PROPORTION,
88     N,
89     NTILES,
90     SAVAGE,
91     n_RANK_FUNCS
92   };
93
94 static const struct fmt_spec dest_format[n_RANK_FUNCS] = {
95   {FMT_F, 9, 3}, /* rank */
96   {FMT_F, 6, 4}, /* normal */
97   {FMT_F, 6, 2}, /* percent */
98   {FMT_F, 6, 4}, /* rfraction */
99   {FMT_F, 6, 4}, /* proportion */
100   {FMT_F, 6, 0}, /* n */
101   {FMT_F, 3, 0}, /* ntiles */
102   {FMT_F, 8, 4}  /* savage */
103 };
104
105 static const char * const function_name[n_RANK_FUNCS] = {
106   "RANK",
107   "NORMAL",
108   "PERCENT",
109   "RFRACTION",
110   "PROPORTION",
111   "N",
112   "NTILES",
113   "SAVAGE"
114 };
115
116 static const rank_function_t rank_func[n_RANK_FUNCS] = {
117   rank_rank,
118   rank_normal,
119   rank_percent,
120   rank_rfraction,
121   rank_proportion,
122   rank_n,
123   rank_ntiles,
124   rank_savage
125 };
126
127
128 enum ties
129   {
130     TIES_LOW,
131     TIES_HIGH,
132     TIES_MEAN,
133     TIES_CONDENSE
134   };
135
136 enum fraction
137   {
138     FRAC_BLOM,
139     FRAC_RANKIT,
140     FRAC_TUKEY,
141     FRAC_VW
142   };
143
144 struct rank_spec
145 {
146   enum rank_func rfunc;
147   const char **dest_names;
148   const char **dest_labels;
149 };
150
151 /* If NEW_NAME exists in DICT or NEW_NAMES, returns NULL without changing
152    anything.  Otherwise, inserts NEW_NAME in NEW_NAMES and returns the copy of
153    NEW_NAME now in NEW_NAMES. */
154 static const char *
155 try_new_name (const char *new_name,
156               const struct dictionary *dict, struct string_set *new_names)
157 {
158   return (!dict_lookup_var (dict, new_name)
159           && string_set_insert (new_names, new_name)
160           ? string_set_find_node (new_names, new_name)->string
161           : NULL);
162 }
163
164 /* Returns a variable name for storing ranks of a variable named SRC_NAME
165    accoring to the rank function F.  The name chosen will not be one already in
166    DICT or NEW_NAMES.
167
168    If successful, adds the new name to NEW_NAMES and returns the name added.
169    If no name can be generated, returns NULL. */
170 static const char *
171 rank_choose_dest_name (struct dictionary *dict, struct string_set *new_names,
172                        enum rank_func f, const char *src_name)
173 {
174   char *src_name_7;
175   char name[128];
176   const char *s;
177   int i;
178
179   /* Try the first character of the ranking function followed by the first 7
180      bytes of the srcinal variable name. */
181   src_name_7 = utf8_encoding_trunc (src_name, dict_get_encoding (dict), 7);
182   snprintf (name, sizeof name, "%c%s", function_name[f][0], src_name_7);
183   free (src_name_7);
184   s = try_new_name (name, dict, new_names);
185   if (s != NULL)
186     return s;
187
188   /* Try "fun###". */
189   for (i = 1; i <= 999; i++)
190     {
191       sprintf (name, "%.3s%03d", function_name[f], i);
192       s = try_new_name (name, dict, new_names);
193       if (s != NULL)
194         return s;
195     }
196
197   /* Try "RNKfn##". */
198   for (i = 1; i <= 99; i++)
199     {
200       sprintf (name, "RNK%.2s%02d", function_name[f], i);
201       s = try_new_name (name, dict, new_names);
202       if (s != NULL)
203         return s;
204     }
205
206   msg (ME, _("Cannot generate variable name for ranking %s with %s.  "
207              "All candidates in use."),
208        src_name, function_name[f]);
209   return NULL;
210 }
211
212 struct rank
213 {
214   struct dictionary *dict;
215
216   struct subcase sc;
217
218   const struct variable **vars;
219   size_t n_vars;
220
221   const struct variable **group_vars;
222   size_t n_group_vars;
223
224
225   enum mv_class exclude;
226
227   struct rank_spec *rs;
228   size_t n_rs;
229
230   enum ties ties;
231
232   enum fraction fraction;
233   int k_ntiles;
234
235   bool print;
236
237   /* Pool on which cell functions may allocate data */
238   struct pool *pool;
239 };
240
241
242 static void
243 destroy_rank (struct rank *rank)
244 {
245  free (rank->vars);
246  free (rank->group_vars);
247  subcase_destroy (&rank->sc);
248  pool_destroy (rank->pool);
249 }
250
251 static bool
252 parse_into (struct lexer *lexer, struct rank *cmd,
253             struct string_set *new_names)
254 {
255   int var_count = 0;
256   struct rank_spec *rs = NULL;
257
258   cmd->rs = pool_realloc (cmd->pool, cmd->rs, sizeof (*cmd->rs) * (cmd->n_rs + 1));
259   rs = &cmd->rs[cmd->n_rs];
260       
261   if (lex_match_id (lexer, "RANK"))
262     {
263       rs->rfunc = RANK;
264     }
265   else if (lex_match_id (lexer, "NORMAL"))
266     {
267       rs->rfunc = NORMAL;
268     }
269   else if (lex_match_id (lexer, "RFRACTION"))
270     {
271       rs->rfunc = RFRACTION;
272     }
273   else if (lex_match_id (lexer, "N"))
274     {
275       rs->rfunc = N;
276     }
277   else if (lex_match_id (lexer, "SAVAGE"))
278     {
279       rs->rfunc = SAVAGE;
280     }
281   else if (lex_match_id (lexer, "PERCENT"))
282     {
283       rs->rfunc = PERCENT;
284     }
285   else if (lex_match_id (lexer, "PROPORTION"))
286     {
287       rs->rfunc = PROPORTION;
288     }
289   else if (lex_match_id (lexer, "NTILES"))
290     {
291       if ( !lex_force_match (lexer, T_LPAREN))
292         return false;
293       
294       if (! lex_force_int (lexer) )
295         return false;
296       
297       cmd->k_ntiles = lex_integer (lexer);
298       lex_get (lexer);
299       
300       if ( !lex_force_match (lexer, T_RPAREN))
301         return false;
302
303       rs->rfunc = NTILES;
304     }
305   else
306     {
307       return false;
308     }
309
310   cmd->n_rs++;
311   rs->dest_names = pool_calloc (cmd->pool, cmd->n_vars,
312                                 sizeof *rs->dest_names);
313
314   if (lex_match_id (lexer, "INTO"))
315     {
316       while( lex_token (lexer) == T_ID )
317         {
318           const char *name = lex_tokcstr (lexer);
319
320           if ( var_count >= subcase_get_n_fields (&cmd->sc) )
321             msg (SE, _("Too many variables in INTO clause."));
322           else if ( dict_lookup_var (cmd->dict, name) != NULL )
323             msg (SE, _("Variable %s already exists."), name);
324           else if (string_set_contains (new_names, name))
325             msg (SE, _("Duplicate variable name %s."), name);
326           else
327             {
328               string_set_insert (new_names, name);
329               rs->dest_names[var_count++] = pool_strdup (cmd->pool, name);
330               lex_get (lexer);
331               continue;
332             }
333
334           /* Error path. */
335           return false;
336         }
337     }
338
339   return true;
340 }
341
342 /* Hardly a rank function !! */
343 static double
344 rank_n (const struct rank *cmd UNUSED, double c UNUSED, double cc UNUSED, double cc_1 UNUSED,
345         int i UNUSED, double w)
346 {
347   return w;
348 }
349
350
351 static double
352 rank_rank (const struct rank *cmd, double c, double cc, double cc_1,
353            int i, double w UNUSED)
354 {
355   double rank;
356
357   if ( c >= 1.0 )
358     {
359       switch (cmd->ties)
360         {
361         case TIES_LOW:
362           rank = cc_1 + 1;
363           break;
364         case TIES_HIGH:
365           rank = cc;
366           break;
367         case TIES_MEAN:
368           rank = cc_1 + (c + 1.0)/ 2.0;
369           break;
370         case TIES_CONDENSE:
371           rank = i;
372           break;
373         default:
374           NOT_REACHED ();
375         }
376     }
377   else
378     {
379       switch (cmd->ties)
380         {
381         case TIES_LOW:
382           rank = cc_1;
383           break;
384         case TIES_HIGH:
385           rank = cc;
386           break;
387         case TIES_MEAN:
388           rank = cc_1 + c / 2.0 ;
389           break;
390         case TIES_CONDENSE:
391           rank = i;
392           break;
393         default:
394           NOT_REACHED ();
395         }
396     }
397
398   return rank;
399 }
400
401
402 static double
403 rank_rfraction (const struct rank *cmd, double c, double cc, double cc_1,
404                 int i, double w)
405 {
406   return rank_rank (cmd, c, cc, cc_1, i, w) / w ;
407 }
408
409
410 static double
411 rank_percent (const struct rank *cmd, double c, double cc, double cc_1,
412               int i, double w)
413 {
414   return rank_rank (cmd, c, cc, cc_1, i, w) * 100.0 / w ;
415 }
416
417
418 static double
419 rank_proportion (const struct rank *cmd, double c, double cc, double cc_1,
420                  int i, double w)
421 {
422   const double r =  rank_rank (cmd, c, cc, cc_1, i, w) ;
423
424   double f;
425
426   switch ( cmd->fraction )
427     {
428     case FRAC_BLOM:
429       f =  (r - 3.0/8.0) / (w + 0.25);
430       break;
431     case FRAC_RANKIT:
432       f = (r - 0.5) / w ;
433       break;
434     case FRAC_TUKEY:
435       f = (r - 1.0/3.0) / (w + 1.0/3.0);
436       break;
437     case FRAC_VW:
438       f = r / ( w + 1.0);
439       break;
440     default:
441       NOT_REACHED ();
442     }
443
444
445   return (f > 0) ? f : SYSMIS;
446 }
447
448 static double
449 rank_normal (const struct rank *cmd, double c, double cc, double cc_1,
450              int i, double w)
451 {
452   double f = rank_proportion (cmd, c, cc, cc_1, i, w);
453
454   return gsl_cdf_ugaussian_Pinv (f);
455 }
456
457 static double
458 rank_ntiles (const struct rank *cmd, double c, double cc, double cc_1,
459              int i, double w)
460 {
461   double r = rank_rank (cmd, c, cc, cc_1, i, w);
462
463
464   return ( floor (( r * cmd->k_ntiles) / ( w + 1) ) + 1);
465 }
466
467 /* Expected value of the order statistics from an exponential distribution */
468 static double
469 ee (int j, double w_star)
470 {
471   int k;
472   double sum = 0.0;
473
474   for (k = 1 ; k <= j; k++)
475     sum += 1.0 / ( w_star + 1 - k );
476
477   return sum;
478 }
479
480
481 static double
482 rank_savage (const struct rank *cmd UNUSED, double c, double cc, double cc_1,
483              int i UNUSED, double w)
484 {
485   double int_part;
486   const int i_1 = floor (cc_1);
487   const int i_2 = floor (cc);
488
489   const double w_star = (modf (w, &int_part) == 0 ) ? w : floor (w) + 1;
490
491   const double g_1 = cc_1 - i_1;
492   const double g_2 = cc - i_2;
493
494   /* The second factor is infinite, when the first is zero.
495      Therefore, evaluate the second, only when the first is non-zero */
496   const double expr1 =  (1 - g_1) ? (1 - g_1) * ee(i_1+1, w_star) : ( 1 - g_1);
497   const double expr2 =  g_2 ? g_2 * ee (i_2+1, w_star) : g_2 ;
498
499   if ( i_1 == i_2 )
500     return ee (i_1 + 1, w_star) - 1;
501
502   if ( i_1 + 1 == i_2 )
503     return ( ( expr1 + expr2 )/c ) - 1;
504
505   if ( i_1 + 2 <= i_2 )
506     {
507       int j;
508       double sigma = 0.0;
509       for (j = i_1 + 2 ; j <= i_2; ++j )
510         sigma += ee (j, w_star);
511       return ( (expr1 + expr2 + sigma) / c) -1;
512     }
513
514   NOT_REACHED();
515 }
516
517 static double
518 sum_weights (const struct casereader *input, int weight_idx)
519 {
520   if (weight_idx == -1)
521     return casereader_count_cases (input);
522   else
523     {
524       struct casereader *pass;
525       struct ccase *c;
526       double w;
527
528       w = 0.0;
529       pass = casereader_clone (input);
530       for (; (c = casereader_read (pass)) != NULL; case_unref (c))
531         w += case_num_idx (c, weight_idx);
532       casereader_destroy (pass);
533
534       return w;
535     }
536 }
537
538 static void
539 rank_sorted_file (struct casereader *input,
540                   struct casewriter *output,
541                   int weight_idx,
542                   const struct rank *cmd)
543 {
544   struct casegrouper *tie_grouper;
545   struct casereader *tied_cases;
546   struct subcase input_var;
547   int tie_group = 1;
548   struct ccase *c;
549   double cc = 0.0;
550   double w;
551
552   /* Get total group weight. */
553   w = sum_weights (input, weight_idx);
554
555   /* Do ranking. */
556   subcase_init (&input_var, 0, 0, SC_ASCEND);
557   tie_grouper = casegrouper_create_subcase (input, &input_var);
558   subcase_destroy (&input_var);
559   for (; casegrouper_get_next_group (tie_grouper, &tied_cases);
560        casereader_destroy (tied_cases))
561     {
562       double tw = sum_weights (tied_cases, weight_idx);
563       double cc_1 = cc;
564       cc += tw;
565
566       taint_propagate (casereader_get_taint (tied_cases),
567                        casewriter_get_taint (output));
568
569       /* Rank tied cases. */
570       for (; (c = casereader_read (tied_cases)) != NULL; case_unref (c))
571         {
572           struct ccase *out_case;
573           size_t i;
574
575           out_case = case_create (casewriter_get_proto (output));
576           case_data_rw_idx (out_case, 0)->f = case_num_idx (c, 1);
577           for (i = 0; i < cmd->n_rs; ++i)
578             {
579               rank_function_t func = rank_func[cmd->rs[i].rfunc];
580               double rank = func (cmd, tw, cc, cc_1, tie_group, w);
581               case_data_rw_idx (out_case, i + 1)->f = rank;
582             }
583
584           casewriter_write (output, out_case);
585         }
586       tie_group++;
587     }
588   casegrouper_destroy (tie_grouper);
589 }
590
591
592 static bool
593 rank_cmd (struct dataset *ds,  const struct rank *cmd);
594
595 static const char *
596 fraction_name (const struct rank *cmd)
597 {
598   switch (cmd->fraction )
599     {
600     case FRAC_BLOM:   return "BLOM";
601     case FRAC_RANKIT: return "RANKIT";
602     case FRAC_TUKEY:  return "TUKEY";
603     case FRAC_VW:     return "VW";
604     default:          NOT_REACHED ();
605     }
606 }
607
608 /* Returns a label for a variable derived from SRC_VAR with function F. */
609 static const char *
610 create_var_label (struct rank *cmd, const struct variable *src_var,
611                   enum rank_func f)
612 {
613   struct string label;
614   const char *pool_label;
615
616   ds_init_empty (&label);
617
618   if ( cmd->n_group_vars > 0 )
619     {
620       struct string group_var_str;
621       int g;
622
623       ds_init_empty (&group_var_str);
624
625       for (g = 0 ; g < cmd->n_group_vars ; ++g )
626         {
627           if ( g > 0 ) ds_put_cstr (&group_var_str, " ");
628           ds_put_cstr (&group_var_str, var_get_name (cmd->group_vars[g]));
629         }
630
631       ds_put_format (&label, _("%s of %s by %s"), function_name[f],
632                      var_get_name (src_var), ds_cstr (&group_var_str));
633       ds_destroy (&group_var_str);
634     }
635   else
636     ds_put_format (&label, _("%s of %s"),
637                    function_name[f], var_get_name (src_var));
638
639   pool_label = pool_strdup (cmd->pool, ds_cstr (&label));
640   
641   ds_destroy (&label);
642
643   return pool_label;
644 }
645
646 int
647 cmd_rank (struct lexer *lexer, struct dataset *ds)
648 {
649   struct string_set new_names;
650   struct rank rank;
651   struct rank_spec *rs;
652   int i;
653
654   subcase_init_empty (&rank.sc);
655
656   rank.rs = NULL;
657   rank.n_rs = 0;
658   rank.exclude = MV_ANY;
659   rank.n_group_vars = 0;
660   rank.group_vars = NULL;
661   rank.dict = dataset_dict (ds);
662   rank.ties = TIES_MEAN;
663   rank.fraction = FRAC_BLOM;
664   rank.print = true;
665   rank.pool = pool_create ();
666
667   string_set_init (&new_names);
668
669   if (lex_match_id (lexer, "VARIABLES"))
670     lex_force_match (lexer, T_EQUALS);
671
672   if (!parse_sort_criteria (lexer, rank.dict,
673                             &rank.sc,
674                             &rank.vars, NULL))
675     goto error;
676
677   rank.n_vars = rank.sc.n_fields;
678
679   if (lex_match (lexer, T_BY) )
680     {
681       if ( ! parse_variables_const (lexer, rank.dict,
682                                     &rank.group_vars, &rank.n_group_vars,
683                                     PV_NO_DUPLICATE | PV_NO_SCRATCH))
684         goto error;
685     }
686
687
688   while (lex_token (lexer) != T_ENDCMD )
689     {
690       lex_force_match (lexer, T_SLASH);
691       if (lex_match_id (lexer, "TIES"))
692         {
693           lex_force_match (lexer, T_EQUALS);
694           if (lex_match_id (lexer, "MEAN"))
695             {
696               rank.ties = TIES_MEAN;
697             }
698           else if (lex_match_id (lexer, "LOW"))
699             {
700               rank.ties = TIES_LOW;
701             }
702           else if (lex_match_id (lexer, "HIGH"))
703             {
704               rank.ties = TIES_HIGH;
705             }
706           else if (lex_match_id (lexer, "CONDENSE"))
707             {
708               rank.ties = TIES_CONDENSE;
709             }
710           else
711             {
712               lex_error (lexer, NULL);
713               goto error;
714             }
715         }
716       else if (lex_match_id (lexer, "FRACTION"))
717         {
718           lex_force_match (lexer, T_EQUALS);
719           if (lex_match_id (lexer, "BLOM"))
720             {
721               rank.fraction = FRAC_BLOM;
722             }
723           else if (lex_match_id (lexer, "TUKEY"))
724             {
725               rank.fraction = FRAC_TUKEY;
726             }
727           else if (lex_match_id (lexer, "VW"))
728             {
729               rank.fraction = FRAC_VW;
730             }
731           else if (lex_match_id (lexer, "RANKIT"))
732             {
733               rank.fraction = FRAC_RANKIT;
734             }
735           else
736             {
737               lex_error (lexer, NULL);
738               goto error;
739             }
740         }
741       else if (lex_match_id (lexer, "PRINT"))
742         {
743           lex_force_match (lexer, T_EQUALS);
744           if (lex_match_id (lexer, "YES"))
745             {
746               rank.print = true;
747             }
748           else if (lex_match_id (lexer, "NO"))
749             {
750               rank.print = false;
751             }
752           else
753             {
754               lex_error (lexer, NULL);
755               goto error;
756             }
757         }
758       else if (lex_match_id (lexer, "MISSING"))
759         {
760           lex_force_match (lexer, T_EQUALS);
761           if (lex_match_id (lexer, "INCLUDE"))
762             {
763               rank.exclude = MV_SYSTEM;
764             }
765           else if (lex_match_id (lexer, "EXCLUDE"))
766             {
767               rank.exclude = MV_ANY;
768             }
769           else
770             {
771               lex_error (lexer, NULL);
772               goto error;
773             }
774         }
775       else if (! parse_into (lexer, &rank, &new_names))
776         goto error;
777     }
778
779
780   /* If no rank specs are given, then apply a default */
781   if ( rank.n_rs == 0)
782     {
783       struct rank_spec *rs;
784
785       rs = pool_calloc (rank.pool, 1, sizeof *rs);
786       rs->rfunc = RANK;
787       rs->dest_names = pool_calloc (rank.pool, rank.n_vars,
788                                     sizeof *rs->dest_names);
789
790       rank.rs = rs;
791       rank.n_rs = 1;
792     }
793
794   /* Choose variable names for all rank destinations which haven't already been
795      created with INTO. */
796   for (rs = rank.rs; rs < &rank.rs[rank.n_rs]; rs++)
797     {
798       int v;
799
800       rs->dest_labels = pool_calloc (rank.pool, rank.n_vars,
801                                      sizeof *rs->dest_labels);
802       for ( v = 0 ; v < rank.n_vars ;  v ++ )
803         {
804           const char **dst_name = &rs->dest_names[v];
805           if ( *dst_name == NULL )
806             {
807               *dst_name = rank_choose_dest_name (rank.dict, &new_names,
808                                                  rs->rfunc,
809                                                  var_get_name (rank.vars[v]));
810               if (*dst_name == NULL)
811                 goto error;
812             }
813
814           rs->dest_labels[v] = create_var_label (&rank, rank.vars[v],
815                                                  rs->rfunc);
816         }
817     }
818
819   if ( rank.print )
820     {
821       int v;
822
823       tab_output_text (0, _("Variables Created By RANK"));
824       tab_output_text (0, "");
825
826       for (i = 0 ; i <  rank.n_rs ; ++i )
827         {
828           for ( v = 0 ; v < rank.n_vars ;  v ++ )
829             {
830               if ( rank.n_group_vars > 0 )
831                 {
832                   struct string varlist;
833                   int g;
834
835                   ds_init_empty (&varlist);
836                   for ( g = 0 ; g < rank.n_group_vars ; ++g )
837                     {
838                       ds_put_cstr (&varlist, var_get_name (rank.group_vars[g]));
839
840                       if ( g < rank.n_group_vars - 1)
841                         ds_put_cstr (&varlist, " ");
842                     }
843
844                   if ( rank.rs[i].rfunc == NORMAL ||
845                        rank.rs[i].rfunc == PROPORTION )
846                     tab_output_text_format (0,
847                                             _("%s into %s(%s of %s using %s BY %s)"),
848                                             var_get_name (rank.vars[v]),
849                                             rank.rs[i].dest_names[v],
850                                             function_name[rank.rs[i].rfunc],
851                                             var_get_name (rank.vars[v]),
852                                             fraction_name (&rank),
853                                             ds_cstr (&varlist));
854
855                   else
856                     tab_output_text_format (0,
857                                             _("%s into %s(%s of %s BY %s)"),
858                                             var_get_name (rank.vars[v]),
859                                             rank.rs[i].dest_names[v],
860                                             function_name[rank.rs[i].rfunc],
861                                             var_get_name (rank.vars[v]),
862                                             ds_cstr (&varlist));
863                   ds_destroy (&varlist);
864                 }
865               else
866                 {
867                   if ( rank.rs[i].rfunc == NORMAL ||
868                        rank.rs[i].rfunc == PROPORTION )
869                     tab_output_text_format (0,
870                                             _("%s into %s(%s of %s using %s)"),
871                                             var_get_name (rank.vars[v]),
872                                             rank.rs[i].dest_names[v],
873                                             function_name[rank.rs[i].rfunc],
874                                             var_get_name (rank.vars[v]),
875                                             fraction_name (&rank));
876
877                   else
878                     tab_output_text_format (0,
879                                             _("%s into %s(%s of %s)"),
880                                             var_get_name (rank.vars[v]),
881                                             rank.rs[i].dest_names[v],
882                                             function_name[rank.rs[i].rfunc],
883                                             var_get_name (rank.vars[v]));
884                 }
885             }
886         }
887     }
888
889   /* Do the ranking */
890   rank_cmd (ds, &rank);
891
892   destroy_rank (&rank);
893   string_set_destroy (&new_names);
894   return CMD_SUCCESS;
895
896  error:
897
898   destroy_rank (&rank);
899   string_set_destroy (&new_names);
900   return CMD_FAILURE;
901 }
902
903 /* RANK transformation. */
904 struct rank_trns
905   {
906     int order_case_idx;
907
908     struct rank_trns_input_var *input_vars;
909     size_t n_input_vars;
910
911     size_t n_funcs;
912   };
913
914 struct rank_trns_input_var
915   {
916     struct casereader *input;
917     struct ccase *current;
918
919     struct variable **output_vars;
920   };
921
922 static void
923 advance_ranking (struct rank_trns_input_var *iv)
924 {
925   case_unref (iv->current);
926   iv->current = casereader_read (iv->input);
927 }
928
929 static int
930 rank_trns_proc (void *trns_, struct ccase **c, casenumber case_idx UNUSED)
931 {
932   struct rank_trns *trns = trns_;
933   double order = case_num_idx (*c, trns->order_case_idx);
934   struct rank_trns_input_var *iv;
935
936   *c = case_unshare (*c);
937   for (iv = trns->input_vars; iv < &trns->input_vars[trns->n_input_vars]; iv++)
938     while (iv->current != NULL)
939       {
940         double iv_order = case_num_idx (iv->current, 0);
941         if (iv_order == order)
942           {
943             size_t i;
944
945             for (i = 0; i < trns->n_funcs; i++)
946               case_data_rw (*c, iv->output_vars[i])->f
947                 = case_num_idx (iv->current, i + 1);
948             advance_ranking (iv);
949             break;
950           }
951         else if (iv_order > order)
952           break;
953         else
954           advance_ranking (iv);
955       }
956   return TRNS_CONTINUE;
957 }
958
959 static bool
960 rank_trns_free (void *trns_)
961 {
962   struct rank_trns *trns = trns_;
963   struct rank_trns_input_var *iv;
964
965   for (iv = trns->input_vars; iv < &trns->input_vars[trns->n_input_vars]; iv++)
966     {
967       casereader_destroy (iv->input);
968       case_unref (iv->current);
969
970       free (iv->output_vars);
971     }
972   free (trns->input_vars);
973   free (trns);
974
975   return true;
976 }
977
978 static bool
979 rank_cmd (struct dataset *ds, const struct rank *cmd)
980 {
981   struct dictionary *d = dataset_dict (ds);
982   struct variable *weight_var = dict_get_weight (d);
983   struct casewriter **outputs;
984   struct variable *order_var;
985   struct casereader *input;
986   struct rank_trns *trns;
987   bool ok = true;
988   int i;
989
990   order_var = add_permanent_ordering_transformation (ds);
991
992   /* Create output files. */
993   {
994     struct caseproto *output_proto;
995     struct subcase by_order;
996
997     output_proto = caseproto_create ();
998     for (i = 0; i < cmd->n_rs + 1; i++)
999       output_proto = caseproto_add_width (output_proto, 0);
1000
1001     subcase_init (&by_order, 0, 0, SC_ASCEND);
1002
1003     outputs = xnmalloc (cmd->n_vars, sizeof *outputs);
1004     for (i = 0; i < cmd->n_vars; i++)
1005       outputs[i] = sort_create_writer (&by_order, output_proto);
1006
1007     subcase_destroy (&by_order);
1008     caseproto_unref (output_proto);
1009   }
1010
1011   /* Open the active file and make one pass per input variable. */
1012   input = proc_open (ds);
1013   input = casereader_create_filter_weight (input, d, NULL, NULL);
1014   for (i = 0 ; i < cmd->n_vars ; ++i )
1015     {
1016       const struct variable *input_var = cmd->vars[i];
1017       struct casereader *input_pass;
1018       struct casegrouper *split_grouper;
1019       struct casereader *split_group;
1020       struct subcase rank_ordering;
1021       struct subcase projection;
1022       struct subcase split_vars;
1023       struct subcase group_vars;
1024       int weight_idx;
1025       int j;
1026
1027       /* Discard cases that have missing values of input variable. */
1028       input_pass = i == cmd->n_vars - 1 ? input : casereader_clone (input);
1029       input_pass = casereader_create_filter_missing (input_pass, &input_var, 1,
1030                                                      cmd->exclude, NULL, NULL);
1031
1032       /* Keep only the columns we really need, to save time and space when we
1033          sort them just below.
1034
1035          After this projection, the input_pass case indexes look like:
1036
1037            - 0: input_var.
1038            - 1: order_var.
1039            - 2 and up: cmd->n_group_vars group variables
1040            - 2 + cmd->n_group_vars and up: split variables
1041            - 2 + cmd->n_group_vars + n_split_vars: weight var
1042       */
1043       subcase_init_empty (&projection);
1044       subcase_add_var_always (&projection, input_var, SC_ASCEND);
1045       subcase_add_var_always (&projection, order_var, SC_ASCEND);
1046       subcase_add_vars_always (&projection,
1047                                cmd->group_vars, cmd->n_group_vars);
1048       subcase_add_vars_always (&projection, dict_get_split_vars (d),
1049                                dict_get_split_cnt (d));
1050       if (weight_var != NULL)
1051         {
1052           subcase_add_var_always (&projection, weight_var, SC_ASCEND);
1053           weight_idx = 2 + cmd->n_group_vars + dict_get_split_cnt (d);
1054         }
1055       else
1056         weight_idx = -1;
1057       input_pass = casereader_project (input_pass, &projection);
1058       subcase_destroy (&projection);
1059
1060       /* Prepare 'group_vars' as the set of grouping variables. */
1061       subcase_init_empty (&group_vars);
1062       for (j = 0; j < cmd->n_group_vars; j++)
1063         subcase_add_always (&group_vars,
1064                             j + 2, var_get_width (cmd->group_vars[j]),
1065                             SC_ASCEND);
1066
1067       /* Prepare 'rank_ordering' for sorting with the group variables as
1068          primary key and the input variable as secondary key. */
1069       subcase_clone (&rank_ordering, &group_vars);
1070       subcase_add (&rank_ordering, 0, 0, subcase_get_direction (&cmd->sc, i));
1071
1072       /* Group by split variables */
1073       subcase_init_empty (&split_vars);
1074       for (j = 0; j < dict_get_split_cnt (d); j++)
1075         subcase_add_always (&split_vars, 2 + j + cmd->n_group_vars,
1076                             var_get_width (dict_get_split_vars (d)[j]),
1077                             SC_ASCEND);
1078       split_grouper = casegrouper_create_subcase (input_pass, &split_vars);
1079       subcase_destroy (&split_vars);
1080       while (casegrouper_get_next_group (split_grouper, &split_group))
1081         {
1082           struct casereader *ordered;
1083           struct casegrouper *by_grouper;
1084           struct casereader *by_group;
1085
1086           ordered = sort_execute (split_group, &rank_ordering);
1087           by_grouper = casegrouper_create_subcase (ordered, &group_vars);
1088           while (casegrouper_get_next_group (by_grouper, &by_group))
1089             rank_sorted_file (by_group, outputs[i], weight_idx, cmd);
1090           ok = casegrouper_destroy (by_grouper) && ok;
1091         }
1092       subcase_destroy (&group_vars);
1093       subcase_destroy (&rank_ordering);
1094
1095       ok = casegrouper_destroy (split_grouper) && ok;
1096     }
1097   ok = proc_commit (ds) && ok;
1098
1099   /* Re-fetch the dictionary and order variable, because if TEMPORARY was in
1100      effect then there's a new dictionary. */
1101   d = dataset_dict (ds);
1102   order_var = dict_lookup_var_assert (d, "$ORDER");
1103
1104   /* Merge the original data set with the ranks (which we already sorted on
1105      $ORDER). */
1106   trns = xmalloc (sizeof *trns);
1107   trns->order_case_idx = var_get_case_index (order_var);
1108   trns->input_vars = xnmalloc (cmd->n_vars, sizeof *trns->input_vars);
1109   trns->n_input_vars = cmd->n_vars;
1110   trns->n_funcs = cmd->n_rs;
1111   for (i = 0; i < trns->n_input_vars; i++)
1112     {
1113       struct rank_trns_input_var *iv = &trns->input_vars[i];
1114       int j;
1115
1116       iv->input = casewriter_make_reader (outputs[i]);
1117       iv->current = casereader_read (iv->input);
1118       iv->output_vars = xnmalloc (trns->n_funcs, sizeof *iv->output_vars);
1119       for (j = 0; j < trns->n_funcs; j++)
1120         {
1121           struct rank_spec *rs = &cmd->rs[j];
1122           struct variable *var;
1123
1124           var = dict_create_var_assert (d, rs->dest_names[i], 0);
1125           var_set_both_formats (var, &dest_format[rs->rfunc]);
1126           var_set_label (var, rs->dest_labels[i], false);
1127
1128           iv->output_vars[j] = var;
1129         }
1130     }
1131   free (outputs);
1132
1133   add_transformation (ds, rank_trns_proc, rank_trns_free, trns);
1134
1135   /* Delete our sort key, which we don't need anymore. */
1136   dict_delete_var (d, order_var);
1137
1138   return ok;
1139 }