Matrix reader, fix comment to reflect reality.
[pspp] / src / language / data-io / matrix-reader.c
1 /* PSPP - a program for statistical analysis.
2    Copyright (C) 2017 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 "matrix-reader.h"
20
21 #include <stdbool.h>
22
23 #include <libpspp/message.h>
24 #include <data/casegrouper.h>
25 #include <data/casereader.h>
26 #include <data/dictionary.h>
27 #include <data/variable.h>
28
29 #include "gettext.h"
30 #define _(msgid) gettext (msgid)
31 #define N_(msgid) msgid
32
33
34 /*
35 This module interprets a "data matrix", typically generated by the command
36 MATRIX DATA.  The dictionary of such a matrix takes the form:
37
38  s_0, s_1, ... s_m, ROWTYPE_, VARNAME_, v_0, v_1, .... v_n
39
40 where s_0, s_1 ... s_m are the variables defining the splits, and
41 v_0, v_1 ... v_n are the continuous variables.
42
43 m >= 0; n >= 0
44
45 The ROWTYPE_ variable is of type A8.
46 The VARNAME_ variable is a string type whose width is not predetermined.
47 The variables s_x are of type F4.0 (although this reader accepts any type),
48 and v_x are of any numeric type.
49
50 The values of the ROWTYPE_ variable are in the set {MEAN, STDDEV, N, CORR, COV}
51 and determine the purpose of that case.
52 The values of the VARNAME_ variable must correspond to the names of the varibles
53 in {v_0, v_1 ... v_n} and indicate the rows of the correlation or covariance
54 matrices.
55
56
57
58 A typical example is as follows:
59
60 s_0 ROWTYPE_   VARNAME_   v_0         v_1         v_2
61
62 0   MEAN                5.0000       4.0000       3.0000
63 0   STDDEV              1.0000       2.0000       3.0000
64 0   N                   9.0000       9.0000       9.0000
65 0   CORR       V1       1.0000        .6000        .7000
66 0   CORR       V2        .6000       1.0000        .8000
67 0   CORR       V3        .7000        .8000       1.0000
68 1   MEAN                9.0000       8.0000       7.0000
69 1   STDDEV              5.0000       6.0000       7.0000
70 1   N                   9.0000       9.0000       9.0000
71 1   CORR       V1       1.0000        .4000        .3000
72 1   CORR       V2        .4000       1.0000        .2000
73 1   CORR       V3        .3000        .2000       1.0000
74
75 */
76
77 struct matrix_reader
78 {
79   const struct dictionary *dict;
80   const struct variable *varname;
81   const struct variable *rowtype;
82   struct casegrouper *grouper;
83
84   gsl_matrix *n_vectors;
85   gsl_matrix *mean_vectors;
86   gsl_matrix *var_vectors;
87
88   gsl_matrix *correlation;
89   gsl_matrix *covariance;
90 };
91
92 struct matrix_reader *
93 create_matrix_reader_from_case_reader (const struct dictionary *dict, struct casereader *in_reader,
94                                        const struct variable ***vars, size_t *n_vars)
95 {
96   struct matrix_reader *mr = xzalloc (sizeof *mr);
97
98   mr->dict = dict;
99   mr->varname = dict_lookup_var (dict, "varname_");
100   if (mr->varname == NULL)
101     {
102       msg (ME, _("Matrix dataset lacks a variable called %s."), "VARNAME_");
103       free (mr);
104       return NULL;
105     }
106
107   mr->rowtype = dict_lookup_var (dict, "rowtype_");
108   if (mr->rowtype == NULL)
109     {
110       msg (ME, _("Matrix dataset lacks a variable called %s."), "ROWTYPE_");
111       free (mr);
112       return NULL;
113     }
114
115   size_t dvarcnt;
116   const struct variable **dvars = NULL;
117   dict_get_vars (dict, &dvars, &dvarcnt, DC_SCRATCH);
118
119   if (n_vars)
120     *n_vars = dvarcnt - var_get_dict_index (mr->varname) - 1;
121
122   if (vars)
123     {
124       int i;
125       *vars = xcalloc (sizeof (struct variable **), *n_vars);
126
127       for (i = 0; i < *n_vars; ++i)
128         {
129           (*vars)[i] = dvars[i + var_get_dict_index (mr->varname) + 1];
130         }
131     }
132
133   /* All the variables before ROWTYPE_ (if any) are split variables */
134   mr->grouper = casegrouper_create_vars (in_reader, dvars, var_get_dict_index (mr->rowtype));
135
136   free (dvars);
137
138   return mr;
139 }
140
141 bool
142 destroy_matrix_reader (struct matrix_reader *mr)
143 {
144   if (mr == NULL)
145     return false;
146   bool ret = casegrouper_destroy (mr->grouper);
147   free (mr);
148   return ret;
149 }
150
151
152 bool
153 next_matrix_from_reader (struct matrix_material *mm,
154                          struct matrix_reader *mr,
155                          const struct variable **vars, int n_vars)
156 {
157   struct casereader *group;
158
159   gsl_matrix_free (mr->n_vectors);
160   gsl_matrix_free (mr->mean_vectors);
161   gsl_matrix_free (mr->var_vectors);
162   gsl_matrix_free (mr->correlation);
163   gsl_matrix_free (mr->covariance);
164
165   if (!casegrouper_get_next_group (mr->grouper, &group))
166     return false;
167
168   mr->n_vectors    = gsl_matrix_alloc (n_vars, n_vars);
169   mr->mean_vectors = gsl_matrix_alloc (n_vars, n_vars);
170   mr->var_vectors  = gsl_matrix_alloc (n_vars, n_vars);
171
172   mm->n = mr->n_vectors;
173   mm->mean_matrix = mr->mean_vectors;
174   mm->var_matrix = mr->var_vectors;
175
176   mr->correlation  = NULL;
177   mr->covariance   = NULL;
178
179   struct ccase *c;
180   int crow = 0;
181   for ( ; (c = casereader_read (group) ); case_unref (c))
182     {
183       const union value *uv  = case_data (c, mr->rowtype);
184       int col, row;
185       for (col = 0; col < n_vars; ++col)
186         {
187           const struct variable *cv
188             = vars ? vars[col] : dict_get_var (mr->dict, var_get_dict_index (mr->varname) + 1 + col);
189           double x = case_data (c, cv)->f;
190           if (0 == strncasecmp ((char *)value_str (uv, 8), "N       ", 8))
191             for (row = 0; row < n_vars; ++row)
192               gsl_matrix_set (mr->n_vectors, row, col, x);
193           else if (0 == strncasecmp ((char *) value_str (uv, 8), "MEAN    ", 8))
194             for (row = 0; row < n_vars; ++row)
195               gsl_matrix_set (mr->mean_vectors, row, col, x);
196           else if (0 == strncasecmp ((char *) value_str (uv, 8), "STDDEV  ", 8))
197             for (row = 0; row < n_vars; ++row)
198               gsl_matrix_set (mr->var_vectors, row, col, x * x);
199         }
200       if (0 == strncasecmp ((char *) value_str (uv, 8), "CORR    ", 8))
201         {
202           if (mr->correlation == NULL)
203             mr->correlation  = gsl_matrix_alloc (n_vars, n_vars);
204           for (col = 0; col < n_vars; ++col)
205             {
206               const struct variable *cv
207                 = vars ? vars[col] : dict_get_var (mr->dict, var_get_dict_index (mr->varname) + 1 + col);
208               double x = case_data (c, cv)->f;
209               gsl_matrix_set (mr->correlation, crow, col, x);
210             }
211           crow++;
212         }
213       else if (0 == strncasecmp ((char *) value_str (uv, 8), "COV     ", 8))
214         {
215           if (mr->covariance == NULL)
216             mr->covariance   = gsl_matrix_alloc (n_vars, n_vars);
217           for (col = 0; col < n_vars; ++col)
218             {
219               const struct variable *cv
220                 = vars ? vars[col] : dict_get_var (mr->dict, var_get_dict_index (mr->varname) + 1 + col);
221               double x = case_data (c, cv)->f;
222               gsl_matrix_set (mr->covariance, crow, col, x);
223             }
224           crow++;
225         }
226     }
227
228   casereader_destroy (group);
229
230   mm->cov = mr->covariance;
231   mm->corr = mr->correlation;
232
233   return true;
234 }