Allocate omxData with Calloc
[openmx:openmx.git] / src / omxData.c
1 /*
2  *  Copyright 2007-2013 The OpenMx Project
3  *
4  *  Licensed under the Apache License, Version 2.0 (the "License");
5  *  you may not use this file except in compliance with the License.
6  *  You may obtain a copy of the License at
7  *
8  *       http://www.apache.org/licenses/LICENSE-2.0
9  *
10  *  Unless required by applicable law or agreed to in writing, software
11  *  distributed under the License is distributed on an "AS IS" BASIS,
12  *  WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13  *  See the License for the specific language governing permissions and
14  *  limitations under the License.
15  */
16
17 /***********************************************************
18  *
19  *  omxData.cc
20  *
21  *  Created: Timothy R. Brick   Date: 2009-07-15
22  *
23  *      Contains code for the omxData class
24  *   omxData objects hold data in whatever form it takes
25  *
26  **********************************************************/
27 #include "omxData.h"
28 #include "npsolWrap.h"
29
30 omxData* omxInitData(omxData* od, omxState* os) {
31
32         if(od == NULL) {
33                 od = Calloc(1, omxData);
34         }
35
36         od->currentState = os;
37         
38         if(OMX_DEBUG) {Rprintf("Data's state object is at 0x%x.\n", od->currentState);}
39
40         return od;
41
42 }
43
44 omxData* omxDataLookupFromState(SEXP dataObject, omxState* state) {
45         int dataIdx = INTEGER(dataObject)[0];
46
47         return state->dataList[dataIdx];
48 }
49
50 omxData* omxNewDataFromMxData(omxData* data, SEXP dataObject, omxState* state) {
51         if(OMX_DEBUG) {Rprintf("Initializing data Element.\n");}
52         if(dataObject == NULL) {
53                 error("Null Data Object detected.  This is an internal error, and should be reported on the forums.\n");
54         }
55
56         omxData* od = data;
57         if(od == NULL) {
58                 od = omxInitData(data, state);
59         }
60
61         SEXP dataLoc, dataVal;
62         int numCols;
63         int numInts=0;
64         int numReals=0;
65
66         // PARSE MxData Structure
67         if(OMX_DEBUG) {Rprintf("Processing Data Type.\n");}
68         PROTECT(dataLoc = GET_SLOT(dataObject, install("type")));
69         if(dataLoc == NULL) { error("Data has no type.  Sorry.\nThis is an internal error, and should be reported on the forums.\n");}
70         PROTECT(dataVal = STRING_ELT(dataLoc,0));
71         od->_type = CHAR(dataVal);
72         if(OMX_DEBUG) {Rprintf("Element is type %s.\n", od->_type);}
73
74         PROTECT(dataLoc = GET_SLOT(dataObject, install("observed")));
75         if(OMX_DEBUG) {Rprintf("Processing Data Elements.\n");}
76         if(isFrame(dataLoc)) {
77                 if(OMX_DEBUG) {Rprintf("Data is a frame.\n");}
78                 // Process Data Frame into Columns
79                 od->cols = length(dataLoc);
80                 if(OMX_DEBUG) {Rprintf("Data has %d columns.\n", od->cols);}
81                 numCols = od->cols;
82                 od->columns = (SEXP*) R_alloc(numCols, sizeof(SEXP));
83                 od->realData = (double**) R_alloc(numCols, sizeof(double*));
84                 od->intData = (int**) R_alloc(numCols, sizeof(int*));
85                 od->location = (int*) R_alloc(numCols, sizeof(int));
86                 for(int j = 0; j < numCols; j++) {
87                         PROTECT(od->columns[j] = VECTOR_ELT(dataLoc, j));
88                         if(isFactor(od->columns[j])) {
89                                 if(OMX_DEBUG) {Rprintf("Column %d is a factor.\n", j);}
90                                 od->intData[numInts] = INTEGER(od->columns[j]);
91                                 od->location[j] = ~(numInts++);
92                                 od->numFactor++;
93                         } else if (isInteger(od->columns[j])) {
94                                 error("Internal error: Column %d is in integer format.", j);
95                         } else {
96                                 if(OMX_DEBUG) {Rprintf("Column %d is a numeric.\n", j);}
97                                 od->realData[numReals] = REAL(od->columns[j]);
98                                 od->location[j] = (numReals++);
99                                 od->numNumeric++;
100                         }
101                 }
102                 od->rows = length(VECTOR_ELT(dataLoc, 0));
103                 if(OMX_DEBUG) {Rprintf("And %d rows.\n", od->rows);}
104         } else {
105                 if(OMX_DEBUG) {Rprintf("Data contains a matrix.\n");}
106                 od->dataMat = omxNewMatrixFromRPrimitive(dataLoc, od->currentState, 0, 0);
107                 
108                 if (od->dataMat->colMajor && strncmp(od->_type, "raw", 3) == 0) { 
109                         omxToggleRowColumnMajor(od->dataMat);
110                 }
111                 od->cols = od->dataMat->cols;
112                 od->rows = od->dataMat->rows;
113                 od->numNumeric = od->cols;
114         }
115
116         if(OMX_DEBUG) {Rprintf("Processing Means Matrix.\n");}
117         PROTECT(dataLoc = GET_SLOT(dataObject, install("means")));
118         od->meansMat = omxNewMatrixFromRPrimitive(dataLoc, od->currentState, 0, 0);
119         if(od->meansMat->rows == 1 && od->meansMat->cols == 1 && 
120            (!R_finite(omxMatrixElement(od->meansMat, 0, 0)) ||
121             !isfinite(omxMatrixElement(od->meansMat, 0, 0)))) {
122                 omxFreeMatrixData(od->meansMat); // Clear just-allocated memory.
123                 od->meansMat = NULL;  // 1-by-1 matrix of NAs is a null means matrix.
124                 // FIXME: The above check may cause problems for dynamic data if the means
125                 //          originally is a 1x1 that has not yet been calculated.  This should be
126                 //          adjusted.
127         }
128         
129         if(OMX_DEBUG) {
130                 if(od->meansMat == NULL) {Rprintf("No means found.\n");}
131                 else {omxPrint(od->meansMat, "Means Matrix is:");}
132         }
133
134         if(strncmp(od->_type, "raw", 3) != 0) {
135                 if(OMX_DEBUG) {Rprintf("Processing Observation Count.\n");}
136                 PROTECT(dataLoc = GET_SLOT(dataObject, install("numObs")));
137                 od->numObs = REAL(dataLoc)[0];
138         } else {
139                 od->numObs = od->rows;
140                 if(OMX_DEBUG) {Rprintf("Processing presort metadata.\n");}
141                 /* For raw data, process sorting metadata. */
142                 // Process unsorted indices:  // TODO: Generate reverse lookup table
143                 PROTECT(dataLoc = GET_SLOT(dataObject, install("indexVector")));
144                 od->indexVector = INTEGER(dataLoc);
145                 if(od->indexVector[0] == R_NaInt) od->indexVector = NULL;
146                 // Process pre-computed identicality checks
147                 if(OMX_DEBUG) {Rprintf("Processing definition variable identicality.\n");}
148                 PROTECT(dataLoc = GET_SLOT(dataObject, install("identicalDefVars")));
149                 od->identicalDefs = INTEGER(dataLoc);
150                 if(od->identicalDefs[0] == R_NaInt) od->identicalDefs = NULL;
151                 if(OMX_DEBUG) {Rprintf("Processing missingness identicality.\n");}
152                 PROTECT(dataLoc = GET_SLOT(dataObject, install("identicalMissingness")));
153                 od->identicalMissingness = INTEGER(dataLoc);
154                 if(od->identicalMissingness[0] == R_NaInt) od->identicalMissingness = NULL;
155                 if(OMX_DEBUG) {Rprintf("Processing row identicality.\n");}
156                 PROTECT(dataLoc = GET_SLOT(dataObject, install("identicalRows")));
157                 od->identicalRows = INTEGER(dataLoc);
158                 if(od->identicalRows[0] == R_NaInt) od->identicalRows = NULL;
159         }
160
161         return od;
162 }
163
164 void resetDefinitionVariables(double *oldDefs, int numDefs) {
165         int nextDef;
166
167         for(nextDef = 0; nextDef < numDefs; nextDef++) {
168                 oldDefs[nextDef] = NA_REAL;                                     // Def Vars default to NA
169         }
170
171 }
172
173 void omxFreeData(omxData* od) {
174         omxFreeAllMatrixData(od->dataMat);
175         omxFreeAllMatrixData(od->meansMat);
176         Free(od);
177 }
178
179 double omxDoubleDataElement(omxData *od, int row, int col) {
180         if(od->dataMat != NULL) {
181                 return omxMatrixElement(od->dataMat, row, col);
182         }
183         int location = od->location[col];
184         if(location < 0) {
185                 return (double)(od->intData[~location][row]);
186         } else {
187                 return od->realData[location][row];
188         }
189 }
190
191 int omxIntDataElement(omxData *od, int row, int col) {
192         if(od->dataMat != NULL) {
193                 error("Use a data frame for factor data");
194         }
195
196         int location = od->location[col];
197         if(location < 0) {
198                 return (od->intData[~location][row]);
199         } else {
200                 return (int)(od->realData[location][row]);
201         }
202 }
203
204 omxMatrix* omxDataMatrix(omxData *od, omxMatrix* om) {
205         double dataElement;
206
207         if(od->dataMat != NULL) {               // Data was entered as a matrix.
208                 if(om != NULL) {                        // It stays as such
209                         omxCopyMatrix(om, od->dataMat);
210                         return om;
211                 }
212                 return od->dataMat;
213         }
214         // Otherwise, we must construct the matrix.
215         int numRows = od->rows, numCols = od->cols;
216
217         if(om == NULL) {
218                 om = omxInitMatrix(om, numRows, numCols, TRUE, od->currentState);
219         }
220
221         if(om->rows != numRows || om->cols != numCols) {
222                 omxResizeMatrix(om, numRows, numCols, FALSE);
223         }
224
225         for(int j = 0; j < numCols; j++) {
226                 for(int k = 0; k < numRows; k++) {
227                         int location = od->location[j];
228                         if(location < 0) {
229                                 dataElement = (double) od->intData[~location][k];
230                         } else {
231                                 dataElement = od->realData[location][k];
232                         }
233                         omxSetMatrixElement(om, k, j, dataElement);
234                 }
235         }
236         return om;
237 }
238
239 unsigned short int omxDataColumnIsFactor(omxData *od, int col) {
240         if(od->dataMat != NULL) return FALSE;
241         if(col <= od->cols) return (od->location[col] < 0);
242
243         error("Attempted to access column %d of a %d-column data object", col, od->cols);
244         return 0; // not reached
245 }
246
247 omxMatrix* omxDataMeans(omxData *od, omxMatrix* colList, omxMatrix* om) {
248
249         if(od->meansMat == NULL) return NULL;
250
251         if(colList == NULL) {
252                 if(om == NULL) return od->meansMat;
253                 omxCopyMatrix(om, od->meansMat);
254                 return om;
255         }
256
257         int cols = colList->cols;
258
259         if(colList == NULL || cols == 0 || cols > od->cols) {
260                 cols = od->cols;
261                 if(om == NULL) return od->meansMat;
262                 omxCopyMatrix(om, od->meansMat);
263                 return om;
264         }
265
266         if(om == NULL) {
267                 om = omxInitMatrix(om, 1, cols, TRUE, od->currentState);
268         }
269
270         for(int i = 0; i < cols; i++) {
271                 omxSetMatrixElement(om, 1, i, omxVectorElement(od->meansMat, omxVectorElement(colList, i)));
272         }
273
274         return om;
275 }
276
277
278 void omxSetContiguousDataColumns(omxContiguousData* contiguous, omxData* data, omxMatrix* colList) {
279
280         contiguous->isContiguous = FALSE;   // Assume not contiguous
281
282         if (data->dataMat == NULL) return; // Data has no matrix elements, so skip.
283
284         omxMatrix* dataMat = data->dataMat;
285         if (dataMat->colMajor) return;      // If data matrix is column-major, there's no continuity
286         
287         int colListLength = colList->cols;              // # of columns in the cov matrix
288         double start = omxVectorElement(colList, 0);    // Data col of first column of the covariance
289         contiguous->start = (int) start;                // That's our starting point.
290         contiguous->length = colListLength;             // And the length is ncol(cov)
291         
292         for(int i = 1; i < colListLength; i++) {        // Make sure that the col list is 
293                 double next = omxVectorElement(colList, i); // contiguously increasing in column number
294                 if (next != (start + i)) return;            // If it isn't, it's not contiguous data
295         }
296         
297         contiguous->isContiguous = TRUE;    // Passed.  This is contiguous.
298 }
299
300 omxMatrix* omxContiguousDataRow(omxData *od, int row, int start, int length, omxMatrix* om) {
301         // TODO: Might be better to combine this with omxDataRow to make a single accessor omxDataRow with a second signature that accepts an omxContiguousData argument.
302         if(row > od->rows) return NULL; // Sanity check
303
304         if(om == NULL) {
305                 om = omxInitMatrix(om, 1, od->cols, TRUE, od->currentState);
306         }
307         
308         int numcols = od->cols;
309         omxMatrix* dataMat = od->dataMat;
310         double *dest = om->data;
311         double *source = dataMat->data + row * numcols + start;
312         memcpy(dest, source, sizeof(double) * length);
313         return(om);
314 }
315
316 omxMatrix* omxDataRow(omxData *od, int row, omxMatrix* colList, omxMatrix* om) {
317
318         if(colList == NULL || row > od->rows) return NULL;      // Sanity check
319
320         if(om == NULL) {
321                 om = omxInitMatrix(om, 1, od->cols, TRUE, od->currentState);
322         }
323
324         int numcols = om->cols;
325         if(od->dataMat != NULL) { // Matrix Object
326                 omxMatrix* dataMat = od->dataMat;
327                 for(int j = 0; j < numcols; j++) {
328                         omxSetMatrixElement(om, 0, j, omxMatrixElement(dataMat, row, 
329                                                                        omxVectorElement(colList, j)));
330                 }
331         } else {                // Data Frame object
332                 double dataElement;
333                 int* locations = od->location;
334                 int** intDataColumns = od->intData;
335                 double **realDataColumns = od->realData;
336                 for(int j = 0; j < numcols; j++) {
337                         int location = locations[(int)omxVectorElement(colList, j)];
338                         if(location < 0) {
339                                 dataElement = (double) intDataColumns[~location][row];
340                         } else {
341                                 dataElement = realDataColumns[location][row];
342                         }
343                         omxSetMatrixElement(om, 0, j, dataElement);
344                 }
345         }
346         return om;
347 }
348
349 int omxDataIndex(omxData *od, int row) {
350         if(od->indexVector != NULL)
351                 return od->indexVector[row];
352         else return row;
353 }
354
355 int omxDataNumIdenticalRows(omxData *od, int row) {
356         if(od->identicalRows != NULL)
357                 return od->identicalRows[row];
358         else return 1;
359 }
360 int omxDataNumIdenticalMissingness(omxData *od, int row) {
361         if(od->identicalMissingness != NULL)
362                 return od->identicalMissingness[row];
363         else return 1;
364 }
365
366 int omxDataNumIdenticalDefs(omxData *od, int row){
367         if(od->identicalDefs != NULL)
368                 return od->identicalDefs[row];
369         else return 1;
370 }
371
372 int omxDataNumIdenticalContinuousRows(omxData *od, int row) {
373         if(od->numNumeric <= 0) {
374                 return od->rows;
375         }
376         return omxDataNumIdenticalRows(od, row);
377 }
378
379 int omxDataNumIdenticalContinuousMissingness(omxData *od, int row) {
380         if(od->numNumeric <= 0) {
381                 return od->rows;
382         }
383         return omxDataNumIdenticalMissingness(od, row);
384 }
385
386 int omxDataNumIdenticalOrdinalRows(omxData *od, int row) {
387         if(od->numFactor <= 0) {
388                 return od->rows;
389         }
390         return omxDataNumIdenticalRows(od, row);
391 }
392
393 int omxDataNumIdenticalOrdinalMissingness(omxData *od, int row) {
394         if(od->numFactor <= 0) {
395                 return od->rows;
396         }
397         return omxDataNumIdenticalMissingness(od, row);
398 }
399
400
401 double omxDataNumObs(omxData *od) {
402         return od->numObs;
403 }
404
405 int omxDataNumFactor(omxData *od) {
406         return od->numFactor;
407 }
408
409 int omxDataNumNumeric(omxData *od) {
410         return od->numNumeric;
411 }
412
413
414 const char *omxDataType(omxData *od) {
415         return od->_type;
416 }
417
418 int elementEqualsDataframe(SEXP column, int offset1, int offset2) {
419         switch (TYPEOF(column)) {
420         case REALSXP:
421                 if(ISNA(REAL(column)[offset1])) return ISNA(REAL(column)[offset2]);
422                 if(ISNA(REAL(column)[offset2])) return ISNA(REAL(column)[offset1]);
423                 return(REAL(column)[offset1] == REAL(column)[offset2]);
424         case LGLSXP:
425         case INTSXP:
426                 return(INTEGER(column)[offset1] == INTEGER(column)[offset2]);           
427         }
428         return(0);
429 }
430
431 int testRowDataframe(SEXP data, int numrow, int numcol, int i, int *row, int base) {
432         SEXP column;
433         int j, equal = TRUE;
434
435         if (i == numrow) {
436                 equal = FALSE;
437         } else {
438                 for(j = 0; j < numcol && equal; j++) {
439                         column = VECTOR_ELT(data, j);
440                         equal = elementEqualsDataframe(column, base, i);
441                 }
442         }
443
444         if (!equal) {
445                 int gap = i - base;
446                 for(j = 0; j < gap; j++) {
447                         row[base + j] = gap - j;
448                 }
449                 base = i;
450         }
451         return(base);
452 }
453
454 int elementEqualsMatrix(SEXP data, int row1, int row2, int numrow, int col) {
455         int coloffset = col * numrow;
456         switch (TYPEOF(data)) {
457         case REALSXP:
458                 if(ISNA(REAL(data)[row1 + coloffset])) return ISNA(REAL(data)[row2 + coloffset]);
459                 if(ISNA(REAL(data)[row2 + coloffset])) return ISNA(REAL(data)[row1 + coloffset]);
460                 return(REAL(data)[row1 + coloffset] == REAL(data)[row2 + coloffset]);
461         case LGLSXP:
462         case INTSXP:
463                 return(INTEGER(data)[row1 + coloffset] == INTEGER(data)[row2 + coloffset]);
464         }
465         return(0);
466 }
467
468 int testRowMatrix(SEXP data, int numrow, int numcol, int i, int *row, int base) {
469         int j, equal = TRUE;
470
471         if (i == numrow) {
472                 equal = FALSE;
473         } else {
474                 for(j = 0; j < numcol && equal; j++) {
475                         equal = elementEqualsMatrix(data, i, base, numrow, j);
476                 }
477         }
478
479         if (!equal) {
480                 int gap = i - base;
481                 for(j = 0; j < gap; j++) {
482                         row[base + j] = gap - j;
483                 }
484                 base = i;
485         }
486         return(base);
487 }
488
489 SEXP findIdenticalMatrix(SEXP data, SEXP missing, SEXP defvars,
490                          SEXP skipMissingExp, SEXP skipDefvarsExp) {
491
492         SEXP retval, identicalRows, identicalMissing, identicalDefvars;
493         int i, numrow, numcol, defvarcol;
494         int *irows, *imissing, *idefvars;
495         int baserows, basemissing, basedefvars;
496         int skipMissing, skipDefvars;
497
498         skipMissing = LOGICAL(skipMissingExp)[0];
499         skipDefvars = LOGICAL(skipDefvarsExp)[0];
500         numrow = nrows(data);
501         numcol = ncols(data);
502         defvarcol = ncols(defvars);
503         PROTECT(retval = allocVector(VECSXP, 3));
504         PROTECT(identicalRows = allocVector(INTSXP, numrow));
505         PROTECT(identicalMissing = allocVector(INTSXP, numrow));
506         PROTECT(identicalDefvars = allocVector(INTSXP, numrow));
507         irows = INTEGER(identicalRows);
508         imissing = INTEGER(identicalMissing);
509         idefvars = INTEGER(identicalDefvars);
510         if (skipMissing) {
511                 for(i = 0; i < numrow; i++) {
512                         imissing[i] = numrow - i;
513                 }
514         }
515         if (skipDefvars) {
516                 for(i = 0; i < numrow; i++) {
517                         idefvars[i] = numrow - i;
518                 }
519         }
520         baserows = 0;
521         basemissing = 0;
522         basedefvars = 0;
523         for(i = 1; i <= numrow; i++) {
524                 baserows = testRowMatrix(data, numrow, numcol, i, irows, baserows); 
525                 if (!skipMissing) {
526                         basemissing = testRowMatrix(missing, numrow, numcol, i, imissing, basemissing); 
527                 }
528                 if (!skipDefvars) {
529                         basedefvars = testRowMatrix(defvars, numrow, defvarcol, i, idefvars, basedefvars);
530                 }
531         }
532         SET_VECTOR_ELT(retval, 0, identicalRows);
533         SET_VECTOR_ELT(retval, 1, identicalMissing);
534         SET_VECTOR_ELT(retval, 2, identicalDefvars);
535         UNPROTECT(4); // retval, identicalRows, identicalMissing, identicalDefvars
536         return retval;
537 }
538
539 SEXP findIdenticalDataFrame(SEXP data, SEXP missing, SEXP defvars,
540                             SEXP skipMissingExp, SEXP skipDefvarsExp) {
541
542         SEXP retval, identicalRows, identicalMissing, identicalDefvars;
543         int i, numrow, numcol, defvarcol;
544         int *irows, *imissing, *idefvars;
545         int baserows, basemissing, basedefvars;
546         int skipMissing, skipDefvars;
547
548         skipMissing = LOGICAL(skipMissingExp)[0];
549         skipDefvars = LOGICAL(skipDefvarsExp)[0];
550         numrow = length(VECTOR_ELT(data, 0));
551         numcol = length(data);
552         defvarcol = length(defvars);
553         PROTECT(retval = allocVector(VECSXP, 3));
554         PROTECT(identicalRows = allocVector(INTSXP, numrow));
555         PROTECT(identicalMissing = allocVector(INTSXP, numrow));
556         PROTECT(identicalDefvars = allocVector(INTSXP, numrow));
557         irows = INTEGER(identicalRows);
558         imissing = INTEGER(identicalMissing);
559         idefvars = INTEGER(identicalDefvars);
560         if (skipMissing) {
561                 for(i = 0; i < numrow; i++) {
562                         imissing[i] = numrow - i;
563                 }
564         }
565         if (skipDefvars) {
566                 for(i = 0; i < numrow; i++) {
567                         idefvars[i] = numrow - i;
568                 }
569         }
570         baserows = 0;
571         basemissing = 0;
572         basedefvars = 0;
573         for(i = 1; i <= numrow; i++) {
574                 baserows = testRowDataframe(data, numrow, numcol, i, irows, baserows); 
575                 if (!skipMissing) {
576                         basemissing = testRowMatrix(missing, numrow, numcol, i, imissing, basemissing);
577                 }
578                 if (!skipDefvars) {
579                         basedefvars = testRowDataframe(defvars, numrow, defvarcol, i, idefvars, basedefvars);
580                 }
581         }
582         SET_VECTOR_ELT(retval, 0, identicalRows);
583         SET_VECTOR_ELT(retval, 1, identicalMissing);
584         SET_VECTOR_ELT(retval, 2, identicalDefvars);
585         UNPROTECT(4); // retval, identicalRows, identicalMissing, identicalDefvars
586         return retval;
587 }
588
589 SEXP findIdenticalRowsData(SEXP data, SEXP missing, SEXP defvars,
590                            SEXP skipMissingness, SEXP skipDefvars) {
591         if (isMatrix(data)) {
592                 return(findIdenticalMatrix(data, missing, defvars, skipMissingness, skipDefvars));
593         } else {
594                 return(findIdenticalDataFrame(data, missing, defvars, skipMissingness, skipDefvars));
595         }
596 }
597
598 void omxPrintData(omxData *od, const char *header) {
599         if (!header) error("omxPrintData: header is NULL");
600
601         if (!od) {
602                 Rprintf("%s: NULL\n", header);
603                 return;
604         }
605
606         Rprintf("%s(%s): %f observations %d x %d\n", header, od->_type, od->numObs,
607                 od->rows, od->cols);
608         Rprintf("numNumeric %d numFactor %d\n", od->numNumeric, od->numFactor);
609
610         if (od->columns) {
611                 for(int j = 0; j < od->cols; j++) {
612                         PrintValue(od->columns[j]);
613                 }
614         }
615
616         if (od->location) {
617                 Rprintf("Location: ");
618                 for(int j = 0; j < od->cols; j++) {
619                         Rprintf("%d ", od->location[j]);
620                 }
621                 Rprintf("\n");
622         }
623
624         if (od->identicalRows) {
625                 Rprintf("\trow\tmissing\tdefvars\n");
626                 for(int j = 0; j < od->rows; j++) {
627                         Rprintf("[%d]\t%d\t%d\t%d\n", j, od->identicalRows[j],
628                                 od->identicalMissingness[j], od->identicalDefs[j]);
629                 }
630         }
631
632         if (od->dataMat) omxPrintMatrix(od->dataMat, "dataMat");
633         if (od->meansMat) omxPrintMatrix(od->meansMat, "meansMat");
634 }
635