IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 26748


Ignore:
Timestamp:
Jan 31, 2010, 5:01:05 PM (16 years ago)
Author:
eugene
Message:

updates relative to 20091201, fixes for all psphot variants

Location:
branches/eam_branches/psphot.stack.20100120
Files:
1 deleted
31 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/psphot.stack.20100120

  • branches/eam_branches/psphot.stack.20100120/doc/footprints.txt

    r17710 r26748  
     1
     22010.01.18
     3
     4  * we are having some problems with footprints being freed to early in the second pass
     5
     6  * footprint memory history:
     7
     8    * entering psphotFindFootprints, on the first pass, we have an
     9      array of detected peak in detections->peaks
     10
     11    * in pmFootprintsFind, footprints array owns allocated footprints
     12      footprints each own their spans
     13
     14    * from pmFootprintAssignPeaks:
     15      * the footprints have an array of peaks with valid memory references
     16      * the peaks have just the pointer to the corresponding footprint (not owned)
     17      * all peaks get assigned to a footprint
     18      * at end of pmFootprintAssignPeaks, only a single copy of each peak is/should be available from the collection of footprints
    119
    220pmFootprintCullPeaks is very expensive (15 - 30 msec per object) and

  • branches/eam_branches/psphot.stack.20100120/doc/stack.txt

    r26688 r26748  
     1
     220100126:
     3
     4  * watch out for psphotSetMomentsWindow & MOMENTS_SX_MAX,etc
     5  * watch out for psphotSignificanceImage.c:,psphotEfficiency using the FWHM_MAJ from psphotChoosePSF
     6    * ppSimDetections.c : SIGMA_SMOOTH
     7ppSmooth/src/ppSmoothReadout.c:    psMetadataAddF32(recipe, PS_LIST_TAIL, "EFFECTIVE_AREA", PS_META_REPLACE, "Effective Area", effArea);
     8ppSmooth/src/ppSmoothReadout.c:    psMetadataAddF32(recipe, PS_LIST_TAIL, "SIGNIFICANCE_SCALE_FACTOR", PS_META_REPLACE, "Signicance scale factor", factor);
     9
    110
    21120100120 : more stack processing mods:

  • branches/eam_branches/psphot.stack.20100120/src/Makefile.am

    r26691 r26748  
    2525libpsphot_la_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
    2626
    27 bin_PROGRAMS = psphot
    28 # bin_PROGRAMS = psphotPetrosianStudy psphotForced psphotMakePSF psphotTest psphotMomentsStudy
     27bin_PROGRAMS = psphot psphotForced psphotMakePSF
     28# bin_PROGRAMS = psphotPetrosianStudy psphotTest psphotMomentsStudy
    2929
    3030psphot_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
    … …  
    3232psphot_LDADD = libpsphot.la
    3333
    34 # psphotForced_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
    35 # psphotForced_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
    36 # psphotForced_LDADD = libpsphot.la
    37 #
    38 # psphotMakePSF_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
    39 # psphotMakePSF_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
    40 # psphotMakePSF_LDADD = libpsphot.la
    41 #
     34psphotForced_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
     35psphotForced_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
     36psphotForced_LDADD = libpsphot.la
     37
     38psphotMakePSF_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
     39psphotMakePSF_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
     40psphotMakePSF_LDADD = libpsphot.la
     41
     42# psphotMomentsStudy_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
     43# psphotMomentsStudy_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
     44# psphotMomentsStudy_LDADD = libpsphot.la
     45
     46# psphotPetrosianStudy_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
     47# psphotPetrosianStudy_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
     48# psphotPetrosianStudy_LDADD = libpsphot.la
     49
    4250# psphotTest_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
    4351# psphotTest_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
    4452# psphotTest_LDADD = libpsphot.la
    45 #
    46 # psphotMomentsStudy_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
    47 # psphotMomentsStudy_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
    48 # psphotMomentsStudy_LDADD = libpsphot.la
    49 
    50 # psphotPetrosianStudy_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
    51 # psphotPetrosianStudy_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
    52 # psphotPetrosianStudy_LDADD = libpsphot.la
    5353
    5454# standard psphot for generic photometry
    … …  
    6161        psphotCleanup.c
    6262
    63 # # forced photometry of specified positions given a specified psf
    64 # psphotForced_SOURCES = \
    65 #         psphotForced.c             \
    66 #       psphotForcedArguments.c    \
    67 #       psphotForcedImageLoop.c    \
    68 #       psphotForcedReadout.c      \
    69 #       psphotParseCamera.c        \
    70 #       psphotMosaicChip.c         \
    71 #       psphotCleanup.c
    72 #
    73 # # a psphot-variant that simply generates the PSF model
    74 # psphotMakePSF_SOURCES = \
    75 #         psphotMakePSF.c            \
    76 #       psphotMakePSFArguments.c   \
    77 #       psphotMakePSFImageLoop.c   \
    78 #       psphotMakePSFReadout.c     \
    79 #       psphotParseCamera.c        \
    80 #       psphotMosaicChip.c         \
    81 #       psphotCleanup.c
    82 #
    83 #
     63# forced photometry of specified positions given a specified psf
     64psphotForced_SOURCES = \
     65        psphotForced.c             \
     66        psphotForcedArguments.c    \
     67        psphotForcedImageLoop.c    \
     68        psphotForcedReadout.c      \
     69        psphotParseCamera.c        \
     70        psphotMosaicChip.c         \
     71        psphotCleanup.c
     72
     73# a psphot-variant that simply generates the PSF model
     74psphotMakePSF_SOURCES = \
     75        psphotMakePSF.c            \
     76        psphotMakePSFArguments.c   \
     77        psphotMakePSFImageLoop.c   \
     78        psphotMakePSFReadout.c     \
     79        psphotParseCamera.c        \
     80        psphotMosaicChip.c         \
     81        psphotCleanup.c
     82
    8483# # psphot analysis of the detectability of specified positions
    8584# psphotDetect_SOURCES =            \
    … …  
    9190#       psphotMosaicChip.c        \
    9291#       psphotCleanup.c
    93 #
     92
    9493# psphotTest_SOURCES = \
    9594#         psphotTest.c
    … …  
    9796# psphotMomentsStudy_SOURCES = \
    9897#         psphotMomentsStudy.c
    99 
     98#
    10099# psphotPetrosianStudy_SOURCES = \
    101100#         psphotPetrosianStudy.c
    … …  
    112111        psphotDefineFiles.c            \
    113112        psphotReadout.c                \
     113        psphotReadoutFindPSF.c         \
     114        psphotReadoutKnownSources.c    \
     115        psphotReadoutMinimal.c         \
    114116        psphotModelBackground.c        \
    115117        psphotSubtractBackground.c     \
    … …  
    170172        psphotEfficiency.c
    171173
    172 # XXX need to fix these for the new apis
    173 #       psphotReadoutFindPSF.c         
    174 #       psphotReadoutKnownSources.c   
    175 #       psphotReadoutMinimal.c         
     174# XXX need to fix this for the new apis
    176175#       psphotModelTest.c             
    177176

  • branches/eam_branches/psphot.stack.20100120/src/psphot.h

    r26691 r26748  
    6363bool            psphotSubtractBackgroundReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index, psMetadata *recipe);
    6464
    65 bool            psphotFindDetections (pmConfig *config, const pmFPAview *view);
    66 bool            psphotFindDetectionsReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index, psMetadata *recipe);
     65bool            psphotFindDetections (pmConfig *config, const pmFPAview *view, bool firstPass);
     66bool            psphotFindDetectionsReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index, psMetadata *recipe, bool firstPass);
    6767
    6868bool            psphotSourceStats (pmConfig *config, const pmFPAview *view, bool setWindow);
    … …  
    140140
    141141// in psphotChoosePSF.c:
    142 bool            psphotPSFstats (pmReadout *readout, psMetadata *recipe, pmPSF *psf);
    143 bool            psphotMomentsStats (pmReadout *readout, psMetadata *recipe, psArray *sources);
     142bool            psphotPSFstats (pmReadout *readout, pmPSF *psf);
     143bool            psphotMomentsStats (pmReadout *readout, psArray *sources);
    144144
    145145// in psphotGuessModel.c
    … …  
    147147
    148148// in psphotMergeSources.c:
    149 bool            psphotLoadExtSources (pmConfig *config, const pmFPAview *view, psArray *sources);
     149bool            psphotLoadExtSources (pmConfig *config, const pmFPAview *view);
    150150psArray        *psphotLoadPSFSources (pmConfig *config, const pmFPAview *view);
     151bool            psphotRepairLoadedSources (pmConfig *config, const pmFPAview *view);
     152bool            psphotCheckExtSources (pmConfig *config, const pmFPAview *view);
    151153
    152154// generate the detection structure for the supplied array of sources
    153 pmDetections   *psphotDetectionsFromSources (pmConfig *config, psArray *sources);
     155bool            psphotDetectionsFromSources (pmConfig *config, const pmFPAview *view, psArray *sources);
    154156
    155157// generate the detection structure for the supplied array of sources
    … …  
    159161// Create a background model for a readout, without saving the result as a pmFPAfile on config->files.  Otherwise identical to psphotModelBackgroundFileIndex.
    160162psImage        *psphotModelBackgroundReadoutNoFile(pmReadout *readout, const pmConfig *config);
    161 psImageBinning *psphotBackgroundBinning(const psImage *, const pmConfig *);
     163psImageBinning *psphotBackgroundBinning(const psImage *image, const pmConfig *config);
    162164
    163165// in psphotReplaceUnfit.c:
    … …  
    186188
    187189// functions to set the correct source pixels
    188 bool            psphotInitRadiusPSF (const psMetadata *recipe, const pmModelType type);
     190bool            psphotInitRadiusPSF(const psMetadata *recipe, const psMetadata *analysis, const pmModelType type);
     191
    189192bool            psphotCheckRadiusPSF (pmReadout *readout, pmSource *source, pmModel *model, psImageMaskType markVal);
    190193bool            psphotCheckRadiusPSFBlend (pmReadout *readout, pmSource *source, pmModel *model, psImageMaskType markVal, float dR);
    … …  
    224227
    225228bool            psphotLoadPSF (pmConfig *config, const pmFPAview *view);
    226 bool            psphotLoadPSFReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index, psMetadata *recipe);
     229bool            psphotLoadPSFReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index);
    227230
    228231bool            psphotSetHeaderNstars (psMetadata *recipe, psArray *sources);

  • branches/eam_branches/psphot.stack.20100120/src/psphotApResid.c

    r26691 r26748  
    6767    if (!measureAptrend) {
    6868        // save nan values since these were not calculated
    69         psMetadataAdd (recipe, PS_LIST_TAIL, "APMIFIT",  PS_DATA_F32 | PS_META_REPLACE, "aperture residual",   NAN);
    70         psMetadataAdd (recipe, PS_LIST_TAIL, "DAPMIFIT", PS_DATA_F32 | PS_META_REPLACE, "ap residual scatter", NAN);
    71         psMetadataAdd (recipe, PS_LIST_TAIL, "NAPMIFIT", PS_DATA_S32 | PS_META_REPLACE, "number of apresid stars", 0);
    72         psMetadataAdd (recipe, PS_LIST_TAIL, "APLOSS",   PS_DATA_F32 | PS_META_REPLACE, "aperture loss (mag)", NAN);
     69        psMetadataAdd (readout->analysis, PS_LIST_TAIL, "APMIFIT",  PS_DATA_F32 | PS_META_REPLACE, "aperture residual",   NAN);
     70        psMetadataAdd (readout->analysis, PS_LIST_TAIL, "DAPMIFIT", PS_DATA_F32 | PS_META_REPLACE, "ap residual scatter", NAN);
     71        psMetadataAdd (readout->analysis, PS_LIST_TAIL, "NAPMIFIT", PS_DATA_S32 | PS_META_REPLACE, "number of apresid stars", 0);
     72        psMetadataAdd (readout->analysis, PS_LIST_TAIL, "APLOSS",   PS_DATA_F32 | PS_META_REPLACE, "aperture loss (mag)", NAN);
    7373        return true;
    7474    }
    … …  
    322322
    323323    // save results for later output
    324     psMetadataAdd (recipe, PS_LIST_TAIL, "APMIFIT",  PS_DATA_F32 | PS_META_REPLACE, "aperture residual",   psf->ApResid);
    325     psMetadataAdd (recipe, PS_LIST_TAIL, "DAPMIFIT", PS_DATA_F32 | PS_META_REPLACE, "ap residual scatter", psf->dApResid);
    326     psMetadataAdd (recipe, PS_LIST_TAIL, "NAPMIFIT", PS_DATA_S32 | PS_META_REPLACE, "number of apresid stars", psf->nApResid);
    327     psMetadataAdd (recipe, PS_LIST_TAIL, "APLOSS",   PS_DATA_F32 | PS_META_REPLACE, "aperture loss (mag)", psf->growth->apLoss);
     324    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "APMIFIT",  PS_DATA_F32 | PS_META_REPLACE, "aperture residual",   psf->ApResid);
     325    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "DAPMIFIT", PS_DATA_F32 | PS_META_REPLACE, "ap residual scatter", psf->dApResid);
     326    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "NAPMIFIT", PS_DATA_S32 | PS_META_REPLACE, "number of apresid stars", psf->nApResid);
     327    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "APLOSS",   PS_DATA_F32 | PS_META_REPLACE, "aperture loss (mag)", psf->growth->apLoss);
    328328
    329329    psLogMsg ("psphot.apresid", PS_LOG_DETAIL, "aperture residual: %f +/- %f\n", psf->ApResid, psf->dApResid);
    … …  
    342342escape:
    343343    // save nan values since these were not calculated
    344     psMetadataAdd (recipe, PS_LIST_TAIL, "APMIFIT",  PS_DATA_F32 | PS_META_REPLACE, "aperture residual",   NAN);
    345     psMetadataAdd (recipe, PS_LIST_TAIL, "DAPMIFIT", PS_DATA_F32 | PS_META_REPLACE, "ap residual scatter", NAN);
    346     psMetadataAdd (recipe, PS_LIST_TAIL, "NAPMIFIT", PS_DATA_S32 | PS_META_REPLACE, "number of apresid stars", 0);
    347     psMetadataAdd (recipe, PS_LIST_TAIL, "APLOSS",   PS_DATA_F32 | PS_META_REPLACE, "aperture loss (mag)", NAN);
     344    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "APMIFIT",  PS_DATA_F32 | PS_META_REPLACE, "aperture residual",   NAN);
     345    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "DAPMIFIT", PS_DATA_F32 | PS_META_REPLACE, "ap residual scatter", NAN);
     346    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "NAPMIFIT", PS_DATA_S32 | PS_META_REPLACE, "number of apresid stars", 0);
     347    psMetadataAdd (readout->analysis, PS_LIST_TAIL, "APLOSS",   PS_DATA_F32 | PS_META_REPLACE, "aperture loss (mag)", NAN);
    348348
    349349    psFree (xPos);

  • branches/eam_branches/psphot.stack.20100120/src/psphotBlendFit.c

    r26691 r26748  
    7979    psphotInitLimitsPSF (recipe, readout);
    8080    psphotInitLimitsEXT (recipe);
    81     psphotInitRadiusPSF (recipe, psf->type);
     81    psphotInitRadiusPSF (recipe, readout->analysis, psf->type);
    8282
    8383    // starts the timer, sets up the array of fitSets

  • branches/eam_branches/psphot.stack.20100120/src/psphotChoosePSF.c

    r26691 r26748  
    107107    // assert (status);
    108108
    109     // We have calculated a Gaussian window function, use that for both the PSF fit radius and
    110     // the aperture radius (scaling SIGMA)
    111     float gaussSigma = psMetadataLookupF32(&status, recipe, "MOMENTS_GAUSS_SIGMA");
     109    // values on readout->analysis if we have calculated a Gaussian window function, use that
     110    // for both the PSF fit radius and the aperture radius (scaling SIGMA), otherwise base the value on the recipe value for MOMENTS_GAUSS_SIGMA
     111    float gaussSigma = psMetadataLookupF32 (&status, readout->analysis, "MOMENTS_GAUSS_SIGMA");
     112    if (!status) {
     113        gaussSigma = psMetadataLookupF32 (&status, recipe, "MOMENTS_GAUSS_SIGMA");
     114    }
    112115    float fitScale = psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS_SCALE");
    113116    float apScale = psMetadataLookupF32(&status, recipe, "PSF_APERTURE_SCALE");
    … …  
    116119
    117120    // XXX use the same radii for standard analysis as for the PSF creation
    118     psMetadataAddF32(recipe, PS_LIST_TAIL, "PSF_FIT_RADIUS", PS_META_REPLACE, "fit radius", options->fitRadius);
    119     psMetadataAddF32(recipe, PS_LIST_TAIL, "PSF_APERTURE", PS_META_REPLACE, "psf aperture", options->apRadius);
     121    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "PSF_FIT_RADIUS", PS_META_REPLACE, "fit radius", options->fitRadius);
     122    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "PSF_APERTURE", PS_META_REPLACE, "psf aperture", options->apRadius);
    120123
    121124    // dimensions of the field for which the PSF is defined
    … …  
    204207        bool status = true;
    205208        status &= psphotMakeFluxScale (readout->image, recipe, psf);
    206         status &= psphotPSFstats (readout, recipe, psf);
     209        status &= psphotPSFstats (readout, psf);
    207210        if (!status) {
    208211            psError(PSPHOT_ERR_PSF, false, "Failed to fit any models when choosing PSF");
    … …  
    279282        bool status = true;
    280283        status &= psphotMakeFluxScale (readout->image, recipe, psf);
    281         status &= psphotPSFstats (readout, recipe, psf);
     284        status &= psphotPSFstats (readout, psf);
    282285        if (!status) {
    283286            psError(PSPHOT_ERR_PSF, false, "Failed to fit any models when choosing PSF");
    … …  
    336339    psFree (models);
    337340
    338     if (!psphotPSFstats (readout, recipe, psf)) {
     341    if (!psphotPSFstats (readout, psf)) {
    339342        psError(PSPHOT_ERR_PSF, false, "cannot measure PSF shape terms");
    340343        psFree(options);
    … …  
    363366
    364367// measure average parameters of the PSF model
    365 bool psphotPSFstats (pmReadout *readout, psMetadata *recipe, pmPSF *psf) {
     368bool psphotPSFstats (pmReadout *readout, pmPSF *psf) {
    366369
    367370    psEllipseShape shape;
    … …  
    369372
    370373    PS_ASSERT_PTR_NON_NULL(readout, false);
    371     PS_ASSERT_PTR_NON_NULL(recipe, false);
    372374    PS_ASSERT_PTR_NON_NULL(psf, false);
    373375
    … …  
    419421    }
    420422
    421     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FWHM_MAJ",   PS_META_REPLACE, "PSF FWHM Major axis (mean)", stats->sampleMean);
    422     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MJ_SG",   PS_META_REPLACE, "PSF FWHM Major axis (sigma)", stats->sampleStdev);
    423     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MJ_LQ",   PS_META_REPLACE, "PSF FWHM Major axis (lower quartile)", stats->sampleLQ);
    424     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MJ_UQ",   PS_META_REPLACE, "PSF FWHM Major axis (upper quartile)", stats->sampleUQ);
     423    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FWHM_MAJ",   PS_META_REPLACE, "PSF FWHM Major axis (mean)", stats->sampleMean);
     424    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_SG",   PS_META_REPLACE, "PSF FWHM Major axis (sigma)", stats->sampleStdev);
     425    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_LQ",   PS_META_REPLACE, "PSF FWHM Major axis (lower quartile)", stats->sampleLQ);
     426    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_UQ",   PS_META_REPLACE, "PSF FWHM Major axis (upper quartile)", stats->sampleUQ);
    425427
    426428    if (!psVectorStats (stats, fwhmMinor, NULL, NULL, 0)) {
    … …  
    428430        return false;
    429431    }
    430     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FWHM_MIN",   PS_META_REPLACE, "PSF FWHM Minor axis (mean)", stats->sampleMean);
    431     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MN_SG",   PS_META_REPLACE, "PSF FWHM Minor axis (sigma)", stats->sampleStdev);
    432     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MN_LQ",   PS_META_REPLACE, "PSF FWHM Minor axis (lower quartile)", stats->sampleLQ);
    433     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MN_UQ",   PS_META_REPLACE, "PSF FWHM Minor axis (upper quartile)", stats->sampleUQ);
    434 
    435     psMetadataAddF32 (recipe, PS_LIST_TAIL, "ANGLE",    PS_META_REPLACE, "PSF angle",           axes.theta);
    436     psMetadataAddS32 (recipe, PS_LIST_TAIL, "NPSFSTAR", PS_META_REPLACE, "Number of stars used to make PSF", psf->nPSFstars);
    437     psMetadataAddBool(recipe, PS_LIST_TAIL, "PSFMODEL", PS_META_REPLACE, "Valid PSF Model?", true);
     432    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FWHM_MIN",   PS_META_REPLACE, "PSF FWHM Minor axis (mean)", stats->sampleMean);
     433    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_SG",   PS_META_REPLACE, "PSF FWHM Minor axis (sigma)", stats->sampleStdev);
     434    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_LQ",   PS_META_REPLACE, "PSF FWHM Minor axis (lower quartile)", stats->sampleLQ);
     435    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_UQ",   PS_META_REPLACE, "PSF FWHM Minor axis (upper quartile)", stats->sampleUQ);
     436
     437    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "ANGLE",    PS_META_REPLACE, "PSF angle",           axes.theta);
     438    psMetadataAddS32 (readout->analysis, PS_LIST_TAIL, "NPSFSTAR", PS_META_REPLACE, "Number of stars used to make PSF", psf->nPSFstars);
     439    psMetadataAddBool(readout->analysis, PS_LIST_TAIL, "PSFMODEL", PS_META_REPLACE, "Valid PSF Model?", true);
    438440
    439441    psFree (fwhmMajor);
    … …  
    444446
    445447// determine approximate PSF shape parameters based on the moments
    446 bool psphotMomentsStats (pmReadout *readout, psMetadata *recipe, psArray *sources) {
     448bool psphotMomentsStats (pmReadout *readout, psArray *sources) {
    447449
    448450    // without the PSF model, we can only rely on the source->moments
    … …  
    458460
    459461    PS_ASSERT_PTR_NON_NULL(readout, false);
    460     PS_ASSERT_PTR_NON_NULL(recipe, false);
    461462    PS_ASSERT_PTR_NON_NULL(sources, false);
    462463
    … …  
    486487    FWHM_T /= FWHM_N;
    487488
    488     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FWHM_MAJ",   PS_META_REPLACE, "PSF FWHM Major axis (mean)", FWHM_X);
    489     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MJ_SG",   PS_META_REPLACE, "PSF FWHM Major axis (sigma)", 0);
    490     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MJ_LQ",   PS_META_REPLACE, "PSF FWHM Major axis (lower quartile)", 0);
    491     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MJ_UQ",   PS_META_REPLACE, "PSF FWHM Major axis (upper quartile)", 0);
    492 
    493     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FWHM_MIN",   PS_META_REPLACE, "PSF FWHM Minor axis (mean)", FWHM_Y);
    494     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MN_SG",   PS_META_REPLACE, "PSF FWHM Minor axis (sigma)", 0);
    495     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MN_LQ",   PS_META_REPLACE, "PSF FWHM Minor axis (lower quartile)", 0);
    496     psMetadataAddF32 (recipe, PS_LIST_TAIL, "FW_MN_UQ",   PS_META_REPLACE, "PSF FWHM Minor axis (upper quartile)", 0);
    497 
    498     psMetadataAddF32 (recipe, PS_LIST_TAIL, "ANGLE",    PS_META_REPLACE, "PSF angle",           FWHM_T);
    499     psMetadataAddS32 (recipe, PS_LIST_TAIL, "NPSFSTAR", PS_META_REPLACE, "Number of stars used to make PSF", 0);
    500     psMetadataAddBool(recipe, PS_LIST_TAIL, "PSFMODEL", PS_META_REPLACE, "Valid PSF Model?", false);
     489    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FWHM_MAJ",   PS_META_REPLACE, "PSF FWHM Major axis (mean)", FWHM_X);
     490    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_SG",   PS_META_REPLACE, "PSF FWHM Major axis (sigma)", 0);
     491    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_LQ",   PS_META_REPLACE, "PSF FWHM Major axis (lower quartile)", 0);
     492    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_UQ",   PS_META_REPLACE, "PSF FWHM Major axis (upper quartile)", 0);
     493    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FWHM_MIN",   PS_META_REPLACE, "PSF FWHM Minor axis (mean)", FWHM_Y);
     494    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_SG",   PS_META_REPLACE, "PSF FWHM Minor axis (sigma)", 0);
     495    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_LQ",   PS_META_REPLACE, "PSF FWHM Minor axis (lower quartile)", 0);
     496    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_UQ",   PS_META_REPLACE, "PSF FWHM Minor axis (upper quartile)", 0);
     497    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "ANGLE",    PS_META_REPLACE, "PSF angle",           FWHM_T);
     498    psMetadataAddS32 (readout->analysis, PS_LIST_TAIL, "NPSFSTAR", PS_META_REPLACE, "Number of stars used to make PSF", 0);
     499    psMetadataAddBool(readout->analysis, PS_LIST_TAIL, "PSFMODEL", PS_META_REPLACE, "Valid PSF Model?", false);
    501500
    502501    return true;

  • branches/eam_branches/psphot.stack.20100120/src/psphotEfficiency.c

    r26691 r26748  
    3535
    3636    psKernel *kernel = psImageSmoothKernel(smoothSigma, smoothNsigma); // Kernel used for smoothing
    37     psKernel *newCovar = psImageCovarianceCalculate(kernel, ro->covariance); // Covariance matrix
     37    *covarFactor = psImageCovarianceCalculateFactor(kernel, ro->covariance);
    3838    psFree(kernel);
    39     *covarFactor = psImageCovarianceFactor(newCovar); // Variance factor
    40     psFree(newCovar);
    4139
    4240    psStats *stats = psStatsAlloc(PS_STAT_ROBUST_MEDIAN); // Statistics for variance

  • branches/eam_branches/psphot.stack.20100120/src/psphotFake.c

    r26266 r26748  
    2828    float smoothSigma  = 0.5*(FWHM_X + FWHM_Y) / (2.0*sqrtf(2.0*log(2.0)));
    2929    psKernel *kernel = psImageSmoothKernel(smoothSigma, smoothNsigma); // Kernel used for smoothing
    30     psKernel *newCovar = psImageCovarianceCalculate(kernel, readout->covariance); // Covariance matrix
     30    float factor = psImageCovarianceCalculateFactor(kernel, readout->covariance); // Covariance matrix
    3131    psFree(kernel);
    32     float factor = psImageCovarianceFactor(newCovar); // Variance factor
    33     psFree(newCovar);
    3432
    3533    psStats *stats = psStatsAlloc(PS_STAT_ROBUST_MEDIAN); // Statistics for variance
    … …  
    213211
    214212    // Putting results on recipe because that appears to be the psphot standard, but it's not a good idea
    215     psMetadataAddVector(recipe, PS_LIST_TAIL, "FAKE.EFF", PS_META_REPLACE,
    216                         "Efficiency fractions", frac);
    217     psMetadataAddVector(recipe, PS_LIST_TAIL, "FAKE.MAG", PS_META_REPLACE,
    218                         "Efficiency magnitudes", magOffsets);
    219     psMetadataAddF32(recipe, PS_LIST_TAIL, "FAKE.REF", PS_META_REPLACE,
    220                      "Efficiency reference magnitude", magLim);
     213    psMetadataAddVector(readout->analysis, PS_LIST_TAIL, "FAKE.EFF", PS_META_REPLACE, "Efficiency fractions", frac);
     214    psMetadataAddVector(readout->analysis, PS_LIST_TAIL, "FAKE.MAG", PS_META_REPLACE, "Efficiency magnitudes", magOffsets);
     215    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "FAKE.REF", PS_META_REPLACE, "Efficiency reference magnitude", magLim);
    221216
    222217    psFree(frac);

  • branches/eam_branches/psphot.stack.20100120/src/psphotFindDetections.c

    r26691 r26748  
    44// peaks and new footprints.  any old peaks are saved on oldPeaks.  the resulting footprint set
    55// contains all footprints (old and new)
    6 bool psphotFindDetections (pmConfig *config, const pmFPAview *view)
     6bool psphotFindDetections (pmConfig *config, const pmFPAview *view, bool firstPass)
    77{
    88    bool status = true;
    … …  
    1717    // loop over the available readouts
    1818    for (int i = 0; i < num; i++) {
    19         if (!psphotFindDetectionsReadout (config, view, "PSPHOT.INPUT", i, recipe)) {
     19        if (!psphotFindDetectionsReadout (config, view, "PSPHOT.INPUT", i, recipe, firstPass)) {
    2020            psError (PSPHOT_ERR_CONFIG, false, "failed to find initial detections for PSPHOT.INPUT entry %d", i);
    2121            return false;
    … …  
    2626
    2727// smooth the image, search for peaks, optionally define footprints based on the peaks
    28 bool psphotFindDetectionsReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index, psMetadata *recipe) {
     28bool psphotFindDetectionsReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index, psMetadata *recipe, bool firstPass) {
    2929
    3030    bool status;
    … …  
    5454       
    5555        // save detections on the readout->analysis
    56         if (!psMetadataAddPtr (readout->analysis, PS_LIST_TAIL, "PSPHOT.DETECTIONS", PS_META_REPLACE | PS_DATA_UNKNOWN, "psphot detectinos", detections)) {
     56        if (!psMetadataAddPtr (readout->analysis, PS_LIST_TAIL, "PSPHOT.DETECTIONS", PS_META_REPLACE | PS_DATA_UNKNOWN, "psphot detections", detections)) {
    5757            psError (PSPHOT_ERR_CONFIG, false, "problem saving detections on readout");
    5858            return false;
    5959        }
     60    } else {
     61        psMemIncrRefCounter(detections); // so we can free the detections below
     62    }
    6063
     64    if (firstPass) {
    6165        pass = 1;
    6266        NSIGMA_PEAK = psMetadataLookupF32 (&status, recipe, "PEAKS_NSIGMA_LIMIT"); PS_ASSERT (status, NULL);
    6367        NMAX = psMetadataLookupS32 (&status, recipe, "PEAKS_NMAX"); PS_ASSERT (status, NULL);
    6468    } else {
    65         psMemIncrRefCounter(detections); // so we can free the detections below
    6669        pass = 2;
    6770        NSIGMA_PEAK = psMetadataLookupF32 (&status, recipe, "PEAKS_NSIGMA_LIMIT_2"); PS_ASSERT (status, NULL);
    … …  
    9699    // detect the peaks in the significance image
    97100    detections->peaks = psphotFindPeaks (significance, readout, recipe, threshold, NMAX);
    98     psMetadataAddF32  (recipe, PS_LIST_TAIL, "PEAK_THRESHOLD", PS_META_REPLACE, "Peak Detection Threshold", threshold);
     101    psMetadataAddF32  (readout->analysis, PS_LIST_TAIL, "PEAK_THRESHOLD", PS_META_REPLACE, "Peak Detection Threshold", threshold);
    99102    if (!detections->peaks) {
    100103        // we only get a NULL peaks array due to a programming or config error.
    … …  
    123126// if we use the footprints, the output peaks list contains both old and new peaks,
    124127// otherwise it only contains the new peaks.
    125 

  • branches/eam_branches/psphot.stack.20100120/src/psphotForcedReadout.c

    r26542 r26748  
    2525    }
    2626
    27     // find the currently selected readout
    28     pmReadout  *readout = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT");
    29     PS_ASSERT_PTR_NON_NULL (readout, false);
    30 
    3127    // optional break-point for processing
    3228    char *breakPt = psMetadataLookupStr (NULL, recipe, "BREAK_POINT");
    … …  
    3430
    3531    // Generate the mask and weight images, including the user-defined analysis region of interest
    36     psphotSetMaskAndVariance (config, view, recipe);
     32    psphotSetMaskAndVariance (config, view);
    3733    if (!strcasecmp (breakPt, "NOTHING")) {
    38         return psphotReadoutCleanup(config, readout, recipe, NULL, NULL, NULL);
     34        return psphotReadoutCleanup(config, view);
    3935    }
    40 
    41     // display the image, weight, mask (ch 1,2,3)
    42     psphotVisualShowImage (readout);
    4336
    4437    // generate a background model (median, smoothed image)
    4538    if (!psphotModelBackground (config, view)) {
    46         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     39        return psphotReadoutCleanup (config, view);
    4740    }
    4841    if (!psphotSubtractBackground (config, view)) {
    49         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     42        return psphotReadoutCleanup (config, view);
    5043    }
    5144    if (!strcasecmp (breakPt, "BACKMDL")) {
    52         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     45        return psphotReadoutCleanup (config, view);
    5346    }
    5447
    55     // display the backsub and backgnd images
    56     psphotVisualShowBackground (config, view, readout);
    57 
    58     // load the psf model, if suppled.  FWHM_X,FWHM_Y,etc are saved in the recipe
    59     pmPSF *psf = psphotLoadPSF (config, view, recipe);
    60     if (!psf) {
    61         psError(PSPHOT_ERR_CONFIG, false, "unable to load psf model");
    62         return false;
     48    if (!psphotLoadPSF (config, view)) {
     49        // this only happens if we had a programming error in psphotLoadPSF
     50        psError (PSPHOT_ERR_UNKNOWN, false, "error loading psf model");
     51        return psphotReadoutCleanup (config, view);
    6352    }
    6453
    6554    // include externally-supplied sources
    66     psArray *sources = psArrayAllocEmpty(100);
    67     psphotLoadExtSources (config, view, sources);
     55    psphotLoadExtSources (config, view);
    6856
    6957    // construct an initial model for each object, set the radius to fitRadius, set circular fit mask
    70     psphotGuessModels (config, readout, sources, psf);
     58    psphotGuessModels (config, view);
    7159
    7260    // linear PSF fit to source peaks, subtract the models from the image (in PSF mask)
    73     psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE);
     61    psphotFitSourcesLinear (config, view, false);
    7462
    7563    // identify CRs and extended sources
    … …  
    7967
    8068    // calculate source magnitudes
    81     psphotMagnitudes(config, readout, view, sources, psf);
     69    psphotMagnitudes(config, view);
    8270
    8371    // XXX do I want to do this?
    … …  
    9179
    9280    // drop the references to the image pixels held by each source
    93     psphotSourceFreePixels (sources);
     81    psphotSourceFreePixels (config, view);
    9482
    9583    // create the exported-metadata and free local data
    96     return psphotReadoutCleanup(config, readout, recipe, NULL, psf, sources);
     84    return psphotReadoutCleanup(config, view);
    9785}

  • branches/eam_branches/psphot.stack.20100120/src/psphotGuessModels.c

    r26691 r26748  
    7575
    7676    // setup the PSF fit radius details
    77     psphotInitRadiusPSF (recipe, psf->type);
     77    psphotInitRadiusPSF (recipe, readout->analysis, psf->type);
    7878
    7979    // choose Cx, Cy (see psphotThreadTools.c for overview of the concepts)

  • branches/eam_branches/psphot.stack.20100120/src/psphotImageQuality.c

    r26691 r26748  
    132132    }
    133133
    134     psMetadataAddS32(recipe, PS_LIST_TAIL, "IQ_NSTAR", PS_META_REPLACE,
    135                      "Number of stars used for IQ measurements", M2->n);
     134    psMetadataAddS32(readout->analysis, PS_LIST_TAIL, "IQ_NSTAR", PS_META_REPLACE, "Number of stars used for IQ measurements", M2->n);
    136135
    137136// XXX make this a recipe option
    … …  
    145144
    146145    float fwhm_major = stats->sampleMean;
    147     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW1",   PS_META_REPLACE,
    148                      "FWHM of Major Axis from moments", stats->sampleMean);
    149     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW1_E", PS_META_REPLACE,
    150                      "FWHM scatter (Major) from moments", stats->sampleStdev);
     146    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW1",   PS_META_REPLACE, "FWHM of Major Axis from moments", stats->sampleMean);
     147    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW1_E", PS_META_REPLACE, "FWHM scatter (Major) from moments", stats->sampleStdev);
    151148
    152149    if (!psVectorStats(stats, FWHM_MINOR, NULL, NULL, 0)) {
    … …  
    155152    }
    156153    float fwhm_minor = stats->sampleMean;
    157     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW2",   PS_META_REPLACE,
    158                      "FWHM of Minor Axis from moments", stats->sampleMean);
    159     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW2_E", PS_META_REPLACE,
    160                      "FWHM scatter (Minor) from moments", stats->sampleStdev);
     154    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW2",   PS_META_REPLACE, "FWHM of Minor Axis from moments", stats->sampleMean);
     155    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW2_E", PS_META_REPLACE, "FWHM scatter (Minor) from moments", stats->sampleStdev);
    161156
    162157    if (!psVectorStats(stats, M2, NULL, NULL, 0)) {
    … …  
    166161    float vM2 = stats->sampleMean;
    167162    float dM2 = stats->sampleStdev;
    168     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2",    PS_META_REPLACE,
    169                      "M_2 = sqrt (M_c2^2 + M_s2^2)", vM2);
    170     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_ER", PS_META_REPLACE,
    171                      "Stdev of M_2", dM2);
    172     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_LQ", PS_META_REPLACE,
    173                      "Lower Quartile of M_2", stats->sampleLQ);
    174     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_UQ", PS_META_REPLACE,
    175                      "Upper Quartile of M_2", stats->sampleUQ);
     163    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2",    PS_META_REPLACE, "M_2 = sqrt (M_c2^2 + M_s2^2)", vM2);
     164    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_ER", PS_META_REPLACE, "Stdev of M_2", dM2);
     165    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_LQ", PS_META_REPLACE, "Lower Quartile of M_2", stats->sampleLQ);
     166    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_UQ", PS_META_REPLACE, "Upper Quartile of M_2", stats->sampleUQ);
    176167
    177168    if (!psVectorStats(stats, M2c, NULL, NULL, 0)) {
    … …  
    179170        goto FAIL;
    180171    }
    181     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C",   PS_META_REPLACE,
    182                      "M_2c = sum f r^2 cos(2phi) / sum f", stats->sampleMean);
    183     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_E", PS_META_REPLACE,
    184                      "Stdev of M_2c", stats->sampleStdev);
    185     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_L", PS_META_REPLACE,
    186                      "Lower Quartile of M_2c", stats->sampleLQ);
    187     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_U", PS_META_REPLACE,
    188                      "Upper Quartile of M_2c", stats->sampleUQ);
     172    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C",   PS_META_REPLACE, "M_2c = sum f r^2 cos(2phi) / sum f", stats->sampleMean);
     173    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_E", PS_META_REPLACE, "Stdev of M_2c", stats->sampleStdev);
     174    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_L", PS_META_REPLACE, "Lower Quartile of M_2c", stats->sampleLQ);
     175    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_U", PS_META_REPLACE, "Upper Quartile of M_2c", stats->sampleUQ);
    189176
    190177    if (!psVectorStats(stats, M2s, NULL, NULL, 0)) {
    … …  
    192179        goto FAIL;
    193180    }
    194     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S",   PS_META_REPLACE,
    195                      "M_2s = sum f r^2 cos(2phi) / sum f", stats->sampleMean);
    196     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_E", PS_META_REPLACE,
    197                      "Stdev of M_2s", stats->sampleStdev);
    198     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_L", PS_META_REPLACE,
    199                      "Lower Quartile of M_2s", stats->sampleLQ);
    200     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_U", PS_META_REPLACE,
    201                      "Upper Quartile of M_2s", stats->sampleUQ);
     181    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S",   PS_META_REPLACE, "M_2s = sum f r^2 cos(2phi) / sum f", stats->sampleMean);
     182    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_E", PS_META_REPLACE, "Stdev of M_2s", stats->sampleStdev);
     183    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_L", PS_META_REPLACE, "Lower Quartile of M_2s", stats->sampleLQ);
     184    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_U", PS_META_REPLACE, "Upper Quartile of M_2s", stats->sampleUQ);
    202185
    203186    if (!psVectorStats(stats, M3, NULL, NULL, 0)) {
    … …  
    207190    float vM3 = stats->sampleMean;
    208191    float dM3 = stats->sampleStdev;
    209     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3",   PS_META_REPLACE,
    210                      "M_3 = sqrt (M_c3^2 + M_s3^2)", vM3);
    211     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_ER", PS_META_REPLACE,
    212                      "Stdev of M_3", dM3);
    213     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_LQ", PS_META_REPLACE,
    214                      "Lower Quartile of M_3", stats->sampleLQ);
    215     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_UQ", PS_META_REPLACE,
    216                      "Upper Quartile of M_3", stats->sampleUQ);
     192    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3",   PS_META_REPLACE, "M_3 = sqrt (M_c3^2 + M_s3^2)", vM3);
     193    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_ER", PS_META_REPLACE, "Stdev of M_3", dM3);
     194    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_LQ", PS_META_REPLACE, "Lower Quartile of M_3", stats->sampleLQ);
     195    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_UQ", PS_META_REPLACE, "Upper Quartile of M_3", stats->sampleUQ);
    217196
    218197    if (!psVectorStats(stats, M4, NULL, NULL, 0)) {
    … …  
    222201    float vM4 = stats->sampleMean;
    223202    float dM4 = stats->sampleStdev;
    224     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4",   PS_META_REPLACE,
    225                      "M_4 = sqrt (M_c4^2 + M_s4^2)", vM4);
    226     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_ER", PS_META_REPLACE,
    227                      "Stdev of M_4", dM4);
    228     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_LQ", PS_META_REPLACE,
    229                      "Lower Quartile of M_4", stats->sampleLQ);
    230     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_UQ", PS_META_REPLACE,
    231                      "Upper Quartile of M_4", stats->sampleUQ);
     203    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4",   PS_META_REPLACE, "M_4 = sqrt (M_c4^2 + M_s4^2)", vM4);
     204    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_ER", PS_META_REPLACE, "Stdev of M_4", dM4);
     205    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_LQ", PS_META_REPLACE, "Lower Quartile of M_4", stats->sampleLQ);
     206    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_UQ", PS_META_REPLACE, "Upper Quartile of M_4", stats->sampleUQ);
    232207
    233208#else
    … …  
    239214    }
    240215    float fwhm_major = stats->robustMedian;
    241     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW1",   PS_META_REPLACE,
    242                      "FWHM of Major Axis from moments", stats->robustMedian);
    243     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW1_E", PS_META_REPLACE,
    244                      "FWHM scatter (Major) from moments", stats->robustStdev);
     216    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW1",   PS_META_REPLACE, "FWHM of Major Axis from moments", stats->robustMedian);
     217    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW1_E", PS_META_REPLACE, "FWHM scatter (Major) from moments", stats->robustStdev);
    245218
    246219    if (!psVectorStats(stats, FWHM_MINOR, NULL, NULL, 0)) {
    … …  
    249222    }
    250223    float fwhm_minor = stats->robustMedian;
    251     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW2",   PS_META_REPLACE,
    252                      "FWHM of Minor Axis from moments", stats->robustMedian);
    253     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW2_E", PS_META_REPLACE,
    254                      "FWHM scatter (Minor) from moments", stats->robustStdev);
     224    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW2",   PS_META_REPLACE, "FWHM of Minor Axis from moments", stats->robustMedian);
     225    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW2_E", PS_META_REPLACE, "FWHM scatter (Minor) from moments", stats->robustStdev);
    255226
    256227    if (!psVectorStats(stats, M2, NULL, NULL, 0)) {
    … …  
    260231    float vM2 = stats->robustMedian;
    261232    float dM2 = stats->robustStdev;
    262     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2",    PS_META_REPLACE,
    263                      "M_2 = sqrt (M_c2^2 + M_s2^2)", vM2);
    264     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_ER", PS_META_REPLACE,
    265                      "Stdev of M_2", dM2);
    266     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_LQ", PS_META_REPLACE,
    267                      "Lower Quartile of M_2", stats->robustLQ);
    268     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_UQ", PS_META_REPLACE,
    269                      "Upper Quartile of M_2", stats->robustUQ);
     233    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2",    PS_META_REPLACE, "M_2 = sqrt (M_c2^2 + M_s2^2)", vM2);
     234    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_ER", PS_META_REPLACE, "Stdev of M_2", dM2);
     235    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_LQ", PS_META_REPLACE, "Lower Quartile of M_2", stats->robustLQ);
     236    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_UQ", PS_META_REPLACE, "Upper Quartile of M_2", stats->robustUQ);
    270237
    271238    if (!psVectorStats(stats, M2c, NULL, NULL, 0)) {
    … …  
    273240        goto FAIL;
    274241    }
    275     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C",   PS_META_REPLACE,
    276                      "M_2c = sum f r^2 cos(2phi) / sum f", stats->robustMedian);
    277     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_E", PS_META_REPLACE,
    278                      "Stdev of M_2c", stats->robustStdev);
    279     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_L", PS_META_REPLACE,
    280                      "Lower Quartile of M_2c", stats->robustLQ);
    281     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_U", PS_META_REPLACE,
    282                      "Upper Quartile of M_2c", stats->robustUQ);
     242    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C",   PS_META_REPLACE, "M_2c = sum f r^2 cos(2phi) / sum f", stats->robustMedian);
     243    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_E", PS_META_REPLACE, "Stdev of M_2c", stats->robustStdev);
     244    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_L", PS_META_REPLACE, "Lower Quartile of M_2c", stats->robustLQ);
     245    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_U", PS_META_REPLACE, "Upper Quartile of M_2c", stats->robustUQ);
    283246
    284247    if (!psVectorStats(stats, M2s, NULL, NULL, 0)) {
    … …  
    286249        goto FAIL;
    287250    }
    288     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S",   PS_META_REPLACE,
    289                      "M_2s = sum f r^2 cos(2phi) / sum f", stats->robustMedian);
    290     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_E", PS_META_REPLACE,
    291                      "Stdev of M_2s", stats->robustStdev);
    292     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_L", PS_META_REPLACE,
    293                      "Lower Quartile of M_2s", stats->robustLQ);
    294     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_U", PS_META_REPLACE,
    295                      "Upper Quartile of M_2s", stats->robustUQ);
     251    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S",   PS_META_REPLACE, "M_2s = sum f r^2 cos(2phi) / sum f", stats->robustMedian);
     252    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_E", PS_META_REPLACE, "Stdev of M_2s", stats->robustStdev);
     253    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_L", PS_META_REPLACE, "Lower Quartile of M_2s", stats->robustLQ);
     254    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_U", PS_META_REPLACE, "Upper Quartile of M_2s", stats->robustUQ);
    296255
    297256    if (!psVectorStats(stats, M3, NULL, NULL, 0)) {
    … …  
    301260    float vM3 = stats->robustMedian;
    302261    float dM3 = stats->robustStdev;
    303     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3",   PS_META_REPLACE,
    304                      "M_3 = sqrt (M_c3^2 + M_s3^2)", vM3);
    305     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_ER", PS_META_REPLACE,
    306                      "Stdev of M_3", dM3);
    307     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_LQ", PS_META_REPLACE,
    308                      "Lower Quartile of M_3", stats->robustLQ);
    309     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_UQ", PS_META_REPLACE,
    310                      "Upper Quartile of M_3", stats->robustUQ);
     262    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3",   PS_META_REPLACE, "M_3 = sqrt (M_c3^2 + M_s3^2)", vM3);
     263    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_ER", PS_META_REPLACE, "Stdev of M_3", dM3);
     264    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_LQ", PS_META_REPLACE, "Lower Quartile of M_3", stats->robustLQ);
     265    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_UQ", PS_META_REPLACE, "Upper Quartile of M_3", stats->robustUQ);
    311266
    312267    if (!psVectorStats(stats, M4, NULL, NULL, 0)) {
    … …  
    316271    float vM4 = stats->robustMedian;
    317272    float dM4 = stats->robustStdev;
    318     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4",   PS_META_REPLACE,
    319                      "M_4 = sqrt (M_c4^2 + M_s4^2)", vM4);
    320     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_ER", PS_META_REPLACE,
    321                      "Stdev of M_4", dM4);
    322     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_LQ", PS_META_REPLACE,
    323                      "Lower Quartile of M_4", stats->robustLQ);
    324     psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_UQ", PS_META_REPLACE,
    325                      "Upper Quartile of M_4", stats->robustUQ);
     273    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4",   PS_META_REPLACE, "M_4 = sqrt (M_c4^2 + M_s4^2)", vM4);
     274    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_ER", PS_META_REPLACE, "Stdev of M_4", dM4);
     275    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_LQ", PS_META_REPLACE, "Lower Quartile of M_4", stats->robustLQ);
     276    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_UQ", PS_META_REPLACE, "Upper Quartile of M_4", stats->robustUQ);
    326277#endif
    327278

  • branches/eam_branches/psphot.stack.20100120/src/psphotLoadPSF.c

    r26691 r26748  
    11# include "psphotInternal.h"
    22
     3// NOTE : pmPSF_IO.c functions must load the psf model onto the chip->analysis metadata because
     4// the I/O operation likely occurs before the readout exists.  This implementation assumes that
     5// a single psf model is valid for the entire set of readouts (not valid for a time series of readouts)
     6
     7// XXX for now (2010.01.27), the supporting programs do not define multiple PSPHOT.PSF.LOAD
     8// files to go with multiple PSPHOT.INPUT files.  as a result, the implementation below is
     9// currently going to work for the case of a single input file, but will fail if we try with a
     10// stack of images.
     11
    312// load an externally supplied psf model
    4 bool psphotLoadPSFReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index, psMetadata *recipe) {
     13bool psphotLoadPSFReadout (pmConfig *config, const pmFPAview *view, const char *filename, int index) {
    514
    615    bool status;
    … …  
    2534    }
    2635
    27     if (!psphotPSFstats (readout, recipe, psf)) {
     36    if (!psphotPSFstats (readout, psf)) {
    2837        psAbort("cannot measure PSF shape terms");
    2938    }
    … …  
    4453    bool status = false;
    4554
    46     // select the appropriate recipe information
    47     psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
    48     psAssert (recipe, "missing recipe?");
    49 
    50     // XXX PSPHOT.PSF.LOAD vs PSPHOT.INPUT??
     55    // XXX PSPHOT.PSF.LOAD vs PSPHOT.INPUT -- see note at top
    5156    int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM");
    5257    psAssert (status, "programming error: must define PSPHOT.INPUT.NUM");
    … …  
    5661
    5762        // Generate the mask and weight images, including the user-defined analysis region of interest
    58         if (!psphotLoadPSFReadout (config, view, "PSPHOT.PSF.LOAD", i, recipe)) {
     63        if (!psphotLoadPSFReadout (config, view, "PSPHOT.PSF.LOAD", i)) {
    5964            psError (PSPHOT_ERR_CONFIG, false, "failed to load PSF model for PSPHOT.PSF.LOAD entry %d", i);
    6065            return false;

  • branches/eam_branches/psphot.stack.20100120/src/psphotMakePSFReadout.c

    r26542 r26748  
    2424    }
    2525
    26     // find the currently selected readout
    27     pmReadout  *readout = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT");
    28     PS_ASSERT_PTR_NON_NULL (readout, false);
    29 
    3026    // optional break-point for processing
    3127    char *breakPt = psMetadataLookupStr (NULL, recipe, "BREAK_POINT");
    … …  
    3329
    3430    // Generate the mask and weight images, including the user-defined analysis region of interest
    35     psphotSetMaskAndVariance (config, view, recipe);
     31    psphotSetMaskAndVariance (config, view);
    3632    if (!strcasecmp (breakPt, "NOTHING")) {
    37         return psphotReadoutCleanup(config, readout, recipe, NULL, NULL, NULL);
     33        return psphotReadoutCleanup(config, view);
    3834    }
    39 
    40     // display the image, weight, mask (ch 1,2,3)
    41     psphotVisualShowImage (readout);
    4235
    4336    // generate a background model (median, smoothed image)
    4437    if (!psphotModelBackground (config, view)) {
    45         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     38        return psphotReadoutCleanup (config, view);
    4639    }
    4740    if (!psphotSubtractBackground (config, view)) {
    48         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     41        return psphotReadoutCleanup (config, view);
    4942    }
    5043    if (!strcasecmp (breakPt, "BACKMDL")) {
    51         return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     44        return psphotReadoutCleanup (config, view);
    5245    }
    5346
    54     // display the backsub and backgnd images
    55     psphotVisualShowBackground (config, view, readout);
     47    psphotLoadExtSources (config, view);
    5648
    57     pmDetections *detections = NULL;
    58 
    59     // If sources have been supplied, then these should be used to measure the PSF
    60     // include externally-supplied sources
    61 
    62     // XXX sources loaded from a text file have no valid SN values, but psphotChoosePSF
    63     // selected the top PSF_MAX_NSTARS to generate the PSF, excluding an arbitrary subset.
    64     psArray *sources = psArrayAllocEmpty(100);
    65     psphotLoadExtSources (config, view, sources);
    66     if (sources->n) {
    67         // the user wants to make the psf from these stars; define them as psf stars:
    68         for (int i = 0; i < sources->n; i++) {
    69             pmSource *source = sources->data[i];
    70             source->mode |= PM_SOURCE_MODE_PSFSTAR;
    71         }
    72         // force psphotChoosePSF to use all loaded sources
    73         psMetadataAddS32 (recipe, PS_LIST_TAIL, "PSF_MAX_NSTARS", PS_META_REPLACE, "fit radius", sources->n);
    74 
    75         // measure stats of externally specified sources
    76         if (!psphotSourceStatsUpdate (sources, config, readout)) {
    77             psError(PSPHOT_ERR_CONFIG, false, "failure to measure stats of existing sources");
    78             return false;
    79         }
    80 
    81     } else {
    82         // find the detections (by peak and/or footprint) in the image.
    83         detections = psphotFindDetections (NULL, readout, recipe);
    84         if (!detections) {
    85             psLogMsg ("psphot", 3, "unable to find detections in this image");
    86             return psphotReadoutCleanup (config, readout, recipe, detections, NULL, sources);
    87         }
    88 
    89         // construct sources and measure basic stats
    90         psFree (sources);
    91         sources = psphotSourceStats (config, readout, detections, true);
    92         if (!sources) return false;
    93 
    94         // find blended neighbors of very saturated stars
    95         // XXX merge this with Basic Deblend?
    96         psphotDeblendSatstars (readout, sources, recipe);
    97 
    98         // mark blended peaks PS_SOURCE_BLEND
    99         if (!psphotBasicDeblend (sources, recipe)) {
    100             psLogMsg ("psphot", 3, "failed on deblend analysis");
    101             return psphotReadoutCleanup (config, readout, recipe, detections, NULL, sources);
    102         }
    103 
    104         // classify sources based on moments, brightness
    105         if (!psphotRoughClass (readout, sources, recipe, false)) {
    106             psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
    107             return psphotReadoutCleanup (config, readout, recipe, detections, NULL, sources);
    108         }
     49    // If sources have been supplied, then these should be used to measure the PSF include
     50    // externally-supplied sources; if not, we need to generate a set of possible PSF sources.
     51    // This function updates the SN entries for the loaded sources or generates a set of
     52    // detections from the image, if no external ones have been supplied.  Sources loaded from
     53    // a text file have no valid SN values, but psphotChoosePSF needs to select the top
     54    // PSF_MAX_NSTARS to generate the PSF.
     55    if (!psphotCheckExtSources (config, view)) {
     56        psLogMsg ("psphot", 3, "failure to select possible PSF sources (external or internal)");
     57        return psphotReadoutCleanup (config, view);
    10958    }
    11059
    111     // use bright stellar objects to measure PSF
    112     // XXX if we do not have enough stars to generate the PSF, build one
    113     // from the SEEING guess and model class
    114     pmPSF *psf = psphotChoosePSF (readout, sources, recipe);
    115     if (psf == NULL) {
     60    // Use bright stellar objects to measure PSF. If we do not have enough stars to generate
     61    // the PSF, build one from the SEEING guess and model class
     62    if (!psphotChoosePSF (config, view)) {
    11663        psLogMsg ("psphot", 3, "failure to construct a psf model");
    117         return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
     64        return psphotReadoutCleanup (config, view);
    11865    }
    119     psphotVisualShowPSFModel (readout, psf);
    120     return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
    121 }
    12266
    12367    // measure aperture photometry corrections
    124     // XXX isn't this part of the PSF??
    125     // if (!psphotApResid (config, readout, sources, psf)) {
    126     //     psLogMsg ("psphot", 3, "failed on psphotApResid");
    127     //     return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
    128     // }
     68# if 0
     69    if (!psphotApResid (config, view)) {
     70        psLogMsg ("psphot", 3, "failed on psphotApResid");
     71        return psphotReadoutCleanup (config, view);
     72    }
     73# endif
     74
     75    return psphotReadoutCleanup (config, view);
     76}

  • branches/eam_branches/psphot.stack.20100120/src/psphotMaskReadout.c

    r26691 r26748  
    3636    psAssert (readout, "missing readout?");
    3737
    38     // ** Interpret the mask values:
    39     // XXX drop the write to recipe and move config into psphotRoughClass?
    40     // XXX alternatively, define a function to set the psphot recipe masks
     38    // save maskSat and maskBad on the psphot recipe (mostly for psphotRoughClass)
    4139    psImageMaskType maskSat  = pmConfigMaskGet("SAT", config); // Mask value for saturated pixels
    4240    psMetadataAddImageMask (recipe, PS_LIST_TAIL, "MASK.SAT", PS_META_REPLACE, "user-defined mask", maskSat);
    … …  
    4442    psImageMaskType maskBad  = pmConfigMaskGet("LOW", config); // Mask value for low pixels
    4543    if (!maskBad) {
    46         // XXX: for backward compatability look up old name
     44        // for backward compatability look up old name
    4745        maskBad  = pmConfigMaskGet("BAD", config);
    4846    }

  • branches/eam_branches/psphot.stack.20100120/src/psphotMergeSources.c

    r26691 r26748  
    4141    psAssert (newSources, "missing sources?");
    4242
     43    // XXX TEST:
     44    if (detections->allSources) {
     45        psphotMaskCosmicRayFootprintCheck(detections->allSources);
     46    }
     47    if (detections->newSources) {
     48        psphotMaskCosmicRayFootprintCheck(detections->newSources);
     49    }
     50
    4351    if (!detections->allSources) {
    4452        detections->allSources = psArrayAllocEmpty(newSources->n);
    … …  
    5159    }
    5260
    53     // XXX free detections->newSources?
    54     return true;
    55 }
    56 
    57 // merge the externally supplied sources with the existing sources.  mark them as having
    58 // mode PM_SOURCE_MODE_EXTERNAL
    59 // XXX this function needs to be updated to work with the new context of pshot inputs
    60 bool psphotLoadExtSources (pmConfig *config, const pmFPAview *view, psArray *sources) {
    61 
    62     psArray *extSourcesCMF = NULL;
     61    psFree (detections->newSources);
     62    detections->newSources = NULL;
     63
     64    return true;
     65}
     66
     67// Merge the externally supplied sources with the existing sources.  Mark them as having mode
     68// PM_SOURCE_MODE_EXTERNAL.
     69
     70// XXX This function needs to be updated to loop over set of input files.  At the moment, we
     71// only expect a single entry for PSPHOT.INPUT.CMF and PSPHOT.SOURCES.TEXT, so we can only
     72// associate input sources with a single entry for PSPHOT.INPUT
     73bool psphotLoadExtSources (pmConfig *config, const pmFPAview *view) {
     74
     75    bool status;
     76    pmDetections *extCMF = NULL;
    6377    psArray *extSourcesTXT = NULL;
    6478
     79    // find the currently selected readout
     80    pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", 0); // File of interest
     81    psAssert (file, "missing file?");
     82
     83    pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
     84    psAssert (readout, "missing readout?");
     85
     86    pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
     87    psAssert (detections, "missing detections?");
     88
     89    // XXX allSources of newSources?
     90    psArray *sources = detections->allSources;
     91    psAssert (sources, "missing sources?");
     92
    6593    // load data from input CMF file:
    66     { 
    67         pmReadout *readoutCMF = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT.CMF");
    68         if (!readoutCMF) goto loadTXT;
    69 
    70         extSourcesCMF = psMetadataLookupPtr (NULL, readoutCMF->analysis, "PSPHOT.SOURCES");
    71         if (extSourcesCMF) {
    72             for (int i = 0; i < extSourcesCMF->n; i++) {
    73                 pmSource *source = extSourcesCMF->data[i];
     94    {
     95        pmReadout *readoutCMF = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT.CMF");
     96        if (!readoutCMF) goto loadTXT;
     97
     98        extCMF = psMetadataLookupPtr (NULL, readoutCMF->analysis, "PSPHOT.DETECTIONS");
     99        if (extCMF) {
     100            for (int i = 0; i < extCMF->allSources->n; i++) {
     101                pmSource *source = extCMF->allSources->data[i];
    74102                source->mode |= PM_SOURCE_MODE_EXTERNAL;
    75103
    76                 // the supplied peak flux needs to be re-normalized
    77                 source->peak->flux = 1.0;
    78                 source->peak->value = 1.0;
     104                // the supplied peak flux needs to be re-normalized
     105                source->peak->flux = 1.0;
     106                source->peak->value = 1.0;
    79107
    80108                // drop the loaded source modelPSF
    81109                psFree (source->modelPSF);
    82110                source->modelPSF = NULL;
     111
     112                psArrayAdd (detections->allSources, 100, source);
    83113            }
    84             // XXX psphotMergeSources (sources, extSourcesCMF);
    85114        }
    86115    }
    … …  
    89118
    90119    // load data from input TXT file:
    91     { 
    92         pmChip *chipTXT = pmFPAfileThisChip (config->files, view, "PSPHOT.INPUT");
    93         if (!chipTXT) goto finish;
    94 
    95         extSourcesTXT = psMetadataLookupPtr (NULL, chipTXT->analysis, "PSPHOT.SOURCES.TEXT");
    96         if (extSourcesTXT) {
    97             for (int i = 0; i < extSourcesTXT->n; i++) {
    98                 pmSource *source = extSourcesTXT->data[i];
    99                 source->mode |= PM_SOURCE_MODE_EXTERNAL;
     120    {
     121        pmChip *chipTXT = pmFPAfileThisChip (config->files, view, "PSPHOT.INPUT");
     122        if (!chipTXT) goto finish;
     123
     124        extSourcesTXT = psMetadataLookupPtr (NULL, chipTXT->analysis, "PSPHOT.SOURCES.TEXT");
     125        if (extSourcesTXT) {
     126            for (int i = 0; i < extSourcesTXT->n; i++) {
     127                pmSource *source = extSourcesTXT->data[i];
     128                source->mode |= PM_SOURCE_MODE_EXTERNAL;
    100129
    101130                // drop the loaded source modelPSF
    102131                psFree (source->modelPSF);
    103132                source->modelPSF = NULL;
     133
     134                psArrayAdd (detections->allSources, 100, source);
    104135            }
    105             // XXX psphotMergeSources (sources, extSourcesTXT);
    106136        }
    107137    }
    … …  
    109139finish:
    110140
    111     if (!extSourcesTXT && !extSourcesTXT) {
     141    if (!extCMF && !extSourcesTXT) {
    112142        psLogMsg ("psphot", 3, "no external sources for this readout");
    113143        return true;
    114144    }
    115145
    116     int nCMF = extSourcesCMF ? extSourcesCMF->n : 0;
     146    int nCMF = extCMF        ? extCMF->allSources->n        : 0;
    117147    int nTXT = extSourcesTXT ? extSourcesTXT->n : 0;
    118148
    119     psLogMsg ("psphot", 3, "%d external sources (%d cmf, %d text) merged to yield %ld total sources", 
    120               nCMF + nTXT, nCMF, nTXT, sources->n);
     149    psLogMsg ("psphot", 3, "%d external sources (%d cmf, %d text) merged to yield %ld total sources",
     150              nCMF + nTXT, nCMF, nTXT, sources->n);
    121151    return true;
    122152}
    … …  
    125155// XXX this function needs to be updated to work with the new context of pshot inputs
    126156psArray *psphotLoadPSFSources (pmConfig *config, const pmFPAview *view) {
     157
     158    bool status;
    127159
    128160    // find the currently selected readout
    … …  
    133165    }
    134166
    135     psArray *sources = psMetadataLookupPtr (NULL, readout->analysis, "PSPHOT.SOURCES");
     167    pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
     168    if (!detections) {
     169        psLogMsg ("psphot", 3, "no psf sources for this readout");
     170    }
     171
     172    psArray *sources = detections->allSources;
    136173    if (!sources) {
    137174        psLogMsg ("psphot", 3, "no psf sources for this readout");
    … …  
    141178}
    142179
     180bool psphotRepairLoadedSources (pmConfig *config, const pmFPAview *view) {
     181
     182    // find the currently selected readout
     183    pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", 0); // File of interest
     184    psAssert (file, "missing file?");
     185
     186    pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
     187    psAssert (readout, "missing readout?");
     188
     189    pmDetections *detections = psMetadataLookupPtr (NULL, readout->analysis, "PSPHOT.DETECTIONS");
     190    if (!detections) {
     191        psError(PSPHOT_ERR_CONFIG, false, "missing detections");
     192        return false;
     193    }
     194
     195    psArray *sources = detections->allSources;
     196    psAssert (sources, "missing sources?");
     197
     198    // peak flux is wrong : set based on previous image
     199    // use the peak measured in the moments analysis:
     200    for (int i = 0; i < sources->n; i++) {
     201      pmSource *source = sources->data[i];
     202      source->peak->flux = source->moments->Peak;
     203    }
     204
     205    return true;
     206}
     207
    143208// generate the detection structure for the supplied array of sources
    144 pmDetections *psphotDetectionsFromSources (pmConfig *config, psArray *sources) {
     209// XXX this currently assumes there is a single input file
     210bool psphotDetectionsFromSources (pmConfig *config, const pmFPAview *view, psArray *sources) {
     211
     212    // find the currently selected readout
     213    pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", 0); // File of interest
     214    psAssert (file, "missing file?");
     215
     216    pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
     217    psAssert (readout, "missing readout?");
    145218
    146219    pmDetections *detections = pmDetectionsAlloc();
    … …  
    155228    float snMin = psMetadataLookupF32(NULL, recipe, "MOMENTS_SN_MIN");
    156229    if (!isfinite(snMin)) {
    157         return NULL;
     230        return false;
    158231    }
    159232
    … …  
    166239        }
    167240
    168         // use the existing peak information, otherwise generate a new peak
    169         if (source->peak) {
    170             psArrayAdd (detections->peaks, 100, source->peak);
    171             continue;
    172         }
     241        // use the existing peak information, otherwise generate a new peak
     242        if (source->peak) {
     243            source->peak->assigned = false; // So the moments will be measured
     244            psArrayAdd (detections->peaks, 100, source->peak);
     245            continue;
     246        }
    173247
    174248        float flux = powf(10.0, -0.4 * source->psfMag);
    … …  
    181255        peak->flux = flux; // this are being set wrong, but does it matter?
    182256
    183         if (isfinite (source->errMag) && (source->errMag > 0.0)) {
    184           peak->SN = 1.0 / source->errMag;
    185         } else {
    186           peak->SN = 0.0;
    187         }
     257        if (isfinite (source->errMag) && (source->errMag > 0.0)) {
     258          peak->SN = 1.0 / source->errMag;
     259        } else {
     260          peak->SN = 0.0;
     261        }
    188262
    189263        psArrayAdd (detections->peaks, 100, peak);
    … …  
    193267    psLogMsg ("psphot", 3, "%ld PSF sources loaded", detections->peaks->n);
    194268
    195     return detections;
     269    // save detections on the readout->analysis
     270    if (!psMetadataAddPtr (readout->analysis, PS_LIST_TAIL, "PSPHOT.DETECTIONS", PS_META_REPLACE | PS_DATA_UNKNOWN, "psphot detectinos", detections)) {
     271        psError (PSPHOT_ERR_CONFIG, false, "problem saving detections on readout");
     272        return false;
     273    }
     274    psFree (detections);
     275    return true;
    196276}
    197277
    … …  
    222302        peak->flux = flux; // this are being set wrong, but does it matter?
    223303
    224         if (isfinite (source->errMag) && (source->errMag > 0.0)) {
    225           peak->SN = 1.0 / source->errMag;
    226         } else {
    227           peak->SN = 0.0;
    228         }
     304        if (isfinite (source->errMag) && (source->errMag > 0.0)) {
     305          peak->SN = 1.0 / source->errMag;
     306        } else {
     307          peak->SN = 0.0;
     308        }
    229309
    230310        source->peak = peak;
    … …  
    235315    return true;
    236316}
     317
     318bool psphotCheckExtSources (pmConfig *config, const pmFPAview *view) {
     319
     320    bool status;
     321
     322    psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE);
     323    psAssert (recipe, "missing recipe");
     324
     325    // find the currently selected readout
     326    pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", 0); // File of interest
     327    psAssert (file, "missing file?");
     328
     329    pmReadout *readout = pmFPAviewThisReadout(view, file->fpa);
     330    psAssert (readout, "missing readout?");
     331
     332    pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS");
     333    psAssert (detections, "missing detections?");
     334
     335    // XXX allSources of newSources?
     336    psArray *sources = detections->allSources;
     337    psAssert (sources, "missing sources?");
     338
     339    if (sources->n) {
     340        // the user wants to make the psf from these stars; define them as psf stars:
     341        for (int i = 0; i < sources->n; i++) {
     342            pmSource *source = sources->data[i];
     343            source->mode |= PM_SOURCE_MODE_PSFSTAR;
     344        }
     345        // force psphotChoosePSF to use all loaded sources
     346        psMetadataAddS32 (recipe, PS_LIST_TAIL, "PSF_MAX_NSTARS", PS_META_REPLACE, "max number of sources for PSF model", sources->n);
     347
     348        // measure stats of externally specified sources
     349        if (!psphotSourceStatsUpdate (sources, config, readout)) {
     350            psError(PSPHOT_ERR_CONFIG, false, "failure to measure stats of existing sources");
     351            return false;
     352        }
     353    } else {
     354
     355        // find the detections (by peak and/or footprint) in the image.
     356        if (!psphotFindDetections (config, view, true)) {
     357            psError(PSPHOT_ERR_CONFIG, false, "unable to find detections in this image");
     358            return psphotReadoutCleanup (config, view);
     359        }
     360
     361        // construct sources and measure basic stats
     362        psphotSourceStats (config, view, true);
     363
     364        // find blended neighbors of very saturated stars
     365        psphotDeblendSatstars (config, view);
     366
     367        // mark blended peaks PS_SOURCE_BLEND
     368        if (!psphotBasicDeblend (config, view)) {
     369            psLogMsg ("psphot", 3, "failed on deblend analysis");
     370            return psphotReadoutCleanup (config, view);
     371        }
     372
     373        // classify sources based on moments, brightness
     374        if (!psphotRoughClass (config, view)) {
     375            psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
     376            return psphotReadoutCleanup (config, view);
     377        }
     378    }
     379
     380    return true;
     381}

  • branches/eam_branches/psphot.stack.20100120/src/psphotModelBackground.c

    r26688 r26748  
    3232// generate the median in NxN boxes, clipping heavily
    3333// linear interpolation to generate full-scale model
     34//
     35// NOTE that the 'analysis' metedata pass in here is used to store the binning information.
     36// This may be the analysis for this readout, but it may be the analysis for the pmFPAfile
     37// corresponding to the model.  Other information about the background model is saved on the
     38// readout->analysis
    3439static bool psphotModelBackgroundReadout(psImage *model,  // Model image
    3540                                  psImage *modelStdev, // Model stdev image
    … …  
    140145
    141146    // we save the binning structure for use in psphotMagnitudes
    142     psMetadataAddPtr(analysis, PS_LIST_TAIL, "PSPHOT.BACKGROUND.BINNING",
    143                      PS_DATA_UNKNOWN | PS_META_REPLACE, "Background binning", binning);
     147    psMetadataAddPtr(analysis, PS_LIST_TAIL, "PSPHOT.BACKGROUND.BINNING", PS_DATA_UNKNOWN | PS_META_REPLACE, "Background binning", binning);
    144148
    145149    psF32 **modelData = model->data.F32;
    … …  
    296300    psLogMsg ("psphot", PS_LOG_INFO, "built background image: %f sec\n", psTimerMark ("psphot.background"));
    297301
    298     psMetadataAddF32(recipe, PS_LIST_TAIL, "SKY_MEAN", PS_META_REPLACE, "sky mean", Value);
    299     psMetadataAddF32(recipe, PS_LIST_TAIL, "SKY_STDEV", PS_META_REPLACE, "sky stdev", ValueStdev);
     302    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "SKY_MEAN", PS_META_REPLACE, "sky mean", Value);
     303    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "SKY_STDEV", PS_META_REPLACE, "sky stdev", ValueStdev);
    300304    psLogMsg ("psphot", PS_LOG_INFO, "image sky : mean %f stdev %f", Value, ValueStdev);
    301305
    … …  
    306310                                      PS_STAT_MAX);
    307311    psImageStats (statsBck, model, NULL, 0);
    308     psMetadataAddF32 (recipe, PS_LIST_TAIL, "MSKY_MN",
    309                       PS_META_REPLACE, "sky model mean",          statsBck->sampleMean);
    310     psMetadataAddF32 (recipe, PS_LIST_TAIL, "MSKY_SIG",
    311                       PS_META_REPLACE, "sky model stdev",         statsBck->sampleStdev);
    312     psMetadataAddF32 (recipe, PS_LIST_TAIL, "MSKY_MAX",
    313                       PS_META_REPLACE, "sky model maximum value", statsBck->max);
    314     psMetadataAddF32 (recipe, PS_LIST_TAIL, "MSKY_MIN",
    315                       PS_META_REPLACE, "sky model minimum value", statsBck->min);
    316     psMetadataAddS32 (recipe, PS_LIST_TAIL, "MSKY_NX",
    317                       PS_META_REPLACE, "sky model size (x)",      model->numCols);
    318     psMetadataAddS32 (recipe, PS_LIST_TAIL, "MSKY_NY",
    319                       PS_META_REPLACE, "sky model size (y)",      model->numRows);
     312    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "MSKY_MN", PS_META_REPLACE, "sky model mean",          statsBck->sampleMean);
     313    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "MSKY_SIG", PS_META_REPLACE, "sky model stdev",         statsBck->sampleStdev);
     314    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "MSKY_MAX", PS_META_REPLACE, "sky model maximum value", statsBck->max);
     315    psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "MSKY_MIN", PS_META_REPLACE, "sky model minimum value", statsBck->min);
     316    psMetadataAddS32 (readout->analysis, PS_LIST_TAIL, "MSKY_NX", PS_META_REPLACE, "sky model size (x)",      model->numCols);
     317    psMetadataAddS32 (readout->analysis, PS_LIST_TAIL, "MSKY_NY", PS_META_REPLACE, "sky model size (y)",      model->numRows);
    320318    psLogMsg ("psphot", PS_LOG_INFO, "background sky : min %f mean %f max %f stdev %f",
    321319              statsBck->min, statsBck->sampleMean, statsBck->max, statsBck->sampleStdev);

  • branches/eam_branches/psphot.stack.20100120/src/psphotOutput.c

    r26688 r26748  
    197197
    198198// these values are saved in an output header stub - they are added to either the PHU header
    199 // (CMP) or the MEF table header (CMF) XXX these are currently carried by the recipe -- this
    200 // will not work in a Stack context XXX they should be place in the readout->analysis of the
    201 // given image or directly on the output header
    202 psMetadata *psphotDefineHeader (psMetadata *recipe) {
     199// (CMP) or the MEF table header (CMF).  before the header is created, each readout has these
     200// values stored on readout->analysis
     201psMetadata *psphotDefineHeader (psMetadata *analysis) {
    203202
    204203    bool status = true;
    … …  
    207206
    208207    // write necessary information to output header
    209     psMetadataItemSupplement (&status, header, recipe, "ZERO_PT");
    210     psMetadataItemSupplement (&status, header, recipe, "PHOTCODE");
    211 
    212     psMetadataItemSupplement (&status, header, recipe, "APMIFIT");
    213     psMetadataItemSupplement (&status, header, recipe, "DAPMIFIT");
    214     psMetadataItemSupplement (&status, header, recipe, "NAPMIFIT");
     208    psMetadataItemSupplement (&status, header, analysis, "ZERO_PT");
     209    psMetadataItemSupplement (&status, header, analysis, "PHOTCODE");
     210
     211    psMetadataItemSupplement (&status, header, analysis, "APMIFIT");
     212    psMetadataItemSupplement (&status, header, analysis, "DAPMIFIT");
     213    psMetadataItemSupplement (&status, header, analysis, "NAPMIFIT");
    215214
    216215    // PSF model parameters (shape values for image center)
    217     psMetadataItemSupplement (&status, header, recipe, "NPSFSTAR");
    218     psMetadataItemSupplement (&status, header, recipe, "APLOSS");
    219 
    220     psMetadataItemSupplement (&status, header, recipe, "FWHM_MAJ");
    221     psMetadataItemSupplement (&status, header, recipe, "FW_MJ_SG");
    222     psMetadataItemSupplement (&status, header, recipe, "FW_MJ_LQ");
    223     psMetadataItemSupplement (&status, header, recipe, "FW_MJ_UQ");
    224 
    225     psMetadataItemSupplement (&status, header, recipe, "FWHM_MIN");
    226     psMetadataItemSupplement (&status, header, recipe, "FW_MN_SG");
    227     psMetadataItemSupplement (&status, header, recipe, "FW_MN_LQ");
    228     psMetadataItemSupplement (&status, header, recipe, "FW_MN_UQ");
    229 
    230     psMetadataItemSupplement (&status, header, recipe, "ANGLE");
     216    psMetadataItemSupplement (&status, header, analysis, "NPSFSTAR");
     217    psMetadataItemSupplement (&status, header, analysis, "APLOSS");
     218
     219    psMetadataItemSupplement (&status, header, analysis, "FWHM_MAJ");
     220    psMetadataItemSupplement (&status, header, analysis, "FW_MJ_SG");
     221    psMetadataItemSupplement (&status, header, analysis, "FW_MJ_LQ");
     222    psMetadataItemSupplement (&status, header, analysis, "FW_MJ_UQ");
     223
     224    psMetadataItemSupplement (&status, header, analysis, "FWHM_MIN");
     225    psMetadataItemSupplement (&status, header, analysis, "FW_MN_SG");
     226    psMetadataItemSupplement (&status, header, analysis, "FW_MN_LQ");
     227    psMetadataItemSupplement (&status, header, analysis, "FW_MN_UQ");
     228
     229    psMetadataItemSupplement (&status, header, analysis, "ANGLE");
    231230
    232231    // Image Quality measurements
    233     psMetadataItemSupplement (&status, header, recipe, "IQ_NSTAR");
    234 
    235     psMetadataItemSupplement (&status, header, recipe, "IQ_FW1");
    236     psMetadataItemSupplement (&status, header, recipe, "IQ_FW1_E");
    237     psMetadataItemSupplement (&status, header, recipe, "IQ_FW2");
    238     psMetadataItemSupplement (&status, header, recipe, "IQ_FW2_E");
    239 
    240     psMetadataItemSupplement (&status, header, recipe, "IQ_M2");
    241     psMetadataItemSupplement (&status, header, recipe, "IQ_M2_ER");
    242     psMetadataItemSupplement (&status, header, recipe, "IQ_M2_LQ");
    243     psMetadataItemSupplement (&status, header, recipe, "IQ_M2_UQ");
    244 
    245     psMetadataItemSupplement (&status, header, recipe, "IQ_M2C");
    246     psMetadataItemSupplement (&status, header, recipe, "IQ_M2C_E");
    247     psMetadataItemSupplement (&status, header, recipe, "IQ_M2C_L");
    248     psMetadataItemSupplement (&status, header, recipe, "IQ_M2C_U");
    249 
    250     psMetadataItemSupplement (&status, header, recipe, "IQ_M2S");
    251     psMetadataItemSupplement (&status, header, recipe, "IQ_M2S_E");
    252     psMetadataItemSupplement (&status, header, recipe, "IQ_M2S_L");
    253     psMetadataItemSupplement (&status, header, recipe, "IQ_M2S_U");
    254 
    255     psMetadataItemSupplement (&status, header, recipe, "IQ_M3");
    256     psMetadataItemSupplement (&status, header, recipe, "IQ_M3_ER");
    257     psMetadataItemSupplement (&status, header, recipe, "IQ_M3_LQ");
    258     psMetadataItemSupplement (&status, header, recipe, "IQ_M3_UQ");
    259 
    260     psMetadataItemSupplement (&status, header, recipe, "IQ_M4");
    261     psMetadataItemSupplement (&status, header, recipe, "IQ_M4_ER");
    262     psMetadataItemSupplement (&status, header, recipe, "IQ_M4_LQ");
    263     psMetadataItemSupplement (&status, header, recipe, "IQ_M4_UQ");
     232    psMetadataItemSupplement (&status, header, analysis, "IQ_NSTAR");
     233
     234    psMetadataItemSupplement (&status, header, analysis, "IQ_FW1");
     235    psMetadataItemSupplement (&status, header, analysis, "IQ_FW1_E");
     236    psMetadataItemSupplement (&status, header, analysis, "IQ_FW2");
     237    psMetadataItemSupplement (&status, header, analysis, "IQ_FW2_E");
     238
     239    psMetadataItemSupplement (&status, header, analysis, "IQ_M2");
     240    psMetadataItemSupplement (&status, header, analysis, "IQ_M2_ER");
     241    psMetadataItemSupplement (&status, header, analysis, "IQ_M2_LQ");
     242    psMetadataItemSupplement (&status, header, analysis, "IQ_M2_UQ");
     243
     244    psMetadataItemSupplement (&status, header, analysis, "IQ_M2C");
     245    psMetadataItemSupplement (&status, header, analysis, "IQ_M2C_E");
     246    psMetadataItemSupplement (&status, header, analysis, "IQ_M2C_L");
     247    psMetadataItemSupplement (&status, header, analysis, "IQ_M2C_U");
     248
     249    psMetadataItemSupplement (&status, header, analysis, "IQ_M2S");
     250    psMetadataItemSupplement (&status, header, analysis, "IQ_M2S_E");
     251    psMetadataItemSupplement (&status, header, analysis, "IQ_M2S_L");
     252    psMetadataItemSupplement (&status, header, analysis, "IQ_M2S_U");
     253
     254    psMetadataItemSupplement (&status, header, analysis, "IQ_M3");
     255    psMetadataItemSupplement (&status, header, analysis, "IQ_M3_ER");
     256    psMetadataItemSupplement (&status, header, analysis, "IQ_M3_LQ");
     257    psMetadataItemSupplement (&status, header, analysis, "IQ_M3_UQ");
     258
     259    psMetadataItemSupplement (&status, header, analysis, "IQ_M4");
     260    psMetadataItemSupplement (&status, header, analysis, "IQ_M4_ER");
     261    psMetadataItemSupplement (&status, header, analysis, "IQ_M4_LQ");
     262    psMetadataItemSupplement (&status, header, analysis, "IQ_M4_UQ");
    264263
    265264    // XXX these need to be defined from elsewhere
    266265    psMetadataAdd (header, PS_LIST_TAIL, "FSATUR",   PS_DATA_F32 | PS_META_REPLACE, "SATURATION MAG",      0.0);
    267266    psMetadataAdd (header, PS_LIST_TAIL, "FLIMIT",   PS_DATA_F32 | PS_META_REPLACE, "COMPLETENESS MAG",    0.0);
    268     psMetadataItemSupplement (&status, header, recipe, "NSTARS");
    269 
    270     psMetadataItemSupplement (&status, header, recipe, "NDET_EXT");
    271     psMetadataItemSupplement (&status, header, recipe, "NDET_CR");
     267    psMetadataItemSupplement (&status, header, analysis, "NSTARS");
     268
     269    psMetadataItemSupplement (&status, header, analysis, "NDET_EXT");
     270    psMetadataItemSupplement (&status, header, analysis, "NDET_CR");
    272271
    273272    // sky background model statistics
    274     psMetadataItemSupplement (&status, header, recipe, "MSKY_MN");
    275     psMetadataItemSupplement (&status, header, recipe, "MSKY_SIG");
    276     psMetadataItemSupplement (&status, header, recipe, "MSKY_MIN");
    277     psMetadataItemSupplement (&status, header, recipe, "MSKY_MAX");
    278     psMetadataItemSupplement (&status, header, recipe, "MSKY_NX");
    279     psMetadataItemSupplement (&status, header, recipe, "MSKY_NY");
     273    psMetadataItemSupplement (&status, header, analysis, "MSKY_MN");
     274    psMetadataItemSupplement (&status, header, analysis, "MSKY_SIG");
     275    psMetadataItemSupplement (&status, header, analysis, "MSKY_MIN");
     276    psMetadataItemSupplement (&status, header, analysis, "MSKY_MAX");
     277    psMetadataItemSupplement (&status, header, analysis, "MSKY_NX");
     278    psMetadataItemSupplement (&status, header, analysis, "MSKY_NY");
    280279
    281280    psMetadataAddF32 (header, PS_LIST_TAIL, "DT_PHOT", PS_META_REPLACE, "elapsed psphot time", psTimerMark ("psphotReadout"));

  • branches/eam_branches/psphot.stack.20100120/src/psphotRadiusChecks.c

    r25755 r26748  
    44static float PSF_FIT_NSIGMA;
    55static float PSF_FIT_PADDING;
     6static float PSF_APERTURE = 0;  // radius to use in PSF aperture mags
    67static float PSF_FIT_RADIUS = 0;        // radius to use in fitting (ignored if <= 0,
    78                                        // and a per-object radius is calculated)
    89
    9 static float PSF_APERTURE = 0;  // radius to use in PSF aperture mags
    10 
    11 
    12 bool psphotInitRadiusPSF(const psMetadata *recipe, const pmModelType type) {
     10bool psphotInitRadiusPSF(const psMetadata *recipe, const psMetadata *analysis, const pmModelType type) {
    1311
    1412    bool status = true;
    1513
    16     PSF_FIT_NSIGMA  = psMetadataLookupF32(&status, recipe, "PSF_FIT_NSIGMA");
     14    PSF_FIT_NSIGMA = psMetadataLookupF32(&status, recipe, "PSF_FIT_NSIGMA");
    1715    PSF_FIT_PADDING = psMetadataLookupF32(&status, recipe, "PSF_FIT_PADDING");
    18     PSF_FIT_RADIUS  =  psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS");
    19     PSF_APERTURE    =  psMetadataLookupF32(&status, recipe, "PSF_APERTURE");
     16
     17    PSF_FIT_RADIUS =  psMetadataLookupF32(&status, analysis, "PSF_FIT_RADIUS");
     18    if (!status) {
     19        PSF_FIT_RADIUS = psMetadataLookupF32(&status, recipe, "PSF_FIT_RADIUS");
     20    }
     21
     22    PSF_APERTURE =  psMetadataLookupF32(&status, analysis, "PSF_APERTURE");
     23    if (!status) {
     24        PSF_APERTURE =  psMetadataLookupF32(&status, recipe, "PSF_APERTURE");
     25    }
    2026
    2127    return true;

  • branches/eam_branches/psphot.stack.20100120/src/psphotReadout.c

    r26691 r26748  
    5050    }
    5151
    52     // run a single-model test if desired (exits from here if test is run)
    53     // XXX fix this in the new multi-input context (or drop this and only keep a stand-alone program)
    54     // psphotModelTest (config, view, recipe);
    55 
    56     // load the psf model, if suppled.  FWHM_X,FWHM_Y,etc are saved in the recipe
    57     // XXX the recipe is probably the wrong place to store those results...
    58     // XXX this is using the wrong file... fix in the new context
    59     // if (!psphotLoadPSF (config, view, recipe)) {
    60     //  // this only happens if we had a programming error in psphotLoadPSF
    61     //     psError (PSPHOT_ERR_UNKNOWN, false, "error loading psf model");
    62     //     return psphotReadoutCleanup (config, view);
    63     // }
     52    // load the psf model, if suppled.  FWHM_X,FWHM_Y,etc are determined and saved on
     53    // readout->analysis XXX this function currently only works with a single PSPHOT.INPUT
     54    if (!psphotLoadPSF (config, view)) {
     55        psError (PSPHOT_ERR_UNKNOWN, false, "error loading psf model");
     56        return psphotReadoutCleanup (config, view);
     57    }
    6458       
    6559    // find the detections (by peak and/or footprint) in the image.
    66     if (!psphotFindDetections (config, view)) { // pass 1
     60    if (!psphotFindDetections (config, view, true)) { // pass 1
    6761        // this only happens if we had an error in psphotFindDetections
    6862        psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis");
    … …  
    119113    // include externally-supplied sources
    120114    // XXX fix this in the new multi-input context
    121     // psphotLoadExtSources (config, view, sources); // pass 1
     115    psphotLoadExtSources (config, view); // pass 1
    122116
    123117    // construct an initial model for each object, set the radius to fitRadius, set circular
    … …  
    159153    // find fainter sources
    160154    // NOTE: finds new peaks and new footprints, OLD and FULL set are saved on detections
    161     psphotFindDetections (config, view); // pass 2 (detections->peaks, detections->footprints)
     155    psphotFindDetections (config, view, false); // pass 2 (detections->peaks, detections->footprints)
    162156
    163157    // remove noise for subtracted objects (ie, return to normal noise level)

  • branches/eam_branches/psphot.stack.20100120/src/psphotReadoutCleanup.c

    r26691 r26748  
    5858    // use the psf-model to measure FWHM stats
    5959    if (psf) {
    60         if (!psphotPSFstats (readout, recipe, psf)) {
     60        if (!psphotPSFstats (readout, psf)) {
    6161            psError(PSPHOT_ERR_PROG, false, "Failed to measure PSF shape parameters");
    6262            return false;
    … …  
    6565    // otherwise, use the source moments to measure FWHM stats
    6666    if (!psf && sources) {
    67         if (!psphotMomentsStats (readout, recipe, sources)) {
     67        if (!psphotMomentsStats (readout, sources)) {
    6868            psError(PSPHOT_ERR_PROG, false, "Failed to measure Moment shape parameters");
    6969            return false;
    … …  
    8282    }
    8383
    84     // create an output header with stats results XXX this needs to be reworked :
    85     // psphotImageQuality and others need to write these values to a header.  that needs to be
    86     // stored on the readout->analysis
    87     psMetadata *header = psphotDefineHeader (recipe);
     84    // create an output header with stats results currently saved on readout->analysis
     85    psMetadata *header = psphotDefineHeader (readout->analysis);
    8886
    8987    // write NSTARS to the image header

  • branches/eam_branches/psphot.stack.20100120/src/psphotReadoutFindPSF.c

    r26611 r26748  
    77    psTimerStart ("psphotReadout");
    88
    9     // select the current recipe
    10     psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE);
    11     if (!recipe) {
    12         psError(PSPHOT_ERR_CONFIG, false, "missing recipe %s", PSPHOT_RECIPE);
    13         return false;
    14     }
    15 
    169    // set the photcode for the PSPHOT.INPUT
    1710    if (!psphotAddPhotcode(config, view)) {
    … …  
    2013    }
    2114
    22     // find the currently selected readout
    23     pmReadout  *readout = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT");
    24     PS_ASSERT_PTR_NON_NULL (readout, false);
    25 
    2615    // Generate the mask and variance images, including the user-defined analysis region of interest
    27     psphotSetMaskAndVariance (config, view, recipe);
    28 
    29     // display the image, weight, mask (ch 1,2,3)
    30     psphotVisualShowImage (readout);
     16    psphotSetMaskAndVariance (config, view);
    3117
    3218    // Note that in this implementation, we do NOT model the background and we do not
    … …  
    3420
    3521    // include externally-supplied sources (supplied as PSPHOT.INPUT.CMF)
    36     pmDetections *detections = psphotDetectionsFromSources (config, inSources);
    37     if (!detections || !detections->peaks) {
     22    // XXX we assume a single set of input sources is supplied
     23    if (!psphotDetectionsFromSources (config, view, inSources)) {
    3824        psError(PSPHOT_ERR_ARGUMENTS, true, "Can't find PSF stars");
    39         return psphotReadoutCleanup(config, readout, recipe, detections, NULL, NULL);
     25        return psphotReadoutCleanup(config, view);
    4026    }
    4127
    4228    // construct sources and measure basic stats (moments, local sky)
    43     psArray *sources = psphotSourceStats(config, readout, detections, true);
    44     if (!sources) return false;
     29    if (!psphotSourceStats(config, view, true)) {
     30        psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources");
     31        return false;
     32    }
    4533
    46     // peak flux is wrong : set based on previous image
    47     // use the peak measured in the moments analysis:
    48     for (int i = 0; i < sources->n; i++) {
    49       pmSource *source = sources->data[i];
    50       source->peak->flux = source->moments->Peak;
     34    // peak flux is wrong : use the peak measured in the moments analysis:
     35    if (!psphotRepairLoadedSources(config, view)) {
     36        psError(PSPHOT_ERR_UNKNOWN, false, "failure to repair sources");
     37        return false;
    5138    }
    5239
    5340    // classify sources based on moments, brightness (psf is not known)
    54     if (!psphotRoughClass (readout, sources, recipe, false)) {
    55         psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
    56         return psphotReadoutCleanup (config, readout, recipe, detections, NULL, sources);
     41    if (!psphotRoughClass (config, view)) {
     42        psError (PSPHOT_ERR_UNKNOWN, false, "failed to determine rough source class");
     43        return psphotReadoutCleanup (config, view);
    5744    }
    5845
    59     if (!psphotImageQuality (recipe, sources)) {
    60         psLogMsg("psphot", 3, "failed to measure image quality");
    61         return psphotReadoutCleanup(config, readout, recipe, detections, NULL, sources);
     46    if (!psphotImageQuality (config, view)) {
     47        psError (PSPHOT_ERR_UNKNOWN, false, "failed to measure image quality");
     48        return psphotReadoutCleanup(config, view);
    6249    }
    6350
    64     pmPSF *psf = psphotChoosePSF(readout, sources, recipe);
    65     if (!psf) {
     51    if (!psphotChoosePSF(config, view)) {
    6652        psError(PSPHOT_ERR_PSF, false, "Failed to construct a psf model");
    67         psFree(sources);
    68         return psphotReadoutCleanup(config, readout, recipe, detections, NULL, NULL);
     53        return psphotReadoutCleanup(config, view);
    6954    }
    70     psphotVisualShowPSFModel(readout, psf);
    7155
    7256# if 0
    … …  
    7559    // fits from that analysis, or run the linear PSF fit for all objects currently in hand
    7660    // construct an initial model for each object, set the radius to fitRadius, set circular fit mask
    77     psphotGuessModels (config, readout, sources, psf);
     61    psphotGuessModels (config, view);
    7862
    7963    // measure aperture photometry corrections
    80     if (!psphotApResid (config, readout, sources, psf)) {
     64    if (!psphotApResid (config, view)) {
    8165        psLogMsg ("psphot", 3, "failed on psphotApResid");
    82         return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
     66        return psphotReadoutCleanup (config, view);
    8367    }
    8468# endif
    8569
    8670    // drop the references to the image pixels held by each source
    87     psphotSourceFreePixels(sources);
    88     psFree(sources);
     71    psphotSourceFreePixels(config, view);
    8972
    9073    // create the exported-metadata and free local data
    91     return psphotReadoutCleanup(config, readout, recipe, detections, psf, NULL);
     74    return psphotReadoutCleanup(config, view);
    9275}

  • branches/eam_branches/psphot.stack.20100120/src/psphotReadoutKnownSources.c

    r26542 r26748  
    11# include "psphotInternal.h"
    22
    3 // in this psphotReadout-variant, we are only measuring the photometry for known sources,
    4 // using a PSF generated from this observation those sources
     3// in this psphotReadout-variant, we are only measuring the photometry for known sources, using
     4// a PSF generated for this observation from those sources
    55bool psphotReadoutKnownSources(pmConfig *config, const pmFPAview *view, psArray *inSources) {
    66
    77    psTimerStart ("psphotReadout");
    8 
    9     // select the current recipe
    10     psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE);
    11     if (!recipe) {
    12         psError(PSPHOT_ERR_CONFIG, false, "missing recipe %s", PSPHOT_RECIPE);
    13         return false;
    14     }
    158
    169    // set the photcode for this image
    … …  
    2013    }
    2114
    22     // find the currently selected readout
    23     pmReadout  *readout = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT");
    24     PS_ASSERT_PTR_NON_NULL (readout, false);
    25 
    2615    // Generate the mask and weight images, including the user-defined analysis region of interest
    27     psphotSetMaskAndVariance (config, view, recipe);
    28 
    29     // display the image, weight, mask (ch 1,2,3)
    30     psphotVisualShowImage (readout);
     16    psphotSetMaskAndVariance (config, view);
    3117
    3218    // Note that in this implementation, we do NOT model the background and we do not
    … …  
    3420
    3521    // include externally-supplied sources (supplied as PSPHOT.INPUT.CMF)
    36     pmDetections *detections = psphotDetectionsFromSources (config, inSources);
    37     if (!detections || !detections->peaks) {
     22    if (!psphotDetectionsFromSources (config, view, inSources)) {
    3823        psError(PSPHOT_ERR_ARGUMENTS, true, "Can't find PSF stars");
    39         return psphotReadoutCleanup(config, readout, recipe, detections, NULL, NULL);
     24        return psphotReadoutCleanup(config, view);
    4025    }
    4126
    4227    // construct sources and measure basic stats
    43     psArray *sources = psphotSourceStats (config, readout, detections, true);
    44     if (!sources) return false;
     28    if (!psphotSourceStats (config, view, true)) {
     29        psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources");
     30        return false;
     31    }
    4532
    46     // peak flux is wrong : set based on previous image
    47     // use the peak measured in the moments analysis:
    48     for (int i = 0; i < sources->n; i++) {
    49         pmSource *source = sources->data[i];
    50         source->peak->flux = source->moments->Peak;
     33    // peak flux is wrong : use the peak measured in the moments analysis:
     34    if (!psphotRepairLoadedSources(config, view)) {
     35        psError(PSPHOT_ERR_UNKNOWN, false, "failure to repair sources");
     36        return false;
    5137    }
    5238
    5339    // classify sources based on moments, brightness (psf is not known)
    54     if (!psphotRoughClass (readout, sources, recipe, false)) {
    55         psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
    56         return psphotReadoutCleanup (config, readout, recipe, detections, NULL, sources);
     40    if (!psphotRoughClass (config, view)) {
     41        psError (PSPHOT_ERR_UNKNOWN, false, "failed to determine rough source class");
     42        return psphotReadoutCleanup(config, view);
    5743    }
    5844
    59     pmPSF *psf = psphotChoosePSF (readout, sources, recipe);
    60     if (!psf) {
    61         psLogMsg ("psphot", 3, "failure to construct a psf model");
    62         return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
     45    if (!psphotChoosePSF (config, view)) {
     46        psError(PSPHOT_ERR_PSF, false, "Failed to construct a psf model");
     47        return psphotReadoutCleanup(config, view);
    6348    }
    64     psphotVisualShowPSFModel (readout, psf);
    6549
    6650    // construct an initial model for each object
    67     psphotGuessModels (config, readout, sources, psf);
     51    psphotGuessModels (config, view);
    6852
    6953    // linear PSF fit to source peaks
    70     psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE);
    71 
    72     // We have to place these visualizations here because the models are not realized until
    73     // psphotGuessModels or fitted until psphotFitSourcesLinear.
    74     psphotVisualShowPSFStars (recipe, psf, sources);
    75     psphotVisualShowSatStars (recipe, psf, sources);
     54    psphotFitSourcesLinear (config, view, false);
    7655
    7756    // measure aperture photometry corrections
    78     if (!psphotApResid (config, readout, sources, psf)) {
     57    if (!psphotApResid (config, view)) {
    7958        psLogMsg ("psphot", 3, "failed on psphotApResid");
    80         return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
     59        return psphotReadoutCleanup(config, view);
    8160    }
    8261
    8362    // calculate source magnitudes
    84     psphotMagnitudes(config, readout, view, sources, psf);
     63    psphotMagnitudes(config, view);
    8564
    8665    // drop the references to the image pixels held by each source
    87     psphotSourceFreePixels (sources);
     66    psphotSourceFreePixels (config, view);
    8867
    8968    // create the exported-metadata and free local data
    90     return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
     69    return psphotReadoutCleanup(config, view);
    9170}

  • branches/eam_branches/psphot.stack.20100120/src/psphotReadoutMinimal.c

    r26596 r26748  
    1717    pmModelClassSetLimits(PM_MODEL_LIMITS_LAX);
    1818
    19     // select the current recipe
    20     psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE);
    21     if (!recipe) {
    22         psError(PSPHOT_ERR_CONFIG, false, "missing recipe %s", PSPHOT_RECIPE);
    23         return false;
    24     }
    25 
    2619    // set the photcode for this image
    2720    if (!psphotAddPhotcode(config, view)) {
    … …  
    3023    }
    3124
    32     // find the currently selected readout
    33     pmReadout  *readout = pmFPAfileThisReadout (config->files, view, "PSPHOT.INPUT");
    34     PS_ASSERT_PTR_NON_NULL (readout, false);
     25    // Generate the mask and weight images, including the user-defined analysis region of interest
     26    psphotSetMaskAndVariance (config, view);
    3527
    36     // Generate the mask and weight images, including the user-defined analysis region of interest
    37     psphotSetMaskAndVariance (config, view, recipe);
    38 
    39     // display the image, weight, mask (ch 1,2,3)
    40     psphotVisualShowImage (readout);
    41 
    42     // load the psf model, if suppled.  FWHM_X,FWHM_Y,etc are saved in the recipe
    43     pmPSF *psf = psphotLoadPSF (config, view, recipe);
    44     if (!psf) {
     28    // load the psf model, if suppled.  FWHM_X,FWHM_Y,etc are saved on readout->analysis
     29    if (!psphotLoadPSF (config, view)) {
    4530      psError (PSPHOT_ERR_CONFIG, false, "missing psf model");
    46       return psphotReadoutCleanup (config, readout, recipe, NULL, NULL, NULL);
     31      return psphotReadoutCleanup (config, view);
    4732    }
    4833
    49     // find the detections (by peak and/or footprint) in the image.
    50     pmDetections *detections = pmDetectionsAlloc(); // New detections; allocated to ensure pass=2
    51     detections = psphotFindDetections(detections, readout, recipe);
    52     if (!detections) {
     34    // find the detections (by peak and/or footprint) in the image. (final pass)
     35    if (!psphotFindDetections(config, view, false)) {
    5336        psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis");
    54         return psphotReadoutCleanup (config, readout, recipe, detections, psf, NULL);
    55     }
    56     if (!detections->peaks->n) {
    57         psLogMsg ("psphot", 3, "unable to find detections in this image");
    58         return psphotReadoutCleanup (config, readout, recipe, detections, psf, NULL);
     37        return psphotReadoutCleanup (config, view);
    5938    }
    6039
    61     // construct sources and measure basic stats
    62     psArray *sources = psphotSourceStats (config, readout, detections, false);
    63     if (!sources) return false;
     40    // construct sources and measure basic stats (saved on detections->newSources)
     41    if (!psphotSourceStats (config, view, false)) { // pass 1
     42        psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources");
     43        return psphotReadoutCleanup (config, view);
     44    }
    6445
    6546    // find blended neighbors of very saturated stars
    66     // XXX merge this with Basic Deblend?
    67     psphotDeblendSatstars (readout, sources, recipe);
     47    psphotDeblendSatstars (config, view);
    6848
    6949    // mark blended peaks PS_SOURCE_BLEND
    70     if (!psphotBasicDeblend (sources, recipe)) {
     50    if (!psphotBasicDeblend (config, view)) {
    7151        psLogMsg ("psphot", 3, "failed on deblend analysis");
    72         return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
     52        return psphotReadoutCleanup (config, view);
    7353    }
    7454
    7555    // classify sources based on moments, brightness (use supplied psf shape parameters)
    76     if (!psphotRoughClass (readout, sources, recipe, true)) {
     56    if (!psphotRoughClass (config, view)) {
    7757        psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image");
    78         return psphotReadoutCleanup (config, readout, recipe, detections, psf, sources);
     58        return psphotReadoutCleanup (config, view);
    7959    }
    8060
    8161    // construct an initial model for each object
    82     psphotGuessModels (config, readout, sources, psf);
     62    psphotGuessModels (config, view);
    8363
    8464    // linear PSF fit to source peaks
    85     psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE);
    86 
    87     // We have to place these visualizations here because the models are not realized until
    88     // psphotGuessModels or fitted until psphotFitSourcesLinear.
    89     psphotVisualShowPSFStars (recipe, psf, sources);
    90     psphotVisualShowSatStars (recipe, psf, sources);
     65    psphotFitSourcesLinear (config, view, false);
    9166
    9267// XXX eventually, add the extended source fits here
    9368# if (0)
    9469    // measure source size for the remaining sources
    95     psphotSourceSize (config, readout, sources, recipe, 0);
     70    psphotSourceSize (config, view);
    9671
    97     psphotExtendedSourceAnalysis (readout, sources, recipe);
     72    psphotExtendedSourceAnalysis (config, view);
    9873
    99     psphotExtendedSourceFits (readout, sources, recipe);
     74    psphotExtendedSourceFits (config, view);
    10075# endif
    10176
    10277    // calculate source magnitudes
    103     psphotMagnitudes(config, readout, view, sources, psf);
     78    psphotMagnitudes(config, view);
    10479
    105     if (!psphotEfficiency(config, readout, view, psf, recipe, sources)) {
     80    // XXX ensure this is measured if the analysis succeeds (even if quality is low)
     81    if (!psphotEfficiency(config, view)) {
    10682        psErrorStackPrint(stderr, "Unable to determine detection efficiencies from fake sources");
    10783        psErrorClear();
    … …  
    10985
    11086    // drop the references to the image pixels held by each source
    111     psphotSourceFreePixels (sources);
     87    psphotSourceFreePixels (config, view);
    11288
    11389    // create the exported-metadata and free local data
    114     return psphotReadoutCleanup(config, readout, recipe, detections, psf, sources);
     90    return psphotReadoutCleanup(config, view);
    11591}

  • branches/eam_branches/psphot.stack.20100120/src/psphotRoughClass.c

    r26691 r26748  
    118118        // determine the PSF parameters from the source moment values
    119119        // XXX why not save the psfClump as a PTR?
    120         psfClump = pmSourcePSFClump (region, sources, recipe);
     120
     121        float PSF_SN_LIM = psMetadataLookupF32(&status, recipe, "PSF_SN_LIM"); psAssert (status, "missing PSF_SN_LIM");
     122        float MOMENTS_AR_MAX = psMetadataLookupF32(&status, recipe, "MOMENTS_AR_MAX"); psAssert (status, "missing MOMENTS_AR_MAX");
     123
     124        float PSF_CLUMP_GRID_SCALE = psMetadataLookupF32(&status, analysis, "PSF_CLUMP_GRID_SCALE");
     125        if (!status) {
     126            PSF_CLUMP_GRID_SCALE = psMetadataLookupF32(&status, recipe, "PSF_CLUMP_GRID_SCALE");
     127            psAssert (status, "missing PSF_CLUMP_GRID_SCALE");
     128        }
     129        float MOMENTS_SX_MAX = psMetadataLookupF32(&status, analysis, "MOMENTS_SX_MAX");
     130        if (!status) {
     131            MOMENTS_SX_MAX = psMetadataLookupF32(&status, recipe, "MOMENTS_SX_MAX");
     132            psAssert (status, "missing MOMENTS_SX_MAX");
     133        }
     134        float MOMENTS_SY_MAX = psMetadataLookupF32(&status, analysis, "MOMENTS_SY_MAX");
     135        if (!status) {
     136            MOMENTS_SY_MAX = psMetadataLookupF32(&status, recipe, "MOMENTS_SY_MAX");
     137            psAssert (status, "missing MOMENTS_SY_MAX");
     138        }
     139
     140        psfClump = pmSourcePSFClump (NULL, region, sources, PSF_SN_LIM, PSF_CLUMP_GRID_SCALE, MOMENTS_SX_MAX, MOMENTS_SY_MAX, MOMENTS_AR_MAX);
     141
    121142        psMetadataAddF32 (regionMD, PS_LIST_TAIL, "PSF.CLUMP.X",  PS_META_REPLACE, "psf clump center", psfClump.X);
    122143        psMetadataAddF32 (regionMD, PS_LIST_TAIL, "PSF.CLUMP.Y",  PS_META_REPLACE, "psf clump center", psfClump.Y);
    … …  
    146167    psLogMsg ("psphot", 3, "psf clump DX, DY: %f, %f\n", psfClump.dX, psfClump.dY);
    147168
     169    // get basic parameters, or set defaults
     170    float PSF_SN_LIM = psMetadataLookupF32 (&status, recipe, "PSF_SN_LIM"); psAssert (status, "missing PSF_SN_LIM");
     171    float PSF_CLUMP_NSIGMA = psMetadataLookupF32 (&status, recipe, "PSF_CLUMP_NSIGMA"); psAssert (status, "missing PSF_CLUMP_NSIGMA");
     172
    148173    // group into STAR, COSMIC, EXTENDED, SATURATED, etc.
    149     if (!pmSourceRoughClass (region, sources, recipe, psfClump, maskSat)) {
     174    if (!pmSourceRoughClass (region, sources, PSF_SN_LIM, PSF_CLUMP_NSIGMA, psfClump, maskSat)) {
    150175        psError(PSPHOT_ERR_PROG, false, "programming error calling pmSourceRoughClass");
    151176        return false;

  • branches/eam_branches/psphot.stack.20100120/src/psphotSetThreads.c

    r25755 r26748  
    2525    psFree(task);
    2626
    27     task = psThreadTaskAlloc("PSPHOT_SOURCE_STATS", 5);
     27    task = psThreadTaskAlloc("PSPHOT_SOURCE_STATS", 10);
    2828    task->function = &psphotSourceStats_Threaded;
    2929    psThreadTaskAdd(task);

  • branches/eam_branches/psphot.stack.20100120/src/psphotSignificanceImage.c

    r25383 r26748  
    2222    }
    2323
     24    // if we have already determined the PSF model, then we have a better idea how to smooth this image
    2425    bool statusMajor, statusMinor;
    25     float fwhmMajor = psMetadataLookupF32(&statusMajor, recipe, "FWHM_MAJ");
    26     float fwhmMinor = psMetadataLookupF32(&statusMinor, recipe, "FWHM_MIN");
     26    float fwhmMajor = psMetadataLookupF32(&statusMajor, readout->analysis, "FWHM_MAJ");
     27    float fwhmMinor = psMetadataLookupF32(&statusMinor, readout->analysis, "FWHM_MIN");
    2728    if (statusMajor && statusMinor) {
    2829        // if we know the FHWM, use that to set the smoothing kernel (XXX allow an optional override?)
    … …  
    4344    }
    4445    // record the actual smoothing sigma
    45     psMetadataAddF32(recipe, PS_LIST_TAIL, "SIGMA_SMOOTH", PS_META_REPLACE, "Smoothing sigma for detections",
    46                      SIGMA_SMTH);
     46    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "SIGMA_SMOOTH", PS_META_REPLACE, "Smoothing sigma for detections", SIGMA_SMTH);
    4747
    4848    // smooth the image, applying the mask as we go
    … …  
    5959    // effective per-pixel variance is maintained.
    6060    psImage *smooth_wt = psImageCopy(NULL, readout->variance, PS_TYPE_F32);
    61     psImageSmoothMask_Threaded(smooth_wt, smooth_wt, readout->mask, maskVal, SIGMA_SMTH * M_SQRT1_2,
    62                       NSIGMA_SMTH, minGauss);
     61    psImageSmoothMask_Threaded(smooth_wt, smooth_wt, readout->mask, maskVal, SIGMA_SMTH * M_SQRT1_2, NSIGMA_SMTH, minGauss);
    6362    psLogMsg("psphot", PS_LOG_MINUTIA, "smooth variance: %f sec\n", psTimerMark("psphot.smooth"));
    6463
    … …  
    7776    // Calculate correction factor for the covariance produced by the (potentially multiple) smoothing
    7877    psKernel *kernel = psImageSmoothKernel(SIGMA_SMTH, NSIGMA_SMTH); // Kernel used for smoothing
    79     psKernel *covar = psImageCovarianceCalculate(kernel, readout->covariance); // Covariance matrix
     78    float factor = 1.0 / psImageCovarianceCalculateFactor(kernel, readout->covariance);
    8079    psFree(kernel);
    81     float factor = 1.0 / psImageCovarianceFactor(covar);
    82     psFree(covar);
    8380
    8481    // record the effective area and significance scaling factor
    8582    float effArea = 8.0 * M_PI * PS_SQR(SIGMA_SMTH);
    86     psMetadataAddF32(recipe, PS_LIST_TAIL, "EFFECTIVE_AREA", PS_META_REPLACE, "Effective Area", effArea);
    87     psMetadataAddF32(recipe, PS_LIST_TAIL, "SIGNIFICANCE_SCALE_FACTOR", PS_META_REPLACE,
    88                      "Signicance scale factor", factor);
     83    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "EFFECTIVE_AREA", PS_META_REPLACE, "Effective Area", effArea);
     84    psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "SIGNIFICANCE_SCALE_FACTOR", PS_META_REPLACE, "Signicance scale factor", factor);
    8985
    9086    // XXX multithread this if needed

  • branches/eam_branches/psphot.stack.20100120/src/psphotSourceSize.c

    r26691 r26748  
    179179        float dMag = source->psfMag - apMag;
    180180
    181         psVectorAppend (Ap, 100, dMag);
    182         psVectorAppend (ApErr, 100, source->errMag);
     181        psVectorAppend (Ap, dMag);
     182        psVectorAppend (ApErr, source->errMag);
    183183    }
    184184

  • branches/eam_branches/psphot.stack.20100120/src/psphotSourceStats.c

    r26691 r26748  
    4242    psAssert (!detections->newSources, "new sources already defined?");
    4343
     44    // XXX TEST:
     45    if (detections->allSources) {
     46        psphotMaskCosmicRayFootprintCheck(detections->allSources);
     47    }
     48    if (detections->newSources) {
     49        psphotMaskCosmicRayFootprintCheck(detections->newSources);
     50    }
     51
    4452    // determine the number of allowed threads
    4553    int nThreads = psMetadataLookupS32(&status, config->arguments, "NTHREADS"); // Number of threads
    … …  
    5765    char *breakPt  = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
    5866    psAssert (status, "missing BREAK_POINT?");
     67
     68    float INNER        = psMetadataLookupF32 (&status, recipe, "SKY_INNER_RADIUS"); psAssert (status, "missing SKY_INNER_RADIUS");
     69    float MIN_SN       = psMetadataLookupF32 (&status, recipe, "MOMENTS_SN_MIN"); psAssert (status, "missing MOMENTS_SN_MIN");
     70
     71    // bit-masks to test for good/bad pixels
     72    psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT");
     73    psAssert (maskVal, "missing MASK.PSPHOT");
     74
     75    // bit-mask to mark pixels not used in analysis
     76    psImageMaskType markVal = psMetadataLookupImageMask(&status, recipe, "MARK.PSPHOT");
     77    psAssert (markVal, "missing MARK.PSPHOT");
    5978
    6079    psArray *peaks = detections->peaks;
    … …  
    112131    }
    113132
     133    // if we have measured the window, we will be saving the modified version of these recipe values on readout->analysis
     134    float SIGMA = psMetadataLookupF32 (&status, readout->analysis, "MOMENTS_GAUSS_SIGMA");
     135    if (!status) {
     136        SIGMA = psMetadataLookupF32 (&status, recipe, "MOMENTS_GAUSS_SIGMA");
     137    }
     138    float RADIUS = psMetadataLookupF32 (&status, readout->analysis, "PSF_MOMENTS_RADIUS");
     139    if (!status) {
     140        RADIUS = psMetadataLookupF32 (&status, recipe, "PSF_MOMENTS_RADIUS");
     141    }
     142
    114143    // threaded measurement of the source magnitudes
    115144    int Nfail = 0;
    … …  
    133162
    134163            psArrayAdd(job->args, 1, cells->data[j]); // sources
    135             psArrayAdd(job->args, 1, recipe);
     164
     165            PS_ARRAY_ADD_SCALAR(job->args, INNER,   PS_TYPE_F32);
     166            PS_ARRAY_ADD_SCALAR(job->args, MIN_SN,  PS_TYPE_F32);
     167            PS_ARRAY_ADD_SCALAR(job->args, RADIUS,  PS_TYPE_F32);
     168            PS_ARRAY_ADD_SCALAR(job->args, SIGMA,   PS_TYPE_F32);
     169
     170            PS_ARRAY_ADD_SCALAR(job->args, maskVal, PS_TYPE_IMAGE_MASK);
     171            PS_ARRAY_ADD_SCALAR(job->args, markVal, PS_TYPE_IMAGE_MASK);
     172
    136173            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nmoments
    137174            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfail
    … …  
    152189            psFree(sources);
    153190            return false;
     191        }
     192
     193        // we have only supplied one type of job, so we can assume the types here
     194        psThreadJob *job = NULL;
     195        while ((job = psThreadJobGetDone()) != NULL) {
     196            if (job->args->n < 1) {
     197                fprintf (stderr, "error with job\n");
     198            } else {
     199                psScalar *scalar = NULL;
     200                scalar = job->args->data[7];
     201                Nmoments += scalar->data.S32;
     202                scalar = job->args->data[8];
     203                Nfail += scalar->data.S32;
     204                scalar = job->args->data[9];
     205                Nfaint += scalar->data.S32;
     206            }
     207            psFree(job);
     208        }
     209    }
     210
     211    psFree (cellGroups);
     212
     213    psLogMsg ("psphot", PS_LOG_INFO, "%ld sources, %d moments, %d faint, %d failed: %f sec\n", sources->n, Nmoments, Nfaint, Nfail, psTimerMark ("psphot.stats"));
     214
     215    psphotVisualShowMoments (sources);
     216
     217    // save the new sources on the detection structure:
     218    detections->newSources = sources;
     219
     220
     221    if (detections->allSources) {
     222        psphotMaskCosmicRayFootprintCheck(detections->allSources);
     223    }
     224    if (detections->newSources) {
     225        psphotMaskCosmicRayFootprintCheck(detections->newSources);
     226    }
     227
     228    return true;
     229}
     230
     231// this function is currently only called by psphotCheckExtSources
     232bool psphotSourceStatsUpdate (psArray *sources, pmConfig *config, pmReadout *readout) {
     233
     234    bool status = false;
     235
     236    psTimerStart ("psphot.stats");
     237
     238    // select the appropriate recipe information
     239    psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
     240    assert (recipe);
     241
     242    // determine the number of allowed threads
     243    int nThreads = psMetadataLookupS32(&status, config->arguments, "NTHREADS"); // Number of threads
     244    if (!status) {
     245        nThreads = 0;
     246    }
     247
     248    // determine properties (sky, moments) of initial sources
     249    float OUTER    = psMetadataLookupF32 (&status, recipe, "SKY_OUTER_RADIUS");
     250    if (!status) return false;
     251
     252    OUTER = PS_MAX(OUTER, 20.0); // XXX Guarantee that we can encompass the max moments radius
     253
     254    char *breakPt  = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
     255    if (!status) return false;
     256
     257    for (int i = 0; i < sources->n; i++) {
     258
     259        pmSource *source = sources->data[i];
     260        if (!source->peak) continue; // XXX how can we have a peak-less source?
     261
     262        // allocate space for moments
     263        if (!source->moments) {
     264            source->moments = pmMomentsAlloc();
     265        }
     266
     267        // allocate image, weight, mask arrays for each peak (square of radius OUTER)
     268        pmSourceDefinePixels (source, readout, source->peak->x, source->peak->y, OUTER);
     269        source->peak->assigned = true;
     270    }
     271
     272    if (!strcasecmp (breakPt, "PEAKS")) {
     273        psLogMsg ("psphot", PS_LOG_INFO, "%ld sources : %f sec\n", sources->n, psTimerMark ("psphot.stats"));
     274        psLogMsg ("psphot", PS_LOG_INFO, "break point PEAKS, skipping MOMENTS\n");
     275        psphotVisualShowMoments (sources);
     276        return sources;
     277    }
     278
     279    // XXX how else could we get the window size in?
     280    if (!psphotSetMomentsWindow(recipe, readout->analysis, sources)) {
     281        psError(PS_ERR_UNEXPECTED_NULL, false, "Failed to determine Moments Window!");
     282        return NULL;
     283    }
     284
     285    // threaded measurement of the source magnitudes
     286    int Nfail = 0;
     287    int Nmoments = 0;
     288    int Nfaint = 0;
     289
     290    // choose Cx, Cy (see psphotThreadTools.c for overview of the concepts)
     291    int Cx = 1, Cy = 1;
     292    psphotChooseCellSizes (&Cx, &Cy, readout, nThreads);
     293
     294    psArray *cellGroups = psphotAssignSources (Cx, Cy, sources);
     295
     296    for (int i = 0; i < cellGroups->n; i++) {
     297
     298        psArray *cells = cellGroups->data[i];
     299
     300        for (int j = 0; j < cells->n; j++) {
     301
     302            // allocate a job -- if threads are not defined, this just runs the job
     303            psThreadJob *job = psThreadJobAlloc ("PSPHOT_SOURCE_STATS");
     304
     305            // XXX: this must match the above
     306            psArrayAdd(job->args, 1, cells->data[j]); // sources
     307            psArrayAdd(job->args, 1, recipe);
     308            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nmoments
     309            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfail
     310            PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfaint
     311
     312            if (!psThreadJobAddPending(job)) {
     313                psError(PS_ERR_UNKNOWN, false, "Unable to launch thread job PSPHOT_SOURCE_STATS");
     314                psFree (job);
     315                return NULL;
     316            }
     317            psFree(job);
     318        }
     319
     320        // wait for the threads to finish and manage results
     321        if (!psThreadPoolWait (false)) {
     322            psError(PS_ERR_UNKNOWN, false, "Failure in thread job PSPHOT_SOURCE_STATS");
     323            return NULL;
    154324        }
    155325
    … …  
    178348    psphotVisualShowMoments (sources);
    179349
    180     // save the new sources on the detection structure:
    181     detections->newSources = sources;
    182     return true;
    183 }
    184 
    185 // XXX fix the return status of this function...
    186 bool psphotSourceStatsUpdate (psArray *sources, pmConfig *config, pmReadout *readout) {
    187 
    188     bool status = false;
    189 
    190     psTimerStart ("psphot.stats");
    191 
    192     // select the appropriate recipe information
    193     psMetadata *recipe  = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE);
    194     assert (recipe);
    195 
    196     // determine the number of allowed threads
    197     int nThreads = psMetadataLookupS32(&status, config->arguments, "NTHREADS"); // Number of threads
    198     if (!status) {
    199         nThreads = 0;
    200     }
    201 
    202     // determine properties (sky, moments) of initial sources
    203     float OUTER    = psMetadataLookupF32 (&status, recipe, "SKY_OUTER_RADIUS");
    204     if (!status) return false;
    205 
    206     OUTER = PS_MAX(OUTER, 20.0); // XXX Guarantee that we can encompass the max moments radius
    207 
    208     char *breakPt  = psMetadataLookupStr (&status, recipe, "BREAK_POINT");
    209     if (!status) return false;
    210 
    211     for (int i = 0; i < sources->n; i++) {
    212 
    213         pmSource *source = sources->data[i];
    214         if (!source->peak) continue; // XXX how can we have a peak-less source?
    215 
    216         // allocate space for moments
    217         if (!source->moments) {
    218             source->moments = pmMomentsAlloc();
    219         }
    220 
    221         // allocate image, weight, mask arrays for each peak (square of radius OUTER)
    222         pmSourceDefinePixels (source, readout, source->peak->x, source->peak->y, OUTER);
    223         source->peak->assigned = true;
    224     }
    225 
    226     if (!strcasecmp (breakPt, "PEAKS")) {
    227         psLogMsg ("psphot", PS_LOG_INFO, "%ld sources : %f sec\n", sources->n, psTimerMark ("psphot.stats"));
    228         psLogMsg ("psphot", PS_LOG_INFO, "break point PEAKS, skipping MOMENTS\n");
    229         psphotVisualShowMoments (sources);
    230         return sources;
    231     }
    232 
    233     // XXX how else could we get the window size in?
    234     if (!psphotSetMomentsWindow(recipe, readout->analysis, sources)) {
    235         psError(PS_ERR_UNEXPECTED_NULL, false, "Failed to determine Moments Window!");
    236         return NULL;
    237     }
    238 
    239     // threaded measurement of the source magnitudes
    240     int Nfail = 0;
    241     int Nmoments = 0;
    242     int Nfaint = 0;
    243 
    244     // choose Cx, Cy (see psphotThreadTools.c for overview of the concepts)
    245     int Cx = 1, Cy = 1;
    246     psphotChooseCellSizes (&Cx, &Cy, readout, nThreads);
    247 
    248     psArray *cellGroups = psphotAssignSources (Cx, Cy, sources);
    249 
    250     for (int i = 0; i < cellGroups->n; i++) {
    251 
    252         psArray *cells = cellGroups->data[i];
    253 
    254         for (int j = 0; j < cells->n; j++) {
    255 
    256             // allocate a job -- if threads are not defined, this just runs the job
    257             psThreadJob *job = psThreadJobAlloc ("PSPHOT_SOURCE_STATS");
    258 
    259             psArrayAdd(job->args, 1, cells->data[j]); // sources
    260             psArrayAdd(job->args, 1, recipe);
    261             PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nmoments
    262             PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfail
    263             PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfaint
    264 
    265             if (!psThreadJobAddPending(job)) {
    266                 psError(PS_ERR_UNKNOWN, false, "Unable to launch thread job PSPHOT_SOURCE_STATS");
    267                 psFree (job);
    268                 return NULL;
    269             }
    270             psFree(job);
    271         }
    272 
    273         // wait for the threads to finish and manage results
    274         if (!psThreadPoolWait (false)) {
    275             psError(PS_ERR_UNKNOWN, false, "Failure in thread job PSPHOT_SOURCE_STATS");
    276             return NULL;
    277         }
    278 
    279         // we have only supplied one type of job, so we can assume the types here
    280         psThreadJob *job = NULL;
    281         while ((job = psThreadJobGetDone()) != NULL) {
    282             if (job->args->n < 1) {
    283                 fprintf (stderr, "error with job\n");
    284             } else {
    285                 psScalar *scalar = NULL;
    286                 scalar = job->args->data[2];
    287                 Nmoments += scalar->data.S32;
    288                 scalar = job->args->data[3];
    289                 Nfail += scalar->data.S32;
    290                 scalar = job->args->data[4];
    291                 Nfaint += scalar->data.S32;
    292             }
    293             psFree(job);
    294         }
    295     }
    296 
    297     psFree (cellGroups);
    298 
    299     psLogMsg ("psphot", PS_LOG_INFO, "%ld sources, %d moments, %d faint, %d failed: %f sec\n", sources->n, Nmoments, Nfaint, Nfail, psTimerMark ("psphot.stats"));
    300 
    301     psphotVisualShowMoments (sources);
    302 
    303350    return (sources);
    304351}
    … …  
    311358
    312359    psArray *sources                = job->args->data[0];
    313     psMetadata *recipe              = job->args->data[1];
    314 
    315     float INNER        = psMetadataLookupF32 (&status, recipe, "SKY_INNER_RADIUS");
    316     if (!status) return false;
    317     float MIN_SN       = psMetadataLookupF32 (&status, recipe, "MOMENTS_SN_MIN");
    318     if (!status) return false;
    319     float RADIUS       = psMetadataLookupF32 (&status, recipe, "PSF_MOMENTS_RADIUS");
    320     if (!status) return false;
    321     float SIGMA        = psMetadataLookupF32 (&status, recipe, "MOMENTS_GAUSS_SIGMA");
    322     if (!status) return false;
    323 
    324     // bit-masks to test for good/bad pixels
    325     psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT");
    326     assert (maskVal);
    327 
    328     // bit-mask to mark pixels not used in analysis
    329     psImageMaskType markVal = psMetadataLookupImageMask(&status, recipe, "MARK.PSPHOT");
    330     assert (markVal);
     360
     361    float INNER                     = PS_SCALAR_VALUE(job->args->data[1],F32);
     362    float MIN_SN                    = PS_SCALAR_VALUE(job->args->data[2],F32);
     363    float RADIUS                    = PS_SCALAR_VALUE(job->args->data[3],F32);
     364    float SIGMA                     = PS_SCALAR_VALUE(job->args->data[4],F32);
     365
     366    psImageMaskType maskVal         = PS_SCALAR_VALUE(job->args->data[5],PS_TYPE_IMAGE_MASK_DATA);
     367    psImageMaskType markVal         = PS_SCALAR_VALUE(job->args->data[6],PS_TYPE_IMAGE_MASK_DATA);
    331368
    332369    // maskVal is used to test for rejected pixels, and must include markVal
    … …  
    392429
    393430    // change the value of a scalar on the array (wrap this and put it in psArray.h)
    394     scalar = job->args->data[2];
     431    scalar = job->args->data[7];
    395432    scalar->data.S32 = Nmoments;
    396433
    397     scalar = job->args->data[3];
     434    scalar = job->args->data[8];
    398435    scalar->data.S32 = Nfail;
    399436
    400     scalar = job->args->data[4];
     437    scalar = job->args->data[9];
    401438    scalar->data.S32 = Nfaint;
    402439
    … …  
    410447    bool status;
    411448
    412     float MIN_SN = psMetadataLookupF32 (&status, recipe, "MOMENTS_SN_MIN");
    413     if (!status) return false;
     449    float MIN_SN = psMetadataLookupF32 (&status, recipe, "MOMENTS_SN_MIN"); psAssert (status, "missing MOMENTS_SN_MIN");
     450    float PSF_SN_LIM = psMetadataLookupF32(&status, recipe, "PSF_SN_LIM"); psAssert (status, "missing PSF_SN_LIM");
     451    psF32 MOMENTS_AR_MAX = psMetadataLookupF32(&status, recipe, "MOMENTS_AR_MAX"); psAssert (status, "missing MOMENTS_AR_MAX");
    414452
    415453    // XXX move this to a config file?
    … …  
    441479
    442480        // choose a grid scale that is a fixed fraction of the psf sigma^2
    443         psMetadataAddF32(recipe, PS_LIST_TAIL, "PSF_CLUMP_GRID_SCALE", PS_META_REPLACE, "clump grid", 0.1*PS_SQR(sigma[i]));
    444         psMetadataAddF32(recipe, PS_LIST_TAIL, "MOMENTS_SX_MAX", PS_META_REPLACE, "moments limit", 2.0*PS_SQR(sigma[i]));
    445         psMetadataAddF32(recipe, PS_LIST_TAIL, "MOMENTS_SY_MAX", PS_META_REPLACE, "moments limit", 2.0*PS_SQR(sigma[i]));
     481        float PSF_CLUMP_GRID_SCALE = 0.1*PS_SQR(sigma[i]);
     482        float MOMENTS_SX_MAX = 2.0*PS_SQR(sigma[i]);
     483        float MOMENTS_SY_MAX = 2.0*PS_SQR(sigma[i]);
    446484
    447485        // determine the PSF parameters from the source moment values
    448         pmPSFClump psfClump = pmSourcePSFClump (NULL, sources, recipe);
     486        pmPSFClump psfClump = pmSourcePSFClump (NULL, NULL, sources, PSF_SN_LIM, PSF_CLUMP_GRID_SCALE, MOMENTS_SX_MAX, MOMENTS_SY_MAX, MOMENTS_AR_MAX);
    449487        psLogMsg ("psphot", 3, "radius %.1f, nStars: %d, nSigma: %5.2f, X,  Y: %f, %f (%f, %f)\n", sigma[i], psfClump.nStars, psfClump.nSigma, psfClump.X, psfClump.Y, sqrt(psfClump.X) / sigma[i], sqrt(psfClump.Y) / sigma[i]);
    450488
    … …  
    485523
    486524    // choose a grid scale that is a fixed fraction of the psf sigma^2
    487     psMetadataAddF32(recipe, PS_LIST_TAIL, "PSF_CLUMP_GRID_SCALE", PS_META_REPLACE, "clump grid", 0.1*PS_SQR(Sigma));
    488     psMetadataAddF32(recipe, PS_LIST_TAIL, "MOMENTS_SX_MAX", PS_META_REPLACE, "moments limit", 2.0*PS_SQR(Sigma));
    489     psMetadataAddF32(recipe, PS_LIST_TAIL, "MOMENTS_SY_MAX", PS_META_REPLACE, "moments limit", 2.0*PS_SQR(Sigma));
    490     psMetadataAddF32(recipe, PS_LIST_TAIL, "MOMENTS_GAUSS_SIGMA", PS_META_REPLACE, "moments limit", Sigma);
    491     psMetadataAddF32(recipe, PS_LIST_TAIL, "PSF_MOMENTS_RADIUS", PS_META_REPLACE, "moments limit", 4.0*Sigma);
     525    psMetadataAddF32(analysis, PS_LIST_TAIL, "PSF_CLUMP_GRID_SCALE", PS_META_REPLACE, "clump grid", 0.1*PS_SQR(Sigma));
     526    psMetadataAddF32(analysis, PS_LIST_TAIL, "MOMENTS_SX_MAX", PS_META_REPLACE, "moments limit", 2.0*PS_SQR(Sigma));
     527    psMetadataAddF32(analysis, PS_LIST_TAIL, "MOMENTS_SY_MAX", PS_META_REPLACE, "moments limit", 2.0*PS_SQR(Sigma));
     528    psMetadataAddF32(analysis, PS_LIST_TAIL, "MOMENTS_GAUSS_SIGMA", PS_META_REPLACE, "moments limit", Sigma);
     529    psMetadataAddF32(analysis, PS_LIST_TAIL, "PSF_MOMENTS_RADIUS", PS_META_REPLACE, "moments limit", 4.0*Sigma);
    492530
    493531    psLogMsg ("psphot", 3, "using window function with sigma = %f\n", Sigma);
Note: See TracChangeset for help on using the changeset viewer.