IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 16820


Ignore:
Timestamp:
Mar 4, 2008, 3:10:25 PM (18 years ago)
Author:
eugene
Message:

merge changes from eam_branch_20080229: better flag definitions, cleanup footprint code (move into psphotFindDetections), push mask selection into called functions

Location:
trunk/psphot/src
Files:
1 added
26 edited

Legend:

Unmodified
Added
Removed
  • trunk/psphot/src/Makefile.am

    r15936 r16820  
    2727        psphotModelBackground.c  \
    2828        psphotSubtractBackground.c \
     29        psphotFindDetections.c   \
    2930        psphotFindPeaks.c        \
    3031        psphotSourceStats.c      \
  • trunk/psphot/src/pmFootprint.c

    r16048 r16820  
    12411241   return peaks;
    12421242}
     1243
     1244void pmDetectionsFree (pmDetections *detections) {
     1245
     1246  if (!detections) return;
     1247
     1248  psFree (detections->footprints);
     1249  psFree (detections->peaks);
     1250  psFree (detections->oldPeaks);
     1251  return;
     1252}
     1253
     1254// generate a pmDetections container with empty (allocated) footprints and peaks containers
     1255pmDetections *pmDetectionsAlloc() {
     1256
     1257    pmDetections *detections = (pmDetections *)psAlloc(sizeof(pmDetections));
     1258    psMemSetDeallocator(detections, (psFreeFunc) pmDetectionsFree);
     1259
     1260    detections->footprints = NULL;
     1261    detections->peaks      = NULL;
     1262    detections->oldPeaks   = NULL;
     1263    detections->last       = 0;
     1264
     1265    return (detections);
     1266}
     1267
     1268/************************************************************************************************************/
     1269/*
     1270 * Cull a set of peaks contained in a psArray of pmFootprints
     1271 */
     1272psErrorCode
     1273psphotCullPeaks(const psImage *image,   // the image wherein lives the footprint
     1274                const psImage *weight,  // corresponding variance image
     1275                const psMetadata *recipe,
     1276                psArray *footprints) {  // array of pmFootprints
     1277    bool status = false;
     1278    float nsigma_delta = psMetadataLookupF32(&status, recipe, "FOOTPRINT_CULL_NSIGMA_DELTA");
     1279    if (!status) {
     1280        nsigma_delta = 0; // min.
     1281    }
     1282    float nsigma_min = psMetadataLookupF32(&status, recipe, "FOOTPRINT_CULL_NSIGMA_MIN");
     1283    if (!status) {
     1284        nsigma_min = 0;
     1285    }
     1286    const float skyStdev = psMetadataLookupF32(NULL, recipe, "SKY_STDEV");
     1287
     1288    return pmFootprintArrayCullPeaks(image, weight, footprints,
     1289                                     nsigma_delta, nsigma_min*skyStdev);
     1290}
     1291
  • trunk/psphot/src/pmFootprint.h

    r14655 r16820  
    5252psArray *pmFootprintArrayToPeaks(const psArray *footprints);
    5353
     54typedef struct {
     55  psArray *footprints;        // collection of footprints in the image
     56  psArray *peaks;             // collection of all peaks contained by the footprints
     57  psArray *oldPeaks;          // collection of all peaks previously found
     58  int last;
     59} pmDetections;
     60
     61pmDetections *pmDetectionsAlloc ();
     62
     63
    5464#endif
  • trunk/psphot/src/psphot.h

    r15936 r16820  
    99
    1010#include "psphotErrorCodes.h"
     11#include "pmFootprint.h"
     12
    1113#define PSPHOT_RECIPE "PSPHOT" // Name of the recipe to use
    1214
     
    1618psString        psphotVersionLong(void);
    1719
    18 bool            psphotModelTest (pmConfig *config, const pmFPAview *view, psMetadata *recipe, psMaskType maskVal, psMaskType mark);
     20bool            psphotModelTest (pmConfig *config, const pmFPAview *view, psMetadata *recipe);
    1921bool            psphotReadout (pmConfig *config, const pmFPAview *view);
    20 bool            psphotReadoutCleanup (pmConfig *config, pmReadout *readout, psMetadata *recipe, pmPSF *psf, psArray *sources);
     22bool            psphotReadoutCleanup (pmConfig *config, pmReadout *readout, psMetadata *recipe, pmDetections *detections, pmPSF *psf, psArray *sources);
    2123bool            psphotDefineFiles (pmConfig *config, pmFPAfile *input);
    22 
    2324
    2425// XXX test functions
     
    2627
    2728// psphotReadout functions
    28 bool            psphotModelBackground (pmConfig *config, const pmFPAview *view, const char *filename, psMaskType maskVal);
    29 bool            psphotSubtractBackground (pmConfig *config, const pmFPAview *view, const char *filename, psMaskType maskVal) ;
    30 psArray        *psphotFindPeaks (pmReadout *readout, psMetadata *recipe,
    31                                  const bool returnFootprints, const int pass, psMaskType maskVal);
    32 #include "pmFootprint.h"
    33 psErrorCode     psphotCullPeaks(const psImage *img, const psImage *weight,
    34                                 const psMetadata *recipe, psArray *footprints);
    35 psArray        *psphotSourceStats (pmReadout *readout, psMetadata *recipe, psArray *allpeaks, psMaskType maskVal, psMaskType mark);
    36 bool            psphotRoughClass (psArray *sources, psMetadata *recipe, const bool findPsfClump, psMaskType maskSat);
     29bool            psphotModelBackground (pmConfig *config, const pmFPAview *view, const char *filename);
     30bool            psphotSubtractBackground (pmConfig *config, const pmFPAview *view, const char *filename) ;
     31pmDetections   *psphotFindDetections (pmDetections *detections, pmReadout *readout, psMetadata *recipe);
     32
     33psArray        *psphotSourceStats (pmReadout *readout, psMetadata *recipe, pmDetections *detections);
     34bool            psphotRoughClass (psArray *sources, psMetadata *recipe, const bool findPsfClump);
    3735bool            psphotBasicDeblend (psArray *sources, psMetadata *recipe);
    38 pmPSF          *psphotChoosePSF (pmReadout *readout, psArray *sources, psMetadata *recipe, psMaskType maskVal, psMaskType mark);
     36pmPSF          *psphotChoosePSF (pmReadout *readout, psArray *sources, psMetadata *recipe);
    3937bool            psphotPSFstats (pmReadout *readout, psMetadata *recipe, pmPSF *psf);
    4038bool            psphotMomentsStats (pmReadout *readout, psMetadata *recipe, psArray *sources);
    41 #if NOT_IN_LIBPSPHOT
    42 bool            psphotEnsemblePSF (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, bool final);
    43 #endif
    44 bool            psphotFitSourcesLinear (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, bool final, psMaskType maskVal);
    45 bool            psphotGuessModels (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, psMaskType maskVal);
    46 bool            psphotBlendFit (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, psMaskType maskVal);
     39bool            psphotFitSourcesLinear (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, bool final);
     40bool            psphotGuessModels (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf);
     41bool            psphotBlendFit (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf);
    4742bool            psphotReplaceUnfit (psArray *sources, psMaskType maskVal);
    48 bool            psphotReplaceAll (psArray *sources, psMaskType maskVal);
    49 bool            psphotApResid (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, psMaskType maskVal, psMaskType mark);
     43bool            psphotReplaceAll (psArray *sources, psMetadata *recipe);
     44bool            psphotApResid (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf);
     45bool            psphotMagnitudes (psArray *sources, psMetadata *recipe, pmPSF *psf, pmReadout *background);
     46bool            psphotSkyReplace (pmConfig *config, const pmFPAview *view);
     47bool            psphotExtendedSources (pmReadout *readout, psArray *sources, psMetadata *recipe);
     48
     49// used by psphotFindDetections
     50psArray        *psphotFindPeaks (pmReadout *readout, psMetadata *recipe, const bool returnFootprints, const int pass, psMaskType maskVal);
     51psErrorCode     psphotCullPeaks(const psImage *img, const psImage *weight, const psMetadata *recipe, psArray *footprints);
     52
     53// used by ApResid
    5054bool            psphotMagErrorScale (float *errorScale, float *errorFloor, psVector *dMag, psVector *dap, psVector *mask, int nGroup);
    5155bool            psphotApResidTrend (pmReadout *readout, pmPSF *psf, int Npsf, int scale, float *errorScale, float *errorFloor, psVector *mask, psVector *xPos, psVector *yPos, psVector *apResid, psVector *dMag);
    52 bool            psphotMagnitudes (psArray *sources, psMetadata *recipe, pmPSF *psf, pmReadout *background, psMaskType maskVal, psMaskType mark);
    53 bool            psphotSkyReplace (pmConfig *config, const pmFPAview *view);
    5456
    5557// basic support functions
    5658void            psphotModelClassInit (void);
    5759bool            psphotGrowthCurve (pmReadout *readout, pmPSF *psf, bool ignore, psMaskType maskVal);
    58 bool            psphotMaskReadout (pmReadout *readout, psMetadata *recipe, psMaskType maskVal);
     60bool            psphotSetMaskAndWeight (pmConfig *config, pmReadout *readout, psMetadata *recipe);
    5961void            psphotSourceFreePixels (psArray *sources);
    6062
     
    100102pmPSF          *psphotLoadPSF (pmConfig *config, const pmFPAview *view, psMetadata *recipe);
    101103bool            psphotSetHeaderNstars (psMetadata *recipe, psArray *sources);
    102 bool            psphotAddNoise (pmReadout *readout, psArray *sources, psMetadata *recipe, bool add, psMaskType maskVal);
     104bool            psphotAddNoise (pmReadout *readout, psArray *sources, psMetadata *recipe);
     105bool            psphotSubNoise (pmReadout *readout, psArray *sources, psMetadata *recipe);
     106bool            psphotAddOrSubNoise (pmReadout *readout, psArray *sources, psMetadata *recipe, bool add);
    103107bool            psphotRadialPlot (int *kapa, const char *filename, pmSource *source);
    104 bool            psphotSourcePlots (pmReadout *readout, psArray *sources, psMetadata *recipe, psMaskType maskVal);
     108bool            psphotSourcePlots (pmReadout *readout, psArray *sources, psMetadata *recipe);
    105109bool            psphotMosaicSubimage (psImage *outImage, pmSource *source, int Xo, int Yo, int DX, int DY);
    106110
     
    109113bool            psphotSetState (pmSource *source, bool curState, psMaskType maskVal);
    110114bool            psphotDeblendSatstars (psArray *sources, psMetadata *recipe);
    111 bool            psphotSourceSize (pmReadout *readout, psArray *sources, psMetadata *recipe);
     115bool            psphotSourceSize (pmReadout *readout, psArray *sources, psMetadata *recipe, long first);
    112116
    113117bool            psphotMakeResiduals (psArray *sources, psMetadata *recipe, pmPSF *psf, psMaskType maskVal);
    114118
    115 bool            psphotExtendedSources (pmReadout *readout, psArray *sources, psMetadata *recipe, psMaskType maskVal);
    116119bool            psphotPSFConvModel (pmSource *source, psMetadata *recipe, psMaskType maskVal);
    117120psKernel       *psphotKernelFromPSF (pmSource *source, int nPix);
  • trunk/psphot/src/psphotAddNoise.c

    r14655 r16820  
    11# include "psphotInternal.h"
    22
    3 bool psphotAddNoise (pmReadout *readout, psArray *sources, psMetadata *recipe, bool add, psMaskType maskVal) {
     3bool psphotAddNoise (pmReadout *readout, psArray *sources, psMetadata *recipe) {
     4  return psphotAddOrSubNoise (readout, sources, recipe, true);
     5}
     6
     7bool psphotSubNoise (pmReadout *readout, psArray *sources, psMetadata *recipe) {
     8  return psphotAddOrSubNoise (readout, sources, recipe, false);
     9}
     10
     11bool psphotAddOrSubNoise (pmReadout *readout, psArray *sources, psMetadata *recipe, bool add) {
    412
    513    bool status = false;
     
    1321
    1422    psTimerStart ("psphot");
     23
     24    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     25    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     26    assert (maskVal);
    1527
    1628    // increase variance by factor*(object noise):
  • trunk/psphot/src/psphotApResid.c

    r15143 r16820  
    33# define SKIPSTAR(MSG) { psTrace ("psphot", 3, "invalid : %s", MSG); continue; }
    44// measure the aperture residual statistics and 2D variations
    5 bool psphotApResid (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, psMaskType maskVal, psMaskType mark)
     5
     6bool psphotApResid (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf)
    67{
    78    int Nfail = 0;
     
    1920
    2021    psTimerStart ("psphot");
     22
     23    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     24    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     25    assert (maskVal);
     26
     27    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     28    psMaskType mark = psMetadataLookupU8(&status, recipe, "MASK.MARK"); // Mask value for bad pixels
     29    assert (mark);
    2130
    2231    // S/N limit to perform full non-linear fits
  • trunk/psphot/src/psphotBlendFit.c

    r15800 r16820  
    22
    33// XXX I don't like this name
    4 bool psphotBlendFit (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, psMaskType maskVal) {
     4bool psphotBlendFit (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf) {
    55
    66    int Nfit = 0;
     
    1111
    1212    psTimerStart ("psphot");
     13
     14    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     15    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     16    assert (maskVal);
    1317
    1418    // source analysis is done in S/N order (brightest first)
     
    7377        // try fitting PSFs, then try extended sources
    7478        // these functions subtract the resulting fitted source (XXX and update the modelFlux?)
    75         // XXX consider conditions under which the source has EXT fit
     79        // XXX re-consider conditions under which the source has EXT fit:
     80        // I should try EXT if the source size measurement says it is large
    7681        if (psphotFitBlend (readout, source, psf, maskVal)) {
    7782            psTrace ("psphot", 5, "source at %7.1f, %7.1f is psf", source->moments->x, source->moments->y);
     
    9297        pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
    9398        source->mode |= PM_SOURCE_MODE_SUBTRACTED;
    94         source->mode |= PM_SOURCE_MODE_TEMPSUB;
     99    }
     100
     101    if (psTraceGetLevel("psphot") >= 6) {
     102      psphotSaveImage (NULL, readout->image,  "image.v2.fits");
    95103    }
    96104
  • trunk/psphot/src/psphotChoosePSF.c

    r15023 r16820  
    22
    33// try PSF models and select best option
    4 pmPSF *psphotChoosePSF (pmReadout *readout, psArray *sources, psMetadata *recipe, psMaskType maskVal, psMaskType mark) {
     4pmPSF *psphotChoosePSF (pmReadout *readout, psArray *sources, psMetadata *recipe) {
    55
    66    bool status;
    77
    88    psTimerStart ("psphot");
     9
     10    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     11    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     12    assert (maskVal);
     13
     14    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     15    psMaskType mark = psMetadataLookupU8(&status, recipe, "MASK.MARK"); // Mask value for bad pixels
     16    assert (mark);
    917
    1018    // examine PSF sources in S/N order (brightest first)
  • trunk/psphot/src/psphotEvalPSF.c

    r15060 r16820  
    134134    // assign PM_SOURCE_MODE_GOODSTAR to bright objects within PSF region of dparams[]
    135135    keep = TRUE;
    136     keep &= (fabs(nSx) < PSF_SHAPE_NSIGMA);
    137     keep &= (fabs(nSy) < PSF_SHAPE_NSIGMA);
     136    // keep &= (fabs(nSx) < PSF_SHAPE_NSIGMA);
     137    // keep &= (fabs(nSy) < PSF_SHAPE_NSIGMA);
    138138    keep &= (SN > PSF_MIN_SN);
    139139    keep &= (Chi < PSF_MAX_CHI);
  • trunk/psphot/src/psphotExtendedSources.c

    r15562 r16820  
    11# include "psphot.h"
    22
    3 bool psphotExtendedSources (pmReadout *readout, psArray *sources, psMetadata *recipe, psMaskType maskVal) {
     3bool psphotExtendedSources (pmReadout *readout, psArray *sources, psMetadata *recipe) {
    44
    55    bool status;
    66    int Next = 0;
    77    int Npsf = 0;
     8
     9    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     10    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     11    assert (maskVal);
    812
    913    // S/N limit to perform full non-linear fits
     
    117121        pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
    118122        source->mode |= PM_SOURCE_MODE_SUBTRACTED;
    119         source->mode |= PM_SOURCE_MODE_TEMPSUB;
    120123    }
    121124
  • trunk/psphot/src/psphotFindPeaks.c

    r15962 r16820  
    166166        psFree(peaks);
    167167        peaks = footprints;             // well, you know what I mean
     168
     169        psLogMsg ("psphot", PS_LOG_INFO, "%ld footprints: %f sec\n", peaks->n, psTimerMark ("psphot"));
    168170    }
    169171
     
    174176}
    175177
    176 
    177 /************************************************************************************************************/
    178 /*
    179  * Cull a set of peaks contained in a psArray of pmFootprints
    180  */
    181 psErrorCode
    182 psphotCullPeaks(const psImage *image,   // the image wherein lives the footprint
    183                 const psImage *weight,  // corresponding variance image
    184                 const psMetadata *recipe,
    185                 psArray *footprints) {  // array of pmFootprints
    186     bool status = false;
    187     float nsigma_delta = psMetadataLookupF32(&status, recipe, "FOOTPRINT_CULL_NSIGMA_DELTA");
    188     if (!status) {
    189         nsigma_delta = 0; // min.
    190     }
    191     float nsigma_min = psMetadataLookupF32(&status, recipe, "FOOTPRINT_CULL_NSIGMA_MIN");
    192     if (!status) {
    193         nsigma_min = 0;
    194     }
    195     const float skyStdev = psMetadataLookupF32(NULL, recipe, "SKY_STDEV");
    196 
    197     return pmFootprintArrayCullPeaks(image, weight, footprints,
    198                                      nsigma_delta, nsigma_min*skyStdev);
    199 }
  • trunk/psphot/src/psphotFitSourcesLinear.c

    r15132 r16820  
    1212static bool SetBorderMatrixElements (psSparseBorder *border, pmReadout *readout, psArray *sources, bool constant_weights, int SKY_FIT_ORDER);
    1313
    14 bool psphotFitSourcesLinear (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, bool final, psMaskType maskVal) {
     14bool psphotFitSourcesLinear (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, bool final) {
    1515
    1616    bool status;
     
    2121
    2222    psTimerStart ("psphot");
     23
     24    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     25    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     26    assert (maskVal);
    2327
    2428    // source analysis is done in spatial order
     
    187191        pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
    188192        source->mode |= PM_SOURCE_MODE_SUBTRACTED;
    189         if (!final) source->mode |= PM_SOURCE_MODE_TEMPSUB;
    190         // XXX not sure about the use of TEMPSUB
    191193    }
    192194
  • trunk/psphot/src/psphotGuessModels.c

    r13900 r16820  
    1717
    1818// construct an initial PSF model for each object
    19 bool psphotGuessModels (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, psMaskType maskVal) {
     19bool psphotGuessModels (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf) {
    2020
    21   psTimerStart ("psphot");
     21    bool status;
    2222
    23   // setup the PSF fit radius details
    24   psphotInitRadiusPSF (recipe, psf->type);
     23    psTimerStart ("psphot");
    2524
    26   for (int i = 0; i < sources->n; i++) {
    27     pmSource *source = sources->data[i];
     25    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     26    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     27    assert (maskVal);
    2828
    29     // skip non-astronomical objects (very likely defects)
    30     if (source->type == PM_SOURCE_TYPE_DEFECT) continue;
    31     if (source->type == PM_SOURCE_TYPE_SATURATED) continue;
     29    // setup the PSF fit radius details
     30    psphotInitRadiusPSF (recipe, psf->type);
    3231
    33     // XXX if a source is faint, it will not have moments measured.
    34     // it must be modelled as a PSF.  In this case, we need to use
    35     // the peak centroid to get the coordinates and get the peak flux
    36     // from the image?
    37     pmModel *modelEXT;
    38     if (!source->moments) {
    39         modelEXT = wildGuess(source, psf);
    40     } else {
    41         // use the source moments, etc to guess basic model parameters
    42         modelEXT = pmSourceModelGuess (source, psf->type);
    43         if (!modelEXT) {
    44             modelEXT = wildGuess(source, psf);
    45         }
    46         // these valuse are set in pmSourceModelGuess, should this rule be in there as well?
    47         if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
    48             modelEXT->params->data.F32[PM_PAR_XPOS] = source->moments->x;
    49             modelEXT->params->data.F32[PM_PAR_YPOS] = source->moments->y;
    50         } else {
    51             modelEXT->params->data.F32[PM_PAR_XPOS] = source->peak->xf;
    52             modelEXT->params->data.F32[PM_PAR_YPOS] = source->peak->yf;
    53         }
     32    for (int i = 0; i < sources->n; i++) {
     33        pmSource *source = sources->data[i];
     34
     35        // skip non-astronomical objects (very likely defects)
     36        if (source->type == PM_SOURCE_TYPE_DEFECT) continue;
     37        if (source->type == PM_SOURCE_TYPE_SATURATED) continue;
     38
     39        // XXX if a source is faint, it will not have moments measured.
     40        // it must be modelled as a PSF.  In this case, we need to use
     41        // the peak centroid to get the coordinates and get the peak flux
     42        // from the image?
     43        pmModel *modelEXT;
     44        if (!source->moments) {
     45            modelEXT = wildGuess(source, psf);
     46        } else {
     47            // use the source moments, etc to guess basic model parameters
     48            modelEXT = pmSourceModelGuess (source, psf->type);
     49            if (!modelEXT) {
     50                modelEXT = wildGuess(source, psf);
     51            }
     52            // these valuse are set in pmSourceModelGuess, should this rule be in there as well?
     53            if (source->mode &  PM_SOURCE_MODE_SATSTAR) {
     54                modelEXT->params->data.F32[PM_PAR_XPOS] = source->moments->x;
     55                modelEXT->params->data.F32[PM_PAR_YPOS] = source->moments->y;
     56            } else {
     57                modelEXT->params->data.F32[PM_PAR_XPOS] = source->peak->xf;
     58                modelEXT->params->data.F32[PM_PAR_YPOS] = source->peak->yf;
     59            }
     60        }
     61
     62        // set PSF parameters for this model (apply 2D shape model)
     63        pmModel *modelPSF = pmModelFromPSF (modelEXT, psf);
     64        if (modelPSF == NULL) {
     65            psError(PSPHOT_ERR_PSF, false,
     66                    "Failed to determine PSF model at r,c = (%d,%d); trying centre of image",
     67                    source->peak->y, source->peak->x);
     68            //
     69            // Try the centre of the image
     70            //
     71            modelEXT->params->data.F32[PM_PAR_XPOS] = 0.5*readout->image->numCols;
     72            modelEXT->params->data.F32[PM_PAR_YPOS] = 0.5*readout->image->numRows;
     73            modelPSF = pmModelFromPSF (modelEXT, psf);
     74            if (modelPSF == NULL) {
     75                psError(PSPHOT_ERR_PSF, false,
     76                        "Failed to determine PSF model at centre of image");
     77                psFree(modelEXT);
     78                return false;
     79            }
     80
     81            source->mode |= PM_SOURCE_MODE_BADPSF;
     82        }
     83        psFree (modelEXT);
     84
     85        // XXX need to define the guess flux?
     86        // set the fit radius based on the object flux limit and the model
     87        psphotCheckRadiusPSF (readout, source, modelPSF);
     88
     89        // set the source PSF model
     90        source->modelPSF = modelPSF;
     91        source->modelPSF->residuals = psf->residuals;
     92
     93        pmSourceCacheModel (source, maskVal);
    5494    }
    55 
    56     // set PSF parameters for this model (apply 2D shape model)
    57     pmModel *modelPSF = pmModelFromPSF (modelEXT, psf);
    58     if (modelPSF == NULL) {
    59         psError(PSPHOT_ERR_PSF, false,
    60                 "Failed to determine PSF model at r,c = (%d,%d); trying centre of image",
    61                 source->peak->y, source->peak->x);
    62         //
    63         // Try the centre of the image
    64         //
    65         modelEXT->params->data.F32[PM_PAR_XPOS] = 0.5*readout->image->numCols;
    66         modelEXT->params->data.F32[PM_PAR_YPOS] = 0.5*readout->image->numRows;
    67         modelPSF = pmModelFromPSF (modelEXT, psf);
    68         if (modelPSF == NULL) {
    69             psError(PSPHOT_ERR_PSF, false,
    70                     "Failed to determine PSF model at centre of image");
    71             psFree(modelEXT);
    72             return false;
    73         }
    74 
    75         source->mode |= PM_SOURCE_MODE_BADPSF;
    76     }
    77     psFree (modelEXT);
    78 
    79     // XXX need to define the guess flux?
    80     // set the fit radius based on the object flux limit and the model
    81     psphotCheckRadiusPSF (readout, source, modelPSF);
    82 
    83     // set the source PSF model
    84     source->modelPSF = modelPSF;
    85     source->modelPSF->residuals = psf->residuals;
    86 
    87     pmSourceCacheModel (source, maskVal);
    88   }
    89   psLogMsg ("psphot.models", 4, "built models for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot"));
    90   return true;
     95    psLogMsg ("psphot.models", 4, "built models for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot"));
     96    return true;
    9197}
    9298
  • trunk/psphot/src/psphotMagnitudes.c

    r13900 r16820  
    11# include "psphotInternal.h"
    22
    3 bool psphotMagnitudes(psArray *sources,
    4                       psMetadata *recipe,
    5                       pmPSF *psf,
    6                       pmReadout *background,
    7                       psMaskType maskVal,
    8                       psMaskType mark)
    9 {
     3bool psphotMagnitudes(psArray *sources, psMetadata *recipe, pmPSF *psf, pmReadout *background) {
     4
    105    bool status = false;
    116    int Nap = 0;
    127
    138    psTimerStart ("psphot");
     9
     10    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     11    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     12    assert (maskVal);
     13
     14    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     15    psMaskType mark = psMetadataLookupU8(&status, recipe, "MASK.MARK"); // Mask value for bad pixels
     16    assert (mark);
    1417
    1518    pmSourceMagnitudesInit (recipe);
  • trunk/psphot/src/psphotMaskReadout.c

    r13900 r16820  
    11# include "psphotInternal.h"
    22
    3 bool psphotMaskReadout (pmReadout *readout, psMetadata *recipe, psMaskType maskVal) {
     3// generate mask and weight if not defined, additional mask for restricted subregion
     4bool psphotSetMaskAndWeight (pmConfig *config, pmReadout *readout, psMetadata *recipe) {
    45
    56    bool status;
     7
     8    // ** Interpret the mask values:
     9
     10    // single-bit named masks
     11    psMaskType maskMark = pmConfigMask("MARK",     config); // Mask value for marking
     12    psMetadataAddU8 (recipe, PS_LIST_TAIL, "MASK.MARK", PS_META_REPLACE, "user-defined mask", maskMark);
     13
     14    psMaskType maskSat  = pmConfigMask("SAT",      config); // Mask value for saturated pixels
     15    psMetadataAddU8 (recipe, PS_LIST_TAIL, "MASK.SAT", PS_META_REPLACE, "user-defined mask", maskSat);
     16
     17    psMaskType maskBad  = pmConfigMask("BAD",      config); // Mask value for bad pixels
     18    psMetadataAddU8 (recipe, PS_LIST_TAIL, "MASK.BAD", PS_META_REPLACE, "user-defined mask", maskBad);
     19
     20    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     21    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "PSPHOT"); // Mask value for bad pixels
     22    if (!maskVal) {
     23      const char *maskValStr = psMetadataLookupStr(NULL, recipe, "MASKVAL"); // String with mask names
     24      if (!maskValStr) {
     25        psError(PSPHOT_ERR_CONFIG, false, "Missing recipe item: MASKVAL(STR)");
     26        return false;
     27      }
     28      maskVal = pmConfigMask(maskValStr, config); // Mask values to mask against
     29      psMetadataAddU8 (recipe, PS_LIST_TAIL, "MASK.PSPHOT", PS_META_REPLACE, "user-defined mask", maskVal);
     30      assert (maskVal);
     31    }
     32
     33    // generate mask & weight images if they don't already exit
     34    if (!readout->mask) {
     35        if (!pmReadoutGenerateMask(readout, maskSat, maskBad)) {
     36            psError (PSPHOT_ERR_CONFIG, false, "trouble creating mask");
     37            return false;
     38        }
     39    }
     40    if (!readout->weight) {
     41        if (!pmReadoutGenerateWeight(readout, true)) {
     42            psError (PSPHOT_ERR_CONFIG, false, "trouble creating weight");
     43            return false;
     44        }
     45    }
    646
    747    // mask the excluded outer pixels
     
    1959
    2060    // psImageKeepRegion assumes the region refers to the parent coordinates
    21     psImageKeepRegion (readout->mask, keep, "OR", maskVal);
     61    psImageKeepRegion (readout->mask, keep, "OR", maskBad);
     62
     63    // test output of files at this stage
     64    if (psTraceGetLevel("psphot") >= 5) {
     65        psphotSaveImage (NULL, readout->image,  "image.fits");
     66        psphotSaveImage (NULL, readout->mask,   "mask.fits");
     67        psphotSaveImage (NULL, readout->weight, "weight.fits");
     68    }
    2269
    2370    return true;
  • trunk/psphot/src/psphotModelBackground.c

    r15936 r16820  
    44// generate the median in NxN boxes, clipping heavily
    55// linear interpolation to generate full-scale model
    6 bool psphotModelBackground (pmConfig *config, const pmFPAview *view, const char *filename, psMaskType maskVal)
     6bool psphotModelBackground (pmConfig *config, const pmFPAview *view, const char *filename)
    77{
    88    bool status = true;
     
    2020    // select the appropriate recipe information
    2121    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
     22    assert (recipe);
     23
     24    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     25    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     26    assert (maskVal);
    2227
    2328    // user supplied seed, if available
  • trunk/psphot/src/psphotModelTest.c

    r14799 r16820  
    33
    44// XXX add more test information?
    5 bool psphotModelTest (pmConfig *config, const pmFPAview *view, psMetadata *recipe, psMaskType maskVal, psMaskType mark) {
     5bool psphotModelTest (pmConfig *config, const pmFPAview *view, psMetadata *recipe) {
    66
    77    bool status;
     
    1111    pmPSF *psf = NULL;
    1212    pmSourceFitMode fitMode;
     13
     14    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     15    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     16    psMaskType mark    = psMetadataLookupU8(&status, recipe, "MASK.MARK"); // Mask value for bad pixels
     17    assert (maskVal);
     18    assert (mark);
    1319
    1420    // run model fitting tests on a single source?
  • trunk/psphot/src/psphotReadout.c

    r16250 r16820  
    99    psTimerStart ("psphotReadout");
    1010
    11     bool dump = (psTraceGetLevel("psphot") >= 6);
    12 
    1311    // select the current recipe
    1412    psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE);
     
    1816    }
    1917
    20     // Interpret the mask values
    21     const char *maskValStr = psMetadataLookupStr(NULL, recipe, "MASKVAL"); // String with mask names
    22     if (!maskValStr) {
    23         psError(PSPHOT_ERR_CONFIG, false, "Missing recipe item: MASKVAL(STR)");
    24         return false;
    25     }
    26     psMaskType maskVal = pmConfigMask(maskValStr, config); // Mask values to mask against
    27     psMaskType maskMark = pmConfigMask("MARK", config); // Mask value for marking
    28     psMaskType maskSat = pmConfigMask("SAT", config); // Mask value for saturated pixels
    29     psMaskType maskBad = pmConfigMask("BAD", config); // Mask value for bad pixels
    30 
    31     // find the currently selected readout
    32     pmReadout  *readout = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT");
    33     PS_ASSERT_PTR_NON_NULL (readout, false);
    34 
    35     // optional break-point for processing
    36     char *breakPt = psMetadataLookupStr (NULL, recipe, "BREAK_POINT");
    37     PS_ASSERT_PTR_NON_NULL (breakPt, false);
    38 
    39     // Use the new pmFootprints approach?
    40     const bool useFootprints = psMetadataLookupBool(NULL, recipe, "USE_FOOTPRINTS");
    41 
    42     // generate mask & weight images if they don't already exit
    43     if (!readout->mask) {
    44         if (!pmReadoutGenerateMask(readout, maskSat, maskBad)) {
    45             psError (PSPHOT_ERR_CONFIG, false, "trouble creating mask");
    46             return false;
    47         }
    48     }
    49     if (!readout->weight) {
    50         if (!pmReadoutGenerateWeight(readout, true)) {
    51             psError (PSPHOT_ERR_CONFIG, false, "trouble creating weight");
    52             return false;
    53         }
    54     }
    55 
    5618    // set the photcode for this image
    5719    if (!psphotAddPhotcode (recipe, config, view)) {
     
    6022    }
    6123
    62     // I have a valid mask, now mask in the analysis region of interest
    63     psphotMaskReadout (readout, recipe, maskBad);
    64 
    65     if (psTraceGetLevel("psphot") >= 5) {
    66         psphotSaveImage (NULL, readout->image,  "image.fits");
    67         psphotSaveImage (NULL, readout->mask,   "mask.fits");
    68         psphotSaveImage (NULL, readout->weight, "weight.fits");
    69     }
    70 
     24    // find the currently selected readout
     25    pmReadout  *readout = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT");
     26    PS_ASSERT_PTR_NON_NULL (readout, false);
     27
     28    // optional break-point for processing
     29    char *breakPt = psMetadataLookupStr (NULL, recipe, "BREAK_POINT");
     30    PS_ASSERT_PTR_NON_NULL (breakPt, false);
     31
     32    // Generate the mask and weight images, including the user-defined analysis region of interest
     33    psphotSetMaskAndWeight (config, readout, recipe);
    7134    if (!strcasecmp (breakPt, "NOTHING")) {
    72         return psphotReadoutCleanup(config, readout, recipe, NULL, NULL);
     35        return psphotReadoutCleanup(config, readout, recipe, NULL, NULL, NULL);
    7336    }
    7437
    7538    // generate a background model (median, smoothed image)
    76     if (!psphotModelBackground (config, view, "PSPHOT.INPUT", maskVal)) {
    77         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL);
    78     }
    79     if (!psphotSubtractBackground (config, view, "PSPHOT.INPUT", maskVal)) {
    80         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL);
    81     }
    82 
     39    if (!psphotModelBackground (config, view, "PSPHOT.INPUT")) {
     40        return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     41    }
     42    if (!psphotSubtractBackground (config, view, "PSPHOT.INPUT")) {
     43        return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     44    }
    8345    if (!strcasecmp (breakPt, "BACKMDL")) {
    84         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL);
    85     }
    86 
    87     // run a single-model test if desired
    88     psphotModelTest (config, view, recipe, maskVal, maskMark);
     46        return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     47    }
     48
     49    // run a single-model test if desired (exits from here if test is run)
     50    psphotModelTest (config, view, recipe);
    8951
    9052    // load the psf model, if suppled.  FWHM_X,FWHM_Y,etc are saved in the recipe
     
    9456    bool havePSF = (psf != NULL);
    9557
    96     // find the peaks in the image.
    97 
    98     // XXX clean this up into a single function.  if psf is defined, we should treat this as
    99     // pass "2", not "1".  re-define this boolean to be "have PSF"
    100     psArray *peaks;
    101     psArray *footprints = NULL;
    102     if (useFootprints) {
    103        footprints = psphotFindPeaks (readout, recipe, useFootprints, 1, maskVal);
    104        int growRadius = psMetadataLookupS32(NULL, recipe, "FOOTPRINT_GROW_RADIUS");
    105        if (growRadius > 0) {
    106            psArray *tmp = pmGrowFootprintArray(footprints, growRadius);
    107            psFree(footprints);
    108            footprints = tmp;
    109        }
    110        psphotCullPeaks(readout->image, readout->weight, recipe, footprints);
    111 
    112        peaks = pmFootprintArrayToPeaks(footprints);
    113     } else {
    114        peaks = psphotFindPeaks (readout, recipe, useFootprints, 1, maskVal);
    115     }
    116 
    117     if (!peaks) {
    118         psLogMsg ("psphot", 3, "unable to find peaks in this image");
    119         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL);
     58    // find the detections (by peak and/or footprint) in the image.
     59    pmDetections *detections = psphotFindDetections (NULL, readout, recipe);
     60    if (!detections) {
     61        psLogMsg ("psphot", 3, "unable to find detections in this image");
     62        return psphotReadoutCleanup (config, readout, recipe, detections, psf, NULL);
    12063    }
    12164
    12265    // construct sources and measure basic stats
    123     psArray *sources = psphotSourceStats (readout, recipe, peaks, maskVal, maskMark);
     66    psArray *sources = psphotSourceStats (readout, recipe, detections);
    12467    if (!sources) return false;
    125     psFree (peaks);
    126 
    12768    if (!strcasecmp (breakPt, "PEAKS")) {
    128         return psphotReadoutCleanup(config, readout, recipe, NULL, sources);
    129     }
    130 
     69        return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
     70    }
     71
     72    // find blended neighbors of very saturated stars
     73    // XXX merge this with Basic Deblend?
    13174    psphotDeblendSatstars (sources, recipe);
    13275
     
    13477    if (!psphotBasicDeblend (sources, recipe)) {
    13578        psLogMsg ("psphot", 3, "failed on deblend analysis");
    136         return psphotReadoutCleanup (config, readout, recipe, NULL, sources);
     79        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
    13780    }
    13881
    13982    // classify sources based on moments, brightness
    140     if (!psphotRoughClass (sources, recipe, havePSF, maskSat)) {
     83    if (!psphotRoughClass (sources, recipe, havePSF)) {
    14184        psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
    142         return psphotReadoutCleanup (config, readout, recipe, NULL, sources);
     85        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
    14386    }
    14487    if (!strcasecmp (breakPt, "MOMENTS")) {
    145         return psphotReadoutCleanup(config, readout, recipe, NULL, sources);
     88        return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
    14689    }
    14790
     
    15194        // XXX if we do not have enough stars to generate the PSF, build one
    15295        // from the SEEING guess and model class
    153         psf = psphotChoosePSF (readout, sources, recipe, maskVal, maskMark);
     96        psf = psphotChoosePSF (readout, sources, recipe);
    15497        if (psf == NULL) {
    15598            psLogMsg ("psphot", 3, "failure to construct a psf model");
    156             return psphotReadoutCleanup (config, readout, recipe, psf, sources);
     99            return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
    157100        }
    158101        havePSF = true;
    159102    }
    160 
    161     // XXX Test dump of PSF parameters across image field.
    162 
    163103    if (!strcasecmp (breakPt, "PSFMODEL")) {
    164         return psphotReadoutCleanup (config, readout, recipe, psf, sources);
     104        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
    165105    }
    166106
     
    169109
    170110    // construct an initial model for each object
    171     psphotGuessModels (readout, sources, recipe, psf, maskVal);
     111    psphotGuessModels (readout, sources, recipe, psf);
    172112
    173113    // XXX test output of models
    174114    // psphotTestSourceOutput (readout, sources, recipe, psf);
    175115
    176     if (dump) psphotSaveImage (NULL, readout->image,  "image.v0.fits");
    177 
    178     // linear PSF fit to peaks
    179     psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE, maskVal);
    180     if (dump) psphotSaveImage (NULL, readout->image,  "image.v1.fits");
     116    // linear PSF fit to source peaks
     117    psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE);
    181118
    182119    // XXX test the CR/EXT measurement
    183120    // XXX we need to call this here and after the second pass: option for starting number?
    184     psphotSourceSize (readout, sources, recipe);
    185 
     121    psphotSourceSize (readout, sources, recipe, 0);
    186122    if (!strcasecmp (breakPt, "ENSEMBLE")) {
    187123        goto finish;
     
    189125
    190126    // non-linear PSF and EXT fit to brighter sources
    191     psphotBlendFit (readout, sources, recipe, psf, maskVal);
    192     if (dump) psphotSaveImage (NULL, readout->image,  "image.v2.fits");
    193 
    194     // replace all sources
    195     psphotReplaceAll (sources, maskVal);
    196     if (dump) psphotSaveImage (NULL, readout->image,  "image.v3.fits");
    197 
    198     // linear PSF fit to remaining peaks
    199     psphotFitSourcesLinear (readout, sources, recipe, psf, TRUE, maskVal);
    200     if (dump) psphotSaveImage (NULL, readout->image,  "image.v4.fits");
     127    psphotBlendFit (readout, sources, recipe, psf);
     128
     129    // replace all sources (make this part of psphotFitSourcesLinear?)
     130    psphotReplaceAll (sources, recipe);
     131
     132    // linear fit to include all sources
     133    psphotFitSourcesLinear (readout, sources, recipe, psf, TRUE);
    201134    if (!strcasecmp (breakPt, "PASS1")) {
    202135        goto finish;
     
    210143
    211144    // add noise for subtracted objects
    212     psphotAddNoise (readout, sources, recipe, true, maskVal);
    213 
    214     // find the peaks in the image
    215     psArray *newPeaks;
    216     if (useFootprints) {
    217        psArray *newFootprints = psphotFindPeaks (readout, recipe, useFootprints, 2, maskVal);
    218 
    219        int growRadius = psMetadataLookupS32(NULL, recipe, "FOOTPRINT_GROW_RADIUS_2");
    220        if (growRadius > 0) {
    221            psArray *tmp = pmGrowFootprintArray(newFootprints, growRadius);
    222            psFree(newFootprints);
    223            newFootprints = tmp;
    224        }
    225 
    226        // merge in old peaks
    227        const int includePeaks = 0x0 | 0x2; // i.e. just from newFootprints
    228        psArray *mergedFootprints = pmMergeFootprintArrays(footprints, newFootprints, includePeaks);
    229        psFree(footprints);
    230        psFree(newFootprints);
    231 
    232        psphotCullPeaks(readout->image, readout->weight, recipe, mergedFootprints);
    233 
    234        newPeaks = pmFootprintArrayToPeaks(mergedFootprints);
    235 
    236        psFree(mergedFootprints);
    237     } else {
    238        newPeaks = psphotFindPeaks(readout, recipe, useFootprints, 2, maskVal);
    239     }
     145    psphotAddNoise (readout, sources, recipe);
     146
     147    detections = psphotFindDetections (detections, readout, recipe);
    240148
    241149    // remove noise for subtracted objects
    242     psphotAddNoise (readout, sources, recipe, false, maskVal);
    243 
    244     // define new sources based on the new peaks
    245     psArray *newSources = psphotSourceStats (readout, recipe, newPeaks, maskVal, maskMark);
    246     psFree (newPeaks);
     150    psphotSubNoise (readout, sources, recipe);
     151
     152    // define new sources based on only the new peaks
     153    psArray *newSources = psphotSourceStats (readout, recipe, detections);
    247154
    248155    // set source type
    249     if (!psphotRoughClass (newSources, recipe, havePSF, maskSat)) {
     156    if (!psphotRoughClass (newSources, recipe, havePSF)) {
    250157        psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
    251         return psphotReadoutCleanup (config, readout, recipe, psf, sources);
     158        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
    252159    }
    253160
    254161    // create full input models
    255     psphotGuessModels (readout, newSources, recipe, psf, maskVal);
    256 
    257     // replace all sources
    258     psphotReplaceAll (sources, maskVal);
    259     if (dump) psphotSaveImage (NULL, readout->image,  "image.v5.fits");
    260 
    261     // merge the newly selected peaks into the existing list
     162    psphotGuessModels (readout, newSources, recipe, psf);
     163
     164    // replace all sources so fit below applies to all at once
     165    psphotReplaceAll (sources, recipe);
     166
     167    // merge the newly selected sources into the existing list
    262168    psphotMergeSources (sources, newSources);
    263169    psFree (newSources);
    264170
    265     // linear PSF fit to remaining peaks
    266     psphotFitSourcesLinear (readout, sources, recipe, psf, TRUE, maskVal);
    267     if (dump) psphotSaveImage (NULL, readout->image,  "image.v6.fits");
    268 
    269     psphotExtendedSources (readout, sources, recipe, maskVal);
     171    // linear fit to all sources
     172    psphotFitSourcesLinear (readout, sources, recipe, psf, TRUE);
     173
     174    // measure source size for the remaining sources
     175    psphotSourceSize (readout, sources, recipe, 0);
     176
     177    psphotExtendedSources (readout, sources, recipe);
    270178
    271179finish:
     
    273181    // plot positive sources
    274182    if (!havePSF) {
    275         // psphotSourcePlots (readout, sources, recipe, maskVal);
     183        // psphotSourcePlots (readout, sources, recipe);
    276184    }
    277185
    278186    // measure aperture photometry corrections
    279     if (!havePSF && !psphotApResid (readout, sources, recipe, psf, maskVal, maskMark)) {
    280         psTrace ("psphot", 4, "failure on psphotApResid");
    281         psError(PSPHOT_ERR_PHOTOM, false, "Measure aperture photometry corrections");
    282         return false;
     187    if (!havePSF && !psphotApResid (readout, sources, recipe, psf)) {
     188        psLogMsg ("psphot", 3, "failed on psphotApResid");
     189        return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
    283190    }
    284191
    285192    // calculate source magnitudes
    286193    pmReadout *background = psphotSelectBackground (config, view, false);
    287     psphotMagnitudes(sources, recipe, psf, background, maskVal, maskMark);
     194    psphotMagnitudes(sources, recipe, psf, background);
    288195
    289196    // replace failed sources?
     
    297204
    298205    // create the exported-metadata and free local data
    299     return psphotReadoutCleanup(config, readout, recipe, psf, sources);
     206    return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
    300207}
  • trunk/psphot/src/psphotReadoutCleanup.c

    r14655 r16820  
    11# include "psphotInternal.h"
    22
    3 // psphotReadoutCleanup is called on exit from psphotReadout if the last raised
    4 // error is not a DATA error, then there was a serious problem.  only in this
    5 // case, or if the fail on the stats measurement, do we return false
    6 bool psphotReadoutCleanup (pmConfig *config, pmReadout *readout, psMetadata *recipe, pmPSF *psf, psArray *sources) {
     3// psphotReadoutCleanup is called on exit from psphotReadout.  If the last raised error is
     4// not a DATA error, then there was a serious problem.  Only in this case, or if the fail
     5// on the stats measurement, do we return false
     6bool psphotReadoutCleanup (pmConfig *config, pmReadout *readout, psMetadata *recipe, pmDetections *detections, pmPSF *psf, psArray *sources) {
    77
    88    // remove internal pmFPAfiles, if created
     
    7373    pmKapaClose ();
    7474
     75    psFree (detections);
    7576    psFree (psf);
    7677    psFree (header);
  • trunk/psphot/src/psphotReplaceUnfit.c

    r13900 r16820  
    1818        pmSourceAdd (source, PM_MODEL_OP_FULL, maskVal);
    1919        source->mode &= ~PM_SOURCE_MODE_SUBTRACTED;
    20         source->mode &= ~PM_SOURCE_MODE_TEMPSUB;
    2120    }
    2221    psLogMsg ("psphot.replace", 3, "replace unfitted models: %f sec (%ld objects)\n", psTimerMark ("psphot"), sources->n);
     
    2423}
    2524
    26 bool psphotReplaceAll (psArray *sources, psMaskType maskVal) {
     25bool psphotReplaceAll (psArray *sources, psMetadata *recipe) {
    2726
     27    bool status;
    2828    pmSource *source;
    2929
    3030    psTimerStart ("psphot");
     31
     32    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     33    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     34    assert (maskVal);
    3135
    3236    for (int i = 0; i < sources->n; i++) {
     
    3943      source->mode &= ~PM_SOURCE_MODE_SUBTRACTED;
    4044    }
    41     psLogMsg ("psphot.replace", PS_LOG_INFO, "replace models for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot"));
     45    psLogMsg ("psphot.replace", PS_LOG_INFO, "replaced models for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot"));
    4246    return true;
    4347}
  • trunk/psphot/src/psphotRoughClass.c

    r15059 r16820  
    11# include "psphotInternal.h"
     2
    23# define CHECK_STATUS(S,MSG) { \
    34  if (!status) { \
     
    78
    89// 2006.02.02 : no leaks
    9 bool psphotRoughClass (psArray *sources, psMetadata *recipe, const bool havePSF, psMaskType maskSat) {
     10bool psphotRoughClass (psArray *sources, psMetadata *recipe, const bool havePSF) {
    1011
    1112    bool status;
     
    1314
    1415    psTimerStart ("psphot");
     16
     17    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     18    psMaskType maskSat = psMetadataLookupU8(&status, recipe, "MASK.SAT"); // Mask value for bad pixels
     19    assert (maskSat);
    1520
    1621    if (!havePSF) {
  • trunk/psphot/src/psphotSourceFits.c

    r15058 r16820  
    113113        pmSourceSub (blend, PM_MODEL_OP_FULL, maskVal);
    114114        blend->mode |=  PM_SOURCE_MODE_SUBTRACTED;
    115         blend->mode &= ~PM_SOURCE_MODE_TEMPSUB;
    116115    }
    117116    NfitBlend += modelSet->n;
     
    137136    pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
    138137    source->mode |=  PM_SOURCE_MODE_SUBTRACTED;
    139     source->mode &= ~PM_SOURCE_MODE_TEMPSUB;
    140138    return true;
    141139}
     
    176174    pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
    177175    source->mode |=  PM_SOURCE_MODE_SUBTRACTED;
    178     source->mode &= ~PM_SOURCE_MODE_TEMPSUB;
    179 
    180176    return true;
    181177}
     
    275271    pmSourceCacheModel (source, maskVal);
    276272    pmSourceSub (source, PM_MODEL_OP_FULL, maskVal);
     273    source->mode |=  PM_SOURCE_MODE_SUBTRACTED;
    277274    psTrace ("psphot", 5, "blob as EXT: %f %f\n", EXT->params->data.F32[PM_PAR_XPOS], EXT->params->data.F32[PM_PAR_YPOS]);
    278 
    279     source->mode |=  PM_SOURCE_MODE_SUBTRACTED;
    280     source->mode &= ~PM_SOURCE_MODE_TEMPSUB;
    281275    return true;
    282276
     
    290284    source->mode    |= PM_SOURCE_MODE_SUBTRACTED;
    291285    source->mode    |= PM_SOURCE_MODE_PAIR;
    292     source->mode    &= ~PM_SOURCE_MODE_TEMPSUB;
    293286
    294287    // copy most data from the primary source (modelEXT, blends stay NULL)
  • trunk/psphot/src/psphotSourcePlots.c

    r15000 r16820  
    11# include "psphotInternal.h"
    22
    3 bool psphotSourcePlots (pmReadout *readout, psArray *sources, psMetadata *recipe, psMaskType maskVal) {
     3bool psphotSourcePlots (pmReadout *readout, psArray *sources, psMetadata *recipe) {
    44
    5     // bool status = false;
     5    bool status = false;
    66    psTimerStart ("psphot");
     7
     8    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     9    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     10    assert (maskVal);
    711
    812    // the source images are written to an image 10x the size of a PSF object
  • trunk/psphot/src/psphotSourceSize.c

    r15802 r16820  
    66// compares the observed flux to the model.
    77
    8 bool psphotSourceSize (pmReadout *readout, psArray *sources, psMetadata *recipe) {
     8bool psphotSourceSize (pmReadout *readout, psArray *sources, psMetadata *recipe, long first) {
     9
     10    bool status;
     11
     12    psTimerStart ("psphot");
     13
     14    float CR_NSIGMA_LIMIT = psMetadataLookupF32 (&status, recipe, "PSPHOT.CRNSIGMA.LIMIT");
     15    assert (status);
    916
    1017    // loop over all source
    11     for (int i = 0; i < sources->n; i++) {
     18    for (int i = first; i < sources->n; i++) {
    1219        pmSource *source = sources->data[i];
    1320
     21        // skip source if it was already measured
     22        if (isfinite(source->crNsigma)) continue;
     23
     24        source->crNsigma  = -1.0;
     25        source->extNsigma = -1.0;
     26
    1427        // source must have been subtracted
     28        source->crNsigma  = -3.0;
    1529        if (!(source->mode & PM_SOURCE_MODE_SUBTRACTED)) continue;
    1630
     
    2337
    2438        // skip sources which are too close to a boundary
     39        source->crNsigma  = -4.0;
    2540        if (xPeak < 1) continue;
    2641        if (xPeak > source->pixels->numCols - 2) continue;
     
    2944
    3045        // XXX for now, just skip any sources with masked pixels
     46        source->crNsigma  = -5.0;
    3147        bool keep = true;
    3248        for (int iy = -1; (iy <= +1) && keep; iy++) {
     
    97113        source->extNsigma = (nEXT > 0) ? fabs(fEXT) / nEXT : 0.0;
    98114        // NOTE: abs needed to make the Nsigma value positive
     115
     116        if (!isfinite(source->crNsigma)) {
     117          fprintf (stderr, ".");
     118          source->crNsigma = -6.0;
     119        }
     120
     121        if (source->crNsigma > CR_NSIGMA_LIMIT) {
     122          source->mode |= PM_SOURCE_MODE_CRLIMIT;
     123        }
    99124    }
     125
     126    psLogMsg ("psphot.size", PS_LOG_INFO, "measure source sizes for %ld sources: %f sec\n", sources->n - first, psTimerMark ("psphot"));
     127
    100128    return true;
    101129}
  • trunk/psphot/src/psphotSourceStats.c

    r15405 r16820  
    11# include "psphotInternal.h"
    22
    3 psArray *psphotSourceStats (pmReadout *readout, psMetadata *recipe, psArray *peaks, psMaskType maskVal, psMaskType mark)
     3psArray *psphotSourceStats (pmReadout *readout, psMetadata *recipe, pmDetections *detections)
    44{
    55    bool     status  = false;
     
    88
    99    psTimerStart ("psphot");
     10
     11    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     12    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     13    assert (maskVal);
     14
     15    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     16    psMaskType mark = psMetadataLookupU8(&status, recipe, "MASK.MARK"); // Mask value for bad pixels
     17    assert (mark);
    1018
    1119    // determine properties (sky, moments) of initial sources
     
    2028    char *breakPt  = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
    2129    if (!status) return NULL;
     30
     31    psArray *peaks = detections->peaks;
    2232
    2333    sources = psArrayAllocEmpty (peaks->n);
  • trunk/psphot/src/psphotSubtractBackground.c

    r15936 r16820  
    44// generate the median in NxN boxes, clipping heavily
    55// linear interpolation to generate full-scale model
    6 bool psphotSubtractBackground (pmConfig *config, const pmFPAview *view, const char *filename, psMaskType maskVal)
     6bool psphotSubtractBackground (pmConfig *config, const pmFPAview *view, const char *filename)
    77{
    88    bool status = true;
     
    2525    // select the appropriate recipe information
    2626    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
     27    assert (recipe);
     28
     29    // user-defined masks to test for good/bad pixels (build from recipe list if not yet set)
     30    psMaskType maskVal = psMetadataLookupU8(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels
     31    assert (maskVal);
    2732
    2833    psImageBinning *binning = psMetadataLookupPtr(&status, recipe, "PSPHOT.BACKGROUND.BINNING");
Note: See TracChangeset for help on using the changeset viewer.