IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 15353


Ignore:
Timestamp:
Oct 23, 2007, 10:45:20 AM (19 years ago)
Author:
eugene
Message:

updates and additions to enable reading compressed fits images

Location:
branches/eam_branch_20071015/Ohana/src
Files:
1 added
6 edited

Legend:

Unmodified
Added
Removed
  • branches/eam_branch_20071015/Ohana/src/libfits/Makefile

    r12879 r15353  
    1515
    1616# programs may add their own internal requirements here
    17 FULL_CFLAGS   = $(BASE_CFLAGS) -fPIC
     17FULL_CFLAGS   = $(BASE_CFLAGS) -fPIC -Wall -Werror
    1818FULL_CPPFLAGS = $(BASE_CPPFLAGS)
    1919FULL_LDFLAGS  = $(BASE_LDFLAGS)
     
    3232$(HEAD)/F_copy_H.$(ARCH).o                      $(HEAD)/F_delete.$(ARCH).o    \
    3333$(HEAD)/F_read_XH.$(ARCH).o                     $(HEAD)/F_init_H.$(ARCH).o    \
    34 $(HEAD)/version.$(ARCH).o
     34$(HEAD)/F_convert_H.$(ARCH).o                   $(HEAD)/version.$(ARCH).o
    3535
    3636MATRIX_OBJ = \
     
    4343$(MATR)/F_convert_format.$(ARCH).o              $(MATR)/F_read_segment.$(ARCH).o \
    4444$(MATR)/F_read_portion.$(ARCH).o                $(MATR)/F_load_M.$(ARCH).o       \
    45 $(MATR)/F_matrix.$(ARCH).o
     45$(MATR)/F_matrix.$(ARCH).o                      $(MATR)/F_compress_M.$(ARCH).o
    4646
    4747TABLE_OBJ = \
     
    5151$(TABL)/F_define_column.$(ARCH).o               $(TABL)/F_table_format.$(ARCH).o \
    5252$(TABL)/F_set_column.$(ARCH).o                  $(TABL)/F_get_column.$(ARCH).o   \
    53 $(TABL)/F_table_row.$(ARCH).o                   $(TABL)/F_free_T.$(ARCH).o
     53$(TABL)/F_table_row.$(ARCH).o                   $(TABL)/F_free_T.$(ARCH).o       \
     54$(TABL)/F_table_varlength.$(ARCH).o
    5455
    5556OBJS = $(HEADER_OBJ) $(MATRIX_OBJ) $(TABLE_OBJ)
  • branches/eam_branch_20071015/Ohana/src/libfits/include/gfitsio.h

    r15312 r15353  
    11/* FITS specific macros and structures */
     2
     3# include <assert.h>
    24
    35# ifndef GFITSIO
     
    9193  int  maxlen;
    9294  int  nbytes;
     95  int  Nstart;
     96  int  theap;
    9397  char format;
     98  char cmptype[81];
    9499} VarLengthColumn;
    95100
     
    124129int     gfits_vscan                    PROTO((Header *header, char *field, char *mode, int N, va_list argp));
    125130int     gfits_write_header             PROTO((char *filename, Header *header));
    126 int     gfits_matrix_size              PROTO((Header *header));
     131int     gfits_data_size                PROTO((Header *header));
     132int     gfits_extended_to_primary      PROTO((Header *header, int simple, char *comment));
     133int     gfits_primary_to_extended      PROTO((Header *header, char *exttype, char *comment));
     134int     gfits_modify_extended          PROTO((Header *header, char *exttype, char *comment));
     135
    127136
    128137/******************************* Matrix functions *************/
     
    146155void    gfits_set_matrix_value         PROTO((Matrix *matrix, int x, int y, double value));
    147156int     gfits_write_matrix             PROTO((char *filename, Matrix *matrix));
     157Matrix *gfits_uncompress_image         PROTO((Header *header, FTable *ftable, int primary));
     158int     gfits_uncompress_data          PROTO((char *zdata, int Nzdata, char *cmptype, char **optname, char **optvalue, int Nopt, char *outdata, int *Nout));
     159int     gfits_distribute_data          PROTO((Matrix *matrix, int bitpix, char *data, int Ndata, int *otile, int *ztile, float zscale, float zzero));
     160int     gfits_byteswap_zdata           PROTO((char *zdata, int Nzdata, int bitpix));
    148161
    149162/******************************* Table functions *************/
     
    188201int     gfits_write_Theader            PROTO((char *filename, Header *header));
    189202
     203int     gfits_varlength_column_define  PROTO((FTable *ftable, VarLengthColumn *def, int column));
     204void   *gfits_varlength_column_pointer PROTO((FTable *ftable, VarLengthColumn *column, int row, int *length));
     205
    190206#endif /* FITSIO */
  • branches/eam_branch_20071015/Ohana/src/libfits/matrix/F_compress_M.c

    r15312 r15353  
    11# include <ohana.h>
    22# include <gfitsio.h>
     3# include <zlib.h>
    34
    45# define MOD_KEYWORD(NAME,ZNAME,TYPE,IN,OUT) { \
    56  gfits_scan (header, ZNAME, TYPE, 1, IN); \
    67  gfits_delete_field (header, ZNAME, 1); \
    7   gfits_print (header, NAME, TYPE, 1, OUT);
     8  gfits_print (header, NAME, TYPE, 1, OUT); }
    89
    910# define MOD_KEYWORD_REQUIRED(NAME,ZNAME,TYPE,IN,OUT) { \
    10   if (!gfits_scan (header, ZNAME, TYPE, 1, IN)) return (NULL); \
     11  if (!gfits_scan (header, ZNAME, TYPE, 1, IN)) ESCAPE; \
    1112  gfits_delete_field (header, ZNAME, 1); \
    12   gfits_print (header, NAME, TYPE, 1, OUT);
    13 
    14 Matrix *gfits_uncompress_image (Header *header, FTable *ftable) {
    15 
    16   int status;
    17   int zimage;
     13  gfits_print (header, NAME, TYPE, 1, OUT); }
     14
     15# define ESCAPE { \
     16  if (ztile != NULL) free (ztile); \
     17  if (optname != NULL) { \
     18    for (j = 0; j < Noptions; j++) { \
     19      free (optname[j]); \
     20      free (optvalue[j]); \
     21    } \
     22    free (optname); \
     23    free (optvalue); \
     24  } \
     25  if (out != NULL) free (out); \
     26  if (otile != NULL) free (otile); \
     27  if (ntile != NULL) free (ntile); \
     28  return (NULL); }
     29
     30int gfits_uncompress_image (Header *header, Matrix *matrix, FTable *ftable, int primary) {
     31
     32  int i, j, status, zimage, zcol, Ncur, Nzrows, Nout, Nzdata;
    1833  char cmptype[80];
    19   char azis[10];
     34  char zaxis[10], naxis[10], key[10], word[81], exttype[81], checksum[81], datasum[81];
     35  int Noptions, NOPTIONS;
     36  VarLengthColumn zdef;
     37  float zscale, zzero;
     38
     39  int *ztile = NULL;
     40  int *otile = NULL;
     41  int *ntile = NULL;
     42  char **optname = NULL;
     43  char **optvalue = NULL;
     44  char *out = NULL;
    2045
    2146  // is ZIMAGE present?
    2247  status = gfits_scan (ftable->header, "ZIMAGE", "%t", 1, &zimage);
    23   if (!status || !zimage) return NULL;
     48  if (!status || !zimage) ESCAPE;
    2449
    2550  // copy original header to output header
    26   gfits_copy_header (&ftable->header, header);
     51  gfits_copy_header (ftable->header, header);
    2752
    2853  // extract compression-specific keywords, update header as needed.
    29   if (!gfits_scan (header, "ZCMPTYPE", "%s", 1, cmptype)) return NULL;
     54  if (!gfits_scan (header, "ZCMPTYPE", "%s", 1, cmptype)) ESCAPE;
    3055  gfits_delete_field (header, "ZCMPTYPE", 1);
    3156
     
    3560  for (i = 0; i < header->Naxes; i++) {
    3661    snprintf (zaxis, 10, "ZNAXIS%d", i + 1);
    37     snprintf (axis, 10, "NAXIS%d", i + 1);
    38     MOD_KEYWORD_REQUIRED (zaxis,  axis,  "%d", &header->Naxis[i],  header->Naxis[i]);
     62    snprintf (naxis, 10, "NAXIS%d", i + 1);
     63    MOD_KEYWORD_REQUIRED (zaxis,  naxis,  "%d", &header->Naxis[i],  header->Naxis[i]);
    3964  }   
    4065
     
    4873  } else {
    4974    for (i = 1; i < header->Naxes; i++) {
    50       snprintf (axis, 10, "ZTILE%d", i + 1);
    51       if (!gfits_scan (header, axis, "%d", 1, &ztile[i])) return NULL;
    52       gfits_delete_field (header, axis, 1);
     75      snprintf (key, 10, "ZTILE%d", i + 1);
     76      if (!gfits_scan (header, key, "%d", 1, &ztile[i])) ESCAPE;
     77      gfits_delete_field (header, key, 1);
    5378    }
    5479  }
     
    5984  ALLOCATE (optvalue, char *, NOPTIONS);
    6085  for (Noptions = 0; TRUE; Noptions++) {
    61     snprintf (axis, 10, "ZNAME%d", Noptions + 1);
    62     if (!gfits_scan (header, axis, "%s", 1, word)) break;
    63     gfits_delete_field (header, axis, 1);
     86    snprintf (key, 10, "ZNAME%d", Noptions + 1);
     87    if (!gfits_scan (header, key, "%s", 1, word)) break;
     88    gfits_delete_field (header, key, 1);
    6489    optname[Noptions] = strcreate (word);
    6590
    66     snprintf (axis, 10, "ZVAL%d", Noptions + 1);
    67     if (!gfits_scan (header, axis, "%s", 1, word)) ESCAPE;
    68     gfits_delete_field (header, axis, 1);
     91    snprintf (key, 10, "ZVAL%d", Noptions + 1);
     92    if (!gfits_scan (header, key, "%s", 1, word)) ESCAPE;
     93    gfits_delete_field (header, key, 1);
    6994    optvalue[Noptions] = strcreate (word);
    7095
     
    75100    }
    76101  }
     102
    77103  // XXX is Noptions set to the correct value?
    78 
    79104  // XXX get ZMASKCMP
     105  // XXX can primary be defined on basis of the header?
    80106
    81107  if (primary) {
    82108    // XXX SIMPLE & XTENSION need to be set on the first line
    83     gfits_scan (header, "ZSIMPLE", "%t", 1, &header->simple) return NULL;
     109    if (!gfits_scan (header, "ZSIMPLE", "%t", 1, &header->simple)) ESCAPE;
    84110    gfits_delete_field (header, "ZSIMPLE", 1);
    85     gfits_print (header, "SIMPLE", "%d", 1, header->simple);
     111    gfits_extended_to_primary (header, header->simple, "Image data");
    86112
    87113    MOD_KEYWORD ("ZEXTEND",  "EXTEND",   "%t", &header->extend, header->extend);
     
    89115  } else {
    90116    // XXX XTENSION needs to be set on the first line
    91     gfits_scan (header, "ZTENSION", "%s", 1, exttype) return NULL;
     117    if (!gfits_scan (header, "ZTENSION", "%s", 1, exttype)) ESCAPE;
    92118    gfits_delete_field (header, "ZTENSION", 1);
    93     gfits_print (header, "SIMPLE", "%d", 1, header->simple);
     119    gfits_modify_extended (header, exttype, "Image extension");
    94120
    95121    MOD_KEYWORD ("ZPCOUNT",  "PCOUNT",   "%d", &header->pcount, header->pcount);
     
    99125  MOD_KEYWORD ("ZHECKSUM", "ZHECKSUM", "%s", checksum,        checksum);
    100126  MOD_KEYWORD ("ZDATASUM", "DATASUM",  "%s", datasum,         datasum);
     127
     128  zscale = 1;
     129  gfits_scan (header, "ZSCALE", "%lf", 1, &zscale);
     130
     131  zzero = 0;
     132  gfits_scan (header, "ZZERO", "%lf", 1, &zzero);
    101133
    102134  // find the COMPRESSED_DATA column (format should be 1PB, 1PI, 1PJ)
    103135  for (i = 1; TRUE; i++) {
    104136    snprintf (key, 10, "TTYPE%d", i);
    105     if (!gfits_scan (&ftable->header, key, "%s", 1, ttype)) return NULL;
    106     if (!strcmp (ttype, "COMPRESSED_DATA")) break;
    107   }
    108   compcol = i;
    109 
    110   if (!gfits_varlength_column_define (ftable, compdef, compcol)) return NULL;
     137    if (!gfits_scan (ftable->header, key, "%s", 1, word)) ESCAPE;
     138    if (!strcmp (word, "COMPRESSED_DATA")) break;
     139  }
     140  zcol = i;
     141
     142  if (!gfits_varlength_column_define (ftable, &zdef, zcol)) ESCAPE;
    111143
    112144  // create the output image
    113   gfits_create_matrix (matrix, header);
     145  gfits_create_matrix (header, matrix);
     146
     147  // XXX is maxlen in bytes or in pixels?
     148  ALLOCATE (out, char, zdef.maxlen);
     149
     150  // counters for tile numbers
     151  ALLOCATE (otile, int, matrix->Naxes);
     152  ALLOCATE (ntile, int, matrix->Naxes);
     153  for (i = 0; i < matrix->Naxes; i++) {
     154    otile[i] = 0;
     155    ntile[i] = (matrix->Naxis[i] % ztile[i]) ? (matrix->Naxis[i] / ztile[i] + 1) : (matrix->Naxis[i] / ztile[i]);
     156  }
    114157
    115158  // uncompress the data
    116 
    117   // ZSCALE & ZZERO convert uncompressed values to output values
    118   // data = (uncompressed * ZSCALE) + ZZERO
    119 
     159  Nzrows = ftable->header->Naxis[1];
     160  for (i = 0; i < Nzrows; i++) {
     161    switch (zdef.format) {
     162
     163      // may need to multiple by byte/pix
     164      Nout = zdef.maxlen;
     165
     166      // gfits_varlength returns a pointer to a chunk of Nzdata of data elements starting at zdata.
     167      // uncompress_data uncompresses from zdata to the temporary output buffer which must be allocated
     168      // XXX need to have an API to send the user data
     169
     170      case 'B': {
     171        char *zdata;
     172        zdata = gfits_varlength_column_pointer (ftable, &zdef, i, &Nzdata);
     173        gfits_byteswap_zdata (zdata, Nzdata, 8);
     174        gfits_uncompress_data (zdata, Nzdata, zdef.cmptype, optname, optvalue, Noptions, out, &Nout);
     175        gfits_distribute_data (matrix, 8, out, Nout, otile, ztile, zscale, zzero);      // copy the uncompressed pixels into their correct locations
     176        break;
     177      }
     178
     179      case 'I': {
     180        short *zdata;
     181        zdata = gfits_varlength_column_pointer (ftable, &zdef, i, &Nzdata);
     182        gfits_byteswap_zdata ((char *) zdata, Nzdata, 16);
     183        gfits_uncompress_data ((char *) zdata, Nzdata, zdef.cmptype, optname, optvalue, Noptions, out, &Nout);
     184        gfits_distribute_data (matrix, 16, out, Nout, otile, ztile, zscale, zzero);     // copy the uncompressed pixels into their correct locations
     185        break;
     186      }
     187
     188      case 'J': {
     189        int *zdata;
     190        zdata = gfits_varlength_column_pointer (ftable, &zdef, i, &Nzdata);
     191        gfits_byteswap_zdata ((char *) zdata, Nzdata, 32);
     192        gfits_uncompress_data ((char *) zdata, Nzdata, zdef.cmptype, optname, optvalue, Noptions, out, &Nout);
     193        gfits_distribute_data (matrix, 32, out, Nout, otile, ztile, zscale, zzero);     // copy the uncompressed pixels into their correct locations
     194        break;
     195      }
     196
     197      default:
     198        abort();
     199    }
     200
     201    // update the tile counters, carrying to the next dimension if needed
     202    for (j = 0; j < matrix->Naxes; j++) {
     203      otile[j] ++;
     204      if (otile[j] == ntile[j]) {
     205        otile[j] = 0;
     206      } else {
     207        break;
     208      }
     209    }
     210  }
     211  return (TRUE);
     212}
     213
     214int gfits_uncompress_data (char *zdata, int Nzdata, char *cmptype, char **optname, char **optvalue, int Nopt, char *outdata, int *Nout) {
     215
     216  // XXX allow lowercase value?
     217  if (!strcasecmp(cmptype, "GZIP")) {
     218    unsigned long tNout = *Nout;
     219    uncompress ((Bytef *) outdata, &tNout, (Bytef *) zdata, Nzdata);
     220    *Nout = tNout;
     221    return (TRUE);
     222  }
     223
     224  if (!strcasecmp(cmptype, "RICE")) {
     225  }
    120226 
     227  if (!strcasecmp(cmptype, "PLIO")) {
     228  }
     229
     230  if (!strcasecmp(cmptype, "HCOMPRESS")) {
     231  }
     232  return (TRUE);
     233}
     234
     235// bitpix is the input data size/type
     236int gfits_distribute_data (Matrix *matrix, int bitpix, char *data, int Ndata, int *otile, int *ztile, float zscale, float zzero) {
     237
     238  int i, j, start, offset, coord, Nline;
     239  int *counter = NULL;
     240  int *Ztile = NULL;
     241
     242  // counter for current row in tile to copy
     243  ALLOCATE (counter, int, matrix->Naxes);
     244  memset (counter, 0, matrix->Naxes*sizeof(int));
     245
     246  // true sizes of this tile (in pixels)
     247  ALLOCATE (Ztile, int, matrix->Naxes);
     248  for (i = 0; i < matrix->Naxes; i++) {
     249    Ztile[i] = MIN ((matrix->Naxis[i] - otile[i]*ztile[i]), ztile[i]);
     250  }
     251
     252  // number of lines in the tile (in pixels)
     253  Nline = 1;
     254  for (i = 1; i < matrix->Naxes; i++) {
     255    Nline *= Ztile[i];
     256  }
     257
     258  // double check reported size (needs the Isize)
     259  assert (Ndata == Ztile[0]*Nline);
     260
     261  // set the starting point of the tile:
     262  // start = otile[0]*ztile[0] + otile[1]*ztile[1]*Naxis[0] + otile[2]*ztile[2]*Naxis[0]*Naxis[1] + ...;
     263  // start = otile[0]*ztile[0] + Naxis[0]*(otile[1]*ztile[1] + Naxis[1]*(otile[2]*ztile[2] + ...));
     264  start = otile[matrix->Naxes-1]*ztile[matrix->Naxes-1];
     265  for (i = matrix->Naxes - 2; i >= 0; i--) {
     266    coord = otile[i]*ztile[i];
     267    start = start*matrix->Naxis[i] + coord;
     268  }
     269 
     270  // pixel offset in output array relative to tile start
     271  offset = 0;
     272
     273  // we need to set up switches for all of the possible combinations:
     274  // this macro is used at the inner switch to run the actual loop
     275# define SCALE_AND_DIST(TYPE, SIZE) { \
     276    TYPE *Optr = (TYPE *) &matrix->buffer[SIZE*(offset + start)]; \
     277    for (j = 0; j < Ztile[0]; j++, Iptr++, Optr++) { \
     278      *Optr = *Iptr * zscale + zzero; \
     279    } }
     280
     281  // this macro sets up the outer switch and calls above macro with all output bitpix options
     282# define SETUP_INSIZE(TYPE, SIZE) { \
     283    TYPE *Iptr = (TYPE *) &data[i*SIZE*Ztile[0]]; \
     284    switch (matrix->bitpix) { \
     285      case 8: \
     286        SCALE_AND_DIST (char, 1); \
     287        break; \
     288      case 16: \
     289        SCALE_AND_DIST (short, 2); \
     290        break; \
     291      case 32: \
     292        SCALE_AND_DIST (int, 4); \
     293        break; \
     294      case -32: \
     295        SCALE_AND_DIST (float, 4); \
     296        break; \
     297      case -64: \
     298        SCALE_AND_DIST (double, 8); \
     299        break; \
     300      default: \
     301        abort(); \
     302    } }
     303
     304  // loop over lines in the tile
     305  for (i = 0; i < Nline; i++) {
     306    switch (bitpix) {
     307      case 8:
     308        SETUP_INSIZE (char, 1);
     309        break;
     310      case 16:
     311        SETUP_INSIZE (short, 2);
     312        break;
     313      case 32:
     314        SETUP_INSIZE (int, 4);
     315        break;
     316      default:
     317        abort();
     318    }
     319
     320    // update the counters, carrying to the next dimension if needed
     321    for (j = 1; j < matrix->Naxes; j++) {
     322      counter[j] ++;
     323      if (counter[j] == Ztile[j]) {
     324        counter[j] = 0;
     325      } else {
     326        break;
     327      }
     328    }
     329    if (j == matrix->Naxes) assert (i == Nline - 1); // we should be done here...
     330
     331    // Naxes = 3
     332    // offset = counter[1]*matrix->Naxis[0] + counter[2]*matrix->Naxis[0]*matrix->Naxis[1] +
     333    // offset = matrix->Naxis[0]*(counter[1] + matrix->Naxis[1]*(counter[2] + matrix->Naxis[2]*...))
     334
     335    // determine the offset of the next line relative to the start position
     336    offset = counter[matrix->Naxes - 1];
     337    for (j = 1; j < matrix->Naxes; j++) {
     338      offset = offset*matrix->Naxis[j] + counter[j];
     339    }     
     340  }
     341
     342  free (counter);
     343  free (Ztile);
     344  return (TRUE);
     345}
     346
     347int gfits_byteswap_zdata (char *zdata, int Nzdata, int bitpix) {
     348
     349# ifdef BYTE_SWAP
     350
     351  int i;
     352  char tmp;
     353
     354  switch (bitpix) {
     355    case +8:
     356    case -8:
     357      break;
     358
     359    case +16:
     360    case -16:
     361      for (i = 0; i < Nzdata; i+=2) {
     362        tmp = zdata[i];
     363        zdata[i] = zdata[i+1];
     364        zdata[i+1] = tmp;
     365      }
     366      break;
     367
     368    case +32:
     369    case -32:
     370      for (i = 0; i < Nzdata; i+=4) {
     371        tmp = zdata[i+1];
     372        zdata[i+1] = zdata[i+2];
     373        zdata[i+2] = tmp;
     374        tmp = zdata[i];
     375        zdata[i] = zdata[i+3];
     376        zdata[i+3] = tmp;
     377      }
     378      break;
     379
     380    case +64:
     381    case -64:
     382      for (i = 0; i < Nzdata; i+=8) {
     383        tmp = zdata[i+0];
     384        zdata[i+0] = zdata[i+7];
     385        zdata[i+7] = tmp;
     386        tmp = zdata[i+1];
     387        zdata[i+1] = zdata[i+6];
     388        zdata[i+6] = tmp;
     389        tmp = zdata[i+2];
     390        zdata[i+2] = zdata[i+5];
     391        zdata[i+5] = tmp;
     392        tmp = zdata[i+3];
     393        zdata[i+3] = zdata[i+4];
     394        zdata[i+4] = tmp;
     395      }
     396      break;
     397  }
     398# endif
     399  return (TRUE);
     400}
  • branches/eam_branch_20071015/Ohana/src/libfits/table/F_table_format.c

    r15312 r15353  
    194194  int i, j, n, Nx, Ny, Nfields;
    195195  int off, Nchar, Nval, Nbytes, status;
    196   char *line, format[64], field[16], type[16];
     196  char format[64], field[16], type[16];
    197197  double tzero, tscale;
    198198  char *tmpChar;
     
    233233      for (j = 0; j < Ny; j++) {
    234234        for (n = 0; n < Nval; n++) {
    235           tmpChar = (int *)&ftable[0].buffer[j*Nx + n*Nbytes + off];
     235          tmpChar = (char *)&ftable[0].buffer[j*Nx + n*Nbytes + off];
    236236          *tmpChar = *tmpChar * tscale + tzero;
    237237        }
     
    241241      for (j = 0; j < Ny; j++) {
    242242        for (n = 0; n < Nval; n++) {
    243           tmpShort = (int *)&ftable[0].buffer[j*Nx + n*Nbytes + off];
     243          tmpShort = (short *)&ftable[0].buffer[j*Nx + n*Nbytes + off];
    244244          *tmpShort = *tmpShort * tscale + tzero;
    245245        }
     
    264264  int i, j, n, Nx, Ny, Nfields;
    265265  int off, Nchar, Nval, Nbytes, status;
    266   char *line, format[64], field[16], type[16];
     266  char format[64], field[16], type[16];
    267267  double tzero, tscale;
    268268  char *tmpChar;
     
    303303      for (j = 0; j < Ny; j++) {
    304304        for (n = 0; n < Nval; n++) {
    305           tmpChar = (int *)&ftable[0].buffer[j*Nx + n*Nbytes + off];
     305          tmpChar = (char *)&ftable[0].buffer[j*Nx + n*Nbytes + off];
    306306          *tmpChar = (*tmpChar - tzero) / tscale;
    307307        }
     
    311311      for (j = 0; j < Ny; j++) {
    312312        for (n = 0; n < Nval; n++) {
    313           tmpShort = (int *)&ftable[0].buffer[j*Nx + n*Nbytes + off];
     313          tmpShort = (short *)&ftable[0].buffer[j*Nx + n*Nbytes + off];
    314314          *tmpShort = (*tmpShort - tzero) / tscale;
    315315        }
  • branches/eam_branch_20071015/Ohana/src/libfits/table/F_table_varlength.c

    r15312 r15353  
    44int gfits_varlength_column_define (FTable *ftable, VarLengthColumn *def, int column) {
    55
    6   char field[80];
    7   char format[80];
     6  int i, N, Nv, Nb;
     7  char *p1, *p2, *p3;
     8  char field[81];
     9  char format[81];
     10  char tmpline[81];
    811
    912  // grab the value of TFORMn for this column
    1013  snprintf (field, 80, "TFORM%d", column);
    11   if (!gfits_scan (&ftable[0].header, field, "%s", 1, format)) return (FALSE);
     14  if (!gfits_scan (ftable->header, field, "%s", 1, format)) return (FALSE);
    1215
    1316  // find and remove the max field length element
     
    4952  if (!def[0].nbytes) { return (FALSE); }
    5053 
     54  /* scan columns to find column offset */
     55  def[0].Nstart = 0;
     56  for (i = 1; i < N; i++) {
     57    sprintf (field, "TFORM%d", i);
     58    gfits_scan (ftable->header, field, "%s", 1, format);
     59    gfits_bintable_format (format, tmpline, &Nv, &Nb);
     60    def[0].Nstart += Nv*Nb;
     61  }
     62
     63  if (!gfits_scan (ftable->header, "THEAP", "%d", 1, &def[0].theap)) {
     64    def[0].theap = ftable->header->Naxis[0]*ftable->header->Naxis[1];
     65  }
     66
    5167  return TRUE;
    5268}
    5369
    54 void *gfits_varlength_column_pointer (FTable *ftable, VarLengthColumn *column, int row) {
     70void *gfits_varlength_column_pointer (FTable *ftable, VarLengthColumn *column, int row, int *length) {
     71
     72  int offset, Ny, *ptr;
    5573
    5674  // find the values for the specified row
     75  // the values in the main table for this row and varlength column:
     76  Ny = ftable->header->Naxis[1];
     77  ptr = (int *) &ftable->buffer[row*Ny + column->Nstart];
     78  offset = ptr[0];
     79  *length = ptr[1];
    5780
    58  
    59 
     81  switch (column->format) {
     82    case 'X':
     83    case 'L':
     84    case 'A':
     85    case 'B': {
     86      char *result;
     87      result = (char *) ftable->buffer + column->theap + offset;
     88      return result;
     89    }
     90    case 'I': {
     91      short *result;
     92      result = (short *) ftable->buffer + column->theap + offset;
     93      return result;
     94    }
     95    case 'J': {
     96      int *result;
     97      result = (int *) ftable->buffer + column->theap + offset;
     98      return result;
     99    }
     100    case 'E':
     101    case 'C': {
     102      float *result;
     103      result = (float *) ftable->buffer + column->theap + offset;
     104      return result;
     105    }
     106    case 'D':
     107    case 'M': {
     108      double *result;
     109      result = (double *) ftable->buffer + column->theap + offset;
     110      return result;
     111    }
     112    case 'P':
     113    default:
     114      abort();
     115  }
     116  abort();
    60117}
    61118
  • branches/eam_branch_20071015/Ohana/src/opihi/cmd.data/rd.c

    r13479 r15353  
    8686    }
    8787  }
     88
    8889  /* FITS extension by name */
    8990  if (ccdsel) {
     
    9293    CCDKeyword = get_variable ("CCDKEYWORD");
    9394    if (CCDKeyword == (char *) NULL) {
    94       gprint (GP_ERR, "CCDKEYWORD variable is not set; ");
    95       gprint (GP_ERR, "using EXTNAME as default\n");
     95      // gprint (GP_ERR, "CCDKEYWORD variable is not set; ");
     96      // gprint (GP_ERR, "using EXTNAME as default\n");
    9697      CCDKeyword = strcreate ("EXTNAME");
    9798    }
     
    106107        return (FALSE);
    107108      }
    108       if (!gfits_scan (&buf[0].header, CCDKeyword, "%s", 1, ID)) {
     109
     110      // for compressed data tables, EXTNAME may be duplicated, with the first one containing the
     111      // word 'COMPRESSED_IMAGE'.  in this case, check the second EXTNAME, if CCDKeyword == EXTNAME
     112      // this may have to be a more obscure test specifically for 'imcopy' data...
     113      // need to check each header, since file may contain a mix
     114     
     115      Nword = 1;
     116      IsCompressed = FALSE;
     117      if (gfits_extension_is_compressed_data (&buf[0].header)) {
     118        if (!strcmp (CCDKeyword, "EXTNAME")) Nword = 2;
     119        IsCompressed = TRUE;
     120      }
     121      if (!gfits_scan (&buf[0].header, CCDKeyword, "%s", Nword, ID)) {
    109122        gprint (GP_ERR, "%s not in header\n", CCDKeyword);
    110123        DeleteBuffer (buf);
     
    123136  /* fix up header, if needed */
    124137  if (extend || ccdsel) {
    125     static char simple[] = "SIMPLE  =                    T / Standard FITS";
    126     int Ns, No;
    127 
    128     Ns = strlen (simple);
    129     No = 80 - Ns;
    130     strncpy (buf[0].header.buffer, simple, Ns);
    131     memset (&buf[0].header.buffer[Ns], ' ', No);
     138    if (!IsCompressed) {
     139      gfits_extended_to_primary (buf[0].header, TRUE, "Standard FITS");
     140    }
    132141  } else {
    133142      gfits_fread_header (f, &buf[0].header);
     
    135144
    136145  /* for JustHead, we skip reading the data segment */
     146  // XXX for compressed data, we need to convert the header to the equivalent uncompressed version
    137147  if (JustHead) {
    138     /* should be in CreateMatrix / CreateBuffer function */
     148    // XXX what are we doing here exactly?
    139149    buf[0].header.Naxes = 0;
    140150    ALLOCATE (buf[0].matrix.buffer, char, 1);
     
    164174
    165175  /* load matrix data */
    166   sprintf (region, "-1 -1 -1 -1 %d %d", (plane - 1), plane);
    167   status = gfits_fread_matrix_segment (f, &buf[0].matrix, &buf[0].header, region);
    168   fclose (f);
     176  if (IsCompressed) {
     177    FTable ftable;
     178    Header theader;
     179    ftable.header = &theader;
     180    gfits_copy_header (&buf[0].header, ftable.header);
     181    status = gfits_fread_table (f, &ftable);  // XXX does this do more than read the bytes?
     182    status = gfits_uncompress_image (&buf[0].header, &buf[0].matrix, ftable, ccdsel || extend);
     183    gfits_free_table (&ftable);
     184    // XXX this currently does not work for a cube (we get a cube back, not a specific plane)
     185  } else {
     186    sprintf (region, "-1 -1 -1 -1 %d %d", (plane - 1), plane);
     187    status = gfits_fread_matrix_segment (f, &buf[0].matrix, &buf[0].header, region);
     188    fclose (f);
     189  }
    169190
    170191  if (!status) {
Note: See TracChangeset for help on using the changeset viewer.