IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 27319


Ignore:
Timestamp:
Mar 18, 2010, 12:01:15 PM (16 years ago)
Author:
Paul Price
Message:

Reworking ppStack to not grow the rejected pixels for the unconvolved stack. Also working with convolved and unconvolved stacks separately, so they don't have to be in memory at the same time. This should reduce the memory usage of ppStack.

Location:
trunk
Files:
15 edited

Legend:

Unmodified
Added
Removed

  • trunk/ppStack/src/ppStack.c

    r27160 r27319  
    7070            ppStackFileActivation(config, PPSTACK_FILES_PREPARE, true);
    7171            ppStackFileActivation(config, PPSTACK_FILES_CONVOLVE, true);
    72             ppStackFileActivation(config, PPSTACK_FILES_COMBINE, true);
     72            ppStackFileActivation(config, PPSTACK_FILES_STACK, true);
     73            ppStackFileActivation(config, PPSTACK_FILES_UNCONV, true);
    7374            ppStackFileActivation(config, PPSTACK_FILES_PHOT, true);
    7475            if (!ppStackFilesIterateUp(config)) {

  • trunk/ppStack/src/ppStack.h

    r27108 r27319  
    2727    PPSTACK_FILES_PREPARE,              // Files for preparation
    2828    PPSTACK_FILES_CONVOLVE,             // Files for convolution
    29     PPSTACK_FILES_COMBINE,              // Files for combination
     29    PPSTACK_FILES_STACK,                // Stack files
     30    PPSTACK_FILES_UNCONV,               // Unconvolved stack files
    3031    PPSTACK_FILES_PHOT                  // Files for photometry
    3132} ppStackFileList;

  • trunk/ppStack/src/ppStackCleanup.c

    r27004 r27319  
    66#include <pslib.h>
    77#include <psmodules.h>
     8#include <ppStats.h>
    89
    910#include "ppStack.h"
    … …  
    5960    }
    6061
     62    // Statistics on output
     63    if (options->stats) {
     64        psTrace("ppStack", 1, "Gathering statistics on stacked image....\n");
     65        psString maskBadStr = psMetadataLookupStr(NULL, recipe, "MASK.BAD"); // Name of bits for bad
     66        psImageMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels
     67
     68        pmFPAview *view = pmFPAviewAlloc(0); // View to readout
     69        view->chip = view->cell = view->readout = 0;
     70
     71        ppStatsFPA(options->stats, options->outRO->parent->parent->parent, view, maskBad, config);
     72
     73        psFree(view);
     74    }
     75
    6176    return true;
    6277}

  • trunk/ppStack/src/ppStackCombineFinal.c

    r27004 r27319  
    1515
    1616bool ppStackCombineFinal(pmReadout *target, ppStackThreadData *stack, psArray *covariances,
    17                          ppStackOptions *options, pmConfig *config, bool safe, bool normalise)
     17                         ppStackOptions *options, pmConfig *config, bool safe, bool normalise, bool grow)
    1818{
    1919    psAssert(stack, "Require stack");
    2020    psAssert(options, "Require options");
    2121    psAssert(config, "Require configuration");
     22
     23    int numCols = target->image->numCols, numRows = target->image->numRows; // Size of image
     24
     25    psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSTACK_RECIPE); // ppStack recipe
     26    psAssert(recipe, "We've thrown an error on this before.");
     27    float poorFrac = psMetadataLookupF32(NULL, recipe, "POOR.FRACTION"); // Fraction for "poor"
     28
     29    // Grow the list of rejected pixels, if desired
     30    psArray *reject = psArrayAlloc(options->num); // Pixels rejected for each image
     31    for (int i = 0; i < options->num; i++) {
     32        if (grow) {
     33            reject->data[i] = pmStackRejectGrow(options->rejected->data[i], numCols, numRows, poorFrac,
     34                                                options->regions->data[i], options->kernels->data[i]);
     35            if (!reject->data[i]) {
     36                psError(psErrorCodeLast(), false, "Unable to grow rejected pixels for image %d", i);
     37                psFree(reject);
     38                return false;
     39            }
     40        } else {
     41            reject->data[i] = psMemIncrRefCounter(options->rejected->data[i]);
     42        }
     43    }
    2244
    2345    if (!target->mask) {
    … …  
    3254            // Something went wrong
    3355            psError(psErrorCodeLast(), false, "Unable to read chunk %d", numChunk);
     56            psFree(reject);
    3457            return false;
    3558        }
    … …  
    4366        psArrayAdd(job->args, 1, target);
    4467        psArrayAdd(job->args, 1, thread);
     68        psArrayAdd(job->args, 1, reject);
    4569        psArrayAdd(job->args, 1, options);
    4670        psArrayAdd(job->args, 1, config);
    … …  
    4973        if (!psThreadJobAddPending(job)) {
    5074            psFree(job);
     75            psFree(reject);
    5176            return false;
    5277        }
    … …  
    5681    if (!psThreadPoolWait(true)) {
    5782        psError(psErrorCodeLast(), false, "Unable to do final combination.");
     83        psFree(reject);
    5884        return false;
    5985    }
     86
     87    psFree(reject);
    6088
    6189    // Sum covariance matrices

  • trunk/ppStack/src/ppStackCombinePrepare.c

    r27004 r27319  
    1010#include "ppStackLoop.h"
    1111
    12 bool ppStackCombinePrepare(ppStackThreadData *stack, ppStackOptions *options, pmConfig *config)
     12bool ppStackCombinePrepare(pmReadout **ro, const char *name, ppStackFileList files, ppStackThreadData *stack,
     13                           ppStackOptions *options, pmConfig *config)
    1314{
    1415    psAssert(stack, "Require stack");
    … …  
    2526
    2627    pmFPAfileActivate(config->files, false, NULL);
    27     ppStackFileActivation(config, PPSTACK_FILES_COMBINE, true);
     28    ppStackFileActivation(config, files, true);
    2829    pmFPAview *view = ppStackFilesIterateDown(config); // View to readout
    2930    if (!view) {
    … …  
    3132    }
    3233
    33     pmCell *outCell = pmFPAfileThisCell(config->files, view, "PPSTACK.OUTPUT"); // Output cell
    34     options->outRO = pmReadoutAlloc(outCell); // Output readout
    35 
    36     pmCell *unconvCell = pmFPAfileThisCell(config->files, view, "PPSTACK.UNCONV"); // Unconvolved cell
    37     options->unconvRO = pmReadoutAlloc(unconvCell); // Unconvolved readout
     34    pmCell *cell = pmFPAfileThisCell(config->files, view, name); // Output cell
     35    *ro = pmReadoutAlloc(cell); // Output readout
    3836
    3937    psFree(view);
    … …  
    4442    psImageMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels
    4543
    46     if (!pmReadoutStackDefineOutput(options->outRO, col0, row0, numCols, numRows, true, true, maskBad)) {
     44    if (!pmReadoutStackDefineOutput(*ro, col0, row0, numCols, numRows, true, true, maskBad)) {
    4745        psError(PPSTACK_ERR_ARGUMENTS, false, "Unable to prepare output.");
    4846        return false;
    4947    }
    5048
    51     options->unconvRO->image = psImageCopy(NULL, options->outRO->image, PS_TYPE_F32);
    52 //    options->unconvRO->mask = psImageCopy(NULL, options->outRO->mask, PS_TYPE_IMAGE_MASK);
    53     options->unconvRO->col0 = options->outRO->col0;
    54     options->unconvRO->row0 = options->outRO->row0;
    55 
    56 
    57 
    5849    return true;
    5950}

  • trunk/ppStack/src/ppStackFiles.c

    r27004 r27319  
    2121//                                 "PPSTACK.CONV.KERNEL", NULL };
    2222
    23 /// Output files for the combination
    24 static char *filesCombine[] = { "PPSTACK.OUTPUT", "PPSTACK.OUTPUT.MASK", "PPSTACK.OUTPUT.VARIANCE",
    25                                 "PPSTACK.UNCONV", "PPSTACK.UNCONV.MASK", "PPSTACK.UNCONV.VARIANCE",
    26                                 "PPSTACK.OUTPUT.JPEG1", "PPSTACK.OUTPUT.JPEG2", NULL };
     23/// Regular (convolved) stack files
     24static char *filesStack[] = { "PPSTACK.OUTPUT", "PPSTACK.OUTPUT.MASK", "PPSTACK.OUTPUT.VARIANCE",
     25                              "PPSTACK.OUTPUT.JPEG1", "PPSTACK.OUTPUT.JPEG2", NULL };
     26/// Unconvolved stack files
     27static char *filesUnconv[] = { "PPSTACK.UNCONV", "PPSTACK.UNCONV.MASK", "PPSTACK.UNCONV.VARIANCE", NULL };
    2728
    2829/// Files for photometry
    … …  
    3738      case PPSTACK_FILES_PREPARE:  return filesPrepare;
    3839      case PPSTACK_FILES_CONVOLVE: return filesConvolve;
    39       case PPSTACK_FILES_COMBINE:  return filesCombine;
     40      case PPSTACK_FILES_STACK:    return filesStack;
     41      case PPSTACK_FILES_UNCONV:   return filesUnconv;
    4042      case PPSTACK_FILES_PHOT:     return filesPhot;
    4143      default:

  • trunk/ppStack/src/ppStackFinish.c

    r27189 r27319  
    77#include <pslib.h>
    88#include <psmodules.h>
    9 #include <ppStats.h>
    109#include <psphot.h>
    1110
    … …  
    4948    }
    5049
    51 
    52     // Statistics on output
    53     if (options->stats) {
    54         psTrace("ppStack", 1, "Gathering statistics on stacked image....\n");
    55         psString maskBadStr = psMetadataLookupStr(NULL, recipe, "MASK.BAD"); // Name of bits for bad
    56         psImageMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels
    57 
    58         pmFPAview *view = pmFPAviewAlloc(0); // View to readout
    59         view->chip = view->cell = view->readout = 0;
    60 
    61         ppStatsFPA(options->stats, options->outRO->parent->parent->parent, view, maskBad, config);
    62 
    63         psFree(view);
    64     }
    65 
    66     ppStackMemDump("stats");
    67 
    68     psFree(options->outRO); options->outRO = NULL;
    6950
    7051    return true;

  • trunk/ppStack/src/ppStackLoop.c

    r27309 r27319  
    111111
    112112    // Prepare for combination
    113     if (!ppStackCombinePrepare(stack, options, config)) {
     113    if (!ppStackCombinePrepare(&options->outRO, "PPSTACK.OUTPUT", PPSTACK_FILES_STACK,
     114                               stack, options, config)) {
    114115        psError(psErrorCodeLast(), false, "Unable to prepare for combination.");
    115116        psFree(stack);
    … …  
    148149    // Final combination
    149150    psTrace("ppStack", 2, "Final stack of convolved images....\n");
    150     if (!ppStackCombineFinal(options->outRO, stack, options->convCovars, options, config, false, false)) {
     151    if (!ppStackCombineFinal(options->outRO, stack, options->convCovars, options, config,
     152                             false, false, true)) {
    151153        psError(psErrorCodeLast(), false, "Unable to perform final combination.");
    152154        psFree(stack);
    … …  
    163165        return false;
    164166    }
    165     psLogMsg("ppStack", PS_LOG_INFO, "Stage 6: Cleanup, WCS & JPEGS: %f sec", psTimerClear("PPSTACK_STEPS"));
     167
     168    // Photometry
     169    psTrace("ppStack", 1, "Photometering stacked image....\n");
     170    if (!ppStackPhotometry(options, config)) {
     171        psError(psErrorCodeLast(), false, "Unable to perform photometry.");
     172        return false;
     173    }
     174    psLogMsg("ppStack", PS_LOG_INFO, "Stage 6: Photometry Analysis: %f sec", psTimerClear("PPSTACK_STEPS"));
     175    ppStackMemDump("photometry");
     176
     177    if (!ppStackFilesIterateUp(config)) {
     178        psError(psErrorCodeLast(), false, "Unable to close files.");
     179        return false;
     180    }
     181    ppStackFileActivation(config, PPSTACK_FILES_STACK, false);
     182    ppStackFileActivation(config, PPSTACK_FILES_PHOT, false);
     183    options->outRO->data_exists = false;
     184    options->outRO->parent->data_exists = false;
     185    options->outRO->parent->parent->data_exists = false;
     186    psFree(options->outRO);
     187    options->outRO = NULL;
     188
     189    psLogMsg("ppStack", PS_LOG_INFO, "Stage 7: Cleanup, WCS & JPEGS: %f sec", psTimerClear("PPSTACK_STEPS"));
    166190    ppStackMemDump("cleanup");
    167191
    … …  
    177201            return false;
    178202        }
     203
     204        // Prepare for combination
     205        if (!ppStackCombinePrepare(&options->unconvRO, "PPSTACK.UNCONV", PPSTACK_FILES_UNCONV,
     206                                   stack, options, config)) {
     207            psError(psErrorCodeLast(), false, "Unable to prepare for combination.");
     208            psFree(stack);
     209            return false;
     210        }
     211
    179212        psTrace("ppStack", 2, "Stack of unconvolved images....\n");
    180         if (!ppStackCombineFinal(options->unconvRO, stack, options->origCovars, options, config, false, true)) {
     213        if (!ppStackCombineFinal(options->unconvRO, stack, options->origCovars, options, config,
     214                                 false, true, false)) {
    181215            psError(psErrorCodeLast(), false, "Unable to perform unconvolved combination.");
    182216            psFree(stack);
    183217            return false;
    184218        }
    185         psLogMsg("ppStack", PS_LOG_INFO, "Stage 7: Unconvolved Stack: %f sec", psTimerClear("PPSTACK_STEPS"));
     219        psLogMsg("ppStack", PS_LOG_INFO, "Stage 8: Unconvolved Stack: %f sec", psTimerClear("PPSTACK_STEPS"));
    186220        ppStackMemDump("unconv");
     221
     222        if (!ppStackFilesIterateUp(config)) {
     223            psError(psErrorCodeLast(), false, "Unable to close files.");
     224            psFree(stack);
     225            return false;
     226        }
     227        ppStackFileActivation(config, PPSTACK_FILES_UNCONV, false);
     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;
    187233
    188234        psFree(stack);
    … …  
    190236    psFree(options->cells); options->cells = NULL;
    191237#endif
    192 
    193     // Photometry
    194     psTrace("ppStack", 1, "Photometering stacked image....\n");
    195     if (!ppStackPhotometry(options, config)) {
    196         psError(psErrorCodeLast(), false, "Unable to perform photometry.");
    197         return false;
    198     }
    199     psLogMsg("ppStack", PS_LOG_INFO, "Stage 8: Photometry Analysis: %f sec", psTimerClear("PPSTACK_STEPS"));
    200     ppStackMemDump("photometry");
    201238
    202239    // Finish up

  • trunk/ppStack/src/ppStackLoop.h

    r27160 r27319  
    77#include "ppStackOptions.h"
    88#include "ppStackThread.h"
     9#include "ppStack.h"
    910
    1011
    … …  
    3738// Prepare for combination
    3839bool ppStackCombinePrepare(
     40    pmReadout **readout,                // Readout to set
     41    const char *name,                   // Name of file
     42    ppStackFileList files,              // Files of interest
    3943    ppStackThreadData *stack,           // Stack
    4044    ppStackOptions *options,            // Options
    … …  
    6367    pmConfig *config,                   // Configuration
    6468    bool safe,                          // Allow safe combination?
    65     bool norm                           // Normalise images?
     69    bool norm,                          // Normalise images?
     70    bool grow                           // Grow rejection masks?
    6671    );
    6772

  • trunk/ppStack/src/ppStackPhotometry.c

    r27158 r27319  
    2626    psTimerStart("PPSTACK_PHOT");
    2727
    28     ppStackFileActivation(config, PPSTACK_FILES_COMBINE, false);
     28    pmFPAfileActivate(config->files, false, NULL);
    2929    ppStackFileActivation(config, PPSTACK_FILES_PHOT, true);
    3030    pmFPAview *photView = ppStackFilesIterateDown(config); // View to readout
    31     ppStackFileActivation(config, PPSTACK_FILES_COMBINE, true);
     31    ppStackFileActivation(config, PPSTACK_FILES_STACK, true);
    3232
    3333    pmFPAfile *photFile = psMetadataLookupPtr(NULL, config->files, "PSPHOT.INPUT"); // File for photometry

  • trunk/ppStack/src/ppStackReadout.c

    r27105 r27319  
    3939    pmReadout *target = args->data[0];  // Output readout
    4040    ppStackThread *thread = args->data[1]; // Thread
    41     ppStackOptions *options = args->data[2]; // Options
    42     pmConfig *config = args->data[3];   // Configuration
    43     bool safety = PS_SCALAR_VALUE(args->data[4], U8);    // Safety switch on?
    44     bool normalise = PS_SCALAR_VALUE(args->data[5], U8); // Normalise images?
     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?
    4546
    4647    psVector *mask = options->inputMask; // Mask for inputs
    47     psArray *rejected = options->rejected; // Rejected pixels
    4848    psVector *weightings = options->weightings; // Weightings (1/noise^2) for each image
    4949    psVector *addVariance = options->matchChi2; // Additional variance when rejecting
    5050    psVector *norm = normalise ? options->norm : NULL; // Normalisations to apply to images
    5151
    52     bool status = ppStackReadoutFinal(config, target, thread->readouts, mask, rejected,
     52    bool status = ppStackReadoutFinal(config, target, thread->readouts, mask, reject,
    5353                                      weightings, addVariance, safety, norm); // Status of operation
    5454

  • trunk/ppStack/src/ppStackReject.c

    r27309 r27319  
    3131
    3232    float threshold = psMetadataLookupF32(NULL, recipe, "THRESHOLD.MASK"); // Threshold for mask deconvolution
    33     float poorFrac = psMetadataLookupF32(NULL, recipe, "POOR.FRACTION"); // Fraction for "poor"
    3433    float imageRej = psMetadataLookupF32(NULL, recipe, "IMAGE.REJ"); // Maximum fraction of image to reject
    3534                                                                     // before rejecting entire image
    … …  
    109108
    110109        psPixels *reject = pmStackReject(options->inspect->data[i], options->numCols, options->numRows,
    111                                          threshold, poorFrac, stride, options->regions->data[i],
     110                                         threshold, stride, options->regions->data[i],
    112111                                         options->kernels->data[i]); // Rejected pixels
    113112
    … …  
    130129                options->inputMask->data.PS_TYPE_VECTOR_MASK_DATA[i] |= PPSTACK_MASK_BAD;
    131130                numRejected++;
     131            } else {
     132                // Add to list of pixels already rejected
     133                reject = psPixelsConcatenate(reject, options->rejected->data[i]);
     134                options->rejected->data[i] = psPixelsDuplicates(options->rejected->data[i], reject);
    132135            }
    133         }
    134 
    135         if (reject) {
    136             // Add to list of pixels already rejected
    137             reject = psPixelsConcatenate(reject, options->rejected->data[i]);
    138             options->rejected->data[i] = psPixelsDuplicates(options->rejected->data[i], reject);
    139136        }
    140137
    … …  
    170167
    171168    psFree(options->inspect); options->inspect = NULL;
    172     psFree(options->kernels); options->kernels = NULL;
    173     psFree(options->regions); options->regions = NULL;
    174169
    175170    if (numRejected >= num) {

  • trunk/ppStack/src/ppStackThread.c

    r27004 r27319  
    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/pmStackReject.c

    r27307 r27319  
    1010#include "pmSubtractionThreads.h"
    1111#include "pmSubtractionKernels.h"
     12
     13#include "pmStackReject.h"
    1214
    1315#define PIXEL_LIST_BUFFER 100           // Number of pixels to add to list at a time
    … …  
    115117
    116118
    117 psPixels *pmStackReject(const psPixels *in, int numCols, int numRows, float threshold, float poorFrac,
    118                         int stride, const psArray *subRegions, const psArray *subKernels)
     119psPixels *pmStackReject(const psPixels *in, int numCols, int numRows, float threshold, int stride,
     120                        const psArray *subRegions, const psArray *subKernels)
    119121{
    120122    PS_ASSERT_PIXELS_NON_NULL(in, NULL);
    … …  
    223225    psTrace("psModules.imcombine", 7, "Found %ld bad pixels", bad->n);
    224226
    225     // Now, grow the mask to include everything that touches a bad pixel in the convolution
    226     psImage *source = psPixelsToMask(NULL, bad, psRegionSet(0, numCols - 1, 0, numRows - 1),
     227    return bad;
     228}
     229
     230
     231psPixels *pmStackRejectGrow(const psPixels *in, int numCols, int numRows, float poorFrac,
     232                            const psArray *subRegions, const psArray *subKernels)
     233{
     234    PS_ASSERT_PIXELS_NON_NULL(in, NULL);
     235    PS_ASSERT_ARRAY_NON_NULL(subRegions, NULL);
     236    PS_ASSERT_ARRAY_NON_NULL(subKernels, NULL);
     237    PS_ASSERT_ARRAYS_SIZE_EQUAL(subRegions, subKernels, NULL);
     238
     239    psImage *source = psPixelsToMask(NULL, in, psRegionSet(0, numCols - 1, 0, numRows - 1),
    227240                                     PM_STACK_MASK_BAD); // Mask image to grow
    228241
    … …  
    243256    bool oldThreads = psImageConvolveSetThreads(false); // Old value of threading for psImageColvolve
    244257
    245     psImage *target = psImageRecycle(convolved, numCols, numRows, PS_TYPE_IMAGE_MASK); // Grown image
     258    psImage *target = psImageAlloc(numCols, numRows, PS_TYPE_IMAGE_MASK); // Grown image
    246259    psImageInit(target, 0);
    247260    for (int i = 0; i < subRegions->n; i++) {
    … …  
    326339
    327340    psFree(source);
    328     bad = psPixelsFromMask(bad, target, PM_STACK_MASK_ALL);
     341    psPixels *bad = psPixelsFromMask(NULL, target, PM_STACK_MASK_ALL); // All bad pixels
    329342    psFree(target);
    330343    psTrace("psModules.imcombine", 7, "Total %ld bad pixels", bad->n);

  • trunk/psModules/src/imcombine/pmStackReject.h

    r20568 r27319  
    1212                        int numCols, int numRows, ///< Size of image of interest
    1313                        float threshold, ///< Threshold on convolved image, 0..1
    14                         float poorFrac, ///< Fraction for "poor"
    1514                        int stride,     ///< Size of convolution patches
    1615                        const psArray *regions, ///< Array of image regions for image
    1716                        const psArray *kernels ///< Array of kernel parameters for each region
     17    );
     18
     19/// Given a list of pixels from the convolved image, we grow them by convolution to get the list of all pixels
     20/// which should be rejected.
     21psPixels *pmStackRejectGrow(const psPixels *in, ///< List of pixels in the convolved image
     22                            int numCols, int numRows, ///< Size of image of interest
     23                            float poorFrac, ///< Fraction for "poor"
     24                            const psArray *regions, ///< Array of image regions for image
     25                            const psArray *kernels ///< Array of kernel parameters for each region
    1826    );
    1927
Note: See TracChangeset for help on using the changeset viewer.