IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 27402 for trunk/ppStack


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/ppStack/src
Files:
16 edited

Legend:

Unmodified
Added
Removed
  • 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);
Note: See TracChangeset for help on using the changeset viewer.