Changeset 21365
- Timestamp:
- Feb 5, 2009, 4:47:00 PM (17 years ago)
- Location:
- trunk
- Files:
-
- 16 edited
-
ppMerge/src/ppMergeArguments.c (modified) (4 diffs)
-
ppMerge/src/ppMergeCamera.c (modified) (6 diffs)
-
ppMerge/src/ppMergeFiles.c (modified) (8 diffs)
-
ppMerge/src/ppMergeLoop_Threaded.c (modified) (8 diffs)
-
ppMerge/src/ppMergeMask.c (modified) (13 diffs)
-
ppMerge/src/ppMergeReadChunk.c (modified) (4 diffs)
-
ppMerge/src/ppMergeScaleZero.c (modified) (2 diffs)
-
ppSim/src/ppSimInsertStars.c (modified) (4 diffs)
-
ppSim/src/ppSimLoadForceSources.c (modified) (5 diffs)
-
ppSim/src/ppSimLoop.c (modified) (8 diffs)
-
ppSim/src/ppSimMakeBias.c (modified) (3 diffs)
-
ppSim/src/ppSimMakeDark.c (modified) (2 diffs)
-
ppSim/src/ppSimMakeSky.c (modified) (5 diffs)
-
ppSim/src/ppSimMergeReadouts.c (modified) (3 diffs)
-
ppSim/src/ppSimPhotomReadoutFake.c (modified) (7 diffs)
-
ppSim/src/ppSimPhotomReadoutForce.c (modified) (7 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/ppMerge/src/ppMergeArguments.c
r21244 r21365 6 6 * 7 7 * @author IfA 8 * @version $Revision: 1.1 6$ $Name: not supported by cvs2svn $9 * @date $Date: 2009-02-0 1 21:43:05$8 * @version $Revision: 1.17 $ $Name: not supported by cvs2svn $ 9 * @date $Date: 2009-02-06 02:44:31 $ 10 10 * Copyright 2009 Institute for Astronomy, University of Hawaii 11 11 */ … … 25 25 "\tIMAGE(STR): Image filename\n" 26 26 "\tMASK(STR): Mask filename\n" 27 "\t WEIGHT(STR) Weightmap filename\n"28 "where MASK and WEIGHTare optional.",27 "\tVARIANCE(STR) Variance map filename\n" 28 "where MASK and VARIANCE are optional.", 29 29 program); 30 30 fprintf(stderr, "\n"); … … 155 155 psMetadataAddF32(arguments, PS_LIST_TAIL, "-frachigh", 0, "Fraction of low pixels to discard", NAN); 156 156 psMetadataAddS32(arguments, PS_LIST_TAIL, "-nkeep", 0, "Minimum number of pixels in stack to keep", 0); 157 psMetadataAddBool(arguments, PS_LIST_TAIL, "- weights", 0, "Use image weights in combination?", false);157 psMetadataAddBool(arguments, PS_LIST_TAIL, "-variances", 0, "Use image variances in combination?", false); 158 158 159 159 // XXX EAM : not clear this should be allowed on the command line. … … 300 300 VALUE_ARG_RECIPE_FLOAT("-frachigh", "FRACHIGH", F32); 301 301 VALUE_ARG_RECIPE_INT("-nkeep", "NKEEP", S32, 0); 302 VALUE_ARG_RECIPE_BOOL("- weights", "WEIGHTS");302 VALUE_ARG_RECIPE_BOOL("-variances", "VARIANCES"); 303 303 304 304 // XXX we do not supply this on the command line -
trunk/ppMerge/src/ppMergeCamera.c
r21244 r21365 6 6 * 7 7 * @author IfA 8 * @version $Revision: 1. 8$ $Name: not supported by cvs2svn $9 * @date $Date: 2009-02-0 1 21:43:05$8 * @version $Revision: 1.9 $ $Name: not supported by cvs2svn $ 9 * @date $Date: 2009-02-06 02:44:31 $ 10 10 * Copyright 2009 Institute for Astronomy, University of Hawaii 11 11 */ … … 96 96 { 97 97 bool haveMasks = false; // Do we have masks? 98 bool have Weights = false; // Do we have weightmaps?98 bool haveVariances = false; // Do we have variance maps? 99 99 100 100 ppMergeType type = psMetadataLookupS32(NULL, config->arguments, "TYPE"); ///< Type of frame … … 121 121 122 122 bool mdok; 123 psString mask = psMetadataLookupStr(&mdok, input, "MASK"); // /<Name of mask124 psString weight = psMetadataLookupStr(&mdok, input, "WEIGHT"); ///< Name of weightmap123 psString mask = psMetadataLookupStr(&mdok, input, "MASK"); // Name of mask 124 psString variance = psMetadataLookupStr(&mdok, input, "VARIANCE"); // Name of variance map 125 125 126 126 // Add the image file … … 164 164 } 165 165 166 // Optionally add the weightfile167 if ( weight && strlen(weight) > 0) {168 have Weights = true;169 psArray * weightFiles = psArrayAlloc(1); ///<Array of filenames for this FPA170 weightFiles->data[0] = psMemIncrRefCounter(weight);171 psMetadataAddArray(config->arguments, PS_LIST_TAIL, " WEIGHT.FILENAMES", PS_META_REPLACE,172 "Filenames of weight files", weightFiles);173 psFree( weightFiles);166 // Optionally add the variance file 167 if (variance && strlen(variance) > 0) { 168 haveVariances = true; 169 psArray *varianceFiles = psArrayAlloc(1); // Array of filenames for this FPA 170 varianceFiles->data[0] = psMemIncrRefCounter(variance); 171 psMetadataAddArray(config->arguments, PS_LIST_TAIL, "VARIANCE.FILENAMES", PS_META_REPLACE, 172 "Filenames of variance files", varianceFiles); 173 psFree(varianceFiles); 174 174 175 175 bool status; 176 pmFPAfile * weightFile = pmFPAfileBindFromArgs(&status, imageFile, config, "PPMERGE.INPUT.WEIGHT",177 "WEIGHT.FILENAMES");176 pmFPAfile *varianceFile = pmFPAfileBindFromArgs(&status, imageFile, config, 177 "PPMERGE.INPUT.VARIANCE", "VARIANCE.FILENAMES"); 178 178 if (!status) { 179 psError(PS_ERR_UNKNOWN, false, "Unable to define file from weight %d (%s)", numFiles, weight); 180 return false; 181 } 182 if (weightFile->type != PM_FPA_FILE_WEIGHT) { 183 psError(PS_ERR_IO, true, "PPMERGE.INPUT.WEIGHT is not of type WEIGHT"); 184 return false; 185 } 186 haveWeights = true; 179 psError(PS_ERR_UNKNOWN, false, "Unable to define file from variance %d (%s)", 180 numFiles, variance); 181 return false; 182 } 183 if (varianceFile->type != PM_FPA_FILE_VARIANCE) { 184 psError(PS_ERR_IO, true, "PPMERGE.INPUT.VARIANCE is not of type VARIANCE"); 185 return false; 186 } 187 haveVariances = true; 187 188 } 188 189 … … 194 195 psMetadataRemoveKey(config->arguments, "MASK.FILENAMES"); 195 196 } 196 if (psMetadataLookup(config->arguments, " WEIGHT.FILENAMES")) {197 psMetadataRemoveKey(config->arguments, " WEIGHT.FILENAMES");197 if (psMetadataLookup(config->arguments, "VARIANCE.FILENAMES")) { 198 psMetadataRemoveKey(config->arguments, "VARIANCE.FILENAMES"); 198 199 } 199 200 if (psMetadataLookup(config->arguments, "PSF.FILENAMES")) { … … 203 204 psMetadataAddS32(config->arguments, PS_LIST_TAIL, "INPUTS.NUM", 0, "Number of input files", numFiles); 204 205 psMetadataAddBool(config->arguments, PS_LIST_TAIL, "INPUTS.MASKS", 0, "Got input masks?", haveMasks); 205 psMetadataAddBool(config->arguments, PS_LIST_TAIL, "INPUTS. WEIGHTS", 0,206 "Got input weights?", haveWeights);206 psMetadataAddBool(config->arguments, PS_LIST_TAIL, "INPUTS.VARIANCES", 0, 207 "Got input variances?", haveVariances); 207 208 208 209 -
trunk/ppMerge/src/ppMergeFiles.c
r21244 r21365 6 6 * 7 7 * @author IfA 8 * @version $Revision: 1. 7$ $Name: not supported by cvs2svn $9 * @date $Date: 2009-02-0 1 21:43:05$8 * @version $Revision: 1.8 $ $Name: not supported by cvs2svn $ 9 * @date $Date: 2009-02-06 02:44:31 $ 10 10 * Copyright 2009 Institute for Astronomy, University of Hawaii 11 11 */ … … 13 13 #include "ppMerge.h" 14 14 15 const char *allFiles[] = { "PPMERGE.INPUT", "PPMERGE.INPUT.MASK", "PPMERGE.INPUT. WEIGHT",15 const char *allFiles[] = { "PPMERGE.INPUT", "PPMERGE.INPUT.MASK", "PPMERGE.INPUT.VARIANCE", 16 16 "PPMERGE.OUTPUT", "PPMERGE.OUTPUT.COUNT", "PPMERGE.OUTPUT.SIGMA", 17 NULL }; // /<All files18 const char *inputFiles[] = { "PPMERGE.INPUT", "PPMERGE.INPUT.MASK", "PPMERGE.INPUT. WEIGHT",19 NULL }; // /<Input files17 NULL }; // All files 18 const char *inputFiles[] = { "PPMERGE.INPUT", "PPMERGE.INPUT.MASK", "PPMERGE.INPUT.VARIANCE", 19 NULL }; // Input files 20 20 const char *outputFiles[] = { "PPMERGE.OUTPUT", "PPMERGE.OUTPUT.COUNT", "PPMERGE.OUTPUT.SIGMA", 21 21 NULL }; ///< Output files … … 51 51 } 52 52 } 53 if (psMetadataLookupBool(&mdok, config->arguments, "INPUTS. WEIGHTS")) {54 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPMERGE.INPUT. WEIGHT", num);55 if (!pmReadoutReadChunk Weight(readout, file->fits, 0, rows, 0, config)) {56 psError(PS_ERR_UNKNOWN, false, "Unable to read readout weight.");53 if (psMetadataLookupBool(&mdok, config->arguments, "INPUTS.VARIANCES")) { 54 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPMERGE.INPUT.VARIANCE", num); 55 if (!pmReadoutReadChunkVariance(readout, file->fits, 0, rows, 0, config)) { 56 psError(PS_ERR_UNKNOWN, false, "Unable to read readout variance."); 57 57 return false; 58 58 } … … 91 91 psFree(fileView); 92 92 } 93 if (psMetadataLookupBool(&mdok, config->arguments, "INPUTS.WEIGHTS")) { 94 pmFPAfile *weight = pmFPAfileSelectSingle(config->files, "PPMERGE.INPUT.WEIGHT", num); // Weight file 95 pmFPAview *fileView = pmFPAviewForLevel(weight->fileLevel, view); 96 if (!pmFPAfileOpen(weight, fileView, config)) { 97 psError(PS_ERR_UNKNOWN, false, "Unable to open weight file %d", num); 93 if (psMetadataLookupBool(&mdok, config->arguments, "INPUTS.VARIANCES")) { 94 pmFPAfile *variance = pmFPAfileSelectSingle(config->files, "PPMERGE.INPUT.VARIANCE", 95 num); // Variance file 96 pmFPAview *fileView = pmFPAviewForLevel(variance->fileLevel, view); 97 if (!pmFPAfileOpen(variance, fileView, config)) { 98 psError(PS_ERR_UNKNOWN, false, "Unable to open variance file %d", num); 98 99 psFree(fileView); 99 100 return false; … … 127 128 bool mdok; // Status of MD lookup 128 129 bool haveMasks = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.MASKS"); // Do we have masks? 129 bool have Weights = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.WEIGHTS"); // Got weights?130 bool haveVariances = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.VARIANCES"); // Got variances? 130 131 131 132 const char **fileList = selectFiles(files); // Files to activate … … 134 135 continue; 135 136 } 136 if (!have Weights && strcmp(fileList[i], "PPMERGE.INPUT.WEIGHT") == 0) {137 if (!haveVariances && strcmp(fileList[i], "PPMERGE.INPUT.VARIANCE") == 0) { 137 138 continue; 138 139 } … … 163 164 bool mdok; // Status of MD lookup 164 165 bool haveMasks = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.MASKS"); // Do we have masks? 165 bool have Weights = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.WEIGHTS"); // Got weights?166 bool haveVariances = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.VARIANCES"); // Got variances? 166 167 167 168 psList *list = psListAlloc(NULL); // List of files … … 171 172 continue; 172 173 } 173 if (!have Weights && strcmp(fileList[i], "PPMERGE.INPUT.WEIGHT") == 0) {174 if (!haveVariances && strcmp(fileList[i], "PPMERGE.INPUT.VARIANCE") == 0) { 174 175 continue; 175 176 } -
trunk/ppMerge/src/ppMergeLoop_Threaded.c
r21244 r21365 6 6 * 7 7 * @author IfA 8 * @version $Revision: 1. 8$ $Name: not supported by cvs2svn $9 * @date $Date: 2009-02-0 1 21:43:05$8 * @version $Revision: 1.9 $ $Name: not supported by cvs2svn $ 9 * @date $Date: 2009-02-06 02:44:31 $ 10 10 * Copyright 2009 Institute for Astronomy, University of Hawaii 11 11 */ … … 36 36 bool mdok; ///< Status of MD lookup 37 37 bool haveMasks = psMetadataLookupBool(&mdok, arguments, "INPUTS.MASKS"); // Do we have masks? 38 bool have Weights = psMetadataLookupBool(&mdok, arguments, "INPUTS.WEIGHTS"); ///< Do we have weights?39 40 psArray *inputs = ppMergeFileDataLevel(config, "PPMERGE.INPUT", PM_FPA_LEVEL_READOUT); // /<Input images41 psArray *masks = NULL, * weights = NULL; ///< Input masks and weights38 bool haveVariances = psMetadataLookupBool(&mdok, arguments, "INPUTS.VARIANCES"); // Do we have variances? 39 40 psArray *inputs = ppMergeFileDataLevel(config, "PPMERGE.INPUT", PM_FPA_LEVEL_READOUT); // Input images 41 psArray *masks = NULL, *variances = NULL; // Input masks and variances 42 42 if (haveMasks) { 43 43 masks = ppMergeFileDataLevel(config, "PPMERGE.INPUT.MASK", PM_FPA_LEVEL_READOUT); 44 44 } 45 if (have Weights) {46 weights = ppMergeFileDataLevel(config, "PPMERGE.INPUT.WEIGHT", PM_FPA_LEVEL_READOUT);45 if (haveVariances) { 46 variances = ppMergeFileDataLevel(config, "PPMERGE.INPUT.VARIANCE", PM_FPA_LEVEL_READOUT); 47 47 } 48 48 … … 50 50 if (!mdok) nThreads = 0; 51 51 52 / ** General combination parameters */53 int iter = psMetadataLookupS32(NULL, arguments, "ITER"); // /<Number of rejection iterations54 float rej = psMetadataLookupF32(NULL, arguments, "REJ"); // /<Rejection level55 float fraclow = psMetadataLookupF32(NULL, arguments, "FRACLOW"); // /<Reject fraction of low pixels56 float frachigh = psMetadataLookupF32(NULL, arguments, "FRACHIGH"); // /<Reject fraction of hi pixels57 int nKeep = psMetadataLookupS32(NULL, arguments, "NKEEP"); // /<Minimum number of values to keep58 psStatsOptions combineStat = psMetadataLookupS32(NULL, arguments, "COMBINE"); // /<Combination statistic59 bool use Weights = psMetadataLookupBool(NULL, arguments, "WEIGHTS"); ///< Use weights?60 61 / ** Fringe parameters */62 int fringeNum = psMetadataLookupS32(NULL, arguments, "FRINGE.NUM"); // /<Number of fringe points63 int fringeSize = psMetadataLookupS32(NULL, arguments, "FRINGE.SIZE"); // /<Size of fringe regions64 int fringeSmoothX = psMetadataLookupS32(NULL, arguments, "FRINGE.XSMOOTH"); // /<Smoothing regions in x65 int fringeSmoothY = psMetadataLookupS32(NULL, arguments, "FRINGE.YSMOOTH"); // /<Smoothing regions in y52 // General combination parameters 53 int iter = psMetadataLookupS32(NULL, arguments, "ITER"); // Number of rejection iterations 54 float rej = psMetadataLookupF32(NULL, arguments, "REJ"); // Rejection level 55 float fraclow = psMetadataLookupF32(NULL, arguments, "FRACLOW"); // Reject fraction of low pixels 56 float frachigh = psMetadataLookupF32(NULL, arguments, "FRACHIGH"); // Reject fraction of hi pixels 57 int nKeep = psMetadataLookupS32(NULL, arguments, "NKEEP"); // Minimum number of values to keep 58 psStatsOptions combineStat = psMetadataLookupS32(NULL, arguments, "COMBINE"); // Combination statistic 59 bool useVariances = psMetadataLookupBool(NULL, arguments, "VARIANCES"); // Use variances? 60 61 // Fringe parameters 62 int fringeNum = psMetadataLookupS32(NULL, arguments, "FRINGE.NUM"); // Number of fringe points 63 int fringeSize = psMetadataLookupS32(NULL, arguments, "FRINGE.SIZE"); // Size of fringe regions 64 int fringeSmoothX = psMetadataLookupS32(NULL, arguments, "FRINGE.XSMOOTH"); // Smoothing regions in x 65 int fringeSmoothY = psMetadataLookupS32(NULL, arguments, "FRINGE.YSMOOTH"); // Smoothing regions in y 66 66 67 67 // set the mask and mark bit values based on the named masks 68 psImageMaskType maskVal; 69 psImageMaskType markVal; 70 if (!pmConfigMaskSetBits (&maskVal, &markVal, config)) { 68 psImageMaskType maskVal, markVal; 69 if (!pmConfigMaskSetBits(&maskVal, &markVal, config)) { 71 70 psError (PS_ERR_UNKNOWN, true, "Unable to define the mask bit values"); 72 71 return false; … … 81 80 combination->iter = iter; 82 81 combination->rej = rej; 83 combination-> weights = useWeights;82 combination->variances = useVariances; 84 83 85 84 psMetadata *stats = NULL; ///< Statistics for output … … 344 343 psFree(countsRO); 345 344 346 // XXX EAM 2009.01.18 : sigmaCell and countsCell need to have their concepts copied from outCell.347 // This was causing segfaults for VYSOS5; Why did this ever work for SIMTEST?348 if (!pmConceptsCopyCell(countsCell, outCell)) {349 psError(PS_ERR_UNKNOWN, false, "Unable to copy cell concepts.");350 psFree(outRO);351 goto ERROR;352 }345 // XXX EAM 2009.01.18 : sigmaCell and countsCell need to have their concepts copied from outCell. 346 // This was causing segfaults for VYSOS5; Why did this ever work for SIMTEST? 347 if (!pmConceptsCopyCell(countsCell, outCell)) { 348 psError(PS_ERR_UNKNOWN, false, "Unable to copy cell concepts."); 349 psFree(outRO); 350 goto ERROR; 351 } 353 352 354 353 pmCell *sigmaCell = pmFPAfileThisCell(config->files, view, "PPMERGE.OUTPUT.SIGMA"); … … 360 359 psFree(sigmaRO); 361 360 362 if (!pmConceptsCopyCell(sigmaCell, outCell)) {363 psError(PS_ERR_UNKNOWN, false, "Unable to copy cell concepts.");364 psFree(outRO);365 goto ERROR;366 }361 if (!pmConceptsCopyCell(sigmaCell, outCell)) { 362 psError(PS_ERR_UNKNOWN, false, "Unable to copy cell concepts."); 363 psFree(outRO); 364 goto ERROR; 365 } 367 366 } 368 367 … … 447 446 psFree(inputs); 448 447 psFree(masks); 449 psFree( weights);448 psFree(variances); 450 449 psFree(stats); 451 450 return true; … … 456 455 psFree(inputs); 457 456 psFree(masks); 458 psFree( weights);457 psFree(variances); 459 458 psFree(stats); 460 459 return false; -
trunk/ppMerge/src/ppMergeMask.c
r21244 r21365 6 6 * 7 7 * @author IfA 8 * @version $Revision: 1.1 7$ $Name: not supported by cvs2svn $9 * @date $Date: 2009-02-0 1 21:43:05$8 * @version $Revision: 1.18 $ $Name: not supported by cvs2svn $ 9 * @date $Date: 2009-02-06 02:44:31 $ 10 10 * Copyright 2009 Institute for Astronomy, University of Hawaii 11 11 */ … … 38 38 float smoothScale = psMetadataLookupF32(&mdok, config->arguments, "MASK.SMOOTH.SCALE"); ///< Radius to grow mask 39 39 40 psImageMaskType markVal; 41 psImageMaskType maskValRaw; 42 if (!pmConfigMaskSetBits (&maskValRaw, &markVal, config)) { 43 psError (PS_ERR_UNKNOWN, true, "Unable to define the mask bit values"); 44 return false; 40 psImageMaskType markVal, maskValRaw; 41 if (!pmConfigMaskSetBits(&maskValRaw, &markVal, config)) { 42 psError(PS_ERR_UNKNOWN, true, "Unable to define the mask bit values"); 43 return false; 45 44 } 46 45 … … 48 47 psImageMaskType maskValOut = pmConfigMaskGet (maskOutName, config); 49 48 if (!maskValOut) { 50 psError (PS_ERR_UNKNOWN, true, "Undefined output mask bit value");51 return false;52 } 53 54 psStats *statistics = psStatsAlloc(meanStat | stdevStat); // /<Statistics for background49 psError (PS_ERR_UNKNOWN, true, "Undefined output mask bit value"); 50 return false; 51 } 52 53 psStats *statistics = psStatsAlloc(meanStat | stdevStat); // Statistics for background 55 54 56 55 psString outName = ppMergeOutputFile(config); ///< Name of output file 57 56 pmChip *outChip = pmFPAfileThisChip(config->files, view, outName); ///< Output chip 58 57 psFree(outName); 59 58 60 59 int numCells = 1; 61 60 if (chipStats) { 62 // count the number of active cells for this chip:63 numCells = 0;64 for (int i = 0; i < outChip->cells->n; i++) {65 pmCell *cell = outChip->cells->data[i];66 if (!cell->process) continue;67 numCells ++;68 }61 // count the number of active cells for this chip: 62 numCells = 0; 63 for (int i = 0; i < outChip->cells->n; i++) { 64 pmCell *cell = outChip->cells->data[i]; 65 if (!cell->process) continue; 66 numCells ++; 67 } 69 68 } 70 69 … … 85 84 while ((inCell = pmFPAviewNextCell(inView, inFPA, 1))) { 86 85 87 // the output FPA structure carries the information about which cells to process88 pmCell *outCell = pmFPAfileThisCell(config->files, inView, "PPMERGE.OUTPUT.MASK"); // /<Output cell89 if (!outCell->process) continue;86 // the output FPA structure carries the information about which cells to process 87 pmCell *outCell = pmFPAfileThisCell(config->files, inView, "PPMERGE.OUTPUT.MASK"); // Output cell 88 if (!outCell->process) continue; 90 89 91 90 pmHDU *hdu = pmHDUFromCell(inCell); ///< HDU for cell … … 150 149 int y = pixel / numCols; 151 150 if (mask && (mask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskValRaw)) continue; 152 if (outMask && (outMask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskValOut)) continue;153 if (!isfinite(image->data.F32[y][x])) continue;151 if (outMask && (outMask->data.PS_TYPE_IMAGE_MASK_DATA[y][x] & maskValOut)) continue; 152 if (!isfinite(image->data.F32[y][x])) continue; 154 153 155 154 values->data.F32[valueIndex++] = image->data.F32[y][x]; … … 165 164 } 166 165 167 float mean = psStatsGetValue(statistics, meanStat);168 float stdev = psStatsGetValue(statistics, stdevStat);169 170 // this function increments the count for each suspect pixel in each input plane171 // maskValRaw is used to test for valid input pixels166 float mean = psStatsGetValue(statistics, meanStat); 167 float stdev = psStatsGetValue(statistics, stdevStat); 168 169 // this function increments the count for each suspect pixel in each input plane 170 // maskValRaw is used to test for valid input pixels 172 171 if (!pmMaskFlagSuspectPixels(outRO, readout, mean, stdev, maskSuspect, maskValRaw)) { 173 172 psError(PS_ERR_UNKNOWN, false, "Unable to find suspect values in file %d", i); … … 197 196 inView->cell = -1; 198 197 while ((inCell = pmFPAviewNextCell(inView, inFPA, 1))) { 199 200 // the output FPA structure carries the information about which cells to process 201 pmCell *outCell = pmFPAfileThisCell(config->files, inView, "PPMERGE.OUTPUT.MASK"); ///< Output cell 202 if (!outCell->process) continue; 198 // the output FPA structure carries the information about which cells to process 199 pmCell *outCell = pmFPAfileThisCell(config->files, inView, "PPMERGE.OUTPUT.MASK"); // Output cell 200 if (!outCell->process) continue; 203 201 204 202 pmHDU *hdu = pmHDUFromCell(inCell); ///< HDU for cell … … 210 208 pmReadout *outRO = pmFPAfileThisReadout(config->files, inView, "PPMERGE.OUTPUT.MASK"); 211 209 212 float mean = psStatsGetValue(statistics, meanStat);213 float stdev = psStatsGetValue(statistics, stdevStat);214 215 if (!pmMaskFlagSuspectPixels(outRO, readout, mean, stdev, maskSuspect, maskValRaw)) {210 float mean = psStatsGetValue(statistics, meanStat); 211 float stdev = psStatsGetValue(statistics, stdevStat); 212 213 if (!pmMaskFlagSuspectPixels(outRO, readout, mean, stdev, maskSuspect, maskValRaw)) { 216 214 psError(PS_ERR_UNKNOWN, false, "Unable to find suspect values in file %d", i); 217 215 goto MERGE_MASK_ERROR; … … 243 241 while ((outCell = pmFPAviewNextCell(outView, outFPA, 1))) { 244 242 245 // skip inactive cells246 if (!outCell->process) continue;243 // skip inactive cells 244 if (!outCell->process) continue; 247 245 248 246 pmHDU *hdu = pmHDUFromCell(outCell); ///< HDU for cell … … 254 252 pmReadout *outRO = outCell->readouts->data[0]; ///< Output readout 255 253 256 if (smoothSuspect) {257 // XXX test output of suspect pixel image258 psImage *suspects = psMetadataLookupPtr(NULL, outRO->analysis, PM_MASK_ANALYSIS_SUSPECT); ///<Suspect img259 assert (suspects);260 psImageSmooth (suspects, smoothScale, 3); ///<extend smoothing region to 3-sigma261 }262 263 // set the bad pixels to the value 'maskVal'254 if (smoothSuspect) { 255 // XXX test output of suspect pixel image 256 psImage *suspects = psMetadataLookupPtr(NULL, outRO->analysis, PM_MASK_ANALYSIS_SUSPECT); // Suspect img 257 assert (suspects); 258 psImageSmooth (suspects, smoothScale, 3); // extend smoothing region to 3-sigma 259 } 260 261 // set the bad pixels to the value 'maskVal' 264 262 if (!pmMaskIdentifyBadPixels(outRO, maskValOut, maskBad, maskMode)) { 265 263 psError(PS_ERR_UNKNOWN, false, "Unable to mask bad pixels"); … … 352 350 assert(config); 353 351 354 bool mdok; // /<Status of MD lookup355 int numFiles = psMetadataLookupS32(NULL, config->arguments, "INPUTS.NUM"); // /<Number of inputs356 bool haveMasks = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.MASKS"); // /<Do we have masks?357 bool have Weights = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.WEIGHTS"); ///< Do we have weights?358 int iter = psMetadataLookupS32(NULL, config->arguments, "ITER"); // /<Number of rejection iterations352 bool mdok; // Status of MD lookup 353 int numFiles = psMetadataLookupS32(NULL, config->arguments, "INPUTS.NUM"); // Number of inputs 354 bool haveMasks = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.MASKS"); // Do we have masks? 355 bool haveVariances = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.VARIANCES"); // Got variances? 356 int iter = psMetadataLookupS32(NULL, config->arguments, "ITER"); // Number of rejection iterations 359 357 360 358 PS_ASSERT_INT_POSITIVE(iter, false); … … 375 373 psFree(masks); 376 374 } 377 if (have Weights) {378 psArray * weights = ppMergeFileDataLevel(config, "PPMERGE.INPUT.WEIGHT", PM_FPA_LEVEL_READOUT);379 psFree( weights);375 if (haveVariances) { 376 psArray *variances = ppMergeFileDataLevel(config, "PPMERGE.INPUT.VARIANCE", PM_FPA_LEVEL_READOUT); 377 psFree(variances); 380 378 } 381 379 … … 417 415 } 418 416 419 if (outChip->data_exists) {420 psList *inChips = psListAlloc(NULL);421 for (int i=0; i < numFiles; i++) {422 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPMERGE.INPUT", i); ///<Input file423 pmChip *chip = pmFPAviewThisChip(view, file->fpa);424 psListAdd(inChips, PS_LIST_TAIL, chip);425 }426 427 // XXX I need to call pmConfigMaskWriteHeader for the PHU somewhere, after it is created!428 429 if (!pmConceptsAverageChips(outChip, inChips, true)) {430 psError(PS_ERR_UNKNOWN, false, "Unable to average Chip concepts.");431 psFree(inChips);432 goto PPMERGE_MASK_ERROR;433 }434 psFree(inChips);435 }417 if (outChip->data_exists) { 418 psList *inChips = psListAlloc(NULL); 419 for (int i=0; i < numFiles; i++) { 420 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPMERGE.INPUT", i); // Input file 421 pmChip *chip = pmFPAviewThisChip(view, file->fpa); 422 psListAdd(inChips, PS_LIST_TAIL, chip); 423 } 424 425 // XXX I need to call pmConfigMaskWriteHeader for the PHU somewhere, after it is created! 426 427 if (!pmConceptsAverageChips(outChip, inChips, true)) { 428 psError(PS_ERR_UNKNOWN, false, "Unable to average Chip concepts."); 429 psFree(inChips); 430 goto PPMERGE_MASK_ERROR; 431 } 432 psFree(inChips); 433 } 436 434 if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { 437 435 goto PPMERGE_MASK_ERROR; -
trunk/ppMerge/src/ppMergeReadChunk.c
r21244 r21365 6 6 * 7 7 * @author IfA 8 * @version $Revision: 1. 5$ $Name: not supported by cvs2svn $9 * @date $Date: 2009-02-0 1 21:43:05$8 * @version $Revision: 1.6 $ $Name: not supported by cvs2svn $ 9 * @date $Date: 2009-02-06 02:44:31 $ 10 10 * Copyright 2009 Institute for Astronomy, University of Hawaii 11 11 */ … … 20 20 21 21 bool mdok; 22 bool haveMasks = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.MASKS"); ///< Do we have masks? 23 bool haveWeights = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.WEIGHTS");///< Do we have weights? 24 int rows = psMetadataLookupS32(NULL, config->arguments, "ROWS"); ///< Number of rows to read per chunk 22 bool haveMasks = psMetadataLookupBool(&mdok, config->arguments, "INPUTS.MASKS"); // Do we have masks? 23 bool haveVariances = psMetadataLookupBool(&mdok, config->arguments, 24 "INPUTS.VARIANCES");// Do we have variances? 25 int rows = psMetadataLookupS32(NULL, config->arguments, "ROWS"); // Number of rows to read per chunk 25 26 26 27 // select an available fileGroup … … 52 53 // override the recorded last scan 53 54 inRO->thisImageScan = fileGroup->firstScan; 54 inRO->this WeightScan = fileGroup->firstScan;55 inRO->thisVarianceScan = fileGroup->firstScan; 55 56 inRO->thisMaskScan = fileGroup->firstScan; 56 57 inRO->forceScan = true; … … 81 82 } 82 83 83 if (have Weights && pmReadoutMoreWeight(inRO, file->fits, 0, rows, config)) {84 if (haveVariances && pmReadoutMoreVariance(inRO, file->fits, 0, rows, config)) { 84 85 keepReading = true; 85 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPMERGE.INPUT.WEIGHT", i); 86 if (!pmReadoutReadChunkWeight(inRO, file->fits, 0, rows, 0, config)) { 87 psError(PS_ERR_IO, false, "Unable to read chunk %d for file PPMERGE.INPUT.WEIGHT %d", 86 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPMERGE.INPUT.VARIANCE", i); 87 if (!pmReadoutReadChunkVariance(inRO, file->fits, 0, rows, 0, config)) { 88 psError(PS_ERR_IO, false, 89 "Unable to read chunk %d for file PPMERGE.INPUT.VARIANCE %d", 88 90 numChunk, i); 89 91 *status = false; -
trunk/ppMerge/src/ppMergeScaleZero.c
r21265 r21365 6 6 * 7 7 * @author IfA 8 * @version $Revision: 1.3 0$ $Name: not supported by cvs2svn $9 * @date $Date: 2009-02-0 2 20:26:26$8 * @version $Revision: 1.31 $ $Name: not supported by cvs2svn $ 9 * @date $Date: 2009-02-06 02:44:31 $ 10 10 * Copyright 2009 Institute for Astronomy, University of Hawaii 11 11 */ … … 56 56 for (int i = 0; i < numInputs; i++) { 57 57 pmFPAfileActivate(config->files, false, NULL); 58 psArray *files = ppMergeFileActivateSingle(config, PPMERGE_FILES_INPUT, true, i); // /<Activated files59 pmFPAfile *input = files->data[0]; // /< Representative file; should be the image (not mask or weight)60 pmFPA *fpa = input->fpa; // /<FPA of interest61 view = pmFPAviewAlloc(0); // /<View to component of interest62 int cellNum = 0; // /<Index for cell58 psArray *files = ppMergeFileActivateSingle(config, PPMERGE_FILES_INPUT, true, i); // Activated files 59 pmFPAfile *input = files->data[0]; // Representative file; should be the image (not mask or variance) 60 pmFPA *fpa = input->fpa; // FPA of interest 61 view = pmFPAviewAlloc(0); // View to component of interest 62 int cellNum = 0; // Index for cell 63 63 if (!pmFPAfileIOChecks(config, view, PM_FPA_BEFORE)) { 64 64 goto ERROR; -
trunk/ppSim/src/ppSimInsertStars.c
r20367 r21365 39 39 float skyRate = psMetadataLookupF32(NULL, recipe, "SKY.RATE"); // Sky rate 40 40 if (isnan(skyRate)) { 41 float zp = psMetadataLookupF32(&mdok, recipe, "ZEROPOINT"); assert (mdok);42 float scale = psMetadataLookupF32(&mdok, recipe, "PIXEL.SCALE"); assert (mdok);43 float skyMags = psMetadataLookupF32(&mdok, recipe, "SKY.MAGS"); assert (mdok);41 float zp = psMetadataLookupF32(&mdok, recipe, "ZEROPOINT"); assert (mdok); 42 float scale = psMetadataLookupF32(&mdok, recipe, "PIXEL.SCALE"); assert (mdok); 43 float skyMags = psMetadataLookupF32(&mdok, recipe, "SKY.MAGS"); assert (mdok); 44 44 skyRate = scale * scale * ppSimMagToFlux (skyMags, zp); 45 45 } … … 75 75 FILE *outfile = fopen (outname, "w"); 76 76 77 // add sources to the readout image & weight77 // add sources to the readout image & variance 78 78 psTrace("ppSim", 1, "Inserting %ld stars...\n", stars->n); 79 79 for (long i = 0; i < stars->n; i++) { … … 118 118 psF64 Area = 2.0 * M_PI * axes.major * axes.minor; 119 119 120 // this value is the pure (input) flux, and is saved in the output source cmf files120 // this value is the pure (input) flux, and is saved in the output source cmf files 121 121 source->psfMag = -2.5*log10(star->flux); 122 122 source->errMag = sqrt(Area*PS_SQR(roughNoise) + flux) / flux; … … 132 132 // Blow away the image parts of the source, which makes the memory explode 133 133 RESET(source->pixels); 134 RESET(source-> weight);134 RESET(source->variance); 135 135 RESET(source->maskObj); 136 136 RESET(source->maskView); -
trunk/ppSim/src/ppSimLoadForceSources.c
r18011 r21365 56 56 // Select the spots within range of this readout. Project the spots to this chip 57 57 for (int i = 0; i < spots->n; i++) { 58 pmAstromObj *obj = spots->data[i];58 pmAstromObj *obj = spots->data[i]; 59 59 60 psProject (obj->TP, obj->sky, fpa->toSky);61 psPlaneTransformApply (obj->FP, fpa->fromTPA, obj->TP);62 psPlaneTransformApply (obj->chip, chip->fromFPA, obj->FP);60 psProject (obj->TP, obj->sky, fpa->toSky); 61 psPlaneTransformApply (obj->FP, fpa->fromTPA, obj->TP); 62 psPlaneTransformApply (obj->chip, chip->fromFPA, obj->FP); 63 63 64 // limit the X,Y range of the objs to the selected chip65 if (obj->chip->x < minX) continue;66 if (obj->chip->x > maxX) continue;67 if (obj->chip->y < minY) continue;68 if (obj->chip->y > maxY) continue;64 // limit the X,Y range of the objs to the selected chip 65 if (obj->chip->x < minX) continue; 66 if (obj->chip->x > maxX) continue; 67 if (obj->chip->y < minY) continue; 68 if (obj->chip->y > maxY) continue; 69 69 70 // convert the pmAstromObj to pmSource70 // convert the pmAstromObj to pmSource 71 71 72 72 // instantiate a model for the PSF at this location, Io = 1.0 … … 74 74 75 75 // XXX let the flux limit be a user-defined number of sky sigmas (not just 1.0) 76 // XXX use a fixed radius??76 // XXX use a fixed radius?? 77 77 // float radius = model->modelRadius (model->params, roughNoise); 78 78 // radius = PS_MAX (radius, 1.0); 79 float radius = 5.0;79 float radius = 5.0; 80 80 81 81 // construct a source, with model flux pixels set, based on the model … … 87 87 // psF64 Area = 2.0 * M_PI * axes.major * axes.minor; 88 88 89 // these are not really needed since we will be fitting for them89 // these are not really needed since we will be fitting for them 90 90 source->psfMag = NAN; 91 91 source->errMag = NAN; 92 92 93 93 // insert the source flux in the image 94 // XXX not sure the offset is really 0,094 // XXX not sure the offset is really 0,0 95 95 pmSourceAddWithOffset (source, PM_MODEL_OP_FULL, 0xff, 0.0, 0.0); 96 96 pmSourceAddWithOffset (source, PM_MODEL_OP_FULL | PM_MODEL_OP_NOISE, 0xff, 0.0, 0.0); … … 98 98 // Blow away the image parts of the source, which makes the memory explode 99 99 RESET(source->pixels); 100 RESET(source-> weight);100 RESET(source->variance); 101 101 RESET(source->maskObj); 102 102 RESET(source->maskView); … … 105 105 RESET(source->blends); 106 106 107 psArrayAdd (sources, 100, source);107 psArrayAdd (sources, 100, source); 108 108 psFree(source); // Drop local reference 109 109 } -
trunk/ppSim/src/ppSimLoop.c
r18011 r21365 34 34 psArray *galaxies = psArrayAllocEmpty (1); 35 35 if (type == PPSIM_TYPE_OBJECT) { 36 // Load forced-photometry positions (these are placed on fpa->analysis for use in ppSimPhotomReadout)37 if (!ppSimLoadSpots (fpa, config)) ESCAPE (PS_ERR_UNKNOWN, "failed to load forced-photometry spots");38 39 // Load catalogue stars36 // Load forced-photometry positions (these are placed on fpa->analysis for use in ppSimPhotomReadout) 37 if (!ppSimLoadSpots (fpa, config)) ESCAPE (PS_ERR_UNKNOWN, "failed to load forced-photometry spots"); 38 39 // Load catalogue stars 40 40 if (!ppSimLoadStars (stars, fpa, config)) ESCAPE (PS_ERR_UNKNOWN, "failed to load catalog stars"); 41 41 42 // Add random stars42 // Add random stars 43 43 if (!ppSimMakeStars (stars, fpa, config, rng)) ESCAPE (PS_ERR_UNKNOWN, "failed to make random stars"); 44 44 45 // Add random galaxies45 // Add random galaxies 46 46 if (!ppSimMakeGalaxies (galaxies, fpa, config, rng)) ESCAPE (PS_ERR_UNKNOWN, "failed to make random galaxies"); 47 47 } … … 105 105 // TO DO: Decide if cell is to be windowed, reduce numCols, numRows appropriately 106 106 readout->image = psImageAlloc(numCols, numRows, PS_TYPE_F32); // Signal in pixels 107 readout-> weight= psImageAlloc(numCols, numRows, PS_TYPE_F32); // Noise in pixels107 readout->variance = psImageAlloc(numCols, numRows, PS_TYPE_F32); // Noise in pixels 108 108 109 109 psImageInit (readout->image, 0.0); 110 psImageInit (readout-> weight, 0.0);110 psImageInit (readout->variance, 0.0); 111 111 112 112 psFree(readout); // Drop reference … … 115 115 psVector *biasCols = ppSimMakeBiassec (cell, config); 116 116 117 pmReadout *readout;118 while ((readout = pmFPAviewNextReadout (view, fpa, 1))) {119 120 // if we have not read in a weightor generated a fake image above, we need to117 pmReadout *readout; 118 while ((readout = pmFPAviewNextReadout (view, fpa, 1))) { 119 120 // if we have not read in a variance or generated a fake image above, we need to 121 121 // build one here 122 if (!readout-> weight) {123 if (!pmReadoutGenerate Weight(readout, true)) {124 psError (PS_ERR_UNKNOWN, false, "trouble creating weight");122 if (!readout->variance) { 123 if (!pmReadoutGenerateVariance(readout, true)) { 124 psError (PS_ERR_UNKNOWN, false, "trouble creating variance"); 125 125 return false; 126 126 } … … 133 133 134 134 psVector *biasRows = ppSimMakeBias (&status, readout, config, rng); 135 if (!status) ESCAPE (PS_ERR_UNKNOWN, "problem generating dark structure");135 if (!status) ESCAPE (PS_ERR_UNKNOWN, "problem generating dark structure"); 136 136 if (type == PPSIM_TYPE_BIAS) goto done; 137 137 138 if (!ppSimMakeDark (readout, config)) ESCAPE (PS_ERR_UNKNOWN, "problem generating dark structure");138 if (!ppSimMakeDark (readout, config)) ESCAPE (PS_ERR_UNKNOWN, "problem generating dark structure"); 139 139 if (type == PPSIM_TYPE_DARK) goto done; 140 140 … … 153 153 154 154 done: 155 if (!ppSimAddNoise(readout->image, readout-> weight, cell, config, rng)) ESCAPE (PS_ERR_UNKNOWN, "problem adding noise");155 if (!ppSimAddNoise(readout->image, readout->variance, cell, config, rng)) ESCAPE (PS_ERR_UNKNOWN, "problem adding noise"); 156 156 if (!ppSimSaturate(readout, config)) ESCAPE (PS_ERR_UNKNOWN, "problem setting saturation levels"); 157 157 … … 165 165 readout->parent->parent->data_exists = true; 166 166 167 // if there is an input image, merge it with the simulated image168 if (!ppSimMergeReadouts (config, view)) ESCAPE (PS_ERR_UNKNOWN, "problem merging input image with simulated image");167 // if there is an input image, merge it with the simulated image 168 if (!ppSimMergeReadouts (config, view)) ESCAPE (PS_ERR_UNKNOWN, "problem merging input image with simulated image"); 169 169 } 170 170 psFree(biasCols); … … 173 173 174 174 if (cell->hdu) { 175 // XXX only do this if there is no INPUT image?175 // XXX only do this if there is no INPUT image? 176 176 if (!ppSimInitHeader(config, NULL, NULL, cell)) ESCAPE (PS_ERR_UNKNOWN, "problem setting output header"); 177 177 } … … 186 186 } 187 187 188 // XXX why no UpdateConceptsChip??188 // XXX why no UpdateConceptsChip?? 189 189 190 190 if (chip->hdu) { 191 // XXX only do this if there is no INPUT image191 // XXX only do this if there is no INPUT image 192 192 if (!ppSimInitHeader(config, NULL, chip, NULL)) ESCAPE (PS_ERR_UNKNOWN, "problem setting output header"); 193 193 } 194 194 195 195 // we perform photometry on the readouts of this chip in the output 196 if (!ppSimPhotom (config, view)) ESCAPE (PS_ERR_UNKNOWN, "problem performing photometry");196 if (!ppSimPhotom (config, view)) ESCAPE (PS_ERR_UNKNOWN, "problem performing photometry"); 197 197 198 198 if (!pmFPAfileIOChecks(config, view, PM_FPA_AFTER)) { -
trunk/ppSim/src/ppSimMakeBias.c
r18011 r21365 20 20 float readnoise = psMetadataLookupF32(NULL, cell->concepts, "CELL.READNOISE");// CCD read noise, e 21 21 if (isnan(readnoise)) { 22 psWarning("CELL.READNOISE is not set; reverting to recipe value READNOISE.");23 readnoise = psMetadataLookupF32(&mdok, recipe, "READNOISE");24 if (!mdok) {25 psError(PS_ERR_BAD_PARAMETER_VALUE, false, "Unable to find READNOISE in recipe.");26 *status = false;27 return NULL;28 }22 psWarning("CELL.READNOISE is not set; reverting to recipe value READNOISE."); 23 readnoise = psMetadataLookupF32(&mdok, recipe, "READNOISE"); 24 if (!mdok) { 25 psError(PS_ERR_BAD_PARAMETER_VALUE, false, "Unable to find READNOISE in recipe."); 26 *status = false; 27 return NULL; 28 } 29 29 } 30 30 31 31 psImage *signal = readout->image; 32 psImage *variance = readout-> weight;32 psImage *variance = readout->variance; 33 33 34 34 int numRows = signal->numRows; … … 38 38 psPolynomial1D *biasPoly = psPolynomial1DAlloc(PS_POLYNOMIAL_CHEB, biasOrder); 39 39 for (int j = 0; j < biasOrder + 1; j++) { 40 biasPoly->coeff[j] = biasRange * psRandomGaussian(rng);40 biasPoly->coeff[j] = biasRange * psRandomGaussian(rng); 41 41 } 42 42 … … 45 45 46 46 for (int y = 0; y < numRows; y++) { 47 // Adjust bias level for this row48 biasRows->data.F32[y] = psPolynomial1DEval(biasPoly, (float)(y + biasOffset) /49 (float)numRows - 0.5) + biasLevel;47 // Adjust bias level for this row 48 biasRows->data.F32[y] = psPolynomial1DEval(biasPoly, (float)(y + biasOffset) / 49 (float)numRows - 0.5) + biasLevel; 50 50 51 for (int x = 0; x < numCols; x++) {51 for (int x = 0; x < numCols; x++) { 52 52 53 // Bias level54 signal->data.F32[y][x] += biasRows->data.F32[y];55 variance->data.F32[y][x] += PS_SQR(readnoise);53 // Bias level 54 signal->data.F32[y][x] += biasRows->data.F32[y]; 55 variance->data.F32[y][x] += PS_SQR(readnoise); 56 56 57 }57 } 58 58 } 59 59 psFree(biasPoly); -
trunk/ppSim/src/ppSimMakeDark.c
r18011 r21365 7 7 8 8 psImage *signal = readout->image; 9 psImage *variance = readout-> weight;9 psImage *variance = readout->variance; 10 10 11 11 psMetadata *recipe = psMetadataLookupMetadata(&mdok, config->recipes, PPSIM_RECIPE); // Recipe … … 29 29 30 30 for (int y = 0; y < signal->numRows; y++) { 31 for (int x = 0; x < signal->numCols; x++) {32 33 // Dark current34 float darkCurrent = darkRate * expTime; // Dark current accumulated35 signal->data.F32[y][x] += darkCurrent;36 variance->data.F32[y][x] += darkCurrent;37 }31 for (int x = 0; x < signal->numCols; x++) { 32 33 // Dark current 34 float darkCurrent = darkRate * expTime; // Dark current accumulated 35 signal->data.F32[y][x] += darkCurrent; 36 variance->data.F32[y][x] += darkCurrent; 37 } 38 38 } 39 39 return true; -
trunk/ppSim/src/ppSimMakeSky.c
r20366 r21365 9 9 10 10 psImage *signal = readout->image; 11 psImage *variance = readout-> weight;11 psImage *variance = readout->variance; 12 12 13 13 pmCell *cell = readout->parent; … … 19 19 bool sky = psMetadataLookupBool(&status, recipe, "SKY"); // Generate a SKY flux? 20 20 bool flat = psMetadataLookupBool(&status, recipe, "FLAT"); // Apply flat-field term? 21 21 22 22 float expTime = psMetadataLookupF32(&status, recipe, "EXPTIME"); // Exposure time 23 23 … … 29 29 float skyMags = psMetadataLookupF32(&status, recipe, "SKY.MAGS"); assert (status); 30 30 if (!isnan(skyMags)) { 31 float zp = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); assert (status);32 float scale = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); assert (status);33 skyRate = scale * scale * ppSimMagToFlux (skyMags, zp);31 float zp = psMetadataLookupF32(&status, recipe, "ZEROPOINT"); assert (status); 32 float scale = psMetadataLookupF32(&status, recipe, "PIXEL.SCALE"); assert (status); 33 skyRate = scale * scale * ppSimMagToFlux (skyMags, zp); 34 34 } 35 35 if (type == PPSIM_TYPE_FLAT) { … … 37 37 } 38 38 39 int x0Chip = psMetadataLookupS32(&status, chip->concepts, "CHIP.X0");40 int y0Chip = psMetadataLookupS32(&status, chip->concepts, "CHIP.Y0");39 int x0Chip = psMetadataLookupS32(&status, chip->concepts, "CHIP.X0"); 40 int y0Chip = psMetadataLookupS32(&status, chip->concepts, "CHIP.Y0"); 41 41 int xParityChip = psMetadataLookupS32(&status, chip->concepts, "CHIP.XPARITY"); 42 42 int yParityChip = psMetadataLookupS32(&status, chip->concepts, "CHIP.YPARITY"); 43 43 44 int x0Cell = psMetadataLookupS32(&status, cell->concepts, "CELL.X0");45 int y0Cell = psMetadataLookupS32(&status, cell->concepts, "CELL.Y0");44 int x0Cell = psMetadataLookupS32(&status, cell->concepts, "CELL.X0"); 45 int y0Cell = psMetadataLookupS32(&status, cell->concepts, "CELL.Y0"); 46 46 int xParityCell = psMetadataLookupS32(&status, cell->concepts, "CELL.XPARITY"); 47 47 int yParityCell = psMetadataLookupS32(&status, cell->concepts, "CELL.YPARITY"); … … 72 72 // Gaussian flat-field over the FPA with flatValue = 1.0 at the field center 73 73 float flatValue = 1.0; 74 if (flat) {75 // we make the flat-field have a response of 1.0 at the field center (like a vignetting)76 flatValue = expf(-0.5 / PS_SQR(flatSigma) * (PS_SQR(yFPA) + PS_SQR(xFPA)));77 }74 if (flat) { 75 // we make the flat-field have a response of 1.0 at the field center (like a vignetting) 76 flatValue = expf(-0.5 / PS_SQR(flatSigma) * (PS_SQR(yFPA) + PS_SQR(xFPA))); 77 } 78 78 79 float scatterRate = 0.0;79 float scatterRate = 0.0; 80 80 81 if (sky) {82 // add a scattered light term to the flat-field images83 if (type == PPSIM_TYPE_FLAT) {84 scatterRate = scatterFrac * PS_SQR(xFPA);85 }81 if (sky) { 82 // add a scattered light term to the flat-field images 83 if (type == PPSIM_TYPE_FLAT) { 84 scatterRate = scatterFrac * PS_SQR(xFPA); 85 } 86 86 87 // Sky background88 float skyFlux = (skyRate * (flatValue + scatterRate)) * realExpTime; // Flux from sky89 signal->data.F32[y][x] += skyFlux;90 variance->data.F32[y][x] += skyFlux;91 }87 // Sky background 88 float skyFlux = (skyRate * (flatValue + scatterRate)) * realExpTime; // Flux from sky 89 signal->data.F32[y][x] += skyFlux; 90 variance->data.F32[y][x] += skyFlux; 91 } 92 92 93 // used later to modify the star and galaxy photometry 94 if (expCorr) {95 // exposure correction is (effective exposure time) * (flatValue)96 expCorr->data.F32[y][x] = flatValue * realExpTime / expTime;97 }93 // used later to modify the star and galaxy photometry 94 if (expCorr) { 95 // exposure correction is (effective exposure time) * (flatValue) 96 expCorr->data.F32[y][x] = flatValue * realExpTime / expTime; 97 } 98 98 99 99 // TO DO: Add fringes -
trunk/ppSim/src/ppSimMergeReadouts.c
r18011 r21365 14 14 if (!inReadout) return true; 15 15 16 if (!inReadout-> weight) {17 if (!pmReadoutGenerate Weight(inReadout, true)) {18 psError (PS_ERR_UNKNOWN, false, "trouble creating weight");19 return false;16 if (!inReadout->variance) { 17 if (!pmReadoutGenerateVariance(inReadout, true)) { 18 psError (PS_ERR_UNKNOWN, false, "trouble creating variance"); 19 return false; 20 20 } 21 21 } … … 30 30 31 31 psImage *inSignal = inReadout->image; 32 psImage *inVariance = inReadout-> weight;32 psImage *inVariance = inReadout->variance; 33 33 34 34 psImage *outSignal = outReadout->image; 35 psImage *outVariance = outReadout-> weight;35 psImage *outVariance = outReadout->variance; 36 36 37 37 assert (inSignal->numRows == outSignal->numRows); … … 39 39 40 40 for (int y = 0; y < inSignal->numRows; y++) { 41 for (int x = 0; x < inSignal->numCols; x++) {42 outSignal->data.F32[y][x] += inSignal->data.F32[y][x];43 outVariance->data.F32[y][x] += inVariance->data.F32[y][x];44 }41 for (int x = 0; x < inSignal->numCols; x++) { 42 outSignal->data.F32[y][x] += inSignal->data.F32[y][x]; 43 outVariance->data.F32[y][x] += inVariance->data.F32[y][x]; 44 } 45 45 } 46 46 return true; -
trunk/ppSim/src/ppSimPhotomReadoutFake.c
r21183 r21365 12 12 } 13 13 14 # if 0 14 # if 0 15 15 // set the photcode for this image 16 16 if (!psphotAddPhotcode (recipe, config, view, "PPSIM.FAKE.CHIP")) { … … 22 22 // *** in this section, perform the photometry for real + fake sources on PPSIM.FAKE.CHIP *** 23 23 24 // find the currently selected readout. 24 // find the currently selected readout. 25 25 // we always perform photometry on the mosaiced chip 26 26 pmReadout *readout = pmFPAfileThisReadout (config->files, view, "PPSIM.FAKE.CHIP"); … … 35 35 PS_ASSERT_PTR_NON_NULL (injectedSources, false); 36 36 37 // Generate the mask and weightimages, including the user-defined analysis region of interest38 psphotSetMaskAnd Weight(config, readout, recipe);37 // Generate the mask and variance images, including the user-defined analysis region of interest 38 psphotSetMaskAndVariance (config, readout, recipe); 39 39 40 40 // XXX need to define the source pixels … … 53 53 psArray *realSources = psArrayAlloc (realMeasuredSources->n); 54 54 for (int i = 0; i < realMeasuredSources->n; i++) { 55 realSources->data[i] = pmSourceCopy (realMeasuredSources->data[i]);55 realSources->data[i] = pmSourceCopy (realMeasuredSources->data[i]); 56 56 } 57 57 58 58 // load the psf model, if suppled. FWHM_X,FWHM_Y,etc are saved in the recipe 59 // this function uses PSPHOT.PSF.LOAD as the pmFPAfile 59 // this function uses PSPHOT.PSF.LOAD as the pmFPAfile 60 60 pmPSF *psf = psphotLoadPSF (config, view, recipe); 61 61 assert (psf); 62 62 63 // remove all sources 63 // remove all sources 64 64 psphotRemoveAllSources (realSources, recipe); 65 65 … … 76 76 77 77 // XXX fake sources should measure peak->x,y, force sources should not 78 psImageMaskType maskVal = 0xff; 78 psImageMaskType maskVal = 0xff; 79 79 psImage *significance = psphotSignificanceImage (readout, recipe, 1, maskVal); 80 80 ppSimDetections (significance, recipe, fakeSources); … … 90 90 psphotGuessModels (config, readout, realSources, psf); 91 91 psphotGuessModels (config, readout, fakeSources, psf); 92 92 93 93 // linear fit to real + fake sources 94 94 psArray *sources = ppSimMergeSources (realSources, fakeSources); … … 113 113 pmReadout *fakeReadout = pmFPAfileThisReadout (config->files, view, "PPSIM.FAKE.SOURCES"); 114 114 if (!fakeReadout) { 115 fakeReadout = pmReadoutAlloc (fakeCell);116 psFree (fakeReadout); // there is a copy on 'cell' as well115 fakeReadout = pmReadoutAlloc (fakeCell); 116 psFree (fakeReadout); // there is a copy on 'cell' as well 117 117 } 118 118 psAssert (fakeReadout, "no fakeReadout?"); -
trunk/ppSim/src/ppSimPhotomReadoutForce.c
r21183 r21365 12 12 } 13 13 14 # if 0 14 # if 0 15 15 // set the photcode for this image 16 16 if (!psphotAddPhotcode (recipe, config, view, "PPSIM.FAKE.CHIP")) { … … 22 22 // *** in this section, perform the photometry for real + force sources on PPSIM.FORCE.CHIP *** 23 23 24 // find the currently selected readout. 24 // find the currently selected readout. 25 25 // we always perform photometry on the mosaiced chip 26 26 pmReadout *readout = pmFPAfileThisReadout (config->files, view, "PPSIM.FORCE.CHIP"); … … 35 35 psAssert (forceSources, "failed to load force photometry sources"); 36 36 37 // Generate the mask and weightimages, including the user-defined analysis region of interest38 psphotSetMaskAnd Weight(config, readout, recipe);37 // Generate the mask and variance images, including the user-defined analysis region of interest 38 psphotSetMaskAndVariance (config, readout, recipe); 39 39 40 40 // XXX need to define the source pixels … … 46 46 psArray *realSources = psArrayAlloc (realMeasuredSources->n); 47 47 for (int i = 0; i < realMeasuredSources->n; i++) { 48 realSources->data[i] = pmSourceCopy (realMeasuredSources->data[i]);48 realSources->data[i] = pmSourceCopy (realMeasuredSources->data[i]); 49 49 } 50 50 51 51 // load the psf model, if suppled. FWHM_X,FWHM_Y,etc are saved in the recipe 52 // this function uses PSPHOT.PSF.LOAD as the pmFPAfile 52 // this function uses PSPHOT.PSF.LOAD as the pmFPAfile 53 53 pmPSF *psf = psphotLoadPSF (config, view, recipe); 54 54 assert (psf); 55 55 56 // remove all sources 56 // remove all sources 57 57 psphotRemoveAllSources (realSources, recipe); 58 58 … … 69 69 70 70 // XXX fake sources should measure peak->x,y, force sources should not 71 psImageMaskType maskVal = 0xff; 71 psImageMaskType maskVal = 0xff; 72 72 psImage *significance = psphotSignificanceImage (readout, recipe, 1, maskVal); 73 73 ppSimDetections (significance, recipe, forceSources); … … 83 83 psphotGuessModels (config, readout, realSources, psf); 84 84 psphotGuessModels (config, readout, forceSources, psf); 85 85 86 86 // linear fit to real + force sources 87 87 psArray *sources = ppSimMergeSources (realSources, forceSources); … … 106 106 pmReadout *forceReadout = pmFPAfileThisReadout (config->files, view, "PPSIM.FORCE.SOURCES"); 107 107 if (!forceReadout) { 108 forceReadout = pmReadoutAlloc (forceCell);109 psFree (forceReadout); // there is a copy on 'cell' as well108 forceReadout = pmReadoutAlloc (forceCell); 109 psFree (forceReadout); // there is a copy on 'cell' as well 110 110 } 111 111 psAssert (forceReadout, "no forceReadout?");
Note:
See TracChangeset
for help on using the changeset viewer.
