Changeset 6872

trunk/psModules/configure.ac

-              r6298
+              r6872
 AM_MAINTAINER_MODE
+dnl handle debug building
+AC_ARG_ENABLE(optimize,
+  [AS_HELP_STRING(--enable-optimize,enable compiler optimization)],
+  [AC_MSG_RESULT(compile optimization enabled)
+   CFLAGS="${CFLAGS=} -O2"],
+  [AC_MSG_RESULT([compile optimization disabled])
+   CFLAGS="${CFLAGS=} -O0 -g"]
+)
+dnl handle profiler building
+AC_ARG_ENABLE(profile,
+  [AS_HELP_STRING(--enable-profile,enable compiler profiler information inclusion)],
+  [AC_MSG_RESULT(compile optimization enabled)
+   CFLAGS="${CFLAGS=} -fprofile-arcs -ftest-coverage -pg"]
+)
+dnl otherise AC_PROG_CC will default CFLAGS to "-g -02"
+if test -z ${CFLAGS} ; then
+  CFLAGS=""
+fi
 AC_SUBST([CFLAGS])
 …
 AC_PROG_CC
 dnl FIXME document why GNU_SOURCE is being set
 dnl AC_GNU_SOURCE
+AC_GNU_SOURCE
 AC_PROG_INSTALL
 AM_PROG_LIBTOOL
-AC_SYS_LARGEFILE
 AC_PREFIX_DEFAULT([`pwd`])
 …
 SRCPATH='${top_srcdir}/src'
 SRCDIRS="astrom camera config detrend imcombine imsubtract objects photom"
+SRCDIRS="astrom config detrend imcombine imsubtract objects pslib"
 # escape two escapes at this level so \\ gets passed to the shell and \ to perl
 SRCINC=`echo "${SRCDIRS=}" | ${PERL} -pe "s|(\w+)|-I\\\\${SRCPATH=}/\1|g"`
 SRCINC="-I${SRCPATH=} ${SRCINC=}"
 SRCSUBLIBS=`echo "${SRCDIRS=}" | ${PERL} -pe "s|(\w+)|\1/libpsmodule\1.la|g"`
 AC_SUBST([SRCSUBLIBS],${SRCSUBLIBS=})
 AC_SUBST([SRCINC],${SRCINC=})
+AC_SUBST(SRCSUBLIBS,${SRCSUBLIBS=})
+AC_SUBST(SRCINC,${SRCINC=})
 AC_SUBST([SRCDIRS],${SRCDIRS=})
+dnl handle debug building
+AC_ARG_ENABLE(optimize,
+  [AS_HELP_STRING(--enable-optimize,enable compiler optimization)],
+  [AC_MSG_RESULT(compile optimization enabled)
+   CFLAGS="${CFLAGS=} -O2"],
+  [AC_MSG_RESULT([compile optimization disabled])
+   CFLAGS="${CFLAGS=} -O0 -g"]
+)
+dnl handle profiler building
+AC_ARG_ENABLE(profile,
+  [AS_HELP_STRING(--enable-profile,enable compiler profiler information inclusion)],
+  [AC_MSG_RESULT(compile optimization enabled)
+   CFLAGS="${CFLAGS=} -pg"]
+)
+dnl CFLAGS="${CFLAGS=} -fprofile-arcs -ftest-coverage -pg"]
 dnl doxygen -------------------------------------------------------------------
 …
 AC_CHECK_PROG([doxygen], [doxygen], [true], [false])
 AM_CONDITIONAL(DOXYGEN, test x$doxygen = xtrue)
+dnl libjpeg -------------------------------------------------------------------
+AC_CHECK_LIB(jpeg,jpeg_CreateCompress,[],[AC_MSG_ERROR([jpeg library not found.])])
 dnl pslib ---------------------------------------------------------------------
 …
 if test -z ${PSLIB_CONFIG} ; then
   PKG_CHECK_MODULES([PSLIB], [pslib >= 0.10.0])
+  PKG_CHECK_MODULES([PSLIB], [pslib >= 0.9.0])
 else
   AC_CHECK_FILE($PSLIB_CONFIG,[],
 …
   src/Makefile
   src/astrom/Makefile
-  src/camera/Makefile
   src/config/Makefile
   src/detrend/Makefile
 …
   src/imsubtract/Makefile
   src/objects/Makefile
   src/photom/Makefile
+  src/pslib/Makefile
   test/Makefile
   test/astrom/Makefile
-  test/camera/Makefile
   test/config/Makefile
   test/detrend/Makefile
 …
   test/imsubtract/Makefile
   test/objects/Makefile
   test/photom/Makefile
+  test/pslib/Makefile
   Doxyfile
   psmodule-config

trunk/psModules/src/Makefile.am

r5685	r6872
7	7	libpsmodule_la_DEPENDENCIES = $(SRCSUBLIBS)
8	8
	9	include_HEADERS = \
	10	psmodules.h
	11
9	12	EXTRA_DIST = psErrorCodes.dat

trunk/psModules/src/astrom/Makefile.am

-              r6301
+              r6872
 noinst_LTLIBRARIES = libpsmoduleastrom.la
 libpsmoduleastrom_la_CPPFLAGS = $(SRCINC) $(PSMODULE_CFLAGS)
+libpsmoduleastrom_la_CPPFLAGS = $(SRCINC) $(PSMODULE_CFLAGS) -I../pslib/
 libpsmoduleastrom_la_LDFLAGS  = -release $(PACKAGE_VERSION)
 libpsmoduleastrom_la_SOURCES  = \
+    pmAstrometry.c \
+    pmAstrometryObjects.c
+        pmFPA.c \
+        pmFPAConstruct.c \
+        pmFPACopy.c \
+        pmFPAMaskWeight.c \
+        pmFPARead.c \
+        pmFPAUtils.c \
+        pmFPAWrite.c \
+        pmHDU.c \
+        pmHDUUtils.c \
+        pmReadout.c \
+        pmConcepts.c \
+        pmConceptsRead.c \
+        pmConceptsWrite.c \
+        pmConceptsStandard.c \
+        pmFPA_JPEG.c \
+        pmFPAview.c \
+        pmFPAfile.c
+#       pmFPAAstrometry.c
+#       pmAstrometryObjects.c
+#       pmChipMosaic.c
 psmoduleincludedir = $(includedir)
 psmoduleinclude_HEADERS = \
+    pmAstrometry.h \
+    pmAstrometryObjects.h
+        pmFPA.h \
+        pmFPAConstruct.h \
+        pmFPACopy.h \
+        pmFPAMaskWeight.h \
+        pmFPARead.h \
+        pmFPAUtils.h \
+        pmFPAWrite.h \
+        pmHDU.h \
+        pmHDUUtils.h \
+        pmReadout.h \
+        pmConcepts.h \
+        pmConceptsRead.h \
+        pmConceptsWrite.h \
+        pmConceptsStandard.h \
+        pmFPA_JPEG.h \
+        pmFPAview.h \
+        pmFPAfile.h
+#       pmFPAAstrometry.h
+#       pmAstrometryObjects.h
+#       pmChipMosaic.h

trunk/psModules/src/astrom/pmAstrometry.c

-              r6205
+              r6872
 * XXX: We should review the extent of the warning messages on these functions
 * when the transformations are not successful.
+*
+*
 * XXX: Should we implement non-linear cell->chip transforms?
+*
 *  @version $Revision: 1.12 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2006-01-26 21:10:50 $
+*
+*  @version $Revision: 1.13 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2006-04-17 18:01:04 $
+*
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
 …
 #include <math.h>
 #include "pslib.h"
 #include "pmAstrometry.h"
+#include "pmConcepts.h"
 /*****************************************************************************
 …
         psFree(readout->mask);
         psFree(readout->weight);
-        //
-        // XXX: Not sure if this is the right way to do things.  Currently the psListAdd()
-        // increase the memory reference counter to the list data.  So, we
-        // iterate through the list, and decrement the reference counters.
-        //
-        if (1) {
-            if ((readout->bias != NULL) && (readout->bias->head != NULL)) {
-                psListElem *tmpElem = (psListElem *) readout->bias->head;
-                while (NULL != tmpElem) {
-                    psMemDecrRefCounter((psImage *) tmpElem->data);
-                    tmpElem = tmpElem->next;
+                }
+            }
+        }
-        psFree(readout->bias);
         psFree(readout->analysis);
+        #if 0
         psFree(readout->parent);
+        #endif
+        readout->parent = NULL;
+    }
+}
 …
         psFree(cell->concepts);
         psFree(cell->analysis);
+        psFree(cell->camera);
         //
         // Set the parent to NULL in all cell->readouts before psFree(cell->readouts)
         // in order to avoid memory reference counter problems.
         //
+        #if 0
         for (psS32 i = 0 ; i < cell->readouts->n ; i++) {
             pmReadout *tmpReadout = (pmReadout *) cell->readouts->data[i];
 …
+            }
+        }
+        psFree(cell->parent);
+        #endif
+        cell->parent = NULL;
         psFree(cell->readouts);
-        psFree(cell->parent);
         psFree(cell->hdu);
+    }
+}
 …
         // in order to avoid memory reference counter problems.
         //
+        #if 0
         for (psS32 i = 0 ; i < chip->cells->n ; i++) {
             pmCell *tmpCell = (pmCell *) chip->cells->data[i];
 …
+            }
+        }
+        psFree(chip->parent);
+        #endif
+        chip->parent = NULL;
         psFree(chip->cells);
-        psFree(chip->parent);
         psFree(chip->hdu);
+    }
 …
         // in order to avoid memory reference counter problems.
         //
+        #if 0
         for (psS32 i = 0 ; i < fpa->chips->n ; i++) {
             pmChip *tmpChip = (pmChip *) fpa->chips->data[i];
 …
+            }
+        }
+        #endif
         psFree(fpa->chips);
         psFree(fpa->hdu);
 …
+}
+void p_pmHDUFree(p_pmHDU *hdu)
+{
+    if (hdu) {
+        psFree(hdu->extname);
+        psFree(hdu->header);
+        psFree(hdu->images);
+        psFree(hdu->masks);
+        psFree(hdu->weights);
+    }
+}
 // XXX: Verify these default values for row0, col0, rowBins, colBins
+// PAP: These values may disappear in the future in favour of values in parent->concepts?
 pmReadout *pmReadoutAlloc(pmCell *cell)
+{
     pmReadout *tmpReadout = (pmReadout *) psAlloc(sizeof(pmReadout));
     tmpReadout->col0 = -1;
     tmpReadout->row0 = -1;
     tmpReadout->colBins = -1;
     tmpReadout->rowBins = -1;
+    tmpReadout->col0 = 0;
+    tmpReadout->row0 = 0;
+    tmpReadout->colBins = 0;
+    tmpReadout->rowBins = 0;
     tmpReadout->image = NULL;
     tmpReadout->mask = NULL;
     tmpReadout->weight = NULL;
-    tmpReadout->bias = NULL;
     tmpReadout->analysis = psMetadataAlloc();
     tmpReadout->parent = cell;
 …
 // XXX: Verify these default values for row0, col0.
+// PAP: These values may disappear in the future in favour of values in the "concepts"?
 pmCell *pmCellAlloc(
     pmChip *chip,
     psMetadata *cameradata,
     psString name)
+    psMetadata *cameraData,
+    const char *name)
+{
     pmCell *tmpCell = (pmCell *) psAlloc(sizeof(pmCell));
     tmpCell->col0 = -1;
     tmpCell->row0 = -1;
+    tmpCell->col0 = 0;
+    tmpCell->row0 = 0;
     tmpCell->toChip = NULL;
     tmpCell->toFPA = NULL;
 …
         psLogMsg(__func__, PS_LOG_WARN, "WARNING: Could not add CELL.NAME to metadata.\n");
+    }
     tmpCell->camera = cameradata;
     tmpCell->analysis = NULL;
+    tmpCell->camera = psMemIncrRefCounter(cameraData);
+    tmpCell->analysis = psMetadataAlloc();
     tmpCell->readouts = psArrayAlloc(0);
     tmpCell->parent = chip;
 …
         chip->cells = psArrayAdd(chip->cells, 1, (psPtr) tmpCell);
+    }
+    tmpCell->valid = false;
+    tmpCell->process = true;            // All cells are processed by default
+    tmpCell->exists = false;            // Not yet read in
     tmpCell->hdu = NULL;
+    pmConceptsBlankCell(tmpCell);
     psMemSetDeallocator(tmpCell, (psFreeFunc) cellFree);
 …
 // XXX: Verify these default values for row0, col0.
+// PAP: row0, col0 may disappear in the future in favour of storing values in the "concepts".
 pmChip *pmChipAlloc(
     pmFPA *fpa,
     psString name)
+    const char *name)
+{
     pmChip *tmpChip = (pmChip *) psAlloc(sizeof(pmChip));
     tmpChip->col0 = -1;
     tmpChip->row0 = -1;
+    tmpChip->col0 = 0;
+    tmpChip->row0 = 0;
     tmpChip->toFPA = NULL;
     tmpChip->fromFPA = NULL;
 …
         psLogMsg(__func__, PS_LOG_WARN, "WARNING: Could not add CHIP.NAME to metadata.\n");
+    }
     tmpChip->analysis = NULL;
+    tmpChip->analysis = psMetadataAlloc();
     tmpChip->cells = psArrayAlloc(0);
     tmpChip->parent = fpa;
 …
         fpa->chips = psArrayAdd(fpa->chips, 1, (psPtr) tmpChip);
+    }
+    tmpChip->valid = false;
+    tmpChip->process = true;            // Work on all chips, by default
+    tmpChip->exists = false;            // Not read in yet
     tmpChip->hdu = NULL;
+    pmConceptsBlankChip(tmpChip);
     psMemSetDeallocator(tmpChip, (psFreeFunc) chipFree);
 …
     tmpFPA->concepts = psMetadataAlloc();
     tmpFPA->analysis = NULL;
     tmpFPA->camera = camera;
+    tmpFPA->camera = psMemIncrRefCounter((psPtr)camera);
     tmpFPA->chips = psArrayAlloc(0);
     tmpFPA->hdu = NULL;
     tmpFPA->phu = NULL;
+    pmConceptsBlankFPA(tmpFPA);
     psMemSetDeallocator(tmpFPA, (psFreeFunc) FPAFree);
     return(tmpFPA);
+}
+p_pmHDU *p_pmHDUAlloc(const char *extname)
+{
+    p_pmHDU *hdu = psAlloc(sizeof(p_pmHDU));
+    psMemSetDeallocator(hdu, (psFreeFunc)p_pmHDUFree);
+    hdu->extname = psStringCopy(extname);
+    hdu->header = NULL;
+    hdu->images = NULL;
+    hdu->masks = NULL;
+    hdu->weights = NULL;
+    return hdu;
+}
 …
-/*****************************************************************************/
-/* FUNCTION IMPLEMENTATION - PUBLIC                                          */
-/*****************************************************************************/
-pmCell* pmCellInFPA(
-    const psPlane* fpaCoord,
-    const pmFPA* FPA)
+{
-    PS_ASSERT_PTR_NON_NULL(fpaCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(FPA, NULL);
-    pmChip* tmpChip = NULL;
-    psPlane chipCoord;
-    pmCell* outCell = NULL;
-    // Determine which chip contains the fpaCoords.
-    tmpChip = pmChipInFPA(fpaCoord, FPA);
-    if (tmpChip == NULL) {
-        return(NULL);
+    }
-    // Convert to those chip coordinates.
-    psPlane *rc = pmCoordFPAToChip(&chipCoord, fpaCoord, tmpChip);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine Chip coords.\n");
-        return(NULL);
+    }
-    // Determine which cell contains those chip coordinates.
-    outCell = pmCellInChip(&chipCoord, tmpChip);
-    if (outCell == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the cell.\n");
-        return(NULL);
+    }
-    return (outCell);
+}
-pmChip* pmChipInFPA(
-    const psPlane* fpaCoord,
-    const pmFPA* FPA)
+{
-    PS_ASSERT_PTR_NON_NULL(fpaCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(FPA, NULL);
-    PS_ASSERT_PTR_NON_NULL(FPA->chips, NULL);
-    psArray* chips = FPA->chips;
-    psS32 nChips = chips->n;
-    psPlane chipCoord;
-    pmCell *tmpCell = NULL;
-    //
-    // Loop through every chip in this FPA.  Convert the original FPA
-    // coordinates to chip coordinates for that chip.  Then, determine if any
-    // cells in that chip contain those chip coordinates.
-    // XXX: Depending on the number of chips, and their topology, there may be
-    // a much more efficient way of doing this.
-    //
-    for (psS32 i = 0; i < nChips; i++) {
-        pmChip* tmpChip = chips->data[i];
-        PS_ASSERT_PTR_NON_NULL(tmpChip, NULL);
-        PS_ASSERT_PTR_NON_NULL(tmpChip->fromFPA, NULL);
-        psPlaneTransformApply(&chipCoord, tmpChip->fromFPA, fpaCoord);
-        tmpCell = pmCellInChip(&chipCoord, tmpChip);
-        if (tmpCell != NULL) {
-            return(tmpChip);
+        }
+    }
-    psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the chip.\n");
-    return (NULL);
+}
-pmCell* pmCellInChip(
-    const psPlane* chipCoord,
-    const pmChip* chip)
+{
-    PS_ASSERT_PTR_NON_NULL(chipCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(chip, NULL);
-    psPlane cellCoord;
-    psArray* cells;
-    cells = chip->cells;
-    if (cells == NULL) {
-        return NULL;
+    }
-    //
-    // We loop over each cell in the chip.  We transform the chipCoord into
-    // a cellCoord for that cell and determine if that cellCoord is valid.
-    // If so, then we return that cell.
-    // XXX: Depending on the number of cells, and their topology, there may be
-    // a much more efficient way of doing this.
-    //
-    for (psS32 i = 0; i < cells->n; i++) {
-        pmCell* tmpCell = (pmCell* ) cells->data[i];
-        PS_ASSERT_PTR_NON_NULL(tmpCell, NULL);
-        psPlaneTransform *chipToCell = NULL;
-        if (true ==  p_psIsProjectionLinear(tmpCell->toChip)) {
-            chipToCell = p_psPlaneTransformLinearInvert(tmpCell->toChip);
-        } else {
-            psLogMsg(__func__, PS_LOG_WARN, "WARNING: non-linear cell->chip transforms are not yet implemented.\n");
-            //chipToCell = psPlaneTransformInvert(NULL, tmpCell->toChip, NULL, -1);
-            chipToCell = NULL;
+        }
-        if (chipToCell == NULL) {
-            psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not invert the Cell->toChip transform.\n");
-            return(NULL);
+        }
-        psArray* readouts = tmpCell->readouts;
-        if (readouts != NULL) {
-            for (psS32 j = 0; j < readouts->n; j++) {
-                pmReadout* tmpReadout = readouts->data[j];
-                PS_ASSERT_READOUT_NON_NULL(tmpReadout, NULL);
-                psPlaneTransformApply(&cellCoord,
-                                      chipToCell,
-                                      chipCoord);
-                if (checkValidImageCoords(cellCoord.x,
-                                          cellCoord.y,
-                                          tmpReadout->image)) {
-                    psFree(chipToCell);
-                    return (tmpCell);
+                }
+            }
+        }
-        psFree(chipToCell);
+    }
-    //psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the cell.\n");
-    return (NULL);
+}
-psPlane* pmCoordCellToFPA(
-    psPlane* fpaCoord,
-    const psPlane* cellCoord,
-    const pmCell* cell)
+{
-    PS_ASSERT_PTR_NON_NULL(cellCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(cell, NULL);
-    psPlane *rc = psPlaneTransformApply(fpaCoord, cell->toFPA, cellCoord);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform cell coords to FPA coords.\n");
+    }
-    return(rc);
+}
-psPlane* pmCoordChipToFPA(
-    psPlane* outCoord,
-    const psPlane* inCoord,
-    const pmChip* chip)
+{
-    PS_ASSERT_PTR_NON_NULL(inCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(chip, NULL);
-    psPlane *rc = psPlaneTransformApply(outCoord, chip->toFPA, inCoord);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform chip coords to FPA coords.\n");
+    }
-    return(rc);
+}
-psPlane* pmCoordFPAToChip(
-    psPlane* chipCoord,
-    const psPlane* fpaCoord,
-    const pmChip* chip)
+{
-    PS_ASSERT_PTR_NON_NULL(fpaCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(chip, NULL);
-    PS_ASSERT_PTR_NON_NULL(chip->fromFPA, NULL);
-    psPlane *rc = psPlaneTransformApply(chipCoord, chip->fromFPA, fpaCoord);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform FPA coords to Chip coords.\n");
+    }
-    return(rc);
+}
-psPlane* pmCoordCellToChip(
-    psPlane* outCoord,
-    const psPlane* inCoord,
-    const pmCell* cell)
+{
-    PS_ASSERT_PTR_NON_NULL(inCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(cell, NULL);
-    psPlane *rc = psPlaneTransformApply(outCoord, cell->toChip, inCoord);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform Cell coords to Chip coords.\n");
+    }
-    return(rc);
+}
-psPlane* pmCoordChipToCell(
-    psPlane* cellCoord,
-    const psPlane* chipCoord,
-    const pmCell* cell)
+{
-    PS_ASSERT_PTR_NON_NULL(chipCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(cell, NULL);
-    PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
-    pmCell *tmpCell = pmCellInChip(chipCoord, cell->parent);
-    if (tmpCell == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine the proper cell.\n");
-        return(NULL);
+    }
-    psPlaneTransform *tmpChipToCell = NULL;
-    PS_ASSERT_PTR_NON_NULL(tmpCell->toChip, NULL);
-    if (true ==  p_psIsProjectionLinear(tmpCell->toChip)) {
-        tmpChipToCell = p_psPlaneTransformLinearInvert(tmpCell->toChip);
-    } else {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: non-linear cell->chip transforms are not yet implemented.\n");
-        // XXX: tmpChipToCell = psPlaneTransformInvert(NULL, tmpCell->toChip, NULL, -1);
-        tmpChipToCell = NULL;
+    }
-    if (tmpChipToCell == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not invert the Cell->toChip transform.\n");
-        return(NULL);
+    }
-    psPlane *rc = psPlaneTransformApply(cellCoord, tmpChipToCell, chipCoord);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform Chip coords to Cell coords.\n");
+    }
-    psFree(tmpChipToCell);
-    return(rc);
+}
-psPlane* pmCoordFPAToTP(
-    psPlane* outCoord,
-    const psPlane* inCoord,
-    double color,
-    double magnitude,
-    const pmFPA* fpa)
+{
-    PS_ASSERT_PTR_NON_NULL(inCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(fpa, NULL);
-    psPlane *rc = psPlaneDistortApply(outCoord, fpa->toTangentPlane, inCoord, color, magnitude);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform FPA coords to tangent plane coords.\n");
+    }
-    return(rc);
+}
-psPlane* pmCoordTPToFPA(
-    psPlane* fpaCoord,
-    const psPlane* tpCoord,
-    double color,
-    double magnitude,
-    const pmFPA* fpa)
+{
-    PS_ASSERT_PTR_NON_NULL(tpCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(fpa, NULL);
-    PS_ASSERT_PTR_NON_NULL(fpa->fromTangentPlane, NULL);
-    psPlane *rc = psPlaneDistortApply(fpaCoord, fpa->fromTangentPlane, tpCoord, color, magnitude);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform tangent plane coords to FPA coords.\n");
+    }
-    return(rc);
+}
-/*****************************************************************************
-XXXDeproject(outSphere, coord, projection): This private routine is a wrapper
-for p_psDeproject().  The reason: p_psDeproject() and p_psProject() combined
-do not seem to produce the original coordinates when they even though they
-should.  XXXDeproject() simply negates the ->r and ->d members of the output
-psSphere if the input ->y is larger than 0.0.  I don't know why it works.
-I'm guessing the p_psProject() and p_psDeproject() functions have bugs.
-XXX: It appears that p_psProject() and p_psDeproject() have been fixed.
-Remove this.
- *****************************************************************************/
-psSphere* XXXDeproject(
-    psSphere *outSphere,
-    const psPlane* coord,
-    const psProjection* projection)
+{
-    psSphere *rc = p_psDeproject(outSphere, coord, projection);
-    if (coord->y >= 0.0) {
-        rc->d = -rc->d;
-        rc->r = -rc->r;
+    }
-    return(rc);
+}
-/*****************************************************************************
-  *****************************************************************************/
-psSphere* pmCoordTPToSky(
-    psSphere* outSphere,
-    const psPlane* tpCoord,
-    const psProjection *projection)
+{
-    PS_ASSERT_PTR_NON_NULL(tpCoord, NULL);
-    PS_ASSERT_PTR_NON_NULL(projection, NULL);
-    //    psSphere *rc = XXXDeproject(outSphere, tpCoord, projection);
-    psSphere *rc = p_psDeproject(outSphere, tpCoord, projection);
-    if (rc == NULL) {
-        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform tangent plane coords to sky coords.\n");
+    }
-    return(rc);
+}
 /*****************************************************************************
  *****************************************************************************/
+psPlane* pmCoordSkyToTP(
+    psPlane* tpCoord,
+    const psSphere* in,
+    const psProjection *projection)
+{
+    PS_ASSERT_PTR_NON_NULL(in, NULL);
+    PS_ASSERT_PTR_NON_NULL(projection, NULL);
+    psPlane *rc = p_psProject(tpCoord, in, projection);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform sky to tangent plane coords.\n");
+    }
+    return(rc);
+}
+/*****************************************************************************
+ *****************************************************************************/
+psSphere* pmCoordCellToSky(
+    psSphere* skyCoord,
+    const psPlane* cellCoord,
+    double color,
+    double magnitude,
+    const pmCell* cell)
+{
+    PS_ASSERT_PTR_NON_NULL(cellCoord, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->toFPA, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->parent->parent, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->parent->parent->toTangentPlane, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->parent->parent->projection, NULL);
+    psPlane fpaCoord;
+    psPlane tpCoord;
+    psPlane *rc;
+    pmFPA* parFPA = (cell->parent)->parent;
+    // Convert the input cell coordinates to FPA coordinates.
+    rc = psPlaneTransformApply(&fpaCoord, cell->toFPA, cellCoord);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform cell coords to FPA coords.\n");
+        return(NULL);
+    }
+    // Convert the FPA coordinates to tangent plane Coordinates.
+    rc = psPlaneDistortApply(&tpCoord, parFPA->toTangentPlane, &fpaCoord, color, magnitude);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform FPA coords to tangent plane coords.\n");
+        return(NULL);
+    }
+    // Convert the tangent plane Coordinates to sky coordinates.
+    psSphere *rc2 = pmCoordTPToSky(skyCoord, &tpCoord, parFPA->projection);
+    if (rc2 == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform cell coords to sky coords.\n");
+    }
+    return(rc2);
+}
+/*****************************************************************************
+ *****************************************************************************/
+psPlane* pmCoordSkyToCell(
+    psPlane* cellCoord,
+    const psSphere* skyCoord,
+    float color,
+    float magnitude,
+    const pmCell* cell)
+{
+    PS_ASSERT_PTR_NON_NULL(skyCoord, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->parent, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->parent->parent, NULL);
+    pmChip *parChip = cell->parent;
+    pmFPA *parFPA = parChip->parent;
+    psPlane tpCoord;
+    psPlane fpaCoord;
+    psPlane chipCoord;
+    psPlane *rc;
+    // Convert the skyCoords to tangent plane coords.
+    rc = pmCoordSkyToTP(&tpCoord, skyCoord, parFPA->projection);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine tangent plane coords.\n");
+        return(NULL);
+    }
+    // Convert the tangent plane coords to FPA coords.
+    rc = pmCoordTPToFPA(&fpaCoord, &tpCoord, color, magnitude, parFPA);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine FPA coords.\n");
+        return(NULL);
+    }
+    // Convert the FPA coords to chip coords.
+    rc = pmCoordFPAToChip(&chipCoord, &fpaCoord, parChip);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine chip coords.\n");
+        return(NULL);
+    }
+    // Convert the chip coords to cell coords.
+    rc = pmCoordChipToCell(cellCoord, &chipCoord, cell);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not determine cell coords.\n");
+        return(NULL);
+    }
+    return (cellCoord);
+}
+/*****************************************************************************
+ *****************************************************************************/
+psSphere* pmCoordCellToSkyQuick(
+    psSphere* outSphere,
+    const psPlane* cellCoord,
+    const pmCell* cell)
+{
+    PS_ASSERT_PTR_NON_NULL(cellCoord, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->toSky, NULL);
+    psPlane outPlane;
+    psPlane *rc;
+    rc = psPlaneTransformApply(&outPlane, cell->toSky, cellCoord);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could transform cell coords to sky coords.\n");
+        return(NULL);
+    }
+    psSphere *out = outSphere;
+    if (out == NULL) {
+        out = psSphereAlloc();
+    }
+    out->r = outPlane.y;
+    out->d = outPlane.x;
+    return(out);
+}
+/*****************************************************************************
+ *****************************************************************************/
+psPlane* pmCoordSkyToCellQuick(
+    psPlane* cellCoord,
+    const psSphere* skyCoord,
+    const pmCell* cell)
+{
+    PS_ASSERT_PTR_NON_NULL(skyCoord, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell, NULL);
+    PS_ASSERT_PTR_NON_NULL(cell->toSky, NULL);
+    psPlane skyPlane;
+    skyPlane.y = skyCoord->r;
+    skyPlane.x = skyCoord->d;
+    psPlane *rc = psPlaneTransformApply(cellCoord, cell->toSky, &skyPlane);
+    if (rc == NULL) {
+        psLogMsg(__func__, PS_LOG_WARN, "WARNING: could not transform sky to cell coords.\n");
+    }
+    return(cellCoord);
+// Set cells within a chip to be processed or not
+static bool setCellsProcess(const pmChip *chip, // Chip of interest
+                            bool process  // Process this chip?
+                           )
+{
+    PS_ASSERT_PTR_NON_NULL(chip, false);
+    psArray *cells = chip->cells;       // Component cells
+    if (! cells) {
+        return false;
+    }
+    for (int i = 0; i < cells->n; i++) {
+        pmCell *tmpCell = cells->data[i]; // Cell of interest
+        if (tmpCell) {
+            tmpCell->process = process;
+        }
+    }
+    return true;
+}
 …
+{
     PS_ASSERT_PTR_NON_NULL(fpa, false);
+    if ((fpa->chips == NULL) || (chipNum >= fpa->chips->n)) {
+    psArray *chips = fpa->chips;        // Component chips
+    if ((chips == NULL) || (chipNum >= chips->n)) {
         return(false);
+    }
+    psBool rc = true;
+    for (psS32 i = 0 ; i < fpa->chips->n ; i++) {
+        pmChip *tmpChip = (pmChip *) fpa->chips->data[i];
+    for (int i = 0 ; i < chips->n ; i++) {
+        pmChip *tmpChip = (pmChip *) chips->data[i];
         if (tmpChip == NULL) {
+            rc = false;
+            continue;
+        }
+        if (i == chipNum) {
+            tmpChip->process = true;
+            setCellsProcess(tmpChip, true);
         } else {
+            if (i == chipNum) {
+                tmpChip->valid = true;
+            } else {
+                tmpChip->valid = false;
+            }
+        }
+    }
+    return(rc);
+            tmpChip->process = false;
+            setCellsProcess(tmpChip, false);
+        }
+    }
+    return true;
+}
 …
 /*****************************************************************************
 XXX: The SDRS is ambiguous on a few things:
     Whether or not the other chips should be set valid=true.
     Should we return the number of chip valid=true before or after they're set,
+    Whether or not the other chips should be set process=true. [PAP: No]
+    Should we return the number of chip process=true before or after they're set, [PAP: After]
  *****************************************************************************/
 /**
+ *
  * pmFPAExcludeChip shall set valid to false only for the specified chip
+ *
+ * pmFPAExcludeChip shall set process to false only for the specified chip
  * number (chipNum). In the event that the specified chip number does not exist
  * within the fpa, the function shall generate a warning, and perform no action.
  * The function shall return the number of chips within the fpa that have valid
+ * The function shall return the number of chips within the fpa that have process
  * set to true.
+ *
+ *
  */
 int pmFPAExcludeChip(
 …
     PS_ASSERT_PTR_NON_NULL(fpa, false);
+    if (fpa->chips == NULL) {
+    psArray *chips = fpa->chips;        // Component chips
+    if (chips == NULL) {
         psLogMsg(__func__, PS_LOG_WARN, "WARNING: fpa->chips == NULL\n");
         return(0);
+    }
     if ((chipNum >= fpa->chips->n) || (NULL == (pmChip *) fpa->chips->data[chipNum])) {
+    if ((chipNum >= chips->n) || (NULL == (pmChip *) chips->data[chipNum])) {
         psLogMsg(__func__, PS_LOG_WARN, "WARNING: the specified chip (%d) does not exist.\n", chipNum);
         return(0);
+    }
     psS32 numChips = 0;
     for (psS32 i = 0 ; i < fpa->chips->n ; i++) {
         pmChip *tmpChip = (pmChip *) fpa->chips->data[i];
+    int numChips = 0;                   // Number of chips to be processed
+    for (int i = 0 ; i < chips->n ; i++) {
+        pmChip *tmpChip = (pmChip *) chips->data[i]; // Chip of interest
         if (tmpChip != NULL) {
             if (i == chipNum) {
                 tmpChip->valid = false;
             } else {
                 tmpChip->valid = true;
+                tmpChip->process = false;
+                setCellsProcess(tmpChip, false); // Wipe out the cell as well
+            } else if (tmpChip->process) {
                 numChips++;
+            }
 …
+}
+bool pmCellSetWeights(pmCell *cell // Cell for which to set weights
+                     )
+{
+    float gain = psMetadataLookupF32(NULL, cell->concepts, "CELL.GAIN"); // Cell gain
+    float readnoise = psMetadataLookupF32(NULL, cell->concepts, "CELL.READNOISE"); // Cell read noise
+    psRegion *trimsec = psMetadataLookupPtr(NULL, cell->concepts, "CELL.TRIMSEC"); // Trim section
+    p_pmHDU *hdu = cell->hdu;           // The data unit, containing the weight and mask originals
+    if (!hdu) {
+        pmChip *chip = cell->parent;    // The parent chip
+        if (chip->hdu) {
+            hdu = chip->hdu;
+        } else {
+            pmFPA *fpa = chip->parent;  // The parent FPA
+            hdu = fpa->hdu;
+            if (!hdu) {
+                psError(PS_ERR_UNKNOWN, true, "Unable to find the HDU in the hierarchy!\n");
+                return false;
+            }
+        }
+    }
+    psArray *pixels = hdu->images;      // Array of images
+    psArray *weights = hdu->weights;    // Array of weight images
+    psArray *masks = hdu->masks;        // Array of mask images
+    // Generate the weights and masks if required
+    if (! weights) {
+        weights = psArrayAlloc(pixels->n);
+        for (int i = 0; i < pixels->n; i++) {
+            psImage *image = pixels->data[i];
+            weights->data[i] = psImageAlloc(image->numCols, image->numRows, PS_TYPE_F32);
+            psImageInit(weights->data[i], 0.0);
+        }
+        hdu->weights = weights;
+    }
+    if (! masks) {
+        masks = psArrayAlloc(pixels->n);
+        for (int i = 0; i < pixels->n; i++) {
+            psImage *image = pixels->data[i];
+            masks->data[i] = psImageAlloc(image->numCols, image->numRows, PS_TYPE_U8);
+            psImageInit(masks->data[i], 0);
+        }
+        hdu->masks = masks;
+    }
+    // Set the pixels
+    psArray *readouts = cell->readouts; // Array of readouts
+    for (int i = 0; i < readouts->n; i++) {
+        pmReadout *readout = readouts->data[i]; // The readout of interest
+        if (! readout->weight) {
+            readout->weight = weights->data[i];
+        }
+        if (! readout->mask) {
+            readout->mask = masks->data[i];
+        }
+        // Mask is already set to 0
+        // Set weight image to the variance = g*f + rn^2
+        psImage *image = psImageSubset(readout->image, *trimsec); // The pixels
+        psImage *weight = psImageSubset(readout->weight, *trimsec); // The weight map
+        psBinaryOp(weight, image, "*", psScalarAlloc(gain, PS_TYPE_F32));
+        psBinaryOp(weight, weight, "+", psScalarAlloc(readnoise*readnoise, PS_TYPE_F32));
+    }
+    return true;
+}

trunk/psModules/src/astrom/pmAstrometry.h

-              r6205
+              r6872
 *  @author GLG, MHPCC
+*
 *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2006-01-26 21:10:50 $
+*  @version $Revision: 1.7 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2006-04-17 18:01:04 $
+*
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
 …
+{
     const char *extname;                // Extension name, if it corresponds to this level
-    psArray *pixels;                    // The pixel data, if it corresponds to this level
     psMetadata *header;                 // The FITS header, if it corresponds to this level
+    psArray *images;                    // The pixel data, if it corresponds to this level
+    psArray *masks;                     // The mask data, if it corresponds to this level
+    psArray *weights;                   // The weight data, if it corresponds to this level
+}
 p_pmHDU;
-p_pmHDU *pmHDUAlloc();
-void *pmHDUFree(p_pmHDU *hdu);
 /** Focal plane data structure
+ *
+ *
  *  A focal plane consists of one or more chips (according to the number of
  *  pieces of contiguous silicon). It contains metadata containers for the
 …
  *  transformation which will be derived from numerically inverting the forward
  *  transformation.
+ *
+ *
  */
 typedef struct
 …
 /** Chip data structure
+ *
+ *
  *  A chip consists of one or more cells (according to the number of amplifiers
  *  on the device). The chip contains metadata containers for the concepts and
 …
  *  inverting the forward transformation. A boolean indicates whether the chip is
  *  of interest, allowing it to be excluded from analysis.
+ *
+ *
  */
 typedef struct
 …
     psArray *cells;                     ///< The cells (referred to by name)
     pmFPA *parent;                      ///< Parent FPA
+    bool valid;                         ///< Do we bother about reading and working with this chip?
+    bool process;                       ///< Do we bother about reading and working with this chip?
+    bool exists;                        ///< Does the chip exist --- has it been read in?
     p_pmHDU *hdu;                       ///< FITS data
+}
 …
     psArray *readouts;                  ///< The readouts (referred to by number)
     pmChip *parent;                     ///< Parent chip
+    bool valid;                         ///< Do we bother about reading and working with this cell?
+    bool process;                       ///< Do we bother about reading and working with this cell?
+    bool exists;                        ///< Does the cell exist --- has it been read in?
     p_pmHDU *hdu;                       ///< FITS data
+}
 …
+{
     // Position on the cell
+    // XXX: These may be removed in the future; use parent->concepts instead?
     int col0;                           ///< Offset from the left of chip.
     int row0;                           ///< Offset from the bottom of chip.
 …
     psImage *mask;                      ///< Mask of input image
     psImage *weight;                    ///< Weight of input image
-    psList *bias;                       ///< List of bias section (sub-)images
     psMetadata *analysis;               ///< Readout-level analysis metadata
     pmCell *parent;                     ///< Parent cell
 …
  */
 pmCell *pmCellAlloc(
     pmChip *chip,                        ///< Parent chip
     psMetadata *cameradata,
     psString name
+    pmChip *chip,       ///< Parent chip
+    psMetadata *cameraData, ///< Camera data
+    const char *name    ///< Name of cell
 );
 /** Allocates a pmChip
+ *
+ *
  *  The constructor shall make an empty pmChip. If the parent fpa is not NULL,
  *  the parent link is made and the chip shall be placed in the parent's array
 …
  */
 pmChip *pmChipAlloc(
+    pmFPA *fpa,
+    psString name
+    pmFPA *fpa,                         ///< FPA to which the chip belongs
+    const char *name                    ///< Name of chip
 );
 /** Allocates a pmFPA
+ *
+ *
  *  The constructor shall make an empty pmFPA. The chips array shall be
  *  allocated with a zero size, the camera and db pointers set to the values
  *  provided, and the concepts metadata constructed. All other pointers in the
  *  structure shall be initialized to NULL.
+ *
+ *
  */
 pmFPA *pmFPAAlloc(
 …
 );
+/** Allocates a p_pmHDU
+ *
+ * XXX: More detailed description
+ *
+ * @return p_pmHDU*    newly allocated p_pmHDU
+ */
+p_pmHDU *p_pmHDUAlloc(const char *extname // Extension name
+                     );
 /** Verify parent links.
+ *
+ *
  *  This function checks the validity of the parent links in the FPA hierarchy.
  *  If a parent link is not set (or not set correctly), it is corrected, and the
  *  function shall return false. If all the parent pointers were correct, the
  *  function shall return true.
+ *
+ *
  */
 bool pmFPACheckParents(
 …
 /** FUNC DESC
+ *
+ *
+ *
+ *
+ */
+ *
+ *
+ *
+ *
+ */
+/**
+ *
+ * pmFPASelectChip shall set valid to true for the specified chip number
+ * (chipNum), and all other chips shall have valid set to false. In the event
+ * that the specified chip number does not exist within the fpa, the function
+ * shall return false.
+ *
+ */
+bool pmFPASelectChip(
+    pmFPA *fpa,
+    int chipNum
+);
+/**
+ *
+ * pmFPAExcludeChip shall set valid to false only for the specified chip
+ * number (chipNum). In the event that the specified chip number does not exist
+ * within the fpa, the function shall generate a warning, and perform no action.
+ * The function shall return the number of chips within the fpa that have valid
+ * set to true.
+ *
+ */
+int pmFPAExcludeChip(
+    pmFPA *fpa,
+    int chipNum
+);
 …
 );
-/**
+ *
- * pmFPASelectChip shall set valid to true for the specified chip number
- * (chipNum), and all other chips shall have valid set to false. In the event
- * that the specified chip number does not exist within the fpa, the function
- * shall return false.
+ *
- */
-bool pmFPASelectChip(
-    pmFPA *fpa,
-    int chipNum
-);
-/**
+ *
- * pmFPAExcludeChip shall set valid to false only for the specified chip
- * number (chipNum). In the event that the specified chip number does not exist
- * within the fpa, the function shall generate a warning, and perform no action.
- * The function shall return the number of chips within the fpa that have valid
- * set to true.
+ *
- */
-int pmFPAExcludeChip(
-    pmFPA *fpa,
-    int chipNum
-);
-/**
+ *
- * pmFPAConstruct shall construct a focal plane hierarchy from a camera
- * configuration. A db handle is also provided so that may be set in the pmFPA.
- * The resultant pmFPA and its lower-down components shall be ready for to read a
- * FITS file into it by setting the extname pointers at the appropriate levels to
- * the appropriate FITS extension name.
+ *
- */
-pmFPA *pmFPAConstruct(
-    const psMetadata *camera,
-    psDB *db
-);
 #endif // #ifndef PS_ASTROMETRY_H

trunk/psModules/src/astrom/pmAstrometryObjects.c

-              r6511
+              r6872
 *  @author EAM, IfA
+*
 *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2006-03-04 01:01:33 $
+*  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2006-04-17 18:01:04 $
+*
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
 …
 #include <math.h>
 #include "pslib.h"
 #include "pmAstrometry.h"
+#include "pmFPA.h"
 #include "pmAstrometryObjects.h"
 …
     psArray *rot;
+    pmAstromStats minStat, newStat;
+    pmAstromStats minStat;
+    pmAstromStats newStat = {{0, 0, 0, 0}, {0, 0, 0, 0}, 0, 0, 0, 0};
     psPlane center;

trunk/psModules/src/astrom/pmAstrometryObjects.h

-              r6205
+              r6872
 *  @author EAM, IfA
+*
 *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
 *  @date $Date: 2006-01-26 21:10:50 $
+*  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+*  @date $Date: 2006-04-17 18:01:04 $
+*
 *  Copyright 2004-2005 Maui High Performance Computing Center, University of Hawaii
 …
 # include <unistd.h>   // for unlink
 # include <pslib.h>
 # include <pmAstrometry.h>
 /*
+ *
+# include <pmFPA.h>
+/*
+ *
  * This structure specifies the coordinate of the detection in each of the
  * four necessary coordinate frames: pix defines the position in the psReadout
 …
  * coordinates, while the reference detections will be projected to the other
  * frames from the sky coordinates.
+ *
+ *
  * XXX: There are more members here than in the SDRS.
+ *
+ *
  */
 typedef struct
 …
 /*
+ *
+ *
  * The pmAstromMatch structure defines the cross-correlation between two
  * arrays. A single such data item specifies that item number pmAstromMatch.idx1
  * in the first list corresponds to pmAstromMatch.idx2 in the second list.
+ *
+ *
  */
 typedef struct
 …
 /*
+ *
+ *
  * XXX: Not in SDRS.
+ *
+ *
  */
 typedef struct
 …
 /*
+ *
+ *
  * If the two sets of coordinates are expected to agree very well (ie, the
  * current best-guess astrometric solution is quite close to the radius. The
 …
  * ASTROM.MATCH.RADIUS). The output consists an array of pmAstromMatch values,
  * defined below.
+ *
+ *
  */
 psArray *pmAstromRadiusMatch(
 …
 /*
+ *
+ *
  * This function accepts an array of pmAstromObj objects and rotates them by
  * the given angle about the given center coordinate pCenter,qCenter in the Focal
  * Plane Array coordinates.
+ *
+ *
  * XXX: This differs from the SDRS
+ *
+ *
  */
 /* SDRS
 …
 /*
+ *
+ *
  * If the two sets of coordinates are not known to agree well, but the
  * relative scale and approximate relative rotation is known, then a much faster
 …
  * allowing the procedure to scan over a range of rotations. We define the
  * following function to apply this matching algorithm:
+ *
+ *
  * XXX: In the SDRS, this function is a pointer.
+ *
+ *
  */
 pmAstromStats pmAstromGridMatch(
 …
 /*
+ *
+ *
  * The result of a pmAstromGridMatch may be used to modify the astrometry
  * transformation information for a pmFPA image hierarchy structure. The result
 …
  * the linear terms of the pmFPA.toTangentPlane transformation. These two
  * adjustments are made using the function:
+ *
+ *
  * XXX: This function name is different in the SDRS.
+ *
+ *
  */
 psPlaneTransform *pmAstromGridApply(
 …
 /*
+ *
+ *
  * This function is identical to pmAstromGridMatch, but is valid for only a
  * single relative rotation. The input config information need not contain any of
  * the GRID.*.ANGLE entries (they will be ignored).
+ *
+ *
+ * XXX: This function name is different in the SDRS.
+ *
  */
 /* in pmAstromGrid.c */
 …
 /*
+ *
+ *
  * This function accepts the raw and reference source lists and the list of
  * matched entries. It uses the matched list to determine a polynomial
 …
  * modifications to pmFPA.toTangentPlane incorporate the rotation component of
  * the linear terms and the higher-order terms of the polynomial fits.
+ *
+ *
  * XXX: No prototype code.
+ *
+ *
  */
 bool pmAstromFitFPA(
 …
  *ASTROM.ORDER).  The result of this fit is a set of modifications of the
  *components of the pmChip.toFPA transformation.
+ *
+ *
  * XXX: No prototype code.
+ *
+ *
  */
 bool pmAstromFitChip(
 …
 /*
+ *
+ *
  * The following function determines the position residual, in the tangent
  * plane, as a function of position in the focal plane, for a collection of raw
 …
  * the bin size over which the gradient is measured (keyword: ASTROM.GRAD.BOX).
  * The function returns an array of pmAstromGradient structures, defined below.
+ *
+ *
  * XXX: No prototype code.
+ *
+ *
  */
 psArray pmAstromMeasureGradients(
 …
 /*
+ *
+ *
  * The following data structure carries the information about the residual
  * gradient of source positions in the tangent plane (pmAstromObj.TP) as a
  * function of position in the focal plane (pmAstromObj.FP).
+ *
+ *
  */
 typedef struct
 …
 /*
+ *
+ *
  * The gradient set measured above can be fitted with a pair of 2D
  * polynomials. The resulting fits can then be related back to the implied
 …
  * supplied pmFPA structure. The configuration variable supplies the polynomial
  * order (keyword: ASTROM.DISTORT.ORDER).
+ *
+ *
  * XXX: No prototype code.
+ *
+ *
  */
 psArray pmAstromFitDistortion(
 …
  ******************************************************************************/
 /*
+ *
+ *
+ *
+ *
  */
 /*
+ *
+ *
+ *
+ *
+ *
+ */
+/*
+ *
  * Allocates a pmAstromObj struct.
+ *
+ *
  */
 pmAstromObj *pmAstromObjAlloc (void);
 …
 /*
+ *
+ *
  * Copies a pmAstromObj struct.
+ *
+ *
  */
 pmAstromObj *pmAstromObjCopy(
 …
 /*
+ *
+ *
+ *
+ *
+ *
+ *
  */
 pmAstromMatch *pmAstromMatchAlloc(
 …
 /*
+ *
+ *
+ *
+ *
+ *
+ *
  */
 psPlaneTransform *pmAstromMatchFit(
 …
 /*
+ *
+ *
+ *
+ *
+ *
+ *
  */
 int pmAstromObjSortByFPX(
 …
 /*
+ *
+ *
+ *
+ *
+ *
+ *
  */
 int pmAstromObjSortByMag(

trunk/psModules/src/camera/Makefile.am

-              r5170
+              r6872
 libpsmodulecamera_la_CPPFLAGS = $(SRCINC) $(PSMODULE_CFLAGS)
 libpsmodulecamera_la_LDFLAGS  = -release $(PACKAGE_VERSION)
+libpsmodulecamera_la_SOURCES  =
+libpsmodulecamera_la_SOURCES  = \
+        pmChipMosaic.c \
+        pmFPAConceptsGet.c \
+        pmFPAConceptsSet.c \
+        pmFPAConstruct.c \
+        pmFPAMorph.c \
+        pmFPARead.c \
+        pmFPAWrite.c
 psmoduleincludedir = $(includedir)
+psmoduleinclude_HEADERS =
+psmoduleinclude_HEADERS = \
+        pmChipMosaic.h \
+        pmFPAConceptsGet.h \
+        pmFPAConceptsSet.h \
+        pmFPAConstruct.h \
+        pmFPAMorph.h \
+        pmFPARead.h \
+        pmFPAWrite.h

trunk/psModules/src/config/pmConfig.h

-              r6297
+              r6872
  *  @author PAP, IfA
+ *
  *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2006-02-02 04:51:14 $
+ *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
+// Configuration information
+typedef struct
+{
+    psMetadata *site;                   // Site configuration
+    psMetadata *camera;                 // Camera specification
+    psMetadata *recipes;                // Recipes for processing
+    psMetadata *arguments;              // Command-line arguments
+    psMetadata *files;                  // pmFPAfiles used for analysis
+    //    psMetadata *format;   // camera format for primary image
+    psDB *database;                     // Database handle
+}
+pmConfig;
+pmConfig *pmConfigAlloc(void);
+// Free static variables
+void pmConfigDone(void);
 /** pmConfigRead
+ *
+ *
  * pmConfigRead shall load the site configuration (according to the above rule
  * for determining the source). The camera configuration shall also be loaded if
 …
+ *
  */
+bool pmConfigRead(
+    psMetadata **site,
+    psMetadata **camera,
+    psMetadata **recipe,
+pmConfig *pmConfigRead(
     int *argc,
+    char **argv,
+    const char *recipeName
+    char **argv);
+/** we need this elsewhere; make it public **/
+bool readConfig(
+    psMetadata **config,                // Config to output
+    const char *name,                   // Name of file
+    const char *description             // Description of file
 );
 /** pmConfigValidateCamera
+ *
+ *
  * This function, used by pmConfigCameraFromHeader, shall return true if the
  * FITS header matches the rule contained in the camera configuration (see
  * x2.2.2.3); otherwise it shall return false.
+ *
+ *
  */
+bool pmConfigValidateCamera(
+    const psMetadata *camera,
+    const psMetadata *header
+);
+bool pmConfigValidateCameraFormat(
+    const psMetadata *cameraFormat,
+    const psMetadata *header);
 /** pmConfigCameraFromHeader
+ *
+ *
  * pmConfigCameraFromHeader shall load the camera configuration based on the
  * contents of the FITS header, using the list of known cameras contained in the
  * site configuration. If more than one camera matches the FITS header, a warning
  * shall be generated and the first matching camera returned.
+ *
+ *
  */
 psMetadata *pmConfigCameraFromHeader(
     const psMetadata *site,
     const psMetadata *header
+psMetadata *pmConfigCameraFormatFromHeader(
+    pmConfig *config,                   // The configuration
+    const psMetadata *header           // The FITS header
 );
 /** pmConfigRecipeFromCamera
+ *
+ * pmConfigRecipeFromCamera shall load the recipe configuration based on the
+ * recipeName and the list of known recipes contained in the camera
+ *
+ * pmConfigRecipeFromCamera shall load the recipes from the list of known recipes contained in the camera
  * configuration.
+ *
+ *
  */
+psMetadata *pmConfigRecipeFromCamera(
+    const psMetadata *camera,
+    const char *recipeName
+bool pmConfigReadRecipes(
+    pmConfig *config
 );
 /** pmConfigDB
+ *
+ *
  * pmConfigDB shall use the site configuration data to open a database handle.
  * This is fairly straightforward at the moment, but will change when we beef up
  * security. (TBD)
+ *
+ *
  */
 #ifdef DOMIT_PSDB
 …
 );
+/** pmConfigConformHeader
+ *
+ * Make the supplied header conform to the nominated camera format.
+ */
+bool pmConfigConformHeader(psMetadata *header, // Header to conform
+                           const psMetadata *format // Camera format
+                          );
+psArray *pmConfigFileSets (int *argc, char **argv, char *file, char *list);
+bool pmConfigFileSetsMD (psMetadata *metadata, int *argc, char **argv, char *name, char *file, char *list);
 #endif

trunk/psModules/src/detrend/Makefile.am

-              r5170
+              r6872
 libpsmoduledetrend_la_CPPFLAGS = $(SRCINC) $(PSMODULE_CFLAGS)
 libpsmoduledetrend_la_LDFLAGS  = -release $(PACKAGE_VERSION)
+libpsmoduledetrend_la_SOURCES  = pmFlatField.c \
+    pmMaskBadPixels.c \
+    pmNonLinear.c
+libpsmoduledetrend_la_SOURCES  = \
+        pmFlatField.c \
+        pmFringeStats.c \
+        pmMaskBadPixels.c \
+        pmNonLinear.c
 psmoduleincludedir = $(includedir)
 psmoduleinclude_HEADERS = \
+  pmFlatField.h \
+  pmFlatFieldErrors.h \
+  pmMaskBadPixelsErrors.h \
+  pmMaskBadPixels.h \
+  pmNonLinear.h
+        pmFlatField.h \
+        pmFlatFieldErrors.h \
+        pmFringeStats.h \
+        pmMaskBadPixelsErrors.h \
+        pmMaskBadPixels.h \
+        pmNonLinear.h
 EXTRA_DIST = pmFlatFieldErrors.dat pmMaskBadPixelsErrors.dat

trunk/psModules/src/detrend/pmFlatField.c

-              r5543
+              r6872
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
 /** @file  pmFlatField.c
+ *
 …
  *  @author Ross Harman, MHPCC
+ *
  *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-11-18 19:43:14 $
+ *  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include<stdio.h>
 #include<math.h>
 #include <string.h>
+#include <strings.h>
 #include "pslib.h"
 …
 #include "pmMaskBadPixels.h"
 #include "pmFlatFieldErrors.h"
-#include "pmSubtractBias.h"
-// XXX: This should be removed when the autoconf stuff handles psConstants.h correctly.
-#define PS_WARN_PTR_NON_NULL(NAME) \
-if ((NAME) == NULL) { \
-    psLogMsg(__func__, PS_LOG_WARN, "WARNING: %s is NULL.", #NAME); \
-} \
+bool pmFlatField(
+    pmReadout *in,
+    const pmReadout *flat)
+bool pmFlatField(pmReadout *in, const pmReadout *flat)
+{
-    // XXX: Use the proper image and readout asserts.
-    PS_ASSERT_PTR_NON_NULL(in, true);
-    PS_ASSERT_PTR_NON_NULL(in->image, false);
-    PS_ASSERT_PTR_NON_NULL(in->mask, false);
-    PS_ASSERT_PTR_NON_NULL(flat, false);
-    PS_ASSERT_PTR_NON_NULL(flat->image, false);
-    if (in == NULL)
-        printf("XXX: NULL\n");
-    // XXX: Not sure if this is correct.  Must consult with IfA.
-    PS_ASSERT_PTR_NON_NULL(in->mask, false);
-    PS_WARN_PTR_NON_NULL(in->parent);
-    if (in->parent != NULL) {
-        PS_WARN_PTR_NON_NULL(in->parent->concepts);
+    }
     int i = 0;
     int j = 0;
 …
     psElemType flatType;
     psElemType maskType;
-    psImage *inMask = NULL;
-    psImage *flatImage = NULL;
+    //
+    // Determine trimmed image from metadata.
+    //
+    psImage *trimmedImg = p_psDetermineTrimmedImage(in);
+    flatImage = flat->image;
+    inMask = in->mask;
+    // Check for nulls
+    if (in == NULL) {
+        return true;       // Readout may not have data in it
+    } else if(flat==NULL) {
+        psError( PS_ERR_BAD_PARAMETER_NULL, true,
+                 PS_ERRORTEXT_pmFlatField_NULL_FLAT_READOUT);
+        return false;
+    }
+    psImage *inImage   = in->image;     // Input image
+    psImage *inMask    = in->mask;      // Mask for input image
+    psImage *flatImage = flat->image;   // Flat-field image
+    // Offsets on the chip
+    int x0in = psMetadataLookupS32(NULL, in->parent->concepts, "CELL.X0");
+    int y0in = psMetadataLookupS32(NULL, in->parent->concepts, "CELL.Y0");
+    int x0flat = psMetadataLookupS32(NULL, flat->parent->concepts, "CELL.X0");
+    int y0flat = psMetadataLookupS32(NULL, flat->parent->concepts, "CELL.Y0");
+    if (inImage == NULL) {
+        psError( PS_ERR_BAD_PARAMETER_NULL, true,
+                 PS_ERRORTEXT_pmFlatField_NULL_INPUT_IMAGE);
+        return false;
+    } else if(flatImage == NULL) {
+        psError( PS_ERR_BAD_PARAMETER_NULL, true,
+                 PS_ERRORTEXT_pmFlatField_NULL_FLAT_IMAGE);
+        return false;
+    }
     // Check input image and its mask are not larger than flat image
+    if (trimmedImg == NULL)
+        printf("XXX: 00\n");
+    if (flatImage == NULL)
+        printf("XXX 01\n");
+    if (trimmedImg->numRows>flatImage->numRows || trimmedImg->numCols>flatImage->numCols) {
+    if (inImage->numRows>flatImage->numRows || inImage->numCols>flatImage->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmFlatField_SIZE_INPUT_IMAGE,
                  trimmedImg->numRows, trimmedImg->numCols, flatImage->numRows, flatImage->numCols);
+                 inImage->numRows, inImage->numCols, flatImage->numRows, flatImage->numCols);
         return false;
+    }
     if (inMask->numRows > flatImage->numRows || inMask->numCols > flatImage->numCols) {
+    if (inMask && (inMask->numRows > flatImage->numRows || inMask->numCols > flatImage->numCols)) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmFlatField_SIZE_MASK_IMAGE,
 …
     // Determine total offset based on image offset with chip offset
     totOffCol = trimmedImg->col0 + in->col0;
     totOffRow = trimmedImg->row0 + in->row0;
+    totOffCol = inImage->col0 + y0in - flatImage->col0 - y0flat;
+    totOffRow = inImage->row0 + x0in - flatImage->row0 - x0flat;
     // Check that offsets are within image limits
 …
                  totOffRow, totOffCol, flatImage->numRows, flatImage->numCols);
         return false;
     } else if(totOffRow>=trimmedImg->numRows || totOffCol>=trimmedImg->numCols) {
+    } else if(totOffRow>=inImage->numRows || totOffCol>=inImage->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmFlatField_OFFSET_INPUT_IMAGE,
                  totOffRow, totOffCol, trimmedImg->numRows, trimmedImg->numCols);
+                 totOffRow, totOffCol, inImage->numRows, inImage->numCols);
         return false;
     } else if(totOffRow>=inMask->numRows || totOffCol>=inMask->numCols) {
+    } else if(inMask && (totOffRow>=inMask->numRows || totOffCol>=inMask->numCols)) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmFlatField_OFFSET_MASK_IMAGE,
 …
     // Check for incorrect types
     inType = trimmedImg->type.type;
+    inType = inImage->type.type;
     flatType = flatImage->type.type;
     maskType = inMask->type.type;
 …
                  flatType);
         return false;
     } else if(maskType != PS_TYPE_MASK) {
+    } else if(inMask && inMask->type.type != PS_TYPE_MASK) {
         psError( PS_ERR_BAD_PARAMETER_TYPE, true,
                  PS_ERRORTEXT_pmFlatField_TYPE_MASK_IMAGE,
 …
 case PS_TYPE_##TYPE:                                                                                         \
     /* Per Eugene's request, use two sets of loops: first to fill mask, second to avoid div with bad pix */  \
     for(j = totOffRow; j < trimmedImg->numRows; j++) {                                                          \
         for(i = totOffCol; i < trimmedImg->numCols; i++) {                                                      \
+    for(j = totOffRow; j < inImage->numRows; j++) {                                                          \
+        for(i = totOffCol; i < inImage->numCols; i++) {                                                      \
             if(flatImage->data.TYPE[j][i] <= 0.0) {                                                          \
                 /* Negative or zero flat pixels shall be masked in input image as  PM_MASK_FLAT */           \
+                inMask->data.PS_TYPE_MASK_DATA[j][i] |= PM_MASK_FLAT;                                        \
+                if (inMask) {                                                                                \
+                    inMask->data.PS_TYPE_MASK_DATA[j][i] |= PM_MASK_FLAT;                                    \
+                }                                                                                            \
                 flatImage->data.TYPE[j][i] = 0.0;                                                            \
             }                                                                                                \
         }                                                                                                    \
     }                                                                                                        \
     for(j = totOffRow; j < trimmedImg->numRows; j++) {                                                          \
         for(i = totOffCol; i < trimmedImg->numCols; i++) {                                                      \
             if(!inMask->data.PS_TYPE_MASK_DATA[j][i]) {                                                      \
+    for(j = totOffRow; j < inImage->numRows; j++) {                                                          \
+        for(i = totOffCol; i < inImage->numCols; i++) {                                                      \
+            if(inMask && !inMask->data.PS_TYPE_MASK_DATA[j][i]) {                                            \
                 /* Module shall divide the input image by the flat-fielded image */                          \
                 trimmedImg->data.TYPE[j][i] /= flatImage->data.TYPE[j][i];                                      \
+                inImage->data.TYPE[j][i] /= flatImage->data.TYPE[j][i];                                      \
             }                                                                                                \
         }                                                                                                    \

trunk/psModules/src/detrend/pmFlatField.h

-              r5435
+              r6872
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
 /** @file  pmFlatField.h
+ *
 …
  *  @author Ross Harman, MHPCC
+ *
  *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-10-20 23:06:24 $
+ *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include "pslib.h"
 #include "pmAstrometry.h"
+#include "pmFPA.h"

trunk/psModules/src/detrend/pmMaskBadPixels.c

-              r5543
+              r6872
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
 /** @file  pmMaskBadPixels.c
+ *
 …
  *  @author Ross Harman, MHPCC
+ *
  *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-11-18 19:43:14 $
+ *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include<stdio.h>
 #include<math.h>
 #include<string.h>
+#include<strings.h>
 #include "pmMaskBadPixels.h"
 #include "pmMaskBadPixelsErrors.h"
-#include "pmSubtractBias.h"
+//XXX: REmove, autoconf is broken.
+#define PS_WARN_PTR_NON_NULL(NAME) \
+if ((NAME) == NULL) { \
+    psLogMsg(__func__, PS_LOG_WARN, "WARNING: %s is NULL.", #NAME); \
+} \
+bool pmMaskBadPixels(pmReadout *in, const psImage *mask, unsigned int maskVal, float sat,
+                     unsigned int growVal, int grow)
+{
+    int i = 0;
+    int j = 0;
+    int jj = 0;
+    int ii = 0;
+    int rRound = 0;
+    int rowMin = 0;
+    int rowMax = 0;
+    int colMin = 0;
+    int colMax = 0;
+    int totOffCol = 0;
+    int totOffRow = 0;
+    float r = 0.0f;
+    psElemType inType;
+    psElemType maskType;
+    psImage *inImage = NULL;
+    psImage *inMask = NULL;
+/******************************************************************************
+GrowPixel(inMask, row, col, radius, maskVal): This private routine takes an
+input image mask and a pixel location, then sets (logical or) all pixels with
+parameter radius if that pixel to maskVal parameter.
+*****************************************************************************/
+psBool GrowPixel(
+    psImage *inMask,
+    psS32 col,
+    psS32 row,
+    psS32 radius,
+    psU32 maskVal)
+{
+    psS32 rowMin = PS_MAX(row-radius, 0);
+    psS32 rowMax = PS_MIN(row+radius+1, inMask->numRows);
+    psS32 colMin = PS_MAX(col-radius, 0);
+    psS32 colMax = PS_MIN(col+radius+1, inMask->numCols);
+    psF32 squareRadius = PS_SQR(radius);
+    // Check for nulls
+    if (in == NULL) {
+        return true;   // Readout may not have data in it
+    } else if(mask==NULL) {
+        psError( PS_ERR_BAD_PARAMETER_NULL, true,
+                 PS_ERRORTEXT_pmMaskBadPixels_NULL_MASK_IMAGE);
+        return false;
+    }
+    for(psS32 i=rowMin; i<rowMax; i++) {
+        for(psS32 j=colMin; j<colMax; j++) {
+            // Old code:
+            //            psF32 rRound = 0.5 + sqrtf((psF32) (((col-j)*(col-j)) + ((row-i)*(row-i))));
+            //            if(rRound <= radius) {
+            psF32 squareDis = (psF32) (((col-j)*(col-j)) + ((row-i)*(row-i)));
+            if (squareDis <= squareRadius) {
+                inMask->data.U8[i][j] |= maskVal;
+            }
+        }
+    }
+    return(true);
+}
+/******************************************************************************
+GetRadius(inImg, col, row, sat, growVal, grow): This private routine takes an
+input image and pixel location and determines what radius that pixel must grow
+by.
+//XXX: Inline this or macro it for speed.
+ *****************************************************************************/
+static psS32 GetRadius(
+    psImage *inImg,
+    psS32 col,
+    psS32 row,
+    psF32 sat,
+    psU32 growVal,
+    psS32 grow)
+{
+    psS32 growRadius = 0;
+    if (inImg->type.type == PS_TYPE_F32) {
+        if(inImg->data.F32[row][col] == growVal) {
+            growRadius = grow;
+        }
+        if (inImg->data.F32[row][col] > sat) {
+            growRadius = PS_MAX(growRadius+1, 2);
+        }
+    } else if (inImg->type.type == PS_TYPE_S32) {
+        if(inImg->data.S32[row][col] == growVal) {
+            growRadius = grow;
+        }
+        if (inImg->data.S32[row][col] > sat) {
+            growRadius = PS_MAX(growRadius+1, 2);
+        }
+    } else if (inImg->type.type == PS_TYPE_U16) {
+        if(inImg->data.U16[row][col] == growVal) {
+            growRadius = grow;
+        }
+        if (inImg->data.U16[row][col] > sat) {
+            growRadius = PS_MAX(growRadius+1, 2);
+        }
+    } else {
+        psError( PS_ERR_BAD_PARAMETER_TYPE, true,
+                 PS_ERRORTEXT_pmMaskBadPixels_TYPE_MASK_IMAGE,
+                 inImg->type.type);
+    }
+    return(growRadius);
+}
+bool pmMaskBadPixels(
+    pmReadout *in,
+    const pmReadout *mask,
+    unsigned int maskVal,
+    float sat,
+    unsigned int growVal,
+    int grow)
+{
+    // XXX: Review this code then put all asserts at the top.
+    PS_ASSERT_PTR_NON_NULL(in, true);
+    PS_ASSERT_PTR_NON_NULL(in->image, false);
+    PS_WARN_PTR_NON_NULL(in->parent);
+    if (in->parent != NULL) {
+        PS_WARN_PTR_NON_NULL(in->parent->concepts);
+    }
+    PS_ASSERT_PTR_NON_NULL(mask, false);
+    int i = 0;
+    int j = 0;
+    int totOffCol = 0;
+    int totOffRow = 0;
+    psElemType inType;
+    psElemType maskType;
+    psImage *inMask = NULL;
+    //
+    // Determine trimmed image from metadata.
+    //
+    psImage *trimmedImg = p_psDetermineTrimmedImage(in);
+    if(in->mask == NULL) {
+        in->mask = psImageAlloc(trimmedImg->numCols, trimmedImg->numRows, PS_TYPE_MASK);
+        PS_IMAGE_SET_U8(in->mask, 0);
+    inImage = in->image;
+    if (inImage == NULL) {
+        psError( PS_ERR_BAD_PARAMETER_NULL, true,
+                 PS_ERRORTEXT_pmMaskBadPixels_NULL_INPUT_IMAGE);
+        return false;
+    } else if(in->mask == NULL) {
+        in->mask = psImageAlloc(inImage->numCols, inImage->numRows, PS_TYPE_MASK);
+        memset(in->mask->data.V[0], 0, inImage->numCols*inImage->numRows*PSELEMTYPE_SIZEOF(PS_TYPE_MASK));
+    }
     inMask = in->mask;
     // Check input image and its mask are not larger than mask
     if(trimmedImg->numRows > mask->image->numRows || trimmedImg->numCols > mask->image->numCols) {
+    if(inImage->numRows > mask->numRows || inImage->numCols > mask->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmMaskBadPixels_SIZE_INPUT_IMAGE,
                  trimmedImg->numRows, trimmedImg->numCols, mask->image->numRows, mask->image->numCols);
+                 inImage->numRows, inImage->numCols, mask->numRows, mask->numCols);
         return false;
     } else if(inMask->numRows > mask->image->numRows || inMask->numCols > mask->image->numCols) {
+    } else if(inMask->numRows>mask->numRows || inMask->numCols>mask->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmMaskBadPixels_SIZE_MASK_IMAGE,
                  inMask->numRows, inMask->numCols, mask->image->numRows, mask->image->numCols);
+                 inMask->numRows, inMask->numCols, mask->numRows, mask->numCols);
         return false;
+    }
     // Determine total offset based on image offset with chip offset
     totOffCol = trimmedImg->col0 + in->col0;
     totOffRow = trimmedImg->row0 + in->row0;
+    totOffCol = inImage->col0; // + in->col0;
+    totOffRow = inImage->row0; // + in->row0;
     // Check that offsets are within image limits
     if(totOffRow>=mask->image->numRows || totOffCol>=mask->image->numCols) {
+    if(totOffRow>=mask->numRows || totOffCol>=mask->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmMaskBadPixels_OFFSET_MASK_IMAGE,
                  totOffRow, totOffCol, mask->image->numRows, mask->image->numCols);
+                 totOffRow, totOffCol, mask->numRows, mask->numCols);
         return false;
     } else if(totOffRow>=trimmedImg->numRows || totOffCol>=trimmedImg->numCols) {
+    } else if(totOffRow>=inImage->numRows || totOffCol>=inImage->numCols) {
         psError( PS_ERR_BAD_PARAMETER_SIZE, true,
                  PS_ERRORTEXT_pmMaskBadPixels_OFFSET_INPUT_IMAGE,
                  totOffRow, totOffCol, trimmedImg->numRows, trimmedImg->numCols);
+                 totOffRow, totOffCol, inImage->numRows, inImage->numCols);
         return false;
     } else if(totOffRow>=inMask->numRows || totOffCol>=inMask->numCols) {
 …
     // Check for incorrect types
     inType = trimmedImg->type.type;
     maskType = mask->image->type.type;
+    inType = inImage->type.type;
+    maskType = mask->type.type;
     if(PS_IS_PSELEMTYPE_COMPLEX(inType)) {
         psError( PS_ERR_BAD_PARAMETER_TYPE, true,
 …
+    }
+    //
+    // This is the main loop which examines each pixel and masks pixels if
+    //    A: The mask matches the maskValue
+    //    B: The pixel is larger than the saturation value
+    //    C: The pixel equals the grow value (in which case a circle is masked)
+    //
+    for(i=totOffRow; i<trimmedImg->numRows; i++) {
+        for(j=totOffCol; j<trimmedImg->numCols; j++) {
+            //
+            // A: Pixels which satisfy maskVal shall be masked
+            //
+            if (mask->image->data.U8[i][j] == maskVal) {
+                in->mask->data.U8[i][j] |= maskVal;
+            }
+    // Macro for all PS types
+    #define PM_BAD_PIXELS(TYPE)                                                                              \
+case PS_TYPE_##TYPE:                                                                                         \
+    for(j=totOffRow; j<inImage->numRows; j++) {                                                              \
+        for(i=totOffCol; i<inImage->numCols; i++) {                                                          \
+            \
+            /* Pixels with flux greater than sat shall be masked */                                          \
+            if(inImage->data.TYPE[j][i] > sat) {                                                             \
+                inMask->data.PS_TYPE_MASK_DATA[j][i] |= PM_MASK_SAT;                                         \
+            }                                                                                                \
+            \
+            /* Pixels which satisfy maskVal shall be masked */                                               \
+            inMask->data.PS_TYPE_MASK_DATA[j][i] |= (mask->data.PS_TYPE_MASK_DATA[j][i]&maskVal);            \
+            \
+            /* Pixels which satisfy growVal and within the grow radius shall be masked */                    \
+            if(mask->data.PS_TYPE_MASK_DATA[j][i] & growVal) {                                               \
+                rowMin = MAX(j-grow, 0);                                                                     \
+                rowMax = MIN(j+grow+1, inImage->numRows);                                                    \
+                colMin = MAX(i-grow, 0);                                                                     \
+                colMax = MIN(i+grow+1, inImage->numCols);                                                    \
+                for(jj=rowMin; jj<rowMax; jj++) {                                                            \
+                    for(ii=colMin; ii<colMax; ii++) {                                                        \
+                        r = sqrtf((ii-i)*(ii-i)+(jj-j)*(jj-j));                                              \
+                        rRound = r + 0.5;                                                                    \
+                        if(rRound <= grow) {                                                                 \
+                            inMask->data.PS_TYPE_MASK_DATA[jj][ii] |=                                        \
+                                    (mask->data.PS_TYPE_MASK_DATA[j][i]&growVal);                            \
+                        }                                                                                    \
+                    }                                                                                        \
+                }                                                                                            \
+            }                                                                                                \
+        }                                                                                                    \
+    }                                                                                                        \
+    break;
+            //
+            // We first determine how much we need to grow by and store this in
+            // growRadius.
+            //
+            psS32 growRadius = GetRadius(trimmedImg, j, i, sat, growVal, grow);
+            //
+            // Grow the pixel mask if necessary.
+            //
+            if (growRadius != 0) {
+                GrowPixel(in->mask, j, i, growRadius, maskVal);
+            }
+        }
+    // Switch to call bad pixel masking macro defined above
+    switch(inType) {
+        PM_BAD_PIXELS(U8);
+        PM_BAD_PIXELS(U16);
+        PM_BAD_PIXELS(U32);
+        PM_BAD_PIXELS(U64);
+        PM_BAD_PIXELS(S8);
+        PM_BAD_PIXELS(S16);
+        PM_BAD_PIXELS(S32);
+        PM_BAD_PIXELS(S64);
+        PM_BAD_PIXELS(F32);
+        PM_BAD_PIXELS(F64);
+    default:
+        psError( PS_ERR_BAD_PARAMETER_TYPE, true,
+                 PS_ERRORTEXT_pmMaskBadPixels_TYPE_UNSUPPORTED,
+                 inType);
+    }
     return true;
+    return false;
+}

trunk/psModules/src/detrend/pmMaskBadPixels.h

-              r5516
+              r6872
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
 /** @file  pmMaskBadPixels.h
+ *
 …
  *  @author Ross Harman, MHPCC
+ *
  *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-11-15 20:09:03 $
+ *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include "pslib.h"
 #include "pmAstrometry.h"
+#include "pmFPA.h"
 /** Mask values */
 …
     PM_MASK_BADCOL  = 0x0002,   ///< The pixel is a bad column.
     PM_MASK_SAT     = 0x0004,   ///< The pixel is saturated.
+    PM_MASK_FLAT    = 0x0008    ///< The pixel is non-positive in the flat-field.
+    PM_MASK_BAD     = 0x0008,   ///< The pixel is low
+    PM_MASK_FLAT    = 0x0010    ///< The pixel is non-positive in the flat-field.
 } pmMaskValue;
-// XXX: Use PS_MIN, PS_MAX
 /** Macro to find maximum of two numbers */
 #define MAX(A,B)((A)>=(B)?(A):(B))
 …
 bool pmMaskBadPixels(
     pmReadout *in,          ///< Readout containing input image data.
     const pmReadout *mask,   ///< Mask data to be added to readout mask data.
+    const psImage *mask,    ///< Mask data to be added to readout mask data.
     unsigned int maskVal,   ///< Mask value to determine what to add to input mask.
     float sat,              ///< Saturation limit to mask bad pixels.

trunk/psModules/src/detrend/pmMaskBadPixelsErrors.h

-              r5170
+              r6872
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
 /** @file  pmMaskBadPixelsErrors.h
+ *
 …
  *  @author Ross Harman, MHPCC
+ *
  *  @version $Revision: 1.1 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-09-28 20:43:52 $
+ *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
  *     $2  The error text (rest of the line in pmMaskBadPixelsErrors.h)
  *     $n  The order of the source line in pmMaskBadPixelsErrors.h (comments excluded)
+ *
+ *
  * DO NOT EDIT THE LINES BETWEEN //~Start and //~End!  ANY CHANGES WILL BE OVERWRITTEN.
  */

trunk/psModules/src/detrend/pmNonLinear.c

-              r5675
+              r6872
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
 /** @file  pmNonLinear.c
+ *
 …
  *  @author GLG, MHPCC
+ *
  *  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-12-05 20:49:40 $
+ *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include "pmNonLinear.h"
-#include "pmSubtractBias.h"
-// XXX: Remove, autoconf must be
-#define PS_WARN_PTR_NON_NULL(NAME) \
-if ((NAME) == NULL) { \
-    psLogMsg(__func__, PS_LOG_WARN, "WARNING: %s is NULL.", #NAME); \
-} \
 /******************************************************************************
 pmNonLinearityLookup(): This routine will take an pmReadout image as input
 and a 1-D polynomial.  For each pixel in the input image, the polynomial will
+be evaluated at that pixels value, and the image pixel will then be set
+to
+be evaluated at that pixels value, and the image pixel will then be set to
 that value.
 *****************************************************************************/
+ *****************************************************************************/
+pmReadout *pmNonLinearityPolynomial(
+    pmReadout *inputReadout,
+    const psPolynomial1D *input1DPoly)
+pmReadout *pmNonLinearityPolynomial(pmReadout *inputReadout,
+                                    const psPolynomial1D *input1DPoly)
+{
     PS_ASSERT_PTR_NON_NULL(inputReadout, NULL);
 …
     PS_ASSERT_IMAGE_TYPE(inputReadout->image, PS_TYPE_F32, NULL);
     PS_ASSERT_PTR_NON_NULL(input1DPoly, NULL);
-    PS_WARN_PTR_NON_NULL(inputReadout->parent);
-    if (inputReadout->parent != NULL) {
-        PS_WARN_PTR_NON_NULL(inputReadout->parent->concepts);
+    }
+    //
+    // Determine trimmed image from metadata.
+    //
+    psImage *trimmedImg = p_psDetermineTrimmedImage(inputReadout);
+    psS32 i;
+    psS32 j;
+    for (psS32 i=0;i<trimmedImg->numRows;i++) {
+        for (psS32 j=0;j<trimmedImg->numCols;j++) {
+            trimmedImg->data.F32[i][j] = psPolynomial1DEval(input1DPoly,
+                                         trimmedImg->data.F32[i][j]);
+    for (i=0;i<inputReadout->image->numRows;i++) {
+        for (j=0;j<inputReadout->image->numCols;j++) {
+            inputReadout->image->data.F32[i][j] = psPolynomial1DEval(input1DPoly, inputReadout->image->data.F32[i][j]);
+        }
+    }
 …
 inFluxe, and the corresponding value in outFlux.  The image pixel will then
 be set to the value from outFlux.
-XXX: Must assert that filename exists.  This should probably happen in
-the lookup files.
  *****************************************************************************/
+pmReadout *pmNonLinearityLookup(
+    pmReadout *inputReadout,
+    const char *filename
+)
+pmReadout *pmNonLinearityLookup(pmReadout *inputReadout,
+                                const psVector *inFlux,
+                                const psVector *outFlux)
+{
     PS_ASSERT_PTR_NON_NULL(inputReadout,NULL);
     PS_ASSERT_PTR_NON_NULL(inputReadout->image,NULL);
     PS_ASSERT_IMAGE_TYPE(inputReadout->image, PS_TYPE_F32, NULL);
+    PS_WARN_PTR_NON_NULL(inputReadout->parent);
+    if (inputReadout->parent != NULL) {
+        PS_WARN_PTR_NON_NULL(inputReadout->parent->concepts);
+    PS_ASSERT_PTR_NON_NULL(inFlux,NULL);
+    if (inFlux->n < 2) {
+        psError(PS_ERR_UNKNOWN,true, "pmNonLinearityLookup(): input vector less than 2 elements.  Returning inputReadout image.");
+        return(inputReadout);
+    }
+    //
+    // Determine trimmed image from metadata.
+    //
+    psImage *trimmedImg = p_psDetermineTrimmedImage(inputReadout);
+    PS_ASSERT_PTR_NON_NULL(outFlux,NULL);
+    psS32 tableSize = inFlux->n;
+    if (inFlux->n != outFlux->n) {
+        tableSize = PS_MIN(inFlux->n, outFlux->n);
+        psLogMsg(__func__, PS_LOG_WARN,
+                 "WARNING: pmNonLinear.c: pmNonLinearityLookup(): input vectors have different sizes (%d, %d)\n", inFlux->n, outFlux->n);
+    }
+    PS_ASSERT_VECTOR_TYPE(inFlux, PS_TYPE_F32, NULL);
+    PS_ASSERT_VECTOR_TYPE(outFlux, PS_TYPE_F32, NULL);
+    psLookupTable *tmpLT = psLookupTableAlloc(filename, "%f %f", 0);
+    psS32 numLines = psLookupTableRead(tmpLT);
+    psS32 i;
+    psS32 j;
+    psS32 binNum;
+    psScalar x;
     psS32 numPixels = 0;
+    if (numLines < 2) {
+        psLogMsg(__func__, PS_LOG_WARN,
+                 "WARNING: Lookup Table is too small.  Returning original pmReadout.\n");
+    } else {
+        for (psS32 i=0;i<trimmedImg->numRows;i++) {
+            for (psS32 j=0;j<trimmedImg->numCols;j++) {
+                psF64 tmpD = psLookupTableInterpolate(tmpLT, trimmedImg->data.F32[i][j], 1);
+                if (!isnan(tmpD)) {
+                    trimmedImg->data.F32[i][j] = tmpD;
+                } else {
+                    numPixels++;
+                }
+    psF32 slope;
+    x.type.type = PS_TYPE_F32;
+    for (i=0;i<inputReadout->image->numRows;i++) {
+        for (j=0;j<inputReadout->image->numCols;j++) {
+            x.data.F32 = inputReadout->image->data.F32[i][j];
+            binNum = p_psVectorBinDisect((psVector *)inFlux, &x);
+            if (binNum == -2) {
+                // We get here if x is below the table lookup range.
+                inputReadout->image->data.F32[i][j] = outFlux->data.F32[0];
+                numPixels++;
+            } else if (binNum == -1) {
+                // We get here if x is above the table lookup range.
+                inputReadout->image->data.F32[i][j] = outFlux->data.F32[tableSize-1];
+                numPixels++;
+            } else if (binNum < -2) {
+                // We get here if there was some other problem.
+                psError(PS_ERR_UNKNOWN,true, "pmNonLinearityLookup(): Could not perform p_psVectorBinDisect().  Returning inputReadout image.");
+                return(inputReadout);
+                numPixels++;
+            } else {
+                // Perform linear interpolation.
+                slope = (outFlux->data.F32[binNum+1] - outFlux->data.F32[binNum]) /
+                        (inFlux->data.F32[binNum+1]  - inFlux->data.F32[binNum]);
+                inputReadout->image->data.F32[i][j] = outFlux->data.F32[binNum] +
+                                                      ((x.data.F32 - inFlux->data.F32[binNum]) * slope);
+            }
+        }
-        if (numPixels > 0) {
-            psLogMsg(__func__, PS_LOG_WARN,
-                     "WARNING: pmNonLinear.c: pmNonLinearityLookup(): %d pixels outside table.", numPixels);
+        }
+    }
+    if (numPixels > 0) {
+        psLogMsg(__func__, PS_LOG_WARN,
+                 "WARNING: pmNonLinear.c: pmNonLinearityLookup(): %d pixels outside table.", numPixels);
+    }
     return(inputReadout);
+}

trunk/psModules/src/detrend/pmNonLinear.h

-              r5435
+              r6872
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
 /** @file  pmNonLinear.h
+ *
 …
  *  @author GLG, MHPCC
+ *
  *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-10-20 23:06:24 $
+ *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include "pslib.h"
 #include "pmAstrometry.h"
+#include "pmFPA.h"
+pmReadout *pmNonLinearityPolynomial(
+    pmReadout *in,
+    const psPolynomial1D *coeff
+);
+pmReadout *pmNonLinearityPolynomial(pmReadout *in,
+                                    const psPolynomial1D *coeff);
+pmReadout *pmNonLinearityLookup(
+    pmReadout *in,
+    const char *filename
+);
+pmReadout *pmNonLinearityLookup(pmReadout *in,
+                                const psVector *inFlux,
+                                const psVector *outFlux);
 #endif

trunk/psModules/src/imcombine/Makefile.am

-              r5170
+              r6872
 libpsmoduleimcombine_la_CPPFLAGS = $(SRCINC) $(PSMODULE_CFLAGS)
 libpsmoduleimcombine_la_LDFLAGS  = -release $(PACKAGE_VERSION)
 libpsmoduleimcombine_la_SOURCES  = pmImageCombine.c \
     pmReadoutCombine.c
+libpsmoduleimcombine_la_SOURCES  = pmImageCombine.c
+#    pmReadoutCombine.c
 psmoduleincludedir = $(includedir)
 psmoduleinclude_HEADERS = \
     pmImageCombine.h \
     pmReadoutCombine.h
+    pmImageCombine.h
+#    pmReadoutCombine.h

trunk/psModules/src/imcombine/pmReadoutCombine.h

-              r6205
+              r6872
  *  @author GLG, MHPCC
+ *
  *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2006-01-26 21:10:50 $
+ *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include "pslib.h"
 #include "psConstants.h"
 #include "pmAstrometry.h"
+#include "pmFPA.h"
 typedef struct
 …
     int nKeep;
+}
+psCombineParams;
+pmCombineParams;
+pmCombineParams *pmCombineParamsAlloc (psStatsOptions statsOptions);
 psImage *pmReadoutCombine(psImage *output,
+                          const psList *inputs,
+                          psCombineParams *params,
+                          const psArray *inputs,
                           const psVector *zero,
                           const psVector *scale,
+                          pmCombineParams *params,
                           bool applyZeroScale,
                           float gain,
 …
 pmFringePoint;
-/**
+ *
- * The input array fluxLevels consists of Ni vectors, one per mosaic image.
- * Each vector consists of Nj elements, each a measurement of the input
- * flat-field image flux levels. All of these vectors must be constructed with
- * the same number of elements, or the function will return an error. If a chip
- * is missing from a particular image, that element should be set to NaN. The
- * vector chipGains supplies initial guesses for the chip gains. If the vector
- * contains the values 0.0 or NaN for any of the elements, the gain is set to the
- * mean of the valid values. If the vector length does not match the number of
- * chips, an warning is raised, all chip gain guesses will be set to 1.0, and the
- * vector length modified to match the number of chips defined by the supplied
- * fluxLevels. The sourceFlux input vector must be allocated (not NULL), but the
- * routine will set the vector length to the number of source images regardless
- * of the initial state of the vector. All vectors used by this function must be
- * of type PS_DATA_F64.
+ *
- */
-bool pmFlatNormalization(
-    psVector *sourceFlux,
-    psVector *chipGains,
-    psArray *fluxLevels
-);
 #endif

trunk/psModules/src/imsubtract/Makefile.am

-              r5170
+              r6872
 libpsmoduleimsubtract_la_LDFLAGS  = -release $(PACKAGE_VERSION)
 libpsmoduleimsubtract_la_SOURCES  = pmImageSubtract.c \
     pmSubtractBias.c \
     pmSubtractSky.c
+    pmSubtractBias.c
+#    pmSubtractSky.c
 psmoduleincludedir = $(includedir)
 psmoduleinclude_HEADERS = \
     pmImageSubtract.h \
     pmSubtractBias.h \
     pmSubtractSky.h
+    pmSubtractBias.h
+#    pmSubtractSky.h

trunk/psModules/src/imsubtract/pmSubtractBias.h

-              r5587
+              r6872
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// XXX WARNING: I have completely replaced this file with an OLD VERSION (that works) instead of the
+// one that was being worked on.
+//////////////////////////////////////////////////////////////////////////////////////////////////////////////
 /** @file  pmSubtractBias.h
+ *
 …
  *  @author GLG, MHPCC
+ *
  *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-11-23 23:54:30 $
+ *  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include<math.h>
 #include "pslib.h"
+#include "pmAstrometry.h"
+#include "pmFPA.h"
 typedef enum {
     PM_OVERSCAN_NONE,                         ///< No overscan subtraction
+    PM_OVERSCAN_EDGE,                         ///< Subtract the statistic of pixels along the to-be-determined readout direction
     PM_OVERSCAN_ROWS,                         ///< Subtract rows
     PM_OVERSCAN_COLUMNS,                      ///< Subtract columns
 …
 typedef enum {
+    PM_FIT_NONE,                              ///< No fit
+    PM_FIT_POLYNOMIAL,                        ///< Fit polynomial
+    PM_FIT_SPLINE                             ///< Fit cubic splines
+    PM_FIT_NONE,                        ///< No fit
+    PM_FIT_POLY_ORD,                    ///< Fit ordinary polynomial
+    PM_FIT_POLY_CHEBY,                  ///< Fit Chebyshev polynomial
+    PM_FIT_SPLINE                       ///< Fit cubic splines
 } pmFit;
+pmReadout *pmSubtractBias(
+    pmReadout *in,                      ///< The input pmReadout image
+    void *fitSpec,                      ///< A polynomial or spline, defining the fit type.
+    pmFit fit,                          ///< PM_FIT_SPLINE, PM_FIT_POLYNOMIAL, or PM_FIT_NONE
+    bool overscan,
+    psStats *stat,                      ///< The statistic to be used in combining overscan data
+    int nBin,                           ///< The amount of binning to be done image pixels.
+    const pmReadout *bias,              ///< A possibly NULL bias pmReadout which is to be subtracted
+    const pmReadout *dark               ///< A possibly NULL bias pmReadout which is to be subtracted
+);
+typedef struct
+{
+    // Inputs
+    bool single;                // Reduce all overscan regions to a single value?
+    pmFit fitType;              // Type of fit to overscan
+    unsigned int order;         // Order of polynomial, or number of spline pieces
+    psStats *stat;              // Statistic to use when reducing the minor direction
+    // Outputs
+    psPolynomial1D *poly;       // Result of polynomial fit
+    psSpline1D *spline;         // Result of spline fit
+}
+pmOverscanOptions;
+/******************************************************************************
+DetermineTrimmedImageg() This private routine determines the region of the
+input pmReadout which will be operated on by the various detrend modules.  It
+does a metadata fetch on "CELL.TRIMSEC" for the parent cell of the pmReadout.
+XXX: Consider making a pmUtils.c file and put this routine there.
+ *****************************************************************************/
+psImage *p_psDetermineTrimmedImage(
+    pmReadout *in
+);
+pmOverscanOptions *pmOverscanOptionsAlloc(bool single, pmFit fitType, unsigned int order, psStats *stat);
+pmReadout *pmSubtractBias(pmReadout *in, pmOverscanOptions *overscanOpts,
+                          const pmReadout *bias, const pmReadout *dark);
+#if 0
+pmReadout *pmSubtractBias(pmReadout *in,                ///< The input pmReadout image
+                          void *fitSpec,                ///< A polynomial or spline, defining the fit type.
+                          const psList *overscans,      ///< A psList of overscan images
+                          pmOverscanAxis overScanAxis,  ///< Defines how overscans are applied
+                          psStats *stat,                ///< The statistic to be used in combining overscan data
+                          int nBin,                     ///< The amount of binning to be done image pixels.
+                          pmFit fit,                    ///< PM_FIT_SPLINE, PM_FIT_POLYNOMIAL, or PM_FIT_NONE
+                          const pmReadout *bias);       ///< A possibly NULL bias pmReadout which is to be subtracted
+#endif
 #endif

trunk/psModules/src/imsubtract/pmSubtractSky.h

-              r5170
+              r6872
  *  @author GLG, MHPCC
+ *
  *  @version $Revision: 1.1 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-09-28 20:43:52 $
+ *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include<math.h>
 #include "pslib.h"
 #include "pmAstrometry.h"
+#include "pmFPA.h"
 // XXX: this is pmFit in pmSubtractBias.c, named psFit here.

trunk/psModules/src/objects/Makefile.am

-              r5844
+              r6872
 noinst_LTLIBRARIES = libpsmoduleobjects.la
 libpsmoduleobjects_la_CPPFLAGS = $(SRCINC) $(PSMODULE_CFLAGS)
+libpsmoduleobjects_la_CPPFLAGS = $(SRCINC) $(PSMODULE_CFLAGS) -I../pslib/
 libpsmoduleobjects_la_LDFLAGS  = -release $(PACKAGE_VERSION)
 libpsmoduleobjects_la_SOURCES  = \
+    pmObjects.c \
+    pmPSF.c \
+    pmPSFtry.c \
+    pmModelGroup.c \
+    psEllipse.c
+     pmPeaks.c \
+     pmMoments.c \
+     pmModel.c \
+     pmModelGroup.c \
+     pmSource.c \
+     pmSourceSky.c \
+     pmSourceContour.c \
+     pmSourceFitModel.c \
+     pmSourcePhotometry.c \
+     pmSourceIO.c \
+     pmSourceIO_CMF.c \
+     pmSourceIO_CMP.c \
+     pmSourceIO_OBJ.c \
+     pmSourceIO_SX.c \
+     pmSourceIO_RAW.c \
+     pmPSF.c \
+     pmPSF_IO.c \
+     pmPSFtry.c \
+     pmGrowthCurve.c
 EXTRA_DIST = \
 …
 psmoduleincludedir = $(includedir)
 psmoduleinclude_HEADERS = \
+    pmObjects.h \
+    pmPSF.h \
+    pmPSFtry.h \
+    pmModelGroup.h \
+    psEllipse.h
+     pmPeaks.h \
+     pmMoments.h \
+     pmModel.h \
+     pmModelGroup.h \
+     pmSource.h \
+     pmSourceSky.h \
+     pmSourceContour.h \
+     pmSourceFitModel.h \
+     pmSourcePhotometry.h \
+     pmSourceIO.h \
+     pmPSF.h \
+     pmPSF_IO.h \
+     pmPSFtry.h \
+     pmGrowthCurve.h

trunk/psModules/src/objects/models/pmModel_PGAUSS.c

-              r6511
+              r6872
     if (deriv != NULL) {
         // note difference from a pure gaussian: q = PAR[1]*r
+        psF32 *dPAR = deriv->data.F32;
         psF32 q = PAR[1]*r*r*t;
         deriv->data.F32[0] = +1.0;
         deriv->data.F32[1] = +r;
         deriv->data.F32[2] = q*(2.0*px*PAR[4] + PAR[6]*Y);
         deriv->data.F32[3] = q*(2.0*py*PAR[5] + PAR[6]*X);
         deriv->data.F32[4] = -2.0*q*px*X;
         deriv->data.F32[5] = -2.0*q*py*Y;
         deriv->data.F32[6] = -q*X*Y;
+        dPAR[0] = +1.0;
+        dPAR[1] = +r;
+        dPAR[2] = q*(2.0*px*PAR[4] + PAR[6]*Y);
+        dPAR[3] = q*(2.0*py*PAR[5] + PAR[6]*X);
+        dPAR[4] = -2.0*q*px*X;
+        dPAR[5] = -2.0*q*py*Y;
+        dPAR[6] = -q*X*Y;
+    }
     return(f);
 …
     beta_lim[0][0].data.F32[0] = 1000;
     beta_lim[0][0].data.F32[1] = 10000;
+    beta_lim[0][0].data.F32[1] = 3e6;
     beta_lim[0][0].data.F32[2] = 5;
     beta_lim[0][0].data.F32[3] = 5;
 …
     params_max[0][0].data.F32[0] = 1e5;
     params_max[0][0].data.F32[1] = 1e6;
+    params_max[0][0].data.F32[1] = 1e8;
     params_max[0][0].data.F32[2] = 1e4;  // this should be set by image dimensions!
     params_max[0][0].data.F32[3] = 1e4;  // this should be set by image dimensions!
 …
     params[5] = 1.2 / moments->Sy;
     params[6] = 0.0;
     return(true);
+}

trunk/psModules/src/objects/models/pmModel_QGAUSS.c

-              r6511
+              r6872
     psF32 py = PAR[5]*Y;
     psF32 z  = 0.5*PS_SQR(px) + 0.5*PS_SQR(py) + PAR[6]*X*Y;
+    psF32 r  = 1.0 / (1 + PAR[7]*z + pow(z, 2.25));
+    psF32 f  = PAR[1]*r + PAR[0];
+    psF32 zp = pow(z,1.25);
+    psF32 r  = 1.0 / (1 + PAR[7]*z + z*zp);
+    // test: psF32 r  = 1.0 / (1 + PAR[7]*z + PS_SQR(z));
+    psF32 r1 = PAR[1]*r;
+    psF32 f  = r1 + PAR[0];
     if (deriv != NULL) {
+        psF32 *dPAR = deriv->data.F32;
         // note difference from a pure gaussian: q = params->data.F32[1]*r
+        psF32 t = PAR[1]*r*r;
+        psF32 q = t*(PAR[7] + 2.25*pow(z, 1.25));
+        deriv->data.F32[0] = +1.0;
+        deriv->data.F32[1] = +r;
+        deriv->data.F32[2] = q*(2.0*px*PAR[4] + PAR[6]*Y);
+        deriv->data.F32[3] = q*(2.0*py*PAR[5] + PAR[6]*X);
+        deriv->data.F32[4] = -2.0*q*px*X;
+        deriv->data.F32[5] = -2.0*q*py*Y;
+        deriv->data.F32[6] = -q*X*Y;
+        deriv->data.F32[7] = -t*z;
+        psF32 t = r1*r;
+        psF32 q = t*(PAR[7] + 2.25*zp);
+        // test: psF32 q = t*(PAR[7] + 2*z);
+        dPAR[0] = +1.0;
+        dPAR[1] = +r;
+        dPAR[2] = q*(2.0*px*PAR[4] + PAR[6]*Y);
+        dPAR[3] = q*(2.0*py*PAR[5] + PAR[6]*X);
+        dPAR[4] = -2.0*q*px*X;
+        dPAR[5] = -2.0*q*py*Y;
+        dPAR[6] = -q*X*Y;
+        dPAR[7] = -t*z;
+    }
     return(f);
 …
     beta_lim[0][0].data.F32[0] = 1000;
     beta_lim[0][0].data.F32[1] = 10000;
+    beta_lim[0][0].data.F32[1] = 3e6;
     beta_lim[0][0].data.F32[2] = 5;
     beta_lim[0][0].data.F32[3] = 5;
 …
     params_max[0][0].data.F32[0] = 1e5;
     params_max[0][0].data.F32[1] = 1e6;
+    params_max[0][0].data.F32[1] = 1e8;
     params_max[0][0].data.F32[2] = 1e4;  // this should be set by image dimensions!
     params_max[0][0].data.F32[3] = 1e4;  // this should be set by image dimensions!
 …
     params[6] = 0.0;
     params[7] = 1.0;
     return(true);
+}
 …
     // the area needs to be multiplied by the integral of f(z)
     norm = 0.0;
     for (z = 0.005; z < 50; z += 0.01) {
+    for (z = 0.05; z < 50; z += 0.1) {
         f = 1.0 / (1 + PAR[7]*z + pow(z, 2.25));
+        // test: f = 1.0 / (1 + PAR[7]*z + PS_SQR(z));
         norm += f;
+    }
     norm *= 0.01;
+    norm *= 0.1;
     psF64 Flux = PAR[1] * Area * norm;
 …
     // we can do this much better with intelligent choices here
     for (z = 0.0; z < 20.0; z += dz) {
+    for (z = 0.0; z < 30.0; z += dz) {
         f = 1.0 / (1 + PAR[7]*z + pow(z, 2.25));
+        // test: f = 1.0 / (1 + PAR[7]*z + PS_SQR(z));
         if (f < limit)
             break;
 …
     status &= ((dPAR[1]/PAR[1]) < 0.5);
     if (status)
         return true;
     return false;
+}
+    if (!status)
+        return false;
+    return true;
+}

trunk/psModules/src/objects/models/pmModel_RGAUSS.c

-              r5257
+              r6872
     psVector *params = model->params;
     EllipseAxes axes;
     EllipseShape shape;
     EllipseMoments moments;
+    psEllipseAxes axes;
+    psEllipseShape shape;
+    psEllipseMoments moments;
     moments.x2 = PS_SQR(source->moments->Sx);
 …
     moments.xy = source->moments->Sxy;
     axes = EllipseMomentsToAxes(moments);
     shape = EllipseAxesToShape(axes);
+    axes = psEllipseMomentsToAxes(moments);
+    shape = psEllipseAxesToShape(axes);
     params->data.F32[0] = source->moments->Sky;

trunk/psModules/src/objects/models/pmModel_SGAUSS.c

-              r6511
+              r6872
 # define SQ(A)((A)*(A))
 psF64 psImageEllipseContour (EllipseAxes axes, double xc, double yc, psImage *image);
+psF64 psImageEllipseContour (psEllipseAxes axes, double xc, double yc, psImage *image);
 psF64 p_psImageGetElementF64(psImage *a, int i, int j);
 …
 // measure the flux for the elliptical contour
 psF64 psImageEllipseContour (EllipseAxes axes, double xc, double yc, psImage *image)
+psF64 psImageEllipseContour (psEllipseAxes axes, double xc, double yc, psImage *image)
+{
 …
     psF32     *params   = model->params->data.F32;
     EllipseAxes axes;
     EllipseShape shape;
     EllipseMoments moments;
+    psEllipseAxes axes;
+    psEllipseShape shape;
+    psEllipseMoments moments;
     moments.x2 = PS_SQR(sMoments->Sx);
 …
     // solve the math to go from Moments To Shape
     axes = EllipseMomentsToAxes(moments);
     shape = EllipseAxesToShape(axes);
+    axes = psEllipseMomentsToAxes(moments);
+    shape = psEllipseAxesToShape(axes);
     params[0] = sMoments->Sky;
 …
     float f1, f2;
     EllipseAxes axes;
     EllipseShape shape;
     EllipseMoments moments;
+    psEllipseAxes axes;
+    psEllipseShape shape;
+    psEllipseMoments moments;
     moments.x2 = PS_SQR(sMoments->Sx);
 …
     // solve the math to go from Moments To Shape
     axes = EllipseMomentsToAxes(moments);
     shape = EllipseAxesToShape(axes);
+    axes = psEllipseMomentsToAxes(moments);
+    shape = psEllipseAxesToShape(axes);
     params[0] = sMoments->Sky;
 …
 psF64 pmModelRadius_SGAUSS  (const psVector *params, psF64 flux)
+{
     psF64 r, z, pr, f;
+    psF64 r, z = 0.0, pr, f;
     psF32 *PAR = params->data.F32;
     EllipseAxes axes;
     EllipseShape shape;
+    psEllipseAxes axes;
+    psEllipseShape shape;
     if (flux <= 0)
 …
     shape.sxy = PAR[6];
     axes = EllipseShapeToAxes (shape);
+    axes = psEllipseShapeToAxes (shape);
     psF64 dr = 1.0 / axes.major;
     psF64 limit = flux / PAR[1];
 …
     psF32 dP;
     bool  status;
     EllipseAxes axes;
     EllipseShape shape;
+    psEllipseAxes axes;
+    psEllipseShape shape;
     psF32 *PAR  = model->params->data.F32;
 …
     shape.sxy = PAR[6];
     axes = EllipseShapeToAxes (shape);
+    axes = psEllipseShapeToAxes (shape);
     dP = 0;

trunk/psModules/src/objects/models/pmModel_ZGAUSS.c

-              r5257
+              r6872
     psF32 *PAR = params->data.F32;
     EllipseAxes axes;
     EllipseShape shape;
+    psEllipseAxes axes;
+    psEllipseShape shape;
     if (flux <= 0)
 …
     shape.sxy = PAR[6];
     axes = EllipseShapeToAxes (shape);
+    axes = psEllipseShapeToAxes (shape);
     psF64 dr = 1.0 / axes.major;
     psF64 limit = flux / PAR[1];
 …
     psVector *params = model->params;
     EllipseAxes axes;
     EllipseShape shape;
     EllipseMoments moments;
+    psEllipseAxes axes;
+    psEllipseShape shape;
+    psEllipseMoments moments;
     moments.x2 = PS_SQR(source->moments->Sx);
 …
     moments.xy = source->moments->Sxy;
     axes = EllipseMomentsToAxes(moments);
     shape = EllipseAxesToShape(axes);
+    axes = psEllipseMomentsToAxes(moments);
+    shape = psEllipseAxesToShape(axes);
     params->data.F32[0] = source->moments->Sky;

trunk/psModules/src/objects/pmModelGroup.c

-              r5844
+              r6872
+# include "pmModelGroup.h"
+/** @file  pmModelGroup.c
+ *
+ *  Functions to define and manipulate object model attributes
+ *
+ *  @author GLG, MHPCC
+ *  @author EAM, IfA
+ *
+ *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
+ *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
+ *
+ */
+#include <stdio.h>
+#include <math.h>
+#include <string.h>
+#include "pslib.h"
+#include "psEllipse.h"
+#include "pmHDU.h"
+#include "pmFPA.h"
+#include "pmPeaks.h"
+#include "pmMoments.h"
+#include "pmGrowthCurve.h"
+#include "pmModel.h"
+#include "pmPSF.h"
+#include "pmSource.h"
+#include "pmModelGroup.h"
+// XXX shouldn't these be defined for us in pslib.h ???
 double hypot(double x, double y);
 double sqrt (double x);
-#include "psEllipse.h"
 #include "models/pmModel_GAUSS.c"
 #include "models/pmModel_PGAUSS.c"
 …
     if (modelGroup == NULL)
         return;
-    psFree (modelGroup);
     return;
+}
 …
+    }
     Nmodels = Nnew;
+    return;
+}
+void pmModelGroupCleanup (void)
+{
+    psFree (models);
     return;
+}
 …
     return (models[type].name);
+}
+/******************************************************************************
+    pmSourceModelGuess(source, model): This function allocates a new
+    pmModel structure based on the given modelType specified in the argument list.
+    The corresponding pmModelGuess function is returned, and used to
+    supply the values of the params array in the pmModel structure.
+    XXX: Many parameters are based on the src->moments structure, which is in
+    image, not subImage coords.  Therefore, the calls to the model evaluation
+    functions will be in image, not subImage coords.  Remember this.
+*****************************************************************************/
+pmModel *pmSourceModelGuess(pmSource *source,
+                            pmModelType modelType)
+{
+    psTrace(__func__, 3, "---- %s() begin ----\n", __func__);
+    PS_ASSERT_PTR_NON_NULL(source->moments, false);
+    PS_ASSERT_PTR_NON_NULL(source->peak, false);
+    pmModel *model = pmModelAlloc(modelType);
+    pmModelGuessFunc modelGuessFunc = pmModelGuessFunc_GetFunction(modelType);
+    modelGuessFunc(model, source);
+    psTrace(__func__, 3, "---- %s() end ----\n", __func__);
+    return(model);
+}

trunk/psModules/src/objects/pmModelGroup.h

-              r5844
+              r6872
  *  @author EAM, IfA
+ *
  *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-12-24 01:24:32 $
+ *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
+ *
  */
+#include "pmObjects.h"
+#include "pmPSF.h"
+/**
+ *
+ * This function returns the number of parameters used by the listed function.
+ *
+# ifndef PM_MODEL_GROUP_H
+# define PM_MODEL_GROUP_H
+//  This function is the model chi-square minimization function for this model.
+typedef psMinimizeLMChi2Func pmModelFunc;
+// This function returns the integrated flux for the given model parameters.
+typedef psF64 (*pmModelFlux)(const psVector *params);
+// This function returns the radius at which the given model and parameters
+// achieves the given flux.
+typedef psF64 (*pmModelRadius)(const psVector *params, double flux);
+/*  This function sets the model parameter limits vectors for the given model
+ */
+typedef bool (*pmModelLimits)(psVector **beta_lim, psVector **params_min, psVector **params_max);
+/*  This function provides the model guess parameters based on the details of
+ *   the given source.
+ */
+typedef bool (*pmModelGuessFunc)(pmModel *model, pmSource *source);
+/*  This function constructs the PSF model for the given source based on the
+ *  supplied psf and the EXT model for the object.
+ */
+typedef bool (*pmModelFromPSFFunc)(pmModel *modelPSF, pmModel *modelEXT, pmPSF *psf);
+/*  This function returns the success / failure status of the given model fit
+ */
+typedef bool (*pmModelFitStatusFunc)(pmModel *model);
+/* Every model instance belongs to a class of models, defined by the value of
+ * the pmModelType type entry. Various functions need access to information about
+ * each of the models. Some of this information varies from model to model, and
+ * may depend on the current parameter values or other data quantities. In order
+ * to keep the code from requiring the information about each model to be coded
+ * into the low-level fitting routines, we define a collection of functions which
+ * allow us to abstract this type of model-dependent information. These generic
+ * functions take the model type and return the corresponding function pointer
+ * for the specified model. Each model is defined by creating this collection of
+ * specific functions, and placing them in a single file for each model. We
+ * define the following structure to carry the collection of information about
+ * the models.
+ */
+typedef struct
+{
+    char *name;
+    int nParams;
+    pmModelFunc          modelFunc;
+    pmModelFlux          modelFlux;
+    pmModelRadius        modelRadius;
+    pmModelLimits        modelLimits;
+    pmModelGuessFunc     modelGuessFunc;
+    pmModelFromPSFFunc   modelFromPSFFunc;
+    pmModelFitStatusFunc modelFitStatusFunc;
+}
+pmModelGroup;
+// allocate a pmModelGroup to hold nModels entries
+pmModelGroup *pmModelGroupAlloc (int nModels);
+// initialize the internal (static) model group with the default models
+void pmModelGroupInit (void);
+// free the internal (static) model group
+void pmModelGroupCleanup (void);
+// add a new model to the internal (static) model group
+void pmModelGroupAdd (pmModelGroup *model);
+/* This function returns the number of parameters used by the listed function.
  */
 int pmModelParameterCount(
 …
+/**
+ *
+ * This function returns the user-space model names for the specified model type.
+ *
+/* This function returns the user-space model names for the specified model type.
  */
 char *pmModelGetType(
 …
 );
-#ifndef PM_MODEL_GROUP_H
-#define PM_MODEL_GROUP_H
-/**
+ *
- *  This function is the model chi-square minimization function for this model.
+ *
- */
-typedef psMinimizeLMChi2Func pmModelFunc;
-/**
+ *
- * This function returns the integrated flux for the given model parameters.
- */
-typedef psF64 (*pmModelFlux)(const psVector *params);
-/**
+ *
- *  This function returns the radius at which the given model and parameters
- *  achieves the given flux.
+ *
- */
-typedef psF64 (*pmModelRadius)(const psVector *params, double flux);
-/**
+ *
- *  This function sets the model parameter limits vectors for the given model
+ *
- */
-typedef bool (*pmModelLimits)(psVector **beta_lim, psVector **params_min, psVector **params_max);
-/**
+ *
- *  This function provides the model guess parameters based on the details of
- *   the given source.
+ *
- */
-typedef bool (*pmModelGuessFunc)(pmModel *model, pmSource *source);
-/**
+ *
- *  This function constructs the PSF model for the given source based on the
- *  supplied psf and the FLT model for the object.
+ *
- */
-typedef bool (*pmModelFromPSFFunc)(pmModel *modelPSF, pmModel *modelFLT, pmPSF *psf);
-/**
+ *
- *  This function returns the success / failure status of the given model fit
+ *
- */
-typedef bool (*pmModelFitStatusFunc)(pmModel *model);
 /**
+ *
 …
+ *
  */
 /**
 …
+/** pmSourceModelGuess()
+ *
+ * Convert available data to an initial guess for the given model. This
+ * function allocates a pmModel entry for the pmSource structure based on the
+ * provided model selection. The method of defining the model parameter guesses
+ * are specified for each model below. The guess values are placed in the model
+ * parameters. The function returns TRUE on success or FALSE on failure.
+ *
+ */
+pmModel *pmSourceModelGuess(
+    pmSource *source,   ///< The input pmSource
+    pmModelType model   ///< The type of model to be created.
+);
+/**
+ *
+ * Every model instance belongs to a class of models, defined by the value of
+ * the pmModelType type entry. Various functions need access to information about
+ * each of the models. Some of this information varies from model to model, and
+ * may depend on the current parameter values or other data quantities. In order
+ * to keep the code from requiring the information about each model to be coded
+ * into the low-level fitting routines, we define a collection of functions which
+ * allow us to abstract this type of model-dependent information. These generic
+ * functions take the model type and return the corresponding function pointer
+ * for the specified model. Each model is defined by creating this collection of
+ * specific functions, and placing them in a single file for each model. We
+ * define the following structure to carry the collection of information about
+ * the models.
+ *
+ */
+typedef struct
+{
+    char *name;
+    int nParams;
+    pmModelFunc          modelFunc;
+    pmModelFlux          modelFlux;
+    pmModelRadius        modelRadius;
+    pmModelLimits        modelLimits;
+    pmModelGuessFunc     modelGuessFunc;
+    pmModelFromPSFFunc   modelFromPSFFunc;
+    pmModelFitStatusFunc modelFitStatusFunc;
+}
+pmModelGroup;
+// allocate a pmModelGroup to hold nModels entries
+pmModelGroup *pmModelGroupAlloc (int nModels);
+// initialize the internal (static) model group with the default models
+void pmModelGroupInit (void);
+// add a new model to the internal (static) model group
+void pmModelGroupAdd (pmModelGroup *model);
+# endif
+# endif /* PM_MODEL_GROUP_H */

trunk/psModules/src/objects/pmObjects.h

-              r5844
+              r6872
  *  @author GLG, MHPCC
+ *
  *  @version $Revision: 1.5 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-12-24 01:24:32 $
+ *  @version $Revision: 1.6 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 #include <math.h>
 #include "pslib.h"
-#include "pmAstrometry.h"
-/**
- * In the object analysis process, we will use specific mask values to mark the
- * image pixels. The following structure defines the relevant mask values.
+ *
- * XXX: This is probably a bad solution: we will want to set mask values
- * outside of the PSPHOT code.  Perhaps we can set up a registered set of mask
- * values with specific meanings that other functions can add to or define?
- */
-typedef enum {
-    PSPHOT_MASK_CLEAR     = 0x00,
-    PSPHOT_MASK_INVALID   = 0x01,
-    PSPHOT_MASK_SATURATED = 0x02,
-    PSPHOT_MASK_MARKED    = 0x08,
-} psphotMaskValues;
-/** pmPeakType
+ *
- *  A peak pixel may have several features which may be determined when the
- *  peak is found or measured. These are specified by the pmPeakType enum.
- *  PM_PEAK_LONE represents a single pixel which is higher than its 8 immediate
- *  neighbors.  The PM_PEAK_EDGE represents a peak pixel which touching the image
- *  edge. The PM_PEAK_FLAT represents a peak pixel which has more than a specific
- *  number of neighbors at the same value, within some tolarence:
+ *
- */
-typedef enum {
-    PM_PEAK_LONE,                       ///< Isolated peak.
-    PM_PEAK_EDGE,                       ///< Peak on edge.
-    PM_PEAK_FLAT,                       ///< Peak has equal-value neighbors.
-    PM_PEAK_UNDEF                       ///< Undefined.
-} pmPeakType;
-/** pmPeak data structure
+ *
- *  A source has the capacity for several types of measurements. The
- *  simplest measurement of a source is the location and flux of the peak pixel
- *  associated with the source:
+ *
- */
-typedef struct
+{
-    int x;                              ///< X-coordinate of peak pixel.
-    int y;                              ///< Y-coordinate of peak pixel.
-    float counts;                       ///< Value of peak pixel (above sky?).
-    pmPeakType class;                   ///< Description of peak.
+}
-pmPeak;
-/** pmMoments data structure
+ *
- * One of the simplest measurements which can be made quickly for an object
- * are the object moments. We specify a structure to carry the moment information
- * for a specific source:
+ *
- */
-typedef struct
+{
-    float x;    ///< X-coord of centroid.
-    float y;    ///< Y-coord of centroid.
-    float Sx;    ///< x-second moment.
-    float Sy;    ///< y-second moment.
-    float Sxy;    ///< xy cross moment.
-    float Sum;    ///< Pixel sum above sky (background).
-    float Peak;    ///< Peak counts above sky.
-    float Sky;    ///< Sky level (background).
-    float SN;    ///< approx signal-to-noise
-    int nPixels;   ///< Number of pixels used.
+}
-pmMoments;
-/** pmPSFClump data structure
+ *
- * A collection of object moment measurements can be used to determine
- * approximate object classes. The key to this analysis is the location and
- * statistics (in the second-moment plane,
+ *
- */
-typedef struct
+{
-    float X;
-    float dX;
-    float Y;
-    float dY;
+}
-pmPSFClump;
-// type of model carried by the pmModel structure
-typedef int pmModelType;
-typedef enum {
-    PM_MODEL_UNTRIED,               ///< model fit not yet attempted
-    PM_MODEL_SUCCESS,               ///< model fit succeeded
-    PM_MODEL_NONCONVERGE,           ///< model fit did not converge
-    PM_MODEL_OFFIMAGE,              ///< model fit drove out of range
-    PM_MODEL_BADARGS                ///< model fit called with invalid args
-} pmModelStatus;
-/** pmModel data structure
+ *
- * Every source may have two types of models: a PSF model and a FLT (floating)
- * model. The PSF model represents the best fit of the image PSF to the specific
- * object. In this case, the PSF-dependent parameters are specified for the
- * object by the PSF, not by the fit. The FLT model represents the best fit of
- * the given model to the object, with all parameters floating in the fit.
+ *
- */
-typedef struct
+{
-    pmModelType type;   ///< Model to be used.
-    psVector *params;   ///< Paramater values.
-    psVector *dparams;   ///< Parameter errors.
-    float chisq;   ///< Fit chi-squared.
-    int nDOF;    ///< number of degrees of freedom
-    int nIter;    ///< number of iterations to reach min
-    int status;         ///< fit status
-    float radius;   ///< fit radius actually used
+}
-pmModel;
-/** pmSourceType enumeration
+ *
- * A given source may be identified as most-likely to be one of several source
- * types. The pmSource entry pmSourceType defines the current best-guess for this
- * source.
+ *
- * XXX: The values given below are currently illustrative and will require
- * some modification as the source classification code is developed. (TBD)
+ *
- */
-typedef enum {
-    PM_SOURCE_DEFECT,                   ///< a cosmic-ray
-    PM_SOURCE_SATURATED,                ///< random saturated pixels
-    PM_SOURCE_SATSTAR,                  ///< a saturated star
-    PM_SOURCE_PSFSTAR,                  ///< a PSF star
-    PM_SOURCE_GOODSTAR,                 ///< a good-quality star
-    PM_SOURCE_POOR_FIT_PSF,             ///< poor quality PSF fit
-    PM_SOURCE_FAIL_FIT_PSF,             ///< failed to get a good PSF fit
-    PM_SOURCE_FAINTSTAR,                ///< below S/N cutoff
-    PM_SOURCE_GALAXY,                   ///< an extended object (galaxy)
-    PM_SOURCE_FAINT_GALAXY,             ///< a galaxy below S/N cutoff
-    PM_SOURCE_DROP_GALAXY,              ///< ?
-    PM_SOURCE_FAIL_FIT_GAL,             ///< failed on the galaxy fit
-    PM_SOURCE_POOR_FIT_GAL,             ///< poor quality galaxy fit
-    PM_SOURCE_OTHER,                    ///< unidentified
-} pmSourceType;
-/** pmSource data structure
+ *
- *  This source has the capacity for several types of measurements. The
- *  simplest measurement of a source is the location and flux of the peak pixel
- *  associated with the source:
+ *
- */
-typedef struct
+{
-    pmPeak *peak;   ///< Description of peak pixel.
-    psImage *pixels;   ///< Rectangular region including object pixels.
-    psImage *weight;   ///< Image variance.
-    psImage *mask;   ///< Mask which marks pixels associated with objects.
-    pmMoments *moments;   ///< Basic moments measure for the object.
-    pmModel *modelPSF;   ///< PSF Model fit (parameters and type)
-    pmModel *modelFLT;   ///< FLT (floating) Model fit (parameters and type).
-    pmSourceType type;   ///< Best identification of object.
-    float apMag;
-    float fitMag;
+}
-pmSource;
-/** pmPeakAlloc()
+ *
- *  @return pmPeak*    newly allocated pmPeak with all internal pointers set to NULL
- */
-pmPeak *pmPeakAlloc(
-    int x,    ///< Row-coordinate in image space
-    int y,    ///< Col-coordinate in image space
-    float counts,   ///< The value of the peak pixel
-    pmPeakType class   ///< The type of peak pixel
-);
-/** pmMomentsAlloc()
+ *
- */
-pmMoments *pmMomentsAlloc();
-/** pmModelAlloc()
+ *
- */
-pmModel *pmModelAlloc(pmModelType type);
-/** pmSourceAlloc()
+ *
- */
-pmSource  *pmSourceAlloc();
-/** pmFindVectorPeaks()
+ *
- * Find all local peaks in the given vector above the given threshold. A peak
- * is defined as any element with a value greater than its two neighbors and with
- * a value above the threshold. Two types of special cases must be addressed.
- * Equal value elements: If an element has the same value as the following
- * element, it is not considered a peak. If an element has the same value as the
- * preceding element (but not the following), then it is considered a peak. Note
- * that this rule (arbitrarily) identifies flat regions by their trailing edge.
- * Edge cases: At start of the vector, the element must be higher than its
- * neighbor. At the end of the vector, the element must be higher or equal to its
- * neighbor. These two rules again places the peak associated with a flat region
- * which touches the image edge at the image edge. The result of this function is
- * a vector containing the coordinates (element number) of the detected peaks
- * (type psU32).
+ *
- */
-psVector *pmFindVectorPeaks(
-    const psVector *vector,  ///< The input vector (float)
-    float threshold   ///< Threshold above which to find a peak
-);
-/** pmFindImagePeaks()
+ *
- * Find all local peaks in the given image above the given threshold. This
- * function should find all row peaks using pmFindVectorPeaks, then test each row
- * peak and exclude peaks which are not local peaks. A peak is a local peak if it
- * has a higher value than all 8 neighbors. If the peak has the same value as its
- * +y neighbor or +x neighbor, it is NOT a local peak. If any other neighbors
- * have an equal value, the peak is considered a valid peak. Note two points:
- * first, the +x neighbor condition is already enforced by pmFindVectorPeaks.
- * Second, these rules have the effect of making flat-topped regions have single
- * peaks at the (+x,+y) corner. When selecting the peaks, their type must also be
- * set. The result of this function is an array of pmPeak entries.
+ *
- */
-psArray *pmFindImagePeaks(
-    const psImage *image,  ///< The input image where peaks will be found (float)
-    float threshold   ///< Threshold above which to find a peak
-);
-/** pmCullPeaks()
+ *
- * Eliminate peaks from the psList that have a peak value above the given
- * maximum, or fall outside the valid region.
+ *
- */
-psList *pmCullPeaks(
-    psList *peaks,   ///< The psList of peaks to be culled
-    float maxValue,   ///< Cull peaks above this value
-    const psRegion valid                ///< Cull peaks otside this psRegion
-);
-/** pmPeaksSubset()
+ *
- * Create a new peaks array, removing certain types of peaks from the input
- * array of peaks based on the given criteria. Peaks should be eliminated if they
- * have a peak value above the given maximum value limit or if the fall outside
- * the valid region.  The result of the function is a new array with a reduced
- * number of peaks.
+ *
- */
-psArray *pmPeaksSubset(
-    psArray *peaks,                     ///< Add comment.
-    float maxvalue,                     ///< Add comment.
-    const psRegion valid                ///< Add comment.
-);
-/** pmSourceDefinePixels()
+ *
- * Define psImage subarrays for the source located at coordinates x,y on the
- * image set defined by readout. The pixels defined by this operation consist of
- * a square window (of full width 2Radius+1) centered on the pixel which contains
- * the given coordinate, in the frame of the readout. The window is defined to
- * have limits which are valid within the boundary of the readout image, thus if
- * the radius would fall outside the image pixels, the subimage is truncated to
- * only consist of valid pixels. If readout->mask or readout->weight are not
- * NULL, matching subimages are defined for those images as well. This function
- * fails if no valid pixels can be defined (x or y less than Radius, for
- * example). This function should be used to define a region of interest around a
- * source, including both source and sky pixels.
+ *
- * XXX: must code this.
+ *
- */
-// XXX: Uncommenting the pmReadout causes compile errors.
-bool pmSourceDefinePixels(
-    pmSource *mySource,                 ///< Add comment.
-    pmReadout *readout,                 ///< Add comment.
-    psF32 x,                            ///< Add comment.
-    psF32 y,                            ///< Add comment.
-    psF32 Radius                        ///< Add comment.
-);
-/** pmSourceLocalSky()
+ *
- * Measure the local sky in the vicinity of the given source. The Radius
- * defines the square aperture in which the moments will be measured. This
- * function assumes the source pixels have been defined, and that the value of
- * Radius here is smaller than the value of Radius used to define the pixels. The
- * annular region not contained within the radius defined here is used to measure
- * the local background in the vicinity of the source. The local background
- * measurement uses the specified statistic passed in via the statsOptions entry.
- * This function allocates the pmMoments structure. The resulting sky is used to
- * set the value of the pmMoments.sky element of the provided pmSource structure.
+ *
- */
-bool pmSourceLocalSky(
-    pmSource *source,   ///< The input image (float)
-    psStatsOptions statsOptions, ///< The statistic used in calculating the background sky
-    float Radius   ///< The inner radius of the square annulus to exclude
-);
-/** pmSourceMoments()
+ *
- * Measure source moments for the given source, using the value of
- * source.moments.sky provided as the local background value and the peak
- * coordinates as the initial source location. The resulting moment values are
- * applied to the source.moments entry, and the source is returned. The moments
- * are measured within the given circular radius of the source.peak coordinates.
- * The return value indicates the success (TRUE) of the operation.
+ *
- */
-bool pmSourceMoments(
-    pmSource *source,   ///< The input pmSource for which moments will be computed
-    float radius   ///< Use a circle of pixels around the peak
-);
-/** pmSourcePSFClump()
+ *
- * We use the source moments to make an initial, approximate source
- * classification, and as part of the information needed to build a PSF model for
- * the image. As long as the PSF shape does not vary excessively across the
- * image, the sources which are represented by a PSF (the start) will have very
- * similar second moments. The function pmSourcePSFClump searches a collection of
- * sources with measured moments for a group with moments which are all very
- * similar. The function returns a pmPSFClump structure, representing the
- * centroid and size of the clump in the sigma_x, sigma_y second-moment plane.
+ *
- * The goal is to identify and characterize the stellar clump within the
- * sigma_x, sigma_y second-moment plane.  To do this, an image is constructed to
- * represent this plane.  The units of sigma_x and sigma_y are in image pixels. A
- * pixel in this analysis image represents 0.1 pixels in the input image. The
- * dimensions of the image need only be 10 pixels. The peak pixel in this image
- * (above a threshold of half of the image maximum) is found. The coordinates of
- * this peak pixel represent the 2D mode of the sigma_x, sigma_y distribution.
- * The sources with sigma_x, sigma_y within 0.2 pixels of this value are then
- *  * used to calculate the median and standard deviation of the sigma_x, sigma_y
- * values. These resulting values are returned via the pmPSFClump structure.
+ *
- * The return value indicates the success (TRUE) of the operation.
+ *
- * XXX: Limit the S/N of the candidate sources (part of Metadata)? (TBD).
- * XXX: Save the clump parameters on the Metadata (TBD)
+ *
- */
-pmPSFClump pmSourcePSFClump(
-    psArray *source,   ///< The input pmSource
-    psMetadata *metadata  ///< Contains classification parameters
-);
-/** pmSourceRoughClass()
+ *
- * Based on the specified data values, make a guess at the source
- * classification. The sources are provides as a psArray of pmSource entries.
- * Definable parameters needed to make the classification are provided to the
- * routine with the psMetadata structure. The rules (in SDRS) refer to values which
- * can be extracted from the metadata using the given keywords. Except as noted,
- * the data type for these parameters are psF32.
+ *
- */
-bool pmSourceRoughClass(
-    psArray *source,   ///< The input pmSource
-    psMetadata *metadata,  ///< Contains classification parameters
-    pmPSFClump clump   ///< Statistics about the PSF clump
-);
-/** pmSourceModelGuess()
+ *
- * Convert available data to an initial guess for the given model. This
- * function allocates a pmModel entry for the pmSource structure based on the
- * provided model selection. The method of defining the model parameter guesses
- * are specified for each model below. The guess values are placed in the model
- * parameters. The function returns TRUE on success or FALSE on failure.
+ *
- */
-pmModel *pmSourceModelGuess(
-    pmSource *source,   ///< The input pmSource
-    pmModelType model   ///< The type of model to be created.
-);
-/** pmContourType
+ *
- * Only one type is defined at present.
+ *
- */
-typedef enum {
-    PS_CONTOUR_CRUDE,
-    PS_CONTOUR_UNKNOWN01,
-    PS_CONTOUR_UNKNOWN02
-} pmContourType;
-/** pmSourceContour()
+ *
- * Find points in a contour for the given source at the given level. If type
- * is PM_CONTOUR_CRUDE, the contour is found by starting at the source peak,
- * running along each pixel row until the level is crossed, then interpolating to
- * the level coordinate for that row. This is done for each row, with the
- * starting point determined by the midpoint of the previous row, until the
- * starting point has a value below the contour level. The returned contour
- * consists of two vectors giving the x and y coordinates of the contour levels.
- * This function may be used as part of the model guess inputs.  Other contour
- * types may be specified in the future for more refined contours (TBD)
+ *
- */
-psArray *pmSourceContour(
-    pmSource *source,   ///< The input pmSource
-    const psImage *image,  ///< The input image (float) (this arg should be removed)
-    float level,   ///< The level of the contour
-    pmContourType mode   ///< Currently this must be PS_CONTOUR_CRUDE
-);
-/** pmSourceFitModel()
+ *
- * Fit the requested model to the specified source. The starting guess for the
- * model is given by the input source.model parameter values. The pixels of
- * interest are specified by the source.pixelsand source.maskentries. This
- * function calls psMinimizeLMChi2() on the image data. The function returns TRUE
- * on success or FALSE on failure.
+ *
- */
-bool pmSourceFitModel(
-    pmSource *source,   ///< The input pmSource
-    pmModel *model,   ///< model to be fitted
-    const bool PSF   ///< Treat model as PSF or FLT?
-);
-/** pmModelFitStatus()
+ *
- * This function wraps the call to the model-specific function returned by
- * pmModelFitStatusFunc_GetFunction.  The model-specific function examines the
- * model parameters, parameter errors, Chisq, S/N, and other parameters available
- * from model to decide if the particular fit was successful or not.
+ *
- * XXX: Must code this.
+ *
- */
-bool pmModelFitStatus(
-    pmModel *model                      ///< Add comment.
-);
-/** pmSourceAddModel()
+ *
- * Add the given source model flux to/from the provided image. The boolean
- * option center selects if the source is re-centered to the image center or if
- * it is placed at its centroid location. The boolean option sky selects if the
- * background sky is applied (TRUE) or not. The pixel range in the target image
- * is at most the pixel range specified by the source.pixels image. The success
- * status is returned.
+ *
- */
-bool pmSourceAddModel(
-    psImage *image,   ///< The output image (float)
-    psImage *mask,   ///< The image pixel mask (valid == 0)
-    pmModel *model,   ///< The input pmModel
-    bool center,    ///< A boolean flag that determines whether pixels are centered
-    bool sky        ///< A boolean flag that determines if the sky is subtracted
-);
-/** pmSourceSubModel()
+ *
- * Subtract the given source model flux to/from the provided image. The
- * boolean option center selects if the source is re-centered to the image center
- * or if it is placed at its centroid location. The boolean option sky selects if
- * the background sky is applied (TRUE) or not. The pixel range in the target
- * image is at most the pixel range specified by the source.pixels image. The
- * success status is returned.
+ *
- */
-bool pmSourceSubModel(
-    psImage *image,   ///< The output image (float)
-    psImage *mask,   ///< The image pixel mask (valid == 0)
-    pmModel *model,   ///< The input pmModel
-    bool center,    ///< A boolean flag that determines whether pixels are centered
-    bool sky        ///< A boolean flag that determines if the sky is subtracted
-);
-/**
+ *
- * The function returns both the magnitude of the fit, defined as -2.5log(flux),
- * where the flux is integrated under the model, theoretically from a radius of 0
- * to infinity. In practice, we integrate the model beyond 50sigma.  The aperture magnitude is
- * defined as -2.5log(flux) , where the flux is summed for all pixels which are
- * not excluded by the aperture mask. The model flux is calculated by calling the
- * model-specific function provided by pmModelFlux_GetFunction.
+ *
- * XXX: must code this.
+ *
- */
-bool pmSourcePhotometry(
-    float *fitMag,                      ///< integrated fit magnitude
-    float *obsMag,   ///< aperture flux magnitude
-    pmModel *model,                     ///< model used for photometry
-    psImage *image,                     ///< image pixels to be used
-    psImage *mask                       ///< mask of pixels to ignore
-);
 /**
 …
  * This function converts the source classification into the closest available
  * approximation to the Dophot classification scheme:
+ * XXX EAM : fix this to use current source classification scheme
+ *
  * PM_SOURCE_DEFECT: 8
 …
 );
-/** pmSourceFitModel_v5()
+ *
- * Fit the requested model to the specified source. The starting guess for the
- * model is given by the input source.model parameter values. The pixels of
- * interest are specified by the source.pixelsand source.maskentries. This
- * function calls psMinimizeLMChi2() on the image data. The function returns TRUE
- * on success or FALSE on failure.
+ *
- */
-bool pmSourceFitModel_v5(
-    pmSource *source,   ///< The input pmSource
-    pmModel *model,   ///< model to be fitted
-    const bool PSF   ///< Treat model as PSF or FLT?
-);
-/** pmSourceFitModel_v7()
+ *
- * Fit the requested model to the specified source. The starting guess for the
- * model is given by the input source.model parameter values. The pixels of
- * interest are specified by the source.pixelsand source.maskentries. This
- * function calls psMinimizeLMChi2() on the image data. The function returns TRUE
- * on success or FALSE on failure.
+ *
- */
-bool pmSourceFitModel_v7(
-    pmSource *source,   ///< The input pmSource
-    pmModel *model,   ///< model to be fitted
-    const bool PSF   ///< Treat model as PSF or FLT?
-);
-/** pmSourcePhotometry()
+ *
- * XXX: Need descriptions
+ *
- */
-bool pmSourcePhotometry(
-    float *fitMag,
-    float *obsMag,
-    pmModel *model,
-    psImage *image,
-    psImage *mask
-);
-/** pmModelEval()
+ *
- *  XXX: Need descriptions
+ *
- */
-psF32 pmModelEval(
-    pmModel *model,
-    psImage *image,
-    psS32 col,
-    psS32 row
-);
 #endif

trunk/psModules/src/objects/pmPSF.h

-              r5255
+              r6872
  *  @author EAM, IfA
+ *
  *  @version $Revision: 1.1 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-10-10 19:53:40 $
+ *  @version $Revision: 1.2 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 # ifndef PM_PSF_H
 # define PM_PSF_H
 /** pmPSF data structure
 …
     pmModelType type;   ///< PSF Model in use
     psArray *params;   ///< Model parameters (psPolynomial2D)
+    float chisq;   ///< PSF goodness statistic
+    float ApResid;                      ///< ???
+    float dApResid;                     ///< ???
+    float skyBias;                      ///< ???
+    psPolynomial1D *ChiTrend;  ///< Chisq vs flux fit (correction for systematic errors)
+    psPolynomial4D *ApTrend;  ///< ApResid vs (x,y,rflux) (rflux = ten(0.4*mInst)
+    pmGrowthCurve *growth;  ///< apMag vs Radius
+    float ApResid;   ///< apMag - psfMag (for PSF stars)
+    float dApResid;   ///< scatter of ApResid
+    float skyBias;   ///< implied residual sky offset from ApResid fit
+    float skySat;   ///< roll-over of ApResid fit
+    float chisq;   ///< PSF goodness statistic (unused??)
     int nPSFstars;   ///< number of stars used to measure PSF
+    int nApResid;   ///< number of stars used to measure ApResid
+    bool poissonErrors;
+}
 pmPSF;
+typedef enum {
+    PM_PSF_NONE,
+    PM_PSF_CONSTANT,
+    PM_PSF_SKYBIAS,
+    PM_PSF_SKYSAT,
+    PM_PSF_XY_LIN,
+    PM_PSF_XY_QUAD,
+    PM_PSF_SKY_XY_LIN,
+    PM_PSF_SKYSAT_XY_LIN,
+    PM_PSF_ALL
+} pmPSF_ApTrendOptions;
 /**
 …
  */
 pmPSF *pmPSFAlloc(
+    pmModelType type                    ///< Add comment
+    pmModelType type,   // type of model for PSF
+    bool poissonErrors   ///< use poissonian errors or not?
 );
 …
 );
+bool pmPSF_MaskApTrend (pmPSF *psf, pmPSF_ApTrendOptions option);
 # endif

trunk/psModules/src/objects/pmPSFtry.h

-              r5844
+              r6872
  *  @author EAM, IfA
+ *
  *  @version $Revision: 1.3 $ $Name: not supported by cvs2svn $
  *  @date $Date: 2005-12-24 01:24:32 $
+ *  @version $Revision: 1.4 $ $Name: not supported by cvs2svn $
+ *  @date $Date: 2006-04-17 18:01:05 $
+ *
  *  Copyright 2004 Maui High Performance Computing Center, University of Hawaii
 …
 /**
+ *
+ *
  * This structure contains a pointer to the collection of sources which will
  * be used to test the PSF model form. It lists the pmModelType type of model
  * being tests, and contains an element to store the resulting psf
  * representation. In addition, this structure carries the complete collection of
  * FLT (floating parameter) and PSF (fixed parameter) model fits to each of the
  * sources modelFLT and modelPSF. It also contains a mask which is set by the
+ * EXT (floating parameter) and PSF (fixed parameter) model fits to each of the
+ * sources modelEXT and modelPSF. It also contains a mask which is set by the
  * model fitting and psf fitting steps. For each model, the value of the quality
  * metric is stored in the vector metric and the fitted instrumental magnitude is
 …
  * ultimate metric to intercompare multiple types of PSF models is the value of
  * the aperture correction scatter.
+ *
+ *
  * XXX: There are many more members in the SDRS then in the prototype code.
  * I stuck with the prototype code.
+ *
+ *
+ *
+ *
  */
 typedef struct
 …
     pmPSF      *psf;                    ///< Add comment.
     psArray    *sources;                ///< pointers to the original sources
     psArray    *modelFLT;               ///< model fits, floating parameters
+    psArray    *modelEXT;               ///< model fits, floating parameters
     psArray    *modelPSF;               ///< model fits, PSF parameters
     psVector   *mask;                   ///< Add comment.
 …
 /** pmPSFtryMaskValues
+ *
+ *
  * The following datatype defines the masks used by the pmPSFtry analysis to
  * identify sources which should or should not be included in the analysis.
+ *
+ *
  */
 enum {
+typedef enum {
     PSFTRY_MASK_CLEAR    = 0x00,        ///< Add comment.
     PSFTRY_MASK_OUTLIER  = 0x01,        ///< 1: outlier in psf polynomial fit (provided by psPolynomials)
     PSFTRY_MASK_FLT_FAIL = 0x02,        ///< 2: flt model failed to converge
+    PSFTRY_MASK_EXT_FAIL = 0x02,        ///< 2: ext model failed to converge
     PSFTRY_MASK_PSF_FAIL = 0x04,        ///< 3: psf model failed to converge
     PSFTRY_MASK_BAD_PHOT = 0x08,        ///< 4: invalid source photometry
 …
 /** pmPSFtryAlloc()
+ *
+ *
  * Allocate a pmPSFtry data structure.
+ *
+ *
  */
 pmPSFtry *pmPSFtryAlloc(
     psArray *stars,                     ///< Add comment.
+    char *modelName                     ///< Add comment.
+    char *modelName,                     ///< Add comment.
+    bool poissonErrors   // use poissonian or constant errors?
 );
 /** pmPSFtryModel()
+ *
+ *
  * This function takes the input collection of sources and performs a complete
  * analysis to determine a PSF model of the given type (specified by model name).
  * The result is a pmPSFtry with the results of the analysis.
+ *
+ *
  */
 pmPSFtry *pmPSFtryModel(
     psArray *sources,                   ///< Add comment.
     char *modelName,                    ///< Add comment.
+    float radius                        ///< Add comment.
+    float radius,                     ///< Add comment.
+    bool poissonErrors   // use poissonian or constant errors?
 );
 /** pmPSFtryMetric()
+ *
+ *
  * This function is used to measure the PSF model metric for the set of
  * results contained in the pmPSFtry structure.
+ *
+ *
  */
 bool pmPSFtryMetric(

Context Navigation

Legend:

trunk/psModules/configure.ac

trunk/psModules/src/Makefile.am

trunk/psModules/src/astrom/Makefile.am

trunk/psModules/src/astrom/pmAstrometry.c

trunk/psModules/src/astrom/pmAstrometry.h

trunk/psModules/src/astrom/pmAstrometryObjects.c

trunk/psModules/src/astrom/pmAstrometryObjects.h

trunk/psModules/src/camera/Makefile.am

trunk/psModules/src/config/pmConfig.h

trunk/psModules/src/detrend/Makefile.am

trunk/psModules/src/detrend/pmFlatField.c

trunk/psModules/src/detrend/pmFlatField.h

trunk/psModules/src/detrend/pmMaskBadPixels.c

trunk/psModules/src/detrend/pmMaskBadPixels.h

trunk/psModules/src/detrend/pmMaskBadPixelsErrors.h

trunk/psModules/src/detrend/pmNonLinear.c

trunk/psModules/src/detrend/pmNonLinear.h

trunk/psModules/src/imcombine/Makefile.am

trunk/psModules/src/imcombine/pmReadoutCombine.h

trunk/psModules/src/imsubtract/Makefile.am

trunk/psModules/src/imsubtract/pmSubtractBias.h

trunk/psModules/src/imsubtract/pmSubtractSky.h

trunk/psModules/src/objects/Makefile.am

trunk/psModules/src/objects/models/pmModel_PGAUSS.c

trunk/psModules/src/objects/models/pmModel_QGAUSS.c

trunk/psModules/src/objects/models/pmModel_RGAUSS.c

trunk/psModules/src/objects/models/pmModel_SGAUSS.c

trunk/psModules/src/objects/models/pmModel_ZGAUSS.c

trunk/psModules/src/objects/pmModelGroup.c

trunk/psModules/src/objects/pmModelGroup.h

trunk/psModules/src/objects/pmObjects.h

trunk/psModules/src/objects/pmPSF.h

trunk/psModules/src/objects/pmPSFtry.h

Download in other formats: