IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 30636 for trunk


Ignore:
Timestamp:
Feb 14, 2011, 2:58:01 PM (15 years ago)
Author:
watersc1
Message:

merge from czw_branch of logflux code

Location:
trunk
Files:
14 edited

Legend:

Unmodified
Added
Removed

  • trunk

  • trunk/ippScripts/scripts/ipp_apply_burntool_single.pl

  • trunk/ippScripts/scripts/stack_skycell.pl

    r30598 r30636  
    172172my $photometry = metadataLookupBool($recipe, 'PHOTOMETRY'); # perform photometry?
    173173my $output_nocomp = metadataLookupBool($recipe, 'OUTPUT.NOCOMP'); # change filerules to produced uncompressed output images
    174 
     174my $output_logflux = metadataLookupBool($recipe, 'OUTPUT.LOGFLUX'); # change filerules to produce logflux compressed output images.
     175
     176if ($output_nocomp and $output_logflux) {
     177    &my_die("Unable to not compress and logflux compress simultaneously. Check config.",$stack_id, $PS_EXIT_CONFIG_ERROR);
     178}
    175179
    176180# Generate MDC file with the inputs
    … …  
    252256        $command .= " -F PPSTACK.UNCONV.EXPWT PPSTACK.UNCONV.EXPWT.NOCOMP";
    253257    }
     258    if ($output_logflux) {
     259        $command .= " -R PPSTACK.OUTPUT FITS.TYPE COMP_STACK "; # Just this one output component?
     260        $command .= " -R PPSTACK.OUTPUT.VARIANCE FITS.TYPE COMP_STACK ";
     261        $command .= " -R PPSTACK.OUTPUT.EXPWT FITS.TYPE COMP_STACK ";
     262        $command .= " -R PPSTACK.UNCONV FITS.TYPE COMP_STACK ";
     263        $command .= " -R PPSTACK.UNCONV.VARIANCE FITS.TYPE COMP_STACK ";
     264        $command .= " -R PPSTACK.UNCONV.EXPWT FITS.TYPE COMP_STACK ";
     265    }
    254266    $command .= " -threads $threads" if defined $threads;
    255267    $command .= " -debug-stack" if defined $debug;

  • trunk/ippconfig

  • trunk/ippconfig/gpc1/ppStack.config

    r29766 r30636  
    2020PS1_REFERENCE METADATA
    2121    OUTPUT.NOCOMP   BOOL    TRUE
     22    OUTPUT.LOGFLUX  BOOL    FALSE
    2223END
     24
     25STACK_DEEP    METADATA
     26    OUTPUT.NOCOMP   BOOL    TRUE
     27    OUTPUT.LOGFLUX  BOOL    FALSE
     28END
     29
     30STACK_NIGHTLY  METADATA
     31    OUTPUT.NOCOMP   BOOL    FALSE
     32    OUTPUT.LOGFLUX  BOOL    FALSE
     33END
     34
     35STACK_THREEPI  METADATA
     36    OUTPUT.NOCOMP   BOOL    FALSE
     37    OUTPUT.LOGFLUX  BOOL    TRUE
     38END

  • trunk/ippconfig/recipes/fitstypes.mdc

    r29021 r30636  
    9898        NOISE           S32     8
    9999END
     100# Compressed log flux image for stacks.
     101COMP_STACK     METADATA
     102        BITPIX          S32     16
     103        SCALING         STR     LOG_STDEV_POSITIVE
     104        STDEV.BITS      S32     4
     105        STDEV.NUM       F32     10
     106        COMPRESSION     STR     RICE
     107        TILE.X          S32     0
     108        TILE.Y          S32     1
     109        TILE.Z          S32     1
     110        NOISE           S32     8
     111END
    100112
    101113# Compressed exposure image

  • trunk/ippconfig/recipes/ppStack.config

    r29766 r30636  
    112112        PHOTOMETRY          BOOL  FALSE           # Do basic photometry?
    113113END
     114
     115
     116#
     117
     118STACK_DEEP    METADATA
     119END
     120
     121STACK_NIGHTLY  METADATA
     122END
     123
     124STACK_THREEPI  METADATA
     125END

  • trunk/ippconfig/recipes/reductionClasses.mdc

    r29986 r30636  
    192192        JPEG_BIN1       STR     PPIMAGE_J1
    193193        JPEG_BIN2       STR     PPIMAGE_J2
     194END
     195
     196# reduction classes for different stacks
     197# regular stacks
     198NIGHTLY_STACK         METADATA
     199      STACK_PPSTACK   STR      STACK_NIGHTLY
     200      STACK_PPSUB     STR      STACK
     201      STACK_PSPHOT    STR      STACK
     202END
     203
     204# deep stacks
     205DEEP_STACK            METADATA
     206      STACK_PPSTACK   STR      STACK_DEEP
     207      STACK_PPSUB     STR      STACK
     208      STACK_PSPHOT    STR      STACK
     209END
     210
     211# three pi stacks
     212THREEPI_STACK         METADATA
     213      STACK_PPSTACK   STR      STACK_THREEPI
     214      STACK_PPSUB     STR      STACK
     215      STACK_PSPHOT    STR      STACK
    194216END
    195217

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

    r25383 r30636  
    5050    PS_FITS_SCALE_STDEV_NEGATIVE,       ///< Auto-scale to sample stdev, place mean at upper limit
    5151    PS_FITS_SCALE_STDEV_BOTH,           ///< Auto-scale to sample stdev, place mean at middle
    52     PS_FITS_SCALE_MANUAL                ///< Manual scaling (use specified BSCALE and BZERO)
     52    PS_FITS_SCALE_MANUAL,                ///< Manual scaling (use specified BSCALE and BZERO)
     53    PS_FITS_SCALE_LOG_RANGE,            ///< Take logarithm, Auto-scale to preserve dynamic range
     54    PS_FITS_SCALE_LOG_STDEV_POSITIVE,   ///< Take logarithm, Auto-scale to sample stdev, place mean at lower limit
     55    PS_FITS_SCALE_LOG_STDEV_NEGATIVE,   ///< Take logarithm, Auto-scale to sample stdev, place mean at upper limit
     56    PS_FITS_SCALE_LOG_STDEV_BOTH,       ///< Take logarithm, Auto-scale to sample stdev, place mean at middle
     57    PS_FITS_SCALE_LOG_MANUAL            ///< Manual scaling (use specified BSCALE, BZERO, and BOFFSET)
    5358} psFitsScaling;
    5459
    … …  
    6671    bool fuzz;                          ///< Fuzz the values when quantising floating-point values?
    6772    double bscale, bzero;               ///< Manually specified BSCALE and BZERO (for SCALE_MANUAL)
     73    double boffset;                     ///< Manually specified BOFFSET (for SCALE_MANUAL)
    6874    double mean, stdev;                 ///< Mean and standard deviation of image
    6975    int stdevBits;                      ///< Number of bits to sample a standard deviation (for SCALE_STDEV_*)

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

    r29931 r30636  
    5454    int fitsDatatype;                   // cfitsio data type
    5555    int psDatatype;                     // psLib data type
     56    bool is_logscaled;                  // is this image log scaled using BOFFSET?
    5657} p_psFitsReadInfo;
    5758
    … …  
    128129
    129130    // Check scale and zero
    130     double bscale = 0.0, bzero = 0.0;    // Scale and zero point
     131    double bscale = 0.0, bzero = 0.0, boffset = NAN;    // Scale and zero point
    131132    if (fits_read_key_dbl(fits->fd, "BSCALE", &bscale, NULL, &status) && status != KEY_NO_EXIST) {
    132133        psFitsError(status, true, "Unable to read header.");
    … …  
    137138        psFitsError(status, true, "Unable to read header.");
    138139        goto bad;
     140    }
     141    status = 0;
     142    if (fits_read_key_dbl(fits->fd, "BOFFSET", &boffset, NULL, &status) && status != KEY_NO_EXIST) {
     143        psFitsError(status, true, "Unable to read header.");
     144        goto bad;
     145    }
     146    if (status == KEY_NO_EXIST) {
     147      info->is_logscaled = false;
     148    }
     149    else if (isfinite(boffset)) {
     150      info->is_logscaled = true;
     151    }
     152    else {
     153      info->is_logscaled = false;
    139154    }
    140155    status = 0;
    … …  
    246261static psImage *imageToDiskRepresentation(double *bscale, // Scaling applied
    247262                                          double *bzero, // Zero point applied
     263                                          double *boffset, // Log offset applied
    248264                                          long *blank, // Blank value (integer data)
    249265                                          psFitsFloat *floatType, // Type of custom floating-point
    … …  
    258274    psAssert(bscale, "impossible");
    259275    psAssert(bzero, "impossible");
     276    psAssert(boffset, "impossible");
    260277    psAssert(floatType, "impossible");
    261278    psAssert(fits, "impossible");
    … …  
    305322        if (newScaleZero) {
    306323            // Choose an appropriate BSCALE and BZERO
    307             if (!psFitsScaleDetermine(bscale, bzero, blank, image, mask, maskVal, fits)) {
     324          if (!psFitsScaleDetermine(bscale, bzero, boffset, blank, image, mask, maskVal, fits)) {
    308325                // We can't have the write dying for this reason --- try to save it somehow!
    309326                psWarning("Unable to determine BSCALE and BZERO for image --- refusing to quantise.");
    … …  
    331348        }
    332349
    333         return psFitsScaleForDisk(image, fits, *bscale, *bzero, rng);
     350        return psFitsScaleForDisk(image, fits, *bscale, *bzero, *boffset, rng);
    334351    }
    335352
    … …  
    459476        return NULL;
    460477    }
    461     psFree(info);
    462478
    463479    if (floatType != PS_FITS_FLOAT_NONE) {
    464480        outImage = psFitsFloatImageFromDisk(outImage, inImage, floatType);
    465481    }
     482
     483    // Need to apply BOFFSET if info->is_logscaled is true
     484    if (info->is_logscaled) {
     485      double boffset;
     486      int status;
     487      fits_read_key_dbl(fits->fd, "BOFFSET", &boffset, NULL, &status);
     488      outImage = psFitsScaleFromDisk(outImage,boffset);
     489    }
     490    psFree(info);
     491
    466492    psFree(inImage);
    467493
    … …  
    572598
    573599    double bscale = NAN, bzero = NAN;   // Scale and zero point to put in header
     600    double boffset = NAN;               // Log offset to put into header.
    574601    long blank = 0;                     // Blank (undefined) value for image
    575602    psFitsFloat floatType;              // Custom floating-point convention type
    576     psImage *diskImage = imageToDiskRepresentation(&bscale, &bzero, &blank, &floatType, fits, image,
     603    psImage *diskImage = imageToDiskRepresentation(&bscale, &bzero, &boffset, &blank, &floatType, fits, image,
    577604                                                   mask, maskVal, NULL, true); // Image to write out
    578605    if (!diskImage) {
    … …  
    610637
    611638    psFitsOptions *options = fits->options; // FITS I/O options
     639/*     if (options) { */
     640/*       if (options->scaling == PS_FITS_SCALE_LOG_RANGE) { */
     641/*      fprintf(stderr,"it has the scaling I expect\n"); */
     642/*       } */
     643/*       else { */
     644/*      fprintf(stderr,"it does nto have the scaling I expect\n"); */
     645/*       } */
     646/*     } */
     647/*     else { */
     648/*       fprintf(stderr,"options is null, apparently? \n"); */
     649/*     } */
     650     
    612651    psAssert(!useRequestedScale || !options || bitPix == options->bitpix || options->bitpix == 0,
    613652             "Something's not consistent");
    … …  
    651690    }
    652691
     692    // Remove any BOFFSET values that exist in the header if we are not using that scaling anymore
     693    if (options&&(!((options->scaling == PS_FITS_SCALE_LOG_RANGE)||
     694                    (options->scaling == PS_FITS_SCALE_LOG_MANUAL)||
     695                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_POSITIVE)||
     696                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_NEGATIVE)||
     697                   (options->scaling == PS_FITS_SCALE_LOG_STDEV_BOTH)))) {
     698      if (psMetadataLookup(header,"BOFFSET")) {
     699        psMetadataRemoveKey(header,"BOFFSET");
     700      }
     701    }   
     702
    653703    // write the header, if any.
    654704    if (header && !psFitsWriteHeaderImage(fits, header, createPHU)) {
    … …  
    668718        fits_write_key_dbl(fits->fd, "BSCALE", bscale, 12,
    669719                           "Scaling: TRUE = BZERO + BSCALE * DISK", &status);
     720        if (options&&(((options->scaling == PS_FITS_SCALE_LOG_RANGE)||
     721                       (options->scaling == PS_FITS_SCALE_LOG_MANUAL)||
     722                       (options->scaling == PS_FITS_SCALE_LOG_STDEV_POSITIVE)||
     723                       (options->scaling == PS_FITS_SCALE_LOG_STDEV_NEGATIVE)||
     724                       (options->scaling == PS_FITS_SCALE_LOG_STDEV_BOTH)))) {
     725          fits_write_key_dbl(fits->fd, "BOFFSET", boffset, 12,
     726                             "Scaling: TRUE = BZERO + BSCALE * 10**(DISK) + BOFFSET)", &status);
     727        }       
    670728        if (psFitsError(status, true, "Could not write BSCALE/BZERO headers to file.")) {
    671729            success = false;
    … …  
    769827    bool success = true;                // Successful update?
    770828    double bscale = NAN, bzero = NAN;   // Scale and zero point to put in header
     829    double boffset = NAN;               // Log offset to put in header
    771830    long blank = 0;                     // Blank (undefined) value for image
    772831    psFitsFloat floatType;              // Custom floating-point convention type
    773     psImage *diskImage = imageToDiskRepresentation(&bscale, &bzero, &blank, &floatType, fits, input,
     832    psImage *diskImage = imageToDiskRepresentation(&bscale, &bzero, &boffset, &blank, &floatType, fits, input,
    774833                                                   mask, maskVal, NULL, false); // Image to write out
    775834    if (!diskImage) {

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

    r27417 r30636  
    100100}
    101101
     102   
    102103// Determine appropriate BSCALE and BZERO for an image, mapping the standard deviation to the nominated number
    103104// of bits
    … …  
    212213    switch (options->scaling) {
    213214      case PS_FITS_SCALE_STDEV_POSITIVE:
     215      case PS_FITS_SCALE_LOG_STDEV_POSITIVE:
    214216        // Put (mean - N sigma) at the lowest possible value: predominantly positive images
    215217        imageVal = mean - options->stdevNum * stdev;
    … …  
    217219        break;
    218220      case PS_FITS_SCALE_STDEV_NEGATIVE:
     221      case PS_FITS_SCALE_LOG_STDEV_NEGATIVE:
    219222        // Put (mean + N sigma) at the highest possible value: predominantly negative images
    220223        imageVal = mean + options->stdevNum * stdev;
    … …  
    222225        break;
    223226      case PS_FITS_SCALE_STDEV_BOTH:
     227      case PS_FITS_SCALE_LOG_STDEV_BOTH:
    224228        // Put mean right in the middle: images with an equal abundance of positive and negative values
    225229        imageVal = mean;
    … …  
    237241}
    238242
     243
     244static bool logscaleStdev(double *bscale, // Scaling, to return
     245                          double *bzero, // Zero point, to return
     246                          double *boffset, // Log offset, to return
     247                          const psImage *image, // Image to scale
     248                          const psImage *mask, // Mask image
     249                          psImageMaskType maskVal, // Value to mask
     250                          const psFitsOptions *options // FITS options
     251    )
     252{
     253    psAssert(bscale, "impossible");
     254    psAssert(bzero, "impossible");
     255    psAssert(boffset, "impossible");
     256    psAssert(image, "impossible");
     257    psAssert(options, "impossible");
     258
     259    psTrace("psLib.fits", 3, "Scaling image by logarithm statistics");
     260    int numCols = image->numCols, numRows = image->numRows; // Size of image
     261
     262    psImage *copy;
     263   
     264    *boffset = 99e99;
     265
     266    // Make a copy of the image to pass to get the scaling parameters
     267
     268    switch (image->type.type) {
     269    case PS_TYPE_F32:
     270      copy = psImageCopy(NULL,image,PS_TYPE_F32);
     271      break;
     272    case PS_TYPE_F64:
     273      copy = psImageCopy(NULL,image,PS_TYPE_F64);
     274      break;
     275    default:
     276      psError(PS_ERR_UNKNOWN, true, "Target type is not a float: %d", image->type.type);
     277      return NULL;
     278      break;
     279    }
     280   
     281    // Determine the minimum value on this image.
     282    switch (image->type.type) {
     283    case PS_TYPE_F32:
     284      for (int y = 0; y < numRows; y++) {
     285        for (int x = 0; x < numCols; x++) {
     286          psF32 value = image->data.F32[y][x];
     287          if (isfinite(value)) {
     288            if (value < *boffset) {
     289              *boffset = value;
     290            }
     291          }
     292        }
     293      }
     294      break;
     295    case PS_TYPE_F64:
     296      for (int y = 0; y < numRows; y++) {
     297        for (int x = 0; x < numCols; x++) {
     298          psF64 value = image->data.F64[y][x];
     299          if (isfinite(value)) {
     300            if (value < *boffset) {
     301              *boffset = value;
     302            }
     303          }
     304        }
     305      }
     306      break;
     307    default:
     308      psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Target type is not a float: %d",image->type.type);
     309      return NULL;
     310      break;
     311    }
     312    // We only need to offset images that go negative.
     313    if (*boffset > 0.0) {
     314      *boffset = 0.0;
     315    }
     316    // Write offset to header
     317    // How?
     318    //    psMetadataAddF32(header,PS_LIST_TAIL,"LOGZERO",0,"Flux offset subtracted before taking logarithm.",offset);
     319    // Take the logarithm of the image, applying the offset
     320    switch (image->type.type) {
     321    case PS_TYPE_F32:
     322      for (int y = 0; y < numRows; y++) {
     323        for (int x = 0; x < numCols; x++) {
     324/*        if (x == 2331 && y == 2843) { */
     325/*          fprintf(stderr,"psFS32: %d %d %g %g %g\n",x,y,offset,image->data.F32[y][x],log10(image->data.F32[y][x] - offset)); */
     326/*        } */
     327          copy->data.F32[y][x] = (log10( image->data.F32[y][x] - *boffset));
     328        }
     329      }
     330      break;
     331    case PS_TYPE_F64:
     332        for (int y = 0; y < numRows; y++) {
     333          for (int x = 0; x < numCols; x++) {
     334            //      fprintf(stderr,"psFS64: %d %d %g %g %g\n",x,y,offset,image->data.F64[y][x],log10(image->data.F64[y][x] - offset));
     335            copy->data.F64[y][x] = (log10( image->data.F64[y][x] - *boffset));
     336          }
     337        }
     338      break;
     339    default:
     340      psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Target type is not a float: %d",image->type.type);
     341      return NULL;
     342      break;
     343    }
     344     
     345    // Do regular scaling on the logarithm image
     346    if (!scaleStdev(bscale, bzero, copy, mask, maskVal, options)) {
     347      psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from stdev");
     348      return false;
     349    }
     350    psFree(copy);
     351    return true;
     352}
     353
     354static bool logscaleRange(double *bscale, // Scaling, to return
     355                          double *bzero, // Zero point, to return
     356                          double *boffset, // Log offset, to return
     357                          const psImage *image, // Image to scale
     358                          const psFitsOptions *options // FITS options
     359    )
     360{
     361    psAssert(bscale, "impossible");
     362    psAssert(bzero, "impossible");
     363    psAssert(boffset, "impossible");
     364    psAssert(image, "impossible");
     365    psAssert(options, "impossible");
     366
     367    psTrace("psLib.fits", 3, "Scaling image by logarithm statistics");
     368    int numCols = image->numCols, numRows = image->numRows; // Size of image
     369
     370    psImage *copy;
     371   
     372    *boffset = 99e99;
     373
     374    // Make a copy of the image to pass to get the scaling parameters
     375
     376    switch (image->type.type) {
     377    case PS_TYPE_F32:
     378      copy = psImageCopy(NULL,image,PS_TYPE_F32);
     379      break;
     380    case PS_TYPE_F64:
     381      copy = psImageCopy(NULL,image,PS_TYPE_F64);
     382      break;
     383    default:
     384      psError(PS_ERR_UNKNOWN, true, "Target type is not a float: %d", image->type.type);
     385      return NULL;
     386      break;
     387    }
     388   
     389    // Determine the minimum value on this image.
     390    switch (image->type.type) {
     391    case PS_TYPE_F32:
     392      for (int y = 0; y < numRows; y++) {
     393        for (int x = 0; x < numCols; x++) {
     394          psF32 value = image->data.F32[y][x];
     395          if (!isfinite(value)) {
     396            if (value < *boffset) {
     397              *boffset = value;
     398            }
     399          }
     400        }
     401      }
     402      break;
     403    case PS_TYPE_F64:
     404      for (int y = 0; y < numRows; y++) {
     405        for (int x = 0; x < numCols; x++) {
     406          psF64 value = image->data.F64[y][x];
     407          if (!isfinite(value)) {
     408            if (value < *boffset) {
     409              *boffset = value;
     410            }
     411          }
     412        }
     413      }
     414      break;
     415    default:
     416      psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Target type is not a float: %d",image->type.type);
     417      return NULL;
     418      break;
     419    }
     420    // We only need to offset images that go negative.
     421    if (*boffset > 0.0) {
     422      *boffset = 0.0;
     423    }
     424    // Write offset to header
     425    // How?
     426    //    psMetadataAddF32(header,PS_LIST_TAIL,"LOGZERO",0,"Flux offset subtracted before taking logarithm.",offset);
     427    // Take the logarithm of the image, applying the offset
     428    switch (image->type.type) {
     429    case PS_TYPE_F32:
     430      for (int y = 0; y < numRows; y++) {
     431        for (int x = 0; x < numCols; x++) {
     432          copy->data.F32[y][x] = (log10( image->data.F32[y][x] - *boffset));
     433        }
     434      }
     435      break;
     436    case PS_TYPE_F64:
     437        for (int y = 0; y < numRows; y++) {
     438          for (int x = 0; x < numCols; x++) {
     439            copy->data.F64[y][x] = (log10( image->data.F64[y][x] - *boffset));
     440          }
     441        }
     442      break;
     443    default:
     444      psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Target type is not a float: %d",image->type.type);
     445      return NULL;
     446      break;
     447    }
     448     
     449    // Do regular scaling on the logarithm image
     450    if (!scaleRange(bscale, bzero, copy, options)) {
     451      psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from stdev");
     452      return false;
     453    }
     454    psFree(copy);
     455    return true;
     456}
     457
     458
    239459//////////////////////////////////////////////////////////////////////////////////////////////////////////////
    240460// Public functions
    241461//////////////////////////////////////////////////////////////////////////////////////////////////////////////
    242462
    243 bool psFitsScaleDetermine(double *bscale, double *bzero, long *blank, const psImage *image,
     463bool psFitsScaleDetermine(double *bscale, double *bzero, double *boffset, long *blank, const psImage *image,
    244464                          const psImage *mask, psImageMaskType maskVal, const psFits *fits)
    245465{
    246466    PS_ASSERT_PTR_NON_NULL(bscale, false);
    247467    PS_ASSERT_PTR_NON_NULL(bzero, false);
     468    PS_ASSERT_PTR_NON_NULL(boffset, false);
    248469    PS_ASSERT_PTR_NON_NULL(blank, false);
    249470    PS_ASSERT_IMAGE_NON_NULL(image, false);
    … …  
    256477    *bscale = NAN;
    257478    *bzero = NAN;
     479    *boffset = 0;
    258480    *blank = 0;
    259481
    … …  
    299521        }
    300522        break;
     523    case PS_FITS_SCALE_LOG_RANGE:
     524      if (!logscaleRange(bscale,bzero,boffset,image,options)) {
     525        psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from range");
     526        return false;
     527      }
     528      break;
     529    case PS_FITS_SCALE_LOG_STDEV_POSITIVE:
     530    case PS_FITS_SCALE_LOG_STDEV_NEGATIVE:
     531    case PS_FITS_SCALE_LOG_STDEV_BOTH:
     532      if (!logscaleStdev(bscale, bzero,boffset, image, mask, maskVal, options)) {
     533            psError(PS_ERR_UNKNOWN, false, "Unable to set BSCALE and BZERO from stdev");
     534            return false;
     535        }
     536        break;
    301537      case PS_FITS_SCALE_MANUAL:
    302538        *bscale = options->bscale;
    303539        *bzero = options->bzero;
     540        break;
     541      case PS_FITS_SCALE_LOG_MANUAL:
     542        *bscale = options->bscale;
     543        *bzero = options->bzero;
     544        *boffset = options->boffset;
    304545        break;
    305546      default:
    … …  
    316557
    317558
    318     psTrace("psLib.fits", 3, "BSCALE = %.10lf, BZERO = %.10lf, BLANK = %ld\n", *bscale, *bzero, *blank);
     559    psTrace("psLib.fits", 3, "BSCALE = %.10lf, BZERO = %.10lf, BOFFSET = %.10lf, BLANK = %ld\n",
     560            *bscale, *bzero, *boffset, *blank);
    319561    return true;
    320562}
    321563
    322564
    323 psImage *psFitsScaleForDisk(const psImage *image, const psFits *fits, double bscale, double bzero,
     565psImage *psFitsScaleForDisk(const psImage *image, const psFits *fits, double bscale, double bzero, double boffset,
    324566                            psRandom *rng)
    325567{
    … …  
    377619        for (int y = 0; y < numRows; y++) { \
    378620            for (int x = 0; x < numCols; x++) { \
    379                 ps##INTYPE value = (IN)->data.INTYPE[y][x]; \
     621              ps##INTYPE value; \
     622              if ((options->scaling == PS_FITS_SCALE_LOG_RANGE)||       \
     623                (options->scaling == PS_FITS_SCALE_LOG_MANUAL)||        \
     624                  (options->scaling == PS_FITS_SCALE_LOG_STDEV_POSITIVE)|| \
     625                  (options->scaling == PS_FITS_SCALE_LOG_STDEV_NEGATIVE)|| \
     626                  (options->scaling == PS_FITS_SCALE_LOG_STDEV_BOTH)) { \
     627                value = log10( (IN)->data.INTYPE[y][x] - boffset ); \
     628              }                                                         \
     629                else { \
     630                  value = (IN)->data.INTYPE[y][x];      \
     631                }                                           \
    380632                if (!isfinite(value)) { \
    381633                    /* This choice of "max" for non-finite pixels is mainly cosmetic --- it has to be */ \
    … …  
    423675
    424676
    425 #if 0
     677
    426678// This function to apply BSCALE and BZERO to an image read immediately from disk should not be necessary at
    427679// the present time, since cfitsio should apply the scaling itself in the process of reading.  However, we may
    428680// later desire it (e.g., if we ever make our own FITS implementation).
    429 psImage *psFitsScaleFromDisk(psFits *fits, psImage *image)
     681psImage *psFitsScaleFromDisk(const psImage *image, double boffset)
    430682{
    431683    PS_ASSERT_IMAGE_NON_NULL(image, NULL);
    432 
    433     if (bscale == 0.0) {
    434         // BSCALE = 0 means don't apply anything
    435         return psMemIncrRefCounter(image);
    436     }
    437684
    438685    psElemType inType = image->type.type; // Type for input image
    … …  
    471718      for (int y = 0; y < numRows; y++) { \
    472719          for (int x = 0; x < numCols; x++) { \
    473               out->data.OUTTYPE[y][x] = image->data.INTYPE[y][x] * bscale + bzero; \
     720            out->data.OUTTYPE[y][x] = pow(10,image->data.INTYPE[y][x]) + boffset;; \
    474721          } \
    475722      } \
    … …  
    505752    return out;
    506753}
    507 #endif
    508 
    509 
    510754
    511755psFitsScaling psFitsScalingFromString(const char *string)
    … …  
    517761    if (strcasecmp(string, "STDEV_NEGATIVE") == 0) return PS_FITS_SCALE_STDEV_NEGATIVE;
    518762    if (strcasecmp(string, "STDEV_BOTH") == 0)     return PS_FITS_SCALE_STDEV_BOTH;
     763    if (strcasecmp(string, "LOG_RANGE") == 0)      return PS_FITS_SCALE_LOG_RANGE;
     764    if (strcasecmp(string, "LOG_MANUAL") == 0)      return PS_FITS_SCALE_LOG_MANUAL;
     765    if (strcasecmp(string, "LOG_STDEV_POSITIVE") == 0) return PS_FITS_SCALE_LOG_STDEV_POSITIVE;
     766    if (strcasecmp(string, "LOG_STDEV_NEGATIVE") == 0) return PS_FITS_SCALE_LOG_STDEV_NEGATIVE;
     767    if (strcasecmp(string, "LOG_STDEV_BOTH") == 0) return PS_FITS_SCALE_LOG_STDEV_BOTH;
    519768    if (strcasecmp(string, "MANUAL") == 0)         return PS_FITS_SCALE_MANUAL;
    520769

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

    r21183 r30636  
    1010bool psFitsScaleDetermine(double *bscale, ///< Scaling, to return
    1111                          double *bzero, ///< Zero point, to return
     12                          double *boffset, ///< Log offset, to return
    1213                          long *blank,  ///< Blank value, to return
    1314                          const psImage *image, ///< Image to scale
    … …  
    2728                            double bscale, ///< Scaling
    2829                            double bzero, ///< Zero point
     30                            double boffset, ///< Log offset
    2931                            psRandom *rng ///< Random number generator (for the "fuzz"), or NULL
    3032    );
    31 
     33psImage *psFitsScaleFromDisk(const psImage *image, ///< Image to to unapply BOFFSET
     34                             double boffset        ///< Log offset
     35                             );
    3236/// Interpret a string as a scaling method
    3337psFitsScaling psFitsScalingFromString(const char *string ///< String to interpret

  • trunk/psModules/src/camera/pmFPAfileDefine.c

    r27657 r30636  
    282282              case PS_FITS_SCALE_NONE:
    283283              case PS_FITS_SCALE_RANGE:
     284              case PS_FITS_SCALE_LOG_RANGE:
    284285                // No options required
    285286                break;
    286287              case PS_FITS_SCALE_STDEV_POSITIVE:
    287288              case PS_FITS_SCALE_STDEV_NEGATIVE:
     289              case PS_FITS_SCALE_LOG_STDEV_POSITIVE:
     290              case PS_FITS_SCALE_LOG_STDEV_NEGATIVE:
    288291                options->stdevNum = parseOptionFloat(scheme, "STDEV.NUM", source); // Padding to edge
    289292                if (!isfinite(options->stdevNum)) {
    … …  
    295298                // Flow through
    296299              case PS_FITS_SCALE_STDEV_BOTH:
     300              case PS_FITS_SCALE_LOG_STDEV_BOTH:
    297301                options->stdevBits = parseOptionInt(scheme, "STDEV.BITS", source, 0); // Bits for stdev
    298302                if (options->stdevBits <= 0) {
    … …  
    308312                options->bzero = parseOptionDouble(scheme, "BZERO", source); // Zero point
    309313                break;
     314            case PS_FITS_SCALE_LOG_MANUAL:
     315              options->bscale = parseOptionDouble(scheme, "BSCALE", source); // Scaling
     316              options->bzero = parseOptionDouble(scheme, "BZERO", source); // Zero point
     317              options->boffset = parseOptionDouble(scheme, "BOFFSET", source); // Log offset
     318              break;         
    310319              default:
    311320                psAbort("Should never get here.");

  • trunk/psModules/src/config/pmConfig.c

    r29004 r30636  
    897897            psMetadataAddMetadata(filerules, PS_LIST_TAIL, old, PS_META_REPLACE,
    898898                                  "Original replaced by -F option", newRule);
     899        }
     900        psFree(camerasIter);
     901    }
     902
     903    // Look for command-line options for files to replace
     904    while ((argNum = psArgumentGet(*argc, argv, "-R")) > 0) {
     905        psArgumentRemove(argNum, argc, argv);
     906        if (argNum + 2 >= *argc) {
     907            psError(PM_ERR_CONFIG, true,
     908                    "Filerule element switch (-R) provided without filerule element and value.");
     909            psFree(config);
     910            return NULL;
     911        }
     912
     913        const char *rulename = argv[argNum]; // The filerule, to be modified
     914        psArgumentRemove(argNum, argc, argv);
     915        const char *element  = argv[argNum]; // The element, to be modified
     916        psArgumentRemove(argNum, argc, argv);
     917        const char *value    = argv[argNum]; // The value, to be set
     918        psArgumentRemove(argNum, argc, argv);
     919
     920        psMetadata *cameras = psMetadataLookupMetadata(NULL, config->system, "CAMERAS"); // List of cameras
     921        if (!cameras) {
     922            psError(PM_ERR_CONFIG, false, "Unable to find CAMERAS in the site configuration.\n");
     923            return false;
     924        }
     925
     926        psMetadataIterator *camerasIter = psMetadataIteratorAlloc(cameras, PS_LIST_HEAD, NULL); // Iterator
     927        psMetadataItem *cameraItem;     // Item from iteration
     928        while ((cameraItem = psMetadataGetAndIncrement(camerasIter))) {
     929            // Silently ignore problems --- they will be caught later, because if the user wants the nominated
     930            // file and it's not available for that camera, then they will know.
     931
     932            if (cameraItem->type != PS_DATA_METADATA) {
     933                psTrace("psModules.config", 2,
     934                        "Entry %s in CAMERAS is not of type METADATA --- ignored.", cameraItem->name);
     935                continue;
     936            }
     937            psMetadata *camera = cameraItem->data.md; // Camera configuration
     938
     939            psMetadata *newRule = pmConfigFileRule(config, camera, rulename); // The rule of interest
     940            if (!newRule) {
     941                psTrace("psModules.config", 2,
     942                        "Unable to find filerule %s in camera %s --- ignored.", rulename, cameraItem->name);
     943                continue;
     944            }
     945
     946            // By calling pmConfigFileRule, we've assured that the FILERULES is now a metadata
     947            psMetadata *filerules = psMetadataLookupMetadata(NULL, camera, "FILERULES"); // File rules
     948            if (!filerules) {
     949                psTrace("psModules.config", 2,
     950                        "Can't find FILERULES of type METADATA in camera %s --- ignored.", cameraItem->name);
     951                continue;
     952            }
     953
     954            // Convert newRule to have the element value requested.
     955            if (!psMetadataLookupStr(NULL,newRule,element)) {
     956              psTrace("psModules.config", 2,
     957                      "Unable to find filerule element %s in filerule %s in camera %s --- ignored.",
     958                      element,rulename,cameraItem->name);
     959              continue;
     960            }
     961            psMetadataAddStr(newRule, PS_LIST_TAIL, element, PS_META_REPLACE,
     962                             "Original replaced by -R option", value);
     963           
     964            psMetadataAddMetadata(filerules, PS_LIST_TAIL, rulename, PS_META_REPLACE,
     965                                  "Original replaced by -R option", newRule);
    899966        }
    900967        psFree(camerasIter);
Note: See TracChangeset for help on using the changeset viewer.