Index: trunk/psphot/src/psphotImageQuality.c =================================================================== --- trunk/psphot/src/psphotImageQuality.c (revision 25755) +++ trunk/psphot/src/psphotImageQuality.c (revision 26894) @@ -4,7 +4,54 @@ // XXX Lots of code duplication here +// for now, let's store the detections on the readout->analysis for each readout +bool psphotImageQuality (pmConfig *config, const pmFPAview *view) +{ + bool status = true; + + // select the appropriate recipe information + psMetadata *recipe = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE); + psAssert (recipe, "missing recipe?"); + + int num = psMetadataLookupS32 (&status, config->arguments, "PSPHOT.INPUT.NUM"); + psAssert (status, "programming error: must define PSPHOT.INPUT.NUM"); + + // loop over the available readouts + for (int i = 0; i < num; i++) { + if (!psphotImageQualityReadout (config, view, "PSPHOT.INPUT", i, recipe)) { + psError (PSPHOT_ERR_CONFIG, false, "failed on to measure image quality for PSPHOT.INPUT entry %d", i); + return false; + } + } + return true; +} + // selecting the 'good' stars (likely to be psf stars), measure the M_cn, M_sn terms for n = 2,3,4 -bool psphotImageQuality(psMetadata *recipe, psArray *sources) +bool psphotImageQualityReadout(pmConfig *config, const pmFPAview *view, const char *filename, int index, psMetadata *recipe) { + bool status = true; + + // find the currently selected readout + pmFPAfile *file = pmFPAfileSelectSingle(config->files, filename, index); // File of interest + psAssert (file, "missing file?"); + + pmReadout *readout = pmFPAviewThisReadout(view, file->fpa); + psAssert (readout, "missing readout?"); + + // if we have a PSF, skip this analysis + if (psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.PSF")) { + return true; + } + + pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS"); + psAssert (detections, "missing detections?"); + + psArray *sources = detections->newSources; + psAssert (sources, "missing sources?"); + + if (!sources->n) { + psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping image quality"); + return true; + } + psVector *FWHM_MAJOR = psVectorAllocEmpty(100, PS_TYPE_F32); psVector *FWHM_MINOR = psVectorAllocEmpty(100, PS_TYPE_F32); @@ -81,10 +128,9 @@ if (num == 0) { - psWarning("Unable to find sources from which to measure image quality"); - return false; - } - - psMetadataAddS32(recipe, PS_LIST_TAIL, "IQ_NSTAR", PS_META_REPLACE, - "Number of stars used for IQ measurements", M2->n); + psLogMsg ("psphot", PS_LOG_INFO, "no valid sources for image quality, skipping"); + return true; + } + + psMetadataAddS32(readout->analysis, PS_LIST_TAIL, "IQ_NSTAR", PS_META_REPLACE, "Number of stars used for IQ measurements", M2->n); // XXX make this a recipe option @@ -98,8 +144,6 @@ float fwhm_major = stats->sampleMean; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW1", PS_META_REPLACE, - "FWHM of Major Axis from moments", stats->sampleMean); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW1_E", PS_META_REPLACE, - "FWHM scatter (Major) from moments", stats->sampleStdev); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW1", PS_META_REPLACE, "FWHM of Major Axis from moments", stats->sampleMean); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW1_E", PS_META_REPLACE, "FWHM scatter (Major) from moments", stats->sampleStdev); if (!psVectorStats(stats, FWHM_MINOR, NULL, NULL, 0)) { @@ -108,8 +152,6 @@ } float fwhm_minor = stats->sampleMean; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW2", PS_META_REPLACE, - "FWHM of Minor Axis from moments", stats->sampleMean); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW2_E", PS_META_REPLACE, - "FWHM scatter (Minor) from moments", stats->sampleStdev); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW2", PS_META_REPLACE, "FWHM of Minor Axis from moments", stats->sampleMean); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW2_E", PS_META_REPLACE, "FWHM scatter (Minor) from moments", stats->sampleStdev); if (!psVectorStats(stats, M2, NULL, NULL, 0)) { @@ -119,12 +161,8 @@ float vM2 = stats->sampleMean; float dM2 = stats->sampleStdev; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2", PS_META_REPLACE, - "M_2 = sqrt (M_c2^2 + M_s2^2)", vM2); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_ER", PS_META_REPLACE, - "Stdev of M_2", dM2); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_LQ", PS_META_REPLACE, - "Lower Quartile of M_2", stats->sampleLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_UQ", PS_META_REPLACE, - "Upper Quartile of M_2", stats->sampleUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2", PS_META_REPLACE, "M_2 = sqrt (M_c2^2 + M_s2^2)", vM2); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_ER", PS_META_REPLACE, "Stdev of M_2", dM2); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_LQ", PS_META_REPLACE, "Lower Quartile of M_2", stats->sampleLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_UQ", PS_META_REPLACE, "Upper Quartile of M_2", stats->sampleUQ); if (!psVectorStats(stats, M2c, NULL, NULL, 0)) { @@ -132,12 +170,8 @@ goto FAIL; } - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C", PS_META_REPLACE, - "M_2c = sum f r^2 cos(2phi) / sum f", stats->sampleMean); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_E", PS_META_REPLACE, - "Stdev of M_2c", stats->sampleStdev); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_L", PS_META_REPLACE, - "Lower Quartile of M_2c", stats->sampleLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_U", PS_META_REPLACE, - "Upper Quartile of M_2c", stats->sampleUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C", PS_META_REPLACE, "M_2c = sum f r^2 cos(2phi) / sum f", stats->sampleMean); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_E", PS_META_REPLACE, "Stdev of M_2c", stats->sampleStdev); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_L", PS_META_REPLACE, "Lower Quartile of M_2c", stats->sampleLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_U", PS_META_REPLACE, "Upper Quartile of M_2c", stats->sampleUQ); if (!psVectorStats(stats, M2s, NULL, NULL, 0)) { @@ -145,12 +179,8 @@ goto FAIL; } - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S", PS_META_REPLACE, - "M_2s = sum f r^2 cos(2phi) / sum f", stats->sampleMean); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_E", PS_META_REPLACE, - "Stdev of M_2s", stats->sampleStdev); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_L", PS_META_REPLACE, - "Lower Quartile of M_2s", stats->sampleLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_U", PS_META_REPLACE, - "Upper Quartile of M_2s", stats->sampleUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S", PS_META_REPLACE, "M_2s = sum f r^2 cos(2phi) / sum f", stats->sampleMean); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_E", PS_META_REPLACE, "Stdev of M_2s", stats->sampleStdev); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_L", PS_META_REPLACE, "Lower Quartile of M_2s", stats->sampleLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_U", PS_META_REPLACE, "Upper Quartile of M_2s", stats->sampleUQ); if (!psVectorStats(stats, M3, NULL, NULL, 0)) { @@ -160,12 +190,8 @@ float vM3 = stats->sampleMean; float dM3 = stats->sampleStdev; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3", PS_META_REPLACE, - "M_3 = sqrt (M_c3^2 + M_s3^2)", vM3); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_ER", PS_META_REPLACE, - "Stdev of M_3", dM3); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_LQ", PS_META_REPLACE, - "Lower Quartile of M_3", stats->sampleLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_UQ", PS_META_REPLACE, - "Upper Quartile of M_3", stats->sampleUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3", PS_META_REPLACE, "M_3 = sqrt (M_c3^2 + M_s3^2)", vM3); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_ER", PS_META_REPLACE, "Stdev of M_3", dM3); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_LQ", PS_META_REPLACE, "Lower Quartile of M_3", stats->sampleLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_UQ", PS_META_REPLACE, "Upper Quartile of M_3", stats->sampleUQ); if (!psVectorStats(stats, M4, NULL, NULL, 0)) { @@ -175,12 +201,8 @@ float vM4 = stats->sampleMean; float dM4 = stats->sampleStdev; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4", PS_META_REPLACE, - "M_4 = sqrt (M_c4^2 + M_s4^2)", vM4); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_ER", PS_META_REPLACE, - "Stdev of M_4", dM4); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_LQ", PS_META_REPLACE, - "Lower Quartile of M_4", stats->sampleLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_UQ", PS_META_REPLACE, - "Upper Quartile of M_4", stats->sampleUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4", PS_META_REPLACE, "M_4 = sqrt (M_c4^2 + M_s4^2)", vM4); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_ER", PS_META_REPLACE, "Stdev of M_4", dM4); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_LQ", PS_META_REPLACE, "Lower Quartile of M_4", stats->sampleLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_UQ", PS_META_REPLACE, "Upper Quartile of M_4", stats->sampleUQ); #else @@ -192,8 +214,6 @@ } float fwhm_major = stats->robustMedian; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW1", PS_META_REPLACE, - "FWHM of Major Axis from moments", stats->robustMedian); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW1_E", PS_META_REPLACE, - "FWHM scatter (Major) from moments", stats->robustStdev); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW1", PS_META_REPLACE, "FWHM of Major Axis from moments", stats->robustMedian); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW1_E", PS_META_REPLACE, "FWHM scatter (Major) from moments", stats->robustStdev); if (!psVectorStats(stats, FWHM_MINOR, NULL, NULL, 0)) { @@ -202,8 +222,6 @@ } float fwhm_minor = stats->robustMedian; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW2", PS_META_REPLACE, - "FWHM of Minor Axis from moments", stats->robustMedian); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_FW2_E", PS_META_REPLACE, - "FWHM scatter (Minor) from moments", stats->robustStdev); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW2", PS_META_REPLACE, "FWHM of Minor Axis from moments", stats->robustMedian); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_FW2_E", PS_META_REPLACE, "FWHM scatter (Minor) from moments", stats->robustStdev); if (!psVectorStats(stats, M2, NULL, NULL, 0)) { @@ -213,12 +231,8 @@ float vM2 = stats->robustMedian; float dM2 = stats->robustStdev; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2", PS_META_REPLACE, - "M_2 = sqrt (M_c2^2 + M_s2^2)", vM2); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_ER", PS_META_REPLACE, - "Stdev of M_2", dM2); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_LQ", PS_META_REPLACE, - "Lower Quartile of M_2", stats->robustLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2_UQ", PS_META_REPLACE, - "Upper Quartile of M_2", stats->robustUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2", PS_META_REPLACE, "M_2 = sqrt (M_c2^2 + M_s2^2)", vM2); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_ER", PS_META_REPLACE, "Stdev of M_2", dM2); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_LQ", PS_META_REPLACE, "Lower Quartile of M_2", stats->robustLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2_UQ", PS_META_REPLACE, "Upper Quartile of M_2", stats->robustUQ); if (!psVectorStats(stats, M2c, NULL, NULL, 0)) { @@ -226,12 +240,8 @@ goto FAIL; } - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C", PS_META_REPLACE, - "M_2c = sum f r^2 cos(2phi) / sum f", stats->robustMedian); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_E", PS_META_REPLACE, - "Stdev of M_2c", stats->robustStdev); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_L", PS_META_REPLACE, - "Lower Quartile of M_2c", stats->robustLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2C_U", PS_META_REPLACE, - "Upper Quartile of M_2c", stats->robustUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C", PS_META_REPLACE, "M_2c = sum f r^2 cos(2phi) / sum f", stats->robustMedian); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_E", PS_META_REPLACE, "Stdev of M_2c", stats->robustStdev); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_L", PS_META_REPLACE, "Lower Quartile of M_2c", stats->robustLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2C_U", PS_META_REPLACE, "Upper Quartile of M_2c", stats->robustUQ); if (!psVectorStats(stats, M2s, NULL, NULL, 0)) { @@ -239,12 +249,8 @@ goto FAIL; } - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S", PS_META_REPLACE, - "M_2s = sum f r^2 cos(2phi) / sum f", stats->robustMedian); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_E", PS_META_REPLACE, - "Stdev of M_2s", stats->robustStdev); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_L", PS_META_REPLACE, - "Lower Quartile of M_2s", stats->robustLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M2S_U", PS_META_REPLACE, - "Upper Quartile of M_2s", stats->robustUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S", PS_META_REPLACE, "M_2s = sum f r^2 cos(2phi) / sum f", stats->robustMedian); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_E", PS_META_REPLACE, "Stdev of M_2s", stats->robustStdev); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_L", PS_META_REPLACE, "Lower Quartile of M_2s", stats->robustLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M2S_U", PS_META_REPLACE, "Upper Quartile of M_2s", stats->robustUQ); if (!psVectorStats(stats, M3, NULL, NULL, 0)) { @@ -254,12 +260,8 @@ float vM3 = stats->robustMedian; float dM3 = stats->robustStdev; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3", PS_META_REPLACE, - "M_3 = sqrt (M_c3^2 + M_s3^2)", vM3); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_ER", PS_META_REPLACE, - "Stdev of M_3", dM3); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_LQ", PS_META_REPLACE, - "Lower Quartile of M_3", stats->robustLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M3_UQ", PS_META_REPLACE, - "Upper Quartile of M_3", stats->robustUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3", PS_META_REPLACE, "M_3 = sqrt (M_c3^2 + M_s3^2)", vM3); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_ER", PS_META_REPLACE, "Stdev of M_3", dM3); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_LQ", PS_META_REPLACE, "Lower Quartile of M_3", stats->robustLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M3_UQ", PS_META_REPLACE, "Upper Quartile of M_3", stats->robustUQ); if (!psVectorStats(stats, M4, NULL, NULL, 0)) { @@ -269,12 +271,8 @@ float vM4 = stats->robustMedian; float dM4 = stats->robustStdev; - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4", PS_META_REPLACE, - "M_4 = sqrt (M_c4^2 + M_s4^2)", vM4); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_ER", PS_META_REPLACE, - "Stdev of M_4", dM4); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_LQ", PS_META_REPLACE, - "Lower Quartile of M_4", stats->robustLQ); - psMetadataAddF32(recipe, PS_LIST_TAIL, "IQ_M4_UQ", PS_META_REPLACE, - "Upper Quartile of M_4", stats->robustUQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4", PS_META_REPLACE, "M_4 = sqrt (M_c4^2 + M_s4^2)", vM4); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_ER", PS_META_REPLACE, "Stdev of M_4", dM4); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_LQ", PS_META_REPLACE, "Lower Quartile of M_4", stats->robustLQ); + psMetadataAddF32(readout->analysis, PS_LIST_TAIL, "IQ_M4_UQ", PS_META_REPLACE, "Upper Quartile of M_4", stats->robustUQ); #endif