IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 2810


Ignore:
Timestamp:
Dec 23, 2004, 10:27:51 AM (22 years ago)
Author:
gusciora
Message:

Fixed a few bugs in pmSubtractSky.c. Added tests for various input parameter
error conditions.

Location:
trunk/psModules
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/psModules/src/pmReadoutCombine.c

    r2777 r2810  
    55 *  @author GLG, MHPCC
    66 *
    7  *  @version $Revision: 1.9 $ $Name: not supported by cvs2svn $
    8  *  @date $Date: 2004-12-21 20:41:29 $
     7 *  @version $Revision: 1.10 $ $Name: not supported by cvs2svn $
     8 *  @date $Date: 2004-12-23 20:27:51 $
    99 *
    1010 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii

  • trunk/psModules/src/pmSubtractSky.c

    r2777 r2810  
    66 *  @author GLG, MHPCC
    77 *
    8  *  @version $Revision: 1.10 $ $Name: not supported by cvs2svn $
    9  *  @date $Date: 2004-12-21 20:41:29 $
     8 *  @version $Revision: 1.11 $ $Name: not supported by cvs2svn $
     9 *  @date $Date: 2004-12-23 20:27:51 $
    1010 *
    1111 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    3939
    4040/******************************************************************************
    41 getHighestPriorityStatOption(statOptions): this routine takes as in put a
     41getHighestPriorityStatOption(statOptions): this routine takes as input a
    4242psStats->options with multiple options set and returns one with a single
    4343option set according to the precedence set in the SDRS.
    … …  
    4545psU64 getHighestPriorityStatOption(psU64 statOptions)
    4646{
    47     switch (statOptions) {
    48     case PS_STAT_SAMPLE_MEAN:
     47    if (statOptions & PS_STAT_SAMPLE_MEAN) {
    4948        return(PS_STAT_SAMPLE_MEAN);
    50     case PS_STAT_SAMPLE_MEDIAN:
     49    } else if (statOptions & PS_STAT_SAMPLE_MEDIAN) {
    5150        return(PS_STAT_SAMPLE_MEDIAN);
    52     case PS_STAT_CLIPPED_MEAN:
     51    } else if (statOptions & PS_STAT_CLIPPED_MEAN) {
    5352        return(PS_STAT_CLIPPED_MEAN);
    54     case PS_STAT_ROBUST_MEAN:
     53    } else if (statOptions & PS_STAT_ROBUST_MEAN) {
    5554        return(PS_STAT_ROBUST_MEAN);
    56     case PS_STAT_ROBUST_MEDIAN:
     55    } else if (statOptions & PS_STAT_ROBUST_MEDIAN) {
    5756        return(PS_STAT_ROBUST_MEDIAN);
    58     case PS_STAT_ROBUST_MODE:
     57    } else if (statOptions & PS_STAT_ROBUST_MODE) {
    5958        return(PS_STAT_ROBUST_MODE);
    6059    }
    61     //    psError(XXX);
     60    psError(PS_ERR_UNKNOWN, true, "Unallowable option requested for statistically binning image pixels.\n");
    6261    return(-1);
    6362}
    … …  
    129128/******************************************************************************
    130129CalculatePolyTerms(xOrder, yOrder): this routine will calculate the number of
    131 coefficients (or terms) in a polynomial of order (xOrder, yOrder).
     130coefficients (or terms) in a 2-D polynomial of order (xOrder, yOrder).
    132131 *****************************************************************************/
    133132psS32 CalculatePolyTerms(psS32 xOrder, psS32 yOrder)
    … …  
    187186    }
    188187
    189     // XXX: Get rid of this.
    190     for (i=0; i < localPolyTerms ; i++) {
    191         psTrace(".psModule.pmSubtractSky.buildPolyTerms", 6,
    192                 "x^%d * y^%d\n", polyTerms[i][0], polyTerms[i][1]);
     188    if (psTraceGetLevel(".psModule.pmSubtractSky.buildPolyTerms") >= 10) {
     189        for (i=0; i < localPolyTerms ; i++) {
     190            printf("x^%d * y^%d\n", polyTerms[i][0], polyTerms[i][1]);
     191        }
    193192    }
    194193
    … …  
    198197/******************************************************************************
    199198This procedure calculates various combinations of powers of x and y and stores
    200 them in the data structure sums[][].  After it completes:
     199them in the data structure p_psPolySums[][].  After it completes:
    201200 
    202     sums[i][j] == x^i * y^j
     201    p_psPolySums[i][j] == x^i * y^j
    203202 
    204 XXX: Use a psImage for the sums data structure?
    205 XXX: Check for non-NULL sums argument?
    206 XXX: Check for positive x- and yOrder.
    207 XXX: Use variable size arrays for polynomial sums.
     203XXX: Use a psImage for the p_psPolySums data structure?
     204XXX: Check for positive x- and y-Order.
     205XXX: Use variable size arrays for p_psPolySums[][].
     206XXX" Must initialize p_psPolySums[][]?
    208207 *****************************************************************************/
    209208#define PS_MAX_POLYNOMIAL_ORDER 20
     209psF64 p_psPolySums[PS_MAX_POLYNOMIAL_ORDER+1][PS_MAX_POLYNOMIAL_ORDER+1];
    210210void buildSums(psF64 x,
    211211               psF64 y,
    212                psF64 sums[PS_MAX_POLYNOMIAL_ORDER][PS_MAX_POLYNOMIAL_ORDER],
    213212               psS32 xOrder,
    214213               psS32 yOrder)
    … …  
    224223        ySum = xSum;
    225224        for(j=0;j<=yOrder;j++) {
    226             sums[i][j] = ySum;
     225            p_psPolySums[i][j] = ySum;
    227226            ySum*= y;
    228227        }
    … …  
    232231
    233232/******************************************************************************
    234 ImageFitPolynomial(myPoly, binnedImage, maskImage): this private routine takes
     233ImageFitPolynomial(myPoly, dataImage, maskImage): this private routine takes
    235234an input image along with a mask and fits a polynomial to it.  The degree of
    236235the polynomial is specified by input parameter myPoly, and need not be
    … …  
    242241 *****************************************************************************/
    243242psPolynomial2D *ImageFitPolynomial(psPolynomial2D *myPoly,
    244                                    psImage *binnedImage,
     243                                   psImage *dataImage,
    245244                                   psImage *maskImage)
    246245{
    … …  
    249248    PS_POLY_CHECK_NULL(myPoly, NULL);
    250249    PS_POLY_CHECK_TYPE(myPoly, PS_POLYNOMIAL_ORD, NULL);
    251     PS_IMAGE_CHECK_NULL(binnedImage, NULL);
    252     PS_IMAGE_CHECK_EMPTY(binnedImage, NULL);
    253     PS_IMAGE_CHECK_TYPE(binnedImage, PS_TYPE_F32, NULL);
     250    PS_IMAGE_CHECK_NULL(dataImage, NULL);
     251    PS_IMAGE_CHECK_EMPTY(dataImage, NULL);
     252    PS_IMAGE_CHECK_TYPE(dataImage, PS_TYPE_F32, NULL);
    254253    PS_IMAGE_CHECK_NULL(maskImage, NULL);
    255254    PS_IMAGE_CHECK_EMPTY(maskImage, NULL);
    256255    PS_IMAGE_CHECK_TYPE(maskImage, PS_TYPE_U8, NULL);
    257     PS_IMAGE_CHECK_SIZE_EQUAL(binnedImage, maskImage, NULL);
     256    PS_IMAGE_CHECK_SIZE_EQUAL(dataImage, maskImage, NULL);
    258257    psS32 oldPolyX = -1;
    259258    psS32 oldPolyY = -1;
    … …  
    264263    // polynomial if there are not enough rows/columns in the input image.
    265264
    266     if (myPoly->nX > binnedImage->numRows) {
     265    if (myPoly->nX > dataImage->numRows) {
    267266        psLogMsg(__func__, PS_LOG_WARN,
    268267                 "WARNING: ImageFitPolynomial(): Reducing polynomial complexity in x-dimension.\n");
    269268        oldPolyX = myPoly->nX;
    270         myPoly->nX = binnedImage->numRows;
    271     }
    272     if (myPoly->nY > binnedImage->numCols) {
     269        myPoly->nX = dataImage->numRows;
     270    }
     271    if (myPoly->nY > dataImage->numCols) {
    273272        psLogMsg(__func__, PS_LOG_WARN,
    274273                 "WARNING: ImageFitPolynomial(): Reducing polynomial complexity in y-dimension.\n");
    275274        oldPolyY = myPoly->nY;
    276         myPoly->nY = binnedImage->numCols;
     275        myPoly->nY = dataImage->numCols;
    277276    }
    278277    psS32 i;
    … …  
    282281    psS32 aRow;
    283282    psS32 aCol;
    284     // XXX: Document this.  The myPoly->nX and ->nY terms are actually 1 larger
     283    // The myPoly->nX and ->nY terms are actually 1 larger
    285284    // than the order of the polynomial.
    286285    psS32 **polyTerms = buildPolyTerms(myPoly->nX-1, myPoly->nY-1);
     286    // We determine how many coefficients will be in the polynomial that we
     287    // are fitting to this image.
    287288    psS32 localPolyTerms = CalculatePolyTerms(myPoly->nX-1, myPoly->nY-1);
    288     psF64 sums[PS_MAX_POLYNOMIAL_ORDER][PS_MAX_POLYNOMIAL_ORDER];
    289289    psImage *A = psImageAlloc(localPolyTerms, localPolyTerms, PS_TYPE_F64);
    290290    psImage *Aout = psImageAlloc(localPolyTerms, localPolyTerms, PS_TYPE_F64);
    … …  
    295295    // Initialize A matrix and B vector.
    296296    //
    297     for(i=0;i<A->numRows;i++) {
    298         for(j=0;j<A->numCols;j++) {
    299             A->data.F64[i][j] = 0.0;
    300         }
    301     }
    302     for(i=0;i<B->n;i++) {
    303         B->data.F64[i] = 0.0;
    304     }
    305 
     297    PS_IMAGE_SET_F64(A, 0.0);
     298    PS_VECTOR_SET_F64(B, 0.0);
    306299
    307300    //
    308301    // We build the A matrix and B vector.
    309302    //
    310     for (x=0;x<binnedImage->numRows;x++) {
    311         for (y=0;y<binnedImage->numCols;y++) {
     303    for (x=0;x<dataImage->numRows;x++) {
     304        for (y=0;y<dataImage->numCols;y++) {
    312305            if (maskImage->data.U8[x][y] == 0) {
    313                 buildSums((psF64) x, (psF64) y, sums, myPoly->nX-1, myPoly->nY-1);
     306                buildSums((psF64) x, (psF64) y, myPoly->nX-1, myPoly->nY-1);
    314307
    315308                /************************************************************
     309                This code dervies from equation (7) of the pilot ADD.  However,
     310                it is not exactly the same in that the order of the polynomial
     311                may be different in X And Y.
     312
    316313                Equation (7) from the pilot ADD describes 16 linear equations.
    317314                The i-th equation is simply the partial derivative of the
    … …  
    325322                    for (aCol=0;aCol<localPolyTerms;aCol++) {
    326323                        A->data.F64[aRow][aCol]+=
    327                             (sums[ polyTerms[aCol][0] ][ polyTerms[aCol][1] ] *
    328                              sums[ polyTerms[aRow][0] ][ polyTerms[aRow][1] ]);
     324                            (p_psPolySums[ polyTerms[aCol][0] ][ polyTerms[aCol][1] ] *
     325                             p_psPolySums[ polyTerms[aRow][0] ][ polyTerms[aRow][1] ]);
    329326                    }
    330327                }
    331328                // Build the B[] vector, which is the right-hand side of (7).
    332329                for (i=0;i<localPolyTerms;i++) {
    333                     B->data.F64[i]+= binnedImage->data.F32[x][y] *
    334                                      sums[ polyTerms[i][0] ][ polyTerms[i][1] ];
     330                    B->data.F64[i]+= dataImage->data.F32[x][y] *
     331                                     p_psPolySums[ polyTerms[i][0] ][ polyTerms[i][1] ];
    335332                }
    336333            }
    … …  
    338335    }
    339336
    340     // XXX: Put this loop inside a psTrace conditional, somehow.
    341     for (aRow=0;aRow<localPolyTerms;aRow++) {
    342         for (aCol=0;aCol<localPolyTerms;aCol++) {
    343             psTrace(".psModule.pmSubtractSky.ImageFitPolynomial", 8,
    344                     "A[%d][%d] is %f\n", aRow, aCol, A->data.F64[aRow][aCol]);
    345         }
    346     }
    347     // XXX: Put this loop inside a psTrace conditional, somehow.
    348     for (i=0;i<=localPolyTerms;i++) {
    349         psTrace(".psModule.pmSubtractSky.ImageFitPolynomial", 6,
    350                 "B[%d] is %f\n", i, B->data.F64[i]);
     337    if (psTraceGetLevel(".psModule.pmSubtractSky.ImageFitPolynomial") >= 8) {
     338        for (aRow=0;aRow<localPolyTerms;aRow++) {
     339            for (aCol=0;aCol<localPolyTerms;aCol++) {
     340                printf("A[%d][%d] is %f\n", aRow, aCol,
     341                      A->data.F64[aRow][aCol]);
     342            }
     343        }
     344
     345        for (i=0;i<=localPolyTerms;i++) {
     346            printf("B[%d] is %f\n", i, B->data.F64[i]);
     347        }
    351348    }
    352349
    353350    //
    354351    // Solve the matrix equations for the polynomial coefficients C.
     352    // XXX: How do we know if these matrix operations were successful?
    355353    //
    356354    Aout = psMatrixLUD(Aout, outPerm, A);
     355    PS_IMAGE_CHECK_NULL(Aout, NULL);
     356    PS_IMAGE_CHECK_EMPTY(Aout, NULL);
    357357    psVector *C = psVectorAlloc(localPolyTerms, PS_TYPE_F64);
    358358    psMatrixLUSolve(C, Aout, B, outPerm);
    … …  
    382382    //
    383383    // We restore the original size of the polynomial and set remaining
    384     // coefficients to 0.0.
     384    // coefficients to 0.0, if necessary.
    385385    //
    386386    if (oldPolyX != -1) {
    … …  
    424424    psTrace(".psModule.pmSubtractSky", 4,
    425425            "---- pmSubtractSky() begin ----\n");
     426    psStatsOptions statOptions = 0;
    426427
    427428    // Return the original input readout if the fit specs are poorly defined.
    … …  
    439440    // Determine which statistic to use when binning pixels, if any.
    440441    //
    441     psStatsOptions statOptions = stats->options;
    442     if (1 < p_psDetermineNumBits(statOptions)) {
    443         //XXX  psWarning(PS_ERR_UNKNOWN,true, "Multiple statistical options have been requested.\n");
    444         statOptions = getHighestPriorityStatOption(statOptions);
    445         // XXX: Don't modify input parameter.
    446         oldStatOptions = stats->options;
    447         stats->options = statOptions;
    448     }
    449     if (0 >= p_psDetermineNumBits(statOptions)) {
    450         psLogMsg(__func__, PS_LOG_WARN,
    451                  "WARNING: pmSubtractSky(): no stats->options was requested\n");
    452         return(in);
     442    if (stats != NULL) {
     443        statOptions = stats->options;
     444        if (1 < p_psDetermineNumBits(statOptions)) {
     445            psLogMsg(__func__, PS_LOG_WARN, "WARNING: Multiple statistical options have been requested.\n");
     446            statOptions = getHighestPriorityStatOption(statOptions);
     447
     448            // Save old input "stats" parameter.
     449            oldStatOptions = stats->options;
     450            stats->options = statOptions;
     451        }
     452        if (0 == p_psDetermineNumBits(statOptions)) {
     453            psLogMsg(__func__, PS_LOG_WARN,
     454                     "WARNING: pmSubtractSky(): no stats->options was requested\n");
     455        }
    453456    }
    454457
    … …  
    456459    // Generate required warning messages.
    457460    //
    458     if (0 == p_psDetermineNumBits(statOptions)) {
    459         psLogMsg(__func__, PS_LOG_WARN,
    460                  "WARNING: pmSubtractSky(): no stats->options was requested\n");
    461     }
    462461    if (binFactor <= 0) {
    463462        psLogMsg(__func__, PS_LOG_WARN,
    … …  
    558557        }
    559558    } else {
    560         //XXX: error
     559        // We shouldn't get here since we check this above.
     560        psError(PS_ERR_UNKNOWN, true, "Unallowable fit type.");
    561561    }
    562562

  • trunk/psModules/test/tst_pmSubtractSky.c

    r2777 r2810  
    77 *  @author GLG, MHPCC
    88 *
    9  *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
    10  *  @date $Date: 2004-12-21 20:41:29 $
     9 *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
     10 *  @date $Date: 2004-12-23 20:27:51 $
    1111 *
    1212 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
    … …  
    2222#define OBJECT_INTENSITY 2000.0
    2323static int test00(void);
     24static int test01(void);
    2425testDescription tests[] = {
    2526                              {test00, 000, "pmSubtractSky", 0, false},
     27                              {test01, 000, "pmSubtractSky: warning, error messages", 0, false},
    2628                              {NULL}
    2729                          };
    … …  
    4547    int j;
    4648    int testStatus = 0;
    47     psImage *tmpImage1 = psImageAlloc(numCols, numRows, PS_TYPE_F32);
    48     psReadout *myReadout = psReadoutAlloc(numCols, numRows, tmpImage1);
     49    psImage *tmpImageF32 = psImageAlloc(numCols, numRows, PS_TYPE_F32);
     50    psReadout *myReadout = psReadoutAlloc(numCols, numRows, tmpImageF32);
    4951    psStats *myStats = psStatsAlloc(PS_STAT_SAMPLE_MEAN);
    5052    psPolynomial2D *myPoly = psPolynomial2DAlloc(POLY_X_ORDER, POLY_Y_ORDER, PS_POLYNOMIAL_ORD);
    … …  
    7981    int j;
    8082    int testStatus = 0;
    81     psImage *tmpImage1 = psImageAlloc(numCols, numRows, PS_TYPE_F32);
    82     psReadout *myReadout = psReadoutAlloc(numCols, numRows, tmpImage1);
     83    psImage *tmpImageF32 = psImageAlloc(numCols, numRows, PS_TYPE_F32);
     84    psReadout *myReadout = psReadoutAlloc(numCols, numRows, tmpImageF32);
    8385    psStats *myStats = psStatsAlloc(PS_STAT_SAMPLE_MEAN);
    8486    psPolynomial2D *myPoly = psPolynomial2DAlloc(POLY_X_ORDER, POLY_Y_ORDER, PS_POLYNOMIAL_ORD);
    … …  
    169171    return(testStatus);
    170172}
     173
     174#define NUM_ROWS_SMALL 16
     175#define NUM_COLS_SMALL 16
     176int test01( void )
     177{
     178    int testStatus = 0;
     179    psS32 i;
     180    psS32 j;
     181    psImage *tmpImageF32 = psImageAlloc(NUM_COLS_SMALL, NUM_ROWS_SMALL, PS_TYPE_F32);
     182    psImage *tmpImageF64 = psImageAlloc(NUM_COLS_SMALL, NUM_ROWS_SMALL, PS_TYPE_F64);
     183    psReadout *myReadout = psReadoutAlloc(NUM_COLS_SMALL, NUM_ROWS_SMALL, tmpImageF32);
     184    psReadout *rc = NULL;
     185    psStats *myStats = psStatsAlloc(PS_STAT_SAMPLE_MEAN);
     186    psPolynomial2D *myPoly = psPolynomial2DAlloc(POLY_X_ORDER, POLY_Y_ORDER, PS_POLYNOMIAL_ORD);
     187    for (i=0;i<NUM_ROWS_SMALL;i++) {
     188        for (j=0;j<NUM_COLS_SMALL;j++) {
     189            myReadout->image->data.F32[i][j] = func(i, j);
     190        }
     191    }
     192
     193    printf("----------------------------------------------------------------\n");
     194    printf("Calling pmSubtractSky() with NULL psReadout.  Should error.\n\n");
     195    rc = pmSubtractSky(NULL, (void *) myPoly, PM_FIT_POLYNOMIAL, 1, myStats, 2.0);
     196    if (rc != NULL) {
     197        printf("ERROR: pmSubtractSky() returned a non-NULL psReadout\n");
     198        testStatus = false;
     199    }
     200
     201    printf("----------------------------------------------------------------\n");
     202    printf("Calling pmSubtractSky() with NULL psReadout->image.  Should error.\n\n");
     203    myReadout->image = NULL;
     204    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL, 1, myStats, 2.0);
     205    if (rc != NULL) {
     206        printf("ERROR: pmSubtractSky() returned a non-NULL psReadout\n");
     207        testStatus = false;
     208    }
     209    myReadout->image = tmpImageF32;
     210
     211    printf("----------------------------------------------------------------\n");
     212    printf("Calling pmSubtractSky() with PS_TYPE_F64 psReadout->image.  Should error.\n\n");
     213    myReadout->image = tmpImageF64;
     214    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL, 1, myStats, 2.0);
     215    if (rc != NULL) {
     216        printf("ERROR: pmSubtractSky() returned a non-NULL psReadout\n");
     217        testStatus = false;
     218    }
     219    myReadout->image = tmpImageF32;
     220
     221    printf("----------------------------------------------------------------\n");
     222    printf("Calling pmSubtractSky() with NULL fitSpec.  Should return image, no ERROR, no WARNING.\n\n");
     223    rc = pmSubtractSky(myReadout, NULL, PM_FIT_POLYNOMIAL, 1, myStats, 2.0);
     224    if (rc != myReadout) {
     225        printf("ERROR: pmSubtractSky() returned something other than psReadout\n");
     226        testStatus = false;
     227    }
     228
     229    printf("----------------------------------------------------------------\n");
     230    printf("Calling pmSubtractSky() with PM_FIT_NONE fit.  Should return image, no ERROR, no WARNING.\n\n");
     231    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_NONE, 1, myStats, 2.0);
     232    if (rc != myReadout) {
     233        printf("ERROR: pmSubtractSky() returned something other than psReadout\n");
     234        testStatus = false;
     235    }
     236
     237    printf("----------------------------------------------------------------\n");
     238    printf("Calling pmSubtractSky() with PM_FIT_SPLINE fit.  Should return image, no ERROR, no WARNING.\n\n");
     239    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_SPLINE, 1, myStats, 2.0);
     240    if (rc != myReadout) {
     241        printf("ERROR: pmSubtractSky() returned something other than psReadout\n");
     242        testStatus = false;
     243    }
     244
     245    printf("----------------------------------------------------------------\n");
     246    printf("Calling pmSubtractSky() with NULL myStats.  Should fit entire image, generate WARNING.\n\n");
     247    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL, 1, NULL, 2.0);
     248    for (i=0;i<NUM_ROWS_SMALL;i++) {
     249        for (j=0;j<NUM_COLS_SMALL;j++) {
     250            if (ERROR_TOLERANCE < fabs(rc->image->data.F32[i][j])) {
     251                printf("ERROR: image[%d][%d] is %f, should be 0.0\n", i, j,
     252                       rc->image->data.F32[i][j]);
     253                testStatus = false;
     254            }
     255        }
     256    }
     257
     258    printf("----------------------------------------------------------------\n");
     259    psU64 oldOptions = myStats->options;
     260    myStats->options = 0;
     261    printf("Calling pmSubtractSky() with no myStats->options specified.  Should fit entire image, generate WARNING.\n\n");
     262    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL, 1, myStats, 2.0);
     263    for (i=0;i<NUM_ROWS_SMALL;i++) {
     264        for (j=0;j<NUM_COLS_SMALL;j++) {
     265            if (ERROR_TOLERANCE < fabs(rc->image->data.F32[i][j])) {
     266                printf("ERROR: image[%d][%d] is %f, should be 0.0\n", i, j,
     267                       rc->image->data.F32[i][j]);
     268                testStatus = false;
     269            }
     270        }
     271    }
     272    myStats->options = oldOptions;
     273
     274    printf("----------------------------------------------------------------\n");
     275    oldOptions = myStats->options;
     276    myStats->options = PS_STAT_SAMPLE_MEAN | PS_STAT_SAMPLE_MEDIAN;
     277    printf("Calling pmSubtractSky() with multiple myStats->options specified.  Should fit entire image, generate WARNING.\n\n");
     278    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL, 1, myStats, 2.0);
     279    for (i=0;i<NUM_ROWS_SMALL;i++) {
     280        for (j=0;j<NUM_COLS_SMALL;j++) {
     281            if (ERROR_TOLERANCE < fabs(rc->image->data.F32[i][j])) {
     282                printf("ERROR: image[%d][%d] is %f, should be 0.0\n", i, j,
     283                       rc->image->data.F32[i][j]);
     284                testStatus = false;
     285            }
     286        }
     287    }
     288    myStats->options = oldOptions;
     289
     290    printf("----------------------------------------------------------------\n");
     291    printf("Calling pmSubtractSky() with 0 binFactor.  Should fit entire image, generate WARNING.\n\n");
     292    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL, 0, myStats, 2.0);
     293    for (i=0;i<NUM_ROWS_SMALL;i++) {
     294        for (j=0;j<NUM_COLS_SMALL;j++) {
     295            if (ERROR_TOLERANCE < fabs(rc->image->data.F32[i][j])) {
     296                printf("ERROR: image[%d][%d] is %f, should be 0.0\n", i, j,
     297                       rc->image->data.F32[i][j]);
     298                testStatus = false;
     299            }
     300        }
     301    }
     302
     303    printf("----------------------------------------------------------------\n");
     304    printf("Calling pmSubtractSky() with -1 binFactor.  Should fit entire image, generate WARNING.\n\n");
     305    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL, -1, myStats, 2.0);
     306    for (i=0;i<NUM_ROWS_SMALL;i++) {
     307        for (j=0;j<NUM_COLS_SMALL;j++) {
     308            if (ERROR_TOLERANCE < fabs(rc->image->data.F32[i][j])) {
     309                printf("ERROR: image[%d][%d] is %f, should be 0.0\n", i, j,
     310                       rc->image->data.F32[i][j]);
     311                testStatus = false;
     312            }
     313        }
     314    }
     315
     316    printf("----------------------------------------------------------------\n");
     317    printf("Calling pmSubtractSky() with -1.0 clipSD.  Should fit entire image, generate WARNING.\n\n");
     318    rc = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL, 1, myStats, -1.0);
     319    for (i=0;i<NUM_ROWS_SMALL;i++) {
     320        for (j=0;j<NUM_COLS_SMALL;j++) {
     321            if (ERROR_TOLERANCE < fabs(rc->image->data.F32[i][j])) {
     322                printf("ERROR: image[%d][%d] is %f, should be 0.0\n", i, j,
     323                       rc->image->data.F32[i][j]);
     324                testStatus = false;
     325            }
     326        }
     327    }
     328
     329
     330    //    myReadout = pmSubtractSky(myReadout, (void *) myPoly, PM_FIT_POLYNOMIAL,
     331    //                              1, myStats, 2.0);
     332
     333
     334    printf("----------------------------------------------------------------\n");
     335    psFree(myReadout);
     336    psFree(myStats);
     337    psFree(myPoly);
     338    psFree(tmpImageF64);
     339    return(testStatus);
     340}
     341// This code will
Note: See TracChangeset for help on using the changeset viewer.