IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Aug 22, 2008, 12:45:36 PM (18 years ago)
Author:
Paul Price
Message:

Reworking the masking so as to avoid exploding the number of bad pixels in the convolved image. Introduced the bad/poor distinction. Bad pixels that are convolved now have a large 'poor' halo in the mask, with a smaller 'bad' halo. Got rid of the 'FOOTPRINT' mask value (internal subtraction mask): there's no need to convolve the mask in order to find out if a sprinkling of footprints are bad (not sure this is so efficient when a list of sources isn't supplied), and I needed the mask bits to signify the difference between 'bad' and 'poor'. Tested on a single subtraction, which seems to work.

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/psModules/src/imcombine/pmSubtractionMatch.c

    r18962 r19164  
    121121                        int inner, int ringsOrder, int binning, float penalty,
    122122                        bool optimum, const psVector *optFWHMs, int optOrder, float optThreshold,
    123                         int iter, float rej, psMaskType maskBad, psMaskType maskBlank,
    124                         float badFrac, pmSubtractionMode mode)
     123                        int iter, float rej, psMaskType maskVal, psMaskType maskBad, psMaskType maskPoor,
     124                        float poorFrac, float badFrac, pmSubtractionMode subMode)
    125125{
    126     if (mode != PM_SUBTRACTION_MODE_2) {
     126    if (subMode != PM_SUBTRACTION_MODE_2) {
    127127        PM_ASSERT_READOUT_NON_NULL(conv1, false);
    128128        if (conv1->image) {
    … …  
    139139        }
    140140    }
    141     if (mode != PM_SUBTRACTION_MODE_1) {
     141    if (subMode != PM_SUBTRACTION_MODE_1) {
    142142        PM_ASSERT_READOUT_NON_NULL(conv2, false);
    143143        if (conv2->image) {
    … …  
    190190    PS_ASSERT_INT_POSITIVE(iter, false);
    191191    PS_ASSERT_FLOAT_LARGER_THAN(rej, 0.0, false);
     192    // Don't care about maskVal
    192193    // Don't care about maskBad
    193     // Don't care about maskBlank
     194    // Don't care about maskPoor
     195    PS_ASSERT_FLOAT_LARGER_THAN(poorFrac, 0.0, NULL);
     196    PS_ASSERT_FLOAT_LESS_THAN_OR_EQUAL(poorFrac, 1.0, NULL);
    194197    if (isfinite(badFrac)) {
    195198        PS_ASSERT_FLOAT_LARGER_THAN(badFrac, 0.0, NULL);
    … …  
    230233    memCheck("start");
    231234
    232     subMask = pmSubtractionMask(ro1->mask, ro2 ? ro2->mask : NULL, maskBad, size, footprint,
     235    subMask = pmSubtractionMask(ro1->mask, ro2 ? ro2->mask : NULL, maskVal, size, footprint,
    233236                                badFrac, useFFT);
    234237    if (!subMask) {
    … …  
    267270            if (sources) {
    268271                stamps = pmSubtractionStampsSetFromSources(sources, subMask, region, footprint,
    269                                                            stampSpacing, mode);
     272                                                           stampSpacing, subMode);
    270273            } else if (stampsName && strlen(stampsName) > 0) {
    271274                stamps = pmSubtractionStampsSetFromFile(stampsName, ro1->image, subMask, region, footprint,
    272                                                         stampSpacing, mode);
     275                                                        stampSpacing, subMode);
    273276            }
    274277
    … …  
    276279            // doesn't matter.
    277280            if (!getStamps(&stamps, ro1, ro2, subMask, weight, NULL, threshold, stampSpacing,
    278                            size, footprint, mode)) {
     281                           size, footprint, subMode)) {
    279282                goto MATCH_ERROR;
    280283            }
    281284
    282             if (mode == PM_SUBTRACTION_MODE_UNSURE) {
     285            if (subMode == PM_SUBTRACTION_MODE_UNSURE) {
    283286                // Get backgrounds
    284287                psStats *bgStats = psStatsAlloc(BG_STAT); // Statistics for background
    285288                psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS, 0); // Random number generator
    286289                psVector *buffer = NULL;// Buffer for stats
    287                 if (!psImageBackground(bgStats, &buffer, ro1->image, ro1->mask, maskBad, rng)) {
     290                if (!psImageBackground(bgStats, &buffer, ro1->image, ro1->mask, maskVal, rng)) {
    288291                    psError(PS_ERR_UNKNOWN, false, "Unable to measure background of image 1.");
    289292                    psFree(bgStats);
    … …  
    293296                }
    294297                float bg1 = psStatsGetValue(bgStats, BG_STAT); // Background for image 1
    295                 if (!psImageBackground(bgStats, &buffer, ro2->image, ro2->mask, maskBad, rng)) {
     298                if (!psImageBackground(bgStats, &buffer, ro2->image, ro2->mask, maskVal, rng)) {
    296299                    psError(PS_ERR_UNKNOWN, false, "Unable to measure background of image 2.");
    297300                    psFree(bgStats);
    … …  
    317320                    goto MATCH_ERROR;
    318321                }
    319                 mode = newMode;
     322                subMode = newMode;
    320323            }
    321324
    … …  
    323326            if (optimum && (type == PM_SUBTRACTION_KERNEL_ISIS || type == PM_SUBTRACTION_KERNEL_GUNK)) {
    324327                kernels = pmSubtractionKernelsOptimumISIS(type, size, inner, spatialOrder, optFWHMs, optOrder,
    325                                                           stamps, footprint, optThreshold, penalty, mode);
     328                                                          stamps, footprint, optThreshold, penalty, subMode);
    326329                if (!kernels) {
    327330                    psErrorClear();
    … …  
    332335                // Not an ISIS/GUNK kernel, or the optimum kernel search failed
    333336                kernels = pmSubtractionKernelsGenerate(type, size, spatialOrder, isisWidths, isisOrders,
    334                                                        inner, binning, ringsOrder, penalty, mode);
     337                                                       inner, binning, ringsOrder, penalty, subMode);
    335338            }
    336339
    … …  
    344347                }
    345348
    346                 if (mode == PM_SUBTRACTION_MODE_1 || mode == PM_SUBTRACTION_MODE_DUAL) {
     349                if (subMode == PM_SUBTRACTION_MODE_1 || subMode == PM_SUBTRACTION_MODE_DUAL) {
    347350                    psMetadataAddPtr(conv1->analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_KERNEL,
    348351                                     PS_DATA_UNKNOWN | PS_META_DUPLICATE_OK, "Subtraction kernels", kernels);
    349352                    psMetadataAddS32(conv1->analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_MODE,
    350                                      PS_META_DUPLICATE_OK, "Subtraction kernels", mode);
     353                                     PS_META_DUPLICATE_OK, "Subtraction kernels", subMode);
    351354                    psMetadataAddPtr(conv1->analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_REGION,
    352355                                     PS_DATA_REGION | PS_META_DUPLICATE_OK,
    353356                                     "Region over which subtraction was performed", subRegion);
    354357                }
    355                 if (mode == PM_SUBTRACTION_MODE_2 || mode == PM_SUBTRACTION_MODE_DUAL) {
     358                if (subMode == PM_SUBTRACTION_MODE_2 || subMode == PM_SUBTRACTION_MODE_DUAL) {
    356359                    psMetadataAddPtr(conv2->analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_KERNEL,
    357360                                     PS_DATA_UNKNOWN | PS_META_DUPLICATE_OK, "Subtraction kernels", kernels);
    358361                    psMetadataAddS32(conv2->analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_MODE,
    359                                      PS_META_DUPLICATE_OK, "Subtraction kernels", mode);
     362                                     PS_META_DUPLICATE_OK, "Subtraction kernels", subMode);
    360363                    psMetadataAddPtr(conv2->analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_REGION,
    361364                                     PS_DATA_REGION | PS_META_DUPLICATE_OK,
    … …  
    374377
    375378                if (!getStamps(&stamps, ro1, ro2, subMask, weight, region, threshold, stampSpacing,
    376                                size, footprint, mode)) {
     379                               size, footprint, subMode)) {
    377380                    goto MATCH_ERROR;
    378381                }
    … …  
    434437            stamps = NULL;
    435438
    436             if (mode == PM_SUBTRACTION_MODE_1 || mode == PM_SUBTRACTION_MODE_DUAL) {
     439            if (subMode == PM_SUBTRACTION_MODE_1 || subMode == PM_SUBTRACTION_MODE_DUAL) {
    437440                psMetadataAddS32(conv1->analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_STAMPS_NUM,
    438441                                 PS_META_DUPLICATE_OK, "Number of good stamps", numStamps);
    … …  
    440443                                 PS_META_DUPLICATE_OK, "RMS deviation of stamps", rmsStamps);
    441444            }
    442             if (mode == PM_SUBTRACTION_MODE_2 || mode == PM_SUBTRACTION_MODE_DUAL) {
     445            if (subMode == PM_SUBTRACTION_MODE_2 || subMode == PM_SUBTRACTION_MODE_DUAL) {
    443446                psMetadataAddS32(conv2->analysis, PS_LIST_TAIL, PM_SUBTRACTION_ANALYSIS_STAMPS_NUM,
    444447                                 PS_META_DUPLICATE_OK, "Number of good stamps", numStamps);
    … …  
    454457                int fullSize = 2 * size + 1 + 1;    // Full size of kernel
    455458                int imageSize = (2 * KERNEL_MOSAIC + 1) * fullSize;
    456                 psImage *convKernels = psImageAlloc((mode == PM_SUBTRACTION_MODE_DUAL ? 2 : 1) *
     459                psImage *convKernels = psImageAlloc((subMode == PM_SUBTRACTION_MODE_DUAL ? 2 : 1) *
    457460                                                    imageSize - 1 +
    458                                                     (mode == PM_SUBTRACTION_MODE_DUAL ? 4 : 0),
     461                                                    (subMode == PM_SUBTRACTION_MODE_DUAL ? 4 : 0),
    459462                                                    imageSize - 1, PS_TYPE_F32);
    460463                psImageInit(convKernels, NAN);
    … …  
    479482                        psFree(kernel);
    480483
    481                         if (mode == PM_SUBTRACTION_MODE_DUAL) {
     484                        if (subMode == PM_SUBTRACTION_MODE_DUAL) {
    482485                            kernel = pmSubtractionKernelImage(kernels, (float)i / (float)KERNEL_MOSAIC,
    483486                                                              (float)j / (float)KERNEL_MOSAIC,
    … …  
    534537
    535538            psTrace("psModules.imcombine", 2, "Convolving...\n");
    536             if (!pmSubtractionConvolve(conv1, conv2, ro1, ro2, subMask, maskBlank, region, kernels,
    537                                        true, useFFT)) {
     539            if (!pmSubtractionConvolve(conv1, conv2, ro1, ro2, subMask, maskBad, maskPoor, poorFrac,
     540                                       region, kernels, true, useFFT)) {
    538541                psError(PS_ERR_UNKNOWN, false, "Unable to convolve image.");
    539542                goto MATCH_ERROR;
    … …  
    541544
    542545            // Set the variance factors
    543             switch (mode) {
     546            switch (subMode) {
    544547              case PM_SUBTRACTION_MODE_1: {
    545548                  recordVarianceFactor(conv1, pmSubtractionVarianceFactor(kernels, 0.5, 0.5, false),
    … …  
    565568
    566569            // There is data in the readout now
    567             if (mode == PM_SUBTRACTION_MODE_1 || mode == PM_SUBTRACTION_MODE_DUAL) {
     570            if (subMode == PM_SUBTRACTION_MODE_1 || subMode == PM_SUBTRACTION_MODE_DUAL) {
    568571                conv1->data_exists = true;
    569572                if (conv1->parent) {
    … …  
    572575                }
    573576            }
    574             if (mode == PM_SUBTRACTION_MODE_2 || mode == PM_SUBTRACTION_MODE_DUAL) {
     577            if (subMode == PM_SUBTRACTION_MODE_2 || subMode == PM_SUBTRACTION_MODE_DUAL) {
    575578                conv2->data_exists = true;
    576579                if (conv2->parent) {
    … …  
    590593    weight = NULL;
    591594
    592     if (!pmSubtractionBorder(conv1->image, conv1->weight, conv1->mask, size, maskBlank)) {
     595    if (!pmSubtractionBorder(conv1->image, conv1->weight, conv1->mask, size, maskBad)) {
    593596        psError(PS_ERR_UNKNOWN, false, "Unable to set border of convolved image.");
    594597        goto MATCH_ERROR;
    … …  
    600603#ifdef TESTING
    601604    {
    602         if (mode == PM_SUBTRACTION_MODE_1 || mode == PM_SUBTRACTION_MODE_DUAL) {
     605        if (subMode == PM_SUBTRACTION_MODE_1 || subMode == PM_SUBTRACTION_MODE_DUAL) {
    603606            psFits *fits = psFitsOpen("convolved1.fits", "w");
    604607            psFitsWriteImage(fits, NULL, conv1->image, 0, NULL);
    … …  
    606609        }
    607610
    608         if (mode == PM_SUBTRACTION_MODE_2 || mode == PM_SUBTRACTION_MODE_DUAL) {
     611        if (subMode == PM_SUBTRACTION_MODE_2 || subMode == PM_SUBTRACTION_MODE_DUAL) {
    609612            psFits *fits = psFitsOpen("convolved2.fits", "w");
    610613            psFitsWriteImage(fits, NULL, conv2->image, 0, NULL);
Note: See TracChangeset for help on using the changeset viewer.