Changeset 13035
- Timestamp:
- Apr 25, 2007, 3:35:50 PM (19 years ago)
- Location:
- trunk/psphot/src
- Files:
-
- 1 added
- 19 edited
-
Makefile.am (modified) (3 diffs)
-
psphot.h (modified) (1 diff)
-
psphotAddNoise.c (modified) (1 diff)
-
psphotBlendFit.c (modified) (3 diffs)
-
psphotChoosePSF.c (modified) (1 diff)
-
psphotEnsemblePSF.c (modified) (2 diffs)
-
psphotFitSet.c (modified) (4 diffs)
-
psphotFitSourcesLinear.c (added)
-
psphotGrowthCurve.c (modified) (1 diff)
-
psphotGuessModels.c (modified) (1 diff)
-
psphotMakeResiduals.c (modified) (3 diffs)
-
psphotModelTest.c (modified) (5 diffs)
-
psphotOutput.c (modified) (4 diffs)
-
psphotRadialPlot.c (modified) (1 diff)
-
psphotRadiusChecks.c (modified) (5 diffs)
-
psphotReadout.c (modified) (3 diffs)
-
psphotReplaceUnfit.c (modified) (7 diffs)
-
psphotSourceFits.c (modified) (14 diffs)
-
psphotTestPSF.c (modified) (1 diff)
-
psphotWeightBias.c (modified) (3 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/psphot/src/Makefile.am
r13007 r13035 31 31 psphotChoosePSF.c \ 32 32 psphotGuessModels.c \ 33 psphot EnsemblePSF.c\33 psphotFitSourcesLinear.c \ 34 34 psphotBlendFit.c \ 35 35 psphotReplaceUnfit.c \ … … 47 47 psphotModelTest.c \ 48 48 psphotFitSet.c \ 49 psphotWeightBias.c \50 49 psphotSourceFreePixels.c \ 51 50 psphotSummaryPlots.c \ 52 psphotTestPSF.c \53 51 psphotMergeSources.c \ 54 52 psphotReadoutCleanup.c \ … … 58 56 psphotMosaicSubimage.c \ 59 57 psphotMakeResiduals.c \ 60 psphotTestSourceOutput.c \61 58 psphotAddNoise.c 62 59 -
trunk/psphot/src/psphot.h
r13007 r13035 41 41 bool psphotMomentsStats (pmReadout *readout, psMetadata *recipe, psArray *sources); 42 42 bool psphotEnsemblePSF (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, bool final); 43 bool psphotFitSourcesLinear (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf, bool final); 43 44 bool psphotGuessModels (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf); 44 45 bool psphotBlendFit (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf); -
trunk/psphot/src/psphotAddNoise.c
r12922 r13035 67 67 PAR[PM_PAR_SXY] = newshape.sxy; 68 68 69 if (add) { 70 pmModelAdd (source->weight, source->mask, model, false, false); 71 } else { 72 pmModelSub (source->weight, source->mask, model, false, false); 73 } 69 // XXX if we use pmSourceOp, the size (and possibly Io) will not be respected 70 pmSourceOp (source, PM_MODEL_OP_FULL, add); 74 71 75 72 // restore original values -
trunk/psphot/src/psphotBlendFit.c
r13006 r13035 43 43 44 44 // limit selection to some SN limit 45 // XXX this should use peak? 45 46 if (source->moments == NULL) continue; 46 47 if (source->moments->SN < FIT_SN_LIM) continue; 47 48 49 // XXX this should use peak? 48 50 if (source->moments->x < AnalysisRegion.x0) continue; 49 51 if (source->moments->y < AnalysisRegion.y0) continue; … … 65 67 // replace object in image 66 68 if (source->mode & PM_SOURCE_MODE_SUBTRACTED) { 67 pm ModelAdd (source->pixels, source->mask, source->modelPSF, false, false);69 pmSourceAdd (source, PM_MODEL_OP_FULL); 68 70 } 69 71 Nfit ++; 70 72 71 73 // try fitting PSFs, then try extended sources 74 // these functions subtract the resulting fitted source (XXX and update the modelFlux?) 72 75 if (psphotFitBlend (readout, source, psf)) { 73 76 psTrace ("psphot", 5, "source at %7.1f, %7.1f is psf", source->moments->x, source->moments->y); … … 84 87 Nfail ++; 85 88 86 // re-subtract PSF for object, leave local sky 87 pmModelSub (source->pixels, source->mask, source->modelPSF, false, false); 89 // re-subtract the object, leave local sky 90 pmSourceCacheModel (source); 91 pmSourceSub (source, PM_MODEL_OP_FULL); 88 92 source->mode |= PM_SOURCE_MODE_SUBTRACTED; 89 93 source->mode |= PM_SOURCE_MODE_TEMPSUB; -
trunk/psphot/src/psphotChoosePSF.c
r13011 r13035 239 239 240 240 // set the mask and subtract the PSF model 241 psImageKeepCircle (source->mask, x, y, RADIUS, "OR", PM_MASK_MARK); 242 pmModelSub (source->pixels, source->mask, source->modelPSF, false, false); 243 psImageKeepCircle (source->mask, x, y, RADIUS, "AND", PS_NOT_U8(PM_MASK_MARK)); 241 // XXX should we be using maskObj? should we be unsetting the mask? 242 // use pmModelSub because modelFlux has not been generated 243 assert (source->maskObj); 244 psImageKeepCircle (source->maskObj, x, y, RADIUS, "OR", PM_MASK_MARK); 245 pmModelSub (source->pixels, source->maskObj, source->modelPSF, PM_MODEL_OP_FULL); 246 psImageKeepCircle (source->maskObj, x, y, RADIUS, "AND", PS_NOT_U8(PM_MASK_MARK)); 244 247 } 245 248 -
trunk/psphot/src/psphotEnsemblePSF.c
r13006 r13035 91 91 92 92 // fill in the model pixel values 93 // XXX this function should just work with the (normalized) modelFlux pixels 94 // XXX we should not have to create a complete copy 95 // XXX review the use of the object mask: can we use a copied mask? 93 96 psImageInit (fitSource->pixels, 0.0); 94 97 psImageKeepCircle (fitSource->mask, x, y, model->radiusFit, "OR", PM_MASK_MARK); 95 pmModelAdd (fitSource->pixels, fitSource->mask, model, false, false);98 pmModelAdd (fitSource->pixels, fitSource->mask, model, PM_MODEL_OP_FULL); 96 99 97 100 // save source in list … … 218 221 219 222 // subtract object 220 pmModelSub (Ri->pixels, Ri->mask, model, false, false);223 pmModelSub (Ri->pixels, Ri->mask, model, PM_MODEL_OP_FULL); 221 224 Ri->mode |= PM_SOURCE_MODE_SUBTRACTED; 222 225 if (!final) Ri->mode |= PM_SOURCE_MODE_TEMPSUB; -
trunk/psphot/src/psphotFitSet.c
r12792 r13035 1 1 # include "psphotInternal.h" 2 2 3 // XXX this is not used in main psphot code3 // This is not used in main psphot code (only in psphotModelTest.c) 4 4 bool psphotFitSet (pmSource *source, pmModel *oneModel, char *fitset, pmSourceFitMode mode) { 5 5 … … 24 24 } 25 25 26 // XXX pmSourceFitSet must cache the modelFlux? 26 27 pmSourceFitSet (source, modelSet, mode); 27 28 … … 32 33 for (int i = 0; i < modelSet->n; i++) { 33 34 pmModel *model = modelSet->data[i]; 34 pmModelSub (source->pixels, source->mask , model, false, false);35 pmModelSub (source->pixels, source->maskObj, model, PM_MODEL_OP_FULL); 35 36 36 37 fprintf (stderr, "output parameters (obj %d):\n", i); … … 42 43 // write out 43 44 psphotSaveImage (NULL, source->pixels, "resid.fits"); 44 psphotSaveImage (NULL, source->mask , "mask.fits");45 psphotSaveImage (NULL, source->maskObj, "mask.fits"); 45 46 return true; 46 47 } -
trunk/psphot/src/psphotGrowthCurve.c
r12792 r13035 53 53 radius = psf->growth->radius->data.F32[i]; 54 54 55 // NOTE: we use pmModelAdd not pmSourceAdd because we are not working with a normal source 55 56 psImageKeepCircle (mask, xc, yc, radius, "OR", PM_MASK_MARK); 56 pmModelAdd (image, mask, model, false, false);57 pmModelAdd (image, mask, model, PM_MODEL_OP_FULL); 57 58 pmSourcePhotometryAper (&apMag, model, image, mask); 58 59 psImageKeepCircle (mask, xc, yc, radius, "AND", PS_NOT_U8(PM_MASK_MARK)); -
trunk/psphot/src/psphotGuessModels.c
r12950 r13035 56 56 source->modelPSF = modelPSF; 57 57 source->modelPSF->residuals = psf->residuals; 58 59 pmSourceCacheModel (source); 58 60 } 59 61 psLogMsg ("psphot.models", 4, "built models for %ld objects: %f sec\n", sources->n, psTimerMark ("psphot")); -
trunk/psphot/src/psphotMakeResiduals.c
r12950 r13035 76 76 if (model == NULL) continue; // model must be defined 77 77 78 psImage *image = psImageCopy (NULL, source->pixels, PS_TYPE_F32);79 psImage *weight = psImageCopy (NULL, source->weight, PS_TYPE_F32);80 psImage *mask = psImageCopy (NULL, source->mask ,PS_TYPE_U8);81 pmModelSub (image, mask, model, false, false);78 psImage *image = psImageCopy (NULL, source->pixels, PS_TYPE_F32); 79 psImage *weight = psImageCopy (NULL, source->weight, PS_TYPE_F32); 80 psImage *mask = psImageCopy (NULL, source->maskView, PS_TYPE_U8); 81 pmModelSub (image, mask, model, PM_MODEL_OP_FUNC); 82 82 83 83 // re-normalize image and weight … … 253 253 # endif 254 254 255 psLogMsg ("psphot.pspsf", PS_LOG_INFO, "generate residuals for %ld objects: %f sec\n", input->n, psTimerMark ("residuals")); 256 255 257 psFree (xC); 256 258 psFree (yC); … … 272 274 psphotSaveImage (NULL, resid->mask, "resid.mk.fits"); 273 275 274 psLogMsg ("psphot.pspsf", PS_LOG_INFO, "generate residuals for %ld objects: %f sec\n", input->n, psTimerMark ("residuals"));275 276 276 psf->residuals = resid; 277 277 return true; -
trunk/psphot/src/psphotModelTest.c
r12792 r13035 2 2 static char DEFAULT_MODE[] = "EXT"; 3 3 4 // XXX consider this function : add more test information? 4 5 bool psphotModelTest (pmReadout *readout, psMetadata *recipe) { 5 6 … … 158 159 159 160 // define the pixels used for the fit 160 psImageKeepCircle (source->mask , xObj, yObj, RADIUS, "OR", PM_MASK_MARK);161 psImageKeepCircle (source->maskObj, xObj, yObj, RADIUS, "OR", PM_MASK_MARK); 161 162 162 163 char *fitset = psMetadataLookupStr (&status, recipe, "TEST_FIT_SET"); … … 170 171 // measure the source mags 171 172 pmSourcePhotometryModel (&fitMag, model); 172 pmSourcePhotometryAper (&obsMag, model, source->pixels, source->mask );173 pmSourcePhotometryAper (&obsMag, model, source->pixels, source->maskObj); 173 174 fprintf (stderr, "ap: %f, fit: %f, apmifit: %f\n", obsMag, fitMag, obsMag - fitMag); 174 175 … … 177 178 178 179 // subtract object, leave local sky 179 pmModelSub (source->pixels, source->mask , model, false, false);180 pmModelSub (source->pixels, source->maskObj, model, PM_MODEL_OP_FULL); 180 181 181 182 fprintf (stderr, "output parameters: \n"); … … 186 187 // write out 187 188 psphotSaveImage (NULL, source->pixels, "resid.fits"); 188 psphotSaveImage (NULL, source->mask , "mask.fits");189 psphotSaveImage (NULL, source->maskObj, "mask.fits"); 189 190 190 191 exit (0); -
trunk/psphot/src/psphotOutput.c
r13008 r13035 59 59 for (int j = 0; j < source->pixels->numCols; j++) { 60 60 // skip masked points 61 if (source->mask ->data.U8[i][j]) {61 if (source->maskObj->data.U8[i][j]) { 62 62 continue; 63 63 } … … 72 72 source->pixels->data.F32[i][j], 73 73 1.0 / source->weight->data.F32[i][j], 74 source->mask ->data.U8[i][j]);74 source->maskObj->data.U8[i][j]); 75 75 } 76 76 } … … 104 104 for (int i = 0; (sources != NULL) && (i < sources->n); i++) { 105 105 pmSource *source = (pmSource *) sources->data[i]; 106 pmModel *model = pmSource SelectModel (source);106 pmModel *model = pmSourceGetModel (NULL, source); 107 107 if (model == NULL) 108 108 continue; … … 154 154 return header; 155 155 } 156 157 # if (0)158 // output functions: we have several fixed modes (sx, obj, cmp)159 void psphotOutput (pmReadout *readout, psMetadata *arguments) {160 161 bool status;162 char *outputFile = NULL;163 164 psMetadata *header = pmReadoutGetHeader (readout);165 166 pmPSF *psf = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF");167 // sample PSF images??168 if (psfSample != NULL) psphotSamplePSFs (psf, readout->image, psfSample);169 if (psfSample != NULL) psphotSamplePSFs (psf, readout->image, psfSample);170 171 if (psfFile != NULL) {172 psMetadata *psfData = pmPSFtoMD (NULL, psf);173 psMetadataConfigWrite (psfData, psfFile);174 psFree (psfData);175 }176 }177 178 psImage *pmModelPSFatXY (psImage *image, pmModel *modelEXT, pmPSF *psf, int x, int y, int dx, int dy) {179 180 psRegion region = {x - dx, x + dx, y - dy, y + dy};181 psImage *view = psImageSubset (image, region);182 psImage *sample = psImageCopy (NULL, view, PS_TYPE_F32);183 psImageInit (sample, 0);184 modelEXT->params->data.F32[2] = x;185 modelEXT->params->data.F32[3] = y;186 pmModel *modelPSF = pmModelFromPSF (modelEXT, psf);187 pmModelAdd (sample, NULL, modelPSF, false, false);188 psFree (modelPSF);189 return (sample);190 }191 192 bool psphotSamplePSFs (pmPSF *psf, psImage *image, char *output) {193 194 // make sample PSFs for 4 corners and the center195 psImage *sample;196 197 // optional dump of all rough source data198 if (output[0] == 0) return false;199 200 pmModel *modelEXT = pmModelAlloc (psf->type);201 modelEXT->params->data.F32[0] = 0;202 modelEXT->params->data.F32[1] = 1;203 204 psFits *fits = psFitsOpen (output, "w");205 206 // the centers are in parent coordinates; they do not need to correspond to valid pixels...207 sample = pmModelPSFatXY (image, modelEXT, psf, 25, 25, 25, 25);208 psFitsWriteImage (fits, NULL, sample, 0);209 sample = pmModelPSFatXY (image, modelEXT, psf, image->numCols - 25, image->numRows - 25, 25, 25);210 psFitsWriteImage (fits, NULL, sample, 0);211 sample = pmModelPSFatXY (image, modelEXT, psf, image->numCols - 25, 25, 25, 25);212 psFitsWriteImage (fits, NULL, sample, 0);213 sample = pmModelPSFatXY (image, modelEXT, psf, 25, image->numRows - 25, 25, 25);214 psFitsWriteImage (fits, NULL, sample, 0);215 sample = pmModelPSFatXY (image, modelEXT, psf, image->numCols / 2, image->numRows / 2, 25, 25);216 psFitsWriteImage (fits, NULL, sample, 0);217 218 psFitsClose (fits);219 return (TRUE);220 }221 # endif -
trunk/psphot/src/psphotRadialPlot.c
r12900 r13035 60 60 for (int iy = 0; iy < source->pixels->numRows; iy++) { 61 61 for (int ix = 0; ix < source->pixels->numCols; ix++) { 62 if (source->mask ->data.U8[iy][ix]) {62 if (source->maskObj->data.U8[iy][ix]) { 63 63 rb->data.F32[nb] = hypot (ix - Xo, iy - Yo) ; 64 64 fb->data.F32[nb] = log10(source->pixels->data.F32[iy][ix]); -
trunk/psphot/src/psphotRadiusChecks.c
r12792 r13035 25 25 bool psphotCheckRadiusPSF (pmReadout *readout, pmSource *source, pmModel *model) 26 26 { 27 psF32 *PAR = model->params->data.F32; 28 29 // XXX do we have a better value for the sky noise level? not really... 27 30 pmMoments *moments = source->moments; 28 // do we have a better value for the sky noise level?29 // not really...30 31 31 32 // set the fit radius based on the object flux limit and the model … … 46 47 } 47 48 48 bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radiusFit); 49 bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit); 50 51 // set the mask to flag the excluded pixels 52 psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", PM_MASK_MARK); 49 53 return status; 50 54 } 51 55 52 56 bool psphotCheckRadiusPSFBlend (pmReadout *readout, pmSource *source, pmModel *model, float dR) { 57 58 psF32 *PAR = model->params->data.F32; 53 59 54 60 pmMoments *moments = source->moments; … … 63 69 } 64 70 65 bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radiusFit); 71 bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit); 72 73 // set the mask to flag the excluded pixels 74 psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", PM_MASK_MARK); 66 75 return status; 67 76 } … … 86 95 bool psphotCheckRadiusEXT (pmReadout *readout, pmSource *source, pmModel *model) { 87 96 97 psF32 *PAR = model->params->data.F32; 98 88 99 pmMoments *moments = source->moments; 89 100 if (moments == NULL) return false; … … 94 105 95 106 // redefine the pixels if needed 96 bool status = pmSourceRedefinePixels (source, readout, model->params->data.F32[2], model->params->data.F32[3], model->radiusFit); 107 bool status = pmSourceRedefinePixels (source, readout, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit); 108 109 // set the mask to flag the excluded pixels 110 psImageKeepCircle (source->maskObj, PAR[PM_PAR_XPOS], PAR[PM_PAR_YPOS], model->radiusFit, "OR", PM_MASK_MARK); 97 111 return status; 98 112 } -
trunk/psphot/src/psphotReadout.c
r13012 r13035 126 126 127 127 // linear PSF fit to peaks 128 psphot EnsemblePSF(readout, sources, recipe, psf, FALSE);128 psphotFitSourcesLinear (readout, sources, recipe, psf, FALSE); 129 129 if (dump) psphotSaveImage (NULL, readout->image, "image.v1.fits"); 130 130 if (!strcasecmp (breakPt, "ENSEMBLE")) { … … 141 141 142 142 // linear PSF fit to remaining peaks 143 psphot EnsemblePSF(readout, sources, recipe, psf, TRUE);143 psphotFitSourcesLinear (readout, sources, recipe, psf, TRUE); 144 144 if (dump) psphotSaveImage (NULL, readout->image, "image.v4.fits"); 145 145 if (!strcasecmp (breakPt, "PASS1")) { … … 192 192 193 193 // linear PSF fit to remaining peaks 194 psphot EnsemblePSF(readout, sources, recipe, psf, TRUE);194 psphotFitSourcesLinear (readout, sources, recipe, psf, TRUE); 195 195 if (dump) psphotSaveImage (NULL, readout->image, "image.v6.fits"); 196 196 -
trunk/psphot/src/psphotReplaceUnfit.c
r12792 r13035 1 1 # include "psphotInternal.h" 2 2 3 // replace the flux for sources which failed 3 4 bool psphotReplaceUnfit (psArray *sources) { 4 5 … … 15 16 16 17 replace: 17 if (source->modelPSF == NULL) continue; 18 19 psTrace ("psphot", 3, "replacing object at %f,%f\n", 20 source->modelPSF->params->data.F32[PM_PAR_XPOS], source->modelPSF->params->data.F32[PM_PAR_YPOS]); 21 22 pmModelAdd (source->pixels, source->mask, source->modelPSF, false, false); 18 pmSourceAdd (source, PM_MODEL_OP_FULL); 23 19 source->mode &= ~PM_SOURCE_MODE_SUBTRACTED; 24 20 source->mode &= ~PM_SOURCE_MODE_TEMPSUB; … … 40 36 if (!(source->mode & PM_SOURCE_MODE_SUBTRACTED)) continue; 41 37 42 // select appropriate model 43 pmModel *model = pmSourceGetModel (NULL, source); 44 if (model == NULL) continue; // model must be defined 45 46 psTrace ("psphot", 3, "replacing object at %f,%f\n", 47 model->params->data.F32[PM_PAR_XPOS], model->params->data.F32[PM_PAR_YPOS]); 48 49 pmModelAdd (source->pixels, source->mask, model, false, false); 38 pmSourceAdd (source, PM_MODEL_OP_FULL); 50 39 source->mode &= ~PM_SOURCE_MODE_SUBTRACTED; 51 40 } … … 54 43 } 55 44 56 // add or sub source replace or if the source has45 // add source, if the source has been subtracted (or if we ignore the state) 57 46 bool psphotAddWithTest (pmSource *source, bool useState) { 58 47 … … 60 49 bool state = !(source->mode & PM_SOURCE_MODE_SUBTRACTED); 61 50 if (state && useState) return true; 62 63 // select appropriate model64 pmModel *model = pmSourceGetModel (NULL, source);65 if (model == NULL) return false; // model must be defined66 67 psTrace ("psphot", 3, "replacing object at %f,%f\n",68 model->params->data.F32[PM_PAR_XPOS], model->params->data.F32[PM_PAR_YPOS]);69 51 70 52 // replace the model if 1) state says it is missing or 2) useState is false (just do it) 71 53 if (!state || !useState) { 72 pm ModelAdd (source->pixels, source->mask, model, false, false);54 pmSourceAdd (source, PM_MODEL_OP_FULL); 73 55 } 74 56 return true; 75 57 } 76 58 77 // add or sub source replace or if the source has59 // sub source, if the source has been added (or if we ignore the state) 78 60 bool psphotSubWithTest (pmSource *source, bool useState) { 79 61 … … 82 64 if (state && useState) return true; 83 65 84 // select appropriate model85 pmModel *model = pmSourceGetModel (NULL, source);86 if (model == NULL) return false; // model must be defined87 88 psTrace ("psphot", 3, "replacing object at %f,%f\n",89 model->params->data.F32[PM_PAR_XPOS], model->params->data.F32[PM_PAR_YPOS]);90 91 66 // replace the model if 1) state says it is missing or 2) useState is false (just do it) 92 67 if (!state || !useState) { 93 pm ModelSub (source->pixels, source->mask, model, false, false);68 pmSourceSub (source, PM_MODEL_OP_FULL); 94 69 } 95 70 return true; … … 102 77 bool newState = !(source->mode & PM_SOURCE_MODE_SUBTRACTED); 103 78 if (curState == newState) return true; 104 105 // select appropriate model106 pmModel *model = pmSourceGetModel (NULL, source);107 if (model == NULL) return false; // model must be defined108 109 psTrace ("psphot", 3, "replacing object at %f,%f\n",110 model->params->data.F32[PM_PAR_XPOS], model->params->data.F32[PM_PAR_YPOS]);111 79 112 80 if (curState && !newState) { 113 pm ModelSub (source->pixels, source->mask, model, false, false);81 pmSourceSub (source, PM_MODEL_OP_FULL); 114 82 } 115 83 if (newState && !curState) { 116 pm ModelAdd (source->pixels, source->mask, model, false, false);84 pmSourceAdd (source, PM_MODEL_OP_FULL); 117 85 } 118 86 return true; -
trunk/psphot/src/psphotSourceFits.c
r12792 r13035 67 67 // these should never be invalid values 68 68 // XXX drop these tests eventually 69 if (isnan(model->params->data.F32[PM_PAR_I0])) psAbort("nan in blend fit");69 if (isnan(model->params->data.F32[PM_PAR_I0])) psAbort("nan in blend fit"); 70 70 if (isnan(model->params->data.F32[PM_PAR_XPOS])) psAbort("nan in blend fit"); 71 71 if (isnan(model->params->data.F32[PM_PAR_YPOS])) psAbort("nan in blend fit"); … … 83 83 84 84 // fit PSF model (set/unset the pixel mask) 85 psImageKeepCircle (source->mask, x, y, PSF->radiusFit, "OR", PM_MASK_MARK);86 85 pmSourceFitSet (source, modelSet, PM_SOURCE_FIT_PSF); 87 psImageKeepCircle (source->mask, x, y, PSF->radiusFit, "AND", PS_NOT_U8(PM_MASK_MARK));88 86 89 87 // correct model chisq for flux trend … … 110 108 blend->modelPSF = psMemIncrRefCounter (model); 111 109 psTrace ("psphot", 5, "fitted blend as PSF\n"); 112 pmModelSub (source->pixels, source->mask, model, false, false); 110 111 // build cached model and subtract 112 pmSourceCacheModel (blend); 113 pmSourceSub (blend, PM_MODEL_OP_FULL); 113 114 blend->mode |= PM_SOURCE_MODE_SUBTRACTED; 114 115 blend->mode &= ~PM_SOURCE_MODE_TEMPSUB; … … 127 128 psFree (sourceSet); 128 129 129 psTrace ("psphot", 5, "fitted primary as PSF\n"); 130 pmModelSub (source->pixels, source->mask, PSF, false, false); 130 // save the new, successful model 131 131 psFree (source->modelPSF); 132 132 source->modelPSF = PSF; 133 psTrace ("psphot", 5, "fitted primary as PSF\n"); 134 135 // build cached model and subtract 136 pmSourceCacheModel (source); 137 pmSourceSub (source, PM_MODEL_OP_FULL); 133 138 source->mode |= PM_SOURCE_MODE_SUBTRACTED; 134 139 source->mode &= ~PM_SOURCE_MODE_TEMPSUB; … … 138 143 bool psphotFitPSF (pmReadout *readout, pmSource *source, pmPSF *psf) { 139 144 140 float x, y;141 145 double chiTrend; 142 146 … … 150 154 psphotCheckRadiusPSF (readout, source, PSF); 151 155 152 x = PSF->params->data.F32[PM_PAR_XPOS];153 y = PSF->params->data.F32[PM_PAR_YPOS];154 155 156 // fit PSF model (set/unset the pixel mask) 156 psImageKeepCircle (source->mask, x, y, PSF->radiusFit, "OR", PM_MASK_MARK);157 157 pmSourceFitModel (source, PSF, PM_SOURCE_FIT_PSF); 158 psImageKeepCircle (source->mask, x, y, PSF->radiusFit, "AND", PS_NOT_U8(PM_MASK_MARK));159 158 160 159 // correct model chisq for flux trend … … 168 167 } 169 168 170 psTrace ("psphot", 5, "fitted as PSF\n");171 pmModelSub (source->pixels, source->mask, PSF, false, false);172 173 169 // free old model, save new model 174 170 psFree (source->modelPSF); 175 171 source->modelPSF = PSF; 176 172 psTrace ("psphot", 5, "fitted as PSF\n"); 173 174 // build cached model and subtract 175 pmSourceCacheModel (source); 176 pmSourceSub (source, PM_MODEL_OP_FULL); 177 177 source->mode |= PM_SOURCE_MODE_SUBTRACTED; 178 178 source->mode &= ~PM_SOURCE_MODE_TEMPSUB; … … 264 264 // sub EXT 265 265 psFree (DBL); 266 pmModelSub (source->pixels, source->mask, EXT, false, false);267 psTrace ("psphot", 5, "blob as EXT: %f %f\n", EXT->params->data.F32[PM_PAR_XPOS], EXT->params->data.F32[PM_PAR_YPOS]);268 266 269 267 // save new model 270 268 source->modelEXT = EXT; 269 270 // build cached model and subtract 271 pmSourceCacheModel (source); 272 pmSourceSub (source, PM_MODEL_OP_FULL); 273 psTrace ("psphot", 5, "blob as EXT: %f %f\n", EXT->params->data.F32[PM_PAR_XPOS], EXT->params->data.F32[PM_PAR_YPOS]); 274 271 275 source->mode |= PM_SOURCE_MODE_SUBTRACTED; 272 276 source->mode &= ~PM_SOURCE_MODE_TEMPSUB; … … 276 280 // sub DLB 277 281 psFree (EXT); 278 pmModelSub (source->pixels, source->mask, (pmModel *) DBL->data[0], false, false);279 pmModelSub (source->pixels, source->mask, (pmModel *) DBL->data[1], false, false);280 psTrace ("psphot", 5, "blob as DBL: %f %f\n", ONE->params->data.F32[PM_PAR_XPOS], ONE->params->data.F32[PM_PAR_YPOS]);281 282 282 283 // drop old model, save new second model... … … 288 289 289 290 // copy most data from the primary source (modelEXT, blends stay NULL) 290 pmSource *newSrc = pmSourceAlloc (); 291 newSrc->peak = psMemIncrRefCounter (source->peak); 292 newSrc->pixels = psMemIncrRefCounter (source->pixels); 293 newSrc->weight = psMemIncrRefCounter (source->weight); 294 newSrc->mask = psMemIncrRefCounter (source->mask); 295 newSrc->moments = psMemIncrRefCounter (source->moments); 291 // XXX use pmSourceCopy? 292 pmSource *newSrc = pmSourceCopy (source); 296 293 newSrc->modelPSF = psMemIncrRefCounter (DBL->data[1]); 297 newSrc->type = source->type; 298 newSrc->mode = source->mode; 294 295 // build cached models and subtract 296 pmSourceCacheModel (source); 297 pmSourceSub (source, PM_MODEL_OP_FULL); 298 pmSourceCacheModel (newSrc); 299 pmSourceSub (newSrc, PM_MODEL_OP_FULL); 300 psTrace ("psphot", 5, "blob as DBL: %f %f\n", ONE->params->data.F32[PM_PAR_XPOS], ONE->params->data.F32[PM_PAR_YPOS]); 301 299 302 psArrayAdd (sources, 100, newSrc); 300 303 psFree (newSrc); … … 306 309 psArray *psphotFitDBL (pmReadout *readout, pmSource *source) { 307 310 308 float x, y,dx, dy;311 float dx, dy; 309 312 pmModel *DBL; 310 313 pmModel *PSF; … … 344 347 modelSet->data[1] = DBL; 345 348 346 x = source->moments->x;347 y = source->moments->y;348 349 349 // fit PSF model (set/unset the pixel mask) 350 psImageKeepCircle (source->mask, x, y, PSF->radiusFit, "OR", PM_MASK_MARK);351 350 pmSourceFitSet (source, modelSet, PM_SOURCE_FIT_PSF); 352 psImageKeepCircle (source->mask, x, y, PSF->radiusFit, "AND", PS_NOT_U8(PM_MASK_MARK));353 354 351 return (modelSet); 355 352 } 356 353 357 354 pmModel *psphotFitEXT (pmReadout *readout, pmSource *source) { 358 359 float x, y;360 355 361 356 NfitEXT ++; … … 369 364 psphotCheckRadiusEXT (readout, source, EXT); 370 365 371 x = EXT->params->data.F32[PM_PAR_XPOS];372 y = EXT->params->data.F32[PM_PAR_YPOS];373 374 366 if ((source->moments->Sx < 1e-3) || (source->moments->Sx < 1e-3)) { 375 367 psTrace ("psphot", 5, "problem source: moments: %f %f\n", source->moments->Sx, source->moments->Sy); … … 377 369 378 370 // fit EXT (not PSF) model (set/unset the pixel mask) 379 psImageKeepCircle (source->mask, x, y, EXT->radiusFit, "OR", PM_MASK_MARK);380 371 pmSourceFitModel (source, EXT, PM_SOURCE_FIT_EXT); 381 psImageKeepCircle (source->mask, x, y, EXT->radiusFit, "AND", PS_NOT_U8(PM_MASK_MARK));382 383 372 return (EXT); 384 373 } -
trunk/psphot/src/psphotTestPSF.c
r13006 r13035 124 124 125 125 // subtract model flux 126 pmModelSub (source->pixels, source->mask, model, false, false);126 pmModelSub (source->pixels, source->mask, model, PM_MODEL_OP_FULL); 127 127 } 128 128 fclose (f); -
trunk/psphot/src/psphotWeightBias.c
r13006 r13035 5 5 // only allow the normalization to vary 6 6 // XXX eventually, we should be able to do this with linear fitting... 7 // XXX UNUSED 7 8 bool psphotWeightBias (pmReadout *readout, psArray *sources, psMetadata *recipe, pmPSF *psf) { 8 9 … … 59 60 60 61 // replace object in image 61 pmModelAdd (source->pixels, source->mask, source->modelPSF, false, false);62 pmModelAdd (source->pixels, source->mask, source->modelPSF, PM_MODEL_OP_FULL); 62 63 63 64 // make a temporary model (we don't keep the result of this analysis) … … 76 77 77 78 // re-subtract PSF for object, leave local sky 78 pmModelSub (source->pixels, source->mask, source->modelPSF, false, false);79 pmModelSub (source->pixels, source->mask, source->modelPSF, PM_MODEL_OP_FULL); 79 80 80 81 PARp = source->modelPSF->params->data.F32;
Note:
See TracChangeset
for help on using the changeset viewer.
