Changeset 28029
- Timestamp:
- May 19, 2010, 2:08:46 PM (16 years ago)
- Location:
- branches/czw_branch/20100427/psphot
- Files:
-
- 52 edited
- 7 copied
-
. (modified) (1 prop)
-
doc/stack.txt (modified) (1 diff)
-
src/Makefile.am (modified) (3 diffs)
-
src/psphot.h (modified) (5 diffs)
-
src/psphotAddNoise.c (modified) (4 diffs)
-
src/psphotApResid.c (modified) (4 diffs)
-
src/psphotBasicDeblend.c (modified) (4 diffs)
-
src/psphotBlendFit.c (modified) (4 diffs)
-
src/psphotChoosePSF.c (modified) (4 diffs)
-
src/psphotCleanup.c (modified) (1 diff)
-
src/psphotDeblendSatstars.c (modified) (4 diffs)
-
src/psphotEfficiency.c (modified) (4 diffs)
-
src/psphotErrorCodes.dat (modified) (1 diff)
-
src/psphotExtendedSourceAnalysis.c (modified) (4 diffs)
-
src/psphotExtendedSourceAnalysisByObject.c (copied) (copied from trunk/psphot/src/psphotExtendedSourceAnalysisByObject.c )
-
src/psphotExtendedSourceFits.c (modified) (4 diffs)
-
src/psphotFindDetections.c (modified) (4 diffs)
-
src/psphotFitSourcesLinear.c (modified) (3 diffs)
-
src/psphotFitSourcesLinearStack.c (modified) (1 diff)
-
src/psphotForcedReadout.c (modified) (5 diffs)
-
src/psphotGuessModels.c (modified) (4 diffs)
-
src/psphotImageQuality.c (modified) (3 diffs)
-
src/psphotMagnitudes.c (modified) (3 diffs)
-
src/psphotMakePSFReadout.c (modified) (3 diffs)
-
src/psphotMaskReadout.c (modified) (3 diffs)
-
src/psphotMergeSources.c (modified) (12 diffs)
-
src/psphotModelBackground.c (modified) (2 diffs)
-
src/psphotModelTest.c (modified) (2 diffs)
-
src/psphotOutput.c (modified) (2 diffs)
-
src/psphotPetrosianRadialBins.c (modified) (1 diff)
-
src/psphotPetrosianStats.c (modified) (2 diffs)
-
src/psphotRadialApertures.c (copied) (copied from trunk/psphot/src/psphotRadialApertures.c )
-
src/psphotRadialAperturesByObject.c (copied) (copied from trunk/psphot/src/psphotRadialAperturesByObject.c )
-
src/psphotRadialBins.c (modified) (3 diffs)
-
src/psphotReadout.c (modified) (7 diffs)
-
src/psphotReadoutCleanup.c (modified) (3 diffs)
-
src/psphotReadoutFindPSF.c (modified) (6 diffs)
-
src/psphotReadoutKnownSources.c (modified) (5 diffs)
-
src/psphotReadoutMinimal.c (modified) (3 diffs)
-
src/psphotReplaceUnfit.c (modified) (2 diffs)
-
src/psphotRoughClass.c (modified) (4 diffs)
-
src/psphotSetMaskBits.c (modified) (1 diff)
-
src/psphotSkyReplace.c (modified) (3 diffs)
-
src/psphotSourceFreePixels.c (modified) (3 diffs)
-
src/psphotSourceMatch.c (modified) (7 diffs)
-
src/psphotSourceSize.c (modified) (4 diffs)
-
src/psphotSourceStats.c (modified) (4 diffs)
-
src/psphotStackArguments.c (modified) (4 diffs)
-
src/psphotStackChisqImage.c (modified) (7 diffs)
-
src/psphotStackImageLoop.c (modified) (7 diffs)
-
src/psphotStackMatchPSFs.c (modified) (2 diffs)
-
src/psphotStackMatchPSFsPrepare.c (copied) (copied from trunk/psphot/src/psphotStackMatchPSFsPrepare.c )
-
src/psphotStackMatchPSFsUtils.c (modified) (14 diffs)
-
src/psphotStackObjects.c (copied) (copied from trunk/psphot/src/psphotStackObjects.c )
-
src/psphotStackOptions.c (copied) (copied from trunk/psphot/src/psphotStackOptions.c )
-
src/psphotStackPSF.c (copied) (copied from trunk/psphot/src/psphotStackPSF.c )
-
src/psphotStackParseCamera.c (modified) (3 diffs)
-
src/psphotStackReadout.c (modified) (7 diffs)
-
src/psphotSubtractBackground.c (modified) (4 diffs)
Legend:
- Unmodified
- Added
- Removed
-
branches/czw_branch/20100427/psphot
-
Property svn:mergeinfo
set to
/branches/eam_branches/psphot.20100506 merged eligible /branches/pap/psphot merged eligible /trunk/psphot merged eligible
-
Property svn:mergeinfo
set to
-
branches/czw_branch/20100427/psphot/doc/stack.txt
r28017 r28029 1 2 20100506: 3 4 we may load RAW and/or CNV images. 5 6 * options->convolve: 7 * if only one of RAW or CNV exists, convolve it to match target PSF 8 * if both exist, convolve desired one 9 10 * if matching PSF exists, use it (unless told to re-measure) 11 12 * output goes to OUT (which is then used by psphot analysis routines) 13 14 * detect 15 16 * if (RAW) ? RAW : OUT 17 18 * 19 1 20 2 21 20100503: -
branches/czw_branch/20100427/psphot/src/Makefile.am
r28017 r28029 95 95 psphotFitSourcesLinearStack.c \ 96 96 psphotSourceMatch.c \ 97 psphotStackMatchPSFs.c \ 98 psphotStackMatchPSFsUtils.c \ 99 psphotStackMatchPSFsPrepare.c \ 100 psphotStackOptions.c \ 101 psphotStackObjects.c \ 102 psphotStackPSF.c \ 97 103 psphotCleanup.c 98 104 … … 159 165 psphotModelWithPSF.c \ 160 166 psphotExtendedSourceAnalysis.c \ 167 psphotExtendedSourceAnalysisByObject.c \ 161 168 psphotExtendedSourceFits.c \ 162 169 psphotKernelFromPSF.c \ … … 186 193 psphotEllipticalProfile.c \ 187 194 psphotRadialBins.c \ 195 psphotRadialApertures.c \ 196 psphotRadialAperturesByObject.c \ 188 197 psphotPetrosian.c \ 189 198 psphotPetrosianRadialBins.c \ -
branches/czw_branch/20100427/psphot/src/psphot.h
r28017 r28029 37 37 bool psphotReadoutMinimal(pmConfig *config, const pmFPAview *view); 38 38 39 bool psphotReadoutCleanup (pmConfig *config, const pmFPAview *view );40 bool psphotReadoutCleanupReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);39 bool psphotReadoutCleanup (pmConfig *config, const pmFPAview *view, const char *filerule); 40 bool psphotReadoutCleanupReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 41 41 42 42 bool psphotDefineFiles (pmConfig *config, pmFPAfile *input); … … 51 51 52 52 // psphotReadout functions 53 bool psphotAddPhotcode (pmConfig *config, const pmFPAview *view );53 bool psphotAddPhotcode (pmConfig *config, const pmFPAview *view, const char *filerule); 54 54 bool psphotAddPhotcodeReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index); 55 55 56 bool psphotSetMaskAndVariance (pmConfig *config, const pmFPAview *view );57 bool psphotSetMaskAndVarianceReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);58 59 bool psphotModelBackground (pmConfig *config, const pmFPAview *view );60 bool psphotModelBackgroundReadoutFileIndex (pmConfig *config, const pmFPAview *view, const char *file name, int index);61 62 bool psphotSubtractBackground (pmConfig *config, const pmFPAview *view );63 bool psphotSubtractBackgroundReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);64 65 bool psphotFindDetections (pmConfig *config, const pmFPAview *view, bool firstPass);66 bool psphotFindDetectionsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe, bool firstPass);67 68 bool psphotSourceStats (pmConfig *config, const pmFPAview *view, bool setWindow);69 bool psphotSourceStatsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe, bool setWindow);70 71 bool psphotDeblendSatstars (pmConfig *config, const pmFPAview *view );72 bool psphotDeblendSatstarsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index);73 74 bool psphotBasicDeblend (pmConfig *config, const pmFPAview *view );75 bool psphotBasicDeblendReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index);76 77 bool psphotRoughClass (pmConfig *config, const pmFPAview *view );78 bool psphotRoughClassReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);56 bool psphotSetMaskAndVariance (pmConfig *config, const pmFPAview *view, const char *filerule); 57 bool psphotSetMaskAndVarianceReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 58 59 bool psphotModelBackground (pmConfig *config, const pmFPAview *view, const char *filerule); 60 bool psphotModelBackgroundReadoutFileIndex (pmConfig *config, const pmFPAview *view, const char *filerule, int index); 61 62 bool psphotSubtractBackground (pmConfig *config, const pmFPAview *view, const char *filerule); 63 bool psphotSubtractBackgroundReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 64 65 bool psphotFindDetections (pmConfig *config, const pmFPAview *view, const char *filerule, bool firstPass); 66 bool psphotFindDetectionsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool firstPass); 67 68 bool psphotSourceStats (pmConfig *config, const pmFPAview *view, const char *filerule, bool setWindow); 69 bool psphotSourceStatsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool setWindow); 70 71 bool psphotDeblendSatstars (pmConfig *config, const pmFPAview *view, const char *filerule); 72 bool psphotDeblendSatstarsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index); 73 74 bool psphotBasicDeblend (pmConfig *config, const pmFPAview *view, const char *filerule); 75 bool psphotBasicDeblendReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index); 76 77 bool psphotRoughClass (pmConfig *config, const pmFPAview *view, const char *filerule); 78 bool psphotRoughClassReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 79 79 bool psphotRoughClassRegion (int nRegion, psRegion *region, psArray *sources, psMetadata *analysis, psMetadata *recipe, const bool havePSF); 80 80 81 bool psphotImageQuality (pmConfig *config, const pmFPAview *view );82 bool psphotImageQualityReadout(pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);83 84 bool psphotChoosePSF (pmConfig *config, const pmFPAview *view );85 bool psphotChoosePSFReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);86 87 bool psphotGuessModels (pmConfig *config, const pmFPAview *view );88 bool psphotGuessModelsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index);89 90 bool psphotMergeSources (pmConfig *config, const pmFPAview *view );91 bool psphotMergeSourcesReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index);92 93 bool psphotFitSourcesLinear (pmConfig *config, const pmFPAview *view, bool final);81 bool psphotImageQuality (pmConfig *config, const pmFPAview *view, const char *filerule); 82 bool psphotImageQualityReadout(pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 83 84 bool psphotChoosePSF (pmConfig *config, const pmFPAview *view, const char *filerule); 85 bool psphotChoosePSFReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 86 87 bool psphotGuessModels (pmConfig *config, const pmFPAview *view, const char *filerule); 88 bool psphotGuessModelsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index); 89 90 bool psphotMergeSources (pmConfig *config, const pmFPAview *view, const char *filerule); 91 bool psphotMergeSourcesReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index); 92 93 bool psphotFitSourcesLinear (pmConfig *config, const pmFPAview *view, const char *filerule, bool final); 94 94 bool psphotFitSourcesLinearReadout (psMetadata *recipe, pmReadout *readout, psArray *sources, pmPSF *psf, bool final); 95 95 96 bool psphotSourceSize (pmConfig *config, const pmFPAview *view, bool getPSFsize);97 bool psphotSourceSizeReadout(pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe, bool getPSFsize);98 99 bool psphotBlendFit (pmConfig *config, const pmFPAview *view );100 bool psphotBlendFitReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);96 bool psphotSourceSize (pmConfig *config, const pmFPAview *view, const char *filerule, bool getPSFsize); 97 bool psphotSourceSizeReadout(pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool getPSFsize); 98 99 bool psphotBlendFit (pmConfig *config, const pmFPAview *view, const char *filerule); 100 bool psphotBlendFitReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 101 101 bool psphotBlendFit_Threaded (psThreadJob *job); 102 102 103 bool psphotReplaceAllSources (pmConfig *config, const pmFPAview *view );104 bool psphotReplaceAllSourcesReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);105 106 bool psphotAddNoise (pmConfig *config, const pmFPAview *view );107 bool psphotSubNoise (pmConfig *config, const pmFPAview *view );108 bool psphotAddOrSubNoise (pmConfig *config, const pmFPAview *view, bool add);109 bool psphotAddOrSubNoiseReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe, bool add);110 111 bool psphotExtendedSourceAnalysis (pmConfig *config, const pmFPAview *view );112 bool psphotExtendedSourceAnalysisReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);113 114 bool psphotExtendedSourceFits (pmConfig *config, const pmFPAview *view );115 bool psphotExtendedSourceFitsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);116 117 bool psphotApResid (pmConfig *config, const pmFPAview *view );118 bool psphotApResidReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);119 120 bool psphotMagnitudes (pmConfig *config, const pmFPAview *view );103 bool psphotReplaceAllSources (pmConfig *config, const pmFPAview *view, const char *filerule); 104 bool psphotReplaceAllSourcesReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 105 106 bool psphotAddNoise (pmConfig *config, const pmFPAview *view, const char *filerule); 107 bool psphotSubNoise (pmConfig *config, const pmFPAview *view, const char *filerule); 108 bool psphotAddOrSubNoise (pmConfig *config, const pmFPAview *view, const char *filerule, bool add); 109 bool psphotAddOrSubNoiseReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool add); 110 111 bool psphotExtendedSourceAnalysis (pmConfig *config, const pmFPAview *view, const char *filerule); 112 bool psphotExtendedSourceAnalysisReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 113 114 bool psphotExtendedSourceFits (pmConfig *config, const pmFPAview *view, const char *filerule); 115 bool psphotExtendedSourceFitsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 116 117 bool psphotApResid (pmConfig *config, const pmFPAview *view, const char *filerule); 118 bool psphotApResidReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 119 120 bool psphotMagnitudes (pmConfig *config, const pmFPAview *view, const char *filerule); 121 121 bool psphotMagnitudesReadout(pmConfig *config, psMetadata *recipe, const pmFPAview *view, pmReadout *readout, psArray *sources, pmPSF *psf); 122 122 bool psphotMagnitudes_Threaded (psThreadJob *job); 123 123 124 bool psphotEfficiency (pmConfig *config, const pmFPAview *view );125 bool psphotEfficiencyReadout(pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe);124 bool psphotEfficiency (pmConfig *config, const pmFPAview *view, const char *filerule); 125 bool psphotEfficiencyReadout(pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 126 126 127 127 bool psphotPSFWeights(pmConfig *config, pmReadout *readout, const pmFPAview *view, psArray *sources); 128 128 bool psphotPSFWeights_Threaded (psThreadJob *job); 129 129 130 bool psphotSkyReplace (pmConfig *config, const pmFPAview *view );131 bool psphotSkyReplaceReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index);132 133 bool psphotSourceFreePixels (pmConfig *config, const pmFPAview *view );134 bool psphotSourceFreePixelsReadout(pmConfig *config, const pmFPAview *view, const char *file name, int index);130 bool psphotSkyReplace (pmConfig *config, const pmFPAview *view, const char *filerule); 131 bool psphotSkyReplaceReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index); 132 133 bool psphotSourceFreePixels (pmConfig *config, const pmFPAview *view, const char *filerule); 134 bool psphotSourceFreePixelsReadout(pmConfig *config, const pmFPAview *view, const char *filerule, int index); 135 135 136 136 // in psphotSourceStats.c: … … 147 147 148 148 // in psphotMergeSources.c: 149 bool psphotLoadExtSources (pmConfig *config, const pmFPAview *view );149 bool psphotLoadExtSources (pmConfig *config, const pmFPAview *view, const char *filerule); 150 150 psArray *psphotLoadPSFSources (pmConfig *config, const pmFPAview *view); 151 bool psphotRepairLoadedSources (pmConfig *config, const pmFPAview *view );152 bool psphotCheckExtSources (pmConfig *config, const pmFPAview *view );151 bool psphotRepairLoadedSources (pmConfig *config, const pmFPAview *view, const char *filerule); 152 bool psphotCheckExtSources (pmConfig *config, const pmFPAview *view, const char *filerule); 153 153 154 154 // generate the detection structure for the supplied array of sources 155 bool psphotDetectionsFromSources (pmConfig *config, const pmFPAview *view, psArray *sources);155 bool psphotDetectionsFromSources (pmConfig *config, const pmFPAview *view, const char *filerule, psArray *sources); 156 156 157 157 // generate the detection structure for the supplied array of sources … … 349 349 bool psphotStackImageLoop (pmConfig *config); 350 350 bool psphotStackReadout (pmConfig *config, const pmFPAview *view); 351 bool psphotStackChisqImage (pmConfig *config, const pmFPAview *view );351 bool psphotStackChisqImage (pmConfig *config, const pmFPAview *view, const char *ruleDet, const char *ruleCnv); 352 352 bool psphotStackChisqImageAddReadout(const pmConfig *config, // Configuration 353 353 const pmFPAview *view, 354 const char *filename, 354 355 pmReadout **chiReadout, 355 char *filename,356 356 int index); 357 357 358 bool psphotStackRemoveChisqFromInputs (pmConfig *config );358 bool psphotStackRemoveChisqFromInputs (pmConfig *config, const char *filerule); 359 359 bool pmFPAfileRemoveSingle(psMetadata *files, const char *name, int num); 360 360 361 psArray *psphotMatchSources (pmConfig *config, const pmFPAview *view );362 bool psphotMatchSourcesReadout (psArray *objects, pmConfig *config, const pmFPAview *view, c har *filename, int index);361 psArray *psphotMatchSources (pmConfig *config, const pmFPAview *view, const char *filerule); 362 bool psphotMatchSourcesReadout (psArray *objects, pmConfig *config, const pmFPAview *view, const char *filerule, int index); 363 363 bool psphotMatchSourcesToObjects (psArray *objects, psArray *sources, float RADIUS); 364 364 … … 367 367 int pmPhotObjSortByX (const void **a, const void **b); 368 368 369 typedef enum { 370 PSPHOT_CNV_SRC_NONE, 371 PSPHOT_CNV_SRC_AUTO, 372 PSPHOT_CNV_SRC_CNV, 373 PSPHOT_CNV_SRC_RAW, 374 } psphotStackConvolveSource; 375 376 /// Options for stacking process 377 typedef struct { 378 // Setup 379 380 int numCols; // size of image (X) 381 int numRows; // size of image (Y) 382 383 int num; // Number of inputs 384 bool convolve; // Convolve images? 385 psphotStackConvolveSource convolveSource; 386 // psArray *convImages, *convMasks, *convVariances; // Filenames for the temporary convolved images 387 388 // bool matchZPs; // Adjust relative fluxes based on transparency analysis? 389 // bool photometry; // Perform photometry? 390 // psMetadata *stats; // Statistics for output 391 // FILE *statsFile; // File to which to write statistics 392 // psArray *origImages, *origMasks, *origVariances; // Filenames of the original images 393 // psArray *origCovars; // Original covariances matrices 394 // int quality; // Bad data quality flag 395 396 // Prepare 397 pmPSF *psf; // Target PSF 398 psVector *inputSeeing; // Input seeing FWHMs 399 psVector *inputMask; // Mask for inputs 400 401 float targetSeeing; // Target seeing FWHM 402 psArray *sourceLists; // Individual lists of sources for matching 403 psVector *norm; // Normalisation for each image 404 psArray *psfs; 405 406 // psVector *exposures; // Exposure times 407 // float sumExposure; // Sum of exposure times 408 // float zp; // Zero point for output 409 // psVector *inputMask; // Mask for inputs 410 // psArray *sources; // Matched sources 411 412 // Convolve 413 psArray *kernels; // PSF-matching kernels --- required in the stacking 414 psArray *regions; // PSF-matching regions --- required in the stacking 415 psVector *matchChi2; // chi^2 for stamps from matching 416 psVector *weightings; // Combination weightings for images (1/noise^2) 417 // psArray *cells; // Cells for convolved images --- a handle for reading again 418 // int numCols, numRows; // Size of image 419 // psArray *convCovars; // Convolved covariance matrices 420 421 // Combine initial 422 // pmReadout *outRO; // Output readout 423 // pmReadout *expRO; // Exposure readout 424 // psArray *inspect; // Array of arrays of pixels to inspect 425 426 // Rejection 427 // psArray *rejected; // Rejected pixels 428 } psphotStackOptions; 429 430 /*** psphotStackMatchPSF prototypes ***/ 431 bool psphotStackMatchPSFs (pmConfig *config, const pmFPAview *view); 432 bool psphotStackMatchPSFsReadout (pmConfig *config, const pmFPAview *view, psphotStackOptions *options, int index); 433 bool psphotStackMatchPSFsPrepare (pmConfig *config, const pmFPAview *view, psphotStackOptions *options, int index); 434 435 // psphotStackMatchPSFsUtils 436 psVector *SetOptWidths (bool *optimum, psMetadata *recipe); 437 pmReadout *makeFakeReadout(pmConfig *config, pmReadout *raw, psArray *sources, pmPSF *psf, psImageMaskType maskVal, int fullSize); 438 bool rescaleData(pmReadout *readout, pmConfig *config, psphotStackOptions *options, int index); 439 bool renormKernel(pmReadout *readout, psphotStackOptions *options, int index); 440 bool saveMatchData (pmReadout *readout, psphotStackOptions *options, int index); 441 bool matchKernel(pmConfig *config, pmReadout *cnv, pmReadout *raw, psphotStackOptions *options, int index); 442 bool dumpImageDiff(pmReadout *readoutConv, pmReadout *readoutFake, pmReadout *readoutRef, int index, char *rootname); 443 bool dumpImage(pmReadout *readoutOut, pmReadout *readoutRef, int index, char *rootname); 444 bool loadKernel (pmConfig *config, pmReadout *readoutCnv, psphotStackOptions *options, int index); 445 446 pmPSF *psphotStackPSF(const pmConfig *config, int numCols, int numRows, const psArray *psfs, const psVector *inputMask); 447 448 psphotStackOptions *psphotStackOptionsAlloc (int num); 449 psphotStackConvolveSource psphotStackConvolveSourceFromString (const char *string); 450 pmFPAfile *psphotStackGetConvolveSource (pmConfig *config, psphotStackOptions *options, int index); 451 452 bool psphotCopySources (pmConfig *config, const pmFPAview *view, const char *ruleOut, const char *ruleSrc); 453 bool psphotCopySourcesReadout (pmConfig *config, const pmFPAview *view, const char *ruleOut, const char *ruleSrc, int index); 454 455 bool psphotRadialApertures (pmConfig *config, const pmFPAview *view, const char *filerule); 456 bool psphotRadialAperturesReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); 457 bool psphotRadialApertureSource (pmSource *source, psMetadata *recipe, float skynoise, psImageMaskType maskVal, const psVector *radMax); 458 459 bool psphotExtendedSourceAnalysisByObject (pmConfig *config, psArray *objects, const pmFPAview *view, const char *filerule); 460 bool psphotRadialAperturesByObject (pmConfig *config, psArray *objects, const pmFPAview *view, const char *filerule); 461 462 bool psphotStackObjectsUnifyPosition (psArray *objects); 463 369 464 #endif -
branches/czw_branch/20100427/psphot/src/psphotAddNoise.c
r26894 r28029 1 1 # include "psphotInternal.h" 2 2 3 bool psphotAddNoise (pmConfig *config, const pmFPAview *view ) {4 return psphotAddOrSubNoise (config, view, true);3 bool psphotAddNoise (pmConfig *config, const pmFPAview *view, const char *filerule) { 4 return psphotAddOrSubNoise (config, view, filerule, true); 5 5 } 6 6 7 bool psphotSubNoise (pmConfig *config, const pmFPAview *view ) {8 return psphotAddOrSubNoise (config, view, f alse);7 bool psphotSubNoise (pmConfig *config, const pmFPAview *view, const char *filerule) { 8 return psphotAddOrSubNoise (config, view, filerule, false); 9 9 } 10 10 11 11 // for now, let's store the detections on the readout->analysis for each readout 12 bool psphotAddOrSubNoise (pmConfig *config, const pmFPAview *view, bool add)12 bool psphotAddOrSubNoise (pmConfig *config, const pmFPAview *view, const char *filerule, bool add) 13 13 { 14 14 bool status = true; … … 23 23 // loop over the available readouts 24 24 for (int i = 0; i < num; i++) { 25 if (!psphotAddOrSubNoiseReadout (config, view, "PSPHOT.INPUT", i, recipe, add)) {26 psError (PSPHOT_ERR_CONFIG, false, "failed on to modify noise for PSPHOT.INPUT entry %d", i);25 if (!psphotAddOrSubNoiseReadout (config, view, filerule, i, recipe, add)) { 26 psError (PSPHOT_ERR_CONFIG, false, "failed on to modify noise for %s entry %d", filerule, i); 27 27 return false; 28 28 } … … 31 31 } 32 32 33 bool psphotAddOrSubNoiseReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe, bool add) {33 bool psphotAddOrSubNoiseReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool add) { 34 34 35 35 bool status = false; … … 39 39 40 40 // find the currently selected readout 41 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest41 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 42 42 psAssert (file, "missing file?"); 43 43 -
branches/czw_branch/20100427/psphot/src/psphotApResid.c
r27657 r28029 5 5 6 6 // for now, let's store the detections on the readout->analysis for each readout 7 bool psphotApResid (pmConfig *config, const pmFPAview *view )7 bool psphotApResid (pmConfig *config, const pmFPAview *view, const char *filerule) 8 8 { 9 9 bool status = true; … … 23 23 for (int i = 0; i < num; i++) { 24 24 if (i == chisqNum) continue; // skip chisq image 25 if (!psphotApResidReadout (config, view, "PSPHOT.INPUT", i, recipe)) {26 psError (PSPHOT_ERR_CONFIG, false, "failed to measure aperture residual for PSPHOT.INPUT entry %d", i);25 if (!psphotApResidReadout (config, view, filerule, i, recipe)) { 26 psError (PSPHOT_ERR_CONFIG, false, "failed to measure aperture residual for %s entry %d", filerule, i); 27 27 return false; 28 28 } … … 31 31 } 32 32 33 bool psphotApResidReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe)33 bool psphotApResidReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) 34 34 { 35 35 int Nfail = 0; … … 43 43 44 44 // find the currently selected readout 45 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest45 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 46 46 psAssert (file, "missing file?"); 47 47 -
branches/czw_branch/20100427/psphot/src/psphotBasicDeblend.c
r26894 r28029 2 2 3 3 // for now, let's store the detections on the readout->analysis for each readout 4 bool psphotBasicDeblend (pmConfig *config, const pmFPAview *view )4 bool psphotBasicDeblend (pmConfig *config, const pmFPAview *view, const char *filerule) 5 5 { 6 6 bool status = true; … … 11 11 // loop over the available readouts 12 12 for (int i = 0; i < num; i++) { 13 if (!psphotBasicDeblendReadout (config, view, "PSPHOT.INPUT", i)) {14 psError (PSPHOT_ERR_CONFIG, false, "failed on basic deblend analysis for PSPHOT.INPUT entry %d", i);13 if (!psphotBasicDeblendReadout (config, view, filerule, i)) { 14 psError (PSPHOT_ERR_CONFIG, false, "failed on basic deblend analysis for %s entry %d", filerule, i); 15 15 return false; 16 16 } … … 19 19 } 20 20 21 bool psphotBasicDeblendReadout (pmConfig *config, const pmFPAview *view, const char *file name, int fileIndex) {21 bool psphotBasicDeblendReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int fileIndex) { 22 22 23 23 int N; … … 31 31 32 32 // find the currently selected readout 33 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, fileIndex); // File of interest33 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, fileIndex); // File of interest 34 34 psAssert (file, "missing file?"); 35 35 -
branches/czw_branch/20100427/psphot/src/psphotBlendFit.c
r26894 r28029 2 2 3 3 // for now, let's store the detections on the readout->analysis for each readout 4 bool psphotBlendFit (pmConfig *config, const pmFPAview *view )4 bool psphotBlendFit (pmConfig *config, const pmFPAview *view, const char *filerule) 5 5 { 6 6 bool status = true; … … 15 15 // loop over the available readouts 16 16 for (int i = 0; i < num; i++) { 17 if (!psphotBlendFitReadout (config, view, "PSPHOT.INPUT", i, recipe)) {18 psError (PSPHOT_ERR_CONFIG, false, "failed to fit sources (non-linear) for PSPHOT.INPUT entry %d", i);17 if (!psphotBlendFitReadout (config, view, filerule, i, recipe)) { 18 psError (PSPHOT_ERR_CONFIG, false, "failed to fit sources (non-linear) for %s entry %d", filerule, i); 19 19 return false; 20 20 } … … 24 24 25 25 // XXX I don't like this name 26 bool psphotBlendFitReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe) {26 bool psphotBlendFitReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { 27 27 28 28 int Nfit = 0; … … 35 35 36 36 // find the currently selected readout 37 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest37 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 38 38 psAssert (file, "missing file?"); 39 39 -
branches/czw_branch/20100427/psphot/src/psphotChoosePSF.c
r27657 r28029 2 2 3 3 // generate a PSF model for inputs without PSF models already loaded 4 bool psphotChoosePSF (pmConfig *config, const pmFPAview *view )4 bool psphotChoosePSF (pmConfig *config, const pmFPAview *view, const char *filerule) 5 5 { 6 6 bool status = true; … … 20 20 for (int i = 0; i < num; i++) { 21 21 if (i == chisqNum) continue; // skip chisq image 22 if (!psphotChoosePSFReadout (config, view, "PSPHOT.INPUT", i, recipe)) {23 psError (PSPHOT_ERR_CONFIG, false, "failed to choose a psf model for PSPHOT.INPUT entry %d", i);22 if (!psphotChoosePSFReadout (config, view, filerule, i, recipe)) { 23 psError (PSPHOT_ERR_CONFIG, false, "failed to choose a psf model for %s entry %d", filerule, i); 24 24 return false; 25 25 } … … 29 29 30 30 // try PSF models and select best option 31 bool psphotChoosePSFReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe) {31 bool psphotChoosePSFReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { 32 32 33 33 bool status; … … 36 36 37 37 // find the currently selected readout 38 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest38 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 39 39 psAssert (file, "missing file?"); 40 40 -
branches/czw_branch/20100427/psphot/src/psphotCleanup.c
r28017 r28029 19 19 pmConceptsDone (); 20 20 pmConfigDone (); 21 psLibFinalize(); 21 22 // fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, NULL, false), "psphot"); 22 23 fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, stdout, false), "psphot"); -
branches/czw_branch/20100427/psphot/src/psphotDeblendSatstars.c
r26894 r28029 2 2 3 3 // for now, let's store the detections on the readout->analysis for each readout 4 bool psphotDeblendSatstars (pmConfig *config, const pmFPAview *view )4 bool psphotDeblendSatstars (pmConfig *config, const pmFPAview *view, const char *filerule) 5 5 { 6 6 bool status = true; … … 11 11 // loop over the available readouts 12 12 for (int i = 0; i < num; i++) { 13 if (!psphotDeblendSatstarsReadout (config, view, "PSPHOT.INPUT", i)) {14 psError (PSPHOT_ERR_CONFIG, false, "failed on saturated star deblend analysis for PSPHOT.INPUT entry %d", i);13 if (!psphotDeblendSatstarsReadout (config, view, filerule, i)) { 14 psError (PSPHOT_ERR_CONFIG, false, "failed on saturated star deblend analysis for %s entry %d", filerule, i); 15 15 return false; 16 16 } … … 19 19 } 20 20 21 bool psphotDeblendSatstarsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int fileIndex) {21 bool psphotDeblendSatstarsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int fileIndex) { 22 22 23 23 int N; … … 31 31 32 32 // find the currently selected readout 33 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, fileIndex); // File of interest33 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, fileIndex); // File of interest 34 34 psAssert (file, "missing file?"); 35 35 -
branches/czw_branch/20100427/psphot/src/psphotEfficiency.c
r28017 r28029 156 156 } 157 157 158 bool psphotEfficiency (pmConfig *config, const pmFPAview *view )158 bool psphotEfficiency (pmConfig *config, const pmFPAview *view, const char *filerule) 159 159 { 160 160 bool status = true; … … 173 173 // loop over the available readouts 174 174 for (int i = 0; i < num; i++) { 175 if (i == chisqNum) continue; // skip chisq image176 if (!psphotEfficiencyReadout (config, view, "PSPHOT.INPUT", i, recipe)) {177 psError (PSPHOT_ERR_CONFIG, false, "failed to measure detection efficiency for PSPHOT.INPUT entry %d", i);175 if (i == chisqNum) continue; // skip chisq image 176 if (!psphotEfficiencyReadout (config, view, filerule, i, recipe)) { 177 psError (PSPHOT_ERR_CONFIG, false, "failed to measure detection efficiency for %s entry %d", filerule, i); 178 178 return false; 179 179 } … … 183 183 184 184 // Determine detection efficiency 185 bool psphotEfficiencyReadout(pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe)185 bool psphotEfficiencyReadout(pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) 186 186 { 187 187 bool status = true; … … 190 190 191 191 // find the currently selected readout 192 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest192 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 193 193 psAssert (file, "missing file?"); 194 194 -
branches/czw_branch/20100427/psphot/src/psphotErrorCodes.dat
r27002 r28029 11 11 APERTURE Problem with aperture photometry 12 12 SKY Problem in sky determination 13 IO Problem in data I/O 13 14 # these errors correspond to standard exit conditions 14 15 ARGUMENTS Incorrect arguments -
branches/czw_branch/20100427/psphot/src/psphotExtendedSourceAnalysis.c
r28017 r28029 1 1 # include "psphotInternal.h" 2 2 3 // ?? these cannot happen here --> we would need to do this in psphotExtendedSourceAnalysis 4 // XXX option to choose a consistent position 5 // XXX option to choose a consistent elliptical contour 6 // XXX SDSS uses the r-band petrosian radius to measure petrosian fluxes in all bands 7 // XXX consistent choice of extendedness... 8 3 9 // for now, let's store the detections on the readout->analysis for each readout 4 bool psphotExtendedSourceAnalysis (pmConfig *config, const pmFPAview *view )10 bool psphotExtendedSourceAnalysis (pmConfig *config, const pmFPAview *view, const char *filerule) 5 11 { 6 12 bool status = true; … … 21 27 // loop over the available readouts 22 28 for (int i = 0; i < num; i++) { 23 if (!psphotExtendedSourceAnalysisReadout (config, view, "PSPHOT.INPUT", i, recipe)) {24 psError (PSPHOT_ERR_CONFIG, false, "failed on measure extended source aperture-like parameters for PSPHOT.INPUT entry %d", i);29 if (!psphotExtendedSourceAnalysisReadout (config, view, filerule, i, recipe)) { 30 psError (PSPHOT_ERR_CONFIG, false, "failed on measure extended source aperture-like parameters for %s entry %d", filerule, i); 25 31 return false; 26 32 } … … 30 36 31 37 // aperture-like measurements for extended sources 32 bool psphotExtendedSourceAnalysisReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe) {38 bool psphotExtendedSourceAnalysisReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { 33 39 34 40 bool status; … … 42 48 43 49 // find the currently selected readout 44 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest50 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 45 51 psAssert (file, "missing file?"); 46 52 -
branches/czw_branch/20100427/psphot/src/psphotExtendedSourceFits.c
r26894 r28029 2 2 3 3 // for now, let's store the detections on the readout->analysis for each readout 4 bool psphotExtendedSourceFits (pmConfig *config, const pmFPAview *view )4 bool psphotExtendedSourceFits (pmConfig *config, const pmFPAview *view, const char *filerule) 5 5 { 6 6 bool status = true; … … 21 21 // loop over the available readouts 22 22 for (int i = 0; i < num; i++) { 23 if (!psphotExtendedSourceFitsReadout (config, view, "PSPHOT.INPUT", i, recipe)) {24 psError (PSPHOT_ERR_CONFIG, false, "failed on to fit extended sources for PSPHOT.INPUT entry %d", i);23 if (!psphotExtendedSourceFitsReadout (config, view, filerule, i, recipe)) { 24 psError (PSPHOT_ERR_CONFIG, false, "failed on to fit extended sources for %s entry %d", filerule, i); 25 25 return false; 26 26 } … … 30 30 31 31 // non-linear model fitting for extended sources 32 bool psphotExtendedSourceFitsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe) {32 bool psphotExtendedSourceFitsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { 33 33 34 34 bool status; … … 41 41 42 42 // find the currently selected readout 43 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest43 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 44 44 psAssert (file, "missing file?"); 45 45 -
branches/czw_branch/20100427/psphot/src/psphotFindDetections.c
r26894 r28029 4 4 // peaks and new footprints. any old peaks are saved on oldPeaks. the resulting footprint set 5 5 // contains all footprints (old and new) 6 bool psphotFindDetections (pmConfig *config, const pmFPAview *view, bool firstPass)6 bool psphotFindDetections (pmConfig *config, const pmFPAview *view, const char *filerule, bool firstPass) 7 7 { 8 8 bool status = true; … … 17 17 // loop over the available readouts 18 18 for (int i = 0; i < num; i++) { 19 if (!psphotFindDetectionsReadout (config, view, "PSPHOT.INPUT", i, recipe, firstPass)) {20 psError (PSPHOT_ERR_CONFIG, false, "failed to find initial detections for PSPHOT.INPUT entry %d", i);19 if (!psphotFindDetectionsReadout (config, view, filerule, i, recipe, firstPass)) { 20 psError (PSPHOT_ERR_CONFIG, false, "failed to find initial detections for %s entry %d", filerule, i); 21 21 return false; 22 22 } … … 26 26 27 27 // smooth the image, search for peaks, optionally define footprints based on the peaks 28 bool psphotFindDetectionsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe, bool firstPass) {28 bool psphotFindDetectionsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool firstPass) { 29 29 30 30 bool status; … … 34 34 35 35 // find the currently selected readout 36 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest36 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 37 37 psAssert (file, "missing file?"); 38 38 -
branches/czw_branch/20100427/psphot/src/psphotFitSourcesLinear.c
r27532 r28029 13 13 14 14 // for now, let's store the detections on the readout->analysis for each readout 15 bool psphotFitSourcesLinear (pmConfig *config, const pmFPAview *view, bool final)15 bool psphotFitSourcesLinear (pmConfig *config, const pmFPAview *view, const char *filerule, bool final) 16 16 { 17 17 bool status = true; … … 28 28 29 29 // find the currently selected readout 30 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", i); // File of interest30 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, i); // File of interest 31 31 psAssert (file, "missing file?"); 32 32 … … 44 44 45 45 if (!psphotFitSourcesLinearReadout (recipe, readout, sources, psf, final)) { 46 psError (PSPHOT_ERR_CONFIG, false, "failed to fit sources (linear) for PSPHOT.INPUT entry %d", i);46 psError (PSPHOT_ERR_CONFIG, false, "failed to fit sources (linear) for %s entry %d", filerule, i); 47 47 return false; 48 48 } -
branches/czw_branch/20100427/psphot/src/psphotFitSourcesLinearStack.c
r27657 r28029 28 28 29 29 // analysis is done in spatial order (to speed up overlap search) 30 // sort by first element in each source list31 30 objects = psArraySort (objects, pmPhotObjSortByX); 32 31 -
branches/czw_branch/20100427/psphot/src/psphotForcedReadout.c
r27314 r28029 20 20 21 21 // set the photcode for this image 22 if (!psphotAddPhotcode (config, view )) {22 if (!psphotAddPhotcode (config, view, "PSPHOT.INPUT")) { 23 23 psError (PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); 24 24 return false; … … 30 30 31 31 // Generate the mask and weight images, including the user-defined analysis region of interest 32 psphotSetMaskAndVariance (config, view );32 psphotSetMaskAndVariance (config, view, "PSPHOT.INPUT"); 33 33 if (!strcasecmp (breakPt, "NOTHING")) { 34 return psphotReadoutCleanup (config, view);34 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 35 35 } 36 36 37 37 // generate a background model (median, smoothed image) 38 if (!psphotModelBackground (config, view )) {39 return psphotReadoutCleanup (config, view );38 if (!psphotModelBackground (config, view, "PSPHOT.INPUT")) { 39 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 40 40 } 41 if (!psphotSubtractBackground (config, view )) {42 return psphotReadoutCleanup (config, view );41 if (!psphotSubtractBackground (config, view, "PSPHOT.INPUT")) { 42 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 43 43 } 44 44 if (!strcasecmp (breakPt, "BACKMDL")) { 45 return psphotReadoutCleanup (config, view );45 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 46 46 } 47 47 … … 49 49 // this only happens if we had a programming error in psphotLoadPSF 50 50 psError (PSPHOT_ERR_UNKNOWN, false, "error loading psf model"); 51 return psphotReadoutCleanup (config, view );51 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 52 52 } 53 53 54 54 // include externally-supplied sources 55 psphotLoadExtSources (config, view );55 psphotLoadExtSources (config, view, "PSPHOT.INPUT"); 56 56 57 57 // construct an initial model for each object, set the radius to fitRadius, set circular fit mask 58 psphotGuessModels (config, view );58 psphotGuessModels (config, view, "PSPHOT.INPUT"); 59 59 60 60 // merge the newly selected sources into the existing list 61 61 // NOTE: merge OLD and NEW 62 psphotMergeSources (config, view );62 psphotMergeSources (config, view, "PSPHOT.INPUT"); 63 63 64 64 // linear PSF fit to source peaks, subtract the models from the image (in PSF mask) 65 psphotFitSourcesLinear (config, view, false);65 psphotFitSourcesLinear (config, view, "PSPHOT.INPUT", false); 66 66 67 67 // identify CRs and extended sources … … 71 71 72 72 // calculate source magnitudes 73 psphotMagnitudes(config, view );73 psphotMagnitudes(config, view, "PSPHOT.INPUT"); 74 74 75 75 // XXX do I want to do this? … … 80 80 81 81 // replace background in residual image 82 psphotSkyReplace (config, view );82 psphotSkyReplace (config, view, "PSPHOT.INPUT"); 83 83 84 84 // drop the references to the image pixels held by each source 85 psphotSourceFreePixels (config, view );85 psphotSourceFreePixels (config, view, "PSPHOT.INPUT"); 86 86 87 87 // create the exported-metadata and free local data 88 return psphotReadoutCleanup (config, view);88 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 89 89 } -
branches/czw_branch/20100427/psphot/src/psphotGuessModels.c
r27657 r28029 8 8 9 9 // for now, let's store the detections on the readout->analysis for each readout 10 bool psphotGuessModels (pmConfig *config, const pmFPAview *view )10 bool psphotGuessModels (pmConfig *config, const pmFPAview *view, const char *filerule) 11 11 { 12 12 bool status = true; … … 22 22 for (int i = 0; i < num; i++) { 23 23 if (i == chisqNum) continue; // skip chisq image 24 if (!psphotGuessModelsReadout (config, view, "PSPHOT.INPUT", i)) {25 psError (PSPHOT_ERR_CONFIG, false, "failed on to guess models for PSPHOT.INPUT entry %d", i);24 if (!psphotGuessModelsReadout (config, view, filerule, i)) { 25 psError (PSPHOT_ERR_CONFIG, false, "failed on to guess models for %s entry %d", filerule, i); 26 26 return false; 27 27 } … … 31 31 32 32 // construct an initial PSF model for each object (new sources only) 33 bool psphotGuessModelsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index) {33 bool psphotGuessModelsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index) { 34 34 35 35 bool status; … … 38 38 39 39 // find the currently selected readout 40 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest40 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 41 41 psAssert (file, "missing file?"); 42 42 -
branches/czw_branch/20100427/psphot/src/psphotImageQuality.c
r27657 r28029 5 5 6 6 // for now, let's store the detections on the readout->analysis for each readout 7 bool psphotImageQuality (pmConfig *config, const pmFPAview *view )7 bool psphotImageQuality (pmConfig *config, const pmFPAview *view, const char *filerule) 8 8 { 9 9 bool status = true; … … 23 23 for (int i = 0; i < num; i++) { 24 24 if (i == chisqNum) continue; // skip chisq image 25 if (!psphotImageQualityReadout (config, view, "PSPHOT.INPUT", i, recipe)) {26 psError (PSPHOT_ERR_CONFIG, false, "failed on to measure image quality for PSPHOT.INPUT entry %d", i);25 if (!psphotImageQualityReadout (config, view, filerule, i, recipe)) { 26 psError (PSPHOT_ERR_CONFIG, false, "failed on to measure image quality for %s entry %d", filerule, i); 27 27 return false; 28 28 } … … 32 32 33 33 // selecting the 'good' stars (likely to be psf stars), measure the M_cn, M_sn terms for n = 2,3,4 34 bool psphotImageQualityReadout(pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe)34 bool psphotImageQualityReadout(pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) 35 35 { 36 36 bool status = true; 37 37 38 38 // find the currently selected readout 39 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest39 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 40 40 psAssert (file, "missing file?"); 41 41 -
branches/czw_branch/20100427/psphot/src/psphotMagnitudes.c
r27657 r28029 1 1 # include "psphotInternal.h" 2 2 3 bool psphotMagnitudes (pmConfig *config, const pmFPAview *view )3 bool psphotMagnitudes (pmConfig *config, const pmFPAview *view, const char *filerule) 4 4 { 5 5 bool status = true; … … 21 21 22 22 // find the currently selected readout 23 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", i); // File of interest23 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, i); // File of interest 24 24 psAssert (file, "missing file?"); 25 25 … … 37 37 38 38 if (!psphotMagnitudesReadout (config, recipe, view, readout, sources, psf)) { 39 psError (PSPHOT_ERR_CONFIG, false, "failed to measure magnitudes for PSPHOT.INPUT entry %d", i);39 psError (PSPHOT_ERR_CONFIG, false, "failed to measure magnitudes for %s entry %d", filerule, i); 40 40 return false; 41 41 } -
branches/czw_branch/20100427/psphot/src/psphotMakePSFReadout.c
r26894 r28029 19 19 20 20 // set the photcode for this image 21 if (!psphotAddPhotcode (config, view )) {21 if (!psphotAddPhotcode (config, view, "PSPHOT.INPUT")) { 22 22 psError (PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); 23 23 return false; … … 29 29 30 30 // Generate the mask and weight images, including the user-defined analysis region of interest 31 psphotSetMaskAndVariance (config, view );31 psphotSetMaskAndVariance (config, view, "PSPHOT.INPUT"); 32 32 if (!strcasecmp (breakPt, "NOTHING")) { 33 return psphotReadoutCleanup (config, view);33 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 34 34 } 35 35 36 36 // generate a background model (median, smoothed image) 37 if (!psphotModelBackground (config, view )) {38 return psphotReadoutCleanup (config, view );37 if (!psphotModelBackground (config, view, "PSPHOT.INPUT")) { 38 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 39 39 } 40 if (!psphotSubtractBackground (config, view )) {41 return psphotReadoutCleanup (config, view );40 if (!psphotSubtractBackground (config, view, "PSPHOT.INPUT")) { 41 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 42 42 } 43 43 if (!strcasecmp (breakPt, "BACKMDL")) { 44 return psphotReadoutCleanup (config, view );44 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 45 45 } 46 46 47 psphotLoadExtSources (config, view );47 psphotLoadExtSources (config, view, "PSPHOT.INPUT"); 48 48 49 49 // If sources have been supplied, then these should be used to measure the PSF include … … 53 53 // a text file have no valid SN values, but psphotChoosePSF needs to select the top 54 54 // PSF_MAX_NSTARS to generate the PSF. 55 if (!psphotCheckExtSources (config, view )) {55 if (!psphotCheckExtSources (config, view, "PSPHOT.INPUT")) { 56 56 psLogMsg ("psphot", 3, "failure to select possible PSF sources (external or internal)"); 57 return psphotReadoutCleanup (config, view );57 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 58 58 } 59 59 60 60 // Use bright stellar objects to measure PSF. If we do not have enough stars to generate 61 61 // the PSF, build one from the SEEING guess and model class 62 if (!psphotChoosePSF (config, view )) {62 if (!psphotChoosePSF (config, view, "PSPHOT.INPUT")) { 63 63 psLogMsg ("psphot", 3, "failure to construct a psf model"); 64 return psphotReadoutCleanup (config, view );64 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 65 65 } 66 66 67 67 // measure aperture photometry corrections 68 68 # if 0 69 if (!psphotApResid (config, view )) {69 if (!psphotApResid (config, view, "PSPHOT.INPUT")) { 70 70 psLogMsg ("psphot", 3, "failed on psphotApResid"); 71 return psphotReadoutCleanup (config, view );71 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 72 72 } 73 73 # endif 74 74 75 return psphotReadoutCleanup (config, view );75 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 76 76 } -
branches/czw_branch/20100427/psphot/src/psphotMaskReadout.c
r26894 r28029 1 1 # include "psphotInternal.h" 2 2 3 bool psphotSetMaskAndVariance (pmConfig *config, const pmFPAview *view ) {3 bool psphotSetMaskAndVariance (pmConfig *config, const pmFPAview *view, const char *filerule) { 4 4 5 5 bool status = false; … … 16 16 17 17 // Generate the mask and weight images, including the user-defined analysis region of interest 18 if (!psphotSetMaskAndVarianceReadout (config, view, "PSPHOT.INPUT", i, recipe)) {18 if (!psphotSetMaskAndVarianceReadout (config, view, filerule, i, recipe)) { 19 19 psError (PSPHOT_ERR_CONFIG, false, "failed to generate mask for PSPHOT.INPUT entry %d", i); 20 20 return false; … … 25 25 26 26 // generate mask and variance if not defined, additional mask for restricted subregion 27 bool psphotSetMaskAndVarianceReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe) {27 bool psphotSetMaskAndVarianceReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { 28 28 29 29 bool status; 30 30 31 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", index); // File of interest31 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 32 32 psAssert (file, "missing file?"); 33 33 -
branches/czw_branch/20100427/psphot/src/psphotMergeSources.c
r27657 r28029 6 6 7 7 // for now, let's store the detections on the readout->analysis for each readout 8 bool psphotMergeSources (pmConfig *config, const pmFPAview *view )8 bool psphotMergeSources (pmConfig *config, const pmFPAview *view, const char *filerule) 9 9 { 10 10 bool status = true; … … 15 15 // loop over the available readouts 16 16 for (int i = 0; i < num; i++) { 17 if (!psphotMergeSourcesReadout (config, view, "PSPHOT.INPUT", i)) {18 psError (PSPHOT_ERR_CONFIG, false, "failed to merge sources for PSPHOT.INPUT entry %d", i);17 if (!psphotMergeSourcesReadout (config, view, filerule, i)) { 18 psError (PSPHOT_ERR_CONFIG, false, "failed to merge sources for %s entry %d", filerule, i); 19 19 return false; 20 20 } … … 24 24 25 25 // add newly selected sources to the existing list of sources 26 bool psphotMergeSourcesReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index) {26 bool psphotMergeSourcesReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index) { 27 27 28 28 bool status; 29 29 30 30 // find the currently selected readout 31 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest31 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 32 32 psAssert (file, "missing file?"); 33 33 … … 71 71 // only expect a single entry for PSPHOT.INPUT.CMF and PSPHOT.SOURCES.TEXT, so we can only 72 72 // associate input sources with a single entry for PSPHOT.INPUT 73 bool psphotLoadExtSources (pmConfig *config, const pmFPAview *view ) {73 bool psphotLoadExtSources (pmConfig *config, const pmFPAview *view, const char *filerule) { 74 74 75 75 bool status; … … 79 79 80 80 // find the currently selected readout 81 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", index); // File of interest81 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 82 82 psAssert (file, "missing file?"); 83 83 … … 130 130 // load data from input TXT file: 131 131 { 132 pmChip *chipTXT = pmFPAfileThisChip (config->files, view, "PSPHOT.INPUT");132 pmChip *chipTXT = pmFPAfileThisChip (config->files, view, filerule); 133 133 if (!chipTXT) goto finish; 134 134 … … 167 167 168 168 // extract the input sources corresponding to this readout 169 // XXX this function needs to be updated to work with the new context of ps hot inputs169 // XXX this function needs to be updated to work with the new context of psphot inputs 170 170 psArray *psphotLoadPSFSources (pmConfig *config, const pmFPAview *view) { 171 171 … … 197 197 // psphotDetectionsFromSources to psphotSourceStats and are now stored on 198 198 // detections->newSources. 199 bool psphotRepairLoadedSources (pmConfig *config, const pmFPAview *view ) {200 201 // find the currently selected readout 202 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", 0); // File of interest199 bool psphotRepairLoadedSources (pmConfig *config, const pmFPAview *view, const char *filerule) { 200 201 // find the currently selected readout 202 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, 0); // File of interest 203 203 psAssert (file, "missing file?"); 204 204 … … 227 227 // generate the detection structure for the supplied array of sources 228 228 // XXX this currently assumes there is a single input file 229 bool psphotDetectionsFromSources (pmConfig *config, const pmFPAview *view, psArray *sources) {230 231 // find the currently selected readout 232 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", 0); // File of interest229 bool psphotDetectionsFromSources (pmConfig *config, const pmFPAview *view, const char *filerule, psArray *sources) { 230 231 // find the currently selected readout 232 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, 0); // File of interest 233 233 psAssert (file, "missing file?"); 234 234 … … 335 335 } 336 336 337 bool psphotCheckExtSources (pmConfig *config, const pmFPAview *view ) {337 bool psphotCheckExtSources (pmConfig *config, const pmFPAview *view, const char *filerule) { 338 338 339 339 bool status; … … 343 343 344 344 // find the currently selected readout 345 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", 0); // File of interest345 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, 0); // File of interest 346 346 psAssert (file, "missing file?"); 347 347 … … 373 373 374 374 // find the detections (by peak and/or footprint) in the image. 375 if (!psphotFindDetections (config, view, true)) {375 if (!psphotFindDetections (config, view, filerule, true)) { 376 376 psError(PSPHOT_ERR_CONFIG, false, "unable to find detections in this image"); 377 return psphotReadoutCleanup (config, view );377 return psphotReadoutCleanup (config, view, filerule); 378 378 } 379 379 380 380 // construct sources and measure basic stats 381 psphotSourceStats (config, view, true);381 psphotSourceStats (config, view, filerule, true); 382 382 383 383 // find blended neighbors of very saturated stars 384 psphotDeblendSatstars (config, view );384 psphotDeblendSatstars (config, view, filerule); 385 385 386 386 // mark blended peaks PS_SOURCE_BLEND 387 if (!psphotBasicDeblend (config, view )) {387 if (!psphotBasicDeblend (config, view, filerule)) { 388 388 psLogMsg ("psphot", 3, "failed on deblend analysis"); 389 return psphotReadoutCleanup (config, view );389 return psphotReadoutCleanup (config, view, filerule); 390 390 } 391 391 392 392 // classify sources based on moments, brightness 393 if (!psphotRoughClass (config, view )) {393 if (!psphotRoughClass (config, view, filerule)) { 394 394 psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image"); 395 return psphotReadoutCleanup (config, view); 396 } 397 } 398 399 return true; 400 } 395 return psphotReadoutCleanup (config, view, filerule); 396 } 397 } 398 399 return true; 400 } 401 402 // copy the detections from one pmFPAfile to another 403 bool psphotCopySources (pmConfig *config, const pmFPAview *view, const char *ruleOut, const char *ruleSrc) 404 { 405 bool status = true; 406 407 int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM"); 408 psAssert (status, "programming error: must define PSPHOT.INPUT.NUM"); 409 410 // skip the chisq image because it is a duplicate of the detection version 411 int chisqNum = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.CHISQ.NUM"); 412 if (!status) chisqNum = -1; 413 414 // loop over the available readouts 415 for (int i = 0; i < num; i++) { 416 if (i == chisqNum) continue; // skip chisq image 417 if (!psphotCopySourcesReadout (config, view, ruleOut, ruleSrc, i)) { 418 psError (PSPHOT_ERR_CONFIG, false, "failed to copy sources from %s to %s entry %d", ruleSrc, ruleOut, i); 419 return false; 420 } 421 } 422 return true; 423 } 424 425 // add newly selected sources to the existing list of sources 426 bool psphotCopySourcesReadout (pmConfig *config, const pmFPAview *view, const char *ruleOut, const char *ruleSrc, int index) { 427 428 bool status; 429 430 // find the currently selected readout 431 pmFPAfile *fileSrc = pmFPAfileSelectSingle(config->files, ruleSrc, index); // File of interest 432 psAssert (fileSrc, "missing file?"); 433 434 pmReadout *readoutSrc = pmFPAviewThisReadout(view, fileSrc->fpa); 435 psAssert (readoutSrc, "missing readout?"); 436 437 pmDetections *detections = psMetadataLookupPtr (&status, readoutSrc->analysis, "PSPHOT.DETECTIONS"); 438 psAssert (detections, "missing detections?"); 439 440 // find the currently selected readout 441 pmFPAfile *fileOut = pmFPAfileSelectSingle(config->files, ruleOut, index); // File of interest 442 psAssert (fileOut, "missing file?"); 443 444 pmReadout *readoutOut = pmFPAviewThisReadout(view, fileOut->fpa); 445 psAssert (readoutOut, "missing readout?"); 446 447 // save detections on the readout->analysis 448 // XXX this replaced any existing entry; allow this operation to merge? 449 if (!psMetadataAddPtr (readoutOut->analysis, PS_LIST_TAIL, "PSPHOT.DETECTIONS", PS_META_REPLACE | PS_DATA_UNKNOWN, "psphot detections", detections)) { 450 psError (PSPHOT_ERR_CONFIG, false, "problem saving detections on readout"); 451 return false; 452 } 453 454 return true; 455 } 456 -
branches/czw_branch/20100427/psphot/src/psphotModelBackground.c
r27657 r28029 384 384 385 385 // XXX supply filename or keep PSPHOT.INPUT fixed? 386 bool psphotModelBackground (pmConfig *config, const pmFPAview *view )386 bool psphotModelBackground (pmConfig *config, const pmFPAview *view, const char *filerule) 387 387 { 388 388 bool status = false; … … 393 393 // loop over the available readouts 394 394 for (int i = 0; i < num; i++) { 395 if (!psphotModelBackgroundReadoutFileIndex(config, view, "PSPHOT.INPUT", i)) {395 if (!psphotModelBackgroundReadoutFileIndex(config, view, filerule, i)) { 396 396 psError (PSPHOT_ERR_CONFIG, false, "failed to model background for PSPHOT.INPUT entry %d", i); 397 397 return false; -
branches/czw_branch/20100427/psphot/src/psphotModelTest.c
r26894 r28029 3 3 4 4 // XXX add more test information? 5 bool psphotModelTest (pmConfig *config, const pmFPAview *view, psMetadata *recipe) {5 bool psphotModelTest (pmConfig *config, const pmFPAview *view, const char *filerule, psMetadata *recipe) { 6 6 7 7 bool status; … … 33 33 34 34 // use poissonian errors or local-sky errors 35 bool POISSON_ERRORS = psMetadataLookupBool (&status, recipe, "POISSON_ERRORS");35 bool POISSON_ERRORS = psMetadataLookupBool (&status, recipe, filerule); 36 36 if (!status) POISSON_ERRORS = true; 37 37 pmSourceFitModelInit (15, 0.1, 1.0, POISSON_ERRORS); -
branches/czw_branch/20100427/psphot/src/psphotOutput.c
r26894 r28029 126 126 } 127 127 128 bool psphotAddPhotcode (pmConfig *config, const pmFPAview *view ) {128 bool psphotAddPhotcode (pmConfig *config, const pmFPAview *view, const char *filerule) { 129 129 130 130 bool status = false; … … 135 135 // loop over the available readouts 136 136 for (int i = 0; i < num; i++) { 137 if (!psphotAddPhotcodeReadout (config, view, "PSPHOT.INPUT", i)) {137 if (!psphotAddPhotcodeReadout (config, view, filerule, i)) { 138 138 psError (PSPHOT_ERR_CONFIG, false, "failed to add photcode to PSPHOT.INPUT entry %d", i); 139 139 return false; -
branches/czw_branch/20100427/psphot/src/psphotPetrosianRadialBins.c
r28017 r28029 182 182 psFree(values); 183 183 psFree(stats); 184 return false;184 return true; 185 185 } 186 186 -
branches/czw_branch/20100427/psphot/src/psphotPetrosianStats.c
r28017 r28029 15 15 16 16 pmSourceRadialProfile *profile = source->extpars->petProfile; 17 18 if (!profile->binSB) { 19 psLogMsg ("psphot", PS_LOG_DETAIL, "no petrosian profile, skipping source %f, %f", source->peak->xf, source->peak->yf); 20 return true; 21 } 17 22 18 23 psVector *binSB = profile->binSB; … … 190 195 source->extpars->petrosianR90Err = NAN; 191 196 192 fprintf (stderr, "source @ %f,%f\n", source->peak->xf, source->peak->yf);197 // fprintf (stderr, "source @ %f,%f\n", source->peak->xf, source->peak->yf); 193 198 psphotPetrosianVisualStats (binRad, binSB, refRadius, meanSB, petRatio, petRatioErr, fluxSum, petRadius, PETROSIAN_RATIO, petFlux, apRadius); 194 199 -
branches/czw_branch/20100427/psphot/src/psphotRadialBins.c
r28017 r28029 44 44 psVector *radMin = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.LOWER"); 45 45 psVector *radMax = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.UPPER"); 46 if (!radMin || !radMin->n) return false; 47 if (!radMax || !radMax->n) return false; 46 if (!radMin || !radMin->n) { 47 psError (PSPHOT_ERR_CONFIG, true, "error in definition of annular bins (radMin missing or empty)"); 48 return false; 49 } 50 if (!radMax || !radMax->n) { 51 psError (PSPHOT_ERR_CONFIG, true, "error in definition of annular bins (radMax missing or empty)"); 52 return false; 53 } 48 54 49 55 psVector *binSB = psVectorAllocEmpty(radMin->n, PS_TYPE_F32); // surface brightness of radial bin … … 139 145 } 140 146 } 141 binSB->n = binSBstdev->n = bin Rad->n = binArea->n = nOut;147 binSB->n = binSBstdev->n = binSum->n = binRad->n = binArea->n = nOut; 142 148 143 149 // interpolate any bins that were empty (extrapolate to center if needed) … … 155 161 psFree(values); 156 162 psFree(stats); 157 return false;163 return true; 158 164 } 159 165 -
branches/czw_branch/20100427/psphot/src/psphotReadout.c
r27657 r28029 28 28 29 29 // set the photcode for this image 30 if (!psphotAddPhotcode (config, view )) {30 if (!psphotAddPhotcode (config, view, "PSPHOT.INPUT")) { 31 31 psError (PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); 32 32 return false; … … 34 34 35 35 // Generate the mask and weight images, including the user-defined analysis region of interest 36 if (!psphotSetMaskAndVariance (config, view )) {37 return psphotReadoutCleanup(config, view );36 if (!psphotSetMaskAndVariance (config, view, "PSPHOT.INPUT")) { 37 return psphotReadoutCleanup(config, view, "PSPHOT.INPUT"); 38 38 } 39 39 if (!strcasecmp (breakPt, "NOTHING")) { 40 return psphotReadoutCleanup (config, view);40 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 41 41 } 42 42 43 43 // generate a background model (median, smoothed image) 44 if (!psphotModelBackground (config, view )) {45 return psphotReadoutCleanup (config, view );46 } 47 if (!psphotSubtractBackground (config, view )) {48 return psphotReadoutCleanup (config, view );44 if (!psphotModelBackground (config, view, "PSPHOT.INPUT")) { 45 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 46 } 47 if (!psphotSubtractBackground (config, view, "PSPHOT.INPUT")) { 48 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 49 49 } 50 50 if (!strcasecmp (breakPt, "BACKMDL")) { 51 return psphotReadoutCleanup (config, view );51 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 52 52 } 53 53 54 54 // load the psf model, if suppled. FWHM_X,FWHM_Y,etc are determined and saved on 55 55 // readout->analysis XXX this function currently only works with a single PSPHOT.INPUT 56 if (!psphotLoadPSF (config, view)) { 56 if (!psphotLoadPSF (config, view)) { // ??? need to supply 2 ? 57 57 psError (PSPHOT_ERR_UNKNOWN, false, "error loading psf model"); 58 return psphotReadoutCleanup (config, view );58 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 59 59 } 60 60 61 61 // find the detections (by peak and/or footprint) in the image. 62 if (!psphotFindDetections (config, view, true)) { // pass 162 if (!psphotFindDetections (config, view, "PSPHOT.INPUT", true)) { // pass 1 63 63 // this only happens if we had an error in psphotFindDetections 64 64 psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis"); 65 return psphotReadoutCleanup (config, view );65 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 66 66 } 67 67 68 68 // construct sources and measure basic stats (saved on detections->newSources) 69 if (!psphotSourceStats (config, view, true)) { // pass 169 if (!psphotSourceStats (config, view, "PSPHOT.INPUT", true)) { // pass 1 70 70 psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources"); 71 return psphotReadoutCleanup (config, view );71 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 72 72 } 73 73 if (!strcasecmp (breakPt, "PEAKS")) { 74 return psphotReadoutCleanup (config, view);74 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 75 75 } 76 76 77 77 // find blended neighbors of very saturated stars (detections->newSources) 78 if (!psphotDeblendSatstars (config, view )) {78 if (!psphotDeblendSatstars (config, view, "PSPHOT.INPUT")) { 79 79 psError (PSPHOT_ERR_UNKNOWN, false, "failed on satstar deblend analysis"); 80 return psphotReadoutCleanup (config, view );80 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 81 81 } 82 82 83 83 // mark blended peaks PS_SOURCE_BLEND (detections->newSources) 84 if (!psphotBasicDeblend (config, view )) {84 if (!psphotBasicDeblend (config, view, "PSPHOT.INPUT")) { 85 85 psError (PSPHOT_ERR_UNKNOWN, false, "failed on deblend analysis"); 86 return psphotReadoutCleanup (config, view );86 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 87 87 } 88 88 89 89 // classify sources based on moments, brightness. if a PSF model has been loaded, the PSF 90 90 // clump defined for it is used not measured (detections->newSources) 91 if (!psphotRoughClass (config, view )) { // pass 191 if (!psphotRoughClass (config, view, "PSPHOT.INPUT")) { // pass 1 92 92 psError (PSPHOT_ERR_UNKNOWN, false, "failed to determine rough classifications"); 93 return psphotReadoutCleanup (config, view );93 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 94 94 } 95 95 // if we were not supplied a PSF model, determine the IQ stats here (detections->newSources) 96 if (!psphotImageQuality (config, view )) { // pass 196 if (!psphotImageQuality (config, view, "PSPHOT.INPUT")) { // pass 1 97 97 psError (PSPHOT_ERR_UNKNOWN, false, "failed to measure image quality"); 98 return psphotReadoutCleanup (config, view);98 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 99 99 } 100 100 if (!strcasecmp (breakPt, "MOMENTS")) { 101 return psphotReadoutCleanup (config, view);101 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 102 102 } 103 103 104 104 // use bright stellar objects to measure PSF if we were supplied a PSF for any input file, 105 105 // this step is skipped 106 if (!psphotChoosePSF (config, view )) { // pass 1106 if (!psphotChoosePSF (config, view, "PSPHOT.INPUT")) { // pass 1 107 107 psLogMsg ("psphot", 3, "failure to construct a psf model"); 108 return psphotReadoutCleanup (config, view );108 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 109 109 } 110 110 if (!strcasecmp (breakPt, "PSFMODEL")) { 111 return psphotReadoutCleanup (config, view );111 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 112 112 } 113 113 114 114 // include externally-supplied sources 115 115 // XXX fix this in the new multi-input context 116 // psphotLoadExtSources (config, view ); // pass 1116 // psphotLoadExtSources (config, view, "PSPHOT.INPUT"); // pass 1 117 117 118 118 // construct an initial model for each object, set the radius to fitRadius, set circular 119 119 // fit mask (detections->newSources) 120 psphotGuessModels (config, view ); // pass 1120 psphotGuessModels (config, view, "PSPHOT.INPUT"); // pass 1 121 121 122 122 // merge the newly selected sources into the existing list 123 123 // NOTE: merge OLD and NEW 124 psphotMergeSources (config, view );124 psphotMergeSources (config, view, "PSPHOT.INPUT"); 125 125 126 126 // linear PSF fit to source peaks, subtract the models from the image (in PSF mask) 127 psphotFitSourcesLinear (config, view, false); // pass 1 (detections->allSources)127 psphotFitSourcesLinear (config, view, "PSPHOT.INPUT", false); // pass 1 (detections->allSources) 128 128 129 129 // identify CRs and extended sources (only unmeasured sources are measured) 130 psphotSourceSize (config, view, true); // pass 1 (detections->allSources)130 psphotSourceSize (config, view, "PSPHOT.INPUT", true); // pass 1 (detections->allSources) 131 131 if (!strcasecmp (breakPt, "ENSEMBLE")) { 132 132 goto finish; … … 135 135 // non-linear PSF and EXT fit to brighter sources 136 136 // replace model flux, adjust mask as needed, fit, subtract the models (full stamp) 137 psphotBlendFit (config, view ); // pass 1 (detections->allSources)137 psphotBlendFit (config, view, "PSPHOT.INPUT"); // pass 1 (detections->allSources) 138 138 139 139 // replace all sources 140 psphotReplaceAllSources (config, view ); // pass 1 (detections->allSources)140 psphotReplaceAllSources (config, view, "PSPHOT.INPUT"); // pass 1 (detections->allSources) 141 141 142 142 // linear fit to include all sources (subtract again) 143 143 // NOTE : apply to ALL sources (extended + psf) 144 psphotFitSourcesLinear (config, view, true); // pass 2 (detections->allSources)144 psphotFitSourcesLinear (config, view, "PSPHOT.INPUT", true); // pass 2 (detections->allSources) 145 145 146 146 // if we only do one pass, skip to extended source analysis … … 150 150 151 151 // add noise for subtracted objects 152 psphotAddNoise (config, view ); // pass 1 (detections->allSources)152 psphotAddNoise (config, view, "PSPHOT.INPUT"); // pass 1 (detections->allSources) 153 153 154 154 // find fainter sources 155 155 // NOTE: finds new peaks and new footprints, OLD and FULL set are saved on detections 156 psphotFindDetections (config, view, false); // pass 2 (detections->peaks, detections->footprints)156 psphotFindDetections (config, view, "PSPHOT.INPUT", false); // pass 2 (detections->peaks, detections->footprints) 157 157 158 158 // remove noise for subtracted objects (ie, return to normal noise level) 159 159 // NOTE: this needs to operate only on the OLD sources 160 psphotSubNoise (config, view ); // pass 1 (detections->allSources)160 psphotSubNoise (config, view, "PSPHOT.INPUT"); // pass 1 (detections->allSources) 161 161 162 162 // define new sources based on only the new peaks 163 163 // NOTE: new sources are saved on detections->newSources 164 psphotSourceStats (config, view, false); // pass 2 (detections->newSources)164 psphotSourceStats (config, view, "PSPHOT.INPUT", false); // pass 2 (detections->newSources) 165 165 166 166 // set source type 167 167 // NOTE: apply only to detections->newSources 168 if (!psphotRoughClass (config, view )) { // pass 2 (detections->newSources)168 if (!psphotRoughClass (config, view, "PSPHOT.INPUT")) { // pass 2 (detections->newSources) 169 169 psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image"); 170 return psphotReadoutCleanup (config, view );170 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 171 171 } 172 172 173 173 // create full input models, set the radius to fitRadius, set circular fit mask 174 174 // NOTE: apply only to detections->newSources 175 psphotGuessModels (config, view ); // pass 2 (detections->newSources)175 psphotGuessModels (config, view, "PSPHOT.INPUT"); // pass 2 (detections->newSources) 176 176 177 177 // replace all sources so fit below applies to all at once 178 178 // NOTE: apply only to OLD sources (which have been subtracted) 179 psphotReplaceAllSources (config, view ); // pass 2179 psphotReplaceAllSources (config, view, "PSPHOT.INPUT"); // pass 2 180 180 181 181 // merge the newly selected sources into the existing list 182 182 // NOTE: merge OLD and NEW 183 183 // XXX check on free of sources... 184 psphotMergeSources (config, view ); // (detections->newSources + detections->allSources -> detections->allSources)184 psphotMergeSources (config, view, "PSPHOT.INPUT"); // (detections->newSources + detections->allSources -> detections->allSources) 185 185 186 186 // NOTE: apply to ALL sources 187 psphotFitSourcesLinear (config, view, true); // pass 3 (detections->allSources)187 psphotFitSourcesLinear (config, view, "PSPHOT.INPUT", true); // pass 3 (detections->allSources) 188 188 189 189 pass1finish: … … 191 191 // measure source size for the remaining sources 192 192 // NOTE: applies only to NEW (unmeasured) sources 193 psphotSourceSize (config, view, false); // pass 2 (detections->allSources)194 195 psphotExtendedSourceAnalysis (config, view ); // pass 1 (detections->allSources)196 psphotExtendedSourceFits (config, view ); // pass 1 (detections->allSources)193 psphotSourceSize (config, view, "PSPHOT.INPUT", false); // pass 2 (detections->allSources) 194 195 psphotExtendedSourceAnalysis (config, view, "PSPHOT.INPUT"); // pass 1 (detections->allSources) 196 psphotExtendedSourceFits (config, view, "PSPHOT.INPUT"); // pass 1 (detections->allSources) 197 197 198 198 finish: … … 202 202 203 203 // measure aperture photometry corrections 204 if (!psphotApResid (config, view )) { // pass 1 (detections->allSources)204 if (!psphotApResid (config, view, "PSPHOT.INPUT")) { // pass 1 (detections->allSources) 205 205 psLogMsg ("psphot", 3, "failed on psphotApResid"); 206 return psphotReadoutCleanup (config, view );206 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 207 207 } 208 208 209 209 // calculate source magnitudes 210 psphotMagnitudes(config, view ); // pass 1 (detections->allSources)211 212 if (!psphotEfficiency(config, view )) { // pass 1210 psphotMagnitudes(config, view, "PSPHOT.INPUT"); // pass 1 (detections->allSources) 211 212 if (!psphotEfficiency(config, view, "PSPHOT.INPUT")) { // pass 1 213 213 psErrorStackPrint(stderr, "Unable to determine detection efficiencies from fake sources"); 214 214 psErrorClear(); … … 219 219 220 220 // replace background in residual image 221 psphotSkyReplace (config, view ); // pass 1221 psphotSkyReplace (config, view, "PSPHOT.INPUT"); // pass 1 222 222 223 223 // drop the references to the image pixels held by each source 224 psphotSourceFreePixels (config, view ); // pass 1224 psphotSourceFreePixels (config, view, "PSPHOT.INPUT"); // pass 1 225 225 226 226 // create the exported-metadata and free local data 227 return psphotReadoutCleanup(config, view );227 return psphotReadoutCleanup(config, view, "PSPHOT.INPUT"); 228 228 } -
branches/czw_branch/20100427/psphot/src/psphotReadoutCleanup.c
r27657 r28029 2 2 3 3 // for now, let's store the detections on the readout->analysis for each readout 4 bool psphotReadoutCleanup (pmConfig *config, const pmFPAview *view )4 bool psphotReadoutCleanup (pmConfig *config, const pmFPAview *view, const char *filerule) 5 5 { 6 6 bool status = true; … … 24 24 // loop over the available readouts 25 25 for (int i = 0; i < num; i++) { 26 if (!psphotReadoutCleanupReadout (config, view, "PSPHOT.INPUT", i, recipe)) {27 psError (PSPHOT_ERR_CONFIG, false, "failed on psphotReadoutCleanup for PSPHOT.INPUT entry %d", i);26 if (!psphotReadoutCleanupReadout (config, view, filerule, i, recipe)) { 27 psError (PSPHOT_ERR_CONFIG, false, "failed on psphotReadoutCleanup for %s entry %d", filerule, i); 28 28 return false; 29 29 } … … 39 39 // not a DATA error, then there was a serious problem. Only in this case, or if the fail 40 40 // on the stats measurement, do we return false 41 bool psphotReadoutCleanupReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe) {41 bool psphotReadoutCleanupReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { 42 42 43 43 bool status = true; 44 44 45 45 // find the currently selected readout 46 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest46 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 47 47 psAssert (file, "missing file?"); 48 48 -
branches/czw_branch/20100427/psphot/src/psphotReadoutFindPSF.c
r26894 r28029 8 8 9 9 // set the photcode for the PSPHOT.INPUT 10 if (!psphotAddPhotcode(config, view )) {10 if (!psphotAddPhotcode(config, view, "PSPHOT.INPUT")) { 11 11 psError(PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); 12 12 return false; … … 14 14 15 15 // Generate the mask and variance images, including the user-defined analysis region of interest 16 psphotSetMaskAndVariance (config, view );16 psphotSetMaskAndVariance (config, view, "PSPHOT.INPUT"); 17 17 18 18 // Note that in this implementation, we do NOT model the background and we do not … … 21 21 // include externally-supplied sources (supplied as PSPHOT.INPUT.CMF) 22 22 // XXX we assume a single set of input sources is supplied 23 if (!psphotDetectionsFromSources (config, view, inSources)) {23 if (!psphotDetectionsFromSources (config, view, "PSPHOT.INPUT", inSources)) { 24 24 psError(PSPHOT_ERR_ARGUMENTS, true, "Can't find PSF stars"); 25 return psphotReadoutCleanup (config, view);25 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 26 26 } 27 27 28 28 // construct detections->newSources and measure basic stats (moments, local sky) 29 if (!psphotSourceStats(config, view, true)) {29 if (!psphotSourceStats(config, view, "PSPHOT.INPUT", true)) { 30 30 psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources"); 31 31 return false; … … 33 33 34 34 // peak flux is wrong : use the peak measured in the moments analysis: 35 if (!psphotRepairLoadedSources(config, view )) {35 if (!psphotRepairLoadedSources(config, view, "PSPHOT.INPUT")) { 36 36 psError(PSPHOT_ERR_UNKNOWN, false, "failure to repair sources"); 37 37 return false; … … 39 39 40 40 // classify sources based on moments, brightness (psf is not known) 41 if (!psphotRoughClass (config, view )) {41 if (!psphotRoughClass (config, view, "PSPHOT.INPUT")) { 42 42 psError (PSPHOT_ERR_UNKNOWN, false, "failed to determine rough source class"); 43 return psphotReadoutCleanup (config, view );43 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 44 44 } 45 45 46 if (!psphotImageQuality (config, view )) {46 if (!psphotImageQuality (config, view, "PSPHOT.INPUT")) { 47 47 psError (PSPHOT_ERR_UNKNOWN, false, "failed to measure image quality"); 48 return psphotReadoutCleanup (config, view);48 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 49 49 } 50 50 51 if (!psphotChoosePSF(config, view )) {51 if (!psphotChoosePSF(config, view, "PSPHOT.INPUT")) { 52 52 psError(PSPHOT_ERR_PSF, false, "Failed to construct a psf model"); 53 return psphotReadoutCleanup (config, view);53 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 54 54 } 55 55 … … 59 59 // fits from that analysis, or run the linear PSF fit for all objects currently in hand 60 60 // construct an initial model for each object, set the radius to fitRadius, set circular fit mask 61 psphotGuessModels (config, view );61 psphotGuessModels (config, view, "PSPHOT.INPUT"); 62 62 # endif 63 63 64 64 // merge the newly selected sources into the existing list 65 65 // NOTE: merge OLD and NEW 66 psphotMergeSources (config, view );66 psphotMergeSources (config, view, "PSPHOT.INPUT"); 67 67 68 68 # if 0 69 69 // measure aperture photometry corrections 70 if (!psphotApResid (config, view )) {70 if (!psphotApResid (config, view, "PSPHOT.INPUT")) { 71 71 psLogMsg ("psphot", 3, "failed on psphotApResid"); 72 return psphotReadoutCleanup (config, view );72 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 73 73 } 74 74 # endif 75 75 76 76 // drop the references to the image pixels held by each source 77 psphotSourceFreePixels(config, view );77 psphotSourceFreePixels(config, view, "PSPHOT.INPUT"); 78 78 79 79 // create the exported-metadata and free local data 80 return psphotReadoutCleanup (config, view);80 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 81 81 } -
branches/czw_branch/20100427/psphot/src/psphotReadoutKnownSources.c
r26894 r28029 8 8 9 9 // set the photcode for this image 10 if (!psphotAddPhotcode(config, view )) {10 if (!psphotAddPhotcode(config, view, "PSPHOT.INPUT")) { 11 11 psError(PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); 12 12 return false; … … 14 14 15 15 // Generate the mask and weight images, including the user-defined analysis region of interest 16 psphotSetMaskAndVariance (config, view );16 psphotSetMaskAndVariance (config, view, "PSPHOT.INPUT"); 17 17 18 18 // Note that in this implementation, we do NOT model the background and we do not … … 20 20 21 21 // include externally-supplied sources (supplied as PSPHOT.INPUT.CMF) 22 if (!psphotDetectionsFromSources (config, view, inSources)) {22 if (!psphotDetectionsFromSources (config, view, "PSPHOT.INPUT", inSources)) { 23 23 psError(PSPHOT_ERR_ARGUMENTS, true, "Can't find PSF stars"); 24 return psphotReadoutCleanup (config, view);24 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 25 25 } 26 26 27 27 // construct sources and measure basic stats 28 if (!psphotSourceStats (config, view, true)) {28 if (!psphotSourceStats (config, view, "PSPHOT.INPUT", true)) { 29 29 psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources"); 30 30 return false; … … 32 32 33 33 // peak flux is wrong : use the peak measured in the moments analysis: 34 if (!psphotRepairLoadedSources(config, view )) {34 if (!psphotRepairLoadedSources(config, view, "PSPHOT.INPUT")) { 35 35 psError(PSPHOT_ERR_UNKNOWN, false, "failure to repair sources"); 36 36 return false; … … 38 38 39 39 // classify sources based on moments, brightness (psf is not known) 40 if (!psphotRoughClass (config, view )) {40 if (!psphotRoughClass (config, view, "PSPHOT.INPUT")) { 41 41 psError (PSPHOT_ERR_UNKNOWN, false, "failed to determine rough source class"); 42 return psphotReadoutCleanup (config, view);42 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 43 43 } 44 44 45 if (!psphotChoosePSF (config, view )) {45 if (!psphotChoosePSF (config, view, "PSPHOT.INPUT")) { 46 46 psError(PSPHOT_ERR_PSF, false, "Failed to construct a psf model"); 47 return psphotReadoutCleanup (config, view);47 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 48 48 } 49 49 50 50 // construct an initial model for each object 51 psphotGuessModels (config, view );51 psphotGuessModels (config, view, "PSPHOT.INPUT"); 52 52 53 53 // merge the newly selected sources into the existing list 54 54 // NOTE: merge OLD and NEW 55 psphotMergeSources (config, view );55 psphotMergeSources (config, view, "PSPHOT.INPUT"); 56 56 57 57 // linear PSF fit to source peaks 58 psphotFitSourcesLinear (config, view, false);58 psphotFitSourcesLinear (config, view, "PSPHOT.INPUT", false); 59 59 60 60 // measure aperture photometry corrections 61 if (!psphotApResid (config, view )) {61 if (!psphotApResid (config, view, "PSPHOT.INPUT")) { 62 62 psLogMsg ("psphot", 3, "failed on psphotApResid"); 63 return psphotReadoutCleanup (config, view);63 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 64 64 } 65 65 66 66 // calculate source magnitudes 67 psphotMagnitudes(config, view );67 psphotMagnitudes(config, view, "PSPHOT.INPUT"); 68 68 69 69 // drop the references to the image pixels held by each source 70 psphotSourceFreePixels (config, view );70 psphotSourceFreePixels (config, view, "PSPHOT.INPUT"); 71 71 72 72 // create the exported-metadata and free local data 73 return psphotReadoutCleanup (config, view);73 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 74 74 } -
branches/czw_branch/20100427/psphot/src/psphotReadoutMinimal.c
r26894 r28029 18 18 19 19 // set the photcode for this image 20 if (!psphotAddPhotcode(config, view )) {20 if (!psphotAddPhotcode(config, view, "PSPHOT.INPUT")) { 21 21 psError(PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); 22 22 return false; … … 24 24 25 25 // Generate the mask and weight images, including the user-defined analysis region of interest 26 psphotSetMaskAndVariance (config, view );26 psphotSetMaskAndVariance (config, view, "PSPHOT.INPUT"); 27 27 28 28 // load the psf model, if suppled. FWHM_X,FWHM_Y,etc are saved on readout->analysis 29 29 if (!psphotLoadPSF (config, view)) { 30 30 psError (PSPHOT_ERR_CONFIG, false, "missing psf model"); 31 return psphotReadoutCleanup (config, view );31 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 32 32 } 33 33 34 34 // find the detections (by peak and/or footprint) in the image. (final pass) 35 if (!psphotFindDetections(config, view, false)) {35 if (!psphotFindDetections(config, view, "PSPHOT.INPUT", false)) { 36 36 psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis"); 37 return psphotReadoutCleanup (config, view );37 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 38 38 } 39 39 40 40 // construct sources and measure basic stats (saved on detections->newSources) 41 if (!psphotSourceStats (config, view, false)) { // pass 141 if (!psphotSourceStats (config, view, "PSPHOT.INPUT", false)) { // pass 1 42 42 psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources"); 43 return psphotReadoutCleanup (config, view );43 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 44 44 } 45 45 46 46 // find blended neighbors of very saturated stars 47 psphotDeblendSatstars (config, view );47 psphotDeblendSatstars (config, view, "PSPHOT.INPUT"); 48 48 49 49 // mark blended peaks PS_SOURCE_BLEND 50 if (!psphotBasicDeblend (config, view )) {50 if (!psphotBasicDeblend (config, view, "PSPHOT.INPUT")) { 51 51 psLogMsg ("psphot", 3, "failed on deblend analysis"); 52 return psphotReadoutCleanup (config, view );52 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 53 53 } 54 54 55 55 // classify sources based on moments, brightness (use supplied psf shape parameters) 56 if (!psphotRoughClass (config, view )) {56 if (!psphotRoughClass (config, view, "PSPHOT.INPUT")) { 57 57 psLogMsg ("psphot", 3, "failed to find a valid PSF clump for image"); 58 return psphotReadoutCleanup (config, view );58 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 59 59 } 60 60 61 61 // construct an initial model for each object 62 psphotGuessModels (config, view );62 psphotGuessModels (config, view, "PSPHOT.INPUT"); 63 63 64 64 // merge the newly selected sources into the existing list 65 psphotMergeSources (config, view );65 psphotMergeSources (config, view, "PSPHOT.INPUT"); 66 66 67 67 // linear PSF fit to source peaks 68 psphotFitSourcesLinear (config, view, false);68 psphotFitSourcesLinear (config, view, "PSPHOT.INPUT", false); 69 69 70 70 // XXX eventually, add the extended source fits here 71 71 # if (0) 72 72 // measure source size for the remaining sources 73 psphotSourceSize (config, view );73 psphotSourceSize (config, view, "PSPHOT.INPUT"); 74 74 75 psphotExtendedSourceAnalysis (config, view );75 psphotExtendedSourceAnalysis (config, view, "PSPHOT.INPUT"); 76 76 77 psphotExtendedSourceFits (config, view );77 psphotExtendedSourceFits (config, view, "PSPHOT.INPUT"); 78 78 # endif 79 79 80 80 // calculate source magnitudes 81 psphotMagnitudes(config, view );81 psphotMagnitudes(config, view, "PSPHOT.INPUT"); 82 82 83 83 // XXX ensure this is measured if the analysis succeeds (even if quality is low) 84 if (!psphotEfficiency(config, view )) {84 if (!psphotEfficiency(config, view, "PSPHOT.INPUT")) { 85 85 psErrorStackPrint(stderr, "Unable to determine detection efficiencies from fake sources"); 86 86 psErrorClear(); … … 88 88 89 89 // drop the references to the image pixels held by each source 90 psphotSourceFreePixels (config, view );90 psphotSourceFreePixels (config, view, "PSPHOT.INPUT"); 91 91 92 92 // create the exported-metadata and free local data 93 return psphotReadoutCleanup (config, view);93 return psphotReadoutCleanup (config, view, "PSPHOT.INPUT"); 94 94 } -
branches/czw_branch/20100427/psphot/src/psphotReplaceUnfit.c
r26894 r28029 23 23 24 24 // for now, let's store the detections on the readout->analysis for each readout 25 bool psphotReplaceAllSources (pmConfig *config, const pmFPAview *view )25 bool psphotReplaceAllSources (pmConfig *config, const pmFPAview *view, const char *filerule) 26 26 { 27 27 bool status = true; … … 36 36 // loop over the available readouts 37 37 for (int i = 0; i < num; i++) { 38 if (!psphotReplaceAllSourcesReadout (config, view, "PSPHOT.INPUT", i, recipe)) {38 if (!psphotReplaceAllSourcesReadout (config, view, filerule, i, recipe)) { 39 39 psError (PSPHOT_ERR_CONFIG, false, "failed to replace all sources for PSPHOT.INPUT entry %d", i); 40 40 return false; -
branches/czw_branch/20100427/psphot/src/psphotRoughClass.c
r27657 r28029 8 8 9 9 // for now, let's store the detections on the readout->analysis for each readout 10 bool psphotRoughClass (pmConfig *config, const pmFPAview *view )10 bool psphotRoughClass (pmConfig *config, const pmFPAview *view, const char *filerule) 11 11 { 12 12 bool status = true; … … 26 26 for (int i = 0; i < num; i++) { 27 27 if (i == chisqNum) continue; // skip chisq image 28 if (!psphotRoughClassReadout (config, view, "PSPHOT.INPUT", i, recipe)) {29 psError (PSPHOT_ERR_CONFIG, false, "failed on rough classification for PSPHOT.INPUT entry %d", i);28 if (!psphotRoughClassReadout (config, view, filerule, i, recipe)) { 29 psError (PSPHOT_ERR_CONFIG, false, "failed on rough classification for %s entry %d", filerule, i); 30 30 return false; 31 31 } … … 34 34 } 35 35 36 bool psphotRoughClassReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe) {36 bool psphotRoughClassReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { 37 37 38 38 bool status; … … 41 41 42 42 // find the currently selected readout 43 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest43 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 44 44 psAssert (file, "missing file?"); 45 45 -
branches/czw_branch/20100427/psphot/src/psphotSetMaskBits.c
r21254 r28029 37 37 return true; 38 38 } 39 40 // XXX should these be in config->analysis or somewhere else besides 'recipe'? -
branches/czw_branch/20100427/psphot/src/psphotSkyReplace.c
r27657 r28029 1 1 # include "psphotInternal.h" 2 2 3 bool psphotSkyReplace (pmConfig *config, const pmFPAview *view )3 bool psphotSkyReplace (pmConfig *config, const pmFPAview *view, const char *filerule) 4 4 { 5 5 bool status = true; … … 15 15 for (int i = 0; i < num; i++) { 16 16 if (i == chisqNum) continue; // skip chisq image 17 if (!psphotSkyReplaceReadout (config, view, "PSPHOT.INPUT", i)) {18 psError (PSPHOT_ERR_CONFIG, false, "failed to replace sky for PSPHOT.INPUT entry %d", i);17 if (!psphotSkyReplaceReadout (config, view, filerule, i)) { 18 psError (PSPHOT_ERR_CONFIG, false, "failed to replace sky for %s entry %d", filerule, i); 19 19 return false; 20 20 } … … 25 25 // XXX make this an option? 26 26 // in order to successfully replace the sky, we must define a corresponding file... 27 bool psphotSkyReplaceReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index) {27 bool psphotSkyReplaceReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index) { 28 28 29 29 psTimerStart ("psphot.skyreplace"); 30 30 31 31 // find the currently selected readout 32 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest32 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 33 33 psAssert (file, "missing file?"); 34 34 -
branches/czw_branch/20100427/psphot/src/psphotSourceFreePixels.c
r26894 r28029 1 1 # include "psphotInternal.h" 2 2 3 bool psphotSourceFreePixels (pmConfig *config, const pmFPAview *view )3 bool psphotSourceFreePixels (pmConfig *config, const pmFPAview *view, const char *filerule) 4 4 { 5 5 bool status = true; … … 10 10 // loop over the available readouts 11 11 for (int i = 0; i < num; i++) { 12 if (!psphotSourceFreePixelsReadout (config, view, "PSPHOT.INPUT", i)) {13 psError (PSPHOT_ERR_CONFIG, false, "failed to free source pixels for PSPHOT.INPUT entry %d", i);12 if (!psphotSourceFreePixelsReadout (config, view, filerule, i)) { 13 psError (PSPHOT_ERR_CONFIG, false, "failed to free source pixels for %s entry %d", filerule, i); 14 14 return false; 15 15 } … … 18 18 } 19 19 20 bool psphotSourceFreePixelsReadout(pmConfig *config, const pmFPAview *view, const char *file name, int index) {20 bool psphotSourceFreePixelsReadout(pmConfig *config, const pmFPAview *view, const char *filerule, int index) { 21 21 22 22 bool status; 23 23 24 24 // find the currently selected readout 25 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest25 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 26 26 psAssert (file, "missing file?"); 27 27 -
branches/czw_branch/20100427/psphot/src/psphotSourceMatch.c
r27657 r28029 1 1 # include "psphotInternal.h" 2 2 3 bool psphotMatchSourcesGenerate (pmConfig *config, const pmFPAview *view, psArray *objects); 4 5 psArray *psphotMatchSources (pmConfig *config, const pmFPAview *view) 3 bool psphotMatchSourcesAddMissing (pmConfig *config, const pmFPAview *view, const char *filerule, psArray *objects); 4 bool psphotMatchSourcesSetIDs (psArray *objects); 5 6 psArray *psphotMatchSources (pmConfig *config, const pmFPAview *view, const char *filerule) 6 7 { 7 8 bool status = true; … … 14 15 // loop over the available readouts 15 16 for (int i = 0; i < num; i++) { 16 if (!psphotMatchSourcesReadout (objects, config, view, "PSPHOT.INPUT", i)) {17 if (!psphotMatchSourcesReadout (objects, config, view, filerule, i)) { 17 18 psError (PSPHOT_ERR_CONFIG, false, "failed to merge sources for PSPHOT.INPUT entry %d", i); 18 19 psFree (objects); … … 21 22 } 22 23 23 psphotMatchSourcesGenerate (config, view, objects); 24 // create sources for images where an object has been detected in the other images 25 psphotMatchSourcesAddMissing (config, view, filerule, objects); 26 27 // choose a consistent position; set common sequence values 28 psphotMatchSourcesSetIDs (objects); 24 29 25 30 return objects; 26 31 } 27 32 28 bool psphotMatchSourcesReadout (psArray *objects, pmConfig *config, const pmFPAview *view, c har *filename, int index) {33 bool psphotMatchSourcesReadout (psArray *objects, pmConfig *config, const pmFPAview *view, const char *filerule, int index) { 29 34 30 35 bool status = false; … … 38 43 39 44 // find the currently selected readout 40 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest45 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 41 46 psAssert (file, "missing file?"); 42 47 … … 145 150 } 146 151 147 bool psphotMatchSources Generate (pmConfig *config, const pmFPAview *view, psArray *objects) {152 bool psphotMatchSourcesAddMissing (pmConfig *config, const pmFPAview *view, const char *filerule, psArray *objects) { 148 153 149 154 bool status = false; … … 167 172 168 173 // find the currently selected readout 169 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", i); // File of interest174 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, i); // File of interest 170 175 psAssert (file, "missing file?"); 171 176 … … 255 260 return true; 256 261 } 262 263 bool psphotMatchSourcesSetIDs (psArray *objects) { 264 265 for (int i = 0; i < objects->n; i++) { 266 pmPhotObj *obj = objects->data[i]; 267 268 // set the source->seq values 269 for (int j = 0; j < obj->sources->n; j++) { 270 pmSource *src = obj->sources->data[j]; 271 src->seq = i; 272 } 273 } 274 return true; 275 } -
branches/czw_branch/20100427/psphot/src/psphotSourceSize.c
r27695 r28029 33 33 34 34 // for now, let's store the detections on the readout->analysis for each readout 35 bool psphotSourceSize (pmConfig *config, const pmFPAview *view, bool getPSFsize)35 bool psphotSourceSize (pmConfig *config, const pmFPAview *view, const char *filerule, bool getPSFsize) 36 36 { 37 37 bool status = true; … … 51 51 for (int i = 0; i < num; i++) { 52 52 if (i == chisqNum) continue; // skip chisq image 53 if (!psphotSourceSizeReadout (config, view, "PSPHOT.INPUT", i, recipe, getPSFsize)) {53 if (!psphotSourceSizeReadout (config, view, filerule, i, recipe, getPSFsize)) { 54 54 psError (PSPHOT_ERR_CONFIG, false, "failed on source size analysis for PSPHOT.INPUT entry %d", i); 55 55 return false; … … 60 60 61 61 // this function use an internal flag to mark sources which have already been measured 62 bool psphotSourceSizeReadout(pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe, bool getPSFsize)62 bool psphotSourceSizeReadout(pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool getPSFsize) 63 63 { 64 64 bool status; … … 68 68 69 69 // find the currently selected readout 70 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest70 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 71 71 psAssert (file, "missing file?"); 72 72 -
branches/czw_branch/20100427/psphot/src/psphotSourceStats.c
r27657 r28029 5 5 // The new sources are added to any existing sources on detections->newSources. The sources 6 6 // on detections->allSources are ignored. 7 bool psphotSourceStats (pmConfig *config, const pmFPAview *view, bool setWindow)7 bool psphotSourceStats (pmConfig *config, const pmFPAview *view, const char *filerule, bool setWindow) 8 8 { 9 9 bool status = true; … … 16 16 psAssert (status, "programming error: must define PSPHOT.INPUT.NUM"); 17 17 18 // skip the chisq image (optionally?) 19 // int chisqNum = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.CHISQ.NUM"); 20 // if (!status) chisqNum = -1; 21 18 22 // loop over the available readouts 19 23 for (int i = 0; i < num; i++) { 20 if (!psphotSourceStatsReadout (config, view, "PSPHOT.INPUT", i, recipe, setWindow)) { 21 psError (PSPHOT_ERR_CONFIG, false, "failed to find initial detections for PSPHOT.INPUT entry %d", i); 24 // if (i == chisqNum) continue; // skip chisq image 25 if (!psphotSourceStatsReadout (config, view, filerule, i, recipe, setWindow)) { 26 psError (PSPHOT_ERR_CONFIG, false, "failed to find initial detections for %s entry %d", filerule, i); 22 27 return false; 23 28 } … … 26 31 } 27 32 28 bool psphotSourceStatsReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe, bool setWindow) {33 bool psphotSourceStatsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe, bool setWindow) { 29 34 30 35 bool status = false; … … 34 39 35 40 // find the currently selected readout 36 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest41 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 37 42 psAssert (file, "missing file?"); 38 43 -
branches/czw_branch/20100427/psphot/src/psphotStackArguments.c
r27657 r28029 1 1 # include "psphotStandAlone.h" 2 2 3 static void dumpTemplate(void); 3 4 static void usage(pmConfig *config, int exitCode); 4 5 static void writeHelpInfo(FILE* ofile); … … 11 12 if (psArgumentGet(argc, argv, "-help")) writeHelpInfo(stdout); 12 13 if (psArgumentGet(argc, argv, "-h")) writeHelpInfo(stdout); 14 15 if (psArgumentGet(argc, argv, "-template")) dumpTemplate(); 13 16 14 17 // load config data from default locations … … 84 87 "\n" 85 88 "where INPUTS.mdc contains various METADATAs, each with:\n" 86 "\tIMAGE(STR): Image filename\n" 87 "\tMASK(STR): Mask filename\n" 88 "\tVARIANCE(STR): Variance map filename\n" 89 "\tIMAGE : Image filename\n" 90 "\tMASK : Mask filename\n" 91 "\tVARIANCE : Variance map filename\n" 92 "(use -template to generate a sample input.mdc file)\n" 89 93 "OUTROOT is the 'root name' for output files\n" 90 94 "\n" … … 144 148 } 145 149 150 static void dumpTemplate(void) { 151 152 fprintf (stdout, "# this line is required for multiple INPUT blocks to be accepted\n"); 153 fprintf (stdout, "INPUT MULTI\n\n"); 154 155 fprintf (stdout, "# copy and repeat the following block as needed (one per input image set)\n"); 156 fprintf (stdout, "# either RAW (unconvolved) or CNV (convolved) input images are required\n"); 157 fprintf (stdout, "# if both are supplied, by default RAW is used for detection, CNV is convolved (further) to match target PSF\n"); 158 fprintf (stdout, "# if MASK or VARIANCE images are not supplied, they will be generated\n"); 159 fprintf (stdout, "# if MASK or VARIANCE images are not supplied, they will be generated\n"); 160 fprintf (stdout, "# PSF may be supplied for the convolution target\n"); 161 fprintf (stdout, "INPUT METADATA\n"); 162 fprintf (stdout, " RAW:IMAGE STR file.im.fits # signal image filename\n"); 163 fprintf (stdout, " RAW:MASK STR file.mk.fits # mask image filename\n"); 164 fprintf (stdout, " RAW:VARIANCE STR file.wt.fits # variance image filename\n"); 165 fprintf (stdout, " RAW:PSF STR file.psf.fits # psf from input unconvolved image\n"); 166 167 fprintf (stdout, " CNV:IMAGE STR file.im.fits # signal image filename\n"); 168 fprintf (stdout, " CNV:MASK STR file.mk.fits # mask image filename\n"); 169 fprintf (stdout, " CNV:VARIANCE STR file.wt.fits # variance image filename\n"); 170 fprintf (stdout, " CNV:PSF STR file.psf.fits # psf from input convolved image\n"); 171 172 fprintf (stdout, " SOURCES STR file.cmf # measured source positions\n"); 173 fprintf (stdout, "END\n"); 174 175 psLibFinalize(); 176 exit(PS_EXIT_SUCCESS); 177 } -
branches/czw_branch/20100427/psphot/src/psphotStackChisqImage.c
r27657 r28029 6 6 7 7 // XXX supply filename or keep PSPHOT.INPUT fixed? 8 bool psphotStackChisqImage (pmConfig *config, const pmFPAview *view )8 bool psphotStackChisqImage (pmConfig *config, const pmFPAview *view, const char *ruleDet, const char *ruleCnv) 9 9 { 10 10 bool status = false; … … 21 21 22 22 // loop over the available readouts 23 // generate the chisq image from the 'detection' images 23 24 for (int i = 0; i < num; i++) { 24 if (!psphotStackChisqImageAddReadout(config, view, &chiReadout, "PSPHOT.INPUT", i)) {25 psError (PSPHOT_ERR_CONFIG, false, "failed to model background for PSPHOT.INPUT entry %d", i);25 if (!psphotStackChisqImageAddReadout(config, view, ruleDet, &chiReadout, i)) { 26 psError (PSPHOT_ERR_CONFIG, false, "failed to model background for %s entry %d", ruleDet, i); 26 27 return false; 27 28 } … … 35 36 psMetadataAddS32(config->arguments, PS_LIST_TAIL, "PSPHOT.INPUT.NUM", PS_META_REPLACE, "", num); 36 37 37 // need to save the resulting image somewhere (PSPHOT.INPUT?) 38 if (!psMetadataAddPtr(config->files, PS_LIST_TAIL, "PSPHOT.INPUT", PS_DATA_UNKNOWN | PS_META_DUPLICATE_OK, "", chisqFile)) { 38 // save the resulting image in the 'detection' set 39 if (!psMetadataAddPtr(config->files, PS_LIST_TAIL, ruleDet, PS_DATA_UNKNOWN | PS_META_DUPLICATE_OK, "", chisqFile)) { 40 psError(PM_ERR_CONFIG, false, "could not add chisqFPA to config files"); 41 return false; 42 } 43 44 // save the resulting image in the 'convolved' set 45 if (!psMetadataAddPtr(config->files, PS_LIST_TAIL, ruleCnv, PS_DATA_UNKNOWN | PS_META_DUPLICATE_OK, "", chisqFile)) { 39 46 psError(PM_ERR_CONFIG, false, "could not add chisqFPA to config files"); 40 47 return false; … … 48 55 bool psphotStackChisqImageAddReadout(const pmConfig *config, // Configuration 49 56 const pmFPAview *view, 57 const char *filerule, 50 58 pmReadout **chiReadout, 51 char *filename,52 59 int index) 53 60 { … … 55 62 56 63 // find the currently selected readout 57 pmFPAfile *input = pmFPAfileSelectSingle(config->files, file name, index); // File of interest64 pmFPAfile *input = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 58 65 psAssert (input, "missing file?"); 59 66 … … 111 118 } 112 119 113 bool psphotStackRemoveChisqFromInputs (pmConfig *config ) {120 bool psphotStackRemoveChisqFromInputs (pmConfig *config, const char *filerule) { 114 121 115 122 bool status = false; … … 121 128 psAssert (status, "programming error: must define PSPHOT.INPUT.NUM"); 122 129 123 pmFPAfileRemoveSingle (config->files, "PSPHOT.INPUT", chisqNum);130 pmFPAfileRemoveSingle (config->files, filerule, chisqNum); 124 131 125 132 inputNum --; -
branches/czw_branch/20100427/psphot/src/psphotStackImageLoop.c
r28017 r28029 1 1 # include "psphotStandAlone.h" 2 3 # define ESCAPE(MESSAGE) { \ 4 psError(PSPHOT_ERR_DATA, false, MESSAGE); \ 5 psFree (view); \ 6 return false; \ 7 } 2 #define WCS_NONLIN_TOL 0.001 // Non-linear tolerance for header WCS 3 4 # define ESCAPE(MESSAGE) { \ 5 psError(PSPHOT_ERR_DATA, false, MESSAGE); \ 6 psFree (view); \ 7 return false; \ 8 } 9 10 bool GetAstrometryFPA (pmConfig *config, pmFPAview *view); 11 bool SetAstrometryFPA (pmConfig *config, pmFPAview *view, bool bilevelAstrometry); 12 bool GetAstrometryChip (pmConfig *config, pmFPAview *view, bool bilevelAstrometry); 8 13 9 14 bool psphotStackImageLoop (pmConfig *config) { … … 14 19 pmReadout *readout; 15 20 16 pmFPAfile *input = psMetadataLookupPtr (&status, config->files, "PSPHOT.INPUT"); 17 if (!status) { 21 pmFPAfile *inputRaw = psMetadataLookupPtr (&status, config->files, "PSPHOT.STACK.INPUT.RAW"); 22 pmFPAfile *inputCnv = psMetadataLookupPtr (&status, config->files, "PSPHOT.STACK.INPUT.CNV"); 23 pmFPAfile *input = inputRaw ? inputRaw : inputCnv; 24 25 if (!input) { 18 26 psError(PSPHOT_ERR_PROG, false, "Can't find input data!"); 19 27 return false; … … 24 32 // XXX for now, just load the full set of images up front 25 33 if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) ESCAPE ("failed input for fpa in psphot."); 34 35 bool bilevelAstrometry = GetAstrometryFPA(config, view); 26 36 27 37 // for psphot, we force data to be read at the chip level … … 41 51 if (! readout->data_exists) { continue; } 42 52 43 # if (0) 44 // uncomment to generate matched psfs 53 // PSF matching 45 54 if (!psphotStackMatchPSFs (config, view)) { 46 55 psError(psErrorCodeLast(), false, "failure in psphotStackMatchPSFs for chip %d, cell %d, readout %d\n", view->chip, view->cell, view->readout); … … 48 57 return false; 49 58 } 50 # endif51 59 52 60 // XXX for now, we assume there is only a single chip in the PHU: … … 69 77 } 70 78 } 79 80 GetAstrometryChip(config, view, bilevelAstrometry); 81 71 82 // save output which is saved at the chip level 72 83 if (!pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) ESCAPE ("failed output for Chip in psphot."); 73 84 } 85 86 SetAstrometryFPA(config, view, bilevelAstrometry); 87 74 88 // save output which is saved at the fpa level 75 89 if (!pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) ESCAPE ("failed ouput for FPA in psphot."); … … 87 101 */ 88 102 103 # define UPDATE_HEADER 0 104 105 bool GetAstrometryFPA (pmConfig *config, pmFPAview *view) { 106 107 bool status = false; 108 109 bool foundAstrometry = false; 110 bool bilevelAstrometry = false; 111 112 int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM"); 113 psAssert (status, "programming error: must define PSPHOT.INPUT.NUM"); 114 115 // loop over the available readouts 116 for (int i = 0; i < num; i++) { 117 118 // find the currently selected readout 119 pmFPAfile *output = pmFPAfileSelectSingle(config->files, "PSPHOT.STACK.OUTPUT.IMAGE", i); // File of interest 120 psAssert (output, "missing file?"); 121 122 pmFPAfile *inputRaw = pmFPAfileSelectSingle(config->files, "PSPHOT.STACK.INPUT.RAW", i); // File of interest 123 pmFPAfile *inputCnv = pmFPAfileSelectSingle(config->files, "PSPHOT.STACK.INPUT.CNV", i); // File of interest 124 pmFPAfile *input = inputRaw ? inputRaw : inputCnv; 125 psAssert (input, "missing input file"); 126 127 // find the FPA phu 128 pmHDU *phu = pmFPAviewThisPHU(view, input->fpa); 129 if (!phu) { 130 psWarning("Unable to read bilevel mosaic astrometry for input FPA entry %d", i); 131 continue; 132 } 133 134 char *ctype = psMetadataLookupStr(NULL, phu->header, "CTYPE1"); 135 if (!ctype) { 136 psWarning("Error in WCS keywords for input FPA entry %d", i); 137 continue; 138 } 139 140 if (!foundAstrometry) { 141 bilevelAstrometry = !strcmp (&ctype[4], "-DIS"); 142 } else { 143 if (bilevelAstrometry != !strcmp (&ctype[4], "-DIS")) { 144 psAbort("astrometry type mis-match"); 145 } 146 } 147 148 if (bilevelAstrometry) { 149 // update the output structures 150 if (!pmAstromReadBilevelMosaic(output->fpa, phu->header)) { 151 psWarning("Unable to read bilevel mosaic astrometry for input FPA."); 152 } 153 } 154 } 155 return bilevelAstrometry; 156 } 157 158 bool GetAstrometryChip (pmConfig *config, pmFPAview *view, bool bilevelAstrometry) { 159 160 bool status = false; 161 162 int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM"); 163 psAssert (status, "programming error: must define PSPHOT.INPUT.NUM"); 164 165 // loop over the available readouts 166 for (int i = 0; i < num; i++) { 167 168 // find the currently selected readout 169 pmFPAfile *output = pmFPAfileSelectSingle(config->files, "PSPHOT.STACK.OUTPUT.IMAGE", i); // File of interest 170 psAssert (output, "missing file?"); 171 172 pmFPAfile *inputRaw = pmFPAfileSelectSingle(config->files, "PSPHOT.STACK.INPUT.RAW", i); // File of interest 173 pmFPAfile *inputCnv = pmFPAfileSelectSingle(config->files, "PSPHOT.STACK.INPUT.CNV", i); // File of interest 174 pmFPAfile *input = inputRaw ? inputRaw : inputCnv; 175 psAssert (input, "missing input file"); 176 177 // Need to update the output for astrometry. Read WCS data from the corresponding 178 // header and save in the output fpa 179 pmHDU *inHDU = pmFPAviewThisHDU (view, input->fpa); 180 pmChip *outChip = pmFPAviewThisChip(view, output->fpa); ///< Chip in the output 181 182 # if (UPDATE_HEADER) 183 pmHDU *outHDU = pmFPAviewThisHDU (view, output->fpa); 184 if (!outHDU) { 185 pmFPAAddSourceFromView(output->fpa, "name", view, output->format); 186 outHDU = pmFPAviewThisHDU (view, output->fpa); 187 psAssert (outHDU, "failed to make HDU"); 188 } 189 # endif 190 191 if (bilevelAstrometry) { 192 if (!pmAstromReadBilevelChip (outChip, inHDU->header)) { 193 psWarning("Unable to read bilevel chip astrometry for input FPA."); 194 continue; 195 } 196 # if (UPDATE_HEADER) 197 if (!pmAstromWriteBilevelChip(outHDU->header, outChip, WCS_NONLIN_TOL)) { 198 psWarning("Unable to generate WCS header."); 199 continue; 200 } 201 # endif 202 } else { 203 // we use a default FPA pixel scale of 1.0 204 if (!pmAstromReadWCS (output->fpa, outChip, inHDU->header, 1.0)) { 205 psWarning("Unable to read WCS astrometry for input FPA."); 206 continue; 207 } 208 # if (UPDATE_HEADER) 209 if (UPDATE_HEADER && !pmAstromWriteWCS(outHDU->header, output->fpa, outChip, WCS_NONLIN_TOL)) { 210 psWarning("Unable to generate WCS header."); 211 continue; 212 } 213 # endif 214 } 215 } 216 217 return true; 218 } 219 220 bool SetAstrometryFPA (pmConfig *config, pmFPAview *view, bool bilevelAstrometry) { 221 222 bool status = false; 223 224 if (!bilevelAstrometry) return true; 225 226 int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM"); 227 psAssert (status, "programming error: must define PSPHOT.INPUT.NUM"); 228 229 // loop over the available readouts 230 for (int i = 0; i < num; i++) { 231 232 // find the currently selected readout 233 pmFPAfile *output = pmFPAfileSelectSingle(config->files, "PSPHOT.STACK.OUTPUT.IMAGE", i); // File of interest 234 psAssert (output, "missing file?"); 235 236 # if (UPDATE_HEADER) 237 pmHDU *PHU = pmFPAviewThisPHU(view, output->fpa); 238 if (!PHU) { 239 pmFPAAddSourceFromView(output->fpa, "name", view, output->format); 240 PHU = pmFPAviewThisPHU (view, output->fpa); 241 psAssert (PHU, "failed to make PHU"); 242 } 243 244 if (!pmAstromWriteBilevelMosaic(PHU->header, output->fpa, WCS_NONLIN_TOL)) { 245 psWarning("Unable to generate WCS header."); 246 } 247 # endif 248 } 249 250 return true; 251 } -
branches/czw_branch/20100427/psphot/src/psphotStackMatchPSFs.c
r28017 r28029 1 1 # include "psphotInternal.h" 2 2 3 bool psphotStackMatchPSFs (pmConfig *config, const pmFPAview *view , bool firstPass)3 bool psphotStackMatchPSFs (pmConfig *config, const pmFPAview *view) 4 4 { 5 5 bool status = true; 6 7 psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, "PSPHOT"); 8 psAssert(recipe, "We've thrown an error on this before."); 6 9 7 10 int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM"); 8 11 psAssert (status, "programming error: must define PSPHOT.INPUT.NUM"); 9 12 10 // loop over the available readouts 13 // 'options' carries info needed to perform the stack matching 14 psphotStackOptions *options = psphotStackOptionsAlloc(num); 15 16 options->convolve = psMetadataLookupBool (&status, recipe, "PSPHOT.STACK.MATCH.PSF"); 17 psAssert (status, "PSPHOT.STACK.MATCH.PSF not in recipe"); 18 19 if (options->convolve) { 20 char *convolveSource = psMetadataLookupStr (&status, recipe, "PSPHOT.STACK.MATCH.PSF.SOURCE"); 21 options->convolveSource = psphotStackConvolveSourceFromString (convolveSource); 22 if (options->convolveSource == PSPHOT_CNV_SRC_NONE) { 23 psError (PSPHOT_ERR_CONFIG, true, "stack convolution source not defined in recipe"); 24 return false; 25 } 26 } 27 28 // loop over the available readouts (ignore chisq image) 11 29 for (int i = 0; i < num; i++) { 12 if (!psphotStackMatchPSFsReadout (config, view, i)) { 30 if (!psphotStackMatchPSFsPrepare (config, view, options, i)) { 31 psError (PSPHOT_ERR_CONFIG, false, "failed to set PSF matching options for entry %d", i); 32 return false; 33 } 34 } 35 36 // Generate target PSF 37 if (options->convolve) { 38 options->psf = psphotStackPSF(config, options->numCols, options->numRows, options->psfs, options->inputMask); 39 if (!options->psf) { 40 psError(psErrorCodeLast(), false, "Unable to determine output PSF."); 41 return false; 42 } 43 psMetadataAddPtr(config->arguments, PS_LIST_TAIL, "PSF.TARGET", PS_DATA_UNKNOWN, "Target PSF for stack", options->psf); 44 options->targetSeeing = pmPSFtoFWHM(options->psf, 0.5 * options->numCols, 0.5 * options->numRows); // FWHM for target 45 psLogMsg("psphotStack", PS_LOG_INFO, "Target seeing FWHM: %f\n", options->targetSeeing); 46 47 // XXX is this needed to supply the psf to psphot?? 48 // pmChip *outChip = pmFPAfileThisChip(config->files, view, "PPSTACK.TARGET.PSF"); // Output chip 49 // psMetadataAddPtr(outChip->analysis, PS_LIST_TAIL, "PSPHOT.PSF", PS_DATA_UNKNOWN, "Target PSF", options->psf); 50 // outChip->data_exists = true; 51 } 52 53 // loop over the available readouts (ignore chisq image) 54 for (int i = 0; i < num; i++) { 55 if (!psphotStackMatchPSFsReadout (config, view, options, i)) { 13 56 psError (PSPHOT_ERR_CONFIG, false, "failed to find initial detections for PSPHOT.INPUT entry %d", i); 14 57 return false; 15 58 } 16 59 } 60 61 psFree (options); 17 62 return true; 18 63 } 19 64 20 65 // convolve the image to match desired PSF 21 bool psphotStackMatchPSFsReadout (pmConfig *config, const pmFPAview *view, int index) {66 bool psphotStackMatchPSFsReadout (pmConfig *config, const pmFPAview *view, psphotStackOptions *options, int index) { 22 67 23 bool status; 24 int pass; 25 float NSIGMA_PEAK = 25.0; 26 int NMAX = 0; 68 pmFPAfile *fileSrc = psphotStackGetConvolveSource(config, options, index); 69 if (!fileSrc) { 70 psError(PSPHOT_ERR_CONFIG, false, "desired convolution source is missing"); 71 return false; 72 } 27 73 28 // find the currently selected readout 29 pmFPAfile *fileRaw = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT", index); // File of interest 30 psAssert (fileRaw, "missing file?"); 74 pmFPAfile *fileOut = pmFPAfileSelectSingle(config->files, "PSPHOT.STACK.OUTPUT.IMAGE", index); // File of interest 75 psAssert (fileOut, "missing output file?"); 31 76 32 pm FPAfile *fileCnv = pmFPAfileSelectSingle(config->files, "PSPHOT.INPUT.CONV", index); // File of interest33 psAssert ( fileCnv, "missing file?");77 pmReadout *readoutSrc = pmFPAviewThisReadout(view, fileSrc->fpa); 78 psAssert (readoutSrc, "missing readout?"); 34 79 35 pmReadout *readoutRaw = pmFPAviewThisReadout(view, fileRaw->fpa); 36 psAssert (readoutRaw, "missing readout?"); 80 pmReadout *readoutOut = pmFPAviewThisReadout(view, fileOut->fpa); 81 if (readoutOut == NULL) { 82 readoutOut = pmFPAGenerateReadout(config, view, "PSPHOT.STACK.OUTPUT.IMAGE", fileSrc->fpa, NULL, index); 83 psAssert (readoutOut, "missing readout?"); 84 } 37 85 38 pmReadout *readoutCnv = pmFPAviewThisReadout(view, fileCnv->fpa); 39 psAssert (readoutCnv, "missing readout?"); 40 41 // user-defined masks to test for good/bad pixels (build from recipe list if not yet set) 42 psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels 43 psAssert (maskVal, "missing mask value?"); 44 45 /***** set up recipe options *****/ 46 47 psMetadata *stackRecipe = psMetadataLookupMetadata(NULL, config->recipes, PPSTACK_RECIPE); // ppStack recipe 48 psAssert(stackRecipe, "We've thrown an error on this before."); 49 50 // Look up appropriate values from the ppSub recipe 51 psMetadata *subRecipe = psMetadataLookupMetadata(NULL, config->recipes, "PPSUB"); // PPSUB recipe 52 psAssert(subRecipe, "recipe missing"); 53 54 int size = psMetadataLookupS32(NULL, subRecipe, "KERNEL.SIZE"); // Kernel half-size 55 56 float deconvLimit = psMetadataLookupF32(NULL, stackRecipe, "DECONV.LIMIT"); // Limit on deconvolution fraction 57 58 psString maskValStr = psMetadataLookupStr(NULL, subRecipe, "MASK.VAL"); // Name of bits to mask going in 59 psImageMaskType maskVal = pmConfigMaskGet(maskValStr, config); // Bits to mask going in to pmSubtractionMatch 60 psString maskPoorStr = psMetadataLookupStr(NULL, stackRecipe, "MASK.POOR"); // Name of bits to mask for poor 61 psImageMaskType maskPoor = pmConfigMaskGet(maskPoorStr, config); // Bits to mask for poor pixels 62 psString maskBadStr = psMetadataLookupStr(NULL, stackRecipe, "MASK.BAD"); // Name of bits to mask for bad 63 psImageMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels 64 65 bool mdok; // Status of MD lookup 66 float penalty = psMetadataLookupF32(NULL, subRecipe, "PENALTY"); // Penalty for wideness 67 int threads = psMetadataLookupS32(NULL, config->arguments, "-threads"); // Number of threads 68 69 if (!pmReadoutMaskNonfinite(readout, maskVal)) { 86 // set NAN pixels to 'SAT' 87 psImageMaskType maskVal = pmConfigMaskGet("SAT", config); 88 if (!pmReadoutMaskNonfinite(readoutSrc, maskVal)) { 70 89 psError(psErrorCodeLast(), false, "Unable to mask non-finite pixels in readout."); 71 90 return false; … … 74 93 // Image Matching (PSFs or just flux) 75 94 if (options->convolve) { 76 // Full match of PSFs 77 pmReadout *conv = pmReadoutAlloc(NULL); // Conv readout, for holding results temporarily 78 79 // Read previously produced kernel 80 if (psMetadataLookupBool(NULL, config->arguments, "PPSTACK.DEBUG.STACK")) { 81 loadKernel(); 82 } else { 83 matchKernel(); 84 } // !DEBUG.STACK 85 86 saveMatchData(); 87 88 saveChiSquare(); 89 90 renormKernel(); 91 92 // Reject image completely if the maximum deconvolution fraction exceeds the limit 93 float deconv = psMetadataLookupF32(NULL, conv->analysis, PM_SUBTRACTION_ANALYSIS_DECONV_MAX); // Max deconvolution fraction 94 if (deconv > deconvLimit) { 95 psWarning("Maximum deconvolution fraction (%f) exceeds limit (%f) --- rejecting image %d\n", deconv, deconvLimit, index); 96 psFree(conv); 97 return NULL; 98 } 99 100 readout->analysis = psMetadataCopy(readout->analysis, conv->analysis); 101 102 psFree(conv); 95 matchKernel(config, readoutOut, readoutSrc, options, index); 96 saveMatchData(readoutOut, options, index); 97 // renormKernel(readoutCnv, options, index); 103 98 } else { 104 // only match the flux 105 float norm = powf(10.0, -0.4 * options->norm->data.F32[index]); // Normalisation 106 psBinaryOp(readout->image, readout->image, "*", psScalarAlloc(norm, PS_TYPE_F32)); 107 psBinaryOp(readout->variance, readout->variance, "*", psScalarAlloc(PS_SQR(norm), PS_TYPE_F32)); 99 // only match the flux (NO! not for multi-filter, at least!) 100 // XXX do not generate readoutCnv in this case? 101 // float norm = powf(10.0, -0.4 * options->norm->data.F32[index]); // Normalisation 102 // psBinaryOp(readoutRaw->image, readoutRaw->image, "*", psScalarAlloc(norm, PS_TYPE_F32)); 103 // psBinaryOp(readoutRaw->variance, readoutRaw->variance, "*", psScalarAlloc(PS_SQR(norm), PS_TYPE_F32)); 108 104 } 109 105 110 // Ensure the background value is zero 111 psStats *bg = psStatsAlloc(PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV); // Statistics for background 112 psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS); // Random number generator 113 if (!psImageBackground(bg, NULL, readout->image, readout->mask, maskVal | maskBad, rng)) { 114 psWarning("Can't measure background for image."); 115 psErrorClear(); 116 } else { 117 if (!psMetadataLookupBool(NULL, config->arguments, "PPSTACK.SKIP.BG.SUB")) { 118 psLogMsg("ppStack", PS_LOG_INFO, "Correcting convolved image background by %lf (+/- %lf)", 119 psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN), psStatsGetValue(bg, PS_STAT_ROBUST_STDEV)); 120 (void)psBinaryOp(readout->image, readout->image, "-", 121 psScalarAlloc(psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN), PS_TYPE_F32)); 122 } 123 } 106 rescaleData(readoutOut, config, options, index); 124 107 125 if (!stackRenormaliseReadout(config, readout)) { 126 psFree(rng); 127 psFree(bg); 128 return false; 129 } 130 131 // Measure the variance level for the weighting 132 if (psMetadataLookupBool(NULL, stackRecipe, "WEIGHTS")) { 133 if (!psImageBackground(bg, NULL, readout->variance, readout->mask, maskVal | maskBad, rng)) { 134 psError(PPSTACK_ERR_DATA, false, "Can't measure mean variance for image."); 135 psFree(rng); 136 psFree(bg); 137 return false; 138 } 139 options->weightings->data.F32[index] = 1.0 / (psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN) * psImageCovarianceFactor(readout->covariance)); 140 } else { 141 options->weightings->data.F32[index] = 1.0; 142 } 143 psLogMsg("ppStack", PS_LOG_INFO, "Weighting for image %d is %f\n", 144 index, options->weightings->data.F32[index]); 145 146 psFree(rng); 147 psFree(bg); 148 149 dumpImage3(); 108 dumpImage(readoutOut, readoutSrc, index, "convolved"); 150 109 151 110 return true; 152 111 } 112 113 114 # if (0) 115 // Read previously produced kernel 116 if (psMetadataLookupBool(NULL, config->arguments, "PPSTACK.DEBUG.STACK")) { 117 loadKernel(config, readoutCnv, options, index); 118 } else { 119 matchKernel(config, readoutCnv, readoutRaw, options, index); 120 } 121 # endif -
branches/czw_branch/20100427/psphot/src/psphotStackMatchPSFsUtils.c
r28017 r28029 1 /***** defines *****/ 2 3 #define ARRAY_BUFFER 16 // Number to add to array at a time 4 #define MAG_IGNORE 50 // Ignore magnitudes fainter than this --- they're not real! 5 #define FAKE_SIZE 1 // Size of fake convolution kernel 6 #define SOURCE_MASK (PM_SOURCE_MODE_FAIL | PM_SOURCE_MODE_DEFECT | PM_SOURCE_MODE_SATURATED | PM_SOURCE_MODE_CR_LIMIT | PM_SOURCE_MODE_EXT_LIMIT) // Mask to apply to input sources 7 #define NOISE_FRACTION 0.01 // Set minimum flux to this fraction of noise 8 #define COVAR_FRAC 0.01 // Truncation fraction for covariance matrix 1 # include "psphotInternal.h" 2 # define ARRAY_BUFFER 16 // Number to add to array at a time 9 3 10 4 // XXX better name … … 18 12 psFree(resolved); 19 13 if (!fits) { 20 psError(P PSTACK_ERR_IO, false, "Unable to open previously produced image: %s", name);14 psError(PSPHOT_ERR_IO, false, "Unable to open previously produced image: %s", name); 21 15 return false; 22 16 } 23 17 psImage *image = psFitsReadImage(fits, psRegionSet(0,0,0,0), 0); // Image of interest 24 18 if (!image) { 25 psError(P PSTACK_ERR_IO, false, "Unable to read previously produced image: %s", name);19 psError(PSPHOT_ERR_IO, false, "Unable to read previously produced image: %s", name); 26 20 psFitsClose(fits); 27 21 return false; … … 51 45 } 52 46 47 # define SN_MIN 50.0 53 48 psArray *stackSourcesFilter(psArray *sources, // Source list to filter 54 int exclusion // Exclusion zone, pixels49 int exclusion // Exclusion zone, pixels 55 50 ) 56 51 { … … 68 63 continue; 69 64 } 65 if (!source->peak) continue; 66 if (source->peak->SN < SN_MIN) continue; 70 67 coordsFromSource(&x->data.F32[numGood], &y->data.F32[numGood], source); 71 68 numGood++; … … 83 80 continue; 84 81 } 82 if (!source->peak) continue; 83 if (source->peak->SN < SN_MIN) continue; 85 84 float xSource, ySource; // Coordinates of source 86 85 coordsFromSource(&xSource, &ySource, source); … … 90 89 91 90 long numWithin = psTreeWithin(tree, coords, exclusion); // Number within exclusion zone 92 psTrace("p pStack", 9, "Source at %.0lf,%.0lf has %ld sources in exclusion zone",91 psTrace("psphotStack", 9, "Source at %.0lf,%.0lf has %ld sources in exclusion zone", 93 92 coords->data.F64[0], coords->data.F64[1], numWithin); 94 93 if (numWithin == 1) { … … 104 103 psFree(y); 105 104 106 psLogMsg("p pStack", PS_LOG_INFO, "Filtered out %d of %d sources", numFiltered, numGood);105 psLogMsg("psphotStack", PS_LOG_INFO, "Filtered out %d of %d sources", numFiltered, numGood); 107 106 108 107 return filtered; … … 111 110 // Add background into the fake image 112 111 // Based on ppSubBackground() 113 staticpsImage *stackBackgroundModel(pmReadout *ro, // Readout for which to generate background model112 psImage *stackBackgroundModel(pmReadout *ro, // Readout for which to generate background model 114 113 const pmConfig *config // Configuration 115 114 ) … … 121 120 int numCols = image->numCols, numRows = image->numRows; // Size of image 122 121 123 psMetadata *ppStackRecipe = psMetadataLookupPtr(NULL, config->recipes, PPSTACK_RECIPE);122 psMetadata *ppStackRecipe = psMetadataLookupPtr(NULL, config->recipes, "PPSTACK"); 124 123 psAssert(ppStackRecipe, "Need PPSTACK recipe"); 125 psMetadata *psphotRecipe = psMetadataLookupPtr(NULL, config->recipes, PSPHOT_RECIPE);124 psMetadata *psphotRecipe = psMetadataLookupPtr(NULL, config->recipes, "PSPHOT"); 126 125 psAssert(psphotRecipe, "Need PSPHOT recipe"); 127 126 … … 138 137 psImage *unbinned = psImageAlloc(numCols, numRows, PS_TYPE_F32); // Unbinned background model 139 138 if (!psImageUnbin(unbinned, binned, binning)) { 140 psError(P PSTACK_ERR_DATA, false, "Unable to unbin background model");139 psError(PSPHOT_ERR_DATA, false, "Unable to unbin background model"); 141 140 psFree(binned); 142 141 psFree(unbinned); … … 153 152 ) 154 153 { 155 #if 1156 154 bool mdok; // Status of metadata lookups 157 155 158 psMetadata *recipe = psMetadataLookupPtr(NULL, config->recipes, PPSTACK_RECIPE); // Recipe for ppStack156 psMetadata *recipe = psMetadataLookupPtr(NULL, config->recipes, "PPSTACK"); // Recipe for ppStack 159 157 psAssert(recipe, "Need PPSTACK recipe"); 160 158 … … 163 161 int num = psMetadataLookupS32(&mdok, recipe, "RENORM.NUM"); 164 162 if (!mdok) { 165 psError(P PSTACK_ERR_CONFIG, true, "RENORM.NUM is not set in the recipe");163 psError(PSPHOT_ERR_CONFIG, true, "RENORM.NUM is not set in the recipe"); 166 164 return false; 167 165 } 168 166 float minValid = psMetadataLookupF32(&mdok, recipe, "RENORM.MIN"); 169 167 if (!mdok) { 170 psError(P PSTACK_ERR_CONFIG, true, "RENORM.MIN is not set in the recipe");168 psError(PSPHOT_ERR_CONFIG, true, "RENORM.MIN is not set in the recipe"); 171 169 return false; 172 170 } 173 171 float maxValid = psMetadataLookupF32(&mdok, recipe, "RENORM.MAX"); 174 172 if (!mdok) { 175 psError(P PSTACK_ERR_CONFIG, true, "RENORM.MAX is not set in the recipe");173 psError(PSPHOT_ERR_CONFIG, true, "RENORM.MAX is not set in the recipe"); 176 174 return false; 177 175 } … … 181 179 psImageCovarianceTransfer(readout->variance, readout->covariance); 182 180 return pmReadoutVarianceRenormalise(readout, maskBad, num, minValid, maxValid); 183 #else184 return true;185 #endif186 181 } 187 182 … … 190 185 // It implicitly assumes the output root name is the same between invocations. 191 186 192 bool loadKernel () { 193 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.CONV.KERNEL", index); 194 psAssert(file, "Require file"); 195 196 pmFPAview *view = pmFPAviewAlloc(0); // View to readout of interest 197 view->chip = view->cell = view->readout = 0; 198 psString filename = pmFPAfileNameFromRule(file->filerule, file, view); // Filename of interest 199 200 // Read convolution kernel 201 psString resolved = pmConfigConvertFilename(filename, config, false, false); // Resolved filename 202 psFree(filename); 203 psFits *fits = psFitsOpen(resolved, "r"); // FITS file for subtraction kernel 204 psFree(resolved); 205 if (!fits || !pmReadoutReadSubtractionKernels(conv, fits)) { 206 psError(PPSTACK_ERR_IO, false, "Unable to read previously produced kernel"); 207 psFitsClose(fits); 208 return false; 209 } 210 psFitsClose(fits); 211 212 if (!readImage(&readout->image, options->convImages->data[index], config) || 213 !readImage(&readout->mask, options->convMasks->data[index], config) || 214 !readImage(&readout->variance, options->convVariances->data[index], config)) { 215 psError(PPSTACK_ERR_IO, false, "Unable to read previously produced image."); 216 return false; 217 } 218 219 psRegion *region = psMetadataLookupPtr(NULL, conv->analysis, 220 PM_SUBTRACTION_ANALYSIS_REGION); // Convolution region 221 pmSubtractionKernels *kernels = psMetadataLookupPtr(NULL, conv->analysis, 222 PM_SUBTRACTION_ANALYSIS_KERNEL); 223 224 pmSubtractionAnalysis(conv->analysis, NULL, kernels, region, 225 readout->image->numCols, readout->image->numRows); 226 227 psKernel *kernel = pmSubtractionKernel(kernels, 0.0, 0.0, false); // Convolution kernel 228 bool oldThreads = psImageCovarianceSetThreads(true); // Old thread setting 229 psKernel *covar = psImageCovarianceCalculate(kernel, readout->covariance); // Covariance matrix 230 psImageCovarianceSetThreads(oldThreads); 231 psFree(readout->covariance); 232 readout->covariance = covar; 233 psFree(kernel); 234 } 235 236 bool dumpImage() { 237 // XXX should be optional 238 { 239 pmHDU *hdu = pmHDUFromCell(readout->parent); 240 psString name = NULL; 241 psStringAppend(&name, "fake_%03d.fits", index); 242 pmStackVisualPlotTestImage(fake->image, name); 243 psFits *fits = psFitsOpen(name, "w"); 244 psFree(name); 245 psFitsWriteImage(fits, hdu->header, fake->image, 0, NULL); 246 psFitsClose(fits); 247 } 248 { 249 pmHDU *hdu = pmHDUFromCell(readout->parent); 250 psString name = NULL; 251 psStringAppend(&name, "real_%03d.fits", index); 252 pmStackVisualPlotTestImage(readout->image, name); 253 psFits *fits = psFitsOpen(name, "w"); 254 psFree(name); 255 psFitsWriteImage(fits, hdu->header, readout->image, 0, NULL); 256 psFitsClose(fits); 257 } 258 } 259 260 bool dumpImage2() { 261 // XXX should be optional 262 263 { 264 pmHDU *hdu = pmHDUFromCell(readout->parent); 265 psString name = NULL; 266 psStringAppend(&name, "conv_%03d.fits", index); 267 pmStackVisualPlotTestImage(conv->image, name); 268 psFits *fits = psFitsOpen(name, "w"); 269 psFree(name); 270 psFitsWriteImage(fits, hdu->header, conv->image, 0, NULL); 271 psFitsClose(fits); 272 } 273 { 274 pmHDU *hdu = pmHDUFromCell(readout->parent); 275 psString name = NULL; 276 psStringAppend(&name, "diff_%03d.fits", index); 277 pmStackVisualPlotTestImage(fake->image, name); 278 psFits *fits = psFitsOpen(name, "w"); 279 psFree(name); 280 psBinaryOp(fake->image, conv->image, "-", fake->image); 281 psFitsWriteImage(fits, hdu->header, fake->image, 0, NULL); 282 psFitsClose(fits); 283 } 284 } 285 286 bool dumpImage3() 187 # if (0) 188 bool loadKernel (pmConfig *config, pmReadout *readoutCnv, psphotStackOptions *options, int index) { 189 190 // Read the convolution kernel from the saved file 191 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.CONV.KERNEL", index); 192 psAssert(file, "Require file"); 193 194 pmFPAview *view = pmFPAviewAlloc(0); // View to readout of interest 195 view->chip = view->cell = view->readout = 0; 196 psString filename = pmFPAfileNameFromRule(file->filerule, file, view); // Filename of interest 197 198 // Read convolution kernel data 199 psString resolved = pmConfigConvertFilename(filename, config, false, false); // Resolved filename 200 psFree(filename); 201 psFits *fits = psFitsOpen(resolved, "r"); // FITS file for subtraction kernel 202 psFree(resolved); 203 if (!fits || !pmReadoutReadSubtractionKernels(readoutCnv, fits)) { 204 psError(PSPHOT_ERR_IO, false, "Unable to read previously produced kernel"); 205 psFitsClose(fits); 206 return false; 207 } 208 psFitsClose(fits); 209 210 // read the convolved pixels (image, mask, variance) -- names are pre-defined 211 if (!readImage(&readoutCnv->image, options->convImages->data[index], config) || 212 !readImage(&readoutCnv->mask, options->convMasks->data[index], config) || 213 !readImage(&readoutCnv->variance, options->convVariances->data[index], config)) { 214 psError(PSPHOT_ERR_IO, false, "Unable to read previously produced image."); 215 return false; 216 } 217 218 // XXX ??? not sure what is happening here -- consult Paul 219 psRegion *region = psMetadataLookupPtr(NULL, readoutCnv->analysis, PM_SUBTRACTION_ANALYSIS_REGION); // Convolution region 220 pmSubtractionKernels *kernels = psMetadataLookupPtr(NULL, readoutCnv->analysis, PM_SUBTRACTION_ANALYSIS_KERNEL); 221 222 pmSubtractionAnalysis(readoutCnv->analysis, NULL, kernels, region, readoutCnv->image->numCols, readoutCnv->image->numRows); 223 224 psKernel *kernel = pmSubtractionKernel(kernels, 0.0, 0.0, false); // Convolution kernel 225 226 // update the covariance matrix 227 // XXX why is this needed if we have correctly read the saved data? 228 bool oldThreads = psImageCovarianceSetThreads(true); // Old thread setting 229 psKernel *covar = psImageCovarianceCalculate(kernel, readoutCnv->covariance); // Covariance matrix 230 psImageCovarianceSetThreads(oldThreads); 231 psFree(readoutCnv->covariance); 232 readoutCnv->covariance = covar; 233 psFree(kernel); 234 return true; 235 } 236 # endif 237 238 bool dumpImage(pmReadout *readoutOut, pmReadout *readoutRef, int index, char *rootname) { 239 240 pmHDU *hdu = pmHDUFromCell(readoutRef->parent); 241 psString name = NULL; 242 psStringAppend(&name, "%s_%03d.fits", rootname, index); 243 pmStackVisualPlotTestImage(readoutOut->image, name); 244 psFits *fits = psFitsOpen(name, "w"); 245 psFree(name); 246 psFitsWriteImage(fits, hdu->header, readoutOut->image, 0, NULL); 247 psFitsClose(fits); 248 return true; 249 } 250 251 bool dumpImageDiff(pmReadout *readoutConv, pmReadout *readoutFake, pmReadout *readoutRef, int index, char *rootname) { 252 253 pmHDU *hdu = pmHDUFromCell(readoutRef->parent); 254 psString name = NULL; 255 psStringAppend(&name, "%s_%03d.fits", rootname, index); 256 pmStackVisualPlotTestImage(readoutFake->image, name); 257 psFits *fits = psFitsOpen(name, "w"); 258 psFree(name); 259 psBinaryOp(readoutFake->image, readoutConv->image, "-", readoutFake->image); 260 psFitsWriteImage(fits, hdu->header, readoutFake->image, 0, NULL); 261 psFitsClose(fits); 262 return true; 263 } 264 265 // perform the bulk of the PSF-matching 266 bool matchKernel(pmConfig *config, pmReadout *readoutOut, pmReadout *readoutSrc, psphotStackOptions *options, int index) { 267 268 bool mdok; 269 270 psAssert(options->psf, "Require target PSF"); 271 psAssert(options->sourceLists && options->sourceLists->data[index], "Require source list"); 272 273 psMetadata *stackRecipe = psMetadataLookupMetadata(NULL, config->recipes, "PPSTACK"); // ppStack recipe 274 psAssert(stackRecipe, "We've thrown an error on this before."); 275 276 // Look up appropriate values from the ppSub recipe 277 psMetadata *subRecipe = psMetadataLookupMetadata(NULL, config->recipes, "PPSUB"); // PPSUB recipe 278 psAssert(subRecipe, "recipe missing"); 279 280 psString maskValStr = psMetadataLookupStr(NULL, subRecipe, "MASK.VAL"); // Name of bits to mask going in 281 psString maskPoorStr = psMetadataLookupStr(NULL, stackRecipe, "MASK.POOR"); // Name of bits to mask for poor 282 psString maskBadStr = psMetadataLookupStr(NULL, stackRecipe, "MASK.BAD"); // Name of bits to mask for bad 283 284 psImageMaskType maskVal = pmConfigMaskGet(maskValStr, config); // Bits to mask going in to pmSubtractionMatch 285 psImageMaskType maskPoor = pmConfigMaskGet(maskPoorStr, config); // Bits to mask for poor pixels 286 psImageMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels 287 288 float penalty = psMetadataLookupF32(NULL, subRecipe, "PENALTY"); // Penalty for wideness 289 int threads = psMetadataLookupS32(NULL, config->arguments, "-threads"); // Number of threads 290 291 int order = psMetadataLookupS32(NULL, subRecipe, "SPATIAL.ORDER"); // Spatial polynomial order 292 float regionSize = psMetadataLookupF32(NULL, subRecipe, "REGION.SIZE"); // Size of iso-kernel regs 293 float spacing = psMetadataLookupF32(NULL, subRecipe, "STAMP.SPACING"); // Typical stamp spacing 294 295 int footprint = psMetadataLookupS32(NULL, subRecipe, "STAMP.FOOTPRINT"); // Stamp half-size 296 int size = psMetadataLookupS32(NULL, subRecipe, "KERNEL.SIZE"); // Kernel half-size 297 298 float threshold = psMetadataLookupF32(NULL, subRecipe, "STAMP.THRESHOLD"); // Threshold for stmps 299 int stride = psMetadataLookupS32(NULL, subRecipe, "STRIDE"); // Size of convolution patches 300 int iter = psMetadataLookupS32(NULL, subRecipe, "ITER"); // Rejection iterations 301 float rej = psMetadataLookupF32(NULL, subRecipe, "REJ"); // Rejection threshold 302 float kernelError = psMetadataLookupF32(NULL, subRecipe, "KERNEL.ERR"); // Relative systematic error in kernel 303 float normFrac = psMetadataLookupF32(NULL, subRecipe, "NORM.FRAC"); // Fraction of window for normalisn windw 304 float sysError = psMetadataLookupF32(NULL, subRecipe, "SYS.ERR"); // Relative systematic error in images 305 float skyErr = psMetadataLookupF32(NULL, subRecipe, "SKY.ERR"); // Additional error in sky 306 float covarFrac = psMetadataLookupF32(NULL, subRecipe, "COVAR.FRAC"); // Fraction for covariance calculation 307 308 const char *typeStr = psMetadataLookupStr(NULL, subRecipe, "KERNEL.TYPE"); // Kernel type 309 pmSubtractionKernelsType type = pmSubtractionKernelsTypeFromString(typeStr); // Kernel type 310 psVector *widths = psMetadataLookupPtr(NULL, subRecipe, "ISIS.WIDTHS"); // ISIS Gaussian widths 311 psVector *orders = psMetadataLookupPtr(NULL, subRecipe, "ISIS.ORDERS"); // ISIS Polynomial orders 312 int inner = psMetadataLookupS32(NULL, subRecipe, "INNER"); // Inner radius 313 int ringsOrder = psMetadataLookupS32(NULL, subRecipe, "RINGS.ORDER"); // RINGS polynomial order 314 int binning = psMetadataLookupS32(NULL, subRecipe, "SPAM.BINNING"); // Binning for SPAM kernel 315 float badFrac = psMetadataLookupF32(NULL, subRecipe, "BADFRAC"); // Maximum bad fraction 316 float optThresh = psMetadataLookupF32(&mdok, subRecipe, "OPTIMUM.TOL"); // Tolerance for search 317 int optOrder = psMetadataLookupS32(&mdok, subRecipe, "OPTIMUM.ORDER"); // Order for search 318 float poorFrac = psMetadataLookupF32(&mdok, subRecipe, "POOR.FRACTION"); // Fraction for "poor" 319 320 bool scale = psMetadataLookupBool(NULL, subRecipe, "SCALE"); // Scale kernel parameters? 321 float scaleRef = psMetadataLookupF32(NULL, subRecipe, "SCALE.REF"); // Reference for scaling 322 float scaleMin = psMetadataLookupF32(NULL, subRecipe, "SCALE.MIN"); // Minimum for scaling 323 float scaleMax = psMetadataLookupF32(NULL, subRecipe, "SCALE.MAX"); // Maximum for scaling 324 if (!isfinite(scaleRef) || !isfinite(scaleMin) || !isfinite(scaleMax)) { 325 psError(PSPHOT_ERR_CONFIG, false, 326 "Scale parameters (SCALE.REF=%f, SCALE.MIN=%f, SCALE.MAX=%f) not set in PPSUB recipe.", 327 scaleRef, scaleMin, scaleMax); 328 return false; 329 } 330 331 // These values are specified specifically for stacking 332 const char *stampsName = psMetadataLookupStr(&mdok, config->arguments, "STAMPS");// Stamps filename 333 334 psVector *widthsCopy = NULL; 335 psVector *optWidths = NULL; 336 pmReadout *fake = NULL; 337 psArray *stampSources = NULL; 338 339 bool optimum = false; 340 optWidths = SetOptWidths(&optimum, subRecipe); // Vector with FWHMs for optimum search 341 342 // For the sake of stamps, remove nearby sources 343 stampSources = stackSourcesFilter(options->sourceLists->data[index], footprint); // Filtered list of sources 344 345 fake = makeFakeReadout(config, readoutSrc, stampSources, options->psf, maskVal | maskBad, footprint + size); 346 if (!fake) goto escape; 347 348 dumpImage(fake, readoutSrc, index, "fake"); 349 dumpImage(readoutSrc, readoutSrc, index, "real"); 350 351 if (threads) pmSubtractionThreadsInit(); 352 353 // Do the image matching 354 pmSubtractionKernels *kernel = psMetadataLookupPtr(&mdok, readoutSrc->analysis, PM_SUBTRACTION_ANALYSIS_KERNEL); // Conv kernel 355 if (kernel) { 356 if (!pmSubtractionMatchPrecalc(NULL, readoutOut, fake, readoutSrc, readoutSrc->analysis, stride, kernelError, covarFrac, maskVal, maskBad, maskPoor, poorFrac, badFrac)) { 357 psError(psErrorCodeLast(), false, "Unable to convolve images."); 358 goto escape; 359 } 360 } else { 361 // Scale the input parameters 362 widthsCopy = psVectorCopy(NULL, widths, PS_TYPE_F32); // Copy of kernel widths 363 if (scale && !pmSubtractionParamsScale(&size, &footprint, widthsCopy, options->inputSeeing->data.F32[index], options->targetSeeing, scaleRef, scaleMin, scaleMax)) { 364 psError(psErrorCodeLast(), false, "Unable to scale kernel parameters"); 365 goto escape; 366 } 367 368 if (!pmSubtractionMatch(NULL, readoutOut, fake, readoutSrc, footprint, stride, regionSize, spacing, threshold, stampSources, stampsName, type, size, order, widthsCopy, orders, inner, ringsOrder, binning, penalty, optimum, optWidths, optOrder, optThresh, iter, rej, normFrac, sysError, skyErr, kernelError, covarFrac, maskVal, maskBad, maskPoor, poorFrac, badFrac, PM_SUBTRACTION_MODE_2)) { 369 psError(psErrorCodeLast(), false, "Unable to match images."); 370 goto escape; 371 } 372 } 373 374 // Reject image completely if the maximum deconvolution fraction exceeds the limit 375 float deconvLimit = psMetadataLookupF32(NULL, stackRecipe, "DECONV.LIMIT"); // Limit on deconvolution fraction 376 float deconv = psMetadataLookupF32(NULL, readoutOut->analysis, PM_SUBTRACTION_ANALYSIS_DECONV_MAX); // Max deconvolution fraction 377 if (deconv > deconvLimit) { 378 psWarning("Maximum deconvolution fraction (%f) exceeds limit (%f) --- rejecting image %d\n", deconv, deconvLimit, index); 379 goto escape; 380 } 381 382 dumpImage(readoutOut, readoutSrc, index, "conv"); 383 dumpImageDiff(readoutOut, fake, readoutSrc, index, "diff"); 384 385 psFree(fake); 386 psFree(optWidths); 387 psFree(stampSources); 388 psFree(widthsCopy); 389 pmSubtractionThreadsFinalize(); 390 return true; 391 392 escape: 393 psFree(fake); 394 psFree(optWidths); 395 psFree(stampSources); 396 psFree(widthsCopy); 397 pmSubtractionThreadsFinalize(); 398 return false; 399 } 400 401 // Extract the regions and solutions used in the image matching 402 // This stops them from being freed when we iterate back up the FPA 403 // Record the chi-square value 404 // XXX this function may not be needed for psphotStack 405 bool saveMatchData (pmReadout *readout, psphotStackOptions *options, int index) { 406 407 psArray *regions = options->regions->data[index] = psArrayAllocEmpty(ARRAY_BUFFER); // Match regions 287 408 { 288 pmHDU *hdu = pmHDUFromCell(readout->parent); 289 psString name = NULL; 290 psStringAppend(&name, "convolved_%03d.fits", index); 291 pmStackVisualPlotTestImage(readout->image, name); 292 psFits *fits = psFitsOpen(name, "w"); 293 psFree(name); 294 psFitsWriteImage(fits, hdu->header, readout->image, 0, NULL); 295 psFitsClose(fits); 296 } 297 298 bool matchKernel() { 299 // Normal operations here 300 psAssert(options->psf, "Require target PSF"); 301 psAssert(options->sourceLists && options->sourceLists->data[index], "Require source list"); 302 303 int order = psMetadataLookupS32(NULL, subRecipe, "SPATIAL.ORDER"); // Spatial polynomial order 304 float regionSize = psMetadataLookupF32(NULL, subRecipe, "REGION.SIZE"); // Size of iso-kernel regs 305 float spacing = psMetadataLookupF32(NULL, subRecipe, "STAMP.SPACING"); // Typical stamp spacing 306 int footprint = psMetadataLookupS32(NULL, subRecipe, "STAMP.FOOTPRINT"); // Stamp half-size 307 float threshold = psMetadataLookupF32(NULL, subRecipe, "STAMP.THRESHOLD"); // Threshold for stmps 308 int stride = psMetadataLookupS32(NULL, subRecipe, "STRIDE"); // Size of convolution patches 309 int iter = psMetadataLookupS32(NULL, subRecipe, "ITER"); // Rejection iterations 310 float rej = psMetadataLookupF32(NULL, subRecipe, "REJ"); // Rejection threshold 311 float kernelError = psMetadataLookupF32(NULL, subRecipe, "KERNEL.ERR"); // Relative systematic error in kernel 312 float normFrac = psMetadataLookupF32(NULL, subRecipe, "NORM.FRAC"); // Fraction of window for normalisn windw 313 float sysError = psMetadataLookupF32(NULL, subRecipe, "SYS.ERR"); // Relative systematic error in images 314 float skyErr = psMetadataLookupF32(NULL, subRecipe, "SKY.ERR"); // Additional error in sky 315 float covarFrac = psMetadataLookupF32(NULL, subRecipe, "COVAR.FRAC"); // Fraction for covariance calculation 316 317 const char *typeStr = psMetadataLookupStr(NULL, subRecipe, "KERNEL.TYPE"); // Kernel type 318 pmSubtractionKernelsType type = pmSubtractionKernelsTypeFromString(typeStr); // Kernel type 319 psVector *widths = psMetadataLookupPtr(NULL, subRecipe, "ISIS.WIDTHS"); // ISIS Gaussian widths 320 psVector *orders = psMetadataLookupPtr(NULL, subRecipe, "ISIS.ORDERS"); // ISIS Polynomial orders 321 int inner = psMetadataLookupS32(NULL, subRecipe, "INNER"); // Inner radius 322 int ringsOrder = psMetadataLookupS32(NULL, subRecipe, "RINGS.ORDER"); // RINGS polynomial order 323 int binning = psMetadataLookupS32(NULL, subRecipe, "SPAM.BINNING"); // Binning for SPAM kernel 324 float badFrac = psMetadataLookupF32(NULL, subRecipe, "BADFRAC"); // Maximum bad fraction 325 bool optimum = psMetadataLookupBool(&mdok, subRecipe, "OPTIMUM"); // Derive optimum parameters? 326 float optMin = psMetadataLookupF32(&mdok, subRecipe, "OPTIMUM.MIN"); // Minimum width for search 327 float optMax = psMetadataLookupF32(&mdok, subRecipe, "OPTIMUM.MAX"); // Maximum width for search 328 float optStep = psMetadataLookupF32(&mdok, subRecipe, "OPTIMUM.STEP"); // Step for search 329 float optThresh = psMetadataLookupF32(&mdok, subRecipe, "OPTIMUM.TOL"); // Tolerance for search 330 int optOrder = psMetadataLookupS32(&mdok, subRecipe, "OPTIMUM.ORDER"); // Order for search 331 float poorFrac = psMetadataLookupF32(&mdok, subRecipe, "POOR.FRACTION"); // Fraction for "poor" 332 333 bool scale = psMetadataLookupBool(NULL, subRecipe, "SCALE"); // Scale kernel parameters? 334 float scaleRef = psMetadataLookupF32(NULL, subRecipe, "SCALE.REF"); // Reference for scaling 335 float scaleMin = psMetadataLookupF32(NULL, subRecipe, "SCALE.MIN"); // Minimum for scaling 336 float scaleMax = psMetadataLookupF32(NULL, subRecipe, "SCALE.MAX"); // Maximum for scaling 337 if (!isfinite(scaleRef) || !isfinite(scaleMin) || !isfinite(scaleMax)) { 338 psError(PPSTACK_ERR_CONFIG, false, 339 "Scale parameters (SCALE.REF=%f, SCALE.MIN=%f, SCALE.MAX=%f) not set in PPSUB recipe.", 340 scaleRef, scaleMin, scaleMax); 341 return false; 342 } 343 344 345 // These values are specified specifically for stacking 346 const char *stampsName = psMetadataLookupStr(NULL, config->arguments, "STAMPS");// Stamps filename 347 348 psVector *optWidths = NULL; // Vector with FWHMs for optimum search 349 if (optimum) { 350 optWidths = psVectorCreate(optWidths, optMin, optMax, optStep, PS_TYPE_F32); 351 } 352 353 pmReadout *fake = pmReadoutAlloc(NULL); // Fake readout with target PSF 354 355 psStats *bg = psStatsAlloc(PS_STAT_ROBUST_STDEV); // Statistics for background 356 psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS); // Random number generator 357 if (!psImageBackground(bg, NULL, readout->image, readout->mask, maskVal | maskBad, rng)) { 358 psError(PPSTACK_ERR_DATA, false, "Can't measure background for image."); 359 psFree(fake); 360 psFree(optWidths); 361 psFree(conv); 362 psFree(bg); 363 psFree(rng); 364 return false; 365 } 366 float minFlux = NOISE_FRACTION * bg->robustStdev; // Minimum flux level for fake image 409 psString regex = NULL; // Regular expression 410 psStringAppend(®ex, "^%s$", PM_SUBTRACTION_ANALYSIS_REGION); 411 psMetadataIterator *iter = psMetadataIteratorAlloc(readout->analysis, PS_LIST_HEAD, regex); 412 psFree(regex); 413 psMetadataItem *item = NULL;// Item from iteration 414 while ((item = psMetadataGetAndIncrement(iter))) { 415 assert(item->type == PS_DATA_REGION); 416 regions = psArrayAdd(regions, ARRAY_BUFFER, item->data.V); 417 } 418 psFree(iter); 419 } 420 421 psArray *kernels = options->kernels->data[index] = psArrayAllocEmpty(ARRAY_BUFFER); // Match kernels 422 { 423 psString regex = NULL; // Regular expression 424 psStringAppend(®ex, "^%s$", PM_SUBTRACTION_ANALYSIS_KERNEL); 425 psMetadataIterator *iter = psMetadataIteratorAlloc(readout->analysis, PS_LIST_HEAD, regex); 426 psFree(regex); 427 psMetadataItem *item = NULL;// Item from iteration 428 while ((item = psMetadataGetAndIncrement(iter))) { 429 assert(item->type == PS_DATA_UNKNOWN); 430 pmSubtractionKernels *kernel = item->data.V; // Kernel used in subtraction 431 kernels = psArrayAdd(kernels, ARRAY_BUFFER, kernel); 432 } 433 psFree(iter); 434 } 435 psAssert((regions)->n == (kernels)->n, "Number of match regions and kernels should match"); 436 437 // Record chi^2 438 { 439 double sum = 0.0; // Sum of chi^2 440 int num = 0; // Number of measurements of chi^2 441 psString regex = NULL; // Regular expression 442 psStringAppend(®ex, "^%s$", PM_SUBTRACTION_ANALYSIS_KERNEL); 443 psMetadataIterator *iter = psMetadataIteratorAlloc(readout->analysis, PS_LIST_HEAD, regex); 444 psFree(regex); 445 psMetadataItem *item = NULL;// Item from iteration 446 while ((item = psMetadataGetAndIncrement(iter))) { 447 assert(item->type == PS_DATA_UNKNOWN); 448 pmSubtractionKernels *kernels = item->data.V; // Convolution kernels 449 sum += kernels->mean; 450 num++; 451 } 452 psFree(iter); 453 options->matchChi2->data.F32[index] = sum / (psImageCovarianceFactor(readout->covariance) * num); 454 } 455 456 return true; 457 } 458 459 // Kernel normalisation for convolved readout 460 bool renormKernel(pmReadout *readout, psphotStackOptions *options, int index) { 461 462 double sum = 0.0; // Sum of chi^2 463 int num = 0; // Number of measurements of chi^2 464 psString regex = NULL; // Regular expression 465 psStringAppend(®ex, "^%s$", PM_SUBTRACTION_ANALYSIS_NORM); 466 psMetadataIterator *iter = psMetadataIteratorAlloc(readout->analysis, PS_LIST_HEAD, regex); 467 psFree(regex); 468 psMetadataItem *item = NULL;// Item from iteration 469 while ((item = psMetadataGetAndIncrement(iter))) { 470 assert(item->type == PS_TYPE_F32); 471 float norm = item->data.F32; // Normalisation 472 sum += norm; 473 num++; 474 } 475 psFree(iter); 476 float conv = sum/num; // Mean normalisation from convolution 477 float stars = powf(10.0, -0.4 * options->norm->data.F32[index]); // Normalisation from stars 478 float renorm = stars / conv; // Renormalisation to apply 479 psLogMsg("psphotStack", PS_LOG_INFO, "Renormalising image %d by %f (kernel: %f, stars: %f)\n", index, renorm, conv, stars); 480 481 psBinaryOp(readout->image, readout->image, "*", psScalarAlloc(renorm, PS_TYPE_F32)); 482 psBinaryOp(readout->variance, readout->variance, "*", psScalarAlloc(PS_SQR(renorm), PS_TYPE_F32)); 483 return true; 484 } 485 486 // adjust scaling for readout (remove background, ..., determine weighting) 487 bool rescaleData(pmReadout *readout, pmConfig *config, psphotStackOptions *options, int index) { 488 489 psMetadata *stackRecipe = psMetadataLookupMetadata(NULL, config->recipes, "PPSTACK"); // ppStack recipe 490 psAssert(stackRecipe, "We've thrown an error on this before."); 491 492 // Look up appropriate values from the ppSub recipe 493 psMetadata *subRecipe = psMetadataLookupMetadata(NULL, config->recipes, "PPSUB"); // PPSUB recipe 494 psAssert(subRecipe, "recipe missing"); 495 496 psString maskValStr = psMetadataLookupStr(NULL, subRecipe, "MASK.VAL"); // Name of bits to mask going in 497 psString maskBadStr = psMetadataLookupStr(NULL, stackRecipe, "MASK.BAD"); // Name of bits to mask for bad 498 499 psImageMaskType maskVal = pmConfigMaskGet(maskValStr, config); // Bits to mask going in to pmSubtractionMatch 500 psImageMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels 501 502 // Ensure the background value is zero 503 psStats *bg = psStatsAlloc(PS_STAT_ROBUST_MEDIAN | PS_STAT_ROBUST_STDEV); // Statistics for background 504 psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS); // Random number generator 505 506 // XXX why is this in config->arguments and not recipe? 507 if (!psMetadataLookupBool(NULL, config->arguments, "PPSTACK.SKIP.BG.SUB")) { 508 if (!psImageBackground(bg, NULL, readout->image, readout->mask, maskVal | maskBad, rng)) { 509 psAbort("Can't measure background for image."); 510 // XXX we used to clear error: why is this acceptable? psErrorClear(); 511 } 512 513 float value = psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN); 514 float stdev = psStatsGetValue(bg, PS_STAT_ROBUST_STDEV); 515 516 psLogMsg("psphotStack", PS_LOG_INFO, "Correcting convolved image background by %lf (+/- %lf)", value, stdev); 517 psBinaryOp(readout->image, readout->image, "-", psScalarAlloc(value, PS_TYPE_F32)); 518 } 519 520 if (!stackRenormaliseReadout(config, readout)) { 521 psFree(rng); 522 psFree(bg); 523 return false; 524 } 525 526 // Measure the variance level for the weighting 527 if (psMetadataLookupBool(NULL, stackRecipe, "WEIGHTS")) { 528 if (!psImageBackground(bg, NULL, readout->variance, readout->mask, maskVal | maskBad, rng)) { 529 psError(PSPHOT_ERR_DATA, false, "Can't measure mean variance for image."); 367 530 psFree(rng); 368 531 psFree(bg); 369 370 // For the sake of stamps, remove nearby sources 371 psArray *stampSources = stackSourcesFilter(options->sourceLists->data[index], 372 footprint); // Filtered list of sources 373 374 bool oldThreads = pmReadoutFakeThreads(true); // Old threading state 375 if (!pmReadoutFakeFromSources(fake, readout->image->numCols, readout->image->numRows, 376 stampSources, SOURCE_MASK, NULL, NULL, options->psf, 377 minFlux, footprint + size, false, true)) { 378 psError(PPSTACK_ERR_DATA, false, "Unable to generate fake image with target PSF."); 379 psFree(fake); 380 psFree(optWidths); 381 psFree(conv); 382 return false; 383 } 384 pmReadoutFakeThreads(oldThreads); 385 386 fake->mask = psImageCopy(NULL, readout->mask, PS_TYPE_IMAGE_MASK); 387 388 // Add the background into the target image 389 psImage *bgImage = stackBackgroundModel(readout, config); // Image of background 390 psBinaryOp(fake->image, fake->image, "+", bgImage); 391 psFree(bgImage); 392 393 dumpImage(); 394 395 if (threads > 0) { 396 pmSubtractionThreadsInit(); 397 } 398 399 // Do the image matching 400 pmSubtractionKernels *kernel = psMetadataLookupPtr(&mdok, readout->analysis, PM_SUBTRACTION_ANALYSIS_KERNEL); // Conv kernel 401 if (kernel) { 402 if (!pmSubtractionMatchPrecalc(NULL, conv, fake, readout, readout->analysis, 403 stride, kernelError, covarFrac, maskVal, maskBad, maskPoor, 404 poorFrac, badFrac)) { 405 psError(psErrorCodeLast(), false, "Unable to convolve images."); 406 psFree(fake); 407 psFree(optWidths); 408 psFree(stampSources); 409 psFree(conv); 410 if (threads > 0) { 411 pmSubtractionThreadsFinalize(); 412 } 413 return false; 414 } 415 } else { 416 // Scale the input parameters 417 psVector *widthsCopy = psVectorCopy(NULL, widths, PS_TYPE_F32); // Copy of kernel widths 418 if (scale && !pmSubtractionParamsScale(&size, &footprint, widthsCopy, 419 options->inputSeeing->data.F32[index], 420 options->targetSeeing, scaleRef, scaleMin, scaleMax)) { 421 psError(psErrorCodeLast(), false, "Unable to scale kernel parameters"); 422 psFree(fake); 423 psFree(optWidths); 424 psFree(stampSources); 425 psFree(conv); 426 psFree(widthsCopy); 427 if (threads > 0) { 428 pmSubtractionThreadsFinalize(); 429 } 430 return false; 431 } 432 433 if (!pmSubtractionMatch(NULL, conv, fake, readout, footprint, stride, regionSize, spacing, 434 threshold, stampSources, stampsName, type, size, order, widthsCopy, 435 orders, inner, ringsOrder, binning, penalty, 436 optimum, optWidths, optOrder, optThresh, iter, rej, normFrac, 437 sysError, skyErr, kernelError, covarFrac, maskVal, maskBad, maskPoor, 438 poorFrac, badFrac, PM_SUBTRACTION_MODE_2)) { 439 psError(psErrorCodeLast(), false, "Unable to match images."); 440 psFree(fake); 441 psFree(optWidths); 442 psFree(stampSources); 443 psFree(conv); 444 psFree(widthsCopy); 445 if (threads > 0) { 446 pmSubtractionThreadsFinalize(); 447 } 448 return false; 449 } 450 psFree(widthsCopy); 451 } 452 453 dumpImage2(); 454 455 psFree(fake); 456 psFree(optWidths); 457 psFree(stampSources); 458 459 if (threads > 0) { 460 pmSubtractionThreadsFinalize(); 461 } 462 463 // Replace original images with convolved 464 psFree(readout->image); 465 psFree(readout->mask); 466 psFree(readout->variance); 467 psFree(readout->covariance); 468 readout->image = psMemIncrRefCounter(conv->image); 469 readout->mask = psMemIncrRefCounter(conv->mask); 470 readout->variance = psMemIncrRefCounter(conv->variance); 471 readout->covariance = psImageCovarianceTruncate(conv->covariance, COVAR_FRAC); 472 473 } 474 475 bool saveMatchData () { 476 // Extract the regions and solutions used in the image matching 477 // This stops them from being freed when we iterate back up the FPA 478 psArray *regions = options->regions->data[index] = psArrayAllocEmpty(ARRAY_BUFFER); // Match regions 479 { 480 psString regex = NULL; // Regular expression 481 psStringAppend(®ex, "^%s$", PM_SUBTRACTION_ANALYSIS_REGION); 482 psMetadataIterator *iter = psMetadataIteratorAlloc(conv->analysis, PS_LIST_HEAD, regex); 483 psFree(regex); 484 psMetadataItem *item = NULL;// Item from iteration 485 while ((item = psMetadataGetAndIncrement(iter))) { 486 assert(item->type == PS_DATA_REGION); 487 regions = psArrayAdd(regions, ARRAY_BUFFER, item->data.V); 488 } 489 psFree(iter); 532 return false; 490 533 } 491 psArray *kernels = options->kernels->data[index] = psArrayAllocEmpty(ARRAY_BUFFER); // Match kernels 492 { 493 psString regex = NULL; // Regular expression 494 psStringAppend(®ex, "^%s$", PM_SUBTRACTION_ANALYSIS_KERNEL); 495 psMetadataIterator *iter = psMetadataIteratorAlloc(conv->analysis, PS_LIST_HEAD, regex); 496 psFree(regex); 497 psMetadataItem *item = NULL;// Item from iteration 498 while ((item = psMetadataGetAndIncrement(iter))) { 499 assert(item->type == PS_DATA_UNKNOWN); 500 pmSubtractionKernels *kernel = item->data.V; // Kernel used in subtraction 501 kernels = psArrayAdd(kernels, ARRAY_BUFFER, kernel); 502 } 503 psFree(iter); 504 } 505 psAssert((regions)->n == (kernels)->n, "Number of match regions and kernels should match"); 506 } 507 508 bool saveChiSquare() { 509 // Record chi^2 510 { 511 double sum = 0.0; // Sum of chi^2 512 int num = 0; // Number of measurements of chi^2 513 psString regex = NULL; // Regular expression 514 psStringAppend(®ex, "^%s$", PM_SUBTRACTION_ANALYSIS_KERNEL); 515 psMetadataIterator *iter = psMetadataIteratorAlloc(conv->analysis, PS_LIST_HEAD, regex); 516 psFree(regex); 517 psMetadataItem *item = NULL;// Item from iteration 518 while ((item = psMetadataGetAndIncrement(iter))) { 519 assert(item->type == PS_DATA_UNKNOWN); 520 pmSubtractionKernels *kernels = item->data.V; // Convolution kernels 521 sum += kernels->mean; 522 num++; 523 } 524 psFree(iter); 525 options->matchChi2->data.F32[index] = sum / (psImageCovarianceFactor(readout->covariance) * num); 526 } 527 528 } 529 530 bool renormKernel() { 531 // Kernel normalisation 532 { 533 double sum = 0.0; // Sum of chi^2 534 int num = 0; // Number of measurements of chi^2 535 psString regex = NULL; // Regular expression 536 psStringAppend(®ex, "^%s$", PM_SUBTRACTION_ANALYSIS_NORM); 537 psMetadataIterator *iter = psMetadataIteratorAlloc(conv->analysis, PS_LIST_HEAD, regex); 538 psFree(regex); 539 psMetadataItem *item = NULL;// Item from iteration 540 while ((item = psMetadataGetAndIncrement(iter))) { 541 assert(item->type == PS_TYPE_F32); 542 float norm = item->data.F32; // Normalisation 543 sum += norm; 544 num++; 545 } 546 psFree(iter); 547 float conv = sum/num; // Mean normalisation from convolution 548 float stars = powf(10.0, -0.4 * options->norm->data.F32[index]); // Normalisation from stars 549 float renorm = stars / conv; // Renormalisation to apply 550 psLogMsg("ppStack", PS_LOG_INFO, "Renormalising image %d by %f (kernel: %f, stars: %f)\n", 551 index, renorm, conv, stars); 552 psBinaryOp(readout->image, readout->image, "*", psScalarAlloc(renorm, PS_TYPE_F32)); 553 psBinaryOp(readout->variance, readout->variance, "*", psScalarAlloc(PS_SQR(renorm), PS_TYPE_F32)); 554 } 555 556 } 534 options->weightings->data.F32[index] = 1.0 / (psStatsGetValue(bg, PS_STAT_ROBUST_MEDIAN) * psImageCovarianceFactor(readout->covariance)); 535 } else { 536 options->weightings->data.F32[index] = 1.0; 537 } 538 psLogMsg("psphotStack", PS_LOG_INFO, "Weighting for image %d is %f\n", index, options->weightings->data.F32[index]); 539 540 psFree(rng); 541 psFree(bg); 542 return true; 543 } 544 545 # define NOISE_FRACTION 0.01 // Set minimum flux to this fraction of noise 546 # define SOURCE_MASK (PM_SOURCE_MODE_FAIL | PM_SOURCE_MODE_DEFECT | PM_SOURCE_MODE_SATURATED | PM_SOURCE_MODE_CR_LIMIT | PM_SOURCE_MODE_EXT_LIMIT) // Mask to apply to input sources 547 548 // generate a fake readout against which to PSF match 549 pmReadout *makeFakeReadout(pmConfig *config, pmReadout *readoutSrc, psArray *sources, pmPSF *psf, psImageMaskType maskVal, int fullSize) { 550 551 pmReadout *fake = pmReadoutAlloc(NULL); // Fake readout with target PSF 552 553 psStats *bg = psStatsAlloc(PS_STAT_ROBUST_STDEV); // Statistics for background 554 psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS); // Random number generator 555 if (!psImageBackground(bg, NULL, readoutSrc->image, readoutSrc->mask, maskVal, rng)) { 556 psError(PSPHOT_ERR_DATA, false, "Can't measure background for image."); 557 psFree(fake); 558 psFree(bg); 559 psFree(rng); 560 return NULL; 561 } 562 float minFlux = NOISE_FRACTION * bg->robustStdev; // Minimum flux level for fake image 563 psFree(rng); 564 psFree(bg); 565 566 bool oldThreads = pmReadoutFakeThreads(true); // Old threading state 567 if (!pmReadoutFakeFromSources(fake, readoutSrc->image->numCols, readoutSrc->image->numRows, sources, SOURCE_MASK, NULL, NULL, psf, minFlux, fullSize, false, true)) { 568 psError(PSPHOT_ERR_DATA, false, "Unable to generate fake image with target PSF."); 569 psFree(fake); 570 return NULL; 571 } 572 pmReadoutFakeThreads(oldThreads); 573 574 fake->mask = psImageCopy(NULL, readoutSrc->mask, PS_TYPE_IMAGE_MASK); 575 576 // Add the background into the target image 577 psImage *bgImage = stackBackgroundModel(readoutSrc, config); // Image of background 578 psBinaryOp(fake->image, fake->image, "+", bgImage); 579 psFree(bgImage); 580 581 return fake; 582 } 583 584 // set the widths 585 psVector *SetOptWidths (bool *optimum, psMetadata *recipe) { 586 587 bool status; 588 589 *optimum = psMetadataLookupBool(&status, recipe, "OPTIMUM"); // Derive optimum parameters? 590 psAssert (status, "missing recipe value %s", "OPTIMUM"); 591 592 psVector *optWidths = NULL; // Vector with FWHMs for optimum search 593 594 if (*optimum) { 595 float optMin = psMetadataLookupF32(&status, recipe, "OPTIMUM.MIN"); // Minimum width for search 596 psAssert (status, "missing recipe value %s", "OPTIMUM.MIN"); 597 598 float optMax = psMetadataLookupF32(&status, recipe, "OPTIMUM.MAX"); // Maximum width for search 599 psAssert (status, "missing recipe value %s", "OPTIMUM.MAX"); 600 601 float optStep = psMetadataLookupF32(&status, recipe, "OPTIMUM.STEP"); // Step for search 602 psAssert (status, "missing recipe value %s", "OPTIMUM.STEP"); 603 604 optWidths = psVectorCreate(optWidths, optMin, optMax, optStep, PS_TYPE_F32); 605 } 606 607 return optWidths; 608 } -
branches/czw_branch/20100427/psphot/src/psphotStackParseCamera.c
r27657 r28029 25 25 psMetadata *input = item->data.md; // The input metadata of interest 26 26 27 // look for 'IMAGE', 'MASK', 'VARIANCE' in folder (only 'IMAGE' is required) 28 29 psString image = psMetadataLookupStr(NULL, input, "IMAGE"); // Name of image 30 if (!image || strlen(image) == 0) { 31 psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Component %s lacks IMAGE of type STR", item->name); 27 pmFPAfile *rawInputFile = NULL; 28 pmFPAfile *cnvInputFile = NULL; 29 30 // RAW (unconvolved) input data (RAW:IMAGE, RAW:MASK, RAW:VARIANCE, RAW:PSF) 31 psString rawImage = psMetadataLookupStr(&status, input, "RAW:IMAGE"); // Name of image 32 if (rawImage && strlen(rawImage) > 0) { 33 rawInputFile = defineFile(config, NULL, "PSPHOT.STACK.INPUT.RAW", rawImage, PM_FPA_FILE_IMAGE); // File for image 34 if (!rawInputFile) { 35 psError(PS_ERR_UNKNOWN, false, "Unable to define file from image %d (%s)", i, rawImage); 36 return false; 37 } 38 psString mask = psMetadataLookupStr(&status, input, "RAW:MASK"); // Name of mask 39 if (mask && strlen(mask) > 0) { 40 if (!defineFile(config, rawInputFile, "PSPHOT.STACK.MASK.RAW", mask, PM_FPA_FILE_MASK)) { 41 psError(PS_ERR_UNKNOWN, false, "Unable to define file from mask %d (%s)", i, mask); 42 return false; 43 } 44 } 45 psString variance = psMetadataLookupStr(&status, input, "RAW:VARIANCE"); // Name of variance map 46 if (variance && strlen(variance) > 0) { 47 if (!defineFile(config, rawInputFile, "PSPHOT.STACK.VARIANCE.RAW", variance, PM_FPA_FILE_VARIANCE)) { 48 psError(PS_ERR_UNKNOWN, false, "Unable to define file from variance %d (%s)", i, variance); 49 return false; 50 } 51 } 52 psString psf = psMetadataLookupStr(&status, input, "RAW:PSF"); // Name of mask 53 if (psf && strlen(psf) > 0) { 54 if (!defineFile(config, rawInputFile, "PSPHOT.STACK.PSF.RAW", psf, PM_FPA_FILE_PSF)) { 55 psError(PS_ERR_UNKNOWN, false, "Unable to define file from psf %d (%s)", i, psf); 56 return false; 57 } 58 } 59 } 60 61 // CNV (convolved) input data (CNV:IMAGE, CNV:MASK, CNV:VARIANCE, CNV:PSF) 62 psString cnvImage = psMetadataLookupStr(&status, input, "CNV:IMAGE"); // Name of image 63 if (cnvImage && strlen(cnvImage) > 0) { 64 cnvInputFile = defineFile(config, NULL, "PSPHOT.STACK.INPUT.CNV", cnvImage, PM_FPA_FILE_IMAGE); // File for image 65 if (!cnvInputFile) { 66 psError(PS_ERR_UNKNOWN, false, "Unable to define file from image %d (%s)", i, cnvImage); 67 return false; 68 } 69 psString mask = psMetadataLookupStr(&status, input, "CNV:MASK"); // Name of mask 70 if (mask && strlen(mask) > 0) { 71 if (!defineFile(config, cnvInputFile, "PSPHOT.STACK.MASK.CNV", mask, PM_FPA_FILE_MASK)) { 72 psError(PS_ERR_UNKNOWN, false, "Unable to define file from mask %d (%s)", i, mask); 73 return false; 74 } 75 } 76 psString variance = psMetadataLookupStr(&status, input, "CNV:VARIANCE"); // Name of variance map 77 if (variance && strlen(variance) > 0) { 78 if (!defineFile(config, cnvInputFile, "PSPHOT.STACK.VARIANCE.CNV", variance, PM_FPA_FILE_VARIANCE)) { 79 psError(PS_ERR_UNKNOWN, false, "Unable to define file from variance %d (%s)", i, variance); 80 return false; 81 } 82 } 83 psString psf = psMetadataLookupStr(&status, input, "CNV:PSF"); // Name of mask 84 if (psf && strlen(psf) > 0) { 85 if (!defineFile(config, cnvInputFile, "PSPHOT.STACK.PSF.CNV", psf, PM_FPA_FILE_PSF)) { 86 psError(PS_ERR_UNKNOWN, false, "Unable to define file from psf %d (%s)", i, psf); 87 return false; 88 } 89 } 90 } 91 92 if (!rawInputFile && !cnvInputFile) { 93 psError(PS_ERR_BAD_PARAMETER_VALUE, true, "Component %s (%d) lacks IMAGE of type STR", item->name, i); 32 94 return false; 33 95 } 34 pmFPAfile *imageFile = defineFile(config, NULL, "PSPHOT.INPUT", image, PM_FPA_FILE_IMAGE); // File for image 35 if (!imageFile) { 36 psError(PS_ERR_UNKNOWN, false, "Unable to define file from image %d (%s)", i, image); 37 return false; 38 } 39 40 psString mask = psMetadataLookupStr(&status, input, "MASK"); // Name of mask 41 if (mask && strlen(mask) > 0) { 42 if (!defineFile(config, imageFile, "PSPHOT.MASK", mask, PM_FPA_FILE_MASK)) { 43 psError(PS_ERR_UNKNOWN, false, "Unable to define file from mask %d (%s)", i, mask); 44 return false; 45 } 46 } 47 48 psString variance = psMetadataLookupStr(&status, input, "VARIANCE"); // Name of variance map 49 if (variance && strlen(variance) > 0) { 50 if (!defineFile(config, imageFile, "PSPHOT.VARIANCE", variance, PM_FPA_FILE_VARIANCE)) { 51 psError(PS_ERR_UNKNOWN, false, "Unable to define file from variance %d (%s)", i, variance); 52 return false; 53 } 54 } 55 // the output sources are carried on the input->fpa structures 56 pmFPAfile *outsources = pmFPAfileDefineOutputFromFile (config, imageFile, "PSPHOT.STACK.OUTPUT"); 57 if (!outsources) { 58 psError(PSPHOT_ERR_CONFIG, false, "Cannot find a rule for PSPHOT.STACK.OUTPUT"); 59 return false; 60 } 61 outsources->save = true; 62 outsources->fileID = i; // this is used to generate output names 96 97 psString sources = psMetadataLookupStr(&status, input, "SOURCES"); // Name of mask 98 // pmFPAfile *srcInputFile = rawInputFile ? rawInputFile : cnvInputFile; 99 if (sources && strlen(sources) > 0) { 100 if (!defineFile(config, NULL, "PSPHOT.STACK.SOURCES", sources, PM_FPA_FILE_CMF)) { 101 psError(PS_ERR_UNKNOWN, false, "Unable to define file from sources %d (%s)", i, sources); 102 return false; 103 } 104 } 105 106 // generate an pmFPAimage for the output convolved image 107 // XXX output of these files should be optional 108 { 109 pmFPAfile *outputImage = pmFPAfileDefineOutput(config, NULL, "PSPHOT.STACK.OUTPUT.IMAGE"); 110 if (!outputImage) { 111 psError(PSPHOT_ERR_CONFIG, false, "Trouble defining PSPHOT.STACK.OUTPUT.IMAGE"); 112 return false; 113 } 114 outputImage->save = true; 115 outputImage->fileID = i; // this is used to generate output names 116 117 pmFPAfile *outputMask = pmFPAfileDefineOutput(config, outputImage->fpa, "PSPHOT.STACK.OUTPUT.MASK"); 118 if (!outputMask) { 119 psError(PS_ERR_IO, false, _("Unable to generate output file from PSPHOT.STACK.OUTPUT.MASK")); 120 return NULL; 121 } 122 if (outputMask->type != PM_FPA_FILE_MASK) { 123 psError(PS_ERR_IO, true, "PSPHOT.STACK.OUTPUT.MASK is not of type MASK"); 124 return NULL; 125 } 126 outputMask->save = true; 127 outputMask->fileID = i; // this is used to generate output names 128 129 pmFPAfile *outputVariance = pmFPAfileDefineOutput(config, outputImage->fpa, "PSPHOT.STACK.OUTPUT.VARIANCE"); 130 if (!outputVariance) { 131 psError(PS_ERR_IO, false, _("Unable to generate output file from PSPHOT.STACK.OUTPUT.VARIANCE")); 132 return NULL; 133 } 134 if (outputVariance->type != PM_FPA_FILE_VARIANCE) { 135 psError(PS_ERR_IO, true, "PSPHOT.STACK.OUTPUT.VARIANCE is not of type VARIANCE"); 136 return NULL; 137 } 138 outputVariance->save = true; 139 outputVariance->fileID = i; // this is used to generate output names 140 141 // the output sources are carried on the outputImage->fpa structures 142 pmFPAfile *outsources = pmFPAfileDefineOutputFromFile (config, outputImage, "PSPHOT.STACK.OUTPUT"); 143 if (!outsources) { 144 psError(PSPHOT_ERR_CONFIG, false, "Cannot find a rule for PSPHOT.STACK.OUTPUT"); 145 return false; 146 } 147 outsources->save = true; 148 outsources->fileID = i; // this is used to generate output names 149 } 63 150 } 64 151 psMetadataRemoveKey(config->arguments, "FILENAMES"); … … 71 158 72 159 // generate an pmFPAimage for the chisqImage 73 pmFPAfile *chisqImage = pmFPAfileDefineOutput(config, NULL, "PSPHOT.CHISQ.IMAGE"); 74 if (!chisqImage) { 75 psError(PSPHOT_ERR_CONFIG, false, "Trouble defining PSPHOT.CHISQ.IMAGE"); 76 return false; 77 } 78 chisqImage->save = true; 79 80 pmFPAfile *chisqMask = pmFPAfileDefineOutput(config, chisqImage->fpa, "PSPHOT.CHISQ.MASK"); 81 if (!chisqMask) { 82 psError(PS_ERR_IO, false, _("Unable to generate output file from PSPHOT.CHISQ.MASK")); 83 return NULL; 84 } 85 if (chisqMask->type != PM_FPA_FILE_MASK) { 86 psError(PS_ERR_IO, true, "PSPHOT.CHISQ.MASK is not of type MASK"); 87 return NULL; 88 } 89 chisqMask->save = true; 90 91 pmFPAfile *chisqVariance = pmFPAfileDefineOutput(config, chisqImage->fpa, "PSPHOT.CHISQ.VARIANCE"); 92 if (!chisqVariance) { 93 psError(PS_ERR_IO, false, _("Unable to generate output file from PSPHOT.CHISQ.VARIANCE")); 94 return NULL; 95 } 96 if (chisqVariance->type != PM_FPA_FILE_VARIANCE) { 97 psError(PS_ERR_IO, true, "PSPHOT.CHISQ.VARIANCE is not of type VARIANCE"); 98 return NULL; 99 } 100 chisqVariance->save = true; 101 102 # if (0) 103 // define the additional input/output files associated with psphot 104 // XXX figure out which files are needed by psphotStack 105 if (false && !psphotDefineFiles (config, input)) { 106 psError(PSPHOT_ERR_CONFIG, false, "Trouble defining the additional input/output files"); 107 return false; 108 } 109 # endif 160 // XXX output of these files should be optional 161 { 162 pmFPAfile *chisqImage = pmFPAfileDefineOutput(config, NULL, "PSPHOT.CHISQ.IMAGE"); 163 if (!chisqImage) { 164 psError(PSPHOT_ERR_CONFIG, false, "Trouble defining PSPHOT.CHISQ.IMAGE"); 165 return false; 166 } 167 chisqImage->save = true; 168 169 pmFPAfile *chisqMask = pmFPAfileDefineOutput(config, chisqImage->fpa, "PSPHOT.CHISQ.MASK"); 170 if (!chisqMask) { 171 psError(PS_ERR_IO, false, _("Unable to generate output file from PSPHOT.CHISQ.MASK")); 172 return NULL; 173 } 174 if (chisqMask->type != PM_FPA_FILE_MASK) { 175 psError(PS_ERR_IO, true, "PSPHOT.CHISQ.MASK is not of type MASK"); 176 return NULL; 177 } 178 chisqMask->save = true; 179 180 pmFPAfile *chisqVariance = pmFPAfileDefineOutput(config, chisqImage->fpa, "PSPHOT.CHISQ.VARIANCE"); 181 if (!chisqVariance) { 182 psError(PS_ERR_IO, false, _("Unable to generate output file from PSPHOT.CHISQ.VARIANCE")); 183 return NULL; 184 } 185 if (chisqVariance->type != PM_FPA_FILE_VARIANCE) { 186 psError(PS_ERR_IO, true, "PSPHOT.CHISQ.VARIANCE is not of type VARIANCE"); 187 return NULL; 188 } 189 chisqVariance->save = true; 190 } 110 191 111 192 psTrace("psphot", 1, "Done with psphotStackParseCamera...\n"); … … 145 226 } 146 227 228 229 /*** 230 * 231 * psphotStack : 232 233 * * inputs: 234 * * unconvolved images 235 * * raw convolved images 236 * * psfs (unconvolved or convolved?) 237 * * sources 238 239 * optionally convolve the unconvolved or the raw inputs 240 * optionally perform no convolutions 241 * optionally save the psf-matched images 242 243 */ 244 245 246 # if (0) 247 // define the additional input/output files associated with psphot 248 // XXX figure out which files are needed by psphotStack 249 if (false && !psphotDefineFiles (config, input)) { 250 psError(PSPHOT_ERR_CONFIG, false, "Trouble defining the additional input/output files"); 251 return false; 252 } 253 # endif 254 -
branches/czw_branch/20100427/psphot/src/psphotStackReadout.c
r27657 r28029 1 1 # include "psphotInternal.h" 2 3 # define STACK_RAW "PSPHOT.STACK.INPUT.RAW" 4 # define STACK_OUT "PSPHOT.STACK.OUTPUT.IMAGE" 2 5 3 6 bool psphotStackReadout (pmConfig *config, const pmFPAview *view) { … … 17 20 PS_ASSERT_PTR_NON_NULL (breakPt, false); 18 21 22 // we have 3 relevant files: RAW, CNV, OUT 23 19 24 // set the photcode for each image 20 if (!psphotAddPhotcode (config, view )) {25 if (!psphotAddPhotcode (config, view, STACK_OUT)) { 21 26 psError (PSPHOT_ERR_CONFIG, false, "trouble defining the photcode"); 22 27 return false; … … 24 29 25 30 // Generate the mask and weight images 26 if (!psphotSetMaskAndVariance (config, view)) { 27 return psphotReadoutCleanup (config, view); 31 // XXX this should be done before we perform the convolutions 32 if (!psphotSetMaskAndVariance (config, view, STACK_RAW)) { 33 return psphotReadoutCleanup (config, view, STACK_OUT); 28 34 } 29 35 if (!strcasecmp (breakPt, "NOTHING")) { 30 return psphotReadoutCleanup (config, view );36 return psphotReadoutCleanup (config, view, STACK_OUT); 31 37 } 32 38 33 39 // generate a background model (median, smoothed image) 34 40 // XXX I think this is not defined correctly for an array of images. 35 if (!psphotModelBackground (config, view)) { 36 return psphotReadoutCleanup (config, view); 41 // XXX probably need to subtract the model (same model?) for both RAW and OUT 42 if (!psphotModelBackground (config, view, STACK_RAW)) { 43 return psphotReadoutCleanup (config, view, STACK_OUT); 37 44 } 38 if (!psphotSubtractBackground (config, view )) {39 return psphotReadoutCleanup (config, view );45 if (!psphotSubtractBackground (config, view, STACK_RAW)) { 46 return psphotReadoutCleanup (config, view, STACK_OUT); 40 47 } 41 48 if (!strcasecmp (breakPt, "BACKMDL")) { 42 return psphotReadoutCleanup (config, view );49 return psphotReadoutCleanup (config, view, STACK_OUT); 43 50 } 44 51 … … 47 54 if (!psphotLoadPSF (config, view)) { 48 55 psError (PSPHOT_ERR_UNKNOWN, false, "error loading psf model"); 49 return psphotReadoutCleanup (config, view );56 return psphotReadoutCleanup (config, view, STACK_OUT); 50 57 } 51 58 52 if (!psphotStackChisqImage(config, view )) {59 if (!psphotStackChisqImage(config, view, STACK_RAW, STACK_OUT)) { 53 60 psError (PSPHOT_ERR_UNKNOWN, false, "failure to generate chisq image"); 54 return psphotReadoutCleanup (config, view );61 return psphotReadoutCleanup (config, view, STACK_OUT); 55 62 } 56 63 if (!strcasecmp (breakPt, "CHISQ")) { 57 return psphotReadoutCleanup (config, view );64 return psphotReadoutCleanup (config, view, STACK_OUT); 58 65 } 59 66 60 67 // find the detections (by peak and/or footprint) in the image. 61 68 // This finds the detections on Chisq image as well as the individuals 62 if (!psphotFindDetections (config, view, true)) { // pass 169 if (!psphotFindDetections (config, view, STACK_RAW, true)) { // pass 1 63 70 // this only happens if we had an error in psphotFindDetections 64 71 psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis"); 65 return psphotReadoutCleanup (config, view); 72 return psphotReadoutCleanup (config, view, STACK_OUT); 73 } 74 75 // copy the detections from RAW to OUT 76 if (!psphotCopySources (config, view, STACK_OUT, STACK_RAW)) { 77 psError (PSPHOT_ERR_UNKNOWN, false, "failure in peak analysis"); 78 return psphotReadoutCleanup (config, view, STACK_OUT); 66 79 } 67 80 68 81 // construct sources and measure basic stats (saved on detections->newSources) 69 82 // only run this on detections from the input images, not chisq image 70 if (!psphotSourceStats (config, view, true)) { // pass 183 if (!psphotSourceStats (config, view, STACK_OUT, true)) { // pass 1 71 84 psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources"); 72 return psphotReadoutCleanup (config, view );85 return psphotReadoutCleanup (config, view, STACK_OUT); 73 86 } 74 87 75 // *** generate the objects (whichunify the sources from the different images)76 psArray *objects = psphotMatchSources (config, view );88 // generate the objects (object unify the sources from the different images) 89 psArray *objects = psphotMatchSources (config, view, STACK_OUT); 77 90 78 91 // construct sources for the newly-generated sources (from other images) 79 if (!psphotSourceStats (config, view, false)) { // pass 192 if (!psphotSourceStats (config, view, STACK_OUT, false)) { // pass 1 80 93 psError(PSPHOT_ERR_UNKNOWN, false, "failure to generate sources"); 81 return psphotReadoutCleanup (config, view );94 return psphotReadoutCleanup (config, view, STACK_OUT); 82 95 } 83 96 … … 85 98 // if (!psphotDeblendSatstars (config, view)) { 86 99 // psError (PSPHOT_ERR_UNKNOWN, false, "failed on satstar deblend analysis"); 87 // return psphotReadoutCleanup (config, view );100 // return psphotReadoutCleanup (config, view, STACK_OUT); 88 101 // } 89 102 … … 91 104 // if (!psphotBasicDeblend (config, view)) { 92 105 // psError (PSPHOT_ERR_UNKNOWN, false, "failed on deblend analysis"); 93 // return psphotReadoutCleanup (config, view );106 // return psphotReadoutCleanup (config, view, STACK_OUT); 94 107 // } 95 108 96 109 // classify sources based on moments, brightness 97 110 // only run this on detections from the input images, not chisq image 98 if (!psphotRoughClass (config, view )) {111 if (!psphotRoughClass (config, view, STACK_OUT)) { 99 112 psError (PSPHOT_ERR_UNKNOWN, false, "failed to determine rough classifications"); 100 return psphotReadoutCleanup (config, view );113 return psphotReadoutCleanup (config, view, STACK_OUT); 101 114 } 102 115 // if we were not supplied a PSF model, determine the IQ stats here (detections->newSources) 103 116 // only run this on detections from the input images, not chisq image 104 if (!psphotImageQuality (config, view )) { // pass 1117 if (!psphotImageQuality (config, view, STACK_OUT)) { // pass 1 105 118 psError (PSPHOT_ERR_UNKNOWN, false, "failed to measure image quality"); 106 return psphotReadoutCleanup (config, view);119 return psphotReadoutCleanup (config, view, STACK_OUT); 107 120 } 108 121 if (!strcasecmp (breakPt, "MOMENTS")) { 109 return psphotReadoutCleanup (config, view );122 return psphotReadoutCleanup (config, view, STACK_OUT); 110 123 } 111 124 112 125 // use bright stellar objects to measure PSF if we were supplied a PSF for any input file, 113 126 // this step is skipped 114 if (!psphotChoosePSF (config, view )) { // pass 1127 if (!psphotChoosePSF (config, view, STACK_OUT)) { // pass 1 115 128 psLogMsg ("psphot", 3, "failure to construct a psf model"); 116 return psphotReadoutCleanup (config, view );129 return psphotReadoutCleanup (config, view, STACK_OUT); 117 130 } 118 131 if (!strcasecmp (breakPt, "PSFMODEL")) { 119 return psphotReadoutCleanup (config, view );132 return psphotReadoutCleanup (config, view, STACK_OUT); 120 133 } 121 134 122 // include externally-supplied sources123 // XXX fix this in the new multi-input context124 // psphotLoadExtSources (config, view); // pass 1125 126 135 // construct an initial model for each object, set the radius to fitRadius, set circular fit mask 127 psphotGuessModels (config, view );136 psphotGuessModels (config, view, STACK_OUT); 128 137 129 138 // merge the newly selected sources into the existing list 130 139 // NOTE: merge OLD and NEW 131 psphotMergeSources (config, view );140 psphotMergeSources (config, view, STACK_OUT); 132 141 133 142 // linear PSF fit to source peaks, subtract the models from the image (in PSF mask) 134 143 psphotFitSourcesLinearStack (config, objects, FALSE); 135 psFree (objects);136 144 137 145 // identify CRs and extended sources 138 psphotSourceSize (config, view, TRUE);146 psphotSourceSize (config, view, STACK_OUT, TRUE); 139 147 140 148 // measure aperture photometry corrections 141 if (!psphotApResid (config, view)) { 149 if (!psphotApResid (config, view, STACK_OUT)) { 150 psFree (objects); 142 151 psLogMsg ("psphot", 3, "failed on psphotApResid"); 143 return psphotReadoutCleanup (config, view );152 return psphotReadoutCleanup (config, view, STACK_OUT); 144 153 } 145 154 155 psphotStackObjectsUnifyPosition (objects); 156 psphotRadialAperturesByObject (config, objects, view, STACK_OUT); 157 158 psphotExtendedSourceAnalysisByObject (config, objects, view, STACK_OUT); // pass 1 (detections->allSources) 159 psphotExtendedSourceFits (config, view, STACK_OUT); // pass 1 (detections->allSources) 160 146 161 // calculate source magnitudes 147 psphotMagnitudes(config, view );162 psphotMagnitudes(config, view, STACK_OUT); 148 163 149 if (!psphotEfficiency(config, view )) {164 if (!psphotEfficiency(config, view, STACK_OUT)) { 150 165 psErrorStackPrint(stderr, "Unable to determine detection efficiencies from fake sources"); 151 166 psErrorClear(); … … 156 171 157 172 // replace background in residual image 158 psphotSkyReplace (config, view );173 psphotSkyReplace (config, view, STACK_RAW); 159 174 160 175 // drop the references to the image pixels held by each source 161 psphotSourceFreePixels (config, view );176 psphotSourceFreePixels (config, view, STACK_OUT); 162 177 163 178 // remove chisq image from config->file:PSPHOT.INPUT (why?) 164 psphotStackRemoveChisqFromInputs(config); 179 psphotStackRemoveChisqFromInputs(config, STACK_RAW); 180 181 psFree (objects); 165 182 166 183 // create the exported-metadata and free local data 167 return psphotReadoutCleanup (config, view );184 return psphotReadoutCleanup (config, view, STACK_OUT); 168 185 } 169 186 -
branches/czw_branch/20100427/psphot/src/psphotSubtractBackground.c
r27657 r28029 4 4 // generate the median in NxN boxes, clipping heavily 5 5 // linear interpolation to generate full-scale model 6 bool psphotSubtractBackgroundReadout (pmConfig *config, const pmFPAview *view, const char *file name, int index, psMetadata *recipe)6 bool psphotSubtractBackgroundReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) 7 7 { 8 8 bool status = true; … … 13 13 14 14 // find the currently selected readout 15 pmFPAfile *file = pmFPAfileSelectSingle(config->files, file name, index); // File of interest15 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 16 16 17 17 pmFPA *inFPA = file->fpa; … … 124 124 } 125 125 126 bool psphotSubtractBackground (pmConfig *config, const pmFPAview *view )126 bool psphotSubtractBackground (pmConfig *config, const pmFPAview *view, const char *filerule) 127 127 { 128 128 bool status = false; … … 137 137 // loop over the available readouts 138 138 for (int i = 0; i < num; i++) { 139 if (!psphotSubtractBackgroundReadout (config, view, "PSPHOT.INPUT", i, recipe)) {139 if (!psphotSubtractBackgroundReadout (config, view, filerule, i, recipe)) { 140 140 psError (PSPHOT_ERR_CONFIG, false, "failed to subtract background for PSPHOT.INPUT entry %d", i); 141 141 return false;
Note:
See TracChangeset
for help on using the changeset viewer.
