Changeset 31633 for trunk/psLib/src/fits/psFitsScale.c

Timestamp:

Jun 15, 2011, 6:09:25 PM (15 years ago)

Author:

watersc1

Message:

merge from czw branch: ASINH code

File:

• trunk/psLib/src/fits/psFitsScale.c (modified) (14 diffs)

trunk/psLib/src/fits/psFitsScale.c

@@ -214 +214 @@
       case PS_FITS_SCALE_STDEV_POSITIVE:
       case PS_FITS_SCALE_LOG_STDEV_POSITIVE:
+      case PS_FITS_SCALE_ASINH_STDEV_POSITIVE:
         // Put (mean - N sigma) at the lowest possible value: predominantly positive images
         imageVal = mean - options->stdevNum * stdev;
@@ -220 +221 @@
       case PS_FITS_SCALE_STDEV_NEGATIVE:
       case PS_FITS_SCALE_LOG_STDEV_NEGATIVE:
+      case PS_FITS_SCALE_ASINH_STDEV_NEGATIVE:
         // Put (mean + N sigma) at the highest possible value: predominantly negative images
         imageVal = mean + options->stdevNum * stdev;
@@ -226 +228 @@
       case PS_FITS_SCALE_STDEV_BOTH:
       case PS_FITS_SCALE_LOG_STDEV_BOTH:
+      case PS_FITS_SCALE_ASINH_STDEV_BOTH:
         // Put mean right in the middle: images with an equal abundance of positive and negative values
         imageVal = mean;
@@ -241 +244 @@
 }
 
+
+    
 
 static bool logscaleStdev(double *bscale, // Scaling, to return
@@ -456 +461 @@
 }
 
+static bool asinhStdev(double *bscale, // Scaling, to return
+                       double *bzero,  // Zero point, to return
+                       double *boffset, // asinh flux zeropoint, to return
+                       double *bsoften, // asinh softening parameter, to return
+                       const psImage *image, // Image to scale
+                       const psImage *mask,  // Mask image
+                       psImageMaskType maskVal, // value to mask
+                       const psFitsOptions *options // FITS options
+                       )
+{
+  psAssert(bscale, "impossible");
+  psAssert(bzero, "impossible");
+  psAssert(boffset, "impossible");
+  psAssert(bsoften, "impossible");
+  psAssert(image, "impossible");
+  psAssert(options, "impossible");
+
+  psTrace("psLib.fits", 3, "Scaling image by asinh method");
+  int numCols = image->numCols, numRows = image->numRows; // Size of image
+  
+    // Measure the mean and stdev
+    // psImageBackground automatically excludes pixels that are non-finite, so we don't need to bother about a
+    // mask.
+    psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS);
+    psStats *stats = psStatsAlloc(MEAN_STAT | STDEV_STAT); // Statistics object
+    double mean, stdev;                                    // Mean and standard deviation
+    if (!psImageBackground(stats, NULL, image, mask, maskVal, rng)) {
+        // It could be because the image is entirely masked, in which case we don't want to error
+        bool good = false;              // Any good pixels?
+
+
+// Find good pixels in an image, by image type
+#define GOOD_PIXELS_CASE(TYPE) \
+      case PS_TYPE_##TYPE: \
+        for (int y = 0; y < image->numRows && !good; y++) { \
+            for (int x = 0; x < image->numCols && !good; x++) { \
+                if (mask && (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskVal)) { \
+                    continue; \
+                } \
+                if (!isfinite(image->data.TYPE[y][x])) { \
+                    continue; \
+                } \
+                good = true; \
+            } \
+        } \
+        break;
+
+        switch (image->type.type) {
+            GOOD_PIXELS_CASE(F32);
+            GOOD_PIXELS_CASE(F64);
+          default:
+            psAbort("Unsupported case: %x", image->type.type);
+        }
+
+        if (!good) {
+            psLogMsg("psLib.fits", PS_LOG_DETAIL, "Image has no good pixels, setting BSCALE = 1, BZERO = 0");
+            psErrorClear();
+            *bscale = 1.0;
+            *bzero = 0.0;
+            *bsoften = NAN;
+            psFree(rng);
+            psFree(stats);
+            return true;
+        }
+
+        // There are some good pixels in there somewhere; psImageBackground just didn't find them
+        psLogMsg("psLib.fits", PS_LOG_DETAIL,
+                 "Couldn't measure background statistics for image quantisation; retrying.");
+        psErrorClear();
+        // Retry using all the available pixels
+        stats->nSubsample = image->numCols * image->numRows + 1;
+        if (!psImageStats(stats, image, mask, maskVal)) {
+            psLogMsg("psLib.fits", PS_LOG_DETAIL,
+                     "Couldn't measure background statistics for image quantisation (attempt 2); retrying.");
+            psErrorClear();
+            // Retry with desperate statistic
+            stats->options = DESPERATE_MEAN_STAT | DESPERATE_STDEV_STAT;
+            if (!psImageStats(stats, image, mask, maskVal)) {
+                psError(PS_ERR_UNKNOWN, false, "Unable to measure background statistics for image");
+                psFree(rng);
+                psFree(stats);
+                return false;
+            } else {
+                // Desperate retry
+                mean = psStatsGetValue(stats, DESPERATE_MEAN_STAT);
+                stdev = psStatsGetValue(stats, DESPERATE_STDEV_STAT);
+            }
+        } else {
+            // Retry with all available pixels
+            mean = psStatsGetValue(stats, MEAN_STAT);
+            stdev = psStatsGetValue(stats, STDEV_STAT);
+        }
+    } else {
+        // First attempt
+        mean = psStatsGetValue(stats, MEAN_STAT);
+        stdev = psStatsGetValue(stats, STDEV_STAT);
+    }
+    psFree(rng);
+    psFree(stats);
+    if (!isfinite(mean) || !isfinite(stdev)) {
+        psError(PS_ERR_BAD_PARAMETER_VALUE, true,
+                "Mean (%f) or stdev (%f) is non-finite.", mean, stdev);
+        return false;
+    }
+
+    psImage *copy;
+
+    switch (image->type.type) {
+    case PS_TYPE_F32:
+      copy = psImageCopy(NULL,image,PS_TYPE_F32);
+      break;
+    case PS_TYPE_F64:
+      copy = psImageCopy(NULL,image,PS_TYPE_F64);
+      break;
+    default:
+      psError(PS_ERR_UNKNOWN, true, "Target type is not a float: %d", image->type.type);
+      return NULL;
+      break;
+    }
+    // Do scaling
+    float a = 1.0857362; // 2.5 * log(e);
+    *bsoften = sqrt(a) * stdev;
+    *boffset = mean;
+    //    float m0 = 0; // Can I just arbitrarily set this?
+    
+    switch (image->type.type) {
+    case PS_TYPE_F32: 
+      for (int y = 0; y < numRows; y++) {
+        for (int x = 0; x < numCols; x++) {
+/*        if (x == 266 && y == 4584) { */
+/*          fprintf(stderr,"psFS32: %d %d %g %g %g %g %g\n",x,y,*boffset,*bsoften,image->data.F32[y][x],log10(image->data.F32[y][x] - *boffset), */
+/*                  a * asinh((image->data.F32[y][x] - *boffset) / (2 * *bsoften))); */
+/*        } */
+          if (isfinite(image->data.F32[y][x])) {
+            copy->data.F32[y][x] = a * asinh( (image->data.F32[y][x] - *boffset) / (2 * *bsoften));// - 2.5 * log10(b) + m0;
+          }
+          else {
+            copy->data.F32[y][x] = image->data.F32[y][x];
+          }
+        }
+      }
+      break;
+    case PS_TYPE_F64:
+        for (int y = 0; y < numRows; y++) {
+          for (int x = 0; x < numCols; x++) {
+            //      fprintf(stderr,"psFS64: %d %d %g %g %g\n",x,y,offset,image->data.F64[y][x],log10(image->data.F64[y][x] - offset));
+            if (isfinite(image->data.F64[y][x])) {
+              copy->data.F64[y][x] = a * asinh( (image->data.F64[y][x] - *boffset) / (2 * *bsoften));// - 2.5 * log10(b) + m0;
+            }
+            else {
+              copy->data.F64[y][x] = image->data.F64[y][x];
+            }     
+          }
+        }
+      break;
+    default:
+      psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Target type is not a float: %d",image->type.type);
+      return NULL;
+      break;
+    }
+      
+    // Do regular scaling on the logarithm image
+    if (!scaleStdev(bscale, bzero, copy, mask, maskVal, options)) {
+      psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from stdev");
+      return false;
+    }
+    psFree(copy);
+    return true;
+}
+
+static bool asinhRange(double *bscale, // Scaling, to return
+                       double *bzero,  // Zero point, to return
+                       double *boffset, // asinh flux zeropoint, to return
+                       double *bsoften, // asinh softening parameter, to return
+                       const psImage *image, // Image to scale
+                       const psImage *mask,  // Mask image
+                       psImageMaskType maskVal, // value to mask
+                       const psFitsOptions *options // FITS options
+                       )
+{
+  psAssert(bscale, "impossible");
+  psAssert(bzero, "impossible");
+  psAssert(bsoften, "impossible");
+  psAssert(image, "impossible");
+  psAssert(options, "impossible");
+
+  psTrace("psLib.fits", 3, "Scaling image by asinh method");
+  int numCols = image->numCols, numRows = image->numRows; // Size of image
+  
+    // Measure the mean and stdev
+    // psImageBackground automatically excludes pixels that are non-finite, so we don't need to bother about a
+    // mask.
+    psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS);
+    psStats *stats = psStatsAlloc(MEAN_STAT | STDEV_STAT); // Statistics object
+    double mean, stdev;                                    // Mean and standard deviation
+    if (!psImageBackground(stats, NULL, image, mask, maskVal, rng)) {
+        // It could be because the image is entirely masked, in which case we don't want to error
+        bool good = false;              // Any good pixels?
+
+
+// Find good pixels in an image, by image type
+#define GOOD_PIXELS_CASE(TYPE) \
+      case PS_TYPE_##TYPE: \
+        for (int y = 0; y < image->numRows && !good; y++) { \
+            for (int x = 0; x < image->numCols && !good; x++) { \
+                if (mask && (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskVal)) { \
+                    continue; \
+                } \
+                if (!isfinite(image->data.TYPE[y][x])) { \
+                    continue; \
+                } \
+                good = true; \
+            } \
+        } \
+        break;
+
+        switch (image->type.type) {
+            GOOD_PIXELS_CASE(F32);
+            GOOD_PIXELS_CASE(F64);
+          default:
+            psAbort("Unsupported case: %x", image->type.type);
+        }
+
+        if (!good) {
+            psLogMsg("psLib.fits", PS_LOG_DETAIL, "Image has no good pixels, setting BSCALE = 1, BZERO = 0");
+            psErrorClear();
+            *bscale = 1.0;
+            *bzero = 0.0;
+            *bsoften = NAN;
+            psFree(rng);
+            psFree(stats);
+            return true;
+        }
+
+        // There are some good pixels in there somewhere; psImageBackground just didn't find them
+        psLogMsg("psLib.fits", PS_LOG_DETAIL,
+                 "Couldn't measure background statistics for image quantisation; retrying.");
+        psErrorClear();
+        // Retry using all the available pixels
+        stats->nSubsample = image->numCols * image->numRows + 1;
+        if (!psImageStats(stats, image, mask, maskVal)) {
+            psLogMsg("psLib.fits", PS_LOG_DETAIL,
+                     "Couldn't measure background statistics for image quantisation (attempt 2); retrying.");
+            psErrorClear();
+            // Retry with desperate statistic
+            stats->options = DESPERATE_MEAN_STAT | DESPERATE_STDEV_STAT;
+            if (!psImageStats(stats, image, mask, maskVal)) {
+                psError(PS_ERR_UNKNOWN, false, "Unable to measure background statistics for image");
+                psFree(rng);
+                psFree(stats);
+                return false;
+            } else {
+                // Desperate retry
+                mean = psStatsGetValue(stats, DESPERATE_MEAN_STAT);
+                stdev = psStatsGetValue(stats, DESPERATE_STDEV_STAT);
+            }
+        } else {
+            // Retry with all available pixels
+            mean = psStatsGetValue(stats, MEAN_STAT);
+            stdev = psStatsGetValue(stats, STDEV_STAT);
+        }
+    } else {
+        // First attempt
+        mean = psStatsGetValue(stats, MEAN_STAT);
+        stdev = psStatsGetValue(stats, STDEV_STAT);
+    }
+    psFree(rng);
+    psFree(stats);
+    if (!isfinite(mean) || !isfinite(stdev)) {
+        psError(PS_ERR_BAD_PARAMETER_VALUE, true,
+                "Mean (%f) or stdev (%f) is non-finite.", mean, stdev);
+        return false;
+    }
+
+    psImage *copy;
+
+    switch (image->type.type) {
+    case PS_TYPE_F32:
+      copy = psImageCopy(NULL,image,PS_TYPE_F32);
+      break;
+    case PS_TYPE_F64:
+      copy = psImageCopy(NULL,image,PS_TYPE_F64);
+      break;
+    default:
+      psError(PS_ERR_UNKNOWN, true, "Target type is not a float: %d", image->type.type);
+      return NULL;
+      break;
+    }
+    // Do scaling
+    float a = 1.0857362; // 2.5 * log(e);
+    *bsoften = sqrt(a) * stdev;
+    *boffset = mean;
+    // float m0 = 0; // Can I just arbitrarily set this?
+    
+    switch (image->type.type) {
+    case PS_TYPE_F32: 
+      for (int y = 0; y < numRows; y++) {
+        for (int x = 0; x < numCols; x++) {
+/*        if (x == 2331 && y == 2843) { */
+/*          fprintf(stderr,"psFS32: %d %d %g %g %g\n",x,y,offset,image->data.F32[y][x],log10(image->data.F32[y][x] - offset)); */
+/*        } */
+/*        if (x == 266 && y == 4584) { */
+/*          fprintf(stderr,"psFS32: %d %d %g %g %g %g %g\n",x,y,*boffset,*bsoften,image->data.F32[y][x],log10(image->data.F32[y][x] - *boffset), */
+/*                  a * asinh((image->data.F32[y][x] - *boffset) / (2 * *bsoften))); */
+/*        } */
+
+          if (isfinite(image->data.F32[y][x])) {
+            copy->data.F32[y][x] = a * asinh( (image->data.F32[y][x] - *boffset) / (2 * *bsoften));// - 2.5 * log10(b) + m0;
+          }
+          else {
+            copy->data.F32[y][x] = image->data.F32[y][x];
+          }
+        }
+      }
+      break;
+    case PS_TYPE_F64:
+        for (int y = 0; y < numRows; y++) {
+          for (int x = 0; x < numCols; x++) {
+            //      fprintf(stderr,"psFS64: %d %d %g %g %g\n",x,y,offset,image->data.F64[y][x],log10(image->data.F64[y][x] - offset));
+            if (isfinite(image->data.F64[y][x])) {
+              copy->data.F64[y][x] = a * asinh( (image->data.F64[y][x] - *boffset)/ (2 * *bsoften));// - 2.5 * log10(b) + m0;
+            }
+            else {
+              copy->data.F64[y][x] = image->data.F32[y][x];
+            }
+          }
+        }
+      break;
+    default:
+      psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Target type is not a float: %d",image->type.type);
+      return NULL;
+      break;
+    }
+      
+    // Do regular scaling on the asinh image
+    if (!scaleRange(bscale, bzero, copy, options)) {
+      psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from stdev");
+      return false;
+    }
+    psFree(copy);
+    return true;
+}
+
 
 //////////////////////////////////////////////////////////////////////////////////////////////////////////////
@@ -461 +809 @@
 //////////////////////////////////////////////////////////////////////////////////////////////////////////////
 
-bool psFitsScaleDetermine(double *bscale, double *bzero, double *boffset, long *blank, const psImage *image,
+bool psFitsScaleDetermine(double *bscale, double *bzero, double *boffset, double *bsoften,
+                          long *blank, const psImage *image,
                           const psImage *mask, psImageMaskType maskVal, const psFits *fits)
 {
@@ -521 +870 @@
         }
         break;
-    case PS_FITS_SCALE_LOG_RANGE:
-      if (!logscaleRange(bscale,bzero,boffset,image,options)) {
-        psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from range");
-        return false;
-      }
-      break;
-    case PS_FITS_SCALE_LOG_STDEV_POSITIVE:
-    case PS_FITS_SCALE_LOG_STDEV_NEGATIVE:
-    case PS_FITS_SCALE_LOG_STDEV_BOTH:
-      if (!logscaleStdev(bscale, bzero,boffset, image, mask, maskVal, options)) {
-            psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from stdev");
-            return false;
-        }
-        break;
+      case PS_FITS_SCALE_LOG_RANGE:
+        if (!logscaleRange(bscale,bzero,boffset,image,options)) {
+          psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from range");
+          return false;
+        }
+        break;
+      case PS_FITS_SCALE_LOG_STDEV_POSITIVE:
+      case PS_FITS_SCALE_LOG_STDEV_NEGATIVE:
+      case PS_FITS_SCALE_LOG_STDEV_BOTH:
+        if (!logscaleStdev(bscale, bzero,boffset, image, mask, maskVal, options)) {
+          psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from stdev");
+          return false;
+        }
+        break;
+      case PS_FITS_SCALE_ASINH_RANGE:
+        if (!asinhRange(bscale,bzero,boffset,bsoften,image,mask,maskVal,options)) {
+          psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from range");
+          return false;
+        }
+        break;
+      case PS_FITS_SCALE_ASINH_STDEV_POSITIVE:
+      case PS_FITS_SCALE_ASINH_STDEV_NEGATIVE:
+      case PS_FITS_SCALE_ASINH_STDEV_BOTH:
+        if (!asinhStdev(bscale, bzero,boffset,bsoften, image, mask, maskVal, options)) {
+          psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from stdev");
+          return false;
+        }
+        break;
+        
+        
       case PS_FITS_SCALE_MANUAL:
         *bscale = options->bscale;
@@ -557 +922 @@
 
 
-    psTrace("psLib.fits", 3, "BSCALE = %.10lf, BZERO = %.10lf, BOFFSET = %.10lf, BLANK = %ld\n",
-            *bscale, *bzero, *boffset, *blank);
+    psTrace("psLib.fits", 3, "BSCALE = %.10lf, BZERO = %.10lf, BOFFSET = %.10lf, BSOFTEN = %.10lf, BLANK = %ld\n",
+            *bscale, *bzero, *boffset, *bsoften, *blank);
     return true;
 }
 
 
-psImage *psFitsScaleForDisk(const psImage *image, const psFits *fits, double bscale, double bzero, double boffset,
+psImage *psFitsScaleForDisk(const psImage *image, const psFits *fits, double bscale, double bzero, double boffset, double bsoften,
                             psRandom *rng)
 {
@@ -625 +990 @@
                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_NEGATIVE)|| \
                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_BOTH)) { \
-                value = log10( (IN)->data.INTYPE[y][x] - boffset ); \
+                if (isfinite((IN)->data.INTYPE[y][x])) {                \
+                  value = log10( (IN)->data.INTYPE[y][x] - boffset );   \
+                }                                                       \
+                else {                                                  \
+                  value = (IN)->data.INTYPE[y][x];                      \
+                }                                                       \
+              }                                                         \
+              else if ((options->scaling == PS_FITS_SCALE_ASINH_RANGE)||        \
+                (options->scaling == PS_FITS_SCALE_ASINH_MANUAL)||        \
+                  (options->scaling == PS_FITS_SCALE_ASINH_STDEV_POSITIVE)|| \
+                  (options->scaling == PS_FITS_SCALE_ASINH_STDEV_NEGATIVE)|| \
+                  (options->scaling == PS_FITS_SCALE_ASINH_STDEV_BOTH)) {       \
+                if (isfinite((IN)->data.INTYPE[y][x])) {                \
+                  value = 1.0857362 * (asinh( ((IN)->data.INTYPE[y][x] - boffset) / (2.0 * bsoften))); \
+                }                                                       \
+                else {                                                  \
+                  value = (IN)->data.INTYPE[y][x];                      \
+                }                                                       \
+                if ((x == 1000)&&(y == 1000)) {                         \
+                  fprintf(stderr,"ASINH: %f %f %f",                     \
+                          (IN)->data.INTYPE[y][x],value,bsoften);       \
+                }                                                       \
               }                                                         \
                 else { \
@@ -652 +1038 @@
     case PS_TYPE_##INTYPE: { \
         switch (outType) { \
-            SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, U8); \
-            SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S16); \
-            SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S32); \
-            SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S64); \
+          SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, U8);;   \
+          SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S16);;  \
+          SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S32);;  \
+          SCALE_WRITE_OUT_CASE(IN, INTYPE, OUT, S64);;  \
           default: \
             psAbort("Should be unreachable."); \
@@ -664 +1050 @@
     switch (image->type.type) {
         SCALE_WRITE_IN_CASE(image, F32, out);
-        SCALE_WRITE_IN_CASE(image, F64, out);
+        SCALE_WRITE_IN_CASE(image, F64, out); 
       default:
         psAbort("Should be unreachable.");
@@ -679 +1065 @@
 // the present time, since cfitsio should apply the scaling itself in the process of reading.  However, we may
 // later desire it (e.g., if we ever make our own FITS implementation).
-psImage *psFitsScaleFromDisk(const psImage *image, double boffset)
+psImage *psFitsScaleFromDisk(const psImage *image, double boffset, double bsoften)
 {
     PS_ASSERT_IMAGE_NON_NULL(image, NULL);
@@ -718 +1104 @@
       for (int y = 0; y < numRows; y++) { \
           for (int x = 0; x < numCols; x++) { \
-            out->data.OUTTYPE[y][x] = pow(10,image->data.INTYPE[y][x]) + boffset;; \
+            if (bsoften) {                                              \
+              out->data.OUTTYPE[y][x] = 2 * bsoften * sinh(image->data.INTYPE[y][x] / (1.0857362)) + boffset; \
+            }                                                           \
+            else if (boffset) {                                 \
+              out->data.OUTTYPE[y][x] = pow(10,image->data.INTYPE[y][x]) + boffset;; \
+             }                                                                  \
           } \
       } \
@@ -766 +1157 @@
     if (strcasecmp(string, "LOG_STDEV_NEGATIVE") == 0) return PS_FITS_SCALE_LOG_STDEV_NEGATIVE;
     if (strcasecmp(string, "LOG_STDEV_BOTH") == 0) return PS_FITS_SCALE_LOG_STDEV_BOTH;
+    if (strcasecmp(string, "ASINH_RANGE") == 0)      return PS_FITS_SCALE_ASINH_RANGE;
+    if (strcasecmp(string, "ASINH_MANUAL") == 0)      return PS_FITS_SCALE_ASINH_MANUAL;
+    if (strcasecmp(string, "ASINH_STDEV_POSITIVE") == 0) return PS_FITS_SCALE_ASINH_STDEV_POSITIVE;
+    if (strcasecmp(string, "ASINH_STDEV_NEGATIVE") == 0) return PS_FITS_SCALE_ASINH_STDEV_NEGATIVE;
+    if (strcasecmp(string, "ASINH_STDEV_BOTH") == 0) return PS_FITS_SCALE_ASINH_STDEV_BOTH;
     if (strcasecmp(string, "MANUAL") == 0)         return PS_FITS_SCALE_MANUAL;