- Timestamp:
- Jun 22, 2012, 3:33:32 PM (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/eam_branches/ipp-20120601/psphot/src/psphotMaskReadout.c
r29936 r34053 94 94 95 95 // test output of files at this stage 96 if (psTraceGetLevel("psphot ") >= 5) {96 if (psTraceGetLevel("psphot.imsave") >= 5) { 97 97 psphotSaveImage (NULL, readout->image, "image.fits"); 98 98 psphotSaveImage (NULL, readout->mask, "mask.fits"); … … 105 105 return true; 106 106 } 107 108 // XXX this function and support below was created to test the theory that the faint-end 109 // bias results from the Poisson variation of the background pixels. This is NOT the 110 // case. Using the code below maintains the faint-end bias. 111 bool psphotUpdateVariance (pmConfig *config, const pmFPAview *view, const char *filerule) { 112 113 bool status = false; 114 115 // select the appropriate recipe information 116 psMetadata *recipe = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE); 117 psAssert (recipe, "missing recipe?"); 118 119 int num = psphotFileruleCount(config, filerule); 120 121 // loop over the available readouts 122 for (int i = 0; i < num; i++) { 123 124 // Generate the mask and weight images, including the user-defined analysis region of interest 125 if (!psphotUpdateVarianceReadout (config, view, filerule, i, recipe)) { 126 psError (PSPHOT_ERR_CONFIG, false, "failed to generate mask for %s entry %d", filerule, i); 127 return false; 128 } 129 } 130 return true; 131 } 132 133 // determine the mean variance image (equivalent to the background model, but for the variance image) 134 // set the variance image to the MAX(input, mean) 135 bool psphotUpdateVarianceReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { 136 137 bool status; 138 139 pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest 140 psAssert (file, "missing file?"); 141 142 // find the currently selected readout 143 pmReadout *readout = pmFPAviewThisReadout (view, file->fpa); 144 psAssert (readout, "missing readout?"); 145 146 pmSourceFitVarMode varMode = psphotGetFitVarMode (recipe); 147 if (varMode == PM_SOURCE_PHOTFIT_NONE) { 148 psError (PSPHOT_ERR_CONFIG, false, "need valid LINEAR_FIT_VARIANCE_MODE"); 149 return false; 150 } 151 152 // make this an option via the recipe 153 if (varMode != PM_SOURCE_PHOTFIT_MODEL_SKY) return true; 154 155 // create a model variance image (full-scale image to take result of psImageUnbin below) 156 psImage *modelVar = psImageCopy (NULL, readout->variance, PS_TYPE_F32); 157 158 // find the binning information 159 psImageBinning *backBinning = psphotBackgroundBinning (modelVar, config); 160 assert (backBinning); 161 162 psImage *varModel = psImageAlloc(backBinning->nXruff, backBinning->nYruff, PS_TYPE_F32); // Background model 163 psImage *varModelStdev = psImageAlloc(backBinning->nXruff, backBinning->nYruff, PS_TYPE_F32); // Background model 164 165 if (!psphotModelBackgroundReadout(varModel, varModelStdev, NULL, readout, backBinning, config, true)) { 166 psError(PS_ERR_UNKNOWN, false, "Unable to generate background model"); 167 psFree (varModel); 168 psFree (varModelStdev); 169 return false; 170 } 171 172 // linear interpolation to full-scale 173 if (!psImageUnbin (modelVar, varModel, backBinning)) { 174 psError (PSPHOT_ERR_PROG, true, "inconsistent sizes for unbinning"); 175 psFree (varModel); 176 psFree (varModelStdev); 177 return false; 178 } 179 180 // XXX save these? 181 psFree (varModel); 182 psFree (varModelStdev); 183 184 psImage *im = readout->image; 185 psImage *wt = readout->variance; 186 for (int j = 0; j < im->numRows; j++) { 187 for (int i = 0; i < im->numCols; i++) { 188 if (!isfinite(im->data.F32[j][i])) continue; 189 if (!isfinite(wt->data.F32[j][i])) continue; 190 // XXX for a test, make variance constant wt->data.F32[j][i] = PS_MAX(wt->data.F32[j][i], modelVar->data.F32[j][i]); 191 wt->data.F32[j][i] = modelVar->data.F32[j][i]; 192 } 193 } 194 195 // test output of files at this stage 196 if (psTraceGetLevel("psphot.imsave") >= 5) { 197 psphotSaveImage (NULL, readout->image, "image.varsky.fits"); 198 psphotSaveImage (NULL, readout->mask, "mask.varsky.fits"); 199 psphotSaveImage (NULL, readout->variance, "variance.varsky.fits"); 200 } 201 202 psFree (modelVar); 203 204 return true; 205 }
Note:
See TracChangeset
for help on using the changeset viewer.
