IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Sep 10, 2007, 10:10:05 AM (19 years ago)
Author:
Paul Price
Message:

Preparing to take a (long) list of sources as candidate stamps. Changed the list of stamps from a simple array to a more complicated structure so that I can carry around a list of candidate stamps for each region of interest; when one is rejected, it is replaced by the next brightest within the same region. Optionally apply an exclusion zone around stamps, which is important when we're convolving to a specified PSF (as opposed to convolving to match an image).

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/psModules/src/imcombine/pmSubtractionParams.c

    r14738 r14801  
    171171
    172172pmSubtractionKernels *pmSubtractionKernelsOptimumISIS(pmSubtractionKernelsType type, int size, int inner,
    173                                                      int spatialOrder, psVector *fwhms, int maxOrder,
    174                                                      const psArray *stamps, int footprint, float tolerance)
     173                                                      int spatialOrder, psVector *fwhms, int maxOrder,
     174                                                      const pmSubtractionStampList *stamps, int footprint,
     175                                                      float tolerance)
    175176{
    176177    if (type != PM_SUBTRACTION_KERNEL_ISIS && type != PM_SUBTRACTION_KERNEL_GUNK) {
     
    184185    PS_ASSERT_VECTOR_TYPE(fwhms, PS_TYPE_F32, NULL);
    185186    PS_ASSERT_INT_NONNEGATIVE(maxOrder, NULL);
    186     PS_ASSERT_ARRAY_NON_NULL(stamps, NULL);
     187    PM_ASSERT_SUBTRACTION_STAMP_LIST_NON_NULL(stamps, NULL);
    187188    PS_ASSERT_INT_NONNEGATIVE(footprint, NULL);
    188189    PS_ASSERT_FLOAT_LARGER_THAN(tolerance, 0.0, NULL);
     
    208209
    209210    // Need to save the stamp inputs --- we're changing the values!
    210     int numStamps = stamps->n        // Number of stamps
     211    int numStamps = stamps->num;        // Number of stamps
    211212    psArray *inputs = psArrayAlloc(numStamps); // Deep copies of the inputs
     213    psVector *badStamps = psVectorAlloc(numStamps, PS_TYPE_U8); // Mark the bad stamps
     214    psVectorInit(badStamps, 0);
    212215    for (int i = 0; i < numStamps; i++) {
    213         pmSubtractionStamp *stamp = stamps->data[i]; // Stamp of interest
    214         if (stamp->status == PM_SUBTRACTION_STAMP_REJECTED || stamp->status == PM_SUBTRACTION_STAMP_NONE) {
     216        pmSubtractionStamp *stamp = stamps->stamps->data[i]; // Stamp of interest
     217        if (stamp->status != PM_SUBTRACTION_STAMP_CALCULATE && stamp->status != PM_SUBTRACTION_STAMP_USED) {
     218            badStamps->data.U8[i] = 0xff;
    215219            continue;
    216220        }
     
    230234    int numPixels = 0;                  // Number of pixels contributing to chi^2
    231235    for (int i = 0; i < numStamps; i++) {
    232         pmSubtractionStamp *stamp = stamps->data[i]; // Stamp of interest
    233         if (stamp->status == PM_SUBTRACTION_STAMP_REJECTED || stamp->status == PM_SUBTRACTION_STAMP_NONE) {
     236        pmSubtractionStamp *stamp = stamps->stamps->data[i]; // Stamp of interest
     237        if (badStamps->data.U8[i]) {
    234238            continue;
    235239        }
     
    238242            psFree(inputs);
    239243            psFree(kernels);
     244            psFree(badStamps);
    240245            return NULL;
    241246        }
     
    255260            psFree(inputs);
    256261            psFree(kernels);
     262            psFree(badStamps);
    257263            return NULL;
    258264        }
     
    287293
    288294            for (int j = 0; j < numStamps; j++) {
    289                 pmSubtractionStamp *stamp = stamps->data[j]; // Stamp of interest
    290                 if (stamp->status == PM_SUBTRACTION_STAMP_REJECTED ||
    291                     stamp->status == PM_SUBTRACTION_STAMP_NONE) {
     295                if (badStamps->data.U8[j]) {
    292296                    continue;
    293297                }
     298                pmSubtractionStamp *stamp = stamps->stamps->data[j]; // Stamp of interest
    294299                accumulateCross(&sumI, &sumII, &sumIC, stamp, inputs->data[j], i, footprint);
    295300            }
     
    301306            double chi2 = 0.0;          // Chi^2
    302307            for (int j = 0; j < numStamps; j++) {
    303                 pmSubtractionStamp *stamp = stamps->data[j]; // Stamp of interest
    304                 if (stamp->status == PM_SUBTRACTION_STAMP_REJECTED ||
    305                     stamp->status == PM_SUBTRACTION_STAMP_NONE) {
     308                if (badStamps->data.U8[j]) {
    306309                    continue;
    307310                }
     311                pmSubtractionStamp *stamp = stamps->stamps->data[j]; // Stamp of interest
    308312                chi2 += accumulateChi2(inputs->data[j], stamp, i, coeff, bg, footprint);
    309313            }
     
    329333            psFree(ranking);
    330334            psFree(kernels);
     335            psFree(badStamps);
    331336            return NULL;
    332337        }
     
    336341        // Remove its contribution, and don't include it in the future.
    337342        for (int j = 0; j < numStamps; j++) {
    338             pmSubtractionStamp *stamp = stamps->data[j]; // Stamp of interest
    339             if (stamp->status == PM_SUBTRACTION_STAMP_REJECTED ||
    340                 stamp->status == PM_SUBTRACTION_STAMP_NONE) {
     343            if (badStamps->data.U8[j]) {
    341344                continue;
    342345            }
     346            pmSubtractionStamp *stamp = stamps->stamps->data[j]; // Stamp of interest
    343347            subtractConvolution(inputs->data[j], stamp, bestIndex, bestCoeff, bestBG, footprint);
    344348        }
     
    365369        psFree(ranking);
    366370        psFree(kernels);
     371        psFree(badStamps);
    367372        return NULL;
    368373    }
     
    376381    psArray *convNew = psArrayAlloc(numStamps);
    377382    for (int i = 0; i < numStamps; i++) {
     383        if (badStamps->data.U8[i]) {
     384            continue;
     385        }
    378386        convNew->data[i] = psArrayAlloc(newSize);
    379387    }
     
    388396
    389397            for (int j = 0; j < numStamps; j++) {
     398                if (badStamps->data.U8[j]) {
     399                    continue;
     400                }
     401                pmSubtractionStamp *stamp = stamps->stamps->data[j]; // Stamp of interest
    390402                psArray *convolutions = convNew->data[j]; // Convolutions for this stamp
    391                 pmSubtractionStamp *stamp = stamps->data[j]; // Stamp of interest
    392403                convolutions->data[rank] = psMemIncrRefCounter(stamp->convolutions->data[i]);
    393404            }
     
    405416
    406417    for (int i = 0; i < numStamps; i++) {
    407         pmSubtractionStamp *stamp = stamps->data[i]; // Stamp of interest
     418        if (badStamps->data.U8[i]) {
     419            continue;
     420        }
     421        pmSubtractionStamp *stamp = stamps->stamps->data[i]; // Stamp of interest
    408422        psFree(stamp->convolutions);
    409423        stamp->convolutions = convNew->data[i];
    410424    }
    411425
     426    psFree(badStamps);
    412427    psFree(ranking);
    413428
Note: See TracChangeset for help on using the changeset viewer.