IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 38356


Ignore:
Timestamp:
Jun 3, 2015, 7:00:55 AM (11 years ago)
Author:
eugene
Message:

get the input / output / compressed pixels sizes right and get the byteswaps right

Location:
branches/eam_branches/ohana.20150429/src/libfits
Files:
1 added
11 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ohana.20150429/src/libfits/Makefile

    r38335 r38356  
    2626TEST_LDFLAGS  = $(BASE_LDFLAGS) -lFITS -lohana -ltap_ohana
    2727
    28 TESTPROG = imagecomp tablecomp zlib
     28TESTPROG = imagecomp tablecomp zlib ricetest
    2929$(TESTPROG) : % : $(TESTBIN)/%
    3030test: $(TESTPROG)

  • branches/eam_branches/ohana.20150429/src/libfits/include/gfitsio.h

    r38355 r38356  
    212212int     gfits_uncompressed_data_bitpix  PROTO((char *cmptype, int out_bitpix, char **optname, char **optvalue, int Noptions));
    213213int     gfits_vartable_heap_pixsize     PROTO((char format));
     214int     gfits_cmptype_valid             PROTO((char *cmptype));
    214215
    215216int     gfits_varlength_column_add_data PROTO((FTable *table, char *data, off_t Ndata, int row, VarLengthColumn *column));

  • branches/eam_branches/ohana.20150429/src/libfits/matrix/F_compress_M.c

    r38355 r38356  
    22# include <gfitsio.h>
    33# include <zlib.h>
     4# define VERBOSE_DUMP 0
    45
    56// the user needs to specify the various compression options:
    … …  
    3334
    3435  Header *theader = ftable->header;
     36
     37  if (!gfits_cmptype_valid(zcmptype)) ESCAPE;
    3538
    3639  // determine the number of tiles (from ztile[] and image size)
    … …  
    147150  int max_tile_size = gfits_imtile_maxsize (matrix, ztile);
    148151
    149   // size of a pixel in the image tile -- this is probably not needed
    150   int tile_pixsize = abs(header[0].bitpix) / 8;
    151 
    152152  // size of a pixel in the raw pixel buffer (before compression)
    153153  int raw_pixsize = gfits_uncompressed_data_pixsize (zcmptype, header[0].bitpix, NULL, NULL, 0);
    … …  
    155155  int cmp_pixsize = gfits_compressed_data_pixsize (zcmptype, header[0].bitpix, NULL, NULL, 0);
    156156
    157   fprintf (stderr, "raw_pixsize: %d, cmp_pixsize: %d, tile_pixsize: %d, raw_bitpix: %d\n",
    158            raw_pixsize, cmp_pixsize, tile_pixsize, raw_bitpix);
     157  // size of a pixel in the image tile -- this is probably not needed
     158  // int tile_pixsize = abs(header[0].bitpix) / 8;
     159  // fprintf (stderr, "raw_pixsize: %d, cmp_pixsize: %d, tile_pixsize: %d, raw_bitpix: %d\n", raw_pixsize, cmp_pixsize, tile_pixsize, raw_bitpix);
    159160
    160161  // allocate the buffer for compression work
    … …  
    177178    if (!gfits_collect_data (matrix, raw, Nraw, raw_bitpix, otile, oblank, ztile, zblank, zscale, zzero)) ESCAPE;
    178179
    179     if (0 && (i == 0)) {
     180    if (VERBOSE_DUMP && (i == 0)) {
    180181      int k;
    181182      fprintf (stderr, "cmp mat: ");
    … …  
    189190    // gzip compresses bytes which are in BIG-ENDIAN order
    190191    if (!strcasecmp(zcmptype, "GZIP_1")) {
    191       if (!gfits_byteswap_zdata (raw, Nraw, raw_pixsize)) ESCAPE;
    192     }
    193     if (0 && (i == 0)) {
     192      // Nraw is number of pixels
     193      if (!gfits_byteswap_zdata (raw, Nraw * raw_pixsize, raw_pixsize)) ESCAPE;
     194    }
     195    int Nzdata = Nzdata_alloc; // available space, replaced with actual output size on compression
     196
     197    if (VERBOSE_DUMP && (i == 0)) {
    194198      int k;
    195199      fprintf (stderr, "cmp swp: ");
    196200      for (k = 0; k < 32; k++) { fprintf (stderr, "0x%02hhx ", raw[k]); }
    197201      fprintf (stderr, "\n");
     202      fprintf (stderr, "Nzdata: %d -> ", Nzdata);
    198203    }
    199204
    200205    // optname, optvalue = NULL, Noptions = 0 : defaults for RICE_1, HCOMPRESS_1
    201206
    202     int Nzdata = Nzdata_alloc; // available space, replaced with actual output size on compression
    203207    if (!gfits_compress_data (zdata, &Nzdata, zcmptype, NULL, NULL, 0, raw, Nraw, raw_pixsize, ztile[1], ztile[0])) ESCAPE;
    204     if (0 && (i == 0)) {
     208    if (VERBOSE_DUMP && (i == 0)) {
     209      int k;
     210      fprintf (stderr, "%d\n", Nzdata);
     211      fprintf (stderr, "cmp SWP: ");
     212      for (k = 0; k < 32; k++) { fprintf (stderr, "0x%02hhx ", zdata[k]); }
     213      fprintf (stderr, "\n");
     214    }
     215
     216    // all other compression modes require swapping after compression
     217    // (compresssion & decompression operate on native ENDIAN)
     218    if (strcasecmp(zcmptype, "GZIP_1")) {
     219      // Nzdata is number of bytes
     220      if (!gfits_byteswap_zdata (zdata, Nzdata, cmp_pixsize)) ESCAPE;
     221    }
     222
     223    if (VERBOSE_DUMP && (i == 0)) {
    205224      int k;
    206225      fprintf (stderr, "cmp dat: ");
    … …  
    208227      fprintf (stderr, "\n");
    209228    }
    210 
    211229    if (!gfits_varlength_column_add_data (ftable, zdata, Nzdata, i, &zdef)) ESCAPE;
    212 
    213     // all other compression modes require swapping after compression
    214     // (compresssion & decompression operate on native ENDIAN)
    215     if (strcasecmp(zcmptype, "GZIP_1")) {
    216       if (!gfits_byteswap_zdata (zdata, Nzdata, cmp_pixsize)) ESCAPE;
    217     }
    218230
    219231    // update the otile[] counters, carrying to the next dimension if needed
    … …  
    222234
    223235  if (!gfits_varlength_column_finish (ftable, &zdef)) ESCAPE;
    224   if (0) {
     236  if (VERBOSE_DUMP) {
    225237    int k;
    226238    fprintf (stderr, "cmp tbl: ");

  • branches/eam_branches/ohana.20150429/src/libfits/matrix/F_compress_data.c

    r38355 r38356  
    6565  }
    6666
    67   if (!strcasecmp(cmptype, "RICE_1")) {
     67  if (!strcasecmp(cmptype, "RICE_1") || !strcasecmp(cmptype, "RICE_ONE")) {
    6868    int i, blocksize;
    6969    // look for the BLOCKSIZE
    … …  
    8383    // fprintf (stderr, "%d comp bytes; %d uncomp 'pixels', totals: %d %d\n", Nzdata, Npix, Ninsum, Noutsum);
    8484
     85    char buffer[1000];
     86    int k;
     87
    8588    switch (rawpix_size) {
    8689      case 4:
    … …  
    9093
    9194      case 2:
     95        // Nout is the number of bytes in the compressed buffer
    9296        Nout = fits_rcomp_short ((short *) rawdata, Nrawpix, (unsigned char *) zdata, *Nzdata, blocksize);
     97
     98        if (0) {
     99          fprintf (stderr, "Nout: %d, Nrawpix: %d\n", Nout, Nrawpix);
     100          fprintf (stderr, "cmp inp: ");
     101          for (k = 0; k < Nout; k++) { fprintf (stderr, "0x%02hhx ", zdata[k]); }
     102          fprintf (stderr, "\n");
     103        }
     104
     105        status = fits_rdecomp_short ((unsigned char *) zdata, Nout, (unsigned short *) buffer, Nrawpix, blocksize);
     106        myAssert (!status, "failed to rdecomp");
    93107        break;
    94108

  • branches/eam_branches/ohana.20150429/src/libfits/matrix/F_compress_utils.c

    r38355 r38356  
    101101  int i;
    102102
     103  // fprintf (stderr, "swapping %d bytes in pix of size %d bytes...\n", Nzdata, pixsize);
    103104  switch (pixsize) {
    104105    case 1:
    … …  
    106107
    107108    case 2:
    108       for (i = 0; i < 2*Nzdata; i+=2) {
     109      for (i = 0; i < Nzdata; i+=2) {
    109110        DOSWAP (zdata[i+0], zdata[i+1]);
    110111      }
    … …  
    112113
    113114    case 4:
    114       for (i = 0; i < 4*Nzdata; i+=4) {
     115      for (i = 0; i < Nzdata; i+=4) {
    115116        DOSWAP (zdata[i+0], zdata[i+3]);
    116117        DOSWAP (zdata[i+1], zdata[i+2]);
    … …  
    119120
    120121    case 8:
    121       for (i = 0; i < 8*Nzdata; i+=8) {
     122      for (i = 0; i < Nzdata; i+=8) {
    122123        DOSWAP (zdata[i+0], zdata[i+7]);
    123124        DOSWAP (zdata[i+1], zdata[i+6]);
    … …  
    168169
    169170    return (FALSE);
     171}
     172
     173int gfits_cmptype_valid (char *cmptype) {
     174
     175  if (!strcasecmp(cmptype, "NONE")) return TRUE;
     176  if (!strcasecmp(cmptype, "GZIP_1")) return TRUE;
     177  if (!strcasecmp(cmptype, "GZIP_2")) return TRUE;
     178  if (!strcasecmp(cmptype, "PLIO_1")) return TRUE;
     179  if (!strcasecmp(cmptype, "RICE_1")) return TRUE;
     180  if (!strcasecmp(cmptype, "RICE_ONE")) return TRUE;
     181  if (!strcasecmp(cmptype, "HCOMPRESS_1")) return TRUE;
     182  return FALSE;
    170183}
    171184

  • branches/eam_branches/ohana.20150429/src/libfits/matrix/F_uncompress_M.c

    r38355 r38356  
    22# include <gfitsio.h>
    33# include <zlib.h>
     4# define VERBOSE_DUMP 0
    45
    56# define ESCAPE { \
    … …  
    242243  // zdata_pixsize = gfits_vartable_heap_pixsize (zdef.format); XXX
    243244
    244   // size of a pixel in the final image (not needed)
    245   // idata_pixsize = abs(header[0].bitpix) / 8;
    246 
    247245  // size of a pixel in the output from the decompression routine
    248246  int cmp_pixsize = gfits_compressed_data_pixsize (cmptype, header[0].bitpix, optname, optvalue, Noptions);
    … …  
    250248  int raw_pixsize = gfits_uncompressed_data_pixsize (cmptype, header[0].bitpix, optname, optvalue, Noptions);
    251249
     250  // size of a pixel in the final image (not needed)
     251  // int tile_pixsize = abs(header[0].bitpix) / 8;
     252  // fprintf (stderr, "raw_pixsize: %d, cmp_pixsize: %d, tile_pixsize: %d, raw_bitpix: %d\n", raw_pixsize, cmp_pixsize, tile_pixsize, raw_bitpix);
     253
    252254  int Nout_alloc = raw_pixsize*max_tile_size;
    253255  ALLOCATE (out, char, Nout_alloc);
    … …  
    255257  off_t row;
    256258
    257   if (0) {
     259  if (VERBOSE_DUMP) {
    258260    int k;
    259261    fprintf (stderr, "unc tbl: ");
    … …  
    273275    if (!zdata) ESCAPE;
    274276
    275     if (strcasecmp(cmptype, "GZIP_1")) {
    276       if (!gfits_byteswap_zdata (zdata, Nzdata, cmp_pixsize)) ESCAPE;
    277     }
    278 
    279     if (0 && (row == 0)) {
     277    if (VERBOSE_DUMP && (row == 0)) {
    280278      int k;
    281279      fprintf (stderr, "unc dat: ");
    … …  
    284282    }
    285283
     284    if (strcasecmp(cmptype, "GZIP_1")) {
     285      // Nzdata is number of bytes
     286      if (!gfits_byteswap_zdata (zdata, Nzdata, cmp_pixsize)) ESCAPE;
     287    }
     288
     289    if (VERBOSE_DUMP && (row == 0)) {
     290      int k;
     291      fprintf (stderr, "unc SWP: ");
     292      for (k = 0; k < 32; k++) { fprintf (stderr, "0x%02hhx ", zdata[k]); }
     293      fprintf (stderr, "\n");
     294      fprintf (stderr, "Nout  : %d -> ", Nout);
     295    }
     296   
    286297    // gfits_uncompress_data uncompresses from zdata to the temporary output buffer which must be allocated
    287298    // XXX the tile must not be > 2GB
     299    // Nout going in is number of bytes
    288300    if (!gfits_uncompress_data ((char *)zdata, Nzdata, cmptype, optname, optvalue, Noptions, out, &Nout, raw_pixsize)) ESCAPE;
    289301
    290     if (0 && (row == 0)) {
     302    if (VERBOSE_DUMP && (row == 0)) {
    291303      int k;
     304      fprintf (stderr, "%d\n", Nout);
    292305      fprintf (stderr, "unc swp: ");
    293306      for (k = 0; k < 32; k++) { fprintf (stderr, "0x%02hhx ", out[k]); }
    … …  
    296309   
    297310    if (!strcasecmp(cmptype, "GZIP_1")) {
    298       if (!gfits_byteswap_zdata (out, Nout, raw_pixsize)) ESCAPE;
     311      // Nout is number of pixels
     312      if (!gfits_byteswap_zdata (out, Nout * raw_pixsize, raw_pixsize)) ESCAPE;
    299313    }
    300314   
    301     if (0 && (row == 0)) {
     315    if (VERBOSE_DUMP && (row == 0)) {
    302316      int k;
    303317      fprintf (stderr, "unc raw: ");
    … …  
    309323    if (!gfits_distribute_data (matrix, out, Nout, raw_bitpix, otile, oblank, ztile, zblank, zscale, zzero)) ESCAPE;
    310324
    311     if (0 && (row == 0)) {
     325    if (VERBOSE_DUMP && (row == 0)) {
    312326      int k;
    313327      fprintf (stderr, "unc mat: ");

  • branches/eam_branches/ohana.20150429/src/libfits/matrix/F_uncompress_data.c

    r38340 r38356  
    3131
    3232  int status;
    33   static int Ninsum = 0;
    34   static int Noutsum = 0;
    3533
    3634  // do not actually uncompress : this is used for testing
    … …  
    4240
    4341  if (!strcasecmp(cmptype, "GZIP_1")) {
    44     unsigned long tNout = *Nout * out_pixsize;
     42    // unsigned long tNout = *Nout * out_pixsize;
     43    // input value of *Nout is number of BYTES
     44    unsigned long tNout = *Nout;
    4545
    4646    // for GZIP data, I need to check for and remove the header on the first block:
    … …  
    5353    if (status != Z_OK) ESCAPE(FALSE);
    5454    *Nout = tNout / out_pixsize;
     55    // output value of *Nout is number of PIXELS
    5556
    5657    // the resulting uncompressed data is byteswapped
    … …  
    6061  }
    6162
    62   if (!strcasecmp(cmptype, "RICE_1")) {
     63  if (!strcasecmp(cmptype, "RICE_1") || !strcasecmp(cmptype, "RICE_ONE")) {
    6364    int i, blocksize;
    6465    // look for the BLOCKSIZE
    … …  
    7475    }
    7576
    76     int Npix;
    77     // Npix = *Nout * (out_pixsize / 4.0);
    78     Npix = *Nout;
    79 
    80     Ninsum += Nzdata;
    81     Noutsum += *Nout;
    82     // fprintf (stderr, "%d comp bytes; %d uncomp 'pixels', totals: %d %d\n", Nzdata, Npix, Ninsum, Noutsum);
     77    int status = FALSE;
     78    int Npix = *Nout / out_pixsize;
     79    int k;
    8380
    8481    switch (out_pixsize) {
    8582      case 4:
    8683        // rice decompression from the CFITSIO source tree : we need to tell it the expected number of pixels
    87         // is also REQUIRES 4byte output, which is fairly stupid.
    88         if (fits_rdecomp ((unsigned char *) zdata, Nzdata, (unsigned int *) outdata, Npix, blocksize)) {
    89           fprintf (stderr, "error in rice decompression\n");
    90           return (FALSE);
     84        status = fits_rdecomp ((unsigned char *) zdata, Nzdata, (unsigned int *) outdata, Npix, blocksize);
     85        break;
     86      case 2:
     87        if (0) {
     88          fprintf (stderr, "Nout: %d, Nrawpix: %d\n", Nzdata, Npix);
     89          fprintf (stderr, "cmp out: ");
     90          for (k = 0; k < Nzdata; k++) { fprintf (stderr, "0x%02hhx ", zdata[k]); }
     91          fprintf (stderr, "\n");
    9192        }
    92         return (TRUE);
    93 
    94       case 2:
    95         if (fits_rdecomp_short ((unsigned char *) zdata, Nzdata, (unsigned short *) outdata, Npix, blocksize)) {
    96           fprintf (stderr, "error in rice decompression\n");
    97           return (FALSE);
    98         }
    99         return (TRUE);
    100 
     93        status = fits_rdecomp_short ((unsigned char *) zdata, Nzdata, (unsigned short *) outdata, Npix, blocksize);
     94        break;
    10195      case 1:
    102         if (fits_rdecomp_byte ((unsigned char *) zdata, Nzdata, (unsigned char *) outdata, Npix, blocksize)) {
    103           fprintf (stderr, "error in rice decompression\n");
    104           return (FALSE);
    105         }
    106         return (TRUE);
    107        
     96        status = fits_rdecomp_byte ((unsigned char *) zdata, Nzdata, (unsigned char *) outdata, Npix, blocksize);
     97        break;
    10898      default:
    10999        fprintf (stderr, "invalid output pixel size %d\n", out_pixsize);
    110         return (FALSE);
     100        status = FALSE;
     101        break;
    111102    }
    112    
     103    if (status) {
     104      *Nout = 0;
     105      return FALSE;
     106    }
     107    *Nout = Npix;
     108    return TRUE;
    113109  }
    114110 

  • branches/eam_branches/ohana.20150429/src/libfits/table/F_compress_T.c

    r38355 r38356  
    157157      int Nraw = Nrows*fields[j].Nvalues; // number of pixels
    158158      if (!strcmp(fields[j].zctype, "GZIP_1")) {
    159         if (!gfits_byteswap_zdata (raw, Nraw, fields[j].pixsize)) ESCAPE(A);
     159        if (!gfits_byteswap_zdata (raw, Nraw * fields[j].pixsize, fields[j].pixsize)) ESCAPE(A);
    160160      }
    161161

  • branches/eam_branches/ohana.20150429/src/libfits/table/F_uncompress_T.c

    r38348 r38356  
    161161      if (!strcasecmp(fields[i].zctype, "GZIP_1")) {
    162162        myAssert ((fields[i].zdef.format != 'C') && (fields[i].zdef.format != 'M'), "swap is probably wrong for C or M columns");
    163         if (!gfits_byteswap_zdata (raw, Nraw, fields[i].pixsize)) ESCAPE(A);
     163        if (!gfits_byteswap_zdata (raw, Nraw * fields[i].pixsize, fields[i].pixsize)) ESCAPE(A);
    164164      }
    165165

  • branches/eam_branches/ohana.20150429/src/libfits/test/compress.sh

    r38343 r38356  
    22macro test_image_all
    33
    4   foreach mode NONE GZIP_1
     4  foreach mode NONE GZIP_1 RICE_1
    55    foreach bitpix 8 16 32 -32 -64
    66      echo ===== bitpix = $bitpix =====
    … …  
    4040 end
    4141
     42 if ($cmpmode == RICE_1)
     43   echo "*** funpack fails on RICE_1 ****"
     44   return
     45 end
     46
    4247 if ($bitpix == -64)
    4348   echo "*** bitpix $bitpix fails fpack / funpack ***"
    4449 else
     50   echo "run funpack on our version"
    4551   exec funpack -S test.cmp.fits > test.fun.fits
    4652   rd e test.fun.fits

  • branches/eam_branches/ohana.20150429/src/libfits/test/imagecomp.c

    r38340 r38356  
    44
    55int test_compress (int bitpix, char *zcmptype);
     6
     7// char *cmptype[] = {"NONE", "GZIP_1", "GZIP_2", "PLIO_1", "RICE_1", "RICE_ONE", "HCOMPRESS_1", NULL};
     8char *cmptype[] = {"NONE", "GZIP_1", "RICE_1", "RICE_ONE", NULL};
     9int bitpix[] = {8, 16, 32, -32, -64, 0};
     10
     11// char *cmptype[] = {"NONE", "GZIP_1", NULL};
     12// char *cmptype[] = {"PLIO_1", NULL};
     13// int bitpix[] = {8, 16, 32, 0};
     14
     15// char *cmptype[] = {"RICE_1", NULL};
     16// int bitpix[] = {-32, -64, 0};
    617
    718int main (int argc, char **argv) {
    … …  
    1122  diag ("libfits imagecomp.c tests");
    1223
    13   test_compress (  8, "NONE");
    14   test_compress ( 16, "NONE");
    15   test_compress ( 32, "NONE");
    16   test_compress (-32, "NONE");
    17   test_compress (-64, "NONE");
     24  int i, j;
     25  for (i = 0; cmptype[i]; i++) {
     26    for (j = 0; bitpix[j]; j++) {
     27      test_compress (bitpix[j], cmptype[i]);
     28    }
     29  }
    1830
    19   test_compress (  8, "GZIP_1");
    20   test_compress ( 16, "GZIP_1");
    21   test_compress ( 32, "GZIP_1");
    22   test_compress (-32, "GZIP_1");
    23   test_compress (-64, "GZIP_1");
     31  // test_compress ( 16, "NONE");
     32  // test_compress ( 32, "NONE");
     33  // test_compress (-32, "NONE");
     34  // test_compress (-64, "NONE");
     35  //
     36  // test_compress (  8, "GZIP_1");
     37  // test_compress ( 16, "GZIP_1");
     38  // test_compress ( 32, "GZIP_1");
     39  // test_compress (-32, "GZIP_1");
     40  // test_compress (-64, "GZIP_1");
    2441
    2542  exit (0);
    … …  
    2744
    2845# define NX 100
    29 # define NY 100
     46# define NY 10
    3047
    3148int test_compress (int bitpix, char *zcmptype) { // make a table, compress, uncompress, compare : use compression zcmptype
    … …  
    6582  for (ix = 0; ix < NX; ix++) {
    6683    for (iy = 0; iy < NY; iy++) {
     84      int IY = iy + 1;
    6785      switch (bitpix) {
    68         case   8: rawdata_char  [ix + NX*iy] = ix + iy*iy; break;
    69         case  16: rawdata_short [ix + NX*iy] = ix + iy*iy; break;
    70         case  32: rawdata_int   [ix + NX*iy] = ix + iy*iy; break;
    71         case -32: rawdata_float [ix + NX*iy] = ix + iy*iy; break;
    72         case -64: rawdata_double[ix + NX*iy] = ix + iy*iy; break;
     86        case   8: rawdata_char  [ix + NX*iy] = ix + IY*IY; break;
     87        case  16: rawdata_short [ix + NX*iy] = ix + IY*IY; break;
     88        case  32: rawdata_int   [ix + NX*iy] = ix + IY*IY; break;
     89        case -32: rawdata_float [ix + NX*iy] = ix + IY*IY; break;
     90        case -64: rawdata_double[ix + NX*iy] = ix + IY*IY; break;
    7391        default: myAbort ("bad bitpix value");
    7492      }
Note: See TracChangeset for help on using the changeset viewer.