Changeset 27402 for trunk/ppStack/src
- Timestamp:
- Mar 23, 2010, 9:12:53 AM (16 years ago)
- Location:
- trunk/ppStack/src
- Files:
-
- 16 edited
-
ppStack.h (modified) (2 diffs)
-
ppStackCamera.c (modified) (3 diffs)
-
ppStackCleanup.c (modified) (1 diff)
-
ppStackCombineFinal.c (modified) (8 diffs)
-
ppStackCombineInitial.c (modified) (1 diff)
-
ppStackCombinePrepare.c (modified) (3 diffs)
-
ppStackFiles.c (modified) (1 diff)
-
ppStackLoop.c (modified) (5 diffs)
-
ppStackLoop.h (modified) (2 diffs)
-
ppStackOptions.c (modified) (4 diffs)
-
ppStackOptions.h (modified) (2 diffs)
-
ppStackPrepare.c (modified) (6 diffs)
-
ppStackReadout.c (modified) (8 diffs)
-
ppStackReject.c (modified) (1 diff)
-
ppStackSources.c (modified) (7 diffs)
-
ppStackThread.c (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
trunk/ppStack/src/ppStack.h
r27400 r27402 69 69 const psVector *mask, // Mask for input readouts 70 70 const psVector *weightings, // Weighting factors for each image 71 const psVector *exposures, // Exposure time for each image72 71 const psVector *addVariance // Additional variance for rejection 73 72 ); … … 84 83 bool ppStackReadoutFinal(const pmConfig *config, // Configuration 85 84 pmReadout *outRO, // Output readout 86 pmReadout *expRO, // Exposure readout87 85 const psArray *readouts, // Input readouts 88 86 const psVector *mask, // Mask for input readouts 89 87 const psArray *rejected, // Array with pixels rejected in each image 90 88 const psVector *weightings, // Weighting factors for each image 91 const psVector *exposures, // Exposure times for each image92 89 const psVector *addVariance, // Additional variance for rejection 93 90 bool safety, // Enable safety switch? -
trunk/ppStack/src/ppStackCamera.c
r27400 r27402 283 283 } 284 284 285 286 // Exposure image287 pmFPA *expFPA = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain the output288 if (!expFPA) {289 psError(psErrorCodeLast(), false, "Unable to construct an FPA from camera configuration.");290 return false;291 }292 pmFPAfile *exp = pmFPAfileDefineOutput(config, expFPA, "PPSTACK.OUTPUT.EXP");293 psFree(expFPA); // Drop reference294 if (!exp) {295 psError(psErrorCodeLast(), false, _("Unable to generate output file from PPSTACK.OUTPUT.EXP"));296 return false;297 }298 if (exp->type != PM_FPA_FILE_IMAGE) {299 psError(PPSTACK_ERR_CONFIG, true, "PPSTACK.OUTPUT.EXP is not of type IMAGE");300 return false;301 }302 exp->save = true;303 304 if (!pmFPAAddSourceFromFormat(expFPA, "Stack", exp->format)) {305 psError(psErrorCodeLast(), false, "Unable to generate output FPA.");306 return false;307 }308 309 // Exposure numbers310 pmFPAfile *expNum = pmFPAfileDefineOutput(config, exp->fpa, "PPSTACK.OUTPUT.EXPNUM");311 if (!expNum) {312 psError(psErrorCodeLast(), false, _("Unable to generate output file from PPSTACK.OUTPUT.MASK"));313 return false;314 }315 if (expNum->type != PM_FPA_FILE_MASK) {316 psError(PPSTACK_ERR_CONFIG, true, "PPSTACK.OUTPUT.EXPNUM is not of type MASK");317 return false;318 }319 expNum->save = true;320 321 322 323 285 if (havePSFs) { 324 286 pmFPA *psfFPA = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain PSF … … 340 302 } 341 303 304 #if 1 342 305 // Unconvolved stack 343 306 pmFPA *unconvFPA = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain unconvolved output … … 388 351 unconvVariance->save = true; 389 352 } 390 391 392 // Exposure image 393 pmFPA *unconvExpFPA = pmFPAConstruct(config->camera, config->cameraName); // FPA to contain the output 394 if (!unconvExpFPA) { 395 psError(psErrorCodeLast(), false, "Unable to construct an FPA from camera configuration."); 396 return false; 397 } 398 pmFPAfile *unconvExp = pmFPAfileDefineOutput(config, unconvExpFPA, "PPSTACK.UNCONV.EXP"); 399 psFree(unconvExpFPA); // Drop reference 400 if (!unconvExp) { 401 psError(psErrorCodeLast(), false, _("Unable to generate output file from PPSTACK.UNCONV.EXP")); 402 return false; 403 } 404 if (unconvExp->type != PM_FPA_FILE_IMAGE) { 405 psError(PPSTACK_ERR_CONFIG, true, "PPSTACK.UNCONV.EXP is not of type IMAGE"); 406 return false; 407 } 408 unconvExp->save = true; 409 410 if (!pmFPAAddSourceFromFormat(unconvExpFPA, "Stack", unconvExp->format)) { 411 psError(psErrorCodeLast(), false, "Unable to generate output FPA."); 412 return false; 413 } 414 415 // Output mask 416 pmFPAfile *unconvExpNum = pmFPAfileDefineOutput(config, unconvExp->fpa, "PPSTACK.UNCONV.EXPNUM"); 417 if (!unconvExpNum) { 418 psError(psErrorCodeLast(), false, _("Unable to generate output file from PPSTACK.UNCONV.MASK")); 419 return false; 420 } 421 if (unconvExpNum->type != PM_FPA_FILE_MASK) { 422 psError(PPSTACK_ERR_CONFIG, true, "PPSTACK.UNCONV.EXPNUM is not of type MASK"); 423 return false; 424 } 425 unconvExpNum->save = true; 426 353 #endif 427 354 428 355 // Output JPEGs -
trunk/ppStack/src/ppStackCleanup.c
r27400 r27402 86 86 options->outRO = NULL; 87 87 88 options->expRO->data_exists = false;89 options->expRO->parent->data_exists = false;90 options->expRO->parent->parent->data_exists = false;91 psFree(options->expRO);92 options->expRO = NULL;93 94 88 pmFPAview *view = pmFPAviewAlloc(0);// Pointer into FPA hierarchy 95 89 view->chip = view->cell = 0; // pmFPAviewFreeData doesn't want to deal with readouts -
trunk/ppStack/src/ppStackCombineFinal.c
r27400 r27402 14 14 //#define TESTING // Enable test output 15 15 16 bool ppStackCombineFinal(p pStackThreadData *stack, psArray *covariances, ppStackOptions *options,17 p mConfig *config, bool safe, bool normalise, bool grow)16 bool ppStackCombineFinal(pmReadout *target, ppStackThreadData *stack, psArray *covariances, 17 ppStackOptions *options, pmConfig *config, bool safe, bool normalise, bool grow) 18 18 { 19 19 psAssert(stack, "Require stack"); … … 21 21 psAssert(config, "Require configuration"); 22 22 23 pmReadout *outRO = options->outRO; // Output readout 24 pmReadout *expRO = options->expRO; // Exposure readout 25 int numCols = outRO->image->numCols, numRows = outRO->image->numRows; // Size of image 23 int numCols = target->image->numCols, numRows = target->image->numRows; // Size of image 26 24 27 25 psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, PPSTACK_RECIPE); // ppStack recipe … … 45 43 } 46 44 47 if (!outRO->mask) { 48 outRO->mask = psImageAlloc(numCols, numRows, PS_TYPE_IMAGE_MASK); 49 } 50 if (!expRO->mask) { 51 expRO->mask = psImageAlloc(numCols, numRows, PS_TYPE_IMAGE_MASK); 45 if (!target->mask) { 46 target->mask = psImageAlloc(target->image->numCols, target->image->numRows, PS_TYPE_IMAGE_MASK); 52 47 } 53 48 … … 67 62 } 68 63 69 // call: ppStackReadoutFinal(config, outRO, readouts, rejected)64 // call: ppStackReadoutFinal(config, target, readouts, rejected) 70 65 psThreadJob *job = psThreadJobAlloc("PPSTACK_FINAL_COMBINE"); // Job to start 66 psArrayAdd(job->args, 1, target); 71 67 psArrayAdd(job->args, 1, thread); 72 68 psArrayAdd(job->args, 1, reject); … … 109 105 } 110 106 if (sumWeights > 0.0) { 111 outRO->covariance = psImageCovarianceSum(covariances);112 psBinaryOp( outRO->covariance->image, outRO->covariance->image, "/",107 target->covariance = psImageCovarianceSum(covariances); 108 psBinaryOp(target->covariance->image, target->covariance->image, "/", 113 109 psScalarAlloc(sumWeights, PS_TYPE_F32)); 114 psImageCovarianceTransfer( outRO->variance, outRO->covariance);110 psImageCovarianceTransfer(target->variance, target->covariance); 115 111 } 116 112 } else { 117 outRO->covariance = psImageCovarianceNone();113 target->covariance = psImageCovarianceNone(); 118 114 } 119 115 … … 131 127 wcsDone = true; 132 128 pmHDU *inHDU = pmHDUFromCell(inRO->parent); // Template HDU 133 pmHDU *outHDU = pmHDUFromCell( outRO->parent); // Output HDU134 pmChip *outChip = outRO->parent->parent; // Output chip129 pmHDU *outHDU = pmHDUFromCell(target->parent); // Output HDU 130 pmChip *outChip = target->parent->parent; // Output chip 135 131 pmFPA *outFPA = outChip->parent; // Output FPA 136 132 if (!outHDU || !inHDU) { … … 155 151 } 156 152 157 // Set exposure time correctly158 {159 float exptime = 0.0; // Summed exposure time160 for (int i = 0; i < options->num; i++) {161 if (options->inputMask) {162 continue;163 }164 exptime += options->exposures->data.F32[i];165 }166 167 {168 psMetadataItem *item = psMetadataLookup(outRO->parent->concepts, "CELL.EXPOSURE");169 item->data.F32 = exptime;170 }171 {172 psMetadataItem *item = psMetadataLookup(outRO->parent->parent->parent->concepts, "FPA.EXPOSURE");173 item->data.F32 = exptime;174 }175 }176 177 153 // Put version information into the header 178 pmHDU *hdu = pmHDUFromCell( outRO->parent);154 pmHDU *hdu = pmHDUFromCell(target->parent); 179 155 if (!hdu) { 180 156 psError(PPSTACK_ERR_PROG, false, "Unable to find HDU for output."); … … 192 168 psStringAppend(&name, "combined_image_final_%d.fits", pass); 193 169 pass++; 194 ppStackWriteImage(name, NULL, outRO->image, config);170 ppStackWriteImage(name, NULL, target->image, config); 195 171 psStringSubstitute(&name, "mask", "image"); 196 ppStackWriteImage(name, NULL, outRO->mask, config);172 ppStackWriteImage(name, NULL, target->mask, config); 197 173 psStringSubstitute(&name, "variance", "mask"); 198 ppStackWriteImage(name, NULL, outRO->variance, config);174 ppStackWriteImage(name, NULL, target->variance, config); 199 175 psFree(name); 200 176 201 pmStackVisualPlotTestImage( outRO->image, "combined_image_final.fits");177 pmStackVisualPlotTestImage(target->image, "combined_image_final.fits"); 202 178 #endif 203 179 -
trunk/ppStack/src/ppStackCombineInitial.c
r27400 r27402 89 89 } 90 90 91 ppStackMemDump("initial");92 93 91 #ifdef TESTING 94 92 ppStackWriteImage("combined_image_initial.fits", NULL, options->outRO->image, config); -
trunk/ppStack/src/ppStackCombinePrepare.c
r27400 r27402 10 10 #include "ppStackLoop.h" 11 11 12 bool ppStackCombinePrepare(const char *outName, const char *expName, 13 ppStackFileList files, ppStackThreadData *stack, 12 bool ppStackCombinePrepare(pmReadout **ro, const char *name, ppStackFileList files, ppStackThreadData *stack, 14 13 ppStackOptions *options, pmConfig *config) 15 14 { … … 33 32 } 34 33 35 pmCell *cell = pmFPAfileThisCell(config->files, view, outName); // Output cell 36 options->outRO = pmReadoutAlloc(cell); // Output readout 37 38 pmCell *expCell = pmFPAfileThisCell(config->files, view, expName); // Exposure cell 39 options->expRO = pmReadoutAlloc(expCell); // Output readout 34 pmCell *cell = pmFPAfileThisCell(config->files, view, name); // Output cell 35 *ro = pmReadoutAlloc(cell); // Output readout 40 36 41 37 psFree(view); … … 46 42 psImageMaskType maskBad = pmConfigMaskGet(maskBadStr, config); // Bits to mask for bad pixels 47 43 48 if (!pmReadoutStackDefineOutput(options->outRO, col0, row0, numCols, numRows, true, true, maskBad)) { 49 psError(PPSTACK_ERR_ARGUMENTS, false, "Unable to prepare output."); 50 return false; 51 } 52 53 if (!pmReadoutStackDefineOutput(options->expRO, col0, row0, numCols, numRows, true, true, 0)) { 44 if (!pmReadoutStackDefineOutput(*ro, col0, row0, numCols, numRows, true, true, maskBad)) { 54 45 psError(PPSTACK_ERR_ARGUMENTS, false, "Unable to prepare output."); 55 46 return false; -
trunk/ppStack/src/ppStackFiles.c
r27400 r27402 23 23 /// Regular (convolved) stack files 24 24 static char *filesStack[] = { "PPSTACK.OUTPUT", "PPSTACK.OUTPUT.MASK", "PPSTACK.OUTPUT.VARIANCE", 25 "PPSTACK.OUTPUT.EXP", "PPSTACK.OUTPUT.EXPNUM", 26 "PPSTACK.OUTPUT.JPEG1", "PPSTACK.OUTPUT.JPEG2", 27 NULL }; 25 "PPSTACK.OUTPUT.JPEG1", "PPSTACK.OUTPUT.JPEG2", NULL }; 28 26 /// Unconvolved stack files 29 static char *filesUnconv[] = { "PPSTACK.UNCONV", "PPSTACK.UNCONV.MASK", "PPSTACK.UNCONV.VARIANCE", 30 "PPSTACK.UNCONV.EXP", "PPSTACK.UNCONV.EXPNUM", NULL }; 27 static char *filesUnconv[] = { "PPSTACK.UNCONV", "PPSTACK.UNCONV.MASK", "PPSTACK.UNCONV.VARIANCE", NULL }; 31 28 32 29 /// Files for photometry -
trunk/ppStack/src/ppStackLoop.c
r27400 r27402 111 111 112 112 // Prepare for combination 113 if (!ppStackCombinePrepare( "PPSTACK.OUTPUT", "PPSTACK.OUTPUT.EXP", PPSTACK_FILES_STACK,113 if (!ppStackCombinePrepare(&options->outRO, "PPSTACK.OUTPUT", PPSTACK_FILES_STACK, 114 114 stack, options, config)) { 115 115 psError(psErrorCodeLast(), false, "Unable to prepare for combination."); … … 160 160 // Final combination 161 161 psTrace("ppStack", 2, "Final stack of convolved images....\n"); 162 if (!ppStackCombineFinal(stack, options->convCovars, options, config, false, false, true)) { 162 if (!ppStackCombineFinal(options->outRO, stack, options->convCovars, options, config, 163 false, false, true)) { 163 164 psError(psErrorCodeLast(), false, "Unable to perform final combination."); 164 165 psFree(stack); … … 202 203 203 204 // Prepare for combination 204 if (!ppStackCombinePrepare( "PPSTACK.UNCONV", "PPSTACK.UNCONV.EXP", PPSTACK_FILES_UNCONV,205 if (!ppStackCombinePrepare(&options->unconvRO, "PPSTACK.UNCONV", PPSTACK_FILES_UNCONV, 205 206 stack, options, config)) { 206 207 psError(psErrorCodeLast(), false, "Unable to prepare for combination."); … … 210 211 211 212 psTrace("ppStack", 2, "Stack of unconvolved images....\n"); 212 if (!ppStackCombineFinal( stack, options->origCovars, options, config,213 if (!ppStackCombineFinal(options->unconvRO, stack, options->origCovars, options, config, 213 214 false, true, false)) { 214 215 psError(psErrorCodeLast(), false, "Unable to perform unconvolved combination."); … … 225 226 } 226 227 ppStackFileActivation(config, PPSTACK_FILES_UNCONV, false); 227 options->outRO->data_exists = false; 228 options->outRO->parent->data_exists = false; 229 options->outRO->parent->parent->data_exists = false; 230 psFree(options->outRO); 231 options->outRO = NULL; 232 options->expRO->data_exists = false; 233 options->expRO->parent->data_exists = false; 234 options->expRO->parent->parent->data_exists = false; 235 psFree(options->expRO); 236 options->expRO = NULL; 228 options->unconvRO->data_exists = false; 229 options->unconvRO->parent->data_exists = false; 230 options->unconvRO->parent->parent->data_exists = false; 231 psFree(options->unconvRO); 232 options->unconvRO = NULL; 237 233 238 234 for (int i = 0; i < options->num; i++) { -
trunk/ppStack/src/ppStackLoop.h
r27400 r27402 38 38 // Prepare for combination 39 39 bool ppStackCombinePrepare( 40 const char *outName, // Name of output file41 const char * expName, // Name of exposurefile40 pmReadout **readout, // Readout to set 41 const char *name, // Name of file 42 42 ppStackFileList files, // Files of interest 43 43 ppStackThreadData *stack, // Stack … … 61 61 // Final combination 62 62 bool ppStackCombineFinal( 63 pmReadout *target, // Target readout 63 64 ppStackThreadData *stack, // Stack 64 65 psArray *covariances, // Covariances -
trunk/ppStack/src/ppStackOptions.c
r27400 r27402 21 21 psFree(options->psf); 22 22 psFree(options->inputSeeing); 23 psFree(options->exposures);24 23 psFree(options->inputMask); 25 24 psFree(options->sourceLists); … … 33 32 psFree(options->convCovars); 34 33 psFree(options->outRO); 35 psFree(options-> expRO);34 psFree(options->unconvRO); 36 35 psFree(options->inspect); 37 36 psFree(options->rejected); … … 62 61 options->zp = NAN; 63 62 options->inputSeeing = NULL; 64 options->exposures = NULL;65 63 options->targetSeeing = NAN; 66 64 options->inputMask = NULL; … … 77 75 options->convCovars = NULL; 78 76 options->outRO = NULL; 79 options-> expRO = NULL;77 options->unconvRO = NULL; 80 78 options->inspect = NULL; 81 79 options->rejected = NULL; -
trunk/ppStack/src/ppStackOptions.h
r27400 r27402 22 22 psVector *inputSeeing; // Input seeing FWHMs 23 23 float targetSeeing; // Target seeing FWHM 24 psVector *exposures; // Exposure times25 24 float sumExposure; // Sum of exposure times 26 25 float zp; // Zero point for output … … 39 38 // Combine initial 40 39 pmReadout *outRO; // Output readout 41 pmReadout * expRO; // Exposurereadout40 pmReadout *unconvRO; // Unconvolved readout 42 41 psArray *inspect; // Array of arrays of pixels to inspect 43 42 // Rejection -
trunk/ppStack/src/ppStackPrepare.c
r27400 r27402 126 126 options->inputMask = psVectorAlloc(num, PS_TYPE_VECTOR_MASK); // Mask for inputs 127 127 psVectorInit(options->inputMask, 0); 128 options->exposures = psVectorAlloc(options->num, PS_TYPE_F32);129 psVectorInit(options->exposures, NAN);130 128 131 129 pmFPAfileActivate(config->files, false, NULL); … … 136 134 } 137 135 136 psMetadataIterator *fileIter = psMetadataIteratorAlloc(config->files, PS_LIST_HEAD, "^PPSTACK.INPUT$"); 137 psMetadataItem *fileItem; // Item from iteration 138 138 psArray *psfs = psArrayAlloc(num); // PSFs for PSF envelope 139 139 int numCols = 0, numRows = 0; // Size of image 140 options->sumExposure = 0.0; 141 for (int i = 0; i < num; i++) { 142 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.INPUT", i); // File of interest 143 pmCell *cell = pmFPAviewThisCell(view, file->fpa); // Cell of interest 144 145 options->exposures->data.F32[i] = psMetadataLookupF32(NULL, cell->concepts, "CELL.EXPOSURE"); 146 options->sumExposure += options->exposures->data.F32[i]; 140 int index = 0; // Index for file 141 while ((fileItem = psMetadataGetAndIncrement(fileIter))) { 142 assert(fileItem->type == PS_DATA_UNKNOWN); 143 pmFPAfile *inputFile = fileItem->data.V; // An input file 147 144 148 145 // Get list of PSFs, to determine target PSF 149 146 if (options->convolve) { 150 pmChip *chip = pmFPAviewThisChip(view, file->fpa); // The chip: holds the PSF147 pmChip *chip = pmFPAviewThisChip(view, inputFile->fpa); // The chip: holds the PSF 151 148 pmPSF *psf = psMetadataLookupPtr(NULL, chip->analysis, "PSPHOT.PSF"); // PSF 152 149 if (!psf) { 153 150 psError(PPSTACK_ERR_PROG, false, "Unable to find PSF."); 154 151 psFree(view); 152 psFree(fileIter); 155 153 psFree(psfs); 156 154 return false; 157 155 } 158 psfs->data[i ] = psMemIncrRefCounter(psf);156 psfs->data[index] = psMemIncrRefCounter(psf); 159 157 psMetadataRemoveKey(chip->analysis, "PSPHOT.PSF"); 160 158 161 pmCell *cell = pmFPAviewThisCell(view, file->fpa); // Cell of interest159 pmCell *cell = pmFPAviewThisCell(view, inputFile->fpa); // Cell of interest 162 160 pmHDU *hdu = pmHDUFromCell(cell); 163 161 assert(hdu && hdu->header); … … 167 165 psError(PPSTACK_ERR_PROG, false, "Unable to determine size of image from PSF."); 168 166 psFree(view); 167 psFree(fileIter); 169 168 psFree(psfs); 170 169 return false; … … 181 180 pmDetections *detections = NULL; 182 181 if (options->convolve || options->matchZPs || options->photometry || redoPhot) { 183 pmReadout *ro = pmFPAviewThisReadout(view, file->fpa); // Readout with sources182 pmReadout *ro = pmFPAviewThisReadout(view, inputFile->fpa); // Readout with sources 184 183 detections = psMetadataLookupPtr(NULL, ro->analysis, "PSPHOT.DETECTIONS"); // Sources 185 184 if (!detections) { … … 189 188 psAssert (detections->allSources, "missing sources?"); 190 189 191 options->sourceLists->data[i ] = psMemIncrRefCounter(detections->allSources);190 options->sourceLists->data[index] = psMemIncrRefCounter(detections->allSources); 192 191 } 193 192 194 193 // Re-do photometry if we don't trust the source lists 195 194 if (redoPhot) { 196 psTrace("ppStack", 2, "Photometering input %d of %d....\n", i , num);195 psTrace("ppStack", 2, "Photometering input %d of %d....\n", index, num); 197 196 pmFPAfileActivate(config->files, false, NULL); 198 ppStackFileActivationSingle(config, PPSTACK_FILES_CONVOLVE, true, i );197 ppStackFileActivationSingle(config, PPSTACK_FILES_CONVOLVE, true, index); 199 198 if (options->convolve) { 200 199 pmFPAfileActivate(config->files, true, "PPSTACK.CONV.KERNEL"); 201 200 } 202 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.INPUT", i ); // File201 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.INPUT", index); // File 203 202 pmFPAview *photView = ppStackFilesIterateDown(config); 204 203 if (!photView) { … … 224 223 ppStackFileActivation(config, PPSTACK_FILES_PREPARE, true); 225 224 } 226 } 225 226 index++; 227 } 228 psFree(fileIter); 227 229 228 230 psString log = psStringCopy("Input seeing FWHMs:\n"); // Log message -
trunk/ppStack/src/ppStackReadout.c
r27400 r27402 23 23 psVector *mask = options->inputMask; // Mask for inputs 24 24 psVector *weightings = options->weightings; // Weightings (1/noise^2) for each image 25 psVector *exposures = options->exposures; // Exposure times for each image26 25 psVector *addVariance = options->matchChi2; // Additional variance when rejecting 27 26 28 27 job->results = ppStackReadoutInitial(config, outRO, thread->readouts, mask, 29 weightings, exposures,addVariance);28 weightings, addVariance); 30 29 thread->busy = false; 31 30 … … 38 37 39 38 psArray *args = job->args; // Arguments 40 ppStackThread *thread = args->data[0]; // Thread 41 psArray *reject = args->data[1]; // Rejected pixels for each image 42 ppStackOptions *options = args->data[2]; // Options 43 pmConfig *config = args->data[3]; // Configuration 44 bool safety = PS_SCALAR_VALUE(args->data[4], U8); // Safety switch on? 45 bool normalise = PS_SCALAR_VALUE(args->data[5], U8); // Normalise images? 39 pmReadout *target = args->data[0]; // Output readout 40 ppStackThread *thread = args->data[1]; // Thread 41 psArray *reject = args->data[2]; // Rejected pixels for each image 42 ppStackOptions *options = args->data[3]; // Options 43 pmConfig *config = args->data[4]; // Configuration 44 bool safety = PS_SCALAR_VALUE(args->data[5], U8); // Safety switch on? 45 bool normalise = PS_SCALAR_VALUE(args->data[6], U8); // Normalise images? 46 46 47 47 psVector *mask = options->inputMask; // Mask for inputs 48 48 psVector *weightings = options->weightings; // Weightings (1/noise^2) for each image 49 psVector *exposures = options->exposures; // Exposure times for each image50 49 psVector *addVariance = options->matchChi2; // Additional variance when rejecting 51 50 psVector *norm = normalise ? options->norm : NULL; // Normalisations to apply to images 52 51 53 bool status = ppStackReadoutFinal(config, options->outRO, options->expRO, thread->readouts, mask, reject,54 weightings, exposures,addVariance, safety, norm); // Status of operation52 bool status = ppStackReadoutFinal(config, target, thread->readouts, mask, reject, 53 weightings, addVariance, safety, norm); // Status of operation 55 54 56 55 thread->busy = false; 57 58 psAssert(status, "Stacking failed.");59 56 60 57 return status; … … 104 101 105 102 psArray *ppStackReadoutInitial(const pmConfig *config, pmReadout *outRO, const psArray *readouts, 106 const psVector *mask, const psVector *weightings, const psVector *exposures, 107 const psVector *addVariance) 103 const psVector *mask, const psVector *weightings, const psVector *addVariance) 108 104 { 109 105 assert(config); … … 153 149 } 154 150 155 stack->data[i] = pmStackDataAlloc(ro, weightings->data.F32[i], exposures->data.F32[i], 156 addVariance->data.F32[i]); 157 } 158 159 if (!pmStackCombine(outRO, NULL, stack, maskVal | maskBad, maskSuspect, maskBad, kernelSize, iter, 151 stack->data[i] = pmStackDataAlloc(ro, weightings->data.F32[i], addVariance->data.F32[i]); 152 } 153 154 if (!pmStackCombine(outRO, stack, maskVal | maskBad, maskSuspect, maskBad, kernelSize, iter, 160 155 combineRej, combineSys, combineDiscard, useVariance, safe, false)) { 161 156 psError(PS_ERR_UNKNOWN, false, "Unable to combine input readouts with rejection."); … … 192 187 193 188 194 bool ppStackReadoutFinal(const pmConfig *config, pmReadout *outRO, pmReadout *expRO,const psArray *readouts,189 bool ppStackReadoutFinal(const pmConfig *config, pmReadout *outRO, const psArray *readouts, 195 190 const psVector *mask, const psArray *rejected, const psVector *weightings, 196 const psVector *exposures, const psVector *addVariance, bool safety, 197 const psVector *norm) 191 const psVector *addVariance, bool safety, const psVector *norm) 198 192 { 199 193 assert(config); 200 194 assert(outRO); 201 assert(expRO);202 195 assert(readouts); 203 196 assert(!rejected || readouts->n == rejected->n); … … 245 238 } 246 239 247 pmStackData *data = pmStackDataAlloc(ro, weightings->data.F32[i], exposures->data.F32[i],240 pmStackData *data = pmStackDataAlloc(ro, weightings->data.F32[i], 248 241 addVariance ? addVariance->data.F32[i] : NAN); 249 242 data->reject = rejected ? psMemIncrRefCounter(rejected->data[i]) : NULL; … … 257 250 } 258 251 259 if (!pmStackCombine(outRO, expRO,stack, maskVal | maskBad, maskSuspect, maskBad, 0, iter, combineRej,252 if (!pmStackCombine(outRO, stack, maskVal | maskBad, maskSuspect, maskBad, 0, iter, combineRej, 260 253 combineSys, combineDiscard, useVariance, safe, rejected)) { 261 254 psError(PS_ERR_UNKNOWN, false, "Unable to combine input readouts."); … … 266 259 pmCell *outCell = outRO->parent; // Output cell 267 260 pmChip *outChip = outCell->parent; // Output chip 261 268 262 outRO->data_exists = true; 269 263 outCell->data_exists = true; 270 264 outChip->data_exists = true; 271 265 272 pmCell *expCell = expRO->parent; // Exposure cell273 pmChip *expChip = expCell->parent; // Exposure chip274 expRO->data_exists = true;275 expCell->data_exists = true;276 expChip->data_exists = true;277 278 266 psFree(stack); 279 267 -
trunk/ppStack/src/ppStackReject.c
r27400 r27402 162 162 psLogMsg("ppStack", PS_LOG_INFO, "Time to perform rejection on image %d: %f sec", i, 163 163 psTimerClear("PPSTACK_REJECT")); 164 165 ppStackMemDump("reject"); 164 166 } 165 167 -
trunk/ppStack/src/ppStackSources.c
r27400 r27402 64 64 65 65 if (!options->matchZPs && !options->photometry) { 66 options->norm = psVectorAlloc(options->num, PS_TYPE_F32); 67 psVectorInit(options->norm, 0.0); 66 int num = psMetadataLookupS32(NULL, config->arguments, "INPUTS.NUM"); // Number of inputs 67 options->norm = psVectorAlloc(num, PS_TYPE_F32); 68 psVectorInit (options->norm, 0.0); 69 70 // XXX do I need to set this? 71 // options->sumExposure = sumExpTime; 72 68 73 return true; 69 74 } … … 132 137 } 133 138 134 int num = options->num;// Number of inputs139 int num = psMetadataLookupS32(NULL, config->arguments, "INPUTS.NUM"); // Number of inputs 135 140 psAssert(num == sourceLists->n, "Wrong number of source lists: %ld\n", sourceLists->n); 136 141 … … 141 146 float airmassTerm = NAN; // Airmass term 142 147 float zpTarget = NAN; // Target zero point 148 float sumExpTime = 0.0; // Sum of the exposure time 143 149 int numGoodImages = 0; // Number of good images 144 150 for (int i = 0; i < num; i++) { … … 154 160 155 161 pmFPAfile *file = pmFPAfileSelectSingle(config->files, "PPSTACK.INPUT", i); // File of interest 162 pmCell *cell = pmFPAviewThisCell(view, file->fpa); // Cell of interest 156 163 157 164 #if defined(TESTING) && 0 … … 170 177 #endif 171 178 172 float exptime = options->exposures->data.F32[i]; // Exposure time179 float exptime = psMetadataLookupF32(NULL, cell->concepts, "CELL.EXPOSURE"); // Exposure time 173 180 float airmass = psMetadataLookupF32(NULL, file->fpa->concepts, "FPA.AIRMASS"); // Airmass 174 181 const char *expFilter = psMetadataLookupStr(NULL, file->fpa->concepts, "FPA.FILTER"); // Filter name … … 214 221 215 222 zp->data.F32[i] = airmassTerm * airmass + 2.5 * log10(exptime); 216 } 223 sumExpTime += exptime; 224 225 } 226 227 options->sumExposure = sumExpTime; 217 228 218 229 if (numGoodImages == 0) { … … 280 291 } 281 292 psArray *sources = sourceLists->data[i]; // Sources of interest 282 float magCorr = zp->data.F32[i] + trans->data.F32[i] - 2.5*log10( options->sumExposure);293 float magCorr = zp->data.F32[i] + trans->data.F32[i] - 2.5*log10(sumExpTime); 283 294 if (zpExpNum == numGoodImages) { 284 295 // Using measured zero points, so attempt to set target zero point -
trunk/ppStack/src/ppStackThread.c
r27400 r27402 284 284 285 285 { 286 psThreadTask *task = psThreadTaskAlloc("PPSTACK_FINAL_COMBINE", 6);286 psThreadTask *task = psThreadTaskAlloc("PPSTACK_FINAL_COMBINE", 7); 287 287 task->function = &ppStackReadoutFinalThread; 288 288 psThreadTaskAdd(task);
Note:
See TracChangeset
for help on using the changeset viewer.
