/** @file tst_pmFindObjects.c
 *
 *  @brief Contains the tests for pmSubtractSky.c:
 *
 * test00: This code will ...
 *
 *  @author GLG, MHPCC
 *
 * XXX: Must test
 *       pmSourceRoughClass
 *
 * XXX: Must test output results for many other functions.
 *
 * XXX: There are many cases where row/col can be switched in the code.
 * We must test that here by using non-square images.  All tests
 * in this file should be run with non-square images.
 *
 * XXX: Memory leaks are not being caught.  If I allocated a psVector in these functions
 * abd never deallocate, no error is generated.
 *
 *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2005-04-01 21:39:06 $
 *
 *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 */
#include "psTest.h"
#include "pslib.h"
#include "pmObjects.h"
#define NUM_ROWS 10
#define NUM_COLS 10
#define ERROR_TOLERANCE 1.0
static int test00(void);
static int test01(void);
static int test02(void);
static int test03(void);
static int test04(void);
static int test05(void);
static int test06(void);
static int test07(void);
static int test08(void);
static int test09(void);
static int test15(void);
static int test16(void);
static int test20(void);
testDescription tests[] = {
                              {test00, 000, "pmObjects: structure allocators and deallocators", true, false},
                              {test01, 001, "pmObjects: psFindVectorPeaks()", true, false},
                              {test02, 001, "pmObjects: psFindImagePeaks()", true, false},
                              {test03, 001, "pmObjects: pmCullPeaks()", true, false},
                              {test04, 001, "pmObjects: pmSourceLocalSky()", true, false},
                              {test06, 001, "pmObjects: pmSourceSetPixelCircle()", true, false},
                              {test05, 001, "pmObjects: pmSourceMoments()", true, false},
                              {test07, 001, "pmObjects: pmMin()", true, false},
                              {test08, 001, "pmObjects: pmSourceModelGuess()", true, false},
                              {test09, 001, "pmObjects: pmSourceContour()", true, false},
                              {test15, 001, "pmObjects: pmSourceAddModel()", true, false},
                              {test16, 001, "pmObjects: pmSourceSubModel()", true, false},
                              {test20, 001, "pmObjects: pmSourceSubModel()", true, false},
                              {NULL}
                          };

int main(int argc, char* argv[])
{
    return !runTestSuite(stderr, "Test Point Driver", tests, argc, argv);
    test02();
}

/******************************************************************************
test00(): Test the various allocators and deallocators.
 *****************************************************************************/
int test00( void )
{
    bool testStatus = true;
    psTraceSetLevel(".", 0);

    psPeak *tmpPeak = pmPeakAlloc(0, 0, 0.0, PM_PEAK_UNDEF);
    if (tmpPeak == NULL) {
        printf("TEST ERROR: pmPeakAlloc() returned a NULL psPeak\n");
        testStatus = false;
    }
    // XXX: Test that values were set properly.
    psFree(tmpPeak);

    psMoments *tmpMoments = pmMomentsAlloc();
    if (tmpMoments == NULL) {
        printf("TEST ERROR: pmMomentsAlloc() returned a NULL psMoments\n");
        testStatus = false;
    }
    psFree(tmpMoments);

    // XXX: Must test this with variety of models, and verify param length.
    psModel *tmpModel = pmModelAlloc(PS_MODEL_GAUSS);
    if (tmpModel == NULL) {
        printf("TEST ERROR: pmModelAlloc() returned a NULL psModel\n");
        testStatus = false;
    }
    // XXX: Test that values were set properly.
    psFree(tmpModel);

    // XXX: Must test this with variety of models, and verify param length.
    tmpModel = pmModelAlloc(PS_MODEL_PGAUSS);
    if (tmpModel == NULL) {
        printf("TEST ERROR: pmModelAlloc() returned a NULL psModel\n");
        testStatus = false;
    }
    // XXX: Test that values were set properly.
    psFree(tmpModel);

    // XXX: Must test this with variety of models, and verify param length.
    tmpModel = pmModelAlloc(PS_MODEL_TWIST_GAUSS);
    if (tmpModel == NULL) {
        printf("TEST ERROR: pmModelAlloc() returned a NULL psModel\n");
        testStatus = false;
    }
    // XXX: Test that values were set properly.
    psFree(tmpModel);

    // XXX: Must test this with variety of models, and verify param length.
    tmpModel = pmModelAlloc(PS_MODEL_WAUSS);
    if (tmpModel == NULL) {
        printf("TEST ERROR: pmModelAlloc() returned a NULL psModel\n");
        testStatus = false;
    }
    // XXX: Test that values were set properly.
    psFree(tmpModel);

    // XXX: Must test this with variety of models, and verify param length.
    tmpModel = pmModelAlloc(PS_MODEL_SERSIC);
    if (tmpModel == NULL) {
        printf("TEST ERROR: pmModelAlloc() returned a NULL psModel\n");
        testStatus = false;
    }
    // XXX: Test that values were set properly.
    psFree(tmpModel);

    // XXX: Must test this with variety of models, and verify param length.
    tmpModel = pmModelAlloc(PS_MODEL_SERSIC_CORE);
    if (tmpModel == NULL) {
        printf("TEST ERROR: pmModelAlloc() returned a NULL psModel\n");
        testStatus = false;
    }
    // XXX: Test that values were set properly.
    psFree(tmpModel);

    psSource *tmpSource = pmSourceAlloc();
    if (tmpSource == NULL) {
        printf("TEST ERROR: pmSourceAlloc() returned a NULL psSource\n");
        testStatus = false;
    }
    psFree(tmpSource);

    return(testStatus);
}

#define VECTOR_SIZE 10
bool test_pmFindVectorPeaks(int n)
{
    bool testStatus = true;
    psVector *inData = psVectorAlloc(n, PS_TYPE_F32);
    psVector *outData = NULL;

    printf("-------------- Calling test_pmFindVectorPeaks on an %d size vector. --------------\n", n);
    //
    // Test first pixel peak.
    //
    printf("Test pmFindVectorPeaks() with a first-element peak.\n");
    for (psS32 i = 0 ; i < n ; i++) {
        inData->data.F32[i] = (float) (n-i);
    }
    inData->data.F32[0] = (float) n;

    outData= pmFindVectorPeaks(inData, 0.0);
    if (outData == NULL) {
        printf("TEST ERROR: pmFindVectorPeaks returned a NULL psVector.\n");
        testStatus = false;
    } else {

        if (outData->n != 1) {
            printf("TEST ERROR: outData->n is %d\n", outData->n);
            testStatus = false;
        }
        if (outData->data.U32[0] != 0) {
            printf("TEST ERROR: Did not find peak at element 0.\n");
            testStatus = false;
        }
        psFree(outData);
    }

    //
    // Test first pixel peak, large threshold
    //
    printf("Test pmFindVectorPeaks() with a first-element peak, large threshold.\n");
    for (psS32 i = 0 ; i < n ; i++) {
        inData->data.F32[i] = (float) (n-i);
    }
    inData->data.F32[0] = (float) n;

    outData= pmFindVectorPeaks(inData, (float) (n*n));
    if (outData == NULL) {
        printf("TEST ERROR: pmFindVectorPeaks returned a NULL psVector.\n");
        testStatus = false;
    } else {

        if (outData->n != 0) {
            printf("TEST ERROR: outData->n is %d\n", outData->n);
            testStatus = false;
        }
        psFree(outData);

        // Skip remaining tests if the input vector has length 1.
        if (n == 1) {
            psFree(inData);
            return(testStatus);
        }
    }

    //
    // Test last pixel peak.
    //
    printf("Test pmFindVectorPeaks() with a last-element peak.\n");
    for (psS32 i = 0 ; i < n ; i++) {
        inData->data.F32[i] = (float) (i);
    }
    inData->data.F32[n-1] = (float) n;

    outData= pmFindVectorPeaks(inData, 0.0);
    if (outData == NULL) {
        printf("TEST ERROR: pmFindVectorPeaks returned a NULL psVector.\n");
        testStatus = false;
    } else {

        if (outData->n != 1) {
            printf("TEST ERROR: outData->n is %d\n", outData->n);
            testStatus = false;
        }
        if (outData->data.U32[0] != n-1) {
            printf("TEST ERROR: Did not find peak at element %d.\n", n-1);
            testStatus = false;
        }
        psFree(outData);
    }

    //
    // Test last pixel peak, large threshold.
    //
    printf("Test pmFindVectorPeaks() with a last-element peak, large threshold.\n");
    for (psS32 i = 0 ; i < n ; i++) {
        inData->data.F32[i] = (float) (i);
    }
    inData->data.F32[n-1] = (float) n;

    outData= pmFindVectorPeaks(inData, (float) (n*n));
    if (outData == NULL) {
        printf("TEST ERROR: pmFindVectorPeaks returned a NULL psVector.\n");
        testStatus = false;
    } else {

        if (outData->n != 0) {
            printf("TEST ERROR: outData->n is %d\n", outData->n);
            testStatus = false;
        }
        psFree(outData);
    }

    //
    // Test interior peaks.
    // Set all even number elements to be peaks.
    //
    printf("Test pmFindVectorPeaks() with all even-numbered elements peak.\n");
    for (psS32 i = 0 ; i < n ; i++) {
        if (0 == i%2) {
            inData->data.F32[i] = (float) (2 * i);
        } else {
            inData->data.F32[i] = (float) (i);
        }
    }
    inData->data.F32[0] = (float) n;


    outData= pmFindVectorPeaks(inData, 0.0);
    if (outData == NULL) {
        printf("TEST ERROR: pmFindVectorPeaks returned a NULL psVector.\n");
        testStatus = false;
    } else {

        if (outData->n != n/2) {
            printf("TEST ERROR: outData->n is %d\n", outData->n);
            testStatus = false;
        }

        for (psS32 i = 0 ; i < outData->n ; i++) {
            if (outData->data.U32[i] != (2 * i)) {
                printf("TEST ERROR: the %d-th peak is element number %d\n", i, outData->data.U32[i]);
                testStatus = false;
            }
        }
        psFree(outData);
    }

    //
    // Test interior peaks, with threshold = n*n.
    // Should generate an empty output psVector.
    //
    printf("Test pmFindVectorPeaks() with all even-numbered elements peak, large threshold.\n");
    outData= pmFindVectorPeaks(inData, (float) (n*n));
    if (outData == NULL) {
        printf("TEST ERROR: pmFindVectorPeaks returned a NULL psVector.\n");
        testStatus = false;
    } else {

        if (outData->n != 0) {
            printf("TEST ERROR: outData->n is %d\n", outData->n);
            testStatus = false;
        }
        psFree(outData);
    }


    psFree(inData);
    return(testStatus);
}


int test01( void )
{
    bool testStatus = true;
    psVector *tmpVec = NULL;
    psVector *tmpVecF64 = psVectorAlloc(10, PS_TYPE_F64);
    psVector *tmpVecEmpty = psVectorAlloc(0, PS_TYPE_F32);

    psTraceSetLevel(".", 0);

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmFindVectorPeaks with NULL psVector.  Should generate error and return NULL.\n");
    tmpVec = pmFindVectorPeaks(NULL, 0.0);
    if (tmpVec != NULL) {
        printf("TEST ERROR: pmFindVectorPeaks() returned a non-NULL psVector.\n");
        testStatus = false;
        psFree(tmpVec);
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmFindVectorPeaks with empty psVector.  Should generate error and return NULL.\n");
    tmpVec = pmFindVectorPeaks(tmpVecEmpty, 0.0);
    if (tmpVec != NULL) {
        printf("TEST ERROR: pmFindVectorPeaks() returned a non-NULL psVector.\n");
        testStatus = false;
        psFree(tmpVec);
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmFindVectorPeaks with PS_TYPE_F64 psVector.  Should generate error and return NULL.\n");
    tmpVec = pmFindVectorPeaks(tmpVecF64, 0.0);
    if (tmpVec != NULL) {
        printf("TEST ERROR: pmFindVectorPeaks() returned a non-NULL psVector.\n");
        testStatus = false;
        psFree(tmpVec);
    }
    testStatus&= test_pmFindVectorPeaks(1);
    testStatus&= test_pmFindVectorPeaks(VECTOR_SIZE);

    psFree(tmpVecF64);
    psFree(tmpVecEmpty);
    return(testStatus);
}

// XXX: Add tests for the PEAK_TYPE.
// XXX: Add interior peaks, edge peaks, and flat peaks.
bool test_pmFindImagePeaks(int numRows, int numCols)
{
    printf("-------------- Calling test_pmFindVectorPeaks on an %d-by-%d image. --------------\n", numRows, numCols);
    bool testStatus = true;
    psImage *inData = psImageAlloc(numCols, numRows, PS_TYPE_F32);
    psList *outData = NULL;

    //
    // Test first pixel peak.
    //
    for (psS32 i = 0 ; i < numRows ; i++) {
        for (psS32 j = 0 ; j < numCols ; j++) {
            inData->data.F32[i][j] = PS_SQR(i - numRows/2) + PS_SQR(j-numCols/2);
        }
    }
    inData->data.F32[0][0] = PS_SQR(numRows) + PS_SQR(numCols);
    inData->data.F32[0][numCols-1] = PS_SQR(numRows) + PS_SQR(numCols);
    inData->data.F32[numRows-1][0] = PS_SQR(numRows) + PS_SQR(numCols);
    inData->data.F32[numRows-1][numCols-1] = PS_SQR(numRows) + PS_SQR(numCols);
    for (psS32 i = 0 ; i < numRows ; i++) {
        for (psS32 j = 0 ; j < numCols ; j++) {
            printf("(%.1f) ", inData->data.F32[i][j]);
        }
        printf("\n");
    }

    outData = pmFindImagePeaks(inData, 0.0);

    if (outData == NULL) {
        printf("TEST ERROR: pmFindImagePeaks returned a NULL psList.\n");
        testStatus = false;
    } else {
        if (outData->size != 4) {
            printf("TEST ERROR: pmFindImagePeaks found only %d peaks (should be 4)\n", outData->size);
            testStatus = false;
        }

        psListElem *tmpPeakLE = (psListElem *) outData->head;
        while (tmpPeakLE != NULL) {
            psPeak *tmpPeak = (psPeak *) tmpPeakLE->data;
            if (((tmpPeak->x == 0) || (tmpPeak->x == numRows-1)) &&
                ((tmpPeak->y == 0) || (tmpPeak->y == numCols-1))) {}
            else {
                printf("TEST ERROR: Peak at (%d, %d) (%f)\n", tmpPeak->x, tmpPeak->y, tmpPeak->counts);
                testStatus = false;
            }
            tmpPeakLE = tmpPeakLE->next;
        }
    }

    psFree(inData);
    psFree(outData);
    return(testStatus);
}


int test02( void )
{
    bool testStatus = true;
    psList *tmpList = NULL;
    psImage *tmpImageF64 = psImageAlloc(10, 10, PS_TYPE_F64);
    psImage *tmpImageEmpty = psImageAlloc(0, 0, PS_TYPE_F32);

    psTraceSetLevel(".", 0);

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmFindImagePeaks with NULL psImage.  Should generate error and return NULL.\n");
    tmpList = pmFindImagePeaks(NULL, 0.0);
    if (tmpList != NULL) {
        printf("TEST ERROR: pmFindImagePeaks() returned a non-NULL psImage.\n");
        testStatus = false;
        psFree(tmpList);
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmFindImagePeaks with empty psImage.  Should generate error and return NULL.\n");
    tmpList = pmFindImagePeaks(tmpImageEmpty, 0.0);
    if (tmpList != NULL) {
        printf("TEST ERROR: pmFindImagePeaks() returned a non-NULL psImage.\n");
        testStatus = false;
        psFree(tmpList);
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmFindImagePeaks with PS_TYPE_F64 psImage.  Should generate error and return NULL.\n");
    tmpList = pmFindImagePeaks(tmpImageF64, 0.0);
    if (tmpList != NULL) {
        printf("TEST ERROR: pmFindImagePeaks() returned a non-NULL psImage.\n");
        testStatus = false;
        psFree(tmpList);
    }
    printf("----------------------------------------------------------------------------------\n");
    //    testStatus&= test_pmFindImagePeaks(1, 1);
    //    testStatus&= test_pmFindImagePeaks(2, 5);
    //    testStatus&= test_pmFindImagePeaks(5, 2);
    testStatus&= test_pmFindImagePeaks(10, 10);


    psFree(tmpImageF64);
    psFree(tmpImageEmpty);
    return(testStatus);
}

/******************************************************************************
test03(): We first test pmCullPeaks() with various NULL and unallowable input
parameters.  Then we generate a list of peaks and test that pmCullPeaks()
removes them correctly.
 *****************************************************************************/
int test03( void )
{
    bool testStatus = true;
    psImage *imgData = psImageAlloc(NUM_COLS, NUM_ROWS, PS_TYPE_F32);
    psList *outData = NULL;

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmCullPeaks with NULL psList.  Should generate error and return NULL.\n");
    outData = pmCullPeaks(NULL, 0.0, NULL);
    if (outData != NULL) {
        printf("TEST ERROR: pmCulPeaks() returned a non-NULL psList.\n");
        testStatus = false;
    }

    //
    // Set peaks in input image.  All even-column and even-row pixels are
    // set non-zero, all other pixels are set to zero.
    //
    psS32 numPeaksOrig = 0;
    for (psS32 i = 0 ; i < NUM_ROWS ; i++) {
        for (psS32 j = 0 ; j < NUM_COLS ; j++) {
            if ((0 == i%2) && (0 == j%2)) {
                imgData->data.F32[i][j] = (float) (i + 10);
                numPeaksOrig++;
            } else {
                imgData->data.F32[i][j] = 0.0;
            }
        }
    }
    for (psS32 i = 0 ; i < NUM_ROWS ; i++) {
        for (psS32 j = 0 ; j < NUM_COLS ; j++) {
            printf("(%.1f) ", imgData->data.F32[i][j]);
        }
        printf("\n");
    }
    printf("Set %d peaks\n", numPeaksOrig);

    //
    // Call pmCullPeaks() with HUGE maxValue and NULL psRegion.  Should not
    // remove any peaks.
    //
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmCullPeaks with large maxValue and NULL psRegion.\n");
    outData = pmFindImagePeaks(imgData, 0.0);
    outData = pmCullPeaks(outData, PS_MAX_F32, NULL);

    if (outData == NULL) {
        printf("TEST ERROR: pmCullPeaks() returned a non-NULL psList.\n");
        testStatus = false;
        return(testStatus);
    }
    if (outData->size != numPeaksOrig) {
        printf("TEST ERROR (0): pmCullPeaks incorrectly removed peaks\n");
        printf("The pmCullPeaks() output had %d peaks, should have had %d peaks.\n", outData->size, numPeaksOrig);
        testStatus = false;
    }
    psFree(outData);

    //
    // Call pmCullPeaks() with TINY maxValue and NULL psRegion.  Should
    // remove all peaks.
    //
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmCullPeaks with tiny maxValue and NULL psRegion.\n");
    outData = pmFindImagePeaks(imgData, 0.0);
    printf("pmFindImagePeaks found %d peaks\n", outData->size);
    outData = pmCullPeaks(outData, 0.0, NULL);

    if (outData == NULL) {
        printf("TEST ERROR: pmCullPeaks() returned a non-NULL psList.\n");
        testStatus = false;
        return(testStatus);
    }
    if (outData->size != 0) {
        printf("TEST ERROR (1): pmCullPeaks incorrectly removed peaks\n");
        printf("The pmCullPeaks() output had %d peaks, should have had %d peaks.\n", outData->size, 0);
        testStatus = false;
    }
    psFree(outData);

    //
    // Call pmCullPeaks() with HUGE maxValue and disjoint psRegion.  Should
    // not remove any peaks.
    //
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmCullPeaks with large maxValue and disjoint psRegion.\n");
    outData = pmFindImagePeaks(imgData, 0.0);
    printf("pmFindImagePeaks found %d peaks\n", outData->size);
    psRegion *tmpRegion = psRegionAlloc(10000.0, 20000.0, 10000.0, 20000.0);
    outData = pmCullPeaks(outData, PS_MAX_F32, tmpRegion);

    if (outData == NULL) {
        printf("TEST ERROR: pmCullPeaks() returned a non-NULL psList.\n");
        testStatus = false;
        return(testStatus);
    }
    if (outData->size != numPeaksOrig) {
        printf("TEST ERROR (2): pmCullPeaks incorrectly removed peaks\n");
        printf("The pmCullPeaks() output had %d peaks, should have had %d peaks.\n", outData->size, numPeaksOrig);
        testStatus = false;
    }
    psFree(outData);
    psFree(tmpRegion);

    //
    // Call pmCullPeaks() with HUGE maxValue and non-disjoint psRegion.  Should
    // remove all peaks.
    //
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmCullPeaks with large maxValue and non-disjoint psRegion.\n");
    outData = pmFindImagePeaks(imgData, 0.0);
    printf("pmFindImagePeaks found %d peaks\n", outData->size);
    tmpRegion = psRegionAlloc(-PS_MAX_F32, PS_MAX_F32, -PS_MAX_F32, PS_MAX_F32);
    outData = pmCullPeaks(outData, PS_MAX_F32, tmpRegion);

    if (outData == NULL) {
        printf("TEST ERROR: pmCullPeaks() returned a non-NULL psList.\n");
        testStatus = false;
        return(testStatus);
    }
    if (outData->size != 0) {
        printf("TEST ERROR (3): pmCullPeaks incorrectly removed peaks\n");
        printf("The pmCullPeaks() output had %d peaks, should have had %d peaks.\n", outData->size, 0);
        testStatus = false;
    }
    psFree(outData);
    psFree(tmpRegion);

    printf("----------------------------------------------------------------------------------\n");
    psFree(imgData);
    return(testStatus);
}

#define TST04_NUM_ROWS 100
#define TST04_NUM_COLS 100
#define TST04_SKY 20.0
#define TST04_INNER_RADIUS 3
#define TST04_OUTER_RADIUS 5
/******************************************************************************
test04(): We first test pmSourceLocalSky() with various NULL and unallowable
input parameters.
 
XXX: Should we produce tests with boundary numbers for the inner/outer radius?
 
XXX: Call this with varying sizes for the image.
 *****************************************************************************/
int test04( void )
{
    bool testStatus = true;
    psImage *imgData = psImageAlloc(TST04_NUM_COLS, TST04_NUM_ROWS, PS_TYPE_F32);
    for (psS32 i = 0 ; i < imgData->numRows; i++) {
        for (psS32 j = 0 ; j < imgData->numCols; j++) {
            imgData->data.F32[i][j] = TST04_SKY;
        }
    }
    psImage *imgDataF64 = psImageAlloc(TST04_NUM_COLS, TST04_NUM_ROWS, PS_TYPE_F64);
    for (psS32 i = 0 ; i < imgDataF64->numRows; i++) {
        for (psS32 j = 0 ; j < imgDataF64->numCols; j++) {
            imgDataF64->data.F64[i][j] = 0.0;
        }
    }
    psSource *rc = NULL;
    psPeak *tmpPeak = pmPeakAlloc((psF32) (TST04_NUM_ROWS / 2),
                                  (psF32) (TST04_NUM_COLS / 2),
                                  200.0,
                                  PM_PEAK_LONE);

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with NULL psImage.  Should generate error and return NULL.\n");
    rc = pmSourceLocalSky(NULL, tmpPeak, PS_STAT_SAMPLE_MEAN, 10.0, 20.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with wrong-type psImage.  Should generate error and return NULL.\n");
    rc = pmSourceLocalSky(imgDataF64, tmpPeak, PS_STAT_SAMPLE_MEAN, 10.0, 20.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with NULL psPeak.  Should generate error and return NULL.\n");
    rc = pmSourceLocalSky(imgData, NULL, PS_STAT_SAMPLE_MEAN, 10.0, 20.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with innerRadius<0.0.  Should generate error and return NULL.\n");
    rc = pmSourceLocalSky(imgData, tmpPeak, PS_STAT_SAMPLE_MEAN, -10.0, 20.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with innerRadius>outerRadius.  Should generate error and return NULL.\n");
    rc = pmSourceLocalSky(imgData, tmpPeak, PS_STAT_SAMPLE_MEAN, 10.0, 5.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with subImage startRow < 0.  Should generate error and return NULL.\n");
    tmpPeak->x = (psF32) (TST04_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST04_NUM_COLS / 2);
    tmpPeak->y = 1;
    rc = pmSourceLocalSky(imgData, tmpPeak, PS_STAT_SAMPLE_MEAN, 10.0, 20.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with subImage endRow > numRows.  Should generate error and return NULL.\n");
    tmpPeak->x = (psF32) (TST04_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST04_NUM_COLS / 2);
    tmpPeak->y = TST04_NUM_ROWS;
    rc = pmSourceLocalSky(imgData, tmpPeak, PS_STAT_SAMPLE_MEAN, 10.0, 20.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with subImage startCol < 0.  Should generate error and return NULL.\n");
    tmpPeak->x = (psF32) (TST04_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST04_NUM_COLS / 2);
    tmpPeak->x = 1;
    rc = pmSourceLocalSky(imgData, tmpPeak, PS_STAT_SAMPLE_MEAN, 10.0, 20.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with subImage endCol > numCols.  Should generate error and return NULL.\n");
    tmpPeak->x = (psF32) (TST04_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST04_NUM_COLS / 2);
    tmpPeak->x = TST04_NUM_COLS;
    rc = pmSourceLocalSky(imgData, tmpPeak, PS_STAT_SAMPLE_MEAN, (psF32) TST04_INNER_RADIUS, (psF32) TST04_OUTER_RADIUS);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a non-NULL psSource.\n");
        psFree(rc);
        testStatus = false;
    }


    //
    // XXX: The following code should be a separate function, and we should call it
    // with a variety of image sizes, peaks.
    //
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceLocalSky with valid data.\n");
    tmpPeak->x = (psF32) (TST04_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST04_NUM_COLS / 2);
    rc = pmSourceLocalSky(imgData,
                          tmpPeak,
                          PS_STAT_SAMPLE_MEAN,
                          (psF32) TST04_INNER_RADIUS,
                          (psF32) TST04_OUTER_RADIUS);

    if (rc == NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource.\n");
        testStatus = false;
    } else {
        if (rc->peak == NULL) {
            printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource->peak.\n");
            testStatus = false;
        } else {
            if (rc->peak->x != tmpPeak->x) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->peak->x was %d, should have been %d.\n", rc->peak->x, tmpPeak->x);
                testStatus = false;
            }

            if (rc->peak->y != tmpPeak->y) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->peak->y was %d, should have been %d.\n", rc->peak->y, tmpPeak->y);
                testStatus = false;
            }

            if (rc->peak->counts != tmpPeak->counts) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->peak->counts was %f, should have been %f.\n", rc->peak->counts, tmpPeak->counts);
                testStatus = false;
            }

            if (rc->peak->class != tmpPeak->class) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->peak->class was %d, should have been %d.\n", rc->peak->class, tmpPeak->class);
                testStatus = false;
            }
        }

        if (rc->pixels == NULL) {
            printf("TEST ERROR: pmSourceLocalSky() psSource->pixels was NULL.\n");
            testStatus = false;
        } else {
            if (rc->pixels->numRows != (2 * TST04_OUTER_RADIUS)) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->pixels->numRows was %d, should have been %d.\n",
                       rc->pixels->numRows, (2 * TST04_OUTER_RADIUS));
                testStatus = false;
            }

            if (rc->pixels->numCols != (2 * TST04_OUTER_RADIUS)) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->pixels->numCols was %d, should have been %d.\n",
                       rc->pixels->numCols, (2 * TST04_OUTER_RADIUS));
                testStatus = false;
            }

            if (rc->pixels->col0 != (tmpPeak->x - TST04_OUTER_RADIUS)) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->pixels->col0 was %d, should have been %d.\n",
                       rc->pixels->col0, (tmpPeak->x - TST04_OUTER_RADIUS));
                testStatus = false;
            }

            if (rc->pixels->row0 != (tmpPeak->y - TST04_OUTER_RADIUS)) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->pixels->row0 was %d, should have been %d.\n",
                       rc->pixels->row0, (tmpPeak->y - TST04_OUTER_RADIUS));
                testStatus = false;
            }

            if (rc->pixels->type.type != PS_TYPE_F32) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->pixels->type was %d, should have been %d.\n",
                       rc->pixels->type.type, PS_TYPE_F32);
                testStatus = false;
            }

            // XXX: Test the rc->pixels-> row/col offsets.
            // XXX: Test that the pixels corresponds to the source image pixels.
        }

        if (rc->mask == NULL) {
            printf("TEST ERROR: pmSourceLocalSky() psSource->mask was NULL.\n");
            testStatus = false;
        } else {
            if (rc->mask->numRows != (2 * TST04_OUTER_RADIUS)) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->mask->numRows was %d, should have been %d.\n",
                       rc->mask->numRows, (2 * TST04_OUTER_RADIUS));
                testStatus = false;
            }

            if (rc->mask->numCols != (2 * TST04_OUTER_RADIUS)) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->mask->numCols was %d, should have been %d.\n",
                       rc->mask->numCols, (2 * TST04_OUTER_RADIUS));
                testStatus = false;
            }

            if (rc->mask->type.type != PS_TYPE_U8) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->mask->type was %d, should have been %d.\n",
                       rc->mask->type.type, PS_TYPE_U8);
                testStatus = false;
            }

            // XXX: Test the rc->mask-> row/col offsets.
            // XXX: Test that the correct pixels were masked, not merely the number of masked pixels.
            psS32 unmaskedPixels = 0;
            psS32 maskedPixels = 0;

            for (psS32 row = 0 ; row < rc->mask->numRows; row++) {
                for (psS32 col = 0 ; col < rc->mask->numCols; col++) {
                    if (rc->mask->data.U8[row][col] == 0) {
                        unmaskedPixels++;
                    } else {
                        maskedPixels++;
                    }
                }
            }
            if (maskedPixels != PS_SQR(2*(TST04_INNER_RADIUS-1))) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->mask had %d masked pixels, should have been %d.\n",
                       maskedPixels, PS_SQR(2*(TST04_INNER_RADIUS-1)));
                testStatus = false;
            }
            if (unmaskedPixels != (PS_SQR(2*TST04_OUTER_RADIUS) - PS_SQR(2*(TST04_INNER_RADIUS-1)))) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->mask had %d masked pixels, should have been %d.\n",
                       unmaskedPixels, (PS_SQR(2*TST04_OUTER_RADIUS) - PS_SQR(2*(TST04_INNER_RADIUS-1))));
                testStatus = false;
            }
        }

        if (rc->moments == NULL) {
            printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource->moments.\n");
            testStatus = false;
        } else {
            if (rc->moments->Sky != TST04_SKY) {
                printf("TEST ERROR: pmSourceLocalSky() psSource->moments->Sky was %f, should have been %f.\n", rc->moments->Sky, TST04_SKY);
                testStatus = false;
            }
        }
    }

    printf("----------------------------------------------------------------------------------\n");
    psFree(rc);
    psFree(imgData);
    psFree(imgDataF64);
    return(testStatus);
}

#define TST06_NUM_ROWS 100
#define TST06_NUM_COLS 100
#define TST06_SKY 20.0
#define TST06_INNER_RADIUS 3
#define TST06_OUTER_RADIUS 5
/******************************************************************************
test06(): We first test pmSourceLocalSky() with various NULL and unallowable
input parameters.
 
XXX: Should we produce tests with boundary numbers for the inner/outer radius?
 
XXX: Call this with varying sizes for the image.
 *****************************************************************************/
int test06( void )
{
    bool testStatus = true;
    psSource *tmpSource = NULL;
    bool rc = false;
    // Create the image used in this test.
    psImage *imgData = psImageAlloc(TST06_NUM_COLS, TST06_NUM_ROWS, PS_TYPE_F32);
    for (psS32 i = 0 ; i < imgData->numRows; i++) {
        for (psS32 j = 0 ; j < imgData->numCols; j++) {
            imgData->data.F32[i][j] = TST06_SKY;
        }
    }
    psImage *imgDataF64 = psImageAlloc(TST06_NUM_COLS, TST06_NUM_ROWS, PS_TYPE_F64);
    for (psS32 i = 0 ; i < imgDataF64->numRows; i++) {
        for (psS32 j = 0 ; j < imgDataF64->numCols; j++) {
            imgDataF64->data.F64[i][j] = 0.0;
        }
    }

    //
    // Create a psPeak with the center pixel set to the peak.
    //
    psPeak *tmpPeak = pmPeakAlloc((psF32) (TST06_NUM_ROWS / 2),
                                  (psF32) (TST06_NUM_COLS / 2),
                                  200.0,
                                  PM_PEAK_LONE);

    tmpSource = pmSourceLocalSky(imgData,
                                 tmpPeak,
                                 PS_STAT_SAMPLE_MEAN,
                                 (psF32) TST06_INNER_RADIUS,
                                 (psF32) TST06_OUTER_RADIUS);
    if (tmpSource == NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource.\n");
        psFree(imgData);
        psFree(imgDataF64);
        psFree(tmpSource);
        testStatus = false;
        return(testStatus);
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with NULL psSource.  Should generate error and return NULL.\n");
    rc = pmSourceSetPixelCircle(NULL, imgData, 10.0);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned TRUE.\n");
        testStatus = false;
    }
    // XXX: test with psSource->peaks NULL

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with NULL psImage.  Should generate error and return NULL.\n");
    rc = pmSourceSetPixelCircle(tmpSource, NULL, 10.0);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with wrong type psImage.  Should generate error and return NULL.\n");
    rc = pmSourceSetPixelCircle(tmpSource, imgDataF64, 10.0);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with radius < 0.0.  Should generate error and return NULL.\n");
    rc = pmSourceSetPixelCircle(tmpSource, imgData, -10.0);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with subImage startCol < 0.  Should generate error and return NULL.\n");
    tmpSource->peak->x = (psF32) (TST06_NUM_ROWS / 2);
    tmpSource->peak->y = (psF32) (TST06_NUM_COLS / 2);
    tmpSource->peak->x = 1;
    rc = pmSourceSetPixelCircle(tmpSource, imgData, 10.0);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with subImage endCol > numCols.  Should generate error and return NULL.\n");
    tmpSource->peak->x = (psF32) (TST06_NUM_ROWS / 2);
    tmpSource->peak->y = (psF32) (TST06_NUM_COLS / 2);
    tmpSource->peak->x = TST06_NUM_COLS;
    rc = pmSourceSetPixelCircle(tmpSource, imgData, 10.0);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with subImage startRow < 0.  Should generate error and return NULL.\n");
    tmpSource->peak->x = (psF32) (TST06_NUM_ROWS / 2);
    tmpSource->peak->y = (psF32) (TST06_NUM_COLS / 2);
    tmpSource->peak->y = 1;
    rc = pmSourceSetPixelCircle(tmpSource, imgData, 10.0);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with subImage endRow > numRows.  Should generate error and return NULL.\n");
    tmpSource->peak->x = (psF32) (TST06_NUM_ROWS / 2);
    tmpSource->peak->y = (psF32) (TST06_NUM_COLS / 2);
    tmpSource->peak->y = TST06_NUM_ROWS;
    rc = pmSourceSetPixelCircle(tmpSource, imgData, 10.0);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSetPixelCircle with valid data.\n");
    tmpSource->peak->x = (psF32) (TST06_NUM_ROWS / 2);
    tmpSource->peak->y = (psF32) (TST06_NUM_COLS / 2);
    rc = pmSourceSetPixelCircle(tmpSource, imgData, 10.0);

    if (rc == false) {
        printf("TEST ERROR: pmSourceSetPixelCircle() returned FALSE.\n");
        testStatus = false;
    } else {
        // XXX: Test correctness of the various psSetPixelCircle members.
    }

    printf("----------------------------------------------------------------------------------\n");
    psFree(imgData);
    psFree(imgDataF64);
    psFree(tmpSource);
    return(testStatus);

}

#define TST05_NUM_ROWS 100
#define TST05_NUM_COLS 100
#define TST05_SKY 20.0
#define TST05_INNER_RADIUS 3
#define TST05_OUTER_RADIUS 5
/******************************************************************************
test05(): We first test pmSourceMoments() with various NULL and unallowable
input parameters.
 
XXX: Should we produce tests with boundary numbers for the inner/outer radius?
 
XXX: Call this with varying sizes for the image.
 
XXX: The actual values of the moments are not tested.
 *****************************************************************************/
int test05( void )
{
    bool testStatus = true;
    psSource *tmpSource = NULL;
    psSource *rc = NULL;
    // Create the image used in this test.
    psImage *imgData = psImageAlloc(TST05_NUM_COLS, TST05_NUM_ROWS, PS_TYPE_F32);
    for (psS32 i = 0 ; i < imgData->numRows; i++) {
        for (psS32 j = 0 ; j < imgData->numCols; j++) {
            imgData->data.F32[i][j] = TST05_SKY;
        }
    }

    //
    // Create a psPeak with the center pixel set to the peak.
    //
    psPeak *tmpPeak = pmPeakAlloc((psF32) (TST05_NUM_ROWS / 2),
                                  (psF32) (TST05_NUM_COLS / 2),
                                  200.0,
                                  PM_PEAK_LONE);

    tmpSource = pmSourceLocalSky(imgData,
                                 tmpPeak,
                                 PS_STAT_SAMPLE_MEAN,
                                 (psF32) TST05_INNER_RADIUS,
                                 (psF32) TST05_OUTER_RADIUS);
    if (tmpSource == NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource.\n");
        psFree(imgData);
        psFree(tmpSource);
        testStatus = false;
        return(testStatus);
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceMoments with NULL psSource.  Should generate error and return NULL.\n");
    rc = pmSourceMoments(NULL, 10.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceMoments() returned a non-NULL psSource.\n");
        testStatus = false;
    }
    // XXX: test with psSource->peaks NULL
    // XXX: test with psSource->pixels NULL
    // XXX: test with psSource->mask NULL

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceMoments with radius < 0.0.  Should generate error and return NULL.\n");
    rc = pmSourceMoments(tmpSource, -10.0);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceMoments() returned a non-NULL psSource.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    psFree(imgData);
    psFree(tmpSource);
    return(testStatus);

}

/******************************************************************************
test07(): We first test the various psMinLM_... routines with various NULL and
unallowable input parameters.
 
XXX: We don't verify the numbers.  Must do this.
 *****************************************************************************/
int test07( void )
{
    bool testStatus = true;
    psF32 rc;
    psVector *deriv = psVectorAlloc(7, PS_TYPE_F32);
    psVector *params = psVectorAlloc(7, PS_TYPE_F32);
    psVector *x = psVectorAlloc(2, PS_TYPE_F32);
    for (psS32 i = 0 ; i < 7 ; i++) {
        deriv->data.F32[i] = 1.0;
        params->data.F32[i] = 1.0;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Gauss2D with NULL deriv vector.  Should not generate error.\n");
    rc = pmMinLM_Gauss2D(NULL, params, x);
    if (isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Gauss2D() returned a NAN psF32.\n");
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Gauss2D with NULL params vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_Gauss2D(deriv, NULL, x);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Gauss2D() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Gauss2D with NULL x vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_Gauss2D(deriv, params, NULL);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Gauss2D() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    psFree(deriv);
    psFree(params);

    deriv = psVectorAlloc(7, PS_TYPE_F32);
    params = psVectorAlloc(7, PS_TYPE_F32);
    for (psS32 i = 0 ; i < 7 ; i++) {
        deriv->data.F32[i] = 1.0;
        params->data.F32[i] = 1.0;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_PsuedoGauss2D with NULL deriv vector.  Should not generate error.\n");
    rc = pmMinLM_PsuedoGauss2D(NULL, params, x);
    if (isnan(rc)) {
        printf("TEST ERROR: pmMinLM_PsuedoGauss2D() returned a NAN psF32.\n");
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_PsuedoGauss2D with NULL params vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_PsuedoGauss2D(deriv, NULL, x);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_PsuedoGauss2D() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_PsuedoGauss2D with NULL x vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_PsuedoGauss2D(deriv, params, NULL);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_PsuedoGauss2D() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    psFree(deriv);
    psFree(params);


    deriv = psVectorAlloc(9, PS_TYPE_F32);
    params = psVectorAlloc(9, PS_TYPE_F32);
    for (psS32 i = 0 ; i < 9 ; i++) {
        deriv->data.F32[i] = 1.0;
        params->data.F32[i] = 1.0;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Wauss2D with NULL deriv vector.  Should not generate error.\n");
    rc = pmMinLM_Wauss2D(NULL, params, x);
    if (isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Wauss2D() returned a NAN psF32.\n");
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Wauss2D with NULL params vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_Wauss2D(deriv, NULL, x);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Wauss2D() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Wauss2D with NULL x vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_Wauss2D(deriv, params, NULL);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Wauss2D() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    psFree(deriv);
    psFree(params);

    deriv = psVectorAlloc(11, PS_TYPE_F32);
    params = psVectorAlloc(11, PS_TYPE_F32);
    for (psS32 i = 0 ; i < 11 ; i++) {
        deriv->data.F32[i] = 1.0;
        params->data.F32[i] = 1.0;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_TwistGauss2D with NULL deriv vector.  Should not generate error.\n");
    rc = pmMinLM_TwistGauss2D(NULL, params, x);
    if (isnan(rc)) {
        printf("TEST ERROR: pmMinLM_TwistGauss2D() returned a NAN psF32.\n");
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_TwistGauss2D with NULL params vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_TwistGauss2D(deriv, NULL, x);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_TwistGauss2D() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_TwistGauss2D with NULL x vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_TwistGauss2D(deriv, params, NULL);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_TwistGauss2D() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    psFree(deriv);
    psFree(params);

    deriv = psVectorAlloc(8, PS_TYPE_F32);
    params = psVectorAlloc(8, PS_TYPE_F32);
    for (psS32 i = 0 ; i < 8 ; i++) {
        deriv->data.F32[i] = 1.0;
        params->data.F32[i] = 1.0;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Sersic with NULL deriv vector.  Should not generate error.\n");
    rc = pmMinLM_Sersic(NULL, params, x);
    if (isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Sersic() returned a NAN psF32.\n");
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Sersic with NULL params vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_Sersic(deriv, NULL, x);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Sersic() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_Sersic with NULL x vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_Sersic(deriv, params, NULL);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_Sersic() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    psFree(deriv);
    psFree(params);

    deriv = psVectorAlloc(12, PS_TYPE_F32);
    params = psVectorAlloc(12, PS_TYPE_F32);
    for (psS32 i = 0 ; i < 12 ; i++) {
        deriv->data.F32[i] = 1.0;
        params->data.F32[i] = 1.0;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_SersicCore with NULL deriv vector.  Should not generate error.\n");
    rc = pmMinLM_SersicCore(NULL, params, x);
    if (isnan(rc)) {
        printf("TEST ERROR: pmMinLM_SersicCore() returned a NAN psF32.\n");
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_SersicCore with NULL params vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_SersicCore(deriv, NULL, x);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_SersicCore() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }
    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmMinLM_SersicCore with NULL x vector.  Should generate error and return NAN.\n");
    rc = pmMinLM_SersicCore(deriv, params, NULL);
    if (!isnan(rc)) {
        printf("TEST ERROR: pmMinLM_SersicCore() returned a non-NAN psF32 (%f).\n", rc);
        testStatus = false;
    }

    psFree(deriv);
    psFree(params);
    psFree(x);
    return(testStatus);
}

/******************************************************************************
test08(): We first test pmSourceModelGuess() with various NULL and unallowable
input parameters.
 
XXX: We don't verify the numbers.
 *****************************************************************************/
int test08( void )
{
    bool testStatus = true;
    psImage *imgData = psImageAlloc(TST04_NUM_COLS, TST04_NUM_ROWS, PS_TYPE_F32);
    for (psS32 i = 0 ; i < imgData->numRows; i++) {
        for (psS32 j = 0 ; j < imgData->numCols; j++) {
            imgData->data.F32[i][j] = TST04_SKY;
        }
    }
    psSource *mySrc = NULL;
    bool rc = false;
    psPeak *tmpPeak = pmPeakAlloc((psF32) (TST04_NUM_ROWS / 2),
                                  (psF32) (TST04_NUM_COLS / 2),
                                  200.0,
                                  PM_PEAK_LONE);

    printf("Calling pmSourceLocalSky with valid data.\n");
    tmpPeak->x = (psF32) (TST04_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST04_NUM_COLS / 2);
    mySrc = pmSourceLocalSky(imgData,
                             tmpPeak,
                             PS_STAT_SAMPLE_MEAN,
                             (psF32) TST04_INNER_RADIUS,
                             (psF32) TST04_OUTER_RADIUS);

    if (mySrc == NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with NULL psSource.  Should generate error, return FALSE.\n");
    rc = pmSourceModelGuess(NULL, imgData, PS_MODEL_GAUSS);
    if (rc == true) {
        printf("TEST ERROR: pmSourceModelGuess() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with NULL psImage.  Should generate error, return FALSE.\n");
    rc = pmSourceModelGuess(mySrc, NULL, PS_MODEL_GAUSS);
    if (rc == true) {
        printf("TEST ERROR: pmSourceModelGuess() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with bad model type.  Should generate error, return FALSE.\n");
    rc = pmSourceModelGuess(mySrc, imgData, PS_MODEL_UNDEFINED);
    if (rc == true) {
        printf("TEST ERROR: pmSourceModelGuess() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with PS_MODEL_GAUSS\n");
    rc = pmSourceModelGuess(mySrc, imgData, PS_MODEL_GAUSS);
    if (rc != true) {
        printf("TEST ERROR: pmSourceModelGuess() returned FALSE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with PS_MODEL_PGAUSS\n");
    rc = pmSourceModelGuess(mySrc, imgData, PS_MODEL_PGAUSS);
    if (rc != true) {
        printf("TEST ERROR: pmSourceModelGuess() returned FALSE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with PS_MODEL_TWIST_GAUSS\n");
    rc = pmSourceModelGuess(mySrc, imgData, PS_MODEL_TWIST_GAUSS);
    if (rc != true) {
        printf("TEST ERROR: pmSourceModelGuess() returned FALSE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with PS_MODEL_WAUSS\n");
    rc = pmSourceModelGuess(mySrc, imgData, PS_MODEL_WAUSS);
    if (rc != true) {
        printf("TEST ERROR: pmSourceModelGuess() returned FALSE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with PS_MODEL_SERSIC\n");
    rc = pmSourceModelGuess(mySrc, imgData, PS_MODEL_SERSIC);
    if (rc != true) {
        printf("TEST ERROR: pmSourceModelGuess() returned FALSE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceModelGuess with PS_MODEL_SERSIC_CORE\n");
    rc = pmSourceModelGuess(mySrc, imgData, PS_MODEL_SERSIC_CORE);
    if (rc != true) {
        printf("TEST ERROR: pmSourceModelGuess() returned FALSE.\n");
        testStatus = false;
    }

    psFree(mySrc);
    //    psFree(tmpPeak);
    psFree(imgData);
    return(testStatus);
}


#define TST09_NUM_ROWS 70
#define TST09_NUM_COLS 70
#define TST09_SKY 5.0
#define TST09_INNER_RADIUS 3
#define TST09_OUTER_RADIUS 10
#define LEVEL (TST09_SKY + 10.0)
/******************************************************************************
test09(): We first test pmSourceContour() with various NULL and unallowable
input parameters.
 
XXX: We don't verify the numbers.
 *****************************************************************************/
int test09( void )
{
    bool testStatus = true;
    psImage *imgData = psImageAlloc(TST09_NUM_COLS, TST09_NUM_ROWS, PS_TYPE_F32);
    for (psS32 i = 0 ; i < imgData->numRows; i++) {
        for (psS32 j = 0 ; j < imgData->numCols; j++) {
            imgData->data.F32[i][j] = TST09_SKY;
        }
    }
    psSource *mySrc = NULL;
    psArray *rc = NULL;

    psPeak *tmpPeak = pmPeakAlloc((psF32) (TST09_NUM_ROWS / 2),
                                  (psF32) (TST09_NUM_COLS / 2),
                                  200.0,
                                  PM_PEAK_LONE);

    printf("Calling pmSourceLocalSky with valid data.\n");
    tmpPeak->x = (psF32) (TST09_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST09_NUM_COLS / 2);
    mySrc = pmSourceLocalSky(imgData,
                             tmpPeak,
                             PS_STAT_SAMPLE_MEAN,
                             (psF32) TST09_INNER_RADIUS,
                             (psF32) TST09_OUTER_RADIUS);

    if (mySrc == NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource.\n");
        testStatus = false;
    }

    bool rcBool = pmSourceModelGuess(mySrc, imgData, PS_MODEL_GAUSS);
    if (rcBool != true) {
        printf("TEST ERROR: pmSourceModelGuess() returned FALSE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceContour with NULL psSource .  Should generate error, return NULL.\n");
    rc = pmSourceContour(NULL, imgData, LEVEL, PS_CONTOUR_CRUDE);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceContour() returned non-NULL.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceContour with NULL psImage .  Should generate error, return NULL.\n");
    rc = pmSourceContour(mySrc, NULL, LEVEL, PS_CONTOUR_CRUDE);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceContour() returned non-NULL.\n");
        testStatus = false;
    }


    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceContour with acceptable data.\n");
    printf("NOTE: must figure out the parameters for this test to be meaningful.\n");
    mySrc->models->params[0] = TST09_SKY;
    mySrc->models->params[1] = 15.0;
    mySrc->models->params[2] = (psF32) (TST09_NUM_ROWS / 2);
    mySrc->models->params[3] = (psF32) (TST09_NUM_COLS / 2);
    mySrc->models->params[4] = 2.0;
    mySrc->models->params[5] = 2.0;
    mySrc->models->params[6] = 2.0;
    rc = pmSourceContour(mySrc, imgData, LEVEL, PS_CONTOUR_CRUDE);
    if (rc != NULL) {
        printf("TEST ERROR: pmSourceContour() returned NULL.\n");
        testStatus = false;
    }

    psFree(rc);
    psFree(mySrc);
    //    psFree(tmpPeak);
    psFree(imgData);
    return(testStatus);
}

#define TST15_NUM_ROWS 100
#define TST15_NUM_COLS 100
#define TST15_SKY 10.0
#define TST15_INNER_RADIUS 3
#define TST15_OUTER_RADIUS 5
/******************************************************************************
test15(): We first test pmSourceAddModel() with various NULL and unallowable
input parameters.
 
XXX: We don't verify the numbers.
 *****************************************************************************/
int test15( void )
{
    bool testStatus = true;
    psImage *imgData = psImageAlloc(TST15_NUM_COLS, TST15_NUM_ROWS, PS_TYPE_F32);
    for (psS32 i = 0 ; i < imgData->numRows; i++) {
        for (psS32 j = 0 ; j < imgData->numCols; j++) {
            imgData->data.F32[i][j] = TST15_SKY;
        }
    }
    psSource *mySrc = NULL;
    psBool rc = false;

    psPeak *tmpPeak = pmPeakAlloc((psF32) (TST15_NUM_ROWS / 2),
                                  (psF32) (TST15_NUM_COLS / 2),
                                  200.0,
                                  PM_PEAK_LONE);

    printf("Calling pmSourceLocalSky with valid data.\n");
    tmpPeak->x = (psF32) (TST15_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST15_NUM_COLS / 2);
    mySrc = pmSourceLocalSky(imgData,
                             tmpPeak,
                             PS_STAT_SAMPLE_MEAN,
                             (psF32) TST15_INNER_RADIUS,
                             (psF32) TST15_OUTER_RADIUS);

    if (mySrc == NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource.\n");
        testStatus = false;
    }

    mySrc->models = pmModelAlloc(PS_MODEL_GAUSS);
    mySrc->models->params[0] = 5.0;
    mySrc->models->params[1] = 70.0;
    mySrc->models->params[2] = (psF32) (TST15_NUM_ROWS / 2);
    mySrc->models->params[3] = (psF32) (TST15_NUM_COLS / 2);
    mySrc->models->params[4] = 1.0;
    mySrc->models->params[5] = 1.0;
    mySrc->models->params[6] = 2.0;

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceAddModel with NULL psImage.  Should generate error, return FALSE.\n");
    rc = pmSourceAddModel(NULL, mySrc, true);
    if (rc == true) {
        printf("TEST ERROR: pmSourceAddModel() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceAddModel with NULL psSrc.  Should generate error, return FALSE.\n");
    rc = pmSourceAddModel(imgData, NULL, true);
    if (rc == true) {
        printf("TEST ERROR: pmSourceAddModel() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceAddModel with acceptable data.\n");
    rc = pmSourceAddModel(imgData, mySrc, true);
    if (rc != true) {
        printf("TEST ERROR: pmSourceAddModel() returned FALSE.\n");
        testStatus = false;
    }

    psFree(mySrc);
    psFree(imgData);
    return(testStatus);
}

#define TST16_NUM_ROWS 100
#define TST16_NUM_COLS 100
#define TST16_SKY 10.0
#define TST16_INNER_RADIUS 3
#define TST16_OUTER_RADIUS 5
/******************************************************************************
test16(): We first test pmSourceSubModel() with various NULL and unallowable
input parameters.
 
XXX: We don't verify the numbers.
 *****************************************************************************/
int test16( void )
{
    bool testStatus = true;
    psImage *imgData = psImageAlloc(TST16_NUM_COLS, TST16_NUM_ROWS, PS_TYPE_F32);
    for (psS32 i = 0 ; i < imgData->numRows; i++) {
        for (psS32 j = 0 ; j < imgData->numCols; j++) {
            imgData->data.F32[i][j] = TST16_SKY;
        }
    }
    psSource *mySrc = NULL;
    psBool rc = false;

    psPeak *tmpPeak = pmPeakAlloc((psF32) (TST16_NUM_ROWS / 2),
                                  (psF32) (TST16_NUM_COLS / 2),
                                  200.0,
                                  PM_PEAK_LONE);

    printf("Calling pmSourceLocalSky with valid data.\n");
    tmpPeak->x = (psF32) (TST16_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST16_NUM_COLS / 2);
    mySrc = pmSourceLocalSky(imgData,
                             tmpPeak,
                             PS_STAT_SAMPLE_MEAN,
                             (psF32) TST16_INNER_RADIUS,
                             (psF32) TST16_OUTER_RADIUS);

    if (mySrc == NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource.\n");
        testStatus = false;
    }

    mySrc->models = pmModelAlloc(PS_MODEL_GAUSS);
    mySrc->models->params[0] = 5.0;
    mySrc->models->params[1] = 70.0;
    mySrc->models->params[2] = (psF32) (TST16_NUM_ROWS / 2);
    mySrc->models->params[3] = (psF32) (TST16_NUM_COLS / 2);
    mySrc->models->params[4] = 1.0;
    mySrc->models->params[5] = 1.0;
    mySrc->models->params[6] = 2.0;

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSubModel with NULL psImage.  Should generate error, return FALSE.\n");
    rc = pmSourceSubModel(NULL, mySrc, true);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSubModel() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSubModel with NULL psSrc.  Should generate error, return FALSE.\n");
    rc = pmSourceSubModel(imgData, NULL, true);
    if (rc == true) {
        printf("TEST ERROR: pmSourceSubModel() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceSubModel with acceptable data.\n");
    rc = pmSourceSubModel(imgData, mySrc, true);
    if (rc != true) {
        printf("TEST ERROR: pmSourceSubModel() returned FALSE.\n");
        testStatus = false;
    }

    psFree(mySrc);
    psFree(imgData);
    return(testStatus);
}

#define TST20_NUM_ROWS 100
#define TST20_NUM_COLS 100
#define TST20_SKY 10.0
#define TST20_INNER_RADIUS 3
#define TST20_OUTER_RADIUS 5
/******************************************************************************
test20(): We first test pmSourceSubModel() with various NULL and unallowable
input parameters.
 
XXX: We don't verify the numbers.
 *****************************************************************************/
int test20( void )
{
    bool testStatus = true;
    psImage *imgData = psImageAlloc(TST20_NUM_COLS, TST20_NUM_ROWS, PS_TYPE_F32);
    for (psS32 i = 0 ; i < imgData->numRows; i++) {
        for (psS32 j = 0 ; j < imgData->numCols; j++) {
            imgData->data.F32[i][j] = TST20_SKY;
        }
    }
    psSource *mySrc = NULL;
    psBool rc = false;

    psPeak *tmpPeak = pmPeakAlloc((psF32) (TST20_NUM_ROWS / 2),
                                  (psF32) (TST20_NUM_COLS / 2),
                                  200.0,
                                  PM_PEAK_LONE);

    printf("Calling pmSourceLocalSky with valid data.\n");
    tmpPeak->x = (psF32) (TST20_NUM_ROWS / 2);
    tmpPeak->y = (psF32) (TST20_NUM_COLS / 2);
    mySrc = pmSourceLocalSky(imgData,
                             tmpPeak,
                             PS_STAT_SAMPLE_MEAN,
                             (psF32) TST20_INNER_RADIUS,
                             (psF32) TST20_OUTER_RADIUS);

    if (mySrc == NULL) {
        printf("TEST ERROR: pmSourceLocalSky() returned a NULL psSource.\n");
        testStatus = false;
    }

    mySrc->models = pmModelAlloc(PS_MODEL_GAUSS);


    mySrc->models->params[0] = 5.0;
    mySrc->models->params[1] = 70.0;
    mySrc->models->params[2] = (psF32) (TST20_NUM_ROWS / 2);
    mySrc->models->params[3] = (psF32) (TST20_NUM_COLS / 2);
    mySrc->models->params[4] = 1.0;
    mySrc->models->params[5] = 1.0;
    mySrc->models->params[6] = 2.0;

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceFitModel with NULL psImage.  Should generate error, return FALSE.\n");
    rc = pmSourceFitModel(mySrc, NULL);
    if (rc == true) {
        printf("TEST ERROR: pmSourceFitModel() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceFitModel with NULL psSource.  Should generate error, return FALSE.\n");
    rc = pmSourceFitModel(NULL, imgData);
    if (rc == true) {
        printf("TEST ERROR: pmSourceFitModel() returned TRUE.\n");
        testStatus = false;
    }

    printf("----------------------------------------------------------------------------------\n");
    printf("Calling pmSourceFitModel with acceptable data.\n");
    rc = pmSourceFitModel(mySrc, imgData);
    printf("pmSourceFitModel returned %d\n", rc);

    // XXX: Memory leaks are not being tested
    psVector *junk = psVectorAlloc(10, PS_TYPE_F32);
    junk->data.F32[0] = 0.0;

    psFree(mySrc);
    psFree(imgData);
    return(testStatus);
}


// this code will



