Changeset 26897 for trunk/psastro/src/psastroMaskUpdates.c
- Timestamp:
- Feb 10, 2010, 7:39:43 PM (16 years ago)
- File:
-
- 1 edited
-
trunk/psastro/src/psastroMaskUpdates.c (modified) (12 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/psastro/src/psastroMaskUpdates.c
r26259 r26897 1 1 /** @file psastroMaskUpdates.c 2 2 * 3 * @brief 3 * @brief 4 4 * 5 5 * @ingroup libpsastro … … 35 35 psImageMaskType starMaskValue = pmConfigMaskGet("STARCORE", config); // Mask value for ghost pixels 36 36 psImageMaskType crosstalkMaskValue = pmConfigMaskGet("GHOST", config); // Mask value for crosstalk ghosts 37 37 38 38 // psImageMaskType maskBlank = pmConfigMaskGet("BLANK", config); // Mask value for blank pixels 39 39 … … 55 55 return(false); 56 56 } 57 57 58 58 // convert star positions to ghost positions and add to the readout->analysis data 59 59 if (!psastroLoadGhosts (config)) { … … 81 81 double REFSTAR_MASK_SATSPIKE_MAG_MAX = psMetadataLookupF32 (&status, recipe, "REFSTAR_MASK_SATSPIKE_MAG_MAX"); 82 82 double REFSTAR_MASK_SATSPIKE_WIDTH = psMetadataLookupF32 (&status, recipe, "REFSTAR_MASK_SATSPIKE_WIDTH"); 83 double REFSTAR_MASK_SATSPIKE_OFFSET = psMetadataLookupF32 (&status, recipe, "REFSTAR_MASK_SATSPIKE_OFFSET"); 83 84 double REFSTAR_MASK_BLEED_MAG_MAX = psMetadataLookupF32 (&status, recipe, "REFSTAR_MASK_BLEED_MAG_MAX"); 84 85 double REFSTAR_MASK_BLEED_MAG_SLOPE = psMetadataLookupF32 (&status, recipe, "REFSTAR_MASK_BLEED_MAG_SLOPE"); … … 86 87 double REFSTAR_MASK_CROSSTALK_MAG_MAX = psMetadataLookupF32 (&status, recipe, "REFSTAR_MASK_CROSSTALK_MAG_MAX"); 87 88 double REFSTAR_MASK_CROSSTALK_MAG_SLOPE= psMetadataLookupF32 (&status, recipe, "REFSTAR_MASK_CROSSTALK_MAG_SLOPE"); 88 89 89 90 // select the input data sources 90 91 pmFPAfile *input = psMetadataLookupPtr (NULL, config->files, "PSASTRO.INPUT"); … … 188 189 if (! readout->data_exists) { continue; } 189 190 190 // XXX why not do this? 191 pmReadout *readoutMask = pmFPAviewThisReadout (view, outMask->fpa);192 if (!readoutMask) continue;191 // XXX why not do this? 192 pmReadout *readoutMask = pmFPAviewThisReadout (view, outMask->fpa); 193 if (!readoutMask) continue; 193 194 194 195 // select the raw objects for this readout 195 // XXX : note that we place limits on the refstar sample in psastroChooseRefstars.c:196 // 1) on chip and 2) < PSASTRO.MAX.NREF. magnitude limits and clump exclusion are only 197 // applied to the SUBSETs196 // XXX : note that we place limits on the refstar sample in psastroChooseRefstars.c: 197 // 1) on chip and 2) < PSASTRO.MAX.NREF. magnitude limits and clump exclusion are only 198 // applied to the SUBSETs 198 199 psArray *refstars = psMetadataLookupPtr (&status, readout->analysis, "PSASTRO.REFSTARS"); 199 200 if (!refstars) continue; … … 203 204 // 2) diffraction spikes in direction ROT - ROTo 204 205 // 3) bleed trail in the direction of the readout 205 // update: with 'burntool' applied to the data, the bleed trail mask is not needed.206 // update: with 'burntool' applied to the data, the bleed trail mask is not needed. 206 207 207 208 for (int i = 0; i < refstars->n; i++) { 208 209 pmAstromObj *ref = refstars->data[i]; 209 210 if (COUNT_GHOSTS) {211 if (ref->Mag > GHOST_MAX_MAG) {212 nGhosts ++;213 }214 }210 211 if (COUNT_GHOSTS) { 212 if (ref->Mag > GHOST_MAX_MAG) { 213 nGhosts ++; 214 } 215 } 215 216 216 217 if (ref->Mag > REFSTAR_MASK_MAX_MAG) continue; 217 218 218 219 219 // the reference magnitudes have been converted from the instrumental220 // values supplied in the recipe to apparent mags (above)220 // the reference magnitudes have been converted from the instrumental 221 // values supplied in the recipe to apparent mags (above) 221 222 222 223 // CIRCLE around the stars (scaled by magnitude) … … 224 225 225 226 // XXX for now, assume cell binning is 1x1 relative to chip 226 psastroMaskCircle (readoutMask->mask, starMaskValue, ref->chip->x, ref->chip->y, radius, radius);227 psastroMaskCircle (readoutMask->mask, starMaskValue, ref->chip->x, ref->chip->y, radius, radius); 227 228 228 229 for (float theta = 0.0; theta < 2*M_PI; theta += M_PI / 2.0) { … … 230 231 float Theta = theta - ROTANGLE - REFSTAR_MASK_SATSTAR_POS_ZERO; 231 232 232 // LINE for boundaries of the saturation spikes (scaled by magnitude) 233 // float MAG_MAX = (theta == 0.0) || (theta == M_PI) ? REFSTAR_MASK_SATSPIKE_MAG_MAX1 : REFSTAR_MASK_SATSPIKE_MAG_MAX2; 234 235 float MAG_MAX = REFSTAR_MASK_SATSPIKE_MAG_MAX; 236 float spikeLength = REFSTAR_MASK_SATSPIKE_MAG_SLOPE * (MAG_MAX - ref->Mag); 237 float spikeWidth = 0.5*REFSTAR_MASK_SATSPIKE_WIDTH; 238 // XXX we can make the width depend on the spike as well... 239 // The length should also be a function of the image background level 240 241 psastroMaskBox (readoutMask->mask, spikeMaskValue, ref->chip->x, ref->chip->y, spikeLength, spikeWidth, Theta); 242 } 243 244 // This masking option was needed for persistent charge trails in GPC1; it 245 // has since been replaced with 'burntool', which is applied upon readout 246 // by the camera software, and therefore is aware of the image sequence. 247 if (REFSTAR_MASK_BLEED) { 248 // convert x,y chip coordinates to cells in maskChip 249 pmCell *refCell = pmCellInChip (refChip, ref->chip->x, ref->chip->y); 250 251 // LINE for boundaries of the bleed lines 252 if (refCell) { 253 float xCell = 0.0; 254 float yCell = 0.0; 255 float xEnd = 0.0; 256 float yEnd = 0.0; 257 // find coordinate of star on cell 258 pmCellCoordsForChip (&xCell, &yCell, refCell, ref->chip->x, ref->chip->y); 259 // find coordinate of end-point on chip 260 261 int ySize = psMetadataLookupS32(NULL, refCell->concepts, "CELL.YSIZE"); 262 pmChipCoordsForCell (&xEnd, &yEnd, refCell, xCell, ySize); 263 264 float width = REFSTAR_MASK_BLEED_MAG_SLOPE*(REFSTAR_MASK_BLEED_MAG_MAX - ref->Mag); 265 psastroMaskRectangle (readoutMask->mask, spikeMaskValue, (int) ref->chip->x-0.5*width, (int) ref->chip->y, (int) ref->chip->x+0.5*width+1, yEnd); 233 // LINE for boundaries of the saturation spikes (scaled by magnitude) 234 // float MAG_MAX = (theta == 0.0) || (theta == M_PI) ? REFSTAR_MASK_SATSPIKE_MAG_MAX1 : REFSTAR_MASK_SATSPIKE_MAG_MAX2; 235 236 float MAG_MAX = REFSTAR_MASK_SATSPIKE_MAG_MAX; 237 float spikeLength = pow(10,REFSTAR_MASK_SATSPIKE_MAG_SLOPE * (MAG_MAX - ref->Mag)) - REFSTAR_MASK_SATSPIKE_OFFSET; 238 float spikeWidth = 0.5*REFSTAR_MASK_SATSPIKE_WIDTH; 239 if (spikeLength < 0.0) { 240 spikeLength = 0.0; 266 241 } 267 } 242 // XXX we can make the width depend on the spike as well... 243 // The length should also be a function of the image background level 244 245 psastroMaskBox (readoutMask->mask, spikeMaskValue, ref->chip->x, ref->chip->y, spikeLength, spikeWidth, Theta); 246 } 247 248 // This masking option was needed for persistent charge trails in GPC1; it 249 // has since been replaced with 'burntool', which is applied upon readout 250 // by the camera software, and therefore is aware of the image sequence. 251 if (REFSTAR_MASK_BLEED) { 252 // convert x,y chip coordinates to cells in maskChip 253 pmCell *refCell = pmCellInChip (refChip, ref->chip->x, ref->chip->y); 254 255 // LINE for boundaries of the bleed lines 256 if (refCell) { 257 float xCell = 0.0; 258 float yCell = 0.0; 259 float xEnd = 0.0; 260 float yEnd = 0.0; 261 // find coordinate of star on cell 262 pmCellCoordsForChip (&xCell, &yCell, refCell, ref->chip->x, ref->chip->y); 263 // find coordinate of end-point on chip 264 265 int ySize = psMetadataLookupS32(NULL, refCell->concepts, "CELL.YSIZE"); 266 pmChipCoordsForCell (&xEnd, &yEnd, refCell, xCell, ySize); 267 268 float width = REFSTAR_MASK_BLEED_MAG_SLOPE*(REFSTAR_MASK_BLEED_MAG_MAX - ref->Mag); 269 psastroMaskRectangle (readoutMask->mask, spikeMaskValue, (int) ref->chip->x-0.5*width, (int) ref->chip->y, (int) ref->chip->x+0.5*width+1, yEnd); 270 } 271 } 268 272 } 269 273 … … 274 278 // supplied here. 275 279 psArray *ghosts = psMetadataLookupPtr (&status, readout->analysis, "PSASTRO.GHOSTS"); 276 if (ghosts) { 277 // mask the ghosts on this readout278 for (int i = 0; i < ghosts->n; i++) {279 psastroGhost *ghost = ghosts->data[i];280 // XXX bright vs faint ghost bits? (OR with SUSPECT) 281 psastroMaskEllipticalAnnulus (readoutMask->mask, ghostMaskValue, ghost->chip->x, ghost->chip->y, ghost->inner, ghost->outer);282 }283 }280 if (ghosts) { 281 // mask the ghosts on this readout 282 for (int i = 0; i < ghosts->n; i++) { 283 psastroGhost *ghost = ghosts->data[i]; 284 // XXX bright vs faint ghost bits? (OR with SUSPECT) 285 psastroMaskEllipticalAnnulus (readoutMask->mask, ghostMaskValue, ghost->chip->x, ghost->chip->y, ghost->inner, ghost->outer); 286 } 287 } 284 288 285 289 // Select the glint mask regions for this readout (loaded in … … 290 294 psArray *glints = psMetadataLookupPtr (&status, readout->analysis, "PSASTRO.GLINTS"); 291 295 if (glints) { 292 // mask the glints on this readout 293 for (int i = 0; i < glints->n; i++) { 294 psRegion *glint = glints->data[i]; 295 psastroMaskRectangle (readoutMask->mask, glintMaskValue, glint->x0, glint->y0, glint->x1, glint->y1); 296 } 297 } 298 299 // Select the crosstalk artifact regions for this readout, and mask a circular region 300 // corresponding to the source star 301 302 psArray *crosstalks = psMetadataLookupPtr (&status, readout->analysis, "PSASTRO.CROSSTALKS"); 303 if (crosstalks) { 304 for (int i = 0; i < crosstalks->n; i++) { 305 pmAstromObj *ref = crosstalks->data[i]; 306 float radius = REFSTAR_MASK_CROSSTALK_MAG_SLOPE * (REFSTAR_MASK_CROSSTALK_MAG_MAX - ref->Mag); 307 psTrace("psastro.crosstalk",2,"Masking star on Chip %s @ (%f,%f) Magnitude: %f Radius %f\n", 308 psMetadataLookupStr(&status,readout->parent->parent->concepts,"CHIP.NAME"), 309 ref->chip->x,ref->chip->y,ref->Mag,radius); 310 // XXX for now, assume cell binning is 1x1 relative to chip 311 psastroMaskCircle (readoutMask->mask, crosstalkMaskValue, ref->chip->x, ref->chip->y, radius, radius); 312 } 313 } 314 315 // this probably should move into a function of its own: 316 { 317 // select the raw objects for this readout, flag is they fall in a mask 318 psArray *inSources = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.SOURCES"); 319 if (inSources == NULL) continue; 320 321 // create a replacement output array: 322 // psArray *outSources = psAllocArrayEmpty(100); 323 324 // XXX finish this: raise a bit for stars that land on certain types of masks; 325 // others (eg, bright star core) should be ignored. 326 for (int i = 0; i < inSources->n; i++) { 327 pmSource *source = inSources->data[i]; 328 329 int xChip = source->peak->x; 330 int yChip = source->peak->y; 331 332 bool onChip = true; 333 onChip &= (xChip >= 0); 334 onChip &= (xChip < readoutMask->mask->numCols); 335 onChip &= (yChip >= 0); 336 onChip &= (yChip < readoutMask->mask->numRows); 337 if (!onChip) { 338 // if the source is off the edge of the chip, raise a different bit? 339 source->mode |= PM_SOURCE_MODE_OFF_CHIP; 340 continue; 341 } 342 343 psImageMaskType value = readoutMask->mask->data.PS_TYPE_IMAGE_MASK_DATA[yChip][xChip]; 344 if (value & ghostMaskValue) { 345 source->mode |= PM_SOURCE_MODE_ON_GHOST; 346 } 347 // XXX note that for now, glint and ghost are identical 348 pmSourceMode PM_SOURCE_MODE_ON_GLINT = PM_SOURCE_MODE_ON_GHOST; 349 if (value & glintMaskValue) { 350 source->mode |= PM_SOURCE_MODE_ON_GLINT; 351 } 352 if (value & spikeMaskValue) { 353 source->mode |= PM_SOURCE_MODE_ON_SPIKE; 354 } 355 } 356 } 296 // mask the glints on this readout 297 for (int i = 0; i < glints->n; i++) { 298 psRegion *glint = glints->data[i]; 299 psastroMaskRectangle (readoutMask->mask, glintMaskValue, glint->x0, glint->y0, glint->x1, glint->y1); 300 } 301 } 302 303 // Select the crosstalk artifact regions for this readout, and mask a circular region 304 // corresponding to the source star 305 306 psArray *crosstalks = psMetadataLookupPtr (&status, readout->analysis, "PSASTRO.CROSSTALKS"); 307 if (crosstalks) { 308 for (int i = 0; i < crosstalks->n; i++) { 309 pmAstromObj *ref = crosstalks->data[i]; 310 float radius = REFSTAR_MASK_CROSSTALK_MAG_SLOPE * (REFSTAR_MASK_CROSSTALK_MAG_MAX - ref->Mag); 311 psTrace("psastro.crosstalk",2,"Masking star on Chip %s @ (%f,%f) Magnitude: %f Radius %f\n", 312 psMetadataLookupStr(&status,readout->parent->parent->concepts,"CHIP.NAME"), 313 ref->chip->x,ref->chip->y,ref->Mag,radius); 314 // XXX for now, assume cell binning is 1x1 relative to chip 315 psastroMaskCircle (readoutMask->mask, crosstalkMaskValue, ref->chip->x, ref->chip->y, radius, radius); 316 } 317 } 318 319 // this probably should move into a function of its own: 320 { 321 // select the raw objects for this readout, flag is they fall in a mask 322 pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS"); 323 if (!detections) continue; 324 325 psArray *inSources = detections->allSources; 326 psAssert (inSources, "missing sources?"); 327 328 // create a replacement output array: 329 // psArray *outSources = psAllocArrayEmpty(100); 330 331 // XXX finish this: raise a bit for stars that land on certain types of masks; 332 // others (eg, bright star core) should be ignored. 333 for (int i = 0; i < inSources->n; i++) { 334 pmSource *source = inSources->data[i]; 335 336 int xChip = source->peak->x; 337 int yChip = source->peak->y; 338 339 bool onChip = true; 340 onChip &= (xChip >= 0); 341 onChip &= (xChip < readoutMask->mask->numCols); 342 onChip &= (yChip >= 0); 343 onChip &= (yChip < readoutMask->mask->numRows); 344 if (!onChip) { 345 // if the source is off the edge of the chip, raise a different bit? 346 source->mode |= PM_SOURCE_MODE_OFF_CHIP; 347 continue; 348 } 349 350 psImageMaskType value = readoutMask->mask->data.PS_TYPE_IMAGE_MASK_DATA[yChip][xChip]; 351 if (value & ghostMaskValue) { 352 source->mode |= PM_SOURCE_MODE_ON_GHOST; 353 } 354 // XXX note that for now, glint and ghost are identical 355 pmSourceMode PM_SOURCE_MODE_ON_GLINT = PM_SOURCE_MODE_ON_GHOST; 356 if (value & glintMaskValue) { 357 source->mode |= PM_SOURCE_MODE_ON_GLINT; 358 } 359 if (value & spikeMaskValue) { 360 source->mode |= PM_SOURCE_MODE_ON_SPIKE; 361 } 362 } 363 } 357 364 } 358 365 } … … 375 382 376 383 if (COUNT_GHOSTS) { 377 // save nGhosts to update header.378 psMetadata *updates = psMetadataLookupMetadata (&status, fpa->analysis, "PSASTRO.HEADER");379 if (!updates) {380 updates = psMetadataAlloc ();381 psMetadataAddMetadata (fpa->analysis, PS_LIST_TAIL, "PSASTRO.HEADER", PS_META_REPLACE, "psastro header stats", updates);382 psFree (updates);383 }384 psMetadataAddS32 (updates, PS_LIST_TAIL, "NGHOSTS", PS_META_REPLACE, "total expected ghosts", nGhosts);384 // save nGhosts to update header. 385 psMetadata *updates = psMetadataLookupMetadata (&status, fpa->analysis, "PSASTRO.HEADER"); 386 if (!updates) { 387 updates = psMetadataAlloc (); 388 psMetadataAddMetadata (fpa->analysis, PS_LIST_TAIL, "PSASTRO.HEADER", PS_META_REPLACE, "psastro header stats", updates); 389 psFree (updates); 390 } 391 psMetadataAddS32 (updates, PS_LIST_TAIL, "NGHOSTS", PS_META_REPLACE, "total expected ghosts", nGhosts); 385 392 } 386 393
Note:
See TracChangeset
for help on using the changeset viewer.
