IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 27402


Ignore:
Timestamp:
Mar 23, 2010, 9:12:53 AM (16 years ago)
Author:
Paul Price
Message:

Pulling r27400 out of the trunk until it can be developed further (compression, make optional). Development will continue on branches/pap_stack/

Location:
trunk
Files:
23 edited

Legend:

Unmodified
Added
Removed
  • trunk/ippconfig/recipes/filerules-mef.mdc

    r27400 r27402  
    265265PPSTACK.OUTPUT.MASK     OUTPUT {OUTPUT}.mk.fits                  MASK      COMP_MASK  FPA        TRUE      NONE
    266266PPSTACK.OUTPUT.VARIANCE OUTPUT {OUTPUT}.wt.fits                  VARIANCE  COMP_WT    FPA        TRUE      NONE
    267 PPSTACK.OUTPUT.EXP      OUTPUT {OUTPUT}.exp.fits                 IMAGE     NONE       FPA        TRUE      NONE
    268 PPSTACK.OUTPUT.EXPNUM   OUTPUT {OUTPUT}.num.fits                 MASK      NONE       FPA        TRUE      NONE
    269267PPSTACK.UNCONV          OUTPUT {OUTPUT}.unconv.fits              IMAGE     COMP_IMG   FPA        TRUE      NONE
    270268PPSTACK.UNCONV.MASK     OUTPUT {OUTPUT}.unconv.mask.fits         MASK      COMP_MASK  FPA        TRUE      NONE
    271269PPSTACK.UNCONV.VARIANCE OUTPUT {OUTPUT}.unconv.wt.fits           VARIANCE  COMP_WT    FPA        TRUE      NONE
    272 PPSTACK.UNCONV.EXP      OUTPUT {OUTPUT}.unconv.exp.fits          IMAGE     NONE       FPA        TRUE      NONE
    273 PPSTACK.UNCONV.EXPNUM   OUTPUT {OUTPUT}.unconv.num.fits          MASK      NONE       FPA        TRUE      NONE
    274270PPSTACK.TARGET.PSF      OUTPUT {OUTPUT}.target.psf               PSF       NONE       CHIP       TRUE      NONE
    275271PPSTACK.CONV.KERNEL     OUTPUT {OUTPUT}.{FILE.INDEX}.kernel      SUBKERNEL NONE       FPA        TRUE      NONE
  • trunk/ippconfig/recipes/filerules-simple.mdc

    r27400 r27402  
    214214PPSTACK.OUTPUT.MASK     OUTPUT {OUTPUT}.mask.fits    MASK      NONE       FPA        TRUE      NONE
    215215PPSTACK.OUTPUT.VARIANCE OUTPUT {OUTPUT}.weight.fits  VARIANCE  NONE       FPA        TRUE      NONE
    216 PPSTACK.OUTPUT.EXP      OUTPUT {OUTPUT}.exp.fits     IMAGE     NONE       FPA        TRUE      NONE
    217 PPSTACK.OUTPUT.EXPNUM   OUTPUT {OUTPUT}.num.fits     MASK      NONE       FPA        TRUE      NONE
    218216PPSTACK.UNCONV          OUTPUT {OUTPUT}.unconv.fits  IMAGE     COMP_IMG   FPA        TRUE      NONE
    219217PPSTACK.UNCONV.MASK     OUTPUT {OUTPUT}.unconv.mask.fits MASK  COMP_MASK  FPA        TRUE      NONE
    220218PPSTACK.UNCONV.VARIANCE OUTPUT {OUTPUT}.unconv.wt.fits VARIANCE COMP_WT   FPA        TRUE      NONE
    221 PPSTACK.UNCONV.EXP      OUTPUT {OUTPUT}.unconv.exp.fits IMAGE  NONE       FPA        TRUE      NONE
    222 PPSTACK.UNCONV.EXPNUM   OUTPUT {OUTPUT}.unconv.num.fits MASK   NONE       FPA        TRUE      NONE
    223219PPSTACK.TARGET.PSF      OUTPUT {OUTPUT}.target.psf   PSF       NONE       CHIP       TRUE      NONE
    224220PPSTACK.CONV.KERNEL     OUTPUT {OUTPUT}.{FILE.INDEX}.kernel SUBKERNEL NONE FPA       TRUE      NONE
  • trunk/ippconfig/recipes/filerules-split.mdc

    r27400 r27402  
    233233PPSTACK.OUTPUT.MASK     OUTPUT {OUTPUT}.mask.fits                MASK      COMP_MASK  FPA        TRUE      NONE
    234234PPSTACK.OUTPUT.VARIANCE OUTPUT {OUTPUT}.wt.fits                  VARIANCE  COMP_WT    FPA        TRUE      NONE
    235 PPSTACK.OUTPUT.EXP      OUTPUT {OUTPUT}.exp.fits                 IMAGE     NONE       FPA        TRUE      NONE
    236 PPSTACK.OUTPUT.EXPNUM   OUTPUT {OUTPUT}.num.fits                 MASK      NONE       FPA        TRUE      NONE
    237235PPSTACK.UNCONV          OUTPUT {OUTPUT}.unconv.fits              IMAGE     COMP_IMG   FPA        TRUE      NONE
    238236PPSTACK.UNCONV.MASK     OUTPUT {OUTPUT}.unconv.mask.fits         MASK      COMP_MASK  FPA        TRUE      NONE
    239237PPSTACK.UNCONV.VARIANCE OUTPUT {OUTPUT}.unconv.wt.fits           VARIANCE  COMP_WT    FPA        TRUE      NONE
    240 PPSTACK.UNCONV.EXP      OUTPUT {OUTPUT}.unconv.exp.fits          IMAGE     NONE       FPA        TRUE      NONE
    241 PPSTACK.UNCONV.EXPNUM   OUTPUT {OUTPUT}.unconv.num.fits          MASK      NONE       FPA        TRUE      NONE
    242238PPSTACK.TARGET.PSF      OUTPUT {OUTPUT}.target.psf               PSF       NONE       CHIP       TRUE      NONE
    243239PPSTACK.CONV.KERNEL     OUTPUT {OUTPUT}.{FILE.INDEX}.kernel      SUBKERNEL NONE       FPA        TRUE      NONE
  • trunk/ppStack/src/ppStack.h

    r27400 r27402  
    6969                               const psVector *mask, // Mask for input readouts
    7070                               const psVector *weightings, // Weighting factors for each image
    71                                const psVector *exposures,  // Exposure time for each image
    7271                               const psVector *addVariance // Additional variance for rejection
    7372    );
     
    8483bool ppStackReadoutFinal(const pmConfig *config,   // Configuration
    8584                         pmReadout *outRO,   // Output readout
    86                          pmReadout *expRO,   // Exposure readout
    8785                         const psArray *readouts, // Input readouts
    8886                         const psVector *mask, // Mask for input readouts
    8987                         const psArray *rejected, // Array with pixels rejected in each image
    9088                         const psVector *weightings, // Weighting factors for each image
    91                          const psVector *exposures,  // Exposure times for each image
    9289                         const psVector *addVariance, // Additional variance for rejection
    9390                         bool safety,                 // Enable safety switch?
  • trunk/ppStack/src/ppStackCamera.c

    r27400 r27402  
    283283    }
    284284
    285 
    286     // Exposure image
    287     pmFPA *expFPA = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain the output
    288     if (!expFPA) {
    289         psError(psErrorCodeLast(), false, "Unable to construct an FPA from camera configuration.");
    290         return false;
    291     }
    292     pmFPAfile *exp = pmFPAfileDefineOutput(config, expFPA, "PPSTACK.OUTPUT.EXP");
    293     psFree(expFPA);                        // Drop reference
    294     if (!exp) {
    295         psError(psErrorCodeLast(), false, _("Unable to generate output file from PPSTACK.OUTPUT.EXP"));
    296         return false;
    297     }
    298     if (exp->type != PM_FPA_FILE_IMAGE) {
    299         psError(PPSTACK_ERR_CONFIG, true, "PPSTACK.OUTPUT.EXP is not of type IMAGE");
    300         return false;
    301     }
    302     exp->save = true;
    303 
    304     if (!pmFPAAddSourceFromFormat(expFPA, "Stack", exp->format)) {
    305         psError(psErrorCodeLast(), false, "Unable to generate output FPA.");
    306         return false;
    307     }
    308 
    309     // Exposure numbers
    310     pmFPAfile *expNum = pmFPAfileDefineOutput(config, exp->fpa, "PPSTACK.OUTPUT.EXPNUM");
    311     if (!expNum) {
    312         psError(psErrorCodeLast(), false, _("Unable to generate output file from PPSTACK.OUTPUT.MASK"));
    313         return false;
    314     }
    315     if (expNum->type != PM_FPA_FILE_MASK) {
    316         psError(PPSTACK_ERR_CONFIG, true, "PPSTACK.OUTPUT.EXPNUM is not of type MASK");
    317         return false;
    318     }
    319     expNum->save = true;
    320 
    321 
    322 
    323285    if (havePSFs) {
    324286        pmFPA *psfFPA = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain PSF
     
    340302    }
    341303
     304#if 1
    342305    // Unconvolved stack
    343306    pmFPA *unconvFPA = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain unconvolved output
     
    388351        unconvVariance->save = true;
    389352    }
    390 
    391 
    392     // Exposure image
    393     pmFPA *unconvExpFPA = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain the output
    394     if (!unconvExpFPA) {
    395         psError(psErrorCodeLast(), false, "Unable to construct an FPA from camera configuration.");
    396         return false;
    397     }
    398     pmFPAfile *unconvExp = pmFPAfileDefineOutput(config, unconvExpFPA, "PPSTACK.UNCONV.EXP");
    399     psFree(unconvExpFPA);               // Drop reference
    400     if (!unconvExp) {
    401         psError(psErrorCodeLast(), false, _("Unable to generate output file from PPSTACK.UNCONV.EXP"));
    402         return false;
    403     }
    404     if (unconvExp->type != PM_FPA_FILE_IMAGE) {
    405         psError(PPSTACK_ERR_CONFIG, true, "PPSTACK.UNCONV.EXP is not of type IMAGE");
    406         return false;
    407     }
    408     unconvExp->save = true;
    409 
    410     if (!pmFPAAddSourceFromFormat(unconvExpFPA, "Stack", unconvExp->format)) {
    411         psError(psErrorCodeLast(), false, "Unable to generate output FPA.");
    412         return false;
    413     }
    414 
    415     // Output mask
    416     pmFPAfile *unconvExpNum = pmFPAfileDefineOutput(config, unconvExp->fpa, "PPSTACK.UNCONV.EXPNUM");
    417     if (!unconvExpNum) {
    418         psError(psErrorCodeLast(), false, _("Unable to generate output file from PPSTACK.UNCONV.MASK"));
    419         return false;
    420     }
    421     if (unconvExpNum->type != PM_FPA_FILE_MASK) {
    422         psError(PPSTACK_ERR_CONFIG, true, "PPSTACK.UNCONV.EXPNUM is not of type MASK");
    423         return false;
    424     }
    425     unconvExpNum->save = true;
    426 
     353#endif
    427354
    428355    // Output JPEGs
  • trunk/ppStack/src/ppStackCleanup.c

    r27400 r27402  
    8686        options->outRO = NULL;
    8787
    88         options->expRO->data_exists = false;
    89         options->expRO->parent->data_exists = false;
    90         options->expRO->parent->parent->data_exists = false;
    91         psFree(options->expRO);
    92         options->expRO = NULL;
    93 
    9488        pmFPAview *view = pmFPAviewAlloc(0);// Pointer into FPA hierarchy
    9589        view->chip = view->cell = 0;        // pmFPAviewFreeData doesn't want to deal with readouts
  • trunk/ppStack/src/ppStackCombineFinal.c

    r27400 r27402  
    1414//#define TESTING                         // Enable test output
    1515
    16 bool ppStackCombineFinal(ppStackThreadData *stack, psArray *covariances, ppStackOptions *options,
    17                          pmConfig *config, bool safe, bool normalise, bool grow)
     16bool ppStackCombineFinal(pmReadout *target, ppStackThreadData *stack, psArray *covariances,
     17                         ppStackOptions *options, pmConfig *config, bool safe, bool normalise, bool grow)
    1818{
    1919    psAssert(stack, "Require stack");
     
    2121    psAssert(config, "Require configuration");
    2222
    23     pmReadout *outRO = options->outRO;                                      // Output readout
    24     pmReadout *expRO = options->expRO;                                      // Exposure readout
    25     int numCols = outRO->image->numCols, numRows = outRO->image->numRows; // Size of image
     23    int numCols = target->image->numCols, numRows = target->image->numRows; // Size of image
    2624
    2725    psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSTACK_RECIPE); // ppStack recipe
     
    4543    }
    4644
    47     if (!outRO->mask) {
    48         outRO->mask = psImageAlloc(numCols, numRows, PS_TYPE_IMAGE_MASK);
    49     }
    50     if (!expRO->mask) {
    51         expRO->mask = psImageAlloc(numCols, numRows, PS_TYPE_IMAGE_MASK);
     45    if (!target->mask) {
     46        target->mask = psImageAlloc(target->image->numCols, target->image->numRows, PS_TYPE_IMAGE_MASK);
    5247    }
    5348
     
    6762        }
    6863
    69         // call: ppStackReadoutFinal(config, outRO, readouts, rejected)
     64        // call: ppStackReadoutFinal(config, target, readouts, rejected)
    7065        psThreadJob *job = psThreadJobAlloc("PPSTACK_FINAL_COMBINE"); // Job to start
     66        psArrayAdd(job->args, 1, target);
    7167        psArrayAdd(job->args, 1, thread);
    7268        psArrayAdd(job->args, 1, reject);
     
    109105        }
    110106        if (sumWeights > 0.0) {
    111             outRO->covariance = psImageCovarianceSum(covariances);
    112             psBinaryOp(outRO->covariance->image, outRO->covariance->image, "/",
     107            target->covariance = psImageCovarianceSum(covariances);
     108            psBinaryOp(target->covariance->image, target->covariance->image, "/",
    113109                       psScalarAlloc(sumWeights, PS_TYPE_F32));
    114             psImageCovarianceTransfer(outRO->variance, outRO->covariance);
     110            psImageCovarianceTransfer(target->variance, target->covariance);
    115111        }
    116112    } else {
    117         outRO->covariance = psImageCovarianceNone();
     113        target->covariance = psImageCovarianceNone();
    118114    }
    119115
     
    131127            wcsDone = true;
    132128            pmHDU *inHDU = pmHDUFromCell(inRO->parent); // Template HDU
    133             pmHDU *outHDU = pmHDUFromCell(outRO->parent); // Output HDU
    134             pmChip *outChip = outRO->parent->parent; // Output chip
     129            pmHDU *outHDU = pmHDUFromCell(target->parent); // Output HDU
     130            pmChip *outChip = target->parent->parent; // Output chip
    135131            pmFPA *outFPA = outChip->parent; // Output FPA
    136132            if (!outHDU || !inHDU) {
     
    155151    }
    156152
    157     // Set exposure time correctly
    158     {
    159         float exptime = 0.0;            // Summed exposure time
    160         for (int i = 0; i < options->num; i++) {
    161             if (options->inputMask) {
    162                 continue;
    163             }
    164             exptime += options->exposures->data.F32[i];
    165         }
    166 
    167         {
    168             psMetadataItem *item = psMetadataLookup(outRO->parent->concepts, "CELL.EXPOSURE");
    169             item->data.F32 = exptime;
    170         }
    171         {
    172             psMetadataItem *item = psMetadataLookup(outRO->parent->parent->parent->concepts, "FPA.EXPOSURE");
    173             item->data.F32 = exptime;
    174         }
    175     }
    176 
    177153    // Put version information into the header
    178     pmHDU *hdu = pmHDUFromCell(outRO->parent);
     154    pmHDU *hdu = pmHDUFromCell(target->parent);
    179155    if (!hdu) {
    180156        psError(PPSTACK_ERR_PROG, false, "Unable to find HDU for output.");
     
    192168    psStringAppend(&name, "combined_image_final_%d.fits", pass);
    193169    pass++;
    194     ppStackWriteImage(name, NULL, outRO->image, config);
     170    ppStackWriteImage(name, NULL, target->image, config);
    195171    psStringSubstitute(&name, "mask", "image");
    196     ppStackWriteImage(name, NULL, outRO->mask, config);
     172    ppStackWriteImage(name, NULL, target->mask, config);
    197173    psStringSubstitute(&name, "variance", "mask");
    198     ppStackWriteImage(name, NULL, outRO->variance, config);
     174    ppStackWriteImage(name, NULL, target->variance, config);
    199175    psFree(name);
    200176
    201     pmStackVisualPlotTestImage(outRO->image, "combined_image_final.fits");
     177    pmStackVisualPlotTestImage(target->image, "combined_image_final.fits");
    202178#endif
    203179
  • trunk/ppStack/src/ppStackCombineInitial.c

    r27400 r27402  
    8989    }
    9090
    91     ppStackMemDump("initial");
    92 
    9391#ifdef TESTING
    9492    ppStackWriteImage("combined_image_initial.fits", NULL, options->outRO->image, config);
  • trunk/ppStack/src/ppStackCombinePrepare.c

    r27400 r27402  
    1010#include "ppStackLoop.h"
    1111
    12 bool ppStackCombinePrepare(const char *outName, const char *expName,
    13                            ppStackFileList files, ppStackThreadData *stack,
     12bool ppStackCombinePrepare(pmReadout **ro, const char *name, ppStackFileList files, ppStackThreadData *stack,
    1413                           ppStackOptions *options, pmConfig *config)
    1514{
     
    3332    }
    3433
    35     pmCell *cell = pmFPAfileThisCell(config->files, view, outName); // Output cell
    36     options->outRO = pmReadoutAlloc(cell); // Output readout
    37 
    38     pmCell *expCell = pmFPAfileThisCell(config->files, view, expName); // Exposure cell
    39     options->expRO = pmReadoutAlloc(expCell); // Output readout
     34    pmCell *cell = pmFPAfileThisCell(config->files, view, name); // Output cell
     35    *ro = pmReadoutAlloc(cell); // Output readout
    4036
    4137    psFree(view);
     
    4642    psImageMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels
    4743
    48     if (!pmReadoutStackDefineOutput(options->outRO, col0, row0, numCols, numRows, true, true, maskBad)) {
    49         psError(PPSTACK_ERR_ARGUMENTS, false, "Unable to prepare output.");
    50         return false;
    51     }
    52 
    53     if (!pmReadoutStackDefineOutput(options->expRO, col0, row0, numCols, numRows, true, true, 0)) {
     44    if (!pmReadoutStackDefineOutput(*ro, col0, row0, numCols, numRows, true, true, maskBad)) {
    5445        psError(PPSTACK_ERR_ARGUMENTS, false, "Unable to prepare output.");
    5546        return false;
  • trunk/ppStack/src/ppStackFiles.c

    r27400 r27402  
    2323/// Regular (convolved) stack files
    2424static char *filesStack[] = { "PPSTACK.OUTPUT", "PPSTACK.OUTPUT.MASK", "PPSTACK.OUTPUT.VARIANCE",
    25                               "PPSTACK.OUTPUT.EXP", "PPSTACK.OUTPUT.EXPNUM",
    26                               "PPSTACK.OUTPUT.JPEG1", "PPSTACK.OUTPUT.JPEG2",
    27                               NULL };
     25                              "PPSTACK.OUTPUT.JPEG1", "PPSTACK.OUTPUT.JPEG2", NULL };
    2826/// Unconvolved stack files
    29 static char *filesUnconv[] = { "PPSTACK.UNCONV", "PPSTACK.UNCONV.MASK", "PPSTACK.UNCONV.VARIANCE",
    30                                "PPSTACK.UNCONV.EXP", "PPSTACK.UNCONV.EXPNUM", NULL };
     27static char *filesUnconv[] = { "PPSTACK.UNCONV", "PPSTACK.UNCONV.MASK", "PPSTACK.UNCONV.VARIANCE", NULL };
    3128
    3229/// Files for photometry
  • trunk/ppStack/src/ppStackLoop.c

    r27400 r27402  
    111111
    112112    // Prepare for combination
    113     if (!ppStackCombinePrepare("PPSTACK.OUTPUT", "PPSTACK.OUTPUT.EXP", PPSTACK_FILES_STACK,
     113    if (!ppStackCombinePrepare(&options->outRO, "PPSTACK.OUTPUT", PPSTACK_FILES_STACK,
    114114                               stack, options, config)) {
    115115        psError(psErrorCodeLast(), false, "Unable to prepare for combination.");
     
    160160    // Final combination
    161161    psTrace("ppStack", 2, "Final stack of convolved images....\n");
    162     if (!ppStackCombineFinal(stack, options->convCovars, options, config, false, false, true)) {
     162    if (!ppStackCombineFinal(options->outRO, stack, options->convCovars, options, config,
     163                             false, false, true)) {
    163164        psError(psErrorCodeLast(), false, "Unable to perform final combination.");
    164165        psFree(stack);
     
    202203
    203204        // Prepare for combination
    204         if (!ppStackCombinePrepare("PPSTACK.UNCONV", "PPSTACK.UNCONV.EXP", PPSTACK_FILES_UNCONV,
     205        if (!ppStackCombinePrepare(&options->unconvRO, "PPSTACK.UNCONV", PPSTACK_FILES_UNCONV,
    205206                                   stack, options, config)) {
    206207            psError(psErrorCodeLast(), false, "Unable to prepare for combination.");
     
    210211
    211212        psTrace("ppStack", 2, "Stack of unconvolved images....\n");
    212         if (!ppStackCombineFinal(stack, options->origCovars, options, config,
     213        if (!ppStackCombineFinal(options->unconvRO, stack, options->origCovars, options, config,
    213214                                 false, true, false)) {
    214215            psError(psErrorCodeLast(), false, "Unable to perform unconvolved combination.");
     
    225226        }
    226227        ppStackFileActivation(config, PPSTACK_FILES_UNCONV, false);
    227         options->outRO->data_exists = false;
    228         options->outRO->parent->data_exists = false;
    229         options->outRO->parent->parent->data_exists = false;
    230         psFree(options->outRO);
    231         options->outRO = NULL;
    232         options->expRO->data_exists = false;
    233         options->expRO->parent->data_exists = false;
    234         options->expRO->parent->parent->data_exists = false;
    235         psFree(options->expRO);
    236         options->expRO = NULL;
     228        options->unconvRO->data_exists = false;
     229        options->unconvRO->parent->data_exists = false;
     230        options->unconvRO->parent->parent->data_exists = false;
     231        psFree(options->unconvRO);
     232        options->unconvRO = NULL;
    237233
    238234        for (int i = 0; i < options->num; i++) {
  • trunk/ppStack/src/ppStackLoop.h

    r27400 r27402  
    3838// Prepare for combination
    3939bool ppStackCombinePrepare(
    40     const char *outName,                // Name of output file
    41     const char *expName,                // Name of exposure file
     40    pmReadout **readout,                // Readout to set
     41    const char *name,                   // Name of file
    4242    ppStackFileList files,              // Files of interest
    4343    ppStackThreadData *stack,           // Stack
     
    6161// Final combination
    6262bool ppStackCombineFinal(
     63    pmReadout *target,                  // Target readout
    6364    ppStackThreadData *stack,           // Stack
    6465    psArray *covariances,               // Covariances
  • trunk/ppStack/src/ppStackOptions.c

    r27400 r27402  
    2121    psFree(options->psf);
    2222    psFree(options->inputSeeing);
    23     psFree(options->exposures);
    2423    psFree(options->inputMask);
    2524    psFree(options->sourceLists);
     
    3332    psFree(options->convCovars);
    3433    psFree(options->outRO);
    35     psFree(options->expRO);
     34    psFree(options->unconvRO);
    3635    psFree(options->inspect);
    3736    psFree(options->rejected);
     
    6261    options->zp = NAN;
    6362    options->inputSeeing = NULL;
    64     options->exposures = NULL;
    6563    options->targetSeeing = NAN;
    6664    options->inputMask = NULL;
     
    7775    options->convCovars = NULL;
    7876    options->outRO = NULL;
    79     options->expRO = NULL;
     77    options->unconvRO = NULL;
    8078    options->inspect = NULL;
    8179    options->rejected = NULL;
  • trunk/ppStack/src/ppStackOptions.h

    r27400 r27402  
    2222    psVector *inputSeeing;              // Input seeing FWHMs
    2323    float targetSeeing;                 // Target seeing FWHM
    24     psVector *exposures;                // Exposure times
    2524    float sumExposure;                  // Sum of exposure times
    2625    float zp;                           // Zero point for output
     
    3938    // Combine initial
    4039    pmReadout *outRO;                   // Output readout
    41     pmReadout *expRO;                   // Exposure readout
     40    pmReadout *unconvRO;                // Unconvolved readout
    4241    psArray *inspect;                   // Array of arrays of pixels to inspect
    4342    // Rejection
  • trunk/ppStack/src/ppStackPrepare.c

    r27400 r27402  
    126126    options->inputMask = psVectorAlloc(num, PS_TYPE_VECTOR_MASK); // Mask for inputs
    127127    psVectorInit(options->inputMask, 0);
    128     options->exposures = psVectorAlloc(options->num, PS_TYPE_F32);
    129     psVectorInit(options->exposures, NAN);
    130128
    131129    pmFPAfileActivate(config->files, false, NULL);
     
    136134    }
    137135
     136    psMetadataIterator *fileIter = psMetadataIteratorAlloc(config->files, PS_LIST_HEAD, "^PPSTACK.INPUT$");
     137    psMetadataItem *fileItem; // Item from iteration
    138138    psArray *psfs = psArrayAlloc(num); // PSFs for PSF envelope
    139139    int numCols = 0, numRows = 0;   // Size of image
    140     options->sumExposure = 0.0;
    141     for (int i = 0; i < num; i++) {
    142         pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.INPUT", i); // File of interest
    143         pmCell *cell = pmFPAviewThisCell(view, file->fpa); // Cell of interest
    144 
    145         options->exposures->data.F32[i] = psMetadataLookupF32(NULL, cell->concepts, "CELL.EXPOSURE");
    146         options->sumExposure += options->exposures->data.F32[i];
     140    int index = 0;                  // Index for file
     141    while ((fileItem = psMetadataGetAndIncrement(fileIter))) {
     142        assert(fileItem->type == PS_DATA_UNKNOWN);
     143        pmFPAfile *inputFile = fileItem->data.V; // An input file
    147144
    148145        // Get list of PSFs, to determine target PSF
    149146        if (options->convolve) {
    150             pmChip *chip = pmFPAviewThisChip(view, file->fpa); // The chip: holds the PSF
     147            pmChip *chip = pmFPAviewThisChip(view, inputFile->fpa); // The chip: holds the PSF
    151148            pmPSF *psf = psMetadataLookupPtr(NULL, chip->analysis, "PSPHOT.PSF"); // PSF
    152149            if (!psf) {
    153150                psError(PPSTACK_ERR_PROG, false, "Unable to find PSF.");
    154151                psFree(view);
     152                psFree(fileIter);
    155153                psFree(psfs);
    156154                return false;
    157155            }
    158             psfs->data[i] = psMemIncrRefCounter(psf);
     156            psfs->data[index] = psMemIncrRefCounter(psf);
    159157            psMetadataRemoveKey(chip->analysis, "PSPHOT.PSF");
    160158
    161             pmCell *cell = pmFPAviewThisCell(view, file->fpa); // Cell of interest
     159            pmCell *cell = pmFPAviewThisCell(view, inputFile->fpa); // Cell of interest
    162160            pmHDU *hdu = pmHDUFromCell(cell);
    163161            assert(hdu && hdu->header);
     
    167165                psError(PPSTACK_ERR_PROG, false, "Unable to determine size of image from PSF.");
    168166                psFree(view);
     167                psFree(fileIter);
    169168                psFree(psfs);
    170169                return false;
     
    181180        pmDetections *detections = NULL;
    182181        if (options->convolve || options->matchZPs || options->photometry || redoPhot) {
    183             pmReadout *ro = pmFPAviewThisReadout(view, file->fpa); // Readout with sources
     182            pmReadout *ro = pmFPAviewThisReadout(view, inputFile->fpa); // Readout with sources
    184183            detections = psMetadataLookupPtr(NULL, ro->analysis, "PSPHOT.DETECTIONS"); // Sources
    185184            if (!detections) {
     
    189188            psAssert (detections->allSources, "missing sources?");
    190189
    191             options->sourceLists->data[i] = psMemIncrRefCounter(detections->allSources);
     190            options->sourceLists->data[index] = psMemIncrRefCounter(detections->allSources);
    192191        }
    193192
    194193        // Re-do photometry if we don't trust the source lists
    195194        if (redoPhot) {
    196             psTrace("ppStack", 2, "Photometering input %d of %d....\n", i, num);
     195            psTrace("ppStack", 2, "Photometering input %d of %d....\n", index, num);
    197196            pmFPAfileActivate(config->files, false, NULL);
    198             ppStackFileActivationSingle(config, PPSTACK_FILES_CONVOLVE, true, i);
     197            ppStackFileActivationSingle(config, PPSTACK_FILES_CONVOLVE, true, index);
    199198            if (options->convolve) {
    200199                pmFPAfileActivate(config->files, true, "PPSTACK.CONV.KERNEL");
    201200            }
    202             pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.INPUT", i); // File
     201            pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.INPUT", index); // File
    203202            pmFPAview *photView = ppStackFilesIterateDown(config);
    204203            if (!photView) {
     
    224223            ppStackFileActivation(config, PPSTACK_FILES_PREPARE, true);
    225224        }
    226     }
     225
     226        index++;
     227    }
     228    psFree(fileIter);
    227229
    228230    psString log = psStringCopy("Input seeing FWHMs:\n"); // Log message
  • trunk/ppStack/src/ppStackReadout.c

    r27400 r27402  
    2323    psVector *mask = options->inputMask; // Mask for inputs
    2424    psVector *weightings = options->weightings; // Weightings (1/noise^2) for each image
    25     psVector *exposures = options->exposures;   // Exposure times for each image
    2625    psVector *addVariance = options->matchChi2; // Additional variance when rejecting
    2726
    2827    job->results = ppStackReadoutInitial(config, outRO, thread->readouts, mask,
    29                                          weightings, exposures, addVariance);
     28                                         weightings, addVariance);
    3029    thread->busy = false;
    3130
     
    3837
    3938    psArray *args = job->args;          // Arguments
    40     ppStackThread *thread = args->data[0]; // Thread
    41     psArray *reject = args->data[1];    // Rejected pixels for each image
    42     ppStackOptions *options = args->data[2]; // Options
    43     pmConfig *config = args->data[3];   // Configuration
    44     bool safety = PS_SCALAR_VALUE(args->data[4], U8);    // Safety switch on?
    45     bool normalise = PS_SCALAR_VALUE(args->data[5], U8); // Normalise images?
     39    pmReadout *target = args->data[0];  // Output readout
     40    ppStackThread *thread = args->data[1]; // Thread
     41    psArray *reject = args->data[2];    // Rejected pixels for each image
     42    ppStackOptions *options = args->data[3]; // Options
     43    pmConfig *config = args->data[4];   // Configuration
     44    bool safety = PS_SCALAR_VALUE(args->data[5], U8);    // Safety switch on?
     45    bool normalise = PS_SCALAR_VALUE(args->data[6], U8); // Normalise images?
    4646
    4747    psVector *mask = options->inputMask; // Mask for inputs
    4848    psVector *weightings = options->weightings; // Weightings (1/noise^2) for each image
    49     psVector *exposures = options->exposures;   // Exposure times for each image
    5049    psVector *addVariance = options->matchChi2; // Additional variance when rejecting
    5150    psVector *norm = normalise ? options->norm : NULL; // Normalisations to apply to images
    5251
    53     bool status = ppStackReadoutFinal(config, options->outRO, options->expRO, thread->readouts, mask, reject,
    54                                       weightings, exposures, addVariance, safety, norm); // Status of operation
     52    bool status = ppStackReadoutFinal(config, target, thread->readouts, mask, reject,
     53                                      weightings, addVariance, safety, norm); // Status of operation
    5554
    5655    thread->busy = false;
    57 
    58     psAssert(status, "Stacking failed.");
    5956
    6057    return status;
     
    104101
    105102psArray *ppStackReadoutInitial(const pmConfig *config, pmReadout *outRO, const psArray *readouts,
    106                                const psVector *mask, const psVector *weightings, const psVector *exposures,
    107                                const psVector *addVariance)
     103                               const psVector *mask, const psVector *weightings, const psVector *addVariance)
    108104{
    109105    assert(config);
     
    153149        }
    154150
    155         stack->data[i] = pmStackDataAlloc(ro, weightings->data.F32[i], exposures->data.F32[i],
    156                                           addVariance->data.F32[i]);
    157     }
    158 
    159     if (!pmStackCombine(outRO, NULL, stack, maskVal | maskBad, maskSuspect, maskBad, kernelSize, iter,
     151        stack->data[i] = pmStackDataAlloc(ro, weightings->data.F32[i], addVariance->data.F32[i]);
     152    }
     153
     154    if (!pmStackCombine(outRO, stack, maskVal | maskBad, maskSuspect, maskBad, kernelSize, iter,
    160155                        combineRej, combineSys, combineDiscard, useVariance, safe, false)) {
    161156        psError(PS_ERR_UNKNOWN, false, "Unable to combine input readouts with rejection.");
     
    192187
    193188
    194 bool ppStackReadoutFinal(const pmConfig *config, pmReadout *outRO, pmReadout *expRO, const psArray *readouts,
     189bool ppStackReadoutFinal(const pmConfig *config, pmReadout *outRO, const psArray *readouts,
    195190                         const psVector *mask, const psArray *rejected, const psVector *weightings,
    196                          const psVector *exposures, const psVector *addVariance, bool safety,
    197                          const psVector *norm)
     191                         const psVector *addVariance, bool safety, const psVector *norm)
    198192{
    199193    assert(config);
    200194    assert(outRO);
    201     assert(expRO);
    202195    assert(readouts);
    203196    assert(!rejected || readouts->n == rejected->n);
     
    245238        }
    246239
    247         pmStackData *data = pmStackDataAlloc(ro, weightings->data.F32[i], exposures->data.F32[i],
     240        pmStackData *data = pmStackDataAlloc(ro, weightings->data.F32[i],
    248241                                             addVariance ? addVariance->data.F32[i] : NAN);
    249242        data->reject = rejected ? psMemIncrRefCounter(rejected->data[i]) : NULL;
     
    257250    }
    258251
    259     if (!pmStackCombine(outRO, expRO, stack, maskVal | maskBad, maskSuspect, maskBad, 0, iter, combineRej,
     252    if (!pmStackCombine(outRO, stack, maskVal | maskBad, maskSuspect, maskBad, 0, iter, combineRej,
    260253                        combineSys, combineDiscard, useVariance, safe, rejected)) {
    261254        psError(PS_ERR_UNKNOWN, false, "Unable to combine input readouts.");
     
    266259    pmCell *outCell = outRO->parent;    // Output cell
    267260    pmChip *outChip = outCell->parent;  // Output chip
     261
    268262    outRO->data_exists = true;
    269263    outCell->data_exists = true;
    270264    outChip->data_exists = true;
    271265
    272     pmCell *expCell = expRO->parent;    // Exposure cell
    273     pmChip *expChip = expCell->parent;  // Exposure chip
    274     expRO->data_exists = true;
    275     expCell->data_exists = true;
    276     expChip->data_exists = true;
    277 
    278266    psFree(stack);
    279267
  • trunk/ppStack/src/ppStackReject.c

    r27400 r27402  
    162162        psLogMsg("ppStack", PS_LOG_INFO, "Time to perform rejection on image %d: %f sec", i,
    163163                 psTimerClear("PPSTACK_REJECT"));
     164
     165        ppStackMemDump("reject");
    164166    }
    165167
  • trunk/ppStack/src/ppStackSources.c

    r27400 r27402  
    6464
    6565    if (!options->matchZPs && !options->photometry) {
    66         options->norm = psVectorAlloc(options->num, PS_TYPE_F32);
    67         psVectorInit(options->norm, 0.0);
     66        int num = psMetadataLookupS32(NULL, config->arguments, "INPUTS.NUM"); // Number of inputs
     67        options->norm = psVectorAlloc(num, PS_TYPE_F32);
     68        psVectorInit (options->norm, 0.0);
     69
     70        // XXX do I need to set this?
     71        // options->sumExposure = sumExpTime;
     72
    6873        return true;
    6974    }
     
    132137    }
    133138
    134     int num = options->num;            // Number of inputs
     139    int num = psMetadataLookupS32(NULL, config->arguments, "INPUTS.NUM"); // Number of inputs
    135140    psAssert(num == sourceLists->n, "Wrong number of source lists: %ld\n", sourceLists->n);
    136141
     
    141146    float airmassTerm = NAN;            // Airmass term
    142147    float zpTarget = NAN;               // Target zero point
     148    float sumExpTime = 0.0;             // Sum of the exposure time
    143149    int numGoodImages = 0;              // Number of good images
    144150    for (int i = 0; i < num; i++) {
     
    154160
    155161        pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.INPUT", i); // File of interest
     162        pmCell *cell = pmFPAviewThisCell(view, file->fpa); // Cell of interest
    156163
    157164#if defined(TESTING) && 0
     
    170177#endif
    171178
    172         float exptime = options->exposures->data.F32[i]; // Exposure time
     179        float exptime = psMetadataLookupF32(NULL, cell->concepts, "CELL.EXPOSURE"); // Exposure time
    173180        float airmass = psMetadataLookupF32(NULL, file->fpa->concepts, "FPA.AIRMASS"); // Airmass
    174181        const char *expFilter = psMetadataLookupStr(NULL, file->fpa->concepts, "FPA.FILTER"); // Filter name
     
    214221
    215222        zp->data.F32[i] = airmassTerm * airmass + 2.5 * log10(exptime);
    216     }
     223        sumExpTime += exptime;
     224
     225    }
     226
     227    options->sumExposure = sumExpTime;
    217228
    218229    if (numGoodImages == 0) {
     
    280291                }
    281292                psArray *sources = sourceLists->data[i]; // Sources of interest
    282                 float magCorr = zp->data.F32[i] + trans->data.F32[i] - 2.5*log10(options->sumExposure);
     293                float magCorr = zp->data.F32[i] + trans->data.F32[i] - 2.5*log10(sumExpTime);
    283294                if (zpExpNum == numGoodImages) {
    284295                    // Using measured zero points, so attempt to set target zero point
  • trunk/ppStack/src/ppStackThread.c

    r27400 r27402  
    284284
    285285    {
    286         psThreadTask *task = psThreadTaskAlloc("PPSTACK_FINAL_COMBINE", 6);
     286        psThreadTask *task = psThreadTaskAlloc("PPSTACK_FINAL_COMBINE", 7);
    287287        task->function = &ppStackReadoutFinalThread;
    288288        psThreadTaskAdd(task);
  • trunk/psModules/src/imcombine/pmStack.c

    r27400 r27402  
    3535
    3636//#define TESTING                         // Enable test output
    37 //#define TEST_X 843-1                     // x coordinate to examine
    38 //#define TEST_Y 813-1                     // y coordinate to examine
     37//#define TEST_X 2148-1                     // x coordinate to examine
     38//#define TEST_Y 248-1                     // y coordinate to examine
    3939//#define TEST_RADIUS 0                    // Radius to examine
    4040
     
    4646    psVector *variances;                // Pixel variances
    4747    psVector *weights;                  // Pixel weightings
    48     psVector *exps;                     // Pixel exposures
    4948    psVector *sources;                  // Pixel sources (which image did they come from?)
    5049    psVector *limits;                   // Rejection limits
     
    5857    psFree(buffer->variances);
    5958    psFree(buffer->weights);
    60     psFree(buffer->exps);
    6159    psFree(buffer->sources);
    6260    psFree(buffer->limits);
     
    7573    buffer->variances = psVectorAlloc(numImages, PS_TYPE_F32);
    7674    buffer->weights = psVectorAlloc(numImages, PS_TYPE_F32);
    77     buffer->exps = psVectorAlloc(numImages, PS_TYPE_F32);
    7875    buffer->sources = psVectorAlloc(numImages, PS_TYPE_U16);
    7976    buffer->limits = psVectorAlloc(numImages, PS_TYPE_F32);
     
    9895static bool combinationMeanVariance(float *mean, // Mean value, to return
    9996                                    float *var, // Variance value, to return
    100                                     float *exp, // Exposure time, to return
    101                                     float *expWeight,          // Weighted exposure time, to return
    10297                                    const psVector *values, // Values to combine
    10398                                    const psVector *variances, // Pixel variances to combine
    104                                     const psVector *exps,      // Exposure times to combine
    10599                                    const psVector *weights // Weights to apply
    106100                                    )
     
    127121    float sumVarianceWeight = 0.0;     // Sum of the pixel variances multiplied by the global weights
    128122    float sumWeight = 0.0;              // Sum of the image weights
    129     float sumExp = 0.0;                 // Sum of the exposure time
    130     float sumExpWeight = 0.0;           // Sum of the exposure time multiplied by the global weights
    131123    for (int i = 0; i < values->n; i++) {
    132124        sumValueWeight += values->data.F32[i] * weights->data.F32[i];
     
    135127            sumVarianceWeight += variances->data.F32[i] * PS_SQR(weights->data.F32[i]);
    136128        }
    137         if (exps) {
    138             sumExp += exps->data.F32[i];
    139             sumExpWeight += exps->data.F32[i] * weights->data.F32[i];
    140         }
    141129    }
    142130
     
    148136    if (var) {
    149137        *var = sumVarianceWeight / PS_SQR(sumWeight);
    150     }
    151     if (exp) {
    152         *exp = sumExp;
    153     }
    154     if (expWeight) {
    155         *expWeight = sumExpWeight / sumWeight;
    156138    }
    157139    return true;
     
    294276                           const psArray *inputs, // Stack data
    295277                           const psVector *weights, // Global (single value) weights for data, or NULL
    296                            const psVector *exps,    // Exposures for data, or NULL
    297278                           const psVector *addVariance, // Additional variance for data
    298279                           const psVector *reject, // Indices of pixels to reject, or NULL
     
    311292    psVector *pixelVariances = variance ? buffer->variances : NULL; // Variances for the pixel of interest
    312293    psVector *pixelWeights = buffer->weights; // Image weights for the pixel of interest
    313     psVector *pixelExps = buffer->exps;       // Exposure times
    314294    psVector *pixelSources = buffer->sources; // Sources for the pixel of interest
    315295    psVector *pixelLimits = buffer->limits; // Limits for the pixel of interest
     
    351331        }
    352332        pixelWeights->data.F32[numGood] = data->weight;
    353         pixelExps->data.F32[numGood] = data->exp;
    354333        pixelSources->data.U16[numGood] = i;
    355334        numGood++;
     
    368347    if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    369348        for (int i = 0; i < numGood; i++) {
    370             fprintf(stderr, "Input %d, pixel %d,%d (%" PRIu16 "): %f %f (%f) %f %f %d\n",
     349            fprintf(stderr, "Input %d, pixel %d,%d (%" PRIu16 "): %f %f (%f) %f %d\n",
    371350                    i, x, y, pixelSources->data.U16[i], pixelData->data.F32[i], pixelVariances->data.F32[i],
    372                     addVariance->data.F32[i], pixelWeights->data.F32[i], pixelExps->data.F32[i],
    373                     pixelSuspects->data.U8[i]);
     351                    addVariance->data.F32[i], pixelWeights->data.F32[i], pixelSuspects->data.U8[i]);
    374352        }
    375353    }
     
    384362                          psImage *mask, // Combined mask, for output
    385363                          psImage *variance, // Combined variance map, for output
    386                           psImage *exp,   // Exposure map (time), for output
    387                           psImage *expnum,       // Exposure map (number) for output
    388                           psImage *expweight,    // Exposure map (weighted time) for output
    389364                          int num,      // Number of good pixels
    390365                          combineBuffer *buffer, // Buffer with vectors
     
    397372    psVector *pixelVariances = variance ? buffer->variances : NULL; // Variances for the pixel of interest
    398373    psVector *pixelWeights = buffer->weights; // Image weights for the pixel of interest
    399     psVector *pixelExps = buffer->exps;       // Exposure times
    400374
    401375    // Default option is that the pixel is bad
    402376    float imageValue = NAN, varianceValue = NAN; // Value for combined image and variance map
    403377    psImageMaskType maskValue = bad;    // Value for combined mask
    404     float expValue = 0.0, expWeightValue = NAN; // Exposure value (straight, and weighted)
    405 
    406378    switch (num) {
    407379      case 0: {
     
    421393                  varianceValue = pixelVariances->data.F32[0];
    422394              }
    423               if (exp) {
    424                   expValue = pixelExps->data.F32[0];
    425               }
    426395              maskValue = 0;
    427396#ifdef TESTING
     
    442411          // Automatically accept the mean of the pixels only if we're not playing safe
    443412          if (!safe) {
    444               if (combinationMeanVariance(&imageValue, &varianceValue, &expValue, &expWeightValue,
    445                                           pixelData, pixelVariances, pixelExps, pixelWeights)) {
     413              float mean, var;   // Mean and variance from combination
     414              if (combinationMeanVariance(&mean, &var, pixelData, pixelVariances, pixelWeights)) {
     415                  imageValue = mean;
     416                  varianceValue = var;
    446417                  maskValue = 0;
    447418#ifdef TESTING
    448419                  if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    449420                      fprintf(stderr, "Two inputs to combine using unsafe, pixel %d,%d --> %f %f\n",
    450                               x, y, imageValue, varianceValue);
     421                              x, y, mean, var);
    451422                  }
    452423#endif
     
    464435      default: {
    465436          // Can combine without too much worrying
    466           if (!combinationMeanVariance(&imageValue, &varianceValue, &expValue, &expWeightValue,
    467                                        pixelData, pixelVariances, pixelExps, pixelWeights)) {
     437          float mean, var;           // Mean and variance of the combination
     438          if (!combinationMeanVariance(&mean, &var, pixelData, pixelVariances, pixelWeights)) {
    468439              break;
    469440          }
     441          imageValue = mean;
     442          varianceValue = var;
    470443          maskValue = 0;
    471444#ifdef TESTING
    472445          if (PS_SQR(x - TEST_X) + PS_SQR(y - TEST_Y) <= PS_SQR(TEST_RADIUS)) {
    473               fprintf(stderr, "Combined inputs, pixel %d,%d --> %f %f\n", x, y, imageValue, varianceValue);
     446              fprintf(stderr, "Combined inputs, pixel %d,%d --> %f %f\n", x, y, mean, var);
    474447          }
    475448#endif
     
    483456        variance->data.F32[y][x] = varianceValue;
    484457    }
    485     if (exp) {
    486         exp->data.F32[y][x] = expValue;
    487     }
    488     if (expnum) {
    489         expnum->data.PS_TYPE_IMAGE_MASK_DATA[y][x] = num;
    490     }
    491     if (expweight) {
    492         expweight->data.F32[y][x] = expWeightValue;
    493     }
     458
     459#ifdef TESTING
     460    mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] = num;
     461#endif
     462
    494463
    495464    return;
     
    834803                              int *numCols, int *numRows, // Size of (sky) images
    835804                              const psArray *input, // Input array of pmStackData to validate
    836                               const pmReadout *output, // Output readout
    837                               const pmReadout *exp    // Exposure map
     805                              const pmReadout *output // Output readout
    838806    )
    839807{
     
    901869    }
    902870
    903     if (exp) {
    904         PM_ASSERT_READOUT_NON_NULL(exp, false);
    905         if (exp->image) {
    906             PS_ASSERT_IMAGES_SIZE_EQUAL(exp->image, output->image, false);
    907         }
    908         if (exp->mask) {
    909             PS_ASSERT_IMAGES_SIZE_EQUAL(exp->mask, output->image, false);
    910         }
    911     }
    912 
    913871    return true;
    914872}
     
    979937
    980938/// Constructor
    981 pmStackData *pmStackDataAlloc(pmReadout *readout, float weight, float exp, float addVariance)
     939pmStackData *pmStackDataAlloc(pmReadout *readout, float weight, float addVariance)
    982940{
    983941    pmStackData *data = psAlloc(sizeof(pmStackData)); // Stack data, to return
     
    988946    data->inspect = NULL;
    989947    data->weight = weight;
    990     data->exp = exp;
    991948    data->addVariance = addVariance;
    992949
     
    995952
    996953/// Stack input images
    997 bool pmStackCombine(pmReadout *combined, pmReadout *expmaps, psArray *input,
    998                     psImageMaskType maskVal, psImageMaskType maskSuspect,
    999                     psImageMaskType bad, int kernelSize,
    1000                     float iter, float rej, float sys, float olympic,
     954bool pmStackCombine(pmReadout *combined, psArray *input, psImageMaskType maskVal, psImageMaskType maskSuspect,
     955                    psImageMaskType bad, int kernelSize, float iter, float rej, float sys, float olympic,
    1001956                    bool useVariance, bool safe, bool rejection)
    1002957{
     
    1004959    int num;                            // Number of inputs
    1005960    int numCols, numRows;               // Size of (sky) images
    1006     if (!validateInputData(&haveVariances, &num, &numCols, &numRows, input, combined, expmaps)) {
     961    if (!validateInputData(&haveVariances, &num, &numCols, &numRows, input, combined)) {
    1007962        return false;
    1008963    }
     
    1022977    psVector *addVariance = psVectorAlloc(num, PS_TYPE_F32); // Additional variance for each image
    1023978    psVector *weights = psVectorAlloc(num, PS_TYPE_F32); // Relative weighting for each image
    1024     psVector *exps = psVectorAlloc(num, PS_TYPE_F32);    // Exposure times for each image
    1025979    psArray *stack = psArrayAlloc(num); // Stack of readouts
    1026980    for (int i = 0; i < num; i++) {
     
    1028982        if (!data) {
    1029983            weights->data.F32[i] = 0.0;
    1030             exps->data.F32[i] = NAN;
    1031984            continue;
    1032985        }
    1033986        weights->data.F32[i] = data->weight;
    1034         exps->data.F32[i] = data->exp;
    1035987        pmReadout *ro = data->readout;  // Readout of interest
    1036988        stack->data[i] = psMemIncrRefCounter(ro);
     
    10931045        combined->variance = psImageAlloc(numCols, numRows, PS_TYPE_F32);
    10941046        combinedVariance = combined->variance;
    1095     }
    1096 
    1097     psImage *exp = NULL, *expnum = NULL; // Exposure map and exposure number
    1098     if (expmaps) {
    1099         if (!expmaps->image) {
    1100             expmaps->image = psImageAlloc(numCols, numRows, PS_TYPE_F32);
    1101         }
    1102         exp = expmaps->image;
    1103 
    1104         if (!expmaps->mask) {
    1105             expmaps->mask = psImageAlloc(numCols, numRows, PS_TYPE_IMAGE_MASK);
    1106         }
    1107         expnum = expmaps->mask;
    11081047    }
    11091048
     
    11441083            bool suspect;               // Suspect pixels in stack?
    11451084            combineExtract(&num, &suspect, buffer, combinedImage, combinedMask, combinedVariance,
    1146                            input, weights, exps, addVariance, reject, x, y, maskVal, maskSuspect);
    1147             combinePixels(combinedImage, combinedMask, combinedVariance, exp, expnum, NULL,
    1148                           num, buffer, x, y, bad, safe);
     1085                           input, weights, addVariance, reject, x, y, maskVal, maskSuspect);
     1086            combinePixels(combinedImage, combinedMask, combinedVariance, num, buffer, x, y, bad, safe);
    11491087
    11501088            if (iter > 0) {
  • trunk/psModules/src/imcombine/pmStack.h

    r27400 r27402  
    3131    psPixels *inspect;                  ///< Pixels to inspect
    3232    float weight;                       ///< Relative weighting for image
    33     float exp;                          ///< Exposure for image
    3433    float addVariance;                  ///< Additional variance when rejecting
    3534} pmStackData;
     
    3837pmStackData *pmStackDataAlloc(pmReadout *readout, ///< Warped readout (sky cell)
    3938                              float weight, ///< Weight to apply
    40                               float exp,    ///< Exposure time for input
    4139                              float addVariance ///< Additional variance when rejecting
    4240    );
     
    4442/// Stack input images
    4543bool pmStackCombine(pmReadout *combined,///< Combined readout (output)
    46                     pmReadout *expmap,  ///< Exposure map (output)
    4744                    psArray *input,     ///< Input array of pmStackData
    4845                    psImageMaskType maskVal, ///< Mask value of bad pixels
  • trunk/psModules/src/imcombine/pmSubtractionMask.c

    r27400 r27402  
    9292                    continue;
    9393                }
    94                 if (imageData1[y][x] > 50000) {
    95                     maskData1[y][x] = maskVal;
    96                     numBad++;
    97                     continue;
    98                 }
    99                 if (imageData2[y][x] > 50000) {
    100                     maskData2[y][x] = maskVal;
    101                     numBad++;
    102                     continue;
    103                 }
    10494                xMin = PS_MIN(xMin, x);
    10595                xMax = PS_MAX(xMax, x);
  • trunk/psModules/src/imcombine/pmSubtractionStamps.c

    r27400 r27402  
    107107            if ((image1 && image1->data.F32[y][x] < thresh1) ||
    108108                (image2 && image2->data.F32[y][x] < thresh2)) {
    109                 continue;
    110             }
    111             if ((image1 && image1->data.F32[y][x] > 30000) ||
    112                 (image2 && image2->data.F32[y][x] > 30000)) {
    113109                continue;
    114110            }
Note: See TracChangeset for help on using the changeset viewer.