Index: trunk/psLib/src/fits/Makefile.am =================================================================== --- trunk/psLib/src/fits/Makefile.am (revision 32217) +++ trunk/psLib/src/fits/Makefile.am (revision 32228) @@ -9,4 +9,5 @@ psFitsImage.c \ psFitsTable.c \ + psFitsTableNew.c \ psFitsFloat.c \ psFitsFloatFile.c \ @@ -20,4 +21,5 @@ psFitsImage.h \ psFitsTable.h \ + psFitsTableNew.h \ psFitsFloat.h \ psFitsFloatFile.h \ Index: trunk/psLib/src/fits/psFitsTableNew.c =================================================================== --- trunk/psLib/src/fits/psFitsTableNew.c (revision 32217) +++ trunk/psLib/src/fits/psFitsTableNew.c (revision 32228) @@ -96,4 +96,37 @@ } +void breakhere() +{ +} + +void +freeTable(psFitsTable *table) { + for (int col = 0; col < table->numCols; col++) { + psFitsTableColumn *column = &table->columns[col]; + psFree(column->name); + // all of the members in the data union are pointers so just pick S32 + psFree(column->data.S32); + } + psFree(table->columns); + psFree(table->index); +} + +psFitsTable* psFitsTableAlloc(int numCols, int numRows) { + if (numCols <= 0) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, + "Can't allocate a psFitsTable with %d columns.", numCols); + return NULL; + } + + psFitsTable *table = psAlloc(sizeof(psFitsTable)); + table->index = psMetadataAlloc(); + table->numRows = numRows; + table->numCols = numCols; + table->columns = psAlloc(sizeof(psFitsTableColumn) * numCols); + memset(table->columns, sizeof(psFitsTableColumn) * numCols, 0); + psMemSetDeallocator(table, (psFreeFunc) freeTable); + + return table; +} psFitsTable* psFitsReadTableNew(const psFits* fits) @@ -122,20 +155,6 @@ psTrace("psLib.fits",5,"Table size is %ix%li\n",numCols, numRows); - // the row parameter in the proper range? -#ifdef notdef - if (row < 0 || row >= numRows) { - psError(PS_ERR_BAD_PARAMETER_VALUE, true, - _("Specified row, %d, is not valid for current table of %ld rows."), - row, numRows); - return NULL; - } -#endif - - psFitsTable *table = psAlloc(sizeof(psFitsTable)); - table->index = psMetadataAlloc(); \ - table->numRows = numRows; - table->numCols = numCols; - table->columns = psAlloc(sizeof(psFitsTableColumn) * numCols); -// psMetadata* data = psMetadataAlloc(); + + psFitsTable *table = psFitsTableAlloc(numCols, numRows); int hdutype; // Type of HDU: need to distinguish ASCII and binary tables @@ -151,7 +170,9 @@ for (int col = 1; col <= numCols; col++) { // get the column name + double tscal; + double tzero; if (hdutype == BINARY_TBL) { fits_get_bcolparms(fits->fd, col, name, - NULL, NULL, NULL, NULL, NULL, NULL, NULL, &status); + NULL, NULL, NULL, &tscal, &tzero, NULL, NULL, &status); } else { fits_get_acolparms(fits->fd, col, name, @@ -170,6 +191,8 @@ column->data.TYPE = psAlloc(sizeof(ps##TYPE)*table->numRows); \ if (repeat == 1) { \ + column->elementSize = 1; \ NATIVETYPE *values = (NATIVETYPE *) psAlloc(sizeof(NATIVETYPE) * table->numRows); \ int anynul = 0; \ + breakhere(); \ fits_read_col(fits->fd, FITSTYPE, col, 1, \ 1, numRows, NULL, values, &anynul, &status); \ @@ -180,5 +203,7 @@ col, name, typecode, repeat); \ } else { \ + column->elementSize = repeat; \ column->type = PS_DATA_VECTOR; \ + column->vectorType = PS_DATA_##TYPE; \ column->data.vec = psAlloc(sizeof(psVector*) * table->numRows); \ for (int row = 0; row < table->numRows; row++) { \ @@ -192,7 +217,5 @@ } \ column->data.vec[row] = vec; \ - /* psMetadataAdd(data,PS_LIST_TAIL, name, PS_DATA_VECTOR, "", vec); */ \ psFree(values); \ - /* psFree(vec); */ \ } \ \ @@ -201,11 +224,22 @@ } + // Handle cases where the data is "really" unsigned + if (typecode == TLONG && tzero == 2147483648) { + typecode = TULONG; + } else if (typecode == TSHORT && tzero == 32768) { + typecode = TUSHORT; + } + + // Save the column number for this column in the index metadata psMetadataAddS32(table->index, PS_LIST_TAIL, name, 0, NULL, col-1); + psFitsTableColumn *column = &table->columns[col-1]; - // column->name = psStringCopy(name); + // This second copy of the name makes it easily available when writing the table + column->name = psStringCopy(name); switch (typecode) { - // TBYTE and TSHORT fall though to read into psS32 - case TBYTE: - case TSHORT: + READ_TABLE_ROW_CASE(TBYTE, long, S32, S32); + READ_TABLE_ROW_CASE(TSHORT, short, S16, S16); + READ_TABLE_ROW_CASE(TUSHORT, unsigned short, U16, U16); + READ_TABLE_ROW_CASE(TULONG, unsigned long, U32, U32); READ_TABLE_ROW_CASE(TLONG, long, S32, S32); READ_TABLE_ROW_CASE(TLONGLONG, psS64, S64, S64); @@ -216,4 +250,5 @@ column->data.str = psAlloc(sizeof(psString) * table->numRows); column->type = PS_DATA_STRING; + column->elementSize = 0; for (int row = 0; row < table->numRows; row++) { psString value = psStringAlloc(repeat); @@ -227,4 +262,8 @@ } column->data.str[row] = value; + int len = strlen(value); + if (len > column->elementSize) { + column->elementSize = len; + } } break; @@ -245,4 +284,5 @@ } +// return the index for a given column int psFitsTableGetColumnIndex(psFitsTable *table, psString name) { @@ -258,8 +298,10 @@ } +// return F32 value for given cell psF32 psFitsTableGetCellF32(psFitsTable *table, int row, int col) { return table->columns[col].data.F32[row]; } +// return value for a named column for a given row psF32 psFitsTableGetF32(bool *status, psFitsTable *table, psString name, int row) { bool mdok; @@ -295,6 +337,15 @@ psFitsTableColumn *column = &table->columns[col]; switch (column->type) { + case PS_DATA_S16: + column->data.S16[newRow] = column->data.S16[row]; + break; + case PS_DATA_U16: + column->data.U16[newRow] = column->data.U16[row]; + break; case PS_DATA_S32: column->data.S32[newRow] = column->data.S32[row]; + break; + case PS_DATA_U32: + column->data.U32[newRow] = column->data.U32[row]; break; case PS_DATA_S64: @@ -331,6 +382,6 @@ } } + // set the new number of rows - // XXX: should we do this? table->numRows = newRow; @@ -338,6 +389,325 @@ } +// Get the TFORM character, given a PS type +static inline char getTForm(psDataType type) +{ + switch (type) { + case PS_TYPE_U8: + case PS_TYPE_S8: + return 'B'; + case PS_TYPE_S16: + return 'I'; + case PS_TYPE_S32: + return 'J'; + case PS_TYPE_U64: + case PS_TYPE_S64: + return 'K'; + case PS_TYPE_U16: + return 'U'; + case PS_TYPE_U32: + return 'V'; + case PS_TYPE_F32: + return 'E'; + case PS_TYPE_F64: + return 'D'; + case PS_TYPE_BOOL: + return 'L'; + case PS_DATA_STRING: + return 'A'; + default: + psError(PS_ERR_UNKNOWN, true, "Unknown type: %x\n", type); + return '?'; + } +} + +static bool fitsInsertTableNew(psFits* fits, // FITS file + const psMetadata* header, // FITS header to write + psFitsTable *table, // internal representation of the table + const char *extname, // Extension name to give table + bool after, // Write table after current extension? + bool writeData // Write data? + ) +{ + PS_ASSERT_FITS_NON_NULL(fits, false); + PS_ASSERT_FITS_WRITABLE(fits, false); + PS_ASSERT_PTR_NON_NULL(table, false); + + int status = 0; + + long numRows = table->numRows; + if (writeData && numRows < 1) { + // no table data, what can I do? + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + _("Can't create a table without any rows.")); + return false; + } + #ifdef notdef - + // Find the unique items in the array of metadata 'rows', and their sizes + psMetadata *colSpecs = psMetadataAlloc(); // Column specifications + size_t rowSize = 0; // Size (in bytes) of each row + for (long i = 0; i < numRows; i++) { + psMetadata* row = table->data[i]; + + if (!row) { + continue; + } + psMetadataIterator *rowIter = psMetadataIteratorAlloc(row, PS_LIST_HEAD, NULL); // Iterator + psMetadataItem *colItem; // Column item, from iteration + while ((colItem = psMetadataGetAndIncrement(rowIter))) { + if (!(PS_DATA_IS_PRIMITIVE(colItem->type) || colItem->type == PS_DATA_STRING || + colItem->type == PS_DATA_VECTOR)) { + // unsupported type -- treating as an error + psError(PS_ERR_BAD_PARAMETER_TYPE, true, + "Unsupported data type (%d) for Metadata Item '%s' in row %ld.", + colItem->type, colItem->name, i); + psFree(rowIter); + psFree(colSpecs); + return false; + } + + size_t size = columnSize(colItem); // Size for this column + + // Check to see if we know about this one; or update the size if required + psMetadataItem *colSpecItem = psMetadataLookup(colSpecs, colItem->name); + if (!colSpecItem) { + // A new one! + colSpec *spec = psAlloc(sizeof(colSpec)); // Specification for this column + // BOOL type is not a valid vector type, so we translate it to U8 + spec->type = colItem->type == PS_TYPE_BOOL ? PS_TYPE_U8 : colItem->type; + spec->size = size; + if (colItem->type == PS_DATA_VECTOR) { + psVector *vector = colItem->data.V; // The vector + spec->vectorType = vector->type.type; + } + psMetadataAddPtr(colSpecs, PS_LIST_TAIL, colItem->name, PS_DATA_UNKNOWN, "", spec); + psFree(spec); // Drop reference + rowSize += PSELEMTYPE_SIZEOF(spec->type); + } else { + colSpec *spec = colSpecItem->data.V; // The specification + if (size > spec->size) { + spec->size = size; + } + if (colItem->type != spec->type && + colItem->type != PS_TYPE_BOOL && spec->type != PS_TYPE_U8) { + psWarning("Differing type found for column %s: %x vs %x --- using the first found.\n", + colSpecItem->name, colItem->type, spec->type); + } + if (colItem->type == PS_DATA_VECTOR) { + psVector *vector = colItem->data.V; // The vector + if (vector->type.type != spec->vectorType) { + psWarning("Differing vector type found for column %s: %x vs %x " + "--- using the first found.\n", colSpecItem->name, vector->type.type, + spec->vectorType); + } + } + } + } + psFree(rowIter); + } + + long numColumns = colSpecs->list->n;// Number of columns + if (numColumns == 0) { + // No table columns found + psError(PS_ERR_BAD_PARAMETER_SIZE, true, + "Did not find any column data to write to a table."); + psFree(colSpecs); + return false; + } +#endif // notdef end of OLD column spec determination + + // Create array of column names and types. + psArray *columnNames = psArrayAlloc(table->numCols); // Array of column names, for cfitsio + psArray *columnTypes = psArrayAlloc(table->numCols); // Array of column types, for cfitsio + for (long i = 0; i < table->numCols; i++) { + psFitsTableColumn *column = &table->columns[i]; + columnNames->data[i] = psMemIncrRefCounter(column->name); + psString colType = NULL; // The FITS column type + size_t size = column->elementSize; + if (column->type == PS_DATA_VECTOR) { + psStringAppend(&colType, "%zd%c", size, getTForm(column->vectorType)); + } else { + psStringAppend(&colType, "%zd%c", size, getTForm(column->type)); + } + columnTypes->data[i] = colType; + } + + // Create the table HDU + int numHDUs = psFitsGetSize(fits); // Number of HDUs in file + if (numHDUs == 0) { + // We're creating the first extension + fits_create_tbl(fits->fd, + BINARY_TBL, + writeData ? numRows : 0, // number of rows in table + table->numCols, // number of columns in table + (char**)columnNames->data, // names of the columns + (char**)columnTypes->data, // format of the columns + NULL, // physical unit of columns + NULL, // skip extension name: we set the by hand below + &status); + } else { + if (!after) { + if (psFitsGetExtNum(fits) == 0) { + // We're creating a replacement primary HDU. + // Set status to signal fits_insert_img to insert a new primary HDU + status = PREPEND_PRIMARY; + } else { + // Move back one to perform an insert after the previous HDU + psFitsMoveExtNum(fits, -1, true); + } + } + // Insert the table + fits_insert_btbl(fits->fd, + writeData ? numRows : 0, // number of rows in table + table->numCols, // number of columns in table + (char**)columnNames->data, // names of the columns + (char**)columnTypes->data, // format of the columns + NULL, // physical unit of columns + NULL, // skip extension name: we set this by hand below + 0, &status); + } + psFree(columnNames); + psFree(columnTypes); + + if (status != 0) { + psFitsError(status, true, "Unable to create FITS table with %d columns and %d rows", + table->numCols, table->numRows); + return false; + } + + // Write header + if (header && !psFitsWriteHeader(fits, header)) { + psError(PS_ERR_IO, false, "Unable to write FITS header.\n"); + return false; + } + + // write the header, if any. + if (extname && strlen(extname) > 0) { + if (!psFitsSetExtName(fits, extname)) { + psError(PS_ERR_IO, false, "Unable to write FITS header extension name.\n"); + return false; + } + } + + if (writeData) { + // psMetadataIteratorSet(colSpecsIter, PS_LIST_HEAD); + for (long colNum = 1; colNum <= table->numCols; colNum++) { + // Note: colNum is unit-indexed, because it's for cfitsio + psFitsTableColumn *column = &table->columns[colNum - 1]; + // colSpec *spec = colSpecItem->data.V; // The specification + if (PS_DATA_IS_PRIMITIVE(column->type)) { +#ifdef notdef + size_t dataSize = PSELEMTYPE_SIZEOF(column->type); // Size (in bytes) of this type + psVector *columnData = psVectorAlloc(table->n, column->type); // The raw row data, to be written + psVectorInit(columnData, 0); + for (long i = 0; i < table->numRows; i++) { + psMetadata *row = table->data[i]; // The row of interest + psMetadataItem *dataItem = psMetadataLookup(row, colSpecItem->name); // Value of interest + if (dataItem) { + memcpy(&columnData->data.U8[i * dataSize], &dataItem->data, dataSize); + } else { + // this element is missing from this row; insert an appropriate-sized place holder + // XXX this should insert a NAN for float / double and an appropriate blank for int types + // XXX for the moment I am putting in 0.0 + memset(&columnData->data.U8[i * dataSize], 0, dataSize); + } + } +#endif + + int fitsDataType; // Data type for cfitsio + p_psFitsTypeToCfitsio(column->type, NULL, NULL, &fitsDataType); + fits_write_col(fits->fd, + fitsDataType, + colNum, // column number + 1, // first row + 1, // first element + table->numRows, // number of rows + column->data.U8, // the data + &status); + // psFree(columnData); + } else { + switch (column->type) { + case PS_DATA_STRING: { +#ifdef notdef + psArray *strings = psArrayAlloc(table->n); // Array of strings + for (long i = 0; i < table->n; i++) { + psMetadata *row = table->data[i]; // The row of interest + strings->data[i] = psMemIncrRefCounter(psMetadataLookupStr(NULL, row, + colSpecItem->name)); + } + fits_write_col_str(fits->fd, colNum, 1, 1, table->n, (char**)strings->data, &status); +#endif + fits_write_col_str(fits->fd, colNum, 1, 1, table->numRows, (char**)column->data.str, &status); + // psFree(strings); + break; + } + case PS_DATA_VECTOR: { + size_t dataSize = PSELEMTYPE_SIZEOF(column->vectorType); // Size of data, in bytes + psVector *columnData = psVectorAlloc(column->elementSize * table->numRows * dataSize, PS_TYPE_U8); + psVectorInit(columnData, 0); + for (long i = 0; i < table->numRows; i++) { +#ifdef notdef + psMetadata *row = table->data[i]; // The row of interest + psMetadataItem* dataItem = psMetadataLookup(row, colSpecItem->name); + if (dataItem->type != PS_DATA_VECTOR) { + // Just in case --- get a zero instead of some weird result + continue; + } + psVector *vector = dataItem->data.V; +#endif + psVector *vector = column->data.vec[i]; + memcpy(&columnData->data.U8[i * dataSize * column->elementSize], vector->data.U8, + vector->n * dataSize); + } + + int fitsDataType; // Data type for cfitsio + p_psFitsTypeToCfitsio(column->vectorType, NULL, NULL, &fitsDataType); + fits_write_col(fits->fd, fitsDataType, colNum, 1, 1, table->numRows * column->elementSize, + columnData->data.U8, &status); + psFree(columnData); + break; + } + default: + psAbort("Should never get here.\n"); + } + } + + // Check error status from writing column + if (status != 0) { + psFitsError(status, true, "Unable to write column %ld of FITS table", colNum); + return false; + } + } + } + + // This forces a re-scan of the header to ensure everything's kosher. We found this occassionally + // necessary for compressed images, which are tables, so perhaps it helps here too. I guess it can't + // hurt. + ffrdef(fits->fd, &status); + if (psFitsError(status, true, "Could not re-scan HDU.")) { + return false; + } + + return true; +} + +// END OF new insert table + +bool psFitsWriteTableNew(psFits* fits, + const psMetadata* header, + psFitsTable* table, + const char *extname) +{ + PS_ASSERT_FITS_NON_NULL(fits, false); + PS_ASSERT_FITS_WRITABLE(fits, false); + if (!psFitsMoveLast(fits)) { + psError(PS_ERR_UNKNOWN, false, "Unable to move to last extension to write table"); + return false; + } + bool writeData = table->numRows > 0; + return fitsInsertTableNew(fits, header, table, extname, true, writeData); +} + +#ifdef notdef psArray* psFitsReadTableColumn(const psFits* fits, Index: trunk/psLib/src/fits/psFitsTableNew.h =================================================================== --- trunk/psLib/src/fits/psFitsTableNew.h (revision 32217) +++ trunk/psLib/src/fits/psFitsTableNew.h (revision 32228) @@ -14,6 +14,9 @@ typedef struct { - /* psString name; */ + psString name; psDataType type; + psDataType vectorType; // type inside if Vector + int elementSize; // 1 for primitives, strlen for strings + // vector length for vectors // only types supported by fits tables will be implemented union { @@ -49,5 +52,5 @@ psFitsTable *psFitsReadTableNew(const psFits *fits); -bool psFitsWriteTableNew(psFits *fits, psFitsTable *table, const char *extname); +bool psFitsWriteTableNew(psFits *fits, const psMetadata *header, psFitsTable *table, const char *extname); bool psFitsTableCensor(psFitsTable *table, bool *rowMask);