IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 7846


Ignore:
Timestamp:
Jul 7, 2006, 5:06:40 PM (20 years ago)
Author:
Paul Price
Message:
  • Moving common functions (convertFitsToPsType and convertPsTypeToFits) into

psFits as private functions (p_psFitsTypeFromCfitsio and p_psFitsTypeToCfitsio).

  • Reworked psFitsInsertTable: previous version was not writing the

last column of the table; new version writes all columns.

Location:
trunk/psLib/src/fits
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/psLib/src/fits/psFits.c

    r7300 r7846  
    77 *  @author Robert DeSonia, MHPCC
    88 *
    9  *  @version $Revision: 1.55 $ $Name: not supported by cvs2svn $
    10  *  @date $Date: 2006-06-02 21:33:33 $
     9 *  @version $Revision: 1.56 $ $Name: not supported by cvs2svn $
     10 *  @date $Date: 2006-07-08 03:06:40 $
    1111 *
    1212 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    486486    return true;
    487487}
     488
     489
     490psDataType p_psFitsTypeFromCfitsio(int datatype)
     491{
     492    switch (datatype) {
     493    case TBYTE:
     494        return PS_TYPE_U8;
     495    case TSBYTE:
     496        return PS_TYPE_S8;
     497    case TSHORT:
     498        return PS_TYPE_S16;
     499    case TUSHORT:
     500        return PS_TYPE_U16;
     501    case TLONG:
     502        if (sizeof(long) == 8) {
     503            return PS_TYPE_S64;
     504        }
     505        // no break
     506    case TINT:
     507        return PS_TYPE_S32;
     508    case TULONG:
     509        if (sizeof(unsigned long) == 8) {
     510            return PS_TYPE_U64;
     511        }
     512        // no break
     513    case TUINT:
     514        return PS_TYPE_U32;
     515    case TLONGLONG:
     516        return PS_TYPE_S64;
     517    case TFLOAT:
     518        return PS_TYPE_F32;
     519    case TDOUBLE:
     520        return PS_TYPE_F64;
     521    case TCOMPLEX:
     522        return PS_TYPE_C32;
     523    case TDBLCOMPLEX:
     524        return PS_TYPE_C64;
     525    case TLOGICAL:
     526        return PS_TYPE_BOOL;
     527    default:
     528        psError(PS_ERR_IO, true,
     529                "Unknown FITS datatype, %d.",
     530                datatype);
     531        return 0;
     532    }
     533}
     534
     535bool p_psFitsTypeToCfitsio(psDataType type, int* bitPix, double* bZero, int* dataType)
     536{
     537    int bitpix;
     538    int datatype;
     539    double bzero = 0.0;
     540
     541    switch (type) {
     542
     543    case PS_TYPE_U8:
     544        bitpix = BYTE_IMG;
     545        datatype = TBYTE;
     546        break;
     547
     548    case PS_TYPE_S8:
     549        bitpix = BYTE_IMG;
     550        bzero = INT8_MIN;
     551        datatype = TSBYTE;
     552        break;
     553
     554    case PS_TYPE_U16:
     555        bitpix = SHORT_IMG;
     556        bzero = -1.0 * INT16_MIN;
     557        datatype = TUSHORT;
     558        break;
     559
     560    case PS_TYPE_S16:
     561        bitpix = SHORT_IMG;
     562        datatype = TSHORT;
     563        break;
     564
     565    case PS_TYPE_U32:
     566        bitpix = LONG_IMG;
     567        bzero = -1.0 * INT32_MIN;
     568        datatype = TUINT;
     569        break;
     570
     571    case PS_TYPE_S32:
     572        bitpix = LONG_IMG;
     573        datatype = TINT;
     574        break;
     575
     576    case PS_TYPE_F32:
     577        bitpix = FLOAT_IMG;
     578        datatype = TFLOAT;
     579        break;
     580
     581    case PS_TYPE_F64:
     582        bitpix = DOUBLE_IMG;
     583        datatype = TDOUBLE;
     584        break;
     585
     586    case PS_DATA_STRING:
     587        bitpix = BYTE_IMG;
     588        datatype = TSTRING;
     589        break;
     590
     591    case PS_DATA_BOOL:
     592        bitpix = BYTE_IMG;
     593        datatype = TLOGICAL;
     594        break;
     595
     596    default: {
     597            char* typeStr;
     598            PS_TYPE_NAME(typeStr,type);
     599            psError(PS_ERR_BAD_PARAMETER_TYPE, true,
     600                    PS_ERRORTEXT_psFits_TYPE_UNSUPPORTED,
     601                    typeStr);
     602            return false;
     603        }
     604    }
     605
     606    // pass back the requested parameters  (NULL parameters are not set, of course).
     607    if (bitPix != NULL) {
     608        *bitPix = bitpix;
     609    }
     610
     611    if (dataType != NULL) {
     612        *dataType = datatype;
     613    }
     614
     615    if (bZero != NULL) {
     616        *bZero = bzero;
     617    }
     618
     619    return true;
     620}
     621
     622

  • trunk/psLib/src/fits/psFits.h

    r6767 r7846  
    77 *  @author Robert DeSonia, MHPCC
    88 *
    9  *  @version $Revision: 1.24 $ $Name: not supported by cvs2svn $
    10  *  @date $Date: 2006-04-04 19:52:42 $
     9 *  @version $Revision: 1.25 $ $Name: not supported by cvs2svn $
     10 *  @date $Date: 2006-07-08 03:06:40 $
    1111 *
    1212 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    191191);
    192192
     193// Return the psLib type, given a cfitsio data type
     194psDataType p_psFitsTypeFromCfitsio(int datatype // cfitsio data type
     195                                  );
     196
     197// Return the cfitsio data type, given a psLib type
     198bool p_psFitsTypeToCfitsio(psDataType type, // psLib data type
     199                           int* bitPix, // The corresponding BITPIX (returned)
     200                           double* bZero, // The corresponding BZERO (returned)
     201                           int* dataType// The corresponding cfitsio data type (returned)
     202                          );
     203
    193204/// @}
    194205

  • trunk/psLib/src/fits/psFitsImage.c

    r7540 r7846  
    77 *  @author Robert DeSonia, MHPCC
    88 *
    9  *  @version $Revision: 1.10 $ $Name: not supported by cvs2svn $
    10  *  @date $Date: 2006-06-13 22:28:02 $
     9 *  @version $Revision: 1.11 $ $Name: not supported by cvs2svn $
     10 *  @date $Date: 2006-07-08 03:06:40 $
    1111 *
    1212 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    3131#define MAX_STRING_LENGTH 256  // maximum length string for FITS routines
    3232
    33 
    34 static bool convertPsTypeToFits(int type, int* bitPix, double* bZero, int* dataType)
    35 {
    36 
    37     int bitpix;
    38     int datatype;
    39     double bzero = 0.0;
    40 
    41     switch (type) {
    42 
    43     case PS_TYPE_U8:
    44         bitpix = BYTE_IMG;
    45         datatype = TBYTE;
    46         break;
    47 
    48     case PS_TYPE_S8:
    49         bitpix = BYTE_IMG;
    50         bzero = INT8_MIN;
    51         datatype = TSBYTE;
    52         break;
    53 
    54     case PS_TYPE_U16:
    55         bitpix = SHORT_IMG;
    56         bzero = -1.0 * INT16_MIN;
    57         datatype = TUSHORT;
    58         break;
    59 
    60     case PS_TYPE_S16:
    61         bitpix = SHORT_IMG;
    62         datatype = TSHORT;
    63         break;
    64 
    65     case PS_TYPE_U32:
    66         bitpix = LONG_IMG;
    67         bzero = -1.0 * INT32_MIN;
    68         datatype = TUINT;
    69         break;
    70 
    71     case PS_TYPE_S32:
    72         bitpix = LONG_IMG;
    73         datatype = TINT;
    74         break;
    75 
    76     case PS_TYPE_F32:
    77         bitpix = FLOAT_IMG;
    78         datatype = TFLOAT;
    79         break;
    80 
    81     case PS_TYPE_F64:
    82         bitpix = DOUBLE_IMG;
    83         datatype = TDOUBLE;
    84         break;
    85 
    86     case PS_DATA_STRING:
    87         bitpix = BYTE_IMG;
    88         datatype = TSTRING;
    89         break;
    90 
    91     case PS_DATA_BOOL:
    92         bitpix = BYTE_IMG;
    93         datatype = TLOGICAL;
    94         break;
    95 
    96     default: {
    97             char* typeStr;
    98             PS_TYPE_NAME(typeStr,type);
    99             psError(PS_ERR_BAD_PARAMETER_TYPE, true,
    100                     PS_ERRORTEXT_psFits_TYPE_UNSUPPORTED,
    101                     typeStr);
    102             return false;
    103         }
    104     }
    105 
    106     // pass the requested parameters  (NULL parameters are not set, of course).
    107     if (bitPix != NULL) {
    108         *bitPix = bitpix;
    109     }
    110 
    111     if (dataType != NULL) {
    112         *dataType = datatype;
    113     }
    114 
    115     if (bZero != NULL) {
    116         *bZero = bzero;
    117     }
    118 
    119     return true;
    120 }
    121 
    12233psImage* psFitsReadImage(const psFits* fits,    // the psFits object
    12334                         psRegion region, // the region in the FITS image to read
    … …  
    314225    double bZero;                       // Zero offset
    315226    int dataType;                       // cfitsio data type
    316     if (! convertPsTypeToFits(input->type.type, &bitPix, &bZero, &dataType) ) {
     227    if (! p_psFitsTypeToCfitsio(input->type.type, &bitPix, &bZero, &dataType) ) {
    317228        return false;
    318229    }
    … …  
    440351    double bZero;
    441352    int dataType;
    442     if (! convertPsTypeToFits(input->type.type, &bitPix, &bZero, &dataType) ) {
     353    if (! p_psFitsTypeToCfitsio(input->type.type, &bitPix, &bZero, &dataType) ) {
    443354        return false;
    444355    }

  • trunk/psLib/src/fits/psFitsTable.c

    r7538 r7846  
    77 *  @author Robert DeSonia, MHPCC
    88 *
    9  *  @version $Revision: 1.13 $ $Name: not supported by cvs2svn $
    10  *  @date $Date: 2006-06-13 22:05:23 $
     9 *  @version $Revision: 1.14 $ $Name: not supported by cvs2svn $
     10 *  @date $Date: 2006-07-08 03:06:40 $
    1111 *
    1212 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
    … …  
    2626#include "psFitsTable.h"
    2727#include "psFitsHeader.h"
     28#include "psAbort.h"
     29#include "psAssert.h"
    2830
    2931#define MAX_STRING_LENGTH 256  // maximum length string for FITS routines
    3032
    31 static psElemType convertFitsToPsType(int datatype)
    32 {
    33     switch (datatype) {
    34     case TBYTE:
    35         return PS_TYPE_U8;
    36     case TSBYTE:
    37         return PS_TYPE_S8;
    38     case TSHORT:
    39         return PS_TYPE_S16;
    40     case TUSHORT:
    41         return PS_TYPE_U16;
    42     case TLONG:
    43         if (sizeof(long) == 8) {
    44             return PS_TYPE_S64;
    45         }
    46         // no break
    47     case TINT:
    48         return PS_TYPE_S32;
    49     case TULONG:
    50         if (sizeof(unsigned long) == 8) {
    51             return PS_TYPE_U64;
    52         }
    53         // no break
    54     case TUINT:
    55         return PS_TYPE_U32;
    56     case TLONGLONG:
    57         return PS_TYPE_S64;
    58     case TFLOAT:
    59         return PS_TYPE_F32;
    60     case TDOUBLE:
    61         return PS_TYPE_F64;
    62     case TCOMPLEX:
    63         return PS_TYPE_C32;
    64     case TDBLCOMPLEX:
    65         return PS_TYPE_C64;
    66     case TLOGICAL:
    67         return PS_TYPE_BOOL;
     33
     34psMetadata* psFitsReadTableRow(const psFits* fits,
     35                               int row)
     36{
     37    long numRows;
     38    int numCols;
     39    int status = 0;
     40
     41    if (fits == NULL) {
     42        psError(PS_ERR_BAD_PARAMETER_NULL, true,
     43                PS_ERRORTEXT_psFits_NULL);
     44        return NULL;
     45    }
     46
     47    // check to see if we even are positioned on a table HDU
     48    int hdutype;
     49    fits_get_hdu_type(fits->fd,&hdutype, &status);
     50    if ( status != 0) {
     51        char fitsErr[MAX_STRING_LENGTH];
     52        (void)fits_get_errstatus(status, fitsErr);
     53        psError(PS_ERR_IO, true,
     54                PS_ERRORTEXT_psFits_GET_HDU_TYPE_FAILED,
     55                fitsErr);
     56        return NULL;
     57    }
     58    if (hdutype != ASCII_TBL && hdutype != BINARY_TBL) {
     59        psError(PS_ERR_IO, true,
     60                PS_ERRORTEXT_psFits_NOT_TABLE_TYPE);
     61        return NULL;
     62    }
     63
     64    // get the size of the FITS table
     65    fits_get_num_rows(fits->fd, &numRows, &status);
     66    fits_get_num_cols(fits->fd, &numCols, &status);
     67    if ( status != 0) {
     68        char fitsErr[MAX_STRING_LENGTH];
     69        (void)fits_get_errstatus(status, fitsErr);
     70        psError(PS_ERR_IO, true,
     71                PS_ERRORTEXT_psFits_GET_TABLE_SIZE_FAILED,
     72                fitsErr);
     73        return NULL;
     74    }
     75
     76    psTrace(".psFits.psFitsReadTableRow",5,"Table size is %ix%i\n",numCols, numRows);
     77    // the row parameter in the proper range?
     78    if (row < 0 || row >= numRows) {
     79        psError(PS_ERR_BAD_PARAMETER_VALUE, true,
     80                PS_ERRORTEXT_psFits_ROW_INVALID,
     81                row,numRows);
     82        return NULL;
     83    }
     84
     85    psMetadata* data = psMetadataAlloc();
     86
     87    int typecode;
     88    long repeat;
     89    long width;
     90    char name[60];
     91    for (int col = 1; col <= numCols; col++) {
     92        // get the column name
     93        if (hdutype == BINARY_TBL) {
     94            fits_get_bcolparms(fits->fd, col, name,
     95                               NULL, NULL, NULL, NULL, NULL, NULL, NULL, &status);
     96        } else {
     97            fits_get_acolparms(fits->fd, col, name,
     98                               NULL, NULL, NULL, NULL, NULL, NULL, NULL, &status);
     99        }
     100        // get the column type
     101        fits_get_coltype(fits->fd, col, &typecode, &repeat, &width, &status);
     102
     103        if (status == 0) {
     104
     105            #define READ_TABLE_ROW_CASE(FITSTYPE, NATIVETYPE, TYPE, VECTYPE) \
     106        case FITSTYPE: { \
     107                if (repeat == 1) { \
     108                    NATIVETYPE value; \
     109                    int anynul = 0; \
     110                    fits_read_col(fits->fd, FITSTYPE, col,row+1, \
     111                                  1, 1, NULL, &value, &anynul, &status); \
     112                    psTrace(".psFits.psFitsReadTableRow",5,"Column #%i, '%s', is type %i, repeat %i, Value = %g\n", \
     113                            col, name, typecode, repeat, (double)value); \
     114                    psMetadataAdd(data,PS_LIST_TAIL, name, \
     115                                  PS_DATA_##TYPE, \
     116                                  "", (ps##TYPE)value); \
     117                } else { \
     118                    NATIVETYPE* value = psAlloc(sizeof(NATIVETYPE)*repeat); \
     119                    psVector* vec = psVectorAlloc(repeat,PS_TYPE_##VECTYPE); \
     120                    int anynul = 0; \
     121                    fits_read_col(fits->fd, FITSTYPE, col,row+1, \
     122                                  1, repeat, NULL, value, &anynul, &status); \
     123                    for (int lcv = 0; lcv < repeat; lcv++) { \
     124                        vec->data.VECTYPE[lcv] = value[lcv]; \
     125                    } \
     126                    psMetadataAdd(data,PS_LIST_TAIL, name, PS_DATA_VECTOR, "", vec); \
     127                    psFree(value); \
     128                    psFree(vec); \
     129                } \
     130                break; \
     131            }
     132
     133            switch (typecode) {
     134            case TBYTE:
     135            case TSHORT:
     136            case TLONGLONG:
     137                READ_TABLE_ROW_CASE(TLONG, long, S32,S32)
     138                READ_TABLE_ROW_CASE(TFLOAT, float, F32,F32)
     139                READ_TABLE_ROW_CASE(TDOUBLE, double, F64,F64)
     140                READ_TABLE_ROW_CASE(TLOGICAL, bool, BOOL,S8);
     141            case TSTRING: {
     142                    char* value = psAlloc(repeat+1);
     143                    int anynul = 0;
     144                    fits_read_col(fits->fd, TSTRING, col,row+1,
     145                                  1, 1, NULL, &value, &anynul, &status);
     146                    psTrace(".psFits.psFitsReadTableRow",5,"Column #%i, '%s', is type %i, repeat %i, value = %s\n",
     147                            col, name, typecode, repeat, value);
     148                    if (anynul == 0) {
     149                        psMetadataAdd(data,PS_LIST_TAIL, name,
     150                                      PS_DATA_STRING,
     151                                      "", value);
     152                    }
     153                    psFree(value);
     154                    break;
     155                }
     156            default:
     157                psWarning("Data type (%d) not supportted for column %d",
     158                          typecode, col);
     159
     160                psTrace("psFits.psFitsReadTableRow", 2,
     161                        "Column %d or row %d was of a non primitive type, %d",
     162                        col, row, typecode);
     163            }
     164        }
     165
     166        if ( status != 0) {
     167            char fitsErr[MAX_STRING_LENGTH];
     168            (void)fits_get_errstatus(status, fitsErr);
     169            psError(PS_ERR_IO, true,
     170                    PS_ERRORTEXT_psFits_GET_TABLE_ELEMENT,
     171                    col,row,fitsErr);
     172            psFree(data);
     173            return NULL;
     174        }
     175
     176    }
     177
     178    return data;
     179}
     180
     181psArray* psFitsReadTableColumn(const psFits* fits,
     182                               const char* colname)
     183{
     184    int colnum = 0;
     185    int status = 0;
     186
     187    if (fits == NULL) {
     188        psError(PS_ERR_BAD_PARAMETER_NULL, true,
     189                PS_ERRORTEXT_psFits_NULL);
     190        return NULL;
     191    }
     192
     193    // check to see if we even are positioned on a table HDU
     194    int hdutype;
     195    if ( fits_get_hdu_type(fits->fd,&hdutype, &status) != 0) {
     196        char fitsErr[MAX_STRING_LENGTH];
     197        (void)fits_get_errstatus(status, fitsErr);
     198        psError(PS_ERR_IO, true,
     199                PS_ERRORTEXT_psFits_GET_HDU_TYPE_FAILED,
     200                fitsErr);
     201        return NULL;
     202    }
     203    if (hdutype != ASCII_TBL && hdutype != BINARY_TBL) {
     204        psError(PS_ERR_IO, true,
     205                PS_ERRORTEXT_psFits_NOT_TABLE_TYPE);
     206        return NULL;
     207    }
     208
     209    // find the column by name
     210    if ( fits_get_colnum(fits->fd, CASESEN, (char*)colname, &colnum, &status) != 0) {
     211        char fitsErr[MAX_STRING_LENGTH];
     212        (void)fits_get_errstatus(status, fitsErr);
     213        psError(PS_ERR_IO, true,
     214                PS_ERRORTEXT_psFits_FIND_COLUMN,
     215                colname, fitsErr);
     216        return NULL;
     217    }
     218
     219    // get the number of rows
     220    long numRows = 0;
     221    fits_get_num_rows(fits->fd, &numRows, &status);
     222
     223    // get the column length.
     224    int width;
     225    if ( fits_get_col_display_width(fits->fd, colnum, &width, &status) != 0) {
     226        char fitsErr[MAX_STRING_LENGTH];
     227        (void)fits_get_errstatus(status, fitsErr);
     228        psError(PS_ERR_IO, true,
     229                PS_ERRORTEXT_psFits_GET_COLTYPE,
     230                fitsErr);
     231        return NULL;
     232    }
     233
     234    // allocate the buffers
     235    psArray* result = psArrayAlloc(numRows);
     236    for (int row = 0; row < numRows; row++) {
     237        result->data[row] = psAlloc((width+1)*sizeof(char));
     238    }
     239    result->n = numRows;
     240
     241    fits_read_col_str(fits->fd,
     242                      colnum,
     243                      1, // firstrow
     244                      1, // firestelem
     245                      numRows,
     246                      "", // nulstr
     247                      (char**)result->data,
     248                      NULL,
     249                      &status);
     250
     251    if ( status != 0) {
     252        char fitsErr[MAX_STRING_LENGTH];
     253        (void)fits_get_errstatus(status, fitsErr);
     254        psError(PS_ERR_IO, true,
     255                PS_ERRORTEXT_psFits_TABLE_READ_COL,
     256                fitsErr);
     257        return NULL;
     258    }
     259
     260    return result;
     261}
     262
     263psVector* psFitsReadTableColumnNum(const psFits* fits,
     264                                   const char* colname)
     265{
     266    int status = 0;
     267    int colnum = 0;
     268
     269    if (fits == NULL) {
     270        psError(PS_ERR_BAD_PARAMETER_NULL, true,
     271                PS_ERRORTEXT_psFits_NULL);
     272        return NULL;
     273    }
     274
     275    // check to see if we even are positioned on a table HDU
     276    int hdutype;
     277    if ( fits_get_hdu_type(fits->fd,&hdutype, &status) != 0) {
     278        char fitsErr[MAX_STRING_LENGTH];
     279        (void)fits_get_errstatus(status, fitsErr);
     280        psError(PS_ERR_IO, true,
     281                PS_ERRORTEXT_psFits_GET_HDU_TYPE_FAILED,
     282                fitsErr);
     283        return NULL;
     284    }
     285    if (hdutype != ASCII_TBL && hdutype != BINARY_TBL) {
     286        psError(PS_ERR_IO, true,
     287                PS_ERRORTEXT_psFits_NOT_TABLE_TYPE);
     288        return NULL;
     289    }
     290
     291    // find the column by name
     292    if ( fits_get_colnum(fits->fd, CASESEN, (char*)colname, &colnum, &status) != 0) {
     293        char fitsErr[MAX_STRING_LENGTH];
     294        (void)fits_get_errstatus(status, fitsErr);
     295        psError(PS_ERR_IO, true,
     296                PS_ERRORTEXT_psFits_FIND_COLUMN,
     297                colname, fitsErr);
     298        return NULL;
     299    }
     300
     301    // get the number of rows
     302    long numRows = 0;
     303    fits_get_num_rows(fits->fd,
     304                      &numRows,
     305                      &status);
     306
     307    // get the column datatype.
     308    int typecode;
     309    long repeat;
     310    long width;
     311    if ( fits_get_eqcoltype(fits->fd, colnum, &typecode, &repeat, &width, &status) != 0) {
     312        char fitsErr[MAX_STRING_LENGTH];
     313        (void)fits_get_errstatus(status, fitsErr);
     314        psError(PS_ERR_IO, true,
     315                PS_ERRORTEXT_psFits_GET_COLTYPE,
     316                fitsErr);
     317        return NULL;
     318    }
     319
     320    psVector* result = psVectorAlloc(numRows, p_psFitsTypeFromCfitsio(typecode));
     321    result->n = numRows;
     322
     323    fits_read_col(fits->fd,
     324                  typecode,
     325                  colnum,
     326                  1 /* firstrow */,
     327                  1 /* firstelem */,
     328                  numRows,
     329                  NULL,
     330                  (psPtr)(result->data.U8),
     331                  NULL,
     332                  &status);
     333
     334    if ( status != 0) {
     335        char fitsErr[MAX_STRING_LENGTH];
     336        (void)fits_get_errstatus(status, fitsErr);
     337        psError(PS_ERR_IO, true,
     338                PS_ERRORTEXT_psFits_TABLE_READ_COL,
     339                fitsErr);
     340        return NULL;
     341    }
     342
     343    return result;
     344}
     345
     346
     347psArray* psFitsReadTable(const psFits* fits)
     348{
     349    int status = 0;
     350
     351    if (fits == NULL) {
     352        psError(PS_ERR_BAD_PARAMETER_NULL, true,
     353                PS_ERRORTEXT_psFits_NULL);
     354        return NULL;
     355    }
     356
     357    // check to see if we even are positioned on a table HDU
     358    int hdutype;
     359    if ( fits_get_hdu_type(fits->fd,&hdutype, &status) != 0) {
     360        char fitsErr[MAX_STRING_LENGTH];
     361        (void)fits_get_errstatus(status, fitsErr);
     362        psError(PS_ERR_IO, true,
     363                PS_ERRORTEXT_psFits_GET_HDU_TYPE_FAILED,
     364                fitsErr);
     365        return NULL;
     366    }
     367    if (hdutype != ASCII_TBL && hdutype != BINARY_TBL) {
     368        psError(PS_ERR_IO, true,
     369                PS_ERRORTEXT_psFits_NOT_TABLE_TYPE);
     370        return NULL;
     371    }
     372
     373    // get the size of the FITS table
     374    long numRows = 0;
     375    fits_get_num_rows(fits->fd, &numRows, &status);
     376    if ( status != 0) {
     377        char fitsErr[MAX_STRING_LENGTH];
     378        (void)fits_get_errstatus(status, fitsErr);
     379        psError(PS_ERR_IO, true,
     380                PS_ERRORTEXT_psFits_GET_TABLE_SIZE_FAILED,
     381                fitsErr);
     382        return NULL;
     383    }
     384
     385
     386    psArray* table = psArrayAlloc(numRows);
     387
     388    for (int row = 0; row < numRows; row++) {
     389        psTrace(".psFits.psFitsReadTable",5,"Reading row %i of %i\n",row, numRows);
     390        table->data[row] = psFitsReadTableRow(fits,row);
     391        table->n++;
     392    }
     393
     394    return table;
     395}
     396
     397bool psFitsWriteTable(psFits* fits,
     398                      const psMetadata* header,
     399                      const psArray* table,
     400                      const char *extname)
     401{
     402    psFitsMoveLast(fits);
     403    return psFitsInsertTable(fits, header, table, extname, true);
     404}
     405
     406
     407// Return the size of the column
     408static inline size_t columnSize(const psMetadataItem *item // Item for which to get the size
     409                               )
     410{
     411    if (PS_DATA_IS_PRIMITIVE(item->type)) {
     412        return 1;
     413    }
     414    switch (item->type) {
     415    case PS_DATA_STRING:
     416        return strlen(item->data.V);
     417    case PS_DATA_VECTOR: {
     418            psVector *vector = item->data.V;
     419            return vector->n;
     420        }
    68421    default:
    69         psError(PS_ERR_IO, true,
    70                 "Unknown FITS datatype, %d.",
    71                 datatype);
    72         return 0;
    73     }
    74 }
    75 
    76 static bool convertPsTypeToFits(int type, int* bitPix, double* bZero, int* dataType)
    77 {
    78 
    79     int bitpix;
    80     int datatype;
    81     double bzero = 0.0;
    82 
    83     switch (type) {
    84 
    85     case PS_TYPE_U8:
    86         bitpix = BYTE_IMG;
    87         datatype = TBYTE;
    88         break;
    89 
    90     case PS_TYPE_S8:
    91         bitpix = BYTE_IMG;
    92         bzero = INT8_MIN;
    93         datatype = TSBYTE;
    94         break;
    95 
    96     case PS_TYPE_U16:
    97         bitpix = SHORT_IMG;
    98         bzero = -1.0 * INT16_MIN;
    99         datatype = TUSHORT;
    100         break;
    101 
    102     case PS_TYPE_S16:
    103         bitpix = SHORT_IMG;
    104         datatype = TSHORT;
    105         break;
    106 
    107     case PS_TYPE_U32:
    108         bitpix = LONG_IMG;
    109         bzero = -1.0 * INT32_MIN;
    110         datatype = TUINT;
    111         break;
    112 
    113     case PS_TYPE_S32:
    114         bitpix = LONG_IMG;
    115         datatype = TINT;
    116         break;
    117 
    118     case PS_TYPE_F32:
    119         bitpix = FLOAT_IMG;
    120         datatype = TFLOAT;
    121         break;
    122 
    123     case PS_TYPE_F64:
    124         bitpix = DOUBLE_IMG;
    125         datatype = TDOUBLE;
    126         break;
    127 
    128     case PS_DATA_STRING:
    129         bitpix = BYTE_IMG;
    130         datatype = TSTRING;
    131         break;
    132 
    133     case PS_DATA_BOOL:
    134         bitpix = BYTE_IMG;
    135         datatype = TLOGICAL;
    136         break;
    137 
    138     default: {
    139             char* typeStr;
    140             PS_TYPE_NAME(typeStr,type);
    141             psError(PS_ERR_BAD_PARAMETER_TYPE, true,
    142                     PS_ERRORTEXT_psFits_TYPE_UNSUPPORTED,
    143                     typeStr);
    144             return false;
    145         }
    146     }
    147 
    148     // pass the requested parameters  (NULL parameters are not set, of course).
    149     if (bitPix != NULL) {
    150         *bitPix = bitpix;
    151     }
    152 
    153     if (dataType != NULL) {
    154         *dataType = datatype;
    155     }
    156 
    157     if (bZero != NULL) {
    158         *bZero = bzero;
    159     }
    160 
    161     return true;
    162 }
    163 
    164 static char getTForm(int type)
     422        psAbort(__func__, "Shouldn't ever get here.");
     423    }
     424    return 0;
     425}
     426
     427// Get the TFORM character, given a PS type
     428static inline const char getTForm(psDataType type)
    165429{
    166430    switch (type) {
    … …  
    196460}
    197461
    198 static bool getMetadataTForm(psMetadataItem* item, char** fitsType, int repeat)
    199 {
    200     if (item == NULL) {
    201         return false;
    202     }
    203 
    204     *fitsType = NULL;
    205     char tform = getTForm(item->type);
    206 
    207     switch (item->type) {
    208     case PS_TYPE_U8:
    209     case PS_TYPE_S8:
    210     case PS_TYPE_S16:
    211     case PS_TYPE_S32:
    212     case PS_TYPE_U64:
    213     case PS_TYPE_S64:
    214     case PS_TYPE_U16:
    215     case PS_TYPE_U32:
    216     case PS_TYPE_F32:
    217     case PS_TYPE_F64:
    218     case PS_TYPE_C32:
    219     case PS_TYPE_C64:
    220     case PS_TYPE_BOOL:
    221         psStringAppend(fitsType, "%d%c",repeat,tform);
    222         break;
    223 
    224     case PS_DATA_STRING:
    225         psStringAppend(fitsType,"%dA",
    226                        strlen(item->data.V)*repeat);
    227         break;
    228 
    229     case PS_DATA_VECTOR: {
    230             psVector* vec = item->data.V;
    231             tform = getTForm(vec->type.type);
    232             psStringAppend(fitsType,"%d%c",repeat*vec->n,tform);
    233         }
    234         break;
    235 
    236     default:
    237         return false;
    238     }
    239 
    240     return true;
    241 }
    242 
    243 psMetadata* psFitsReadTableRow(const psFits* fits,
    244                                int row)
    245 {
    246     long numRows;
    247     int numCols;
    248     int status = 0;
    249 
    250     if (fits == NULL) {
    251         psError(PS_ERR_BAD_PARAMETER_NULL, true,
    252                 PS_ERRORTEXT_psFits_NULL);
    253         return NULL;
    254     }
    255 
    256     // check to see if we even are positioned on a table HDU
    257     int hdutype;
    258     fits_get_hdu_type(fits->fd,&hdutype, &status);
    259     if ( status != 0) {
    260         char fitsErr[MAX_STRING_LENGTH];
    261         (void)fits_get_errstatus(status, fitsErr);
    262         psError(PS_ERR_IO, true,
    263                 PS_ERRORTEXT_psFits_GET_HDU_TYPE_FAILED,
    264                 fitsErr);
    265         return NULL;
    266     }
    267     if (hdutype != ASCII_TBL && hdutype != BINARY_TBL) {
    268         psError(PS_ERR_IO, true,
    269                 PS_ERRORTEXT_psFits_NOT_TABLE_TYPE);
    270         return NULL;
    271     }
    272 
    273     // get the size of the FITS table
    274     fits_get_num_rows(fits->fd, &numRows, &status);
    275     fits_get_num_cols(fits->fd, &numCols, &status);
    276     if ( status != 0) {
    277         char fitsErr[MAX_STRING_LENGTH];
    278         (void)fits_get_errstatus(status, fitsErr);
    279         psError(PS_ERR_IO, true,
    280                 PS_ERRORTEXT_psFits_GET_TABLE_SIZE_FAILED,
    281                 fitsErr);
    282         return NULL;
    283     }
    284 
    285     psTrace(".psFits.psFitsReadTableRow",5,"Table size is %ix%i\n",numCols, numRows);
    286     // the row parameter in the proper range?
    287     if (row < 0 || row >= numRows) {
    288         psError(PS_ERR_BAD_PARAMETER_VALUE, true,
    289                 PS_ERRORTEXT_psFits_ROW_INVALID,
    290                 row,numRows);
    291         return NULL;
    292     }
    293 
    294     psMetadata* data = psMetadataAlloc();
    295 
    296     int typecode;
    297     long repeat;
    298     long width;
    299     char name[60];
    300     for (int col = 1; col <= numCols; col++) {
    301         // get the column name
    302         if (hdutype == BINARY_TBL) {
    303             fits_get_bcolparms(fits->fd, col, name,
    304                                NULL, NULL, NULL, NULL, NULL, NULL, NULL, &status);
    305         } else {
    306             fits_get_acolparms(fits->fd, col, name,
    307                                NULL, NULL, NULL, NULL, NULL, NULL, NULL, &status);
    308         }
    309         // get the column type
    310         fits_get_coltype(fits->fd, col, &typecode, &repeat, &width, &status);
    311 
    312         if (status == 0) {
    313 
    314             #define READ_TABLE_ROW_CASE(FITSTYPE, NATIVETYPE, TYPE, VECTYPE) \
    315         case FITSTYPE: { \
    316                 if (repeat == 1) { \
    317                     NATIVETYPE value; \
    318                     int anynul = 0; \
    319                     fits_read_col(fits->fd, FITSTYPE, col,row+1, \
    320                                   1, 1, NULL, &value, &anynul, &status); \
    321                     psTrace(".psFits.psFitsReadTableRow",5,"Column #%i, '%s', is type %i, repeat %i, Value = %g\n", \
    322                             col, name, typecode, repeat, (double)value); \
    323                     psMetadataAdd(data,PS_LIST_TAIL, name, \
    324                                   PS_DATA_##TYPE, \
    325                                   "", (ps##TYPE)value); \
    326                 } else { \
    327                     NATIVETYPE* value = psAlloc(sizeof(NATIVETYPE)*repeat); \
    328                     psVector* vec = psVectorAlloc(repeat,PS_TYPE_##VECTYPE); \
    329                     int anynul = 0; \
    330                     fits_read_col(fits->fd, FITSTYPE, col,row+1, \
    331                                   1, repeat, NULL, value, &anynul, &status); \
    332                     for (int lcv = 0; lcv < repeat; lcv++) { \
    333                         vec->data.VECTYPE[lcv] = value[lcv]; \
    334                     } \
    335                     psMetadataAdd(data,PS_LIST_TAIL, name, PS_DATA_VECTOR, "", vec); \
    336                     psFree(value); \
    337                     psFree(vec); \
    338                 } \
    339                 break; \
    340             }
    341 
    342             switch (typecode) {
    343             case TBYTE:
    344             case TSHORT:
    345             case TLONGLONG:
    346                 READ_TABLE_ROW_CASE(TLONG, long, S32,S32)
    347                 READ_TABLE_ROW_CASE(TFLOAT, float, F32,F32)
    348                 READ_TABLE_ROW_CASE(TDOUBLE, double, F64,F64)
    349                 READ_TABLE_ROW_CASE(TLOGICAL, bool, BOOL,S8);
    350             case TSTRING: {
    351                     char* value = psAlloc(repeat+1);
    352                     int anynul = 0;
    353                     fits_read_col(fits->fd, TSTRING, col,row+1,
    354                                   1, 1, NULL, &value, &anynul, &status);
    355                     psTrace(".psFits.psFitsReadTableRow",5,"Column #%i, '%s', is type %i, repeat %i, value = %s\n",
    356                             col, name, typecode, repeat, value);
    357                     if (anynul == 0) {
    358                         psMetadataAdd(data,PS_LIST_TAIL, name,
    359                                       PS_DATA_STRING,
    360                                       "", value);
    361                     }
    362                     psFree(value);
    363                     break;
    364                 }
    365             default:
    366                 psWarning("Data type (%d) not supportted for column %d",
    367                           typecode, col);
    368 
    369                 psTrace("psFits.psFitsReadTableRow", 2,
    370                         "Column %d or row %d was of a non primitive type, %d",
    371                         col, row, typecode);
    372             }
    373         }
    374 
    375         if ( status != 0) {
    376             char fitsErr[MAX_STRING_LENGTH];
    377             (void)fits_get_errstatus(status, fitsErr);
    378             psError(PS_ERR_IO, true,
    379                     PS_ERRORTEXT_psFits_GET_TABLE_ELEMENT,
    380                     col,row,fitsErr);
    381             psFree(data);
    382             return NULL;
    383         }
    384 
    385     }
    386 
    387     return data;
    388 }
    389 
    390 psArray* psFitsReadTableColumn(const psFits* fits,
    391                                const char* colname)
    392 {
    393     int colnum = 0;
    394     int status = 0;
    395 
    396     if (fits == NULL) {
    397         psError(PS_ERR_BAD_PARAMETER_NULL, true,
    398                 PS_ERRORTEXT_psFits_NULL);
    399         return NULL;
    400     }
    401 
    402     // check to see if we even are positioned on a table HDU
    403     int hdutype;
    404     if ( fits_get_hdu_type(fits->fd,&hdutype, &status) != 0) {
    405         char fitsErr[MAX_STRING_LENGTH];
    406         (void)fits_get_errstatus(status, fitsErr);
    407         psError(PS_ERR_IO, true,
    408                 PS_ERRORTEXT_psFits_GET_HDU_TYPE_FAILED,
    409                 fitsErr);
    410         return NULL;
    411     }
    412     if (hdutype != ASCII_TBL && hdutype != BINARY_TBL) {
    413         psError(PS_ERR_IO, true,
    414                 PS_ERRORTEXT_psFits_NOT_TABLE_TYPE);
    415         return NULL;
    416     }
    417 
    418     // find the column by name
    419     if ( fits_get_colnum(fits->fd, CASESEN, (char*)colname, &colnum, &status) != 0) {
    420         char fitsErr[MAX_STRING_LENGTH];
    421         (void)fits_get_errstatus(status, fitsErr);
    422         psError(PS_ERR_IO, true,
    423                 PS_ERRORTEXT_psFits_FIND_COLUMN,
    424                 colname, fitsErr);
    425         return NULL;
    426     }
    427 
    428     // get the number of rows
    429     long numRows = 0;
    430     fits_get_num_rows(fits->fd, &numRows, &status);
    431 
    432     // get the column length.
    433     int width;
    434     if ( fits_get_col_display_width(fits->fd, colnum, &width, &status) != 0) {
    435         char fitsErr[MAX_STRING_LENGTH];
    436         (void)fits_get_errstatus(status, fitsErr);
    437         psError(PS_ERR_IO, true,
    438                 PS_ERRORTEXT_psFits_GET_COLTYPE,
    439                 fitsErr);
    440         return NULL;
    441     }
    442 
    443     // allocate the buffers
    444     psArray* result = psArrayAlloc(numRows);
    445     for (int row = 0; row < numRows; row++) {
    446         result->data[row] = psAlloc((width+1)*sizeof(char));
    447     }
    448     result->n = numRows;
    449 
    450     fits_read_col_str(fits->fd,
    451                       colnum,
    452                       1, // firstrow
    453                       1, // firestelem
    454                       numRows,
    455                       "", // nulstr
    456                       (char**)result->data,
    457                       NULL,
    458                       &status);
    459 
    460     if ( status != 0) {
    461         char fitsErr[MAX_STRING_LENGTH];
    462         (void)fits_get_errstatus(status, fitsErr);
    463         psError(PS_ERR_IO, true,
    464                 PS_ERRORTEXT_psFits_TABLE_READ_COL,
    465                 fitsErr);
    466         return NULL;
    467     }
    468 
    469     return result;
    470 }
    471 
    472 psVector* psFitsReadTableColumnNum(const psFits* fits,
    473                                    const char* colname)
    474 {
    475     int status = 0;
    476     int colnum = 0;
    477 
    478     if (fits == NULL) {
    479         psError(PS_ERR_BAD_PARAMETER_NULL, true,
    480                 PS_ERRORTEXT_psFits_NULL);
    481         return NULL;
    482     }
    483 
    484     // check to see if we even are positioned on a table HDU
    485     int hdutype;
    486     if ( fits_get_hdu_type(fits->fd,&hdutype, &status) != 0) {
    487         char fitsErr[MAX_STRING_LENGTH];
    488         (void)fits_get_errstatus(status, fitsErr);
    489         psError(PS_ERR_IO, true,
    490                 PS_ERRORTEXT_psFits_GET_HDU_TYPE_FAILED,
    491                 fitsErr);
    492         return NULL;
    493     }
    494     if (hdutype != ASCII_TBL && hdutype != BINARY_TBL) {
    495         psError(PS_ERR_IO, true,
    496                 PS_ERRORTEXT_psFits_NOT_TABLE_TYPE);
    497         return NULL;
    498     }
    499 
    500     // find the column by name
    501     if ( fits_get_colnum(fits->fd, CASESEN, (char*)colname, &colnum, &status) != 0) {
    502         char fitsErr[MAX_STRING_LENGTH];
    503         (void)fits_get_errstatus(status, fitsErr);
    504         psError(PS_ERR_IO, true,
    505                 PS_ERRORTEXT_psFits_FIND_COLUMN,
    506                 colname, fitsErr);
    507         return NULL;
    508     }
    509 
    510     // get the number of rows
    511     long numRows = 0;
    512     fits_get_num_rows(fits->fd,
    513                       &numRows,
    514                       &status);
    515 
    516     // get the column datatype.
    517     int typecode;
    518     long repeat;
    519     long width;
    520     if ( fits_get_eqcoltype(fits->fd, colnum, &typecode, &repeat, &width, &status) != 0) {
    521         char fitsErr[MAX_STRING_LENGTH];
    522         (void)fits_get_errstatus(status, fitsErr);
    523         psError(PS_ERR_IO, true,
    524                 PS_ERRORTEXT_psFits_GET_COLTYPE,
    525                 fitsErr);
    526         return NULL;
    527     }
    528 
    529     psVector* result = psVectorAlloc(numRows, convertFitsToPsType(typecode));
    530     result->n = numRows;
    531 
    532     fits_read_col(fits->fd,
    533                   typecode,
    534                   colnum,
    535                   1 /* firstrow */,
    536                   1 /* firstelem */,
    537                   numRows,
    538                   NULL,
    539                   (psPtr)(result->data.U8),
    540                   NULL,
    541                   &status);
    542 
    543     if ( status != 0) {
    544         char fitsErr[MAX_STRING_LENGTH];
    545         (void)fits_get_errstatus(status, fitsErr);
    546         psError(PS_ERR_IO, true,
    547                 PS_ERRORTEXT_psFits_TABLE_READ_COL,
    548                 fitsErr);
    549         return NULL;
    550     }
    551 
    552     return result;
    553 }
    554 
    555 
    556 psArray* psFitsReadTable(const psFits* fits)
    557 {
    558     int status = 0;
    559 
    560     if (fits == NULL) {
    561         psError(PS_ERR_BAD_PARAMETER_NULL, true,
    562                 PS_ERRORTEXT_psFits_NULL);
    563         return NULL;
    564     }
    565 
    566     // check to see if we even are positioned on a table HDU
    567     int hdutype;
    568     if ( fits_get_hdu_type(fits->fd,&hdutype, &status) != 0) {
    569         char fitsErr[MAX_STRING_LENGTH];
    570         (void)fits_get_errstatus(status, fitsErr);
    571         psError(PS_ERR_IO, true,
    572                 PS_ERRORTEXT_psFits_GET_HDU_TYPE_FAILED,
    573                 fitsErr);
    574         return NULL;
    575     }
    576     if (hdutype != ASCII_TBL && hdutype != BINARY_TBL) {
    577         psError(PS_ERR_IO, true,
    578                 PS_ERRORTEXT_psFits_NOT_TABLE_TYPE);
    579         return NULL;
    580     }
    581 
    582     // get the size of the FITS table
    583     long numRows = 0;
    584     fits_get_num_rows(fits->fd, &numRows, &status);
    585     if ( status != 0) {
    586         char fitsErr[MAX_STRING_LENGTH];
    587         (void)fits_get_errstatus(status, fitsErr);
    588         psError(PS_ERR_IO, true,
    589                 PS_ERRORTEXT_psFits_GET_TABLE_SIZE_FAILED,
    590                 fitsErr);
    591         return NULL;
    592     }
    593 
    594 
    595     psArray* table = psArrayAlloc(numRows);
    596 
    597     for (int row = 0; row < numRows; row++) {
    598         psTrace(".psFits.psFitsReadTable",5,"Reading row %i of %i\n",row, numRows);
    599         table->data[row] = psFitsReadTableRow(fits,row);
    600         table->n++;
    601     }
    602 
    603     return table;
    604 }
    605 
    606 bool psFitsWriteTable(psFits* fits,
    607                       const psMetadata* header,
    608                       const psArray* table,
    609                       const char *extname)
    610 {
    611     psFitsMoveLast(fits);
    612     return psFitsInsertTable(fits, header, table, extname, true);
    613 }
     462
     463// Column specification
     464// Included here, because there's no need for the user to have access to it
     465typedef struct
     466{
     467    psDataType type;                    // psLib type (e.g., PS_DATA_STRING or PS_TYPE_F32)
     468    size_t size;                        // Size (number of repeats)
     469    psElemType vectorType;              // psLib type of vectors
     470}
     471colSpec;
     472
    614473
    615474bool psFitsInsertTable(psFits* fits,
    … …  
    621480    int status = 0;
    622481
    623     if (fits == NULL) {
    624         psError(PS_ERR_BAD_PARAMETER_NULL, true,
    625                 PS_ERRORTEXT_psFits_NULL);
    626         return false;
    627     }
    628 
    629     if (table == NULL) {
    630         psError(PS_ERR_BAD_PARAMETER_NULL, true,
    631                 PS_ERRORTEXT_psFits_IMAGE_NULL);
    632         return false;
    633     }
    634 
    635     int rows = table->n;
    636     if (rows < 1) {
     482    PS_ASSERT_PTR_NON_NULL(fits, false);
     483    PS_ASSERT_ARRAY_NON_NULL(table, false);
     484
     485    long numRows = table->n;
     486    if (numRows < 1) {
    637487        // no table data, what can I do?
    638488        psError(PS_ERR_BAD_PARAMETER_SIZE, true,
    … …  
    641491    }
    642492
    643     // find all the columns needed
    644     psArray* columns = psArrayAlloc(((psMetadata*)table->data[0])->list->n);
    645     columns->n=0;
    646 
    647     // find the unique items in the array of metadata 'rows'
    648     psMetadataItem* item;
    649     for (int row=0; row < rows; row++) {
    650         psMetadata* rowMeta = table->data[row];
    651         if (rowMeta != NULL) {
    652             psListIterator* iter = psListIteratorAlloc(rowMeta->list,
    653                                    PS_LIST_HEAD,true);
    654             while ( (item=psListGetAndIncrement(iter)) != NULL) {
    655                 if (PS_DATA_IS_PRIMITIVE(item->type) ||
    656                         item->type == PS_DATA_STRING ||
    657                         item->type == PS_DATA_VECTOR) {
    658                     bool found = false;
    659                     psMetadataItem* fItem = NULL;
    660                     int n;
    661                     for (n=0; n < columns->n && ! found; n++) {
    662                         fItem = (psMetadataItem*)(columns->data[n]);
    663                         if (strcmp(item->name, fItem->name) == 0) {
    664                             found = true;
    665                             break;
    666                         }
     493    // Find the unique items in the array of metadata 'rows', and their sizes
     494    psMetadata *colSpecs = psMetadataAlloc(); // Column specifications
     495    size_t rowSize = 0;                 // Size (in bytes) of each row
     496    for (long i = 0; i < numRows; i++) {
     497        psMetadata* row = table->data[i];
     498        if (!row) {
     499            continue;
     500        }
     501        psMetadataIterator *rowIter = psMetadataIteratorAlloc(row, PS_LIST_HEAD, NULL); // Iterator
     502        psMetadataItem *colItem;        // Column item, from iteration
     503        while ((colItem = psMetadataGetAndIncrement(rowIter))) {
     504            if (!(PS_DATA_IS_PRIMITIVE(colItem->type) || colItem->type == PS_DATA_STRING ||
     505                    colItem->type == PS_DATA_VECTOR)) {
     506                // unsupported type -- treating as an error
     507                psError(PS_ERR_BAD_PARAMETER_TYPE, true,
     508                        "Unsupported data type (%d) for Metadata Item '%s' in row %d.",
     509                        colItem->type, colItem->name, row);
     510                psFree(rowIter);
     511                psFree(colSpecs);
     512                return false;
     513            }
     514
     515            size_t size = columnSize(colItem); // Size for this column
     516
     517            // Check to see if we know about this one; or update the size if required
     518            psMetadataItem *colSpecItem = psMetadataLookup(colSpecs, colItem->name);
     519            if (!colSpecItem) {
     520                // A new one!
     521                colSpec *spec = psAlloc(sizeof(colSpec)); // Specification for this column
     522                spec->type = colItem->type;
     523                spec->size = size;
     524                if (colItem->type == PS_DATA_VECTOR) {
     525                    psVector *vector = colItem->data.V; // The vector
     526                    spec->vectorType = vector->type.type;
     527                }
     528                psMetadataAddPtr(colSpecs, PS_LIST_TAIL, colItem->name, PS_DATA_UNKNOWN, "", spec);
     529                psFree(spec);           // Drop reference
     530                rowSize += PSELEMTYPE_SIZEOF(spec->type);
     531            } else {
     532                colSpec *spec = colSpecItem->data.V; // The specification
     533                if (size > spec->size) {
     534                    spec->size = size;
     535                }
     536                if (colItem->type != spec->type) {
     537                    psLogMsg(__func__, PS_LOG_WARN, "Differing type found for column %s: %x vs %x --- "
     538                             "using the first found.\n", colSpecItem->name, colItem->type, spec->type);
     539                }
     540                if (colItem->type == PS_DATA_VECTOR) {
     541                    psVector *vector = colItem->data.V; // The vector
     542                    if (vector->type.type != spec->vectorType) {
     543                        psLogMsg(__func__, PS_LOG_WARN, "Differing vector type found for column %s: %x vs %x"
     544                                 "--- using the first found.\n", colSpecItem->name, vector->type.type,
     545                                 spec->vectorType);
    667546                    }
    668                     if (! found) {
    669                         psArrayAdd(columns, columns->nalloc, item);
    670                     } else if (item->type == PS_DATA_STRING &&
    671                                strlen(fItem->data.V) < strlen(item->data.V)) {
    672                         // got to keep the longest string value as to know what size to create the table column
    673                         psMemDecrRefCounter(fItem);
    674                         columns->data[n] = psMemIncrRefCounter(item);
    675                         columns->n++;
    676 
    677                     } else if (item->type == PS_DATA_VECTOR &&
    678                                ((psVector*)(fItem->data.V))->n < ((psVector*)(item->data.V))->n) {
    679                         psMemDecrRefCounter(fItem);
    680                         columns->data[n] = psMemIncrRefCounter(item);
    681                         columns->n++;
    682                     }
    683                 } else {
    684                     // unsupported type -- treating as an error
    685                     psError(PS_ERR_BAD_PARAMETER_TYPE, true,
    686                             "Unsupported data type (%d) for Metadata Item '%s' in row %d.",
    687                             item->type, item->name, row);
    688                     psFree(iter);
    689                     psFree(columns);
    690                     return false;
    691547                }
    692548            }
    693             psFree(iter);
    694         }
    695     }
    696 
    697     if (columns->n == 0) { // no table columns found
     549        }
     550        psFree(rowIter);
     551    }
     552
     553    long numColumns = colSpecs->list->n;// Number of columns
     554    if (numColumns == 0) {
     555        // No table columns found
    698556        psError(PS_ERR_BAD_PARAMETER_SIZE, true,
    699557                "Did not find any column data to write to a table.");
    700         psFree(columns);
     558        psFree(colSpecs);
    701559        return false;
    702560    }
    703561
    704     //create list of column names and types.
    705     if (columns->n > 0) {  // If we're going decrement columns->n (why?), at least free the memory
    706         psFree(columns->data[columns->n - 1]);
    707     }
    708     columns->n--;
    709     psArray* columnNames = psArrayAlloc(columns->n);
    710     psArray* columnTypes = psArrayAlloc(columns->n);
    711     for (int n=0; n < columns->n; n++) {
    712         columnNames->data[n] = psMemIncrRefCounter(((psMetadataItem*)columns->data[n])->name);
    713         columnNames->n++;
    714         columnTypes->n++;
    715         if ( ! getMetadataTForm((psMetadataItem*)columns->data[n],
    716                                 (char**)&(columnTypes->data[n]),1) ) {
    717             psError(PS_ERR_UNKNOWN, true,
    718                     "Failed to determine the FITS data type of '%s' (type=%d).",
    719                     columnNames->data[n],
    720                     ((psMetadataItem*)columns->data[n])->type);
    721             psFree(columns);
    722             psFree(columnNames);
    723             psFree(columnTypes);
    724             return false;
    725         }
     562    // Create array of column names and types.
     563    psArray *columnNames = psArrayAlloc(numColumns); // Array of column names, for cfitsio
     564    psArray *columnTypes = psArrayAlloc(numColumns); // Array of column types, for cfitsio
     565    columnNames->n = numColumns;
     566    columnTypes->n = numColumns;
     567    psMetadataIterator *colSpecsIter = psMetadataIteratorAlloc(colSpecs, PS_LIST_HEAD, NULL); // Iterator
     568    psMetadataItem *colSpecItem;        // Column specification item, from iteration
     569    for (long i = 0; (colSpecItem = psMetadataGetAndIncrement(colSpecsIter)); i++) {
     570        colSpec *spec = colSpecItem->data.V; // The specification
     571        columnNames->data[i] = psMemIncrRefCounter(colSpecItem->name);
     572        psString colType = NULL;        // The column type
     573        psStringAppend(&colType, "%d%c", spec->size, getTForm(spec->type));
     574        columnTypes->data[i] = colType;
    726575    }
    727576
    728577    // Create the table HDU
    729     int hdus = psFitsGetSize(fits);     // Number of HDUs in file
    730     if (hdus == 0) {
     578    int numHDUs = psFitsGetSize(fits);  // Number of HDUs in file
     579    if (numHDUs == 0) {
    731580        // We're creating the first extension
    732581        fits_create_tbl(fits->fd,
    733582                        BINARY_TBL,
    734583                        table->n, // number of rows in table
    735                         columns->n, // number of columns in table
     584                        numColumns, // number of columns in table
    736585                        (char**)columnNames->data, // names of the columns
    737586                        (char**)columnTypes->data, // format of the columns
    … …  
    753602        fits_insert_btbl(fits->fd,
    754603                         table->n, // number of rows in table
    755                          columns->n, // number of columns in table
     604                         numColumns, // number of columns in table
    756605                         (char**)columnNames->data, // names of the columns
    757606                         (char**)columnTypes->data, // format of the columns
    758607                         NULL, // physical unit of columns
    759608                         (char*)extname, // extension name; casting away const because cfitsio is horrible
    760                          0,
    761                          &status);
     609                         0, &status);
    762610    }
    763611    psFree(columnNames);
    … …  
    770618    }
    771619
    772     // fill in the table elements with data
    773     for (int n = 0; n < columns->n; n++) {
    774         int row;
    775         item = columns->data[n];
    776         if (PS_DATA_IS_PRIMITIVE(item->type) ) {
    777             psVector* col = psVectorAlloc(table->n, item->type);
    778             int dataSize = PSELEMTYPE_SIZEOF(item->type);
    779             int dataType;
    780             convertPsTypeToFits(item->type, NULL,NULL,&dataType);
    781             for (row = 0; row < table->n; row++) {
    782                 psMetadataItem* dataItem = psMetadataLookup(
    783                                                table->data[row],
    784                                                item->name);
    785                 memcpy(&col->data.U8[row*dataSize],&dataItem->data,dataSize);
    786             }
     620
     621    // cfitsio requires that we write the data by columns --- urgh!
     622    psMetadataIteratorSet(colSpecsIter, PS_LIST_HEAD);
     623    for (long colNum = 1; (colSpecItem = psMetadataGetAndIncrement(colSpecsIter)); colNum++) {
     624        // Note: colNum is unit-indexed, because it's for cfitsio
     625        colSpec *spec = colSpecItem->data.V; // The specification
     626        if (PS_DATA_IS_PRIMITIVE(spec->type)) {
     627            size_t dataSize = PSELEMTYPE_SIZEOF(spec->type); // Size (in bytes) of this type
     628            psVector *columnData = psVectorAlloc(table->n, spec->type); // The raw row data, to be written
     629            psVectorInit(columnData, 0);
     630            for (long i = 0; i < table->n; i++) {
     631                psMetadata *row = table->data[i]; // The row of interest
     632                psMetadataItem *dataItem = psMetadataLookup(row, colSpecItem->name); // The value of interest
     633                memcpy(&columnData->data.U8[i * dataSize], &dataItem->data, dataSize);
     634            }
     635
     636            int fitsDataType;           // Data type for cfitsio
     637            p_psFitsTypeToCfitsio(spec->type, NULL, NULL, &fitsDataType);
    787638            fits_write_col(fits->fd,
    788                            dataType,
    789                            n+1, // column number
    790                            1, // firstrow
    791                            1, // firstelem
    792                            table->n, // nelements
    793                            col->data.U8,
     639                           fitsDataType,
     640                           colNum, // column number
     641                           1, // first row
     642                           1, // first element
     643                           table->n, // number of rows
     644                           columnData->data.U8, // the data
    794645                           &status);
    795             psFree(col);
    796         } else if (item->type == PS_DATA_STRING) {
    797             psArray* colArray = psArrayAlloc(table->n);
    798             for (row = 0; row < table->n; row++) {
    799                 colArray->data[row] = psStringCopy(psMetadataLookupStr(NULL,
    800                                                    table->data[row],
    801                                                    item->name));
    802                 colArray->n++;
    803             }
    804             fits_write_col_str(fits->fd,
    805                                n+1,
    806                                1,
    807                                1,
    808                                table->n,
    809                                (char**)colArray->data,
     646            psFree(columnData);
     647            continue;
     648        }
     649
     650        switch (spec->type) {
     651        case PS_DATA_STRING: {
     652                psArray *strings = psArrayAlloc(table->n); // Array of strings
     653                strings->n = table->n;
     654                for (long i = 0; i < table->n; i++) {
     655                    psMetadata *row = table->data[i]; // The row of interest
     656                    strings->data[i] = psMemIncrRefCounter(psMetadataLookupStr(NULL, row, colSpecItem->name));
     657                }
     658                fits_write_col_str(fits->fd, colNum, 1, 1, table->n, (char**)strings->data, &status);
     659                psFree(strings);
     660                break;
     661            }
     662        case PS_DATA_VECTOR: {
     663                psVector *columnData = psVectorAlloc(spec->size * table->n, spec->vectorType);
     664                psVectorInit(columnData, 0);
     665                size_t dataSize = PSELEMTYPE_SIZEOF(spec->vectorType); // Size of data, in bytes
     666                for (long i = 0; i < table->n; i++) {
     667                    psMetadata *row = table->data[i]; // The row of interest
     668                    psMetadataItem* dataItem = psMetadataLookup(row, colSpecItem->name);
     669                    if (dataItem->type != PS_DATA_VECTOR) {
     670                        // Just in case --- get a zero instead of some weird result
     671                        continue;
     672                    }
     673                    psVector *vector = dataItem->data.V;
     674                    memcpy(&columnData->data.U8[i * dataSize], vector->data.U8, dataSize);
     675                }
     676
     677                int fitsDataType;           // Data type for cfitsio
     678                p_psFitsTypeToCfitsio(spec->vectorType, NULL, NULL, &fitsDataType);
     679                fits_write_col(fits->fd, fitsDataType, colNum, 1, 1, table->n * spec->size, columnData->data.U8,
    810680                               &status);
    811             psFree(colArray);
    812         } else if (item->type == PS_DATA_VECTOR) {
    813             psVector* vec = (psVector*)(item->data.V);
    814             psVector* col = psVectorAlloc(table->n*vec->n, vec->type.type);
    815             int dataSize = PSELEMTYPE_SIZEOF(vec->type.type)*vec->n;
    816             int dataType;
    817             convertPsTypeToFits(vec->type.type, NULL,NULL,&dataType);
    818             for (row = 0; row < table->n; row++) {
    819                 psMetadataItem* dataItem = psMetadataLookup(
    820                                                table->data[row],
    821                                                item->name);
    822                 memcpy(&col->data.U8[row*dataSize],
    823                        ((psVector*)(dataItem->data.V))->data.U8,
    824                        dataSize);
    825             }
    826             fits_write_col(fits->fd,
    827                            dataType,
    828                            n+1, // column number
    829                            1, // firstrow
    830                            1, // firstelem
    831                            table->n*vec->n, // nelements
    832                            col->data.U8,
    833                            &status);
    834             psFree(col);
    835         }
    836         if (status != 0) {
    837             char fitsErr[MAX_STRING_LENGTH];
    838             (void)fits_get_errstatus(status, fitsErr);
    839             psError(PS_ERR_IO, true,
    840                     PS_ERRORTEXT_psFits_WRITE_FAILED,
    841                     fitsErr);
    842             psFree(columns);
    843             return false;
    844         }
    845 
    846     }
    847 
    848     psFree(columns);
     681                psFree(columnData);
     682                break;
     683            }
     684        default:
     685            psAbort(__func__, "Should never get here.\n");
     686        }
     687    }
     688
     689    psFree(colSpecsIter);
     690    psFree(colSpecs);
    849691
    850692    return true;
    … …  
    899741                // cooresponding column found in table
    900742                int dataType;
    901                 convertPsTypeToFits(item->type, NULL, NULL, &dataType);
     743                p_psFitsTypeToCfitsio(item->type, NULL, NULL, &dataType);
    902744
    903745                if (fits_write_col(fits->fd, dataType, colnum, row+1, 1, 1, &item->data,&status) != 0) {
Note: See TracChangeset for help on using the changeset viewer.