Changeset 41413
- Timestamp:
- Sep 10, 2020, 9:40:15 AM (6 years ago)
- Location:
- branches/eam_branches/gDiff.20200202/src
- Files:
-
- 16 added
- 12 edited
-
Makefile.am (modified) (2 diffs)
-
gDiff.c (modified) (4 diffs)
-
gDiff.h (modified) (1 diff)
-
gDiffArguments.c (modified) (2 diffs)
-
gDiffBackground.c (added)
-
gDiffCamera.c (modified) (3 diffs)
-
gDiffData.c (modified) (1 diff)
-
gDiffDefineOutput.c (modified) (1 diff)
-
gDiffErrorCodes.dat (modified) (1 diff)
-
gDiffExit.c (modified) (2 diffs)
-
gDiffExtras.c (added)
-
gDiffFiles.c (modified) (3 diffs)
-
gDiffFlagNeighbors.c (added)
-
gDiffInputDetections.c (added)
-
gDiffLoop.c (modified) (2 diffs)
-
gDiffMakePSF.c (added)
-
gDiffMaskStats.c (added)
-
gDiffMatchPSFs.c (added)
-
gDiffReadoutForcedPhot.c (added)
-
gDiffReadoutInverse.c (added)
-
gDiffReadoutJpeg.c (added)
-
gDiffReadoutPhotometry.c (added)
-
gDiffReadoutStats.c (added)
-
gDiffReadoutSubtract.c (added)
-
gDiffSetMasks.c (added)
-
gDiffThreshold.c (added)
-
gDiffVarianceRescale.c (added)
-
gDiffVersion.c (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
branches/eam_branches/gDiff.20200202/src/Makefile.am
r41260 r41413 8 8 gDiffVersion.c: gDiffVersionDefinitions.h 9 9 gDiffVersionDefinitions.h: gDiffVersionDefinitions.h.in FORCE 10 ../pslib-setsvnversion.pl GDIFF gDiffVersionDefinitions.h.in gDiffVersionDefinitions.h10 pslib-setsvnversion.pl GDIFF gDiffVersionDefinitions.h.in gDiffVersionDefinitions.h 11 11 FORCE: ; 12 12 … … 17 17 gDiff.c \ 18 18 gDiffArguments.c \ 19 gDiff Version.c\19 gDiffBackground.c \ 20 20 gDiffCamera.c \ 21 21 gDiffData.c \ 22 gDiffDefineOutput.c \ 22 23 gDiffErrorCodes.c \ 24 gDiffExit.c \ 25 gDiffExtras.c \ 26 gDiffFiles.c \ 27 gDiffFlagNeighbors.c \ 28 gDiffInputDetections.c \ 23 29 gDiffLoop.c \ 24 gDiffExit.c \ 25 gDiffFiles.c \ 26 gDiffDefineOutput.c 27 28 # gDiffBackground.c \ 29 # gDiffVarianceRescale.c \ 30 # gDiffInputDetections.c \ 31 # gDiffExtras.c \ 32 # gDiffFlagNeighbors.c \ 33 # gDiffMakePSF.c \ 34 # gDiffMatchPSFs.c \ 35 # gDiffMaskStats.c \ 36 # gDiffReadoutInverse.c \ 37 # gDiffReadoutJpeg.c \ 38 # gDiffReadoutPhotometry.c \ 39 # gDiffReadoutForcedPhot.c \ 40 # gDiffReadoutStats.c \ 41 # gDiffReadoutSubtract.c \ 42 # gDiffSetMasks.c \ 43 # gDiffThreshold.c 30 gDiffMakePSF.c \ 31 gDiffMaskStats.c \ 32 gDiffMatchPSFs.c \ 33 gDiffReadoutForcedPhot.c \ 34 gDiffReadoutInverse.c \ 35 gDiffReadoutJpeg.c \ 36 gDiffReadoutPhotometry.c \ 37 gDiffReadoutStats.c \ 38 gDiffReadoutSubtract.c \ 39 gDiffSetMasks.c \ 40 gDiffThreshold.c \ 41 gDiffVarianceRescale.c \ 42 gDiffVersion.c 44 43 45 44 noinst_HEADERS = \ -
branches/eam_branches/gDiff.20200202/src/gDiff.c
r41260 r41413 18 18 psphotErrorRegister(); 19 19 20 psExit exitValue = PS_EXIT_SUCCESS; // default value 21 20 22 gDiffData *data = NULL; // Processing data 21 23 pmConfig *config = pmConfigRead(&argc, argv, GDIFF_RECIPE); // Configuration … … 26 28 psphotInit(); // cannot return false 27 29 30 // the 'data' container carries pointers to all data structures used by gDiff 28 31 data = gDiffDataAlloc(config); 29 32 33 // parse the command-line options 30 34 if (!gDiffArguments(argc, argv, data)) { 31 35 psError(psErrorCodeLast(), false, "Error reading arguments."); … … 33 37 } 34 38 35 // determine the camera and setup the I/O file handles39 // determine the camera and setup the (many) I/O file handles 36 40 if (!gDiffCamera(data)) goto die; 37 41 … … 40 44 41 45 die: 42 { 43 psExit exitValue = gDiffExitCode(PS_EXIT_SUCCESS); // Exit code 46 exitValue = gDiffExitCode(exitValue); // Exit code 44 47 45 if (data && data->stats && data->statsFile) {46 psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_DIFF", 0, "Total time (sec)", psTimerMark("gDiff"));47 psString stats = psMetadataConfigFormat(data->stats); // Statistics to output48 if (!stats || strlen(stats) == 0) {49 psError(GDIFF_ERR_IO, false, "Unable to format statistics file");50 } else if (fprintf(data->statsFile, "%s", stats) != strlen(stats)) {51 psError(GDIFF_ERR_IO, true, "Unable to write statistics file");52 }53 psFree(stats);54 if (fclose(data->statsFile) == EOF) {55 psError(GDIFF_ERR_IO, true, "Unable to close statistics file");56 }57 data->statsFile = NULL;58 pmConfigRunFilenameAddWrite(data->config, "STATS", data->statsName);59 exitValue = gDiffExitCode(exitValue);60 }48 if (data && data->stats && data->statsFile) { 49 psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_DIFF", 0, "Total time (sec)", psTimerMark("gDiff")); 50 psString stats = psMetadataConfigFormat(data->stats); // Statistics to output 51 if (!stats || strlen(stats) == 0) { 52 psError(GDIFF_ERR_IO, false, "Unable to format statistics file"); 53 } else if (fprintf(data->statsFile, "%s", stats) != strlen(stats)) { 54 psError(GDIFF_ERR_IO, true, "Unable to write statistics file"); 55 } 56 psFree(stats); 57 if (fclose(data->statsFile) == EOF) { 58 psError(GDIFF_ERR_IO, true, "Unable to close statistics file"); 59 } 60 data->statsFile = NULL; 61 pmConfigRunFilenameAddWrite(data->config, "STATS", data->statsName); 62 exitValue = gDiffExitCode(exitValue); 63 } 61 64 62 if (config && !gDiffFilesIterateUp(config, GDIFF_FILES_ALL)) {63 psError(psErrorCodeLast(), false, "Unable to close files.");64 exitValue = gDiffExitCode(exitValue);65 pmFPAfileFreeSetStrict(false);66 }65 if (config && !gDiffFilesIterateUp(config, GDIFF_FILES_ALL)) { 66 psError(psErrorCodeLast(), false, "Unable to close files."); 67 exitValue = gDiffExitCode(exitValue); 68 pmFPAfileFreeSetStrict(false); 69 } 67 70 68 if (data) {69 psString dump_file = psMetadataLookupStr(NULL, data->config->arguments, "DUMP_CONFIG");70 if (dump_file) {71 if (!pmConfigDump(data->config, dump_file)) {72 psError(GDIFF_ERR_IO, false, "Unable to dump configuration.");73 exitValue = gDiffExitCode(exitValue);74 }75 }76 psFree(data);77 }71 if (data) { 72 psString dump_file = psMetadataLookupStr(NULL, data->config->arguments, "DUMP_CONFIG"); 73 if (dump_file) { 74 if (!pmConfigDump(data->config, dump_file)) { 75 psError(GDIFF_ERR_IO, false, "Unable to dump configuration."); 76 exitValue = gDiffExitCode(exitValue); 77 } 78 } 79 psFree(data); 80 } 78 81 79 psTrace("gDiff", 1, "Finished at %f sec\n", psTimerMark("gDiff"));80 psLogMsg("gDiff", PS_LOG_INFO, "Complete gDiff run: %f sec\n", psTimerMark("gDiff"));81 psTimerStop();82 psTrace("gDiff", 1, "Finished at %f sec\n", psTimerMark("gDiff")); 83 psLogMsg("gDiff", PS_LOG_INFO, "Complete gDiff run: %f sec\n", psTimerMark("gDiff")); 84 psTimerStop(); 82 85 83 pmVisualClose(); // close plot windows, if -visual is set84 pmModelClassCleanup();85 pmConfigDone();86 pmVisualCleanup ();87 psLibFinalize();86 pmVisualClose(); // close plot windows, if -visual is set 87 pmModelClassCleanup(); 88 pmConfigDone(); 89 pmVisualCleanup (); 90 psLibFinalize(); 88 91 89 fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, stdout, false), "gDiff"); 90 exitValue = gDiffExitCode(exitValue); 91 exit(exitValue); 92 } 92 fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, stdout, false), "gDiff"); 93 exitValue = gDiffExitCode(exitValue); 94 exit(exitValue); 93 95 } -
branches/eam_branches/gDiff.20200202/src/gDiff.h
r41260 r41413 21 21 #include <pslib.h> 22 22 #include <psmodules.h> 23 #include <ppStats.h> 23 24 #include <psphot.h> 24 25 -
branches/eam_branches/gDiff.20200202/src/gDiffArguments.c
r41260 r41413 1 /** @file gDiffArguments.c 2 * @brief: image subtraction based on matched Gaussians profiles 3 * @author Eugene Magnier @ IfA 4 * @version $Revision: 1.12 $ 5 * @date $Date: 2020-01-25 $ 6 * Copyright 2020 Institute for Astronomy, University of Hawaii 7 */ 8 1 9 #include "gDiff.h" 2 10 … … 60 68 61 69 psMetadata *arguments = config->arguments; // Command-line arguments 62 psMetadataAddStr (arguments, PS_LIST_TAIL, "-inimage", 0, "Input image", NULL);63 psMetadataAddStr (arguments, PS_LIST_TAIL, "-inmask", 0, "Input mask image", NULL);64 psMetadataAddStr (arguments, PS_LIST_TAIL, "-invariance", 0, "Input variance image", NULL);65 psMetadataAddStr (arguments, PS_LIST_TAIL, "-insources", 0, "Input source list", NULL);66 psMetadataAddStr (arguments, PS_LIST_TAIL, "-refimage", 0, "Reference image", NULL);67 psMetadataAddStr (arguments, PS_LIST_TAIL, "-refmask", 0, "Reference mask image", NULL);68 psMetadataAddStr (arguments, PS_LIST_TAIL, "-refvariance", 0, "Reference variance image", NULL);69 psMetadataAddStr (arguments, PS_LIST_TAIL, "-refsources", 0, "Reference source list", NULL);70 psMetadataAddStr (arguments, PS_LIST_TAIL, "-kernel", 0, "Precalculated kernel to apply", NULL);71 psMetadataAddStr (arguments, PS_LIST_TAIL, "-stats", 0, "Statistics file", NULL);72 psMetadataAddStr (arguments,PS_LIST_TAIL, "-stamps", 0, "Stamps filename; x,y on each line", NULL);73 psMetadataAddS32 (arguments, PS_LIST_TAIL, "-convolve", 0, "Image to convolve [1 or 2]", 0);70 psMetadataAddStr (arguments, PS_LIST_TAIL, "-inimage", 0, "Input image", NULL); 71 psMetadataAddStr (arguments, PS_LIST_TAIL, "-inmask", 0, "Input mask image", NULL); 72 psMetadataAddStr (arguments, PS_LIST_TAIL, "-invariance", 0, "Input variance image", NULL); 73 psMetadataAddStr (arguments, PS_LIST_TAIL, "-insources", 0, "Input source list", NULL); 74 psMetadataAddStr (arguments, PS_LIST_TAIL, "-refimage", 0, "Reference image", NULL); 75 psMetadataAddStr (arguments, PS_LIST_TAIL, "-refmask", 0, "Reference mask image", NULL); 76 psMetadataAddStr (arguments, PS_LIST_TAIL, "-refvariance", 0, "Reference variance image", NULL); 77 psMetadataAddStr (arguments, PS_LIST_TAIL, "-refsources", 0, "Reference source list", NULL); 78 psMetadataAddStr (arguments, PS_LIST_TAIL, "-kernel", 0, "Precalculated kernel to apply", NULL); 79 psMetadataAddStr (arguments, PS_LIST_TAIL, "-stats", 0, "Statistics file", NULL); 80 psMetadataAddStr (arguments, PS_LIST_TAIL, "-stamps", 0, "Stamps filename; x,y on each line", NULL); 81 psMetadataAddS32 (arguments, PS_LIST_TAIL, "-convolve", 0, "Image to convolve [1 or 2]", 0); 74 82 psMetadataAddBool(arguments, PS_LIST_TAIL, "-photometry", 0, "Perform photometry?", NULL); 75 83 psMetadataAddBool(arguments, PS_LIST_TAIL, "-forced-phot", 0, "Perform forced photometry?", NULL); 76 84 psMetadataAddBool(arguments, PS_LIST_TAIL, "-forced-input1", 0, "Perform forced photometry?", NULL); 77 85 psMetadataAddBool(arguments, PS_LIST_TAIL, "-forced-input2", 0, "Perform forced photometry?", NULL); 78 psMetadataAddF32 (arguments, PS_LIST_TAIL, "-zp", 0, "Zero point for photometry", NAN);86 psMetadataAddF32 (arguments, PS_LIST_TAIL, "-zp", 0, "Zero point for photometry", NAN); 79 87 psMetadataAddBool(arguments, PS_LIST_TAIL, "-inverse", 0, "Generate inverse subtractions?", false); 80 88 psMetadataAddBool(arguments, PS_LIST_TAIL, "-save-inconv", 0, "Save input convolved images?", false); -
branches/eam_branches/gDiff.20200202/src/gDiffCamera.c
r41260 r41413 1 /** @file gDiffCamera.c 2 * @brief: image subtraction based on matched Gaussians profiles 3 * @author Eugene Magnier @ IfA 4 * @version $Revision: 1.12 $ 5 * @date $Date: 2020-01-25 $ 6 * Copyright 2020 Institute for Astronomy, University of Hawaii 7 */ 8 1 9 #include "gDiff.h" 2 10 … … 207 215 } 208 216 217 /* in ppSub, there is a recipe value 'PHOTOMETRY' which is supposed to enable / disable photometry analysis. 218 This is set to FALSE in the recipes by default. 219 220 In the script diff_skycell.pl, there is a global variable 'do_photom' which is TRUE by default, 221 but is set to FALSE if a) the reduction class is NOCONVDIFF or b) the run_state is not new (i.e, 'update' mode) 222 (the script changes the call in an inconsistent way, if do_photom is true, it appends -photometry, but 223 if it is false, it appends -Db PHOTOMETRY FALSE. this could be made consistent. 224 (also note that it is bad practice to have hard-wired reduction classes) 225 226 But both of these situations are ignored because the test in ppSubCamera for data->photometry 227 (set if PHOTOMETRY is true or -photometry is supplied) is ignored and the photometric analysis 228 is always performed. 229 230 */ 231 209 232 // Now that the camera has been determined, we can read the recipe 210 233 psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, GDIFF_RECIPE); // Recipe for ppSim … … 380 403 } 381 404 382 // psPhot input 383 if (data->photometry || 1) { 384 psphotModelClassInit(); // load implementation-specific models 385 386 pmFPAfile *psphot = pmFPAfileDefineFromFPA(config, output->fpa, 1, 1, "PSPHOT.INPUT"); 387 if (!psphot) { 388 psError(psErrorCodeLast(), false, "Failed to build FPA from PSPHOT.INPUT"); 389 return false; 390 } 391 if (psphot->type != PM_FPA_FILE_IMAGE) { 392 psError(psErrorCodeLast(), true, "PSPHOT.INPUT is not of type IMAGE"); 393 return false; 394 } 395 pmFPAfileActivate(config->files, false, "PSPHOT.INPUT"); 396 397 // Internal file for getting the PSF from the minuend 398 pmFPAfile *psf = pmFPAfileDefineFromFPA(config, output->fpa, 1, 1, "PSPHOT.PSF.LOAD"); 399 if (!psf) { 400 psError(psErrorCodeLast(), false, "Failed to build FPA from PSPHOT.PSF.LOAD"); 401 return false; 402 } 403 if (psf->type != PM_FPA_FILE_PSF) { 404 psError(psErrorCodeLast(), true, "PSPHOT.PSF.LOAD is not of type PSF"); 405 return false; 406 } 407 pmFPAfileActivate(config->files, false, "PSPHOT.PSF.LOAD"); 408 409 if (!psphotDefineFiles(config, psphot)) { 410 psError(psErrorCodeLast(), false, "Unable to set up psphot files."); 411 return false; 412 } 413 414 // Deactivate psphot output sources --- we want to define output source files of our own 415 pmFPAfile *psphotOutput = pmFPAfileSelectSingle(config->files, "PSPHOT.OUTPUT", 0); 416 psphotOutput->save = false; // this one should NOT be set 417 418 pmFPAfile *outSources = defineOutputFile(config, output, false, "GDIFF.OUTPUT.SOURCES", 419 PM_FPA_FILE_CMF); 420 if (!outSources) { 421 psError(psErrorCodeLast(), false, "Unable to set up output source file."); 422 return false; 423 } 424 checkFileruleFileSave(outSources, config); 425 426 if (data->inverse) { 427 pmFPAfile *invSources = defineOutputFile(config, inverse, false, "GDIFF.INVERSE.SOURCES", 428 PM_FPA_FILE_CMF); 429 if (!invSources) { 430 psError(psErrorCodeLast(), false, "Unable to set up inverse source file."); 431 return false; 432 } 433 checkFileruleFileSave(invSources, config); 434 } 435 436 // files need to do the forced photometry on the positions of sources in the positive images 437 if (data->forcedPhot1) { 438 // this pmFPAfile is used to carry sources detected in the positive image #1 439 pmFPAfile *posSources1 = defineOutputFile(config, input, true, "GDIFF.POS1.SOURCES", PM_FPA_FILE_CMF); 440 if (!posSources1) { 441 psError(psErrorCodeLast(), false, "Unable to set up forced source file."); 442 return false; 443 } 444 checkFileruleFileSave(posSources1, config); 445 446 // this pmFPAfile is used to carry sources detected in the diff image @ the positions from positive image #1 447 pmFPAfile *frcSources1 = defineOutputFile(config, input, true, "GDIFF.FORCED1.SOURCES", PM_FPA_FILE_CMF); 448 if (!frcSources1) { 449 psError(psErrorCodeLast(), false, "Unable to set up forced source file."); 450 return false; 451 } 452 checkFileruleFileSave(frcSources1, config); 405 /* in ppSub, there is a recipe value 'PHOTOMETRY' which is supposed to enable / disable photometry analysis. 406 This is set to FALSE in the recipes by default. 407 408 In the script diff_skycell.pl, there is a global variable 'do_photom' which is TRUE by default, 409 but is set to FALSE if a) the reduction class is NOCONVDIFF or b) the run_state is not new (i.e, 'update' mode) 410 (the script changes the call in an inconsistent way, if do_photom is true, it appends -photometry, but 411 if it is false, it appends -Db PHOTOMETRY FALSE. this could be made consistent. 412 (also note that it is bad practice to have hard-wired reduction classes) 413 414 But both of these situations are ignored because the test in ppSubCamera for data->photometry 415 (set if PHOTOMETRY is true or -photometry is supplied) is ignored and the photometric analysis 416 is always performed. 417 418 Should we rationalize this? 419 */ 420 421 // define the psphot I/O files 422 psphotModelClassInit(); // load implementation-specific models 423 424 pmFPAfile *psphot = pmFPAfileDefineFromFPA(config, output->fpa, 1, 1, "PSPHOT.INPUT"); 425 if (!psphot) { 426 psError(psErrorCodeLast(), false, "Failed to build FPA from PSPHOT.INPUT"); 427 return false; 428 } 429 if (psphot->type != PM_FPA_FILE_IMAGE) { 430 psError(psErrorCodeLast(), true, "PSPHOT.INPUT is not of type IMAGE"); 431 return false; 432 } 433 pmFPAfileActivate(config->files, false, "PSPHOT.INPUT"); 434 435 // Internal file for getting the PSF from the minuend 436 pmFPAfile *psf = pmFPAfileDefineFromFPA(config, output->fpa, 1, 1, "PSPHOT.PSF.LOAD"); 437 if (!psf) { 438 psError(psErrorCodeLast(), false, "Failed to build FPA from PSPHOT.PSF.LOAD"); 439 return false; 440 } 441 if (psf->type != PM_FPA_FILE_PSF) { 442 psError(psErrorCodeLast(), true, "PSPHOT.PSF.LOAD is not of type PSF"); 443 return false; 444 } 445 pmFPAfileActivate(config->files, false, "PSPHOT.PSF.LOAD"); 446 447 if (!psphotDefineFiles(config, psphot)) { 448 psError(psErrorCodeLast(), false, "Unable to set up psphot files."); 449 return false; 450 } 451 452 // Deactivate psphot output sources --- we want to define output source files of our own 453 pmFPAfile *psphotOutput = pmFPAfileSelectSingle(config->files, "PSPHOT.OUTPUT", 0); 454 psphotOutput->save = false; // this one should NOT be set 455 456 pmFPAfile *outSources = defineOutputFile(config, output, false, "GDIFF.OUTPUT.SOURCES", 457 PM_FPA_FILE_CMF); 458 if (!outSources) { 459 psError(psErrorCodeLast(), false, "Unable to set up output source file."); 460 return false; 461 } 462 checkFileruleFileSave(outSources, config); 463 464 if (data->inverse) { 465 pmFPAfile *invSources = defineOutputFile(config, inverse, false, "GDIFF.INVERSE.SOURCES", 466 PM_FPA_FILE_CMF); 467 if (!invSources) { 468 psError(psErrorCodeLast(), false, "Unable to set up inverse source file."); 469 return false; 453 470 } 454 455 if (data->forcedPhot2) { 456 // this pmFPAfile is used to carry sources detected in the positive image #2 457 pmFPAfile *posSources2 = defineOutputFile(config, ref, true, "GDIFF.POS2.SOURCES", PM_FPA_FILE_CMF); 458 if (!posSources2) { 459 psError(psErrorCodeLast(), false, "Unable to set up forced source file."); 460 return false; 461 } 462 checkFileruleFileSave(posSources2, config); 463 464 // this pmFPAfile is used to carry sources detected in the diff image @ the positions from positive image #2 465 pmFPAfile *frcSources2 = defineOutputFile(config, ref, true, "GDIFF.FORCED2.SOURCES", PM_FPA_FILE_CMF); 466 if (!frcSources2) { 467 psError(psErrorCodeLast(), false, "Unable to set up forced source file."); 468 return false; 469 } 470 checkFileruleFileSave(frcSources2, config); 471 } 471 checkFileruleFileSave(invSources, config); 472 } 473 474 // files need to do the forced photometry on the positions of sources in the positive images 475 if (data->forcedPhot1) { 476 // this pmFPAfile is used to carry sources detected in the positive image #1 477 pmFPAfile *posSources1 = defineOutputFile(config, input, true, "GDIFF.POS1.SOURCES", PM_FPA_FILE_CMF); 478 if (!posSources1) { 479 psError(psErrorCodeLast(), false, "Unable to set up forced source file."); 480 return false; 481 } 482 checkFileruleFileSave(posSources1, config); 483 484 // this pmFPAfile is used to carry sources detected in the diff image @ the positions from positive image #1 485 pmFPAfile *frcSources1 = defineOutputFile(config, input, true, "GDIFF.FORCED1.SOURCES", PM_FPA_FILE_CMF); 486 if (!frcSources1) { 487 psError(psErrorCodeLast(), false, "Unable to set up forced source file."); 488 return false; 489 } 490 checkFileruleFileSave(frcSources1, config); 491 } 492 493 if (data->forcedPhot2) { 494 // this pmFPAfile is used to carry sources detected in the positive image #2 495 pmFPAfile *posSources2 = defineOutputFile(config, ref, true, "GDIFF.POS2.SOURCES", PM_FPA_FILE_CMF); 496 if (!posSources2) { 497 psError(psErrorCodeLast(), false, "Unable to set up forced source file."); 498 return false; 499 } 500 checkFileruleFileSave(posSources2, config); 501 502 // this pmFPAfile is used to carry sources detected in the diff image @ the positions from positive image #2 503 pmFPAfile *frcSources2 = defineOutputFile(config, ref, true, "GDIFF.FORCED2.SOURCES", PM_FPA_FILE_CMF); 504 if (!frcSources2) { 505 psError(psErrorCodeLast(), false, "Unable to set up forced source file."); 506 return false; 507 } 508 checkFileruleFileSave(frcSources2, config); 472 509 } 473 510 -
branches/eam_branches/gDiff.20200202/src/gDiffData.c
r41260 r41413 1 /** @file gDiffData.c 2 * @brief: image subtraction based on matched Gaussians profiles 3 * @author Eugene Magnier @ IfA 4 * @version $Revision: 1.12 $ 5 * @date $Date: 2020-01-25 $ 6 * Copyright 2020 Institute for Astronomy, University of Hawaii 7 */ 8 1 9 #include "gDiff.h" 2 10 -
branches/eam_branches/gDiff.20200202/src/gDiffDefineOutput.c
r41260 r41413 1 1 /** @file gDiffDefineOutput.c 2 * 3 * @brief 4 * 5 * @ingroup gDiff 6 * 7 * @author IfA 8 * @version $Revision: 1.3 $ $Name: not supported by cvs2svn $ 9 * @date $Date: 2009-02-18 00:31:20 $ 10 * Copyright 2009 Institute for Astronomy, University of Hawaii 2 * @brief: image subtraction based on matched Gaussians profiles 3 * @author Eugene Magnier @ IfA 4 * @version $Revision: 1.12 $ 5 * @date $Date: 2020-01-25 $ 6 * Copyright 2020 Institute for Astronomy, University of Hawaii 11 7 */ 12 13 #ifdef HAVE_CONFIG_H14 #include <config.h>15 #endif16 17 #include <stdio.h>18 #include <libgen.h>19 #include <pslib.h>20 #include <psmodules.h>21 8 22 9 #include "gDiff.h" -
branches/eam_branches/gDiff.20200202/src/gDiffErrorCodes.dat
r41261 r41413 3 3 # 4 4 BASE = 14000 First value we use; lower values belong to psLib 5 UNKNOWN Unknown ppSuberror code5 UNKNOWN Unknown gDiff error code 6 6 NOT_IMPLEMENTED Desired feature is not yet implemented 7 7 ARGUMENTS Incorrect arguments -
branches/eam_branches/gDiff.20200202/src/gDiffExit.c
r41260 r41413 1 #include <stdio.h> 2 #include <pslib.h> 3 #include <psphotErrorCodes.h> 1 /** @file gDiffExit.c 2 * @brief: image subtraction based on matched Gaussians profiles 3 * @author Eugene Magnier @ IfA 4 * @version $Revision: 1.12 $ 5 * @date $Date: 2020-01-25 $ 6 * Copyright 2020 Institute for Astronomy, University of Hawaii 7 */ 4 8 5 9 #include "gDiff.h" … … 16 20 // not include the possible psphot values 17 21 18 // for now, to get around this, we just use an int for the switch 22 // thus, we cannot do this: 23 // psErrorCode errorCode = psErrorCodeLast(); // Error code 19 24 20 // psErrorCode errorCode = psErrorCodeLast(); // Error code25 // for now, to get around this, we just use an int for the switch: 21 26 int errorCode = psErrorCodeLast(); // Error code 27 22 28 if (errorCode != PS_ERR_NONE) { 23 29 pmFPAfileFreeSetStrict(false); -
branches/eam_branches/gDiff.20200202/src/gDiffFiles.c
r41260 r41413 1 #ifdef HAVE_CONFIG_H 2 #include <config.h> 3 #endif 4 5 #include <stdio.h> 6 #include <pslib.h> 7 #include <psmodules.h> 8 #include <psphot.h> 1 /** @file gDiffFiles.c 2 * @brief: image subtraction based on matched Gaussians profiles 3 * @author Eugene Magnier @ IfA 4 * @version $Revision: 1.12 $ 5 * @date $Date: 2020-01-25 $ 6 * Copyright 2020 Institute for Astronomy, University of Hawaii 7 */ 9 8 10 9 #include "gDiff.h" … … 12 11 13 12 // Input files 14 static const char *inputFiles[] = { "GDIFF.INPUT", "GDIFF.INPUT.MASK", "GDIFF.INPUT.VARIANCE", 15 "GDIFF.INPUT.SOURCES", "GDIFF.REF", "GDIFF.REF.MASK", 16 "GDIFF.REF.VARIANCE", "GDIFF.REF.SOURCES", NULL }; 13 static const char *inputFiles[] = { "GDIFF.INPUT", 14 "GDIFF.INPUT.MASK", 15 "GDIFF.INPUT.VARIANCE", 16 "GDIFF.INPUT.SOURCES", 17 "GDIFF.REF", 18 "GDIFF.REF.MASK", 19 "GDIFF.REF.VARIANCE", 20 "GDIFF.REF.SOURCES", NULL }; 17 21 18 22 // Convolved files 19 23 static const char *convFiles[] = { "GDIFF.INPUT.CONV", "GDIFF.INPUT.CONV.MASK", "GDIFF.INPUT.CONV.VARIANCE", 20 "GDIFF.REF.CONV", "GDIFF.REF.CONV.MASK","GDIFF.REF.CONV.VARIANCE",24 "GDIFF.REF.CONV", "GDIFF.REF.CONV.MASK", "GDIFF.REF.CONV.VARIANCE", 21 25 NULL }; 22 26 … … 234 238 return true; 235 239 } 240 -
branches/eam_branches/gDiff.20200202/src/gDiffLoop.c
r41260 r41413 1 /** @file gDiffLoop.c 2 * @brief: image subtraction based on matched Gaussians profiles 3 * @author Eugene Magnier @ IfA 4 * @version $Revision: 1.12 $ 5 * @date $Date: 2020-01-25 $ 6 * Copyright 2020 Institute for Astronomy, University of Hawaii 7 */ 8 1 9 #include "gDiff.h" 10 11 # define HDUCHECK(A) \ 12 fprintf (stderr, "%s\n", A); \ 13 listhdu (config, "GDIFF.OUTPUT"); listhdu (config, "PSPHOT.INPUT"); 14 15 // internal functions: 16 bool dumpout(pmConfig *config, char *name); 17 bool listhdu(pmConfig *config, char *name); 18 19 bool gDiffLoop(gDiffData *data) 20 { 21 bool mdok = false; 22 bool success = true; 23 24 psAssert(data, "Require processing data"); 25 pmConfig *config = data->config; // Configuration 26 psAssert(config, "Require configuration."); 27 28 HDUCHECK ("0a: start gDiffLoop"); 29 30 psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, GDIFF_RECIPE); // Recipe for ppSim 31 psAssert(recipe, "We checked this earlier, so it should be here."); 32 bool noConvolve = psMetadataLookupBool(&mdok, recipe, "NOCONVOLVE"); // Do not use convolved images. 33 34 pmConfigCamerasCull(config, NULL); 35 pmConfigRecipesCull(config, "GDIFF,PPSTATS,PSPHOT,PSASTRO,MASKS,JPEG"); 36 37 pmFPAfile *input = psMetadataLookupPtr(NULL, config->files, "GDIFF.INPUT"); 38 pmFPAfile *reference = psMetadataLookupPtr(NULL, config->files, "GDIFF.REF"); 39 pmFPAfile *output = psMetadataLookupPtr(NULL, config->files, "GDIFF.OUTPUT"); 40 psAssert(input && reference && output, "Require files"); 41 42 // load input files : we are required to load a complete set of files that describe the input. 43 // we are not able to process subdivisions of an image (e.g., cells vs chips) 44 if (!gDiffFilesIterateDown(config, GDIFF_FILES_INPUT | GDIFF_FILES_CONV)) { 45 psError(GDIFF_ERR_IO, false, "Unable to load files."); 46 return false; 47 } 48 49 psTimerStart("GDIFF_MATCH"); 50 51 HDUCHECK ("0b: gDiffSetMasks"); 52 53 if (!gDiffSetMasks(config)) { 54 psError(psErrorCodeLast(), false, "Unable to set masks."); 55 return false; 56 } 57 58 HDUCHECK ("0b.1: gDiffInputDetections"); 59 if (data->forcedPhot1) { 60 bool foundDetections = false; 61 if (!gDiffInputDetections(&foundDetections, "GDIFF.POS1.SOURCES", "GDIFF.INPUT", data)) { 62 psError(psErrorCodeLast(), false, "Unable to measure positive detections (1)"); 63 success = false; 64 goto ESCAPE; 65 } 66 // if nothing was found, don't bother doing the forced photometry below 67 if (!foundDetections) { 68 psWarning ("no sources found in positive image 1, skipping forced photometry"); 69 data->forcedPhot1 = false; 70 } 71 } 72 HDUCHECK ("0b.2: gDiffInputDetections"); 73 if (data->forcedPhot2) { 74 // XXX why is this only done for Phot2 and not Phot1 ?? 75 76 // Change the recipe to use a higher nsigma limit and quit after pass1 77 psMetadata *psphotRecipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE); 78 79 psF32 nsigma_peak_save = psMetadataLookupF32 (NULL, psphotRecipe, "PEAKS_NSIGMA_LIMIT"); 80 char *breakPt_save = psMetadataLookupStr (NULL, psphotRecipe, "BREAK_POINT"); 81 82 psF32 pos2_nsigma_peak = psMetadataLookupF32 (&mdok, psphotRecipe, "PEAKS_POS2_NSIGMA_LIMIT"); 83 if (!mdok) { 84 psWarning("PEAKS_POS2_NSIGMA_LIMIT not found in psphotRecipe. Will use 25.\n"); 85 pos2_nsigma_peak = 25.; 86 } 87 psMetadataAddF32(psphotRecipe, PS_LIST_TAIL, "PEAKS_NSIGMA_LIMIT", PS_META_REPLACE, "", pos2_nsigma_peak); 88 psMetadataAddStr(psphotRecipe, PS_LIST_TAIL, "BREAK_POINT", PS_META_REPLACE, "", "PASS1"); 89 90 bool foundDetections = false; 91 if (!gDiffInputDetections(&foundDetections, "GDIFF.POS2.SOURCES", "GDIFF.REF", data)) { 92 psError(psErrorCodeLast(), false, "Unable to measure positive detections (2)"); 93 psMetadataAddF32(psphotRecipe, PS_LIST_TAIL, "PEAKS_NSIGMA_LIMIT", PS_META_REPLACE, "", nsigma_peak_save); 94 psMetadataAddStr(psphotRecipe, PS_LIST_TAIL, "BREAK_POINT", PS_META_REPLACE, "", breakPt_save); 95 success = false; 96 goto ESCAPE; 97 } 98 psMetadataAddF32(psphotRecipe, PS_LIST_TAIL, "PEAKS_NSIGMA_LIMIT", PS_META_REPLACE, "", nsigma_peak_save); 99 psMetadataAddStr(psphotRecipe, PS_LIST_TAIL, "BREAK_POINT", PS_META_REPLACE, "", breakPt_save); 100 // if nothing was found, don't bother doing the forced photometry below 101 if (!foundDetections) { 102 psWarning ("no sources found in positive image 2, skipping forced photometry"); 103 data->forcedPhot2 = false; 104 } 105 } 106 107 // XXX if it exists, use the POS1, POS2 successs for the FWHMs 108 109 /**** The PSF matching process is performed in this function below. **** 110 **** This is the core of the PSF-matched image subtraction code. ****/ 111 112 HDUCHECK ("1: before gDiffMatchPSFs"); 113 114 if (!gDiffMatchPSFs(data)) { 115 psError(psErrorCodeLast(), false, "Unable to match PSFs."); 116 success = false; 117 goto ESCAPE; 118 } 119 if (data->quality) { 120 // Can't do anything at all 121 success = false; 122 goto ESCAPE; 123 } 124 125 // the functions in this block are not critical for initial testing, so block them out for now 126 127 // generate the residual stamp grid for visualization 128 if (!gDiffResidualSampleJpeg(config)) { 129 psError(psErrorCodeLast(), false, "Unable to update."); 130 success = false; 131 goto ESCAPE; 132 } 133 134 // XXX add in a positive image detection step here (if needed) 135 // XXX EAM what does this mean??? 136 137 psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_MATCH", 0, "Time to match PSFs", psTimerClear("GDIFF_MATCH")); 138 139 // Close input files (freeing up space) : for the noConvolve case, we cannot close the 140 // inputs since we will use them for subtraction below. wait until later to do the work 141 if (!noConvolve) { 142 if (!gDiffFilesIterateUp(config, GDIFF_FILES_INPUT)) { 143 psError(GDIFF_ERR_IO, false, "Unable to close input files."); 144 success = false; 145 goto ESCAPE; 146 } 147 } 148 if (!gDiffLowThreshold(data)) { 149 psError(psErrorCodeLast(), false, "Unable to threshold images."); 150 success = false; 151 goto ESCAPE; 152 } 153 154 HDUCHECK ("2: after gDiffLowThreshold"); 155 156 // Set up subtraction files 157 if (!gDiffFilesIterateDown(config, GDIFF_FILES_SUB)) { 158 psError(GDIFF_ERR_IO, false, "Unable to set up subtraction files."); 159 success = false; 160 goto ESCAPE; 161 } 162 163 if (!gDiffDefineOutput("GDIFF.OUTPUT", config)) { 164 psError(psErrorCodeLast(), false, "Unable to define output."); 165 success = false; 166 goto ESCAPE; 167 } 168 169 HDUCHECK ("2a: gDiffMakePSF"); 170 171 // the functions in this block are not critical for initial testing, so block them out for now 172 if (!data->quality && !gDiffMakePSF(data)) { 173 psError(psErrorCodeLast(), false, "Unable to generate PSF."); 174 success = false; 175 goto ESCAPE; 176 } 177 178 // Now we've got a PSF, blow away detections in case they're confused with real output detections 179 { 180 pmFPAview *view = gDiffViewReadout(); // View to readout 181 pmReadout *out = pmFPAfileThisReadout(config->files, view, "GDIFF.OUTPUT"); 182 if (psMetadataLookup(out->analysis, "PSPHOT.DETECTIONS")) { 183 psMetadataRemoveKey(out->analysis, "PSPHOT.DETECTIONS"); 184 } 185 psFree(view); 186 } 187 188 HDUCHECK ("3: before gDiffReadoutSubtract"); 189 190 if (!gDiffReadoutSubtract(config)) { 191 psError(psErrorCodeLast(), false, "Unable to subtract images."); 192 success = false; 193 goto ESCAPE; 194 } 195 196 // Close convolved files 197 if (!gDiffFilesIterateUp(config, GDIFF_FILES_PSF | GDIFF_FILES_CONV)) { 198 psError(GDIFF_ERR_IO, false, "Unable to close input files."); 199 success = false; 200 goto ESCAPE; 201 } 202 // Close unconvolved input files if we did the noConvolve subtraction 203 // (otherwise they are closed above) 204 if (noConvolve && !gDiffFilesIterateUp(config, GDIFF_FILES_INPUT)) { 205 psError(GDIFF_ERR_IO, false, "Unable to close input files."); 206 success = false; 207 goto ESCAPE; 208 } 209 210 HDUCHECK ("4: before gDiffBackground"); 211 212 // Higher order background subtraction using psphot 213 if (!gDiffBackground(config)) { 214 psError(psErrorCodeLast(), false, "Unable to subtract background."); 215 success = false; 216 goto ESCAPE; 217 } 218 // dumpout(config, "diff.2b.fits"); 219 220 HDUCHECK ("5: before gDiffVarianceRescale"); 221 222 // Perform Variance correction (rescale within a modest range) 223 if (!gDiffVarianceRescale(config, data)) { 224 psError(psErrorCodeLast(), false, "Unable to rescale variance."); 225 success = false; 226 goto ESCAPE; 227 } 228 // dumpout(config, "diff.2c.fits"); 229 230 if (data->quality) { 231 // Done all we can do up to this point 232 success = false; 233 goto ESCAPE; 234 } 235 236 if (!gDiffFilesIterateDown(config, GDIFF_FILES_PHOT_SUB)) { 237 psError(GDIFF_ERR_IO, false, "Unable to set up photometry files."); 238 success = false; 239 goto ESCAPE; 240 } 241 242 HDUCHECK ("6: before gDiffReadoutPhotometry"); 243 244 if (!data->quality && !gDiffReadoutPhotometry("GDIFF.OUTPUT", data)) { 245 psError(psErrorCodeLast(), false, "Unable to perform photometry."); 246 success = false; 247 goto ESCAPE; 248 } 249 // dumpout(config, "diff.3.fits"); 250 251 HDUCHECK ("7: before gDiffReadoutForcedPhot (1)"); 252 253 // forced photometry for positive image 1 254 if (data->forcedPhot1 && !data->quality) { 255 if (!gDiffReadoutForcedPhot("GDIFF.FORCED1.SOURCES", "GDIFF.OUTPUT", "GDIFF.POS1.SOURCES", data)) { 256 psError(psErrorCodeLast(), false, "Unable to perform photometry."); 257 success = false; 258 goto ESCAPE; 259 } 260 } 261 262 // forced photometry for positive image 2 263 if (data->forcedPhot2 && !data->quality) { 264 if (!gDiffReadoutForcedPhot("GDIFF.FORCED2.SOURCES", "GDIFF.OUTPUT", "GDIFF.POS2.SOURCES", data)) { 265 psError(psErrorCodeLast(), false, "Unable to perform photometry."); 266 success = false; 267 goto ESCAPE; 268 } 269 } 270 271 if (!gDiffFilesIterateUp(config, GDIFF_FILES_PHOT_SUB)) { 272 psError(GDIFF_ERR_IO, false, "Unable to set up photometry files."); 273 success = false; 274 goto ESCAPE; 275 } 276 277 // Perform statistics on the cell 278 if (!gDiffReadoutStats(data)) { 279 psError(psErrorCodeLast(), false, "Unable to collect statistics"); 280 success = false; 281 goto ESCAPE; 282 } 283 // Do Mask Stats 284 { 285 pmFPAview *view = gDiffViewReadout(); // View to readout 286 if (!gDiffMaskStats(config, view,data->stats)) { 287 psError(psErrorCodeLast(), false, "Unable to generate mask statistics"); 288 success = false; 289 psFree(view); 290 goto ESCAPE; 291 } 292 psFree(view); 293 } 294 // dumpout(config, "diff.4.fits"); 295 296 // generate the binned image used to write the jpeg 297 if (!gDiffReadoutJpeg(config)) { 298 psError(psErrorCodeLast(), false, "Unable to update."); 299 success = false; 300 goto ESCAPE; 301 } 302 303 if (data->inverse) { 304 // Set up inverse subtraction files 305 if (!gDiffFilesIterateDown(config, GDIFF_FILES_INV)) { 306 psError(GDIFF_ERR_IO, false, "Unable to set up inverse files."); 307 success = false; 308 goto ESCAPE; 309 } 310 311 if (data->inverse && !gDiffDefineOutput("GDIFF.INVERSE", config)) { 312 psError(psErrorCodeLast(), false, "Unable to define inverse."); 313 success = false; 314 goto ESCAPE; 315 } 316 317 if (!gDiffReadoutInverse(config)) { 318 psError(psErrorCodeLast(), false, "Unable to invert images."); 319 success = false; 320 goto ESCAPE; 321 } 322 323 //MEH -- need to also clear out detections or inv.cmf corrupted on update 324 pmFPAview *viewinv = gDiffViewReadout(); // View to readout 325 pmReadout *inv = pmFPAfileThisReadout(config->files, viewinv, "GDIFF.INVERSE"); 326 if (psMetadataLookup(inv->analysis, "PSPHOT.DETECTIONS")) { 327 psMetadataRemoveKey(inv->analysis, "PSPHOT.DETECTIONS"); 328 } 329 psFree(viewinv); 330 331 // Close subtraction files and open inverse photometry files 332 if (!gDiffFilesIterateUp(config, GDIFF_FILES_SUB)) { 333 psError(GDIFF_ERR_IO, false, "Unable to close subtraction files."); 334 success = false; 335 goto ESCAPE; 336 } 337 338 if (!gDiffFilesIterateDown(config, GDIFF_FILES_PHOT_INV)) { 339 psError(GDIFF_ERR_IO, false, "Unable to set up inverse files."); 340 success = false; 341 goto ESCAPE; 342 } 343 344 if (!data->quality && !gDiffReadoutPhotometry("GDIFF.INVERSE", data)) { 345 psError(psErrorCodeLast(), false, "Unable to perform photometry."); 346 success = false; 347 goto ESCAPE; 348 } 349 350 // Close inverse subtraction files 351 if (!gDiffFilesIterateUp(config, GDIFF_FILES_INV | GDIFF_FILES_PHOT_INV)) { 352 psError(GDIFF_ERR_IO, false, "Unable to close subtraction files."); 353 success = false; 354 goto ESCAPE; 355 } 356 } else { 357 // Close subtraction files 358 // dumpout(config, "diff.5.fits"); 359 if (!gDiffFilesIterateUp(config, GDIFF_FILES_SUB)) { 360 psError(GDIFF_ERR_IO, false, "Unable to close subtraction files."); 361 success = false; 362 goto ESCAPE; 363 } 364 } 365 366 ESCAPE: 367 pmFPAfileDropInternal(config->files, "PSPHOT.BACKGND"); 368 pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL"); 369 pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV"); 370 371 return success; 372 } 373 374 /**** internal functions ****/ 2 375 3 376 bool dumpout(pmConfig *config, char *name) … … 10 383 } 11 384 12 13 bool gDiffLoop(gDiffData *data) 385 bool listhdu(pmConfig *config, char *name) 14 386 { 15 bool mdok = false; 16 bool success = true; 17 18 psAssert(data, "Require processing data"); 19 pmConfig *config = data->config; // Configuration 20 psAssert(config, "Require configuration."); 21 22 psMetadata *recipe = psMetadataLookupMetadata(NULL, config->recipes, GDIFF_RECIPE); // Recipe for ppSim 23 psAssert(recipe, "We checked this earlier, so it should be here."); 24 bool noConvolve = psMetadataLookupBool(&mdok, recipe, "NOCONVOLVE"); // Do not use convolved images. 25 26 pmConfigCamerasCull(config, NULL); 27 pmConfigRecipesCull(config, "GDIFF,PPSTATS,PSPHOT,PSASTRO,MASKS,JPEG"); 28 29 pmFPAfile *input = psMetadataLookupPtr(NULL, config->files, "GDIFF.INPUT"); 30 pmFPAfile *reference = psMetadataLookupPtr(NULL, config->files, "GDIFF.REF"); 31 pmFPAfile *output = psMetadataLookupPtr(NULL, config->files, "GDIFF.OUTPUT"); 32 psAssert(input && reference && output, "Require files"); 33 34 if (!gDiffFilesIterateDown(config, GDIFF_FILES_INPUT | GDIFF_FILES_CONV)) { 35 psError(GDIFF_ERR_IO, false, "Unable to load files."); 36 return false; 37 } 38 39 psTimerStart("GDIFF_MATCH"); 40 41 if (!gDiffSetMasks(config)) { 42 psError(psErrorCodeLast(), false, "Unable to set masks."); 43 return false; 44 } 45 46 if (data->forcedPhot1) { 47 bool foundDetections = false; 48 if (!gDiffInputDetections(&foundDetections, "GDIFF.POS1.SOURCES", "GDIFF.INPUT", data)) { 49 psError(psErrorCodeLast(), false, "Unable to measure positive detections (1)"); 50 success = false; 51 goto ESCAPE; 52 } 53 // if nothing was found, don't bother doing the forced photometry below 54 if (!foundDetections) { 55 psWarning ("no sources found in positive image 1, skipping forced photometry"); 56 data->forcedPhot1 = false; 57 } 58 } 59 if (data->forcedPhot2) { 60 // Change the recipe to use a higher nsigma limit and quit after pass1 61 psMetadata *psphotRecipe = psMetadataLookupPtr (NULL, config->recipes, PSPHOT_RECIPE); 62 63 psF32 nsigma_peak_save = psMetadataLookupF32 (NULL, psphotRecipe, "PEAKS_NSIGMA_LIMIT"); 64 char *breakPt_save = psMetadataLookupStr (NULL, psphotRecipe, "BREAK_POINT"); 65 66 psF32 pos2_nsigma_peak = psMetadataLookupF32 (&mdok, psphotRecipe, "PEAKS_POS2_NSIGMA_LIMIT"); 67 if (!mdok) { 68 psWarning("PEAKS_POS2_NSIGMA_LIMIT not found in psphotRecipe. Will use 25.\n"); 69 pos2_nsigma_peak = 25.; 70 } 71 psMetadataAddF32(psphotRecipe, PS_LIST_TAIL, "PEAKS_NSIGMA_LIMIT", PS_META_REPLACE, "", pos2_nsigma_peak); 72 psMetadataAddStr(psphotRecipe, PS_LIST_TAIL, "BREAK_POINT", PS_META_REPLACE, "", "PASS1"); 73 74 bool foundDetections = false; 75 if (!gDiffInputDetections(&foundDetections, "GDIFF.POS2.SOURCES", "GDIFF.REF", data)) { 76 psError(psErrorCodeLast(), false, "Unable to measure positive detections (2)"); 77 psMetadataAddF32(psphotRecipe, PS_LIST_TAIL, "PEAKS_NSIGMA_LIMIT", PS_META_REPLACE, "", nsigma_peak_save); 78 psMetadataAddStr(psphotRecipe, PS_LIST_TAIL, "BREAK_POINT", PS_META_REPLACE, "", breakPt_save); 79 success = false; 80 goto ESCAPE; 81 } 82 psMetadataAddF32(psphotRecipe, PS_LIST_TAIL, "PEAKS_NSIGMA_LIMIT", PS_META_REPLACE, "", nsigma_peak_save); 83 psMetadataAddStr(psphotRecipe, PS_LIST_TAIL, "BREAK_POINT", PS_META_REPLACE, "", breakPt_save); 84 // if nothing was found, don't bother doing the forced photometry below 85 if (!foundDetections) { 86 psWarning ("no sources found in positive image 2, skipping forced photometry"); 87 data->forcedPhot2 = false; 88 } 89 } 90 91 // XXX if it exists, use the POS1, POS2 successs for the FWHMs 92 if (!gDiffMatchPSFs(data)) { 93 psError(psErrorCodeLast(), false, "Unable to match PSFs."); 94 success = false; 95 goto ESCAPE; 96 } 97 98 if (data->quality) { 99 // Can't do anything at all 100 success = false; 101 goto ESCAPE; 102 } 103 // generate the residual stamp grid for visualization 104 if (!gDiffResidualSampleJpeg(config)) { 105 psError(psErrorCodeLast(), false, "Unable to update."); 106 success = false; 107 goto ESCAPE; 108 } 109 110 // XXX add in a positive image detection step here (if needed) 111 112 113 psMetadataAddF32(data->stats, PS_LIST_TAIL, "TIME_MATCH", 0, "Time to match PSFs", 114 psTimerClear("GDIFF_MATCH")); 115 116 // Close input files (freeing up space) : for the noConvolve case, we cannot close the 117 // inputs since we will use them for subtraction below. wait until later to do the work 118 if (!noConvolve) { 119 if (!gDiffFilesIterateUp(config, GDIFF_FILES_INPUT)) { 120 psError(GDIFF_ERR_IO, false, "Unable to close input files."); 121 success = false; 122 goto ESCAPE; 123 } 124 } 125 if (!gDiffLowThreshold(data)) { 126 psError(psErrorCodeLast(), false, "Unable to threshold images."); 127 success = false; 128 goto ESCAPE; 129 } 130 131 // Set up subtraction files 132 if (!gDiffFilesIterateDown(config, GDIFF_FILES_SUB)) { 133 psError(GDIFF_ERR_IO, false, "Unable to set up subtraction files."); 134 success = false; 135 goto ESCAPE; 136 } 137 138 if (!gDiffDefineOutput("GDIFF.OUTPUT", config)) { 139 psError(psErrorCodeLast(), false, "Unable to define output."); 140 success = false; 141 goto ESCAPE; 142 } 143 144 if (!data->quality && !gDiffMakePSF(data)) { 145 psError(psErrorCodeLast(), false, "Unable to generate PSF."); 146 success = false; 147 goto ESCAPE; 148 } 149 150 // Now we've got a PSF, blow away detections in case they're confused with real output detections 151 { 152 pmFPAview *view = gDiffViewReadout(); // View to readout 153 pmReadout *out = pmFPAfileThisReadout(config->files, view, "GDIFF.OUTPUT"); 154 if (psMetadataLookup(out->analysis, "PSPHOT.DETECTIONS")) { 155 psMetadataRemoveKey(out->analysis, "PSPHOT.DETECTIONS"); 156 } 157 psFree(view); 158 } 159 160 if (!gDiffReadoutSubtract(config)) { 161 psError(psErrorCodeLast(), false, "Unable to subtract images."); 162 success = false; 163 goto ESCAPE; 164 } 165 // dumpout(config, "diff.1.fits"); 166 167 // Close convolved files 168 if (!gDiffFilesIterateUp(config, GDIFF_FILES_PSF | GDIFF_FILES_CONV)) { 169 psError(GDIFF_ERR_IO, false, "Unable to close input files."); 170 success = false; 171 goto ESCAPE; 172 } 173 // dumpout(config, "diff.2a.fits"); 174 if (noConvolve) { 175 if (!gDiffFilesIterateUp(config, GDIFF_FILES_INPUT)) { 176 psError(GDIFF_ERR_IO, false, "Unable to close input files."); 177 success = false; 178 goto ESCAPE; 179 } 180 } 181 182 // Higher order background subtraction using psphot 183 if (!gDiffBackground(config)) { 184 psError(psErrorCodeLast(), false, "Unable to subtract background."); 185 success = false; 186 goto ESCAPE; 187 } 188 // dumpout(config, "diff.2b.fits"); 189 190 // Perform Variance correction (rescale within a modest range) 191 if (!gDiffVarianceRescale(config, data)) { 192 psError(psErrorCodeLast(), false, "Unable to rescale variance."); 193 success = false; 194 goto ESCAPE; 195 } 196 // dumpout(config, "diff.2c.fits"); 197 198 if (data->quality) { 199 // Done all we can do up to this point 200 success = false; 201 goto ESCAPE; 202 } 203 204 if (!gDiffFilesIterateDown(config, GDIFF_FILES_PHOT_SUB)) { 205 psError(GDIFF_ERR_IO, false, "Unable to set up photometry files."); 206 success = false; 207 goto ESCAPE; 208 } 209 210 if (!data->quality && !gDiffReadoutPhotometry("GDIFF.OUTPUT", data)) { 211 psError(psErrorCodeLast(), false, "Unable to perform photometry."); 212 success = false; 213 goto ESCAPE; 214 } 215 // dumpout(config, "diff.3.fits"); 216 217 // forced photometry for positive image 1 218 if (data->forcedPhot1 && !data->quality) { 219 if (!gDiffReadoutForcedPhot("GDIFF.FORCED1.SOURCES", "GDIFF.OUTPUT", "GDIFF.POS1.SOURCES", data)) { 220 psError(psErrorCodeLast(), false, "Unable to perform photometry."); 221 success = false; 222 goto ESCAPE; 223 } 224 } 225 226 // forced photometry for positive image 2 227 if (data->forcedPhot2 && !data->quality) { 228 if (!gDiffReadoutForcedPhot("GDIFF.FORCED2.SOURCES", "GDIFF.OUTPUT", "GDIFF.POS2.SOURCES", data)) { 229 psError(psErrorCodeLast(), false, "Unable to perform photometry."); 230 success = false; 231 goto ESCAPE; 232 } 233 } 234 235 if (!gDiffFilesIterateUp(config, GDIFF_FILES_PHOT_SUB)) { 236 psError(GDIFF_ERR_IO, false, "Unable to set up photometry files."); 237 success = false; 238 goto ESCAPE; 239 } 240 241 // Perform statistics on the cell 242 if (!gDiffReadoutStats(data)) { 243 psError(psErrorCodeLast(), false, "Unable to collect statistics"); 244 success = false; 245 goto ESCAPE; 246 } 247 // Do Mask Stats 248 { 249 pmFPAview *view = gDiffViewReadout(); // View to readout 250 if (!gDiffMaskStats(config, view,data->stats)) { 251 psError(psErrorCodeLast(), false, "Unable to generate mask statistics"); 252 success = false; 253 psFree(view); 254 goto ESCAPE; 255 } 256 psFree(view); 257 } 258 // dumpout(config, "diff.4.fits"); 259 260 // generate the binned image used to write the jpeg 261 if (!gDiffReadoutJpeg(config)) { 262 psError(psErrorCodeLast(), false, "Unable to update."); 263 success = false; 264 goto ESCAPE; 265 } 266 267 if (data->inverse) { 268 // Set up inverse subtraction files 269 if (!gDiffFilesIterateDown(config, GDIFF_FILES_INV)) { 270 psError(GDIFF_ERR_IO, false, "Unable to set up inverse files."); 271 success = false; 272 goto ESCAPE; 273 } 274 275 if (data->inverse && !gDiffDefineOutput("GDIFF.INVERSE", config)) { 276 psError(psErrorCodeLast(), false, "Unable to define inverse."); 277 success = false; 278 goto ESCAPE; 279 } 280 281 if (!gDiffReadoutInverse(config)) { 282 psError(psErrorCodeLast(), false, "Unable to invert images."); 283 success = false; 284 goto ESCAPE; 285 } 286 287 //MEH -- need to also clear out detections or inv.cmf corrupted on update 288 pmFPAview *viewinv = gDiffViewReadout(); // View to readout 289 pmReadout *inv = pmFPAfileThisReadout(config->files, viewinv, "GDIFF.INVERSE"); 290 if (psMetadataLookup(inv->analysis, "PSPHOT.DETECTIONS")) { 291 psMetadataRemoveKey(inv->analysis, "PSPHOT.DETECTIONS"); 292 } 293 psFree(viewinv); 294 295 // Close subtraction files and open inverse photometry files 296 if (!gDiffFilesIterateUp(config, GDIFF_FILES_SUB)) { 297 psError(GDIFF_ERR_IO, false, "Unable to close subtraction files."); 298 success = false; 299 goto ESCAPE; 300 } 301 302 if (!gDiffFilesIterateDown(config, GDIFF_FILES_PHOT_INV)) { 303 psError(GDIFF_ERR_IO, false, "Unable to set up inverse files."); 304 success = false; 305 goto ESCAPE; 306 } 307 308 if (!data->quality && !gDiffReadoutPhotometry("GDIFF.INVERSE", data)) { 309 psError(psErrorCodeLast(), false, "Unable to perform photometry."); 310 success = false; 311 goto ESCAPE; 312 } 313 314 // Close inverse subtraction files 315 if (!gDiffFilesIterateUp(config, GDIFF_FILES_INV | GDIFF_FILES_PHOT_INV)) { 316 psError(GDIFF_ERR_IO, false, "Unable to close subtraction files."); 317 success = false; 318 goto ESCAPE; 319 } 320 } else { 321 // Close subtraction files 322 // dumpout(config, "diff.5.fits"); 323 if (!gDiffFilesIterateUp(config, GDIFF_FILES_SUB)) { 324 psError(GDIFF_ERR_IO, false, "Unable to close subtraction files."); 325 success = false; 326 goto ESCAPE; 327 } 328 } 329 330 ESCAPE: 331 pmFPAfileDropInternal(config->files, "PSPHOT.BACKGND"); 332 pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL"); 333 pmFPAfileDropInternal(config->files, "PSPHOT.BACKMDL.STDEV"); 334 335 return success; 387 pmFPAfile *file = psMetadataLookupPtr(NULL, config->files, name); 388 fprintf (stderr, "HDU: %s, %llx\n", name, (long long unsigned int) file->fpa->hdu); 389 return true; 336 390 } -
branches/eam_branches/gDiff.20200202/src/gDiffVersion.c
r41260 r41413 1 /** @file gDiffVersion.c 2 * @brief: image subtraction based on matched Gaussians profiles 3 * @author Eugene Magnier @ IfA 4 * @version $Revision: 1.12 $ 5 * @date $Date: 2020-01-25 $ 6 * Copyright 2020 Institute for Astronomy, University of Hawaii 7 */ 8 1 9 #include "gDiff.h" 2 #include <ppStats.h>3 10 #include "gDiffVersionDefinitions.h" 4 11
Note:
See TracChangeset
for help on using the changeset viewer.
