IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Dec 7, 2005, 2:35:14 PM (21 years ago)
Author:
Paul Price
Message:

Updating to use pslib-0.9.0; replaced psFitsAlloc with psFitsOpen

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/pois/src/poisCalculateDeviations.c

    r5717 r5742  
    66#define MAXCHAR 80
    77
    8 psVector *poisCalculateDeviations(psVector *deviations, // Output array of deviations, or NULL
    9                                   psArray *stamps, // Array of stamps
    10                                   psImage *refImage, // Reference image
    11                                   psImage *inImage, // Input image
    12                                   psImage *mask, // Mask image
    13                                   psArray *kernelParams, // Array of kernel parameters
    14                                   psVector *solution, // Solution vector
    15                                   poisConfig *config // Configuration
     8psVector *poisCalculateDeviations(psVector *deviations, // Output array of deviations, or NULL
     9                                  psArray *stamps, // Array of stamps
     10                                  psImage *refImage, // Reference image
     11                                  psImage *inImage, // Input image
     12                                  psImage *mask, // Mask image
     13                                  psArray *kernelParams, // Array of kernel parameters
     14                                  psVector *solution, // Solution vector
     15                                  poisConfig *config // Configuration
    1616    )
    1717{
     
    3232
    3333    if (!deviations) {
    34         deviations = psVectorAlloc(stamps->n, PS_TYPE_F32);
     34        deviations = psVectorAlloc(stamps->n, PS_TYPE_F32);
    3535    }
    3636
    37     int footprint = config->footprint;  // Size of stamp footprint
     37    int footprint = config->footprint;  // Size of stamp footprint
    3838    int xSize = footprint + config->xKernel; // (Half) size of subimage in x
    3939    int ySize = footprint + config->yKernel; // (Half) size of subimage in y
    40     psStats *stats = psStatsAlloc(PS_STAT_SAMPLE_MEAN); // Statistics
     40    psStats *stats = psStatsAlloc(PS_STAT_SAMPLE_MEAN); // Statistics
    4141
    4242    psImage *subStamp = psImageAlloc(2 * xSize, 2 * ySize, PS_TYPE_F32); // Subtraction of stamp
    4343    for (int s = 0; s < stamps->n; s++) {
    44         poisStamp *stamp = stamps->data[s]; // The coordinates of the stamp of interest
    45         int x = stamp->x;               // Stamp x coord
    46         int y = stamp->y;               // Stamp y coord
    47         if (stamp->status == POIS_STAMP_USED) {
    48             psRegion stampRegion = {x - xSize, x + xSize, y - ySize, y + ySize};
    49             psImage *refStamp = psImageSubset(refImage, stampRegion);
    50             psImage *inStamp = psImageSubset(inImage, stampRegion);
    51             psImage *maskStamp = psImageSubset(mask, stampRegion);
    52             psImage *convRefStamp = poisConvolveImage(refStamp, maskStamp, solution, kernelParams, config);
    53             // Calculate chi^2
    54             (void)psBinaryOp(subStamp, inStamp, "-", convRefStamp);
    55             (void)psBinaryOp(subStamp, subStamp, "*", subStamp);
    56             (void)psBinaryOp(subStamp, subStamp, "/", inStamp);
    57             psRegion stampTrim = { config->xKernel, config->xKernel + 2 * footprint, config->yKernel,
    58                                    config->yKernel + 2 * footprint };
    59             psImage *subStampTrim = psImageSubset(subStamp, stampTrim);
    60             psImage *maskStampTrim = psImageSubset(maskStamp, stampTrim);
    61             // Copy image to workaround bug 305
    62             psImage *tempImage = psImageCopy(NULL, subStampTrim, PS_TYPE_F32);
    63             psImage *tempMask = psImageCopy(NULL, maskStampTrim, PS_TYPE_U8);
    64            
    65             (void)psImageStats(stats, tempImage, tempMask, POIS_MASK_BAD | POIS_MASK_NEAR_BAD);
    66            
    67             psFree(tempImage);
    68             psFree(tempMask);
    69            
    70             deviations->data.F32[s] = sqrtf(stats->sampleMean / 2.0);
    71             psTrace("pois.calculateDeviations", 5, "Deviation of stamp %d (%d,%d) is %f\n", s, x, y,
    72                     deviations->data.F32[s]);
    73            
     44        poisStamp *stamp = stamps->data[s]; // The coordinates of the stamp of interest
     45        int x = stamp->x;               // Stamp x coord
     46        int y = stamp->y;               // Stamp y coord
     47        if (stamp->status == POIS_STAMP_USED) {
     48            psRegion stampRegion = {x - xSize, x + xSize, y - ySize, y + ySize};
     49            psImage *refStamp = psImageSubset(refImage, stampRegion);
     50            psImage *inStamp = psImageSubset(inImage, stampRegion);
     51            psImage *maskStamp = psImageSubset(mask, stampRegion);
     52            psImage *convRefStamp = poisConvolveImage(refStamp, maskStamp, solution, kernelParams, config);
     53            // Calculate chi^2
     54            (void)psBinaryOp(subStamp, inStamp, "-", convRefStamp);
     55            (void)psBinaryOp(subStamp, subStamp, "*", subStamp);
     56            (void)psBinaryOp(subStamp, subStamp, "/", inStamp);
     57            psRegion stampTrim = { config->xKernel, config->xKernel + 2 * footprint, config->yKernel,
     58                                   config->yKernel + 2 * footprint };
     59            psImage *subStampTrim = psImageSubset(subStamp, stampTrim);
     60            psImage *maskStampTrim = psImageSubset(maskStamp, stampTrim);
     61            // Copy image to workaround bug 305
     62            psImage *tempImage = psImageCopy(NULL, subStampTrim, PS_TYPE_F32);
     63            psImage *tempMask = psImageCopy(NULL, maskStampTrim, PS_TYPE_U8);
     64
     65            (void)psImageStats(stats, tempImage, tempMask, POIS_MASK_BAD | POIS_MASK_NEAR_BAD);
     66
     67            psFree(tempImage);
     68            psFree(tempMask);
     69
     70            deviations->data.F32[s] = sqrtf(stats->sampleMean / 2.0);
     71            psTrace("pois.calculateDeviations", 5, "Deviation of stamp %d (%d,%d) is %f\n", s, x, y,
     72                    deviations->data.F32[s]);
     73
    7474#ifdef TESTING
    75             char stampName[MAXCHAR];            // File name for stamp
    76             snprintf(stampName, MAXCHAR, "stamp%d.fits", s);
    77             psFits *stampFile = psFitsAlloc(stampName);
    78             if (!psFitsWriteImage(stampFile, NULL, subStampTrim, 0)) {
    79                 psErrorStackPrint(stderr, "Unable to write stamp: %s\n", stampName);
    80             }
    81             psFree(stampFile);
     75            char stampName[MAXCHAR];            // File name for stamp
     76            snprintf(stampName, MAXCHAR, "stamp%d.fits", s);
     77            psFits *stampFile = psFitsOpen(stampName, "w");
     78            if (!psFitsWriteImage(stampFile, NULL, subStampTrim, 0)) {
     79                psErrorStackPrint(stderr, "Unable to write stamp: %s\n", stampName);
     80            }
     81            psFitsClose(stampFile);
    8282#endif
    8383
    8484#if 0
    85             psFree(convRefStamp);
    86             psFree(maskStampTrim);
    87             psFree(subStampTrim);
    88             psFree(maskStamp);
    89             psFree(refStamp);
    90             psFree(inStamp);
     85            psFree(convRefStamp);
     86            psFree(maskStampTrim);
     87            psFree(subStampTrim);
     88            psFree(maskStamp);
     89            psFree(refStamp);
     90            psFree(inStamp);
    9191#endif
    92         }
     92        }
    9393    }
    9494
Note: See TracChangeset for help on using the changeset viewer.