Changeset 30636 for trunk/psLib/src/fits/psFitsImage.c

Timestamp:

Feb 14, 2011, 2:58:01 PM (15 years ago)

Author:

watersc1

Message:

merge from czw_branch of logflux code

Location:

trunk

Files:

• . (modified) (1 prop)
• psLib/src/fits/psFitsImage.c (modified) (13 diffs)

trunk/psLib/src/fits/psFitsImage.c

@@ -54 +54 @@
     int fitsDatatype;                   // cfitsio data type
     int psDatatype;                     // psLib data type
+    bool is_logscaled;                  // is this image log scaled using BOFFSET?
 } p_psFitsReadInfo;
 
@@ -128 +129 @@
 
     // Check scale and zero
-    double bscale = 0.0, bzero = 0.0;    // Scale and zero point
+    double bscale = 0.0, bzero = 0.0, boffset = NAN;    // Scale and zero point
     if (fits_read_key_dbl(fits->fd, "BSCALE", &bscale, NULL, &status) && status != KEY_NO_EXIST) {
         psFitsError(status, true, "Unable to read header.");
@@ -137 +138 @@
         psFitsError(status, true, "Unable to read header.");
         goto bad;
+    }
+    status = 0;
+    if (fits_read_key_dbl(fits->fd, "BOFFSET", &boffset, NULL, &status) && status != KEY_NO_EXIST) {
+        psFitsError(status, true, "Unable to read header.");
+        goto bad;
+    }
+    if (status == KEY_NO_EXIST) {
+      info->is_logscaled = false;
+    }
+    else if (isfinite(boffset)) {
+      info->is_logscaled = true;
+    }
+    else {
+      info->is_logscaled = false;
     }
     status = 0;
@@ -246 +261 @@
 static psImage *imageToDiskRepresentation(double *bscale, // Scaling applied
                                           double *bzero, // Zero point applied
+                                          double *boffset, // Log offset applied
                                           long *blank, // Blank value (integer data)
                                           psFitsFloat *floatType, // Type of custom floating-point
@@ -258 +274 @@
     psAssert(bscale, "impossible");
     psAssert(bzero, "impossible");
+    psAssert(boffset, "impossible");
     psAssert(floatType, "impossible");
     psAssert(fits, "impossible");
@@ -305 +322 @@
         if (newScaleZero) {
             // Choose an appropriate BSCALE and BZERO
-            if (!psFitsScaleDetermine(bscale, bzero, blank, image, mask, maskVal, fits)) {
+          if (!psFitsScaleDetermine(bscale, bzero, boffset, blank, image, mask, maskVal, fits)) {
                 // We can't have the write dying for this reason --- try to save it somehow!
                 psWarning("Unable to determine BSCALE and BZERO for image --- refusing to quantise.");
@@ -331 +348 @@
         }
 
-        return psFitsScaleForDisk(image, fits, *bscale, *bzero, rng);
+        return psFitsScaleForDisk(image, fits, *bscale, *bzero, *boffset, rng);
     }
 
@@ -459 +476 @@
         return NULL;
     }
-    psFree(info);
 
     if (floatType != PS_FITS_FLOAT_NONE) {
         outImage = psFitsFloatImageFromDisk(outImage, inImage, floatType);
     }
+
+    // Need to apply BOFFSET if info->is_logscaled is true
+    if (info->is_logscaled) {
+      double boffset;
+      int status;
+      fits_read_key_dbl(fits->fd, "BOFFSET", &boffset, NULL, &status);
+      outImage = psFitsScaleFromDisk(outImage,boffset);
+    }
+    psFree(info);
+
     psFree(inImage);
 
@@ -572 +598 @@
 
     double bscale = NAN, bzero = NAN;   // Scale and zero point to put in header
+    double boffset = NAN;               // Log offset to put into header.
     long blank = 0;                     // Blank (undefined) value for image
     psFitsFloat floatType;              // Custom floating-point convention type
-    psImage *diskImage = imageToDiskRepresentation(&bscale, &bzero, &blank, &floatType, fits, image,
+    psImage *diskImage = imageToDiskRepresentation(&bscale, &bzero, &boffset, &blank, &floatType, fits, image,
                                                    mask, maskVal, NULL, true); // Image to write out
     if (!diskImage) {
@@ -610 +637 @@
 
     psFitsOptions *options = fits->options; // FITS I/O options
+/*     if (options) { */
+/*       if (options->scaling == PS_FITS_SCALE_LOG_RANGE) { */
+/*      fprintf(stderr,"it has the scaling I expect\n"); */
+/*       } */
+/*       else { */
+/*      fprintf(stderr,"it does nto have the scaling I expect\n"); */
+/*       } */
+/*     } */
+/*     else { */
+/*       fprintf(stderr,"options is null, apparently? \n"); */
+/*     } */
+     
     psAssert(!useRequestedScale || !options || bitPix == options->bitpix || options->bitpix == 0,
              "Something's not consistent");
@@ -651 +690 @@
     }
 
+    // Remove any BOFFSET values that exist in the header if we are not using that scaling anymore
+    if (options&&(!((options->scaling == PS_FITS_SCALE_LOG_RANGE)||
+                    (options->scaling == PS_FITS_SCALE_LOG_MANUAL)||
+                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_POSITIVE)||
+                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_NEGATIVE)||
+                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_BOTH)))) {
+      if (psMetadataLookup(header,"BOFFSET")) {
+        psMetadataRemoveKey(header,"BOFFSET");
+      }
+    }   
+
     // write the header, if any.
     if (header && !psFitsWriteHeaderImage(fits, header, createPHU)) {
@@ -668 +718 @@
         fits_write_key_dbl(fits->fd, "BSCALE", bscale, 12,
                            "Scaling: TRUE = BZERO + BSCALE * DISK", &status);
+        if (options&&(((options->scaling == PS_FITS_SCALE_LOG_RANGE)||
+                       (options->scaling == PS_FITS_SCALE_LOG_MANUAL)||
+                       (options->scaling == PS_FITS_SCALE_LOG_STDEV_POSITIVE)||
+                       (options->scaling == PS_FITS_SCALE_LOG_STDEV_NEGATIVE)||
+                       (options->scaling == PS_FITS_SCALE_LOG_STDEV_BOTH)))) {
+          fits_write_key_dbl(fits->fd, "BOFFSET", boffset, 12,
+                             "Scaling: TRUE = BZERO + BSCALE * 10**(DISK) + BOFFSET)", &status);
+        }       
         if (psFitsError(status, true, "Could not write BSCALE/BZERO headers to file.")) {
             success = false;
@@ -769 +827 @@
     bool success = true;                // Successful update?
     double bscale = NAN, bzero = NAN;   // Scale and zero point to put in header
+    double boffset = NAN;               // Log offset to put in header
     long blank = 0;                     // Blank (undefined) value for image
     psFitsFloat floatType;              // Custom floating-point convention type
-    psImage *diskImage = imageToDiskRepresentation(&bscale, &bzero, &blank, &floatType, fits, input,
+    psImage *diskImage = imageToDiskRepresentation(&bscale, &bzero, &boffset, &blank, &floatType, fits, input,
                                                    mask, maskVal, NULL, false); // Image to write out
     if (!diskImage) {