Changeset 34259
- Timestamp:
- Jul 31, 2012, 3:59:03 PM (14 years ago)
- Location:
- trunk/psModules
- Files:
-
- 3 deleted
- 9 edited
- 2 copied
-
. (modified) (1 prop)
-
src/objects (modified) (1 prop)
-
src/objects/Makefile.am (modified) (3 diffs)
-
src/objects/mksource.pl (modified) (4 diffs)
-
src/objects/models/pmModel_TRAIL.c (copied) (copied from branches/eam_branches/ipp-20120627/psModules/src/objects/models/pmModel_TRAIL.c )
-
src/objects/models/pmModel_TRAIL.h (copied) (copied from branches/eam_branches/ipp-20120627/psModules/src/objects/models/pmModel_TRAIL.h )
-
src/objects/pmModelClass.c (modified) (2 diffs)
-
src/objects/pmModelFuncs.h (modified) (1 diff)
-
src/objects/pmSourceFitModel.c (modified) (1 diff)
-
src/objects/pmSourceFitModel.h (modified) (1 diff)
-
src/objects/pmSourceIO_CMF.c.in (modified) (15 diffs)
-
src/objects/pmSourceIO_CMF_PS1_DV1.c (deleted)
-
src/objects/pmSourceIO_CMF_PS1_DV2.c (deleted)
-
src/objects/pmSourceIO_CMF_PS1_SV1.c (deleted)
Legend:
- Unmodified
- Added
- Removed
-
trunk/psModules
- Property svn:mergeinfo changed
/branches/eam_branches/ipp-20120627/psModules (added) merged: 34133,34141,34144,34181-34183,34185-34188,34192,34246,34252,34255
- Property svn:mergeinfo changed
-
trunk/psModules/src/objects
- Property svn:ignore
-
old new 9 9 pmSourceIO_CMF_PS1_V3.c 10 10 pmSourceIO_CMF_PS1_V4.c 11 pmSourceIO_CMF_PS1_V3.v1.c 12 pmSourceIO_CMF_PS1_V1.v1.c 13 pmSourceIO_CMF_PS1_V2.v1.c 14 11 pmSourceIO_CMF_PS1_SV1.c 12 pmSourceIO_CMF_PS1_DV1.c 13 pmSourceIO_CMF_PS1_DV2.c
-
- Property svn:ignore
-
trunk/psModules/src/objects/Makefile.am
r32633 r34259 80 80 models/pmModel_SERSIC.c \ 81 81 models/pmModel_EXP.c \ 82 models/pmModel_DEV.c 82 models/pmModel_DEV.c \ 83 models/pmModel_TRAIL.c 83 84 84 85 pkginclude_HEADERS = \ … … 126 127 models/pmModel_SERSIC.h \ 127 128 models/pmModel_EXP.h \ 128 models/pmModel_DEV.h 129 models/pmModel_DEV.h \ 130 models/pmModel_TRAIL.h 129 131 130 132 CLEANFILES = *~ 131 133 132 134 # pmSourceID_CMF_* functions use a common framework 133 BUILT_SOURCES = pmSourceIO_CMF_PS1_V1.c pmSourceIO_CMF_PS1_V2.c pmSourceIO_CMF_PS1_V3.c pmSourceIO_CMF_PS1_ V4.c135 BUILT_SOURCES = pmSourceIO_CMF_PS1_V1.c pmSourceIO_CMF_PS1_V2.c pmSourceIO_CMF_PS1_V3.c pmSourceIO_CMF_PS1_DV1.c pmSourceIO_CMF_PS1_DV2.c pmSourceIO_CMF_PS1_SV1.c 134 136 135 137 pmSourceIO_CMF_PS1_V1.c : pmSourceIO_CMF.c.in mksource.pl … … 145 147 mksource.pl pmSourceIO_CMF.c.in PS1_V4 pmSourceIO_CMF_PS1_V4.c 146 148 149 pmSourceIO_CMF_PS1_DV1.c : pmSourceIO_CMF.c.in mksource.pl 150 mksource.pl pmSourceIO_CMF.c.in PS1_DV1 pmSourceIO_CMF_PS1_DV1.c 151 152 pmSourceIO_CMF_PS1_DV2.c : pmSourceIO_CMF.c.in mksource.pl 153 mksource.pl pmSourceIO_CMF.c.in PS1_DV2 pmSourceIO_CMF_PS1_DV2.c 154 155 pmSourceIO_CMF_PS1_SV1.c : pmSourceIO_CMF.c.in mksource.pl 156 mksource.pl pmSourceIO_CMF.c.in PS1_SV1 pmSourceIO_CMF_PS1_SV1.c 157 147 158 # EXTRA_DIST = pmErrorCodes.h.in pmErrorCodes.dat pmErrorCodes.c.in -
trunk/psModules/src/objects/mksource.pl
r32633 r34259 1 1 #!/usr/bin/env perl 2 $DEBUG = 0; 2 3 3 4 # this program takes the pmSourceIO_CMF.in.c template file and generates the .c version based on the given I/O format made … … 17 18 "PS1_V2", 2, 18 19 "PS1_V3", 3, 19 "PS1_V4", 4); 20 21 print "1: $cmfmodes{1}\n"; 22 print "PS1_V1: $cmfmodes{'PS1_V1'}\n"; 20 "PS1_V4", 4, 21 ); 23 22 24 23 open (FILE, "$template") || die "failed to open template $template\n"; … … 35 34 # @<MODE@ : remove and keep if cmfmode > MODE 36 35 36 # XXX need to add features: split @foo,bar,baz@ by commas 37 # treat each chunk as a rule 38 # add the following options: 39 # !MODE -- exclude the given entry (defaults to all, or is ALL required?) 40 # * and ? regexp 41 42 # some examples: 43 # @ALL,!PS1_V1@ 44 # @PS1_DV?@ 45 # @PS1_V?,!PS1_V1@ 46 37 47 foreach $line (@list) { 38 48 … … 40 50 $line =~ s|\@CMFMODE\@|$cmfmode|g; 41 51 42 if ($line =~ m|\@ALL\@|) { 43 $line =~ s|\@ALL\@\s*||; 52 # print and continue if we do not match @RULES@ 53 unless ($line =~ m|\@.*\@|) { 54 print "no rule\n" if $DEBUG; 55 print FILE $line; 56 next; 44 57 } 45 58 46 ($isMode) = ($line =~ m|\@=(\S*)\@|); 47 ($gtMode) = ($line =~ m|\@>(\S*)\@|); 48 ($ltMode) = ($line =~ m|\@<(\S*)\@|); 59 # grab the rules and the rest of the line 60 ($prefix,$rules,$content) = ($line =~ m|(.*)\@(.*)\@\s*(.*)|); 49 61 50 if ($isMode) { 51 if ($isMode ne $cmfmode) { next; } 52 $line =~ s|\@=\S*\@\s*||; 62 # split the rules into separate items 63 @rules = split (",", $rules); 64 65 $keepLine = 0; 66 # does $cmfmode match any of the rules? 67 foreach $rule (@rules) { 68 print "rule: $rule\n" if $DEBUG; 69 70 # special rule "ALL" 71 if ($rule eq "ALL") { 72 print "ALL match\n" if $DEBUG; 73 $keepLine = 1; 74 next; 75 } # look for other rules (esp !foo) 76 77 # pure match 78 if ($cmfmode eq $rule) { 79 print "simple match\n" if $DEBUG; 80 $keepLine = 1; 81 next; 82 } # skip to end? 83 84 # NOT match 85 if ($rule =~ m|^!|) { 86 print "NOT rule: $rule\n" if $DEBUG; 87 ($realrule) = ($rule =~ m|^!(.*)|); 88 if ($cmfmode eq $realrule) { $keepLine = 0; } # skip to end? 89 next; 90 } 91 92 # simple regexp: foo* 93 if ($rule =~ m|\*$|) { 94 print "regexp * rule: $rule\n" if $DEBUG; 95 ($realrule) = ($rule =~ m|(.*)\*$|); 96 if ($cmfmode =~ m|$realrule|) { $keepLine = 1; } # skip to end? 97 next; 98 } 99 100 # simple regexp: foo? 101 if ($rule =~ m|\?$|) { 102 print "regexp ? rule: $rule\n" if $DEBUG; 103 ($realrule) = ($rule =~ m|(.*)\?$|); 104 if (substr($cmfmode,0,-1) eq $realrule) { $keepLine = 1; } # skip to end? 105 next; 106 } 107 108 # rule: =FOO 109 if ($rule =~ m|^=|) { 110 print "= rule: $rule\n" if $DEBUG; 111 $realrule = substr($cmfmode,1); 112 if ($cmfmode eq $realrule) { $keepLine = 1; } # skip to end? 113 next; 114 } 115 116 # only apply the < > <= >= rules if cmfmode is one of cmfmodes 117 # rule: >=FOO 118 if ($rule =~ m|^>=|) { 119 print ">= rule: $rule\n" if $DEBUG; 120 if ($cmfmodes{$cmfmode} == 0) { next; } 121 $realrule = substr($rule,2); 122 $thisLevel = $cmfmodes{$realrule}; 123 $myLevel = $cmfmodes{$cmfmode}; 124 print "levels: $cmfmode, $realrule, $myLevel, $thisLevel\n" if $DEBUG; 125 if ($myLevel >= $thisLevel) { $keepLine = 1; } 126 next; 127 } 128 129 # rule: >FOO 130 if ($rule =~ m|^>|) { 131 print "> rule: $rule\n" if $DEBUG; 132 if ($cmfmodes{$cmfmode} == 0) { next; } 133 $realrule = substr($rule,1); 134 $thisLevel = $cmfmodes{$realrule}; 135 $myLevel = $cmfmodes{$cmfmode}; 136 print "levels: $cmfmode, $realrule, $myLevel, $thisLevel\n" if $DEBUG; 137 if ($myLevel > $thisLevel) { $keepLine = 1; } 138 next; 139 } 140 141 # rule: <=FOO 142 if ($rule =~ m|^<=|) { 143 print "<= rule: $rule\n" if $DEBUG; 144 if ($cmfmodes{$cmfmode} == 0) { next; } 145 $realrule = substr($rule,2); 146 $thisLevel = $cmfmodes{$realrule}; 147 $myLevel = $cmfmodes{$cmfmode}; 148 print "levels: $cmfmode, $realrule, $myLevel, $thisLevel\n" if $DEBUG; 149 if ($myLevel <= $thisLevel) { $keepLine = 1; } 150 next; 151 } 152 153 # rule: <FOO 154 if ($rule =~ m|^<|) { 155 print "< rule: $rule\n" if $DEBUG; 156 if ($cmfmodes{$cmfmode} == 0) { next; } 157 $realrule = substr($rule,1); 158 $thisLevel = $cmfmodes{$realrule}; 159 $myLevel = $cmfmodes{$cmfmode}; 160 print "levels: $cmfmode, $realrule, $myLevel, $thisLevel\n" if $DEBUG; 161 if ($myLevel < $thisLevel) { $keepLine = 1; } 162 next; 163 } 164 53 165 } 166 print "line: $line\n" if $DEBUG; 54 167 55 if ($gtMode) { 56 # print "gtMode : $line\n"; 57 $thisLevel = $cmfmodes{$gtMode}; 58 $myLevel = $cmfmodes{$cmfmode}; 59 print "gtMode : $gtMode vs $cmfmode, $thisLevel, $myLevel\n"; 60 if ($myLevel <= $thisLevel) { next; } 61 $line =~ s|\@>\S*\@\s*||; 168 if ($keepLine) { 169 print FILE "$prefix $content\n"; 62 170 } 63 64 if ($ltMode) {65 # print "ltMode : $line\n";66 $thisLevel = $cmfmodes{$ltMode};67 $myLevel = $cmfmodes{$cmfmode};68 print "ltMode : $ltMode vs $cmfmode, $thisLevel, $myLevel\n";69 if ($myLevel >= $thisLevel) { next; }70 $line =~ s|\@<\S*\@\s*||;71 }72 73 print FILE $line;74 171 } 75 172 -
trunk/psModules/src/objects/pmModelClass.c
r29004 r34259 51 51 # include "models/pmModel_EXP.h" 52 52 # include "models/pmModel_DEV.h" 53 # include "models/pmModel_TRAIL.h" 53 54 54 55 static pmModelClass defaultModels[] = { … … 59 60 {"PS_MODEL_RGAUSS", 8, (pmModelFunc)pmModelFunc_RGAUSS, (pmModelFlux)pmModelFlux_RGAUSS, (pmModelRadius)pmModelRadius_RGAUSS, (pmModelLimits)pmModelLimits_RGAUSS, (pmModelGuessFunc)pmModelGuess_RGAUSS, (pmModelFromPSFFunc)pmModelFromPSF_RGAUSS, (pmModelParamsFromPSF)pmModelParamsFromPSF_RGAUSS, (pmModelFitStatusFunc)pmModelFitStatus_RGAUSS, (pmModelSetLimitsFunc)pmModelSetLimits_RGAUSS }, 60 61 {"PS_MODEL_SERSIC", 8, (pmModelFunc)pmModelFunc_SERSIC, (pmModelFlux)pmModelFlux_SERSIC, (pmModelRadius)pmModelRadius_SERSIC, (pmModelLimits)pmModelLimits_SERSIC, (pmModelGuessFunc)pmModelGuess_SERSIC, (pmModelFromPSFFunc)pmModelFromPSF_SERSIC, (pmModelParamsFromPSF)pmModelParamsFromPSF_SERSIC, (pmModelFitStatusFunc)pmModelFitStatus_SERSIC, (pmModelSetLimitsFunc)pmModelSetLimits_SERSIC }, 61 {"PS_MODEL_EXP", 7, (pmModelFunc)pmModelFunc_EXP, (pmModelFlux)pmModelFlux_EXP, (pmModelRadius)pmModelRadius_EXP, (pmModelLimits)pmModelLimits_EXP, (pmModelGuessFunc)pmModelGuess_EXP, (pmModelFromPSFFunc)pmModelFromPSF_EXP, (pmModelParamsFromPSF)pmModelParamsFromPSF_EXP, (pmModelFitStatusFunc)pmModelFitStatus_EXP, (pmModelSetLimitsFunc)pmModelSetLimits_EXP }, 62 {"PS_MODEL_DEV", 7, (pmModelFunc)pmModelFunc_DEV, (pmModelFlux)pmModelFlux_DEV, (pmModelRadius)pmModelRadius_DEV, (pmModelLimits)pmModelLimits_DEV, (pmModelGuessFunc)pmModelGuess_DEV, (pmModelFromPSFFunc)pmModelFromPSF_DEV, (pmModelParamsFromPSF)pmModelParamsFromPSF_DEV, (pmModelFitStatusFunc)pmModelFitStatus_DEV, (pmModelSetLimitsFunc)pmModelSetLimits_DEV }, 62 {"PS_MODEL_EXP", 7, (pmModelFunc)pmModelFunc_EXP, (pmModelFlux)pmModelFlux_EXP, (pmModelRadius)pmModelRadius_EXP, (pmModelLimits)pmModelLimits_EXP, (pmModelGuessFunc)pmModelGuess_EXP, (pmModelFromPSFFunc)pmModelFromPSF_EXP, (pmModelParamsFromPSF)pmModelParamsFromPSF_EXP, (pmModelFitStatusFunc)pmModelFitStatus_EXP, (pmModelSetLimitsFunc)pmModelSetLimits_EXP }, 63 {"PS_MODEL_DEV", 7, (pmModelFunc)pmModelFunc_DEV, (pmModelFlux)pmModelFlux_DEV, (pmModelRadius)pmModelRadius_DEV, (pmModelLimits)pmModelLimits_DEV, (pmModelGuessFunc)pmModelGuess_DEV, (pmModelFromPSFFunc)pmModelFromPSF_DEV, (pmModelParamsFromPSF)pmModelParamsFromPSF_DEV, (pmModelFitStatusFunc)pmModelFitStatus_DEV, (pmModelSetLimitsFunc)pmModelSetLimits_DEV }, 64 {"PS_MODEL_TRAIL", 7, (pmModelFunc)pmModelFunc_TRAIL, (pmModelFlux)pmModelFlux_TRAIL, (pmModelRadius)pmModelRadius_TRAIL, (pmModelLimits)pmModelLimits_TRAIL, (pmModelGuessFunc)pmModelGuess_TRAIL, (pmModelFromPSFFunc)pmModelFromPSF_TRAIL, (pmModelParamsFromPSF)pmModelParamsFromPSF_TRAIL, (pmModelFitStatusFunc)pmModelFitStatus_TRAIL, (pmModelSetLimitsFunc)pmModelSetLimits_TRAIL }, 63 65 }; 64 66 -
trunk/psModules/src/objects/pmModelFuncs.h
r34085 r34259 81 81 #define PM_PAR_8 8 ///< Model-dependent parameter 82 82 83 // these are used by pmModel_TRAIL, with refers to L and Theta explicitly 84 #define PM_PAR_LENGTH 4 ///< trail length 85 #define PM_PAR_THETA 5 ///< position angle 86 #define PM_PAR_SIGMA 6 ///< position angle 87 83 88 /*** these prototype classes are used to define elements of the pmModelClass structure below ***/ 84 89 -
trunk/psModules/src/objects/pmSourceFitModel.c
r32347 r34259 182 182 constraint->paramMask->data.PS_TYPE_VECTOR_MASK_DATA[PM_PAR_YPOS] = 1; 183 183 constraint->paramMask->data.PS_TYPE_VECTOR_MASK_DATA[PM_PAR_SKY] = 1; 184 break; 185 case PM_SOURCE_FIT_TRAIL: 186 // special mode for pmModel_TRAIL: Io, Xo, Yo, Length, and Theta (not Io or Sigma) 187 nParams = params->n - 3; 188 psVectorInit (constraint->paramMask, 0); 189 constraint->paramMask->data.PS_TYPE_VECTOR_MASK_DATA[PM_PAR_SKY] = 1; 190 constraint->paramMask->data.PS_TYPE_VECTOR_MASK_DATA[PM_PAR_SIGMA] = 1; 184 191 break; 185 192 case PM_SOURCE_FIT_INDEX: -
trunk/psModules/src/objects/pmSourceFitModel.h
r31153 r34259 22 22 PM_SOURCE_FIT_INDEX, 23 23 PM_SOURCE_FIT_NO_INDEX, 24 PM_SOURCE_FIT_TRAIL, 24 25 } pmSourceFitMode; 25 26 -
trunk/psModules/src/objects/pmSourceIO_CMF.c.in
r33693 r34259 89 89 // by the time we call this function, all values should be assigned. let's use asserts to be sure in some cases. 90 90 for (int i = 0; i < sources->n; i++) { 91 pmSource *source = (pmSource *) sources->data[i]; 91 // this is the source associated with this image 92 pmSource *thisSource = sources->data[i]; 93 94 // this is the "real" version of this source 95 pmSource *source = thisSource->parent ? thisSource->parent : thisSource; 92 96 93 97 // If source->seq is -1, source was generated in this analysis. If source->seq is … … 106 110 pmSourceOutputsSetMoments (&moments, source); 107 111 112 @PS1_DV?@ pmSourceDiffStats diffStats; 113 @PS1_DV?@ pmSourceDiffStatsInit(&diffStats); 114 @PS1_DV?@ if (source->diffStats) { 115 @PS1_DV?@ diffStats = *source->diffStats; 116 @PS1_DV?@ } 117 108 118 row = psMetadataAlloc (); 109 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "IPP_IDET", PS_DATA_U32, "IPP detection identifier index", source->seq);110 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "X_PSF", PS_DATA_F32, "PSF x coordinate", outputs.xPos);111 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "Y_PSF", PS_DATA_F32, "PSF y coordinate", outputs.yPos);112 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "X_PSF_SIG", PS_DATA_F32, "Sigma in PSF x coordinate", outputs.xErr);113 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "Y_PSF_SIG", PS_DATA_F32, "Sigma in PSF y coordinate", outputs.yErr);119 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "IPP_IDET", PS_DATA_U32, "IPP detection identifier index", source->seq); 120 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "X_PSF", PS_DATA_F32, "PSF x coordinate", outputs.xPos); 121 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "Y_PSF", PS_DATA_F32, "PSF y coordinate", outputs.yPos); 122 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "X_PSF_SIG", PS_DATA_F32, "Sigma in PSF x coordinate", outputs.xErr); 123 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "Y_PSF_SIG", PS_DATA_F32, "Sigma in PSF y coordinate", outputs.yErr); 114 124 115 125 // NOTE: pre-PS1_V2, we only reported RA & DEC in floats for reference, not precision 116 @=PS1_V1@ psMetadataAdd (row, PS_LIST_TAIL, "RA_PSF", PS_DATA_F32, "PSF RA coordinate (degrees)", outputs.ra); 117 @=PS1_V1@ psMetadataAdd (row, PS_LIST_TAIL, "DEC_PSF", PS_DATA_F32, "PSF DEC coordinate (degrees)", outputs.dec); 118 119 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "POSANGLE", PS_DATA_F32, "position angle at source (degrees)", outputs.posAngle); 120 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PLTSCALE", PS_DATA_F32, "plate scale at source (arcsec/pixel)", outputs.pltScale); 121 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_MAG", PS_DATA_F32, "PSF fit instrumental magnitude", source->psfMag); 122 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_MAG_SIG", PS_DATA_F32, "Sigma of PSF instrumental magnitude", source->psfMagErr); 123 124 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX", PS_DATA_F32, "PSF fit instrumental flux (counts)", source->psfFlux); 125 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX_SIG",PS_DATA_F32, "Sigma of PSF instrumental flux", source->psfFluxErr); 126 127 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG", PS_DATA_F32, "magnitude in standard aperture", source->apMag); 128 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RAW", PS_DATA_F32, "magnitude in reported aperture", source->apMagRaw); 129 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RADIUS", PS_DATA_F32, "radius used for aperture mags", outputs.apRadius); 130 @<PS1_V3@ psMetadataAdd (row, PS_LIST_TAIL, "PEAK_FLUX_AS_MAG", PS_DATA_F32, "Peak flux expressed as magnitude", outputs.peakMag); 131 132 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "CAL_PSF_MAG", PS_DATA_F32, "PSF Magnitude using supplied calibration", outputs.calMag); 133 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "CAL_PSF_MAG_SIG", PS_DATA_F32, "measured scatter of zero point calibration", zeroptErr); 126 @PS1_V1@ psMetadataAdd (row, PS_LIST_TAIL, "RA_PSF", PS_DATA_F32, "PSF RA coordinate (degrees)", outputs.ra); 127 @PS1_V1@ psMetadataAdd (row, PS_LIST_TAIL, "DEC_PSF", PS_DATA_F32, "PSF DEC coordinate (degrees)", outputs.dec); 128 129 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "POSANGLE", PS_DATA_F32, "position angle at source (degrees)", outputs.posAngle); 130 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PLTSCALE", PS_DATA_F32, "plate scale at source (arcsec/pixel)", outputs.pltScale); 131 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_MAG", PS_DATA_F32, "PSF fit instrumental magnitude", source->psfMag); 132 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_MAG_SIG", PS_DATA_F32, "Sigma of PSF instrumental magnitude", source->psfMagErr); 133 134 @ALL,!PS1_V1,!PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX", PS_DATA_F32, "PSF fit instrumental flux (counts)", source->psfFlux); 135 @ALL,!PS1_V1,!PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_INST_FLUX_SIG",PS_DATA_F32, "Sigma of PSF instrumental flux", source->psfFluxErr); 136 137 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG", PS_DATA_F32, "magnitude in standard aperture", source->apMag); 138 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RAW", PS_DATA_F32, "magnitude in reported aperture", source->apMagRaw); 139 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "AP_MAG_RADIUS", PS_DATA_F32, "radius used for aperture mags", outputs.apRadius); 140 @PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "AP_FLUX", PS_DATA_F32, "instrumental flux in standard aperture", source->apFlux); 141 @PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "AP_FLUX_SIG", PS_DATA_F32, "aperture flux error", source->apFluxErr); 142 143 @<PS1_V3,PS1_SV1,PS1_DV?@ psMetadataAdd (row, PS_LIST_TAIL, "PEAK_FLUX_AS_MAG", PS_DATA_F32, "Peak flux expressed as magnitude", outputs.peakMag); 144 145 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "CAL_PSF_MAG", PS_DATA_F32, "PSF Magnitude using supplied calibration", outputs.calMag); 146 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "CAL_PSF_MAG_SIG", PS_DATA_F32, "measured scatter of zero point calibration", zeroptErr); 134 147 135 148 // NOTE: RA & DEC (both double) need to be on an 8-byte boundary... 136 @>PS1_V1@ psMetadataAdd (row, PS_LIST_TAIL, "RA_PSF", PS_DATA_F64, "PSF RA coordinate (degrees)", outputs.ra); 137 @>PS1_V1@ psMetadataAdd (row, PS_LIST_TAIL, "DEC_PSF", PS_DATA_F64, "PSF DEC coordinate (degrees)", outputs.dec); 138 139 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "PEAK_FLUX_AS_MAG", PS_DATA_F32, "Peak flux expressed as magnitude", outputs.peakMag); 140 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "SKY", PS_DATA_F32, "Sky level", source->sky); 141 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "SKY_SIGMA", PS_DATA_F32, "Sigma of sky level", source->skyErr); 142 143 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_CHISQ", PS_DATA_F32, "Chisq of PSF-fit", outputs.chisq); 144 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "CR_NSIGMA", PS_DATA_F32, "Nsigma deviations from PSF to CF", source->crNsigma); 145 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "EXT_NSIGMA", PS_DATA_F32, "Nsigma deviations from PSF to EXT", source->extNsigma); 146 147 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_MAJOR", PS_DATA_F32, "PSF width (major axis)", outputs.psfMajor); 148 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_MINOR", PS_DATA_F32, "PSF width (minor axis)", outputs.psfMinor); 149 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_THETA", PS_DATA_F32, "PSF orientation angle", outputs.psfTheta); 150 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF", PS_DATA_F32, "PSF coverage/quality factor (bad)", source->pixWeightNotBad); 151 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF_PERFECT", PS_DATA_F32, "PSF coverage/quality factor (poor)", source->pixWeightNotPoor); 152 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_NDOF", PS_DATA_S32, "degrees of freedom", outputs.nDOF); 153 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_NPIX", PS_DATA_S32, "number of pixels in fit", outputs.nPix); 154 155 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_XX", PS_DATA_F32, "second moments (X^2)", moments.Mxx); 156 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_XY", PS_DATA_F32, "second moments (X*Y)", moments.Mxy); 157 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_YY", PS_DATA_F32, "second moments (Y*Y)", moments.Myy); 158 159 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M3C", PS_DATA_F32, "third momemt cos theta", moments.M_c3); 160 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M3S", PS_DATA_F32, "third momemt sin theta", moments.M_s3); 161 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M4C", PS_DATA_F32, "fourth momemt cos theta", moments.M_c4); 162 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M4S", PS_DATA_F32, "fourth momemt sin theta", moments.M_s4); 163 164 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_R1", PS_DATA_F32, "first radial moment", moments.Mrf); 165 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_RH", PS_DATA_F32, "half radial moment", moments.Mrh); 166 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX", PS_DATA_F32, "Kron Flux (in 2.5 R1)", moments.Krf); 167 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_ERR", PS_DATA_F32, "Kron Flux Error", moments.dKrf); 168 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_INNER", PS_DATA_F32, "Kron Flux (in 2.5 R1)", moments.Kinner); 169 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_OUTER", PS_DATA_F32, "Kron Flux (in 2.5 R1)", moments.Kouter); 170 171 @>PS1_V3@ psMetadataAdd (row, PS_LIST_TAIL, "SKY_LIMIT_RAD", PS_DATA_F32, "Radius where object hits sky", source->skyRadius); 172 @>PS1_V3@ psMetadataAdd (row, PS_LIST_TAIL, "SKY_LIMIT_FLUX", PS_DATA_F32, "Flux / pix where object hits sky", source->skyFlux); 173 @>PS1_V3@ psMetadataAdd (row, PS_LIST_TAIL, "SKY_LIMIT_SLOPE", PS_DATA_F32, "d(Flux/pix)/dRadius where object hits sky", source->skySlope); 174 175 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "FLAGS", PS_DATA_U32, "psphot analysis flags", source->mode); 176 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "FLAGS2", PS_DATA_U32, "psphot analysis flags", source->mode2); 177 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "PADDING2", PS_DATA_S32, "more padding", 0); 178 179 // XXX not sure how to get this : need to load Nimages with weight? 180 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "N_FRAMES", PS_DATA_U16, "Number of frames overlapping source center", source->nFrames); 181 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PADDING", PS_DATA_S16, "padding", 0); 149 @ALL,!PS1_V1@ psMetadataAdd (row, PS_LIST_TAIL, "RA_PSF", PS_DATA_F64, "PSF RA coordinate (degrees)", outputs.ra); 150 @ALL,!PS1_V1@ psMetadataAdd (row, PS_LIST_TAIL, "DEC_PSF", PS_DATA_F64, "PSF DEC coordinate (degrees)", outputs.dec); 151 152 @>=PS1_V3@ psMetadataAdd (row, PS_LIST_TAIL, "PEAK_FLUX_AS_MAG", PS_DATA_F32, "Peak flux expressed as magnitude", outputs.peakMag); 153 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "SKY", PS_DATA_F32, "Sky level", source->sky); 154 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "SKY_SIGMA", PS_DATA_F32, "Sigma of sky level", source->skyErr); 155 156 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_CHISQ", PS_DATA_F32, "Chisq of PSF-fit", outputs.chisq); 157 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "CR_NSIGMA", PS_DATA_F32, "Nsigma deviations from PSF to CF", source->crNsigma); 158 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "EXT_NSIGMA", PS_DATA_F32, "Nsigma deviations from PSF to EXT", source->extNsigma); 159 160 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_MAJOR", PS_DATA_F32, "PSF width (major axis)", outputs.psfMajor); 161 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_MINOR", PS_DATA_F32, "PSF width (minor axis)", outputs.psfMinor); 162 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_THETA", PS_DATA_F32, "PSF orientation angle", outputs.psfTheta); 163 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF", PS_DATA_F32, "PSF coverage/quality factor (bad)", source->pixWeightNotBad); 164 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_QF_PERFECT", PS_DATA_F32, "PSF coverage/quality factor (poor)", source->pixWeightNotPoor); 165 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_NDOF", PS_DATA_S32, "degrees of freedom", outputs.nDOF); 166 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PSF_NPIX", PS_DATA_S32, "number of pixels in fit", outputs.nPix); 167 168 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_XX", PS_DATA_F32, "second moments (X^2)", moments.Mxx); 169 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_XY", PS_DATA_F32, "second moments (X*Y)", moments.Mxy); 170 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_YY", PS_DATA_F32, "second moments (Y*Y)", moments.Myy); 171 172 @>PS1_V2,PS1_SV1@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M3C", PS_DATA_F32, "third momemt cos theta", moments.M_c3); 173 @>PS1_V2,PS1_SV1@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M3S", PS_DATA_F32, "third momemt sin theta", moments.M_s3); 174 @>PS1_V2,PS1_SV1@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M4C", PS_DATA_F32, "fourth momemt cos theta", moments.M_c4); 175 @>PS1_V2,PS1_SV1@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_M4S", PS_DATA_F32, "fourth momemt sin theta", moments.M_s4); 176 177 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_R1", PS_DATA_F32, "first radial moment", moments.Mrf); 178 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "MOMENTS_RH", PS_DATA_F32, "half radial moment", moments.Mrh); 179 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX", PS_DATA_F32, "Kron Flux (in 2.5 R1)", moments.Krf); 180 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_ERR", PS_DATA_F32, "Kron Flux Error", moments.dKrf); 181 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_INNER", PS_DATA_F32, "Kron Flux (in 1.0 R1)", moments.Kinner); 182 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "KRON_FLUX_OUTER", PS_DATA_F32, "Kron Flux (in 4.0 R1)", moments.Kouter); 183 184 @>PS1_V3@ psMetadataAdd (row, PS_LIST_TAIL, "SKY_LIMIT_RAD", PS_DATA_F32, "Radius where object hits sky", source->skyRadius); 185 @>PS1_V3@ psMetadataAdd (row, PS_LIST_TAIL, "SKY_LIMIT_FLUX", PS_DATA_F32, "Flux / pix where object hits sky", source->skyFlux); 186 @>PS1_V3@ psMetadataAdd (row, PS_LIST_TAIL, "SKY_LIMIT_SLOPE", PS_DATA_F32, "d(Flux/pix)/dRadius where object hits sky", source->skySlope); 187 188 @PS1_DV?@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_NPOS", PS_DATA_S32, "nPos (n pix > 3 sigma)", diffStats.nGood); 189 @PS1_DV?@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_FRATIO", PS_DATA_F32, "fPos / (fPos + fNeg)", diffStats.fRatio); 190 @PS1_DV?@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_NRATIO_BAD", PS_DATA_F32, "nPos / (nPos + nNeg)", diffStats.nRatioBad); 191 @PS1_DV?@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_NRATIO_MASK", PS_DATA_F32, "nPos / (nPos + nMask)", diffStats.nRatioMask); 192 @PS1_DV?@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_NRATIO_ALL", PS_DATA_F32, "nPos / (nGood + nMask + nBad)", diffStats.nRatioAll); 193 194 @PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_R_P", PS_DATA_F32, "distance to positive match source", diffStats.Rp); 195 @PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_SN_P", PS_DATA_F32, "signal-to-noise of pos match src", diffStats.SNp); 196 @PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_R_M", PS_DATA_F32, "distance to negative match source", diffStats.Rm); 197 @PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "DIFF_SN_M", PS_DATA_F32, "signal-to-noise of neg match src", diffStats.SNm); 198 199 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "FLAGS", PS_DATA_U32, "psphot analysis flags", source->mode); 200 @>PS1_V2,PS1_SV1,PS1_DV2@ psMetadataAdd (row, PS_LIST_TAIL, "FLAGS2", PS_DATA_U32, "psphot analysis flags", source->mode2); 201 @>PS1_V2@ psMetadataAdd (row, PS_LIST_TAIL, "PADDING2", PS_DATA_S32, "more padding", 0); 202 @PS1_SV1@ 203 204 // note that this definition is inconsistent with the definition in 205 // Ohana/src/libautocode/def. This version creates a table with data not 206 // properly aligned with the 8-byte boundaries. The structure defined by 207 // libautocode does this, but has a different order of elements (and adds 208 // padding2 to fix things). We have generated may files with PS1_SV1 as is, so 209 // I'll leave it. But in future a PS1_SV2 should be forced to match 210 // libautocode. Note that addstar knows to detect the alternate version of 211 // PS1_SV1 and correctly interpret its fields. 212 213 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "N_FRAMES", PS_DATA_U16, "Number of frames overlapping source center", source->nFrames); 214 @ALL@ psMetadataAdd (row, PS_LIST_TAIL, "PADDING", PS_DATA_S16, "padding", 0); 182 215 183 216 psArrayAdd (table, 100, row); … … 282 315 @ALL@ source->psfMagErr = psMetadataLookupF32 (&status, row, "PSF_INST_MAG_SIG"); 283 316 @ALL@ source->apMag = psMetadataLookupF32 (&status, row, "AP_MAG"); 284 @>PS1_V2 @ source->apMagRaw = psMetadataLookupF32 (&status, row, "AP_MAG_RAW");317 @>PS1_V2,PS1_SV1,PS1_DV2@ source->apMagRaw = psMetadataLookupF32 (&status, row, "AP_MAG_RAW"); 285 318 286 319 // XXX use these to determine PAR[PM_PAR_I0] if they exist? 287 @>PS1_V2@ source->psfFlux = psMetadataLookupF32 (&status, row, "PSF_INST_FLUX"); 288 @>PS1_V2@ source->psfFluxErr= psMetadataLookupF32 (&status, row, "PSF_INST_FLUX_SIG"); 320 // XXX add these to PS1_SV1? 321 @>PS1_V2,PS1_SV1,PS1_DV?@ source->psfFlux = psMetadataLookupF32 (&status, row, "PSF_INST_FLUX"); 322 @>PS1_V2,PS1_SV1,PS1_DV?@ source->psfFluxErr= psMetadataLookupF32 (&status, row, "PSF_INST_FLUX_SIG"); 289 323 290 324 // XXX this scaling is incorrect: does not include the 2 \pi AREA factor … … 307 341 308 342 @ALL@ source->pixWeightNotBad = psMetadataLookupF32 (&status, row, "PSF_QF"); 309 @>PS1_V2@ source->pixWeightNotPoor = psMetadataLookupF32 (&status, row, "PSF_QF_PERFECT");343 @>PS1_V2,PS1_SV1,PS1_DV2@ source->pixWeightNotPoor = psMetadataLookupF32 (&status, row, "PSF_QF_PERFECT"); 310 344 @ALL@ source->crNsigma = psMetadataLookupF32 (&status, row, "CR_NSIGMA"); 311 345 @ALL@ source->extNsigma = psMetadataLookupF32 (&status, row, "EXT_NSIGMA"); … … 325 359 @ALL@ source->moments->Myy = psMetadataLookupF32 (&status, row, "MOMENTS_YY"); 326 360 327 @>PS1_V2@ source->moments->Mrf = psMetadataLookupF32 (&status, row, "MOMENTS_R1"); 328 @>PS1_V2@ source->moments->Mrh = psMetadataLookupF32 (&status, row, "MOMENTS_RH"); 329 @>PS1_V2@ source->moments->KronFlux = psMetadataLookupF32 (&status, row, "KRON_FLUX"); 330 @>PS1_V2@ source->moments->KronFluxErr = psMetadataLookupF32 (&status, row, "KRON_FLUX_ERR"); 331 332 @>PS1_V2@ source->moments->KronFinner = psMetadataLookupF32 (&status, row, "KRON_FLUX_INNER"); 333 @>PS1_V2@ source->moments->KronFouter = psMetadataLookupF32 (&status, row, "KRON_FLUX_OUTER"); 361 // XXX we do not save all of the 3rd and 4th moment parameters. when we load in data, 362 // we are storing enough information so the output will be consistent with the input 363 @>PS1_V2,PS1_SV1@ source->moments->Mxxx = +1.0 * psMetadataLookupF32 (&status, row, "MOMENTS_M3C"); 364 @>PS1_V2,PS1_SV1@ source->moments->Mxxy = 0.0; 365 @>PS1_V2,PS1_SV1@ source->moments->Mxyy = 0.0; 366 @>PS1_V2,PS1_SV1@ source->moments->Myyy = -1.0 * psMetadataLookupF32 (&status, row, "MOMENTS_M3S"); 367 368 @>PS1_V2,PS1_SV1@ source->moments->Mxxxx = +1.00 * psMetadataLookupF32 (&status, row, "MOMENTS_M4C"); 369 @>PS1_V2,PS1_SV1@ source->moments->Mxxxy = 0.0; 370 @>PS1_V2,PS1_SV1@ source->moments->Mxxyy = 0.0; 371 @>PS1_V2,PS1_SV1@ source->moments->Mxyyy = -0.25 * psMetadataLookupF32 (&status, row, "MOMENTS_M4S"); 372 @>PS1_V2,PS1_SV1@ source->moments->Myyyy = 0.0; 373 374 @>PS1_V2,PS1_SV1,PS1_DV2@ source->moments->Mrf = psMetadataLookupF32 (&status, row, "MOMENTS_R1"); 375 @>PS1_V2,PS1_SV1,PS1_DV2@ source->moments->Mrh = psMetadataLookupF32 (&status, row, "MOMENTS_RH"); 376 @>PS1_V2,PS1_SV1,PS1_DV2@ source->moments->KronFlux = psMetadataLookupF32 (&status, row, "KRON_FLUX"); 377 @>PS1_V2,PS1_SV1,PS1_DV2@ source->moments->KronFluxErr = psMetadataLookupF32 (&status, row, "KRON_FLUX_ERR"); 378 @>PS1_V2,PS1_SV1,PS1_DV2@ source->moments->KronFinner = psMetadataLookupF32 (&status, row, "KRON_FLUX_INNER"); 379 @>PS1_V2,PS1_SV1,PS1_DV2@ source->moments->KronFouter = psMetadataLookupF32 (&status, row, "KRON_FLUX_OUTER"); 334 380 335 381 @>PS1_V3@ source->skyRadius = psMetadataLookupF32 (&status, row, "SKY_LIMIT_RAD"); … … 337 383 @>PS1_V3@ source->skySlope = psMetadataLookupF32 (&status, row, "SKY_LIMIT_SLOPE"); 338 384 339 // XXX we do not save all of the 3rd and 4th moment parameters. when we load in data, 340 // we are storing enough information so the output will be consistent with the input 341 @>PS1_V2@ source->moments->Mxxx = +1.0 * psMetadataLookupF32 (&status, row, "MOMENTS_M3C"); 342 @>PS1_V2@ source->moments->Mxxy = 0.0; 343 @>PS1_V2@ source->moments->Mxyy = 0.0; 344 @>PS1_V2@ source->moments->Myyy = -1.0 * psMetadataLookupF32 (&status, row, "MOMENTS_M3S"); 345 346 @>PS1_V2@ source->moments->Mxxxx = +1.00 * psMetadataLookupF32 (&status, row, "MOMENTS_M4C"); 347 @>PS1_V2@ source->moments->Mxxxy = 0.0; 348 @>PS1_V2@ source->moments->Mxxyy = 0.0; 349 @>PS1_V2@ source->moments->Mxyyy = -0.25 * psMetadataLookupF32 (&status, row, "MOMENTS_M4S"); 350 @>PS1_V2@ source->moments->Myyyy = 0.0; 351 352 @ALL@ source->mode = psMetadataLookupU32 (&status, row, "FLAGS"); 353 @>PS1_V2@ source->mode2 = psMetadataLookupU32 (&status, row, "FLAGS2"); 385 @PS1_DV?@ int nPos = psMetadataLookupS32 (&status, row, "DIFF_NPOS"); 386 @PS1_DV?@ if (nPos) { 387 @PS1_DV?@ source->diffStats = pmSourceDiffStatsAlloc(); 388 @PS1_DV?@ source->diffStats->nGood = nPos; 389 @PS1_DV?@ source->diffStats->fRatio = psMetadataLookupF32 (&status, row, "DIFF_FRATIO"); 390 @PS1_DV?@ source->diffStats->nRatioBad = psMetadataLookupF32 (&status, row, "DIFF_NRATIO_BAD"); 391 @PS1_DV?@ source->diffStats->nRatioMask = psMetadataLookupF32 (&status, row, "DIFF_NRATIO_MASK"); 392 @PS1_DV?@ source->diffStats->nRatioAll = psMetadataLookupF32 (&status, row, "DIFF_NRATIO_ALL"); 393 394 @PS1_DV2@ source->diffStats->Rp = psMetadataLookupF32 (&status, row, "DIFF_R_P"); 395 @PS1_DV2@ source->diffStats->SNp = psMetadataLookupF32 (&status, row, "DIFF_SN_P"); 396 @PS1_DV2@ source->diffStats->Rm = psMetadataLookupF32 (&status, row, "DIFF_R_M"); 397 @PS1_DV2@ source->diffStats->SNm = psMetadataLookupF32 (&status, row, "DIFF_SN_M"); 398 @PS1_DV?@ } 399 400 @ALL@ source->mode = psMetadataLookupU32 (&status, row, "FLAGS"); 401 @>PS1_V2,PS1_SV1,PS1_DV2@ source->mode2 = psMetadataLookupU32 (&status, row, "FLAGS2"); 354 402 assert (status); 355 403 … … 409 457 // write the radial profile apertures to header 410 458 for (int i = 0; i < radMax->n; i++) { 411 412 413 414 459 sprintf (keyword1, "RMIN_%02d", i); 460 sprintf (keyword2, "RMAX_%02d", i); 461 psMetadataAddF32 (outhead, PS_LIST_TAIL, keyword1, PS_META_REPLACE, "min radius for SB profile", radMin->data.F32[i]); 462 psMetadataAddF32 (outhead, PS_LIST_TAIL, keyword2, PS_META_REPLACE, "min radius for SB profile", radMax->data.F32[i]); 415 463 } 416 464 417 465 // we write out all sources, regardless of quality. the source flags tell us the state 418 466 for (int i = 0; i < sources->n; i++) { 419 pmSource *source = sources->data[i]; 467 // this is the source associated with this image 468 pmSource *thisSource = sources->data[i]; 469 470 // this is the "real" version of this source 471 pmSource *source = thisSource->parent ? thisSource->parent : thisSource; 420 472 421 473 // skip sources without measurements … … 470 522 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90", PS_DATA_F32, "Petrosian R90 (pix)", extpars->petrosianR90); 471 523 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90_ERR", PS_DATA_F32, "Petrosian R90 Error (pix)", extpars->petrosianR90Err); 524 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_FILL", PS_DATA_F32, "Petrosian Fill Factor", extpars->petrosianFill); 472 525 } else { 473 526 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG", PS_DATA_F32, "Petrosian Magnitude", NAN); … … 479 532 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90", PS_DATA_F32, "Petrosian R90 (pix)", NAN); 480 533 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90_ERR", PS_DATA_F32, "Petrosian R90 Error (pix)",NAN); 534 psMetadataAdd (row, PS_LIST_TAIL, "PETRO_FILL", PS_DATA_F32, "Petrosian Fill Factor", NAN); 481 535 } 482 536 } … … 539 593 psError(psErrorCodeLast(), false, "writing ext data %s\n", extname); 540 594 psFree (outhead); 541 542 595 psFree(table); 596 return false; 543 597 } 544 598 psFree (outhead); 545 599 psFree (table); 546 600 601 return true; 602 } 603 604 bool pmSourcesRead_CMF_@CMFMODE@_XSRC(psFits *fits, psMetadata *hduHeader, psArray *sources, long *sourceIndex) 605 { 606 PS_ASSERT_PTR_NON_NULL(fits, false); 607 PS_ASSERT_PTR_NON_NULL(sources, false); 608 609 bool status; 610 long numSources = psFitsTableSize(fits); // Number of sources in table 611 if (numSources == 0) { 612 psError(psErrorCodeLast(), false, "XSRC Table contains no entries\n"); 613 return false; 614 } 615 616 // petrosian mags are not saved, we need to calculate fluxes. For this we need exptime and zero point 617 float zeropt = psMetadataLookupF32(&status, hduHeader, "FPA.ZP"); 618 float exptime = psMetadataLookupF32(&status, hduHeader, "EXPTIME"); 619 float magOffset = zeropt + 2.5*log10(exptime); 620 621 for (long i = 0; i < numSources; i++) { 622 psMetadata *row = psFitsReadTableRow(fits, i); // Table row 623 if (!row) { 624 psError(psErrorCodeLast(), false, "Unable to read row %ld of sources", i); 625 psFree(row); 626 return false; 627 } 628 // Find the source with this sequence number. 629 // XXX: I am assuming that sources is sorted in order of seq 630 long seq = psMetadataLookupU32 (&status, row, "IPP_IDET"); 631 pmSource *source = NULL; 632 #ifndef ASSUME_SORTED 633 long j = seq < sources->n ? seq : sources->n - 1; 634 for (; j >= 0; j--) { 635 source = sources->data[j]; 636 if (source->seq == seq) { 637 break; 638 } 639 } 640 #else 641 long j = sourceIndex[seq]; 642 psAssert(j >= 0 && j < sources->n, "invalid sourceIndex"); 643 source = sources->data[j]; 644 #endif 645 if (!source) { 646 psError(PS_ERR_UNKNOWN, false, "Failed to find source for row %ld sequence number %ld\n", i, seq); 647 psFree(row); 648 return false; 649 } 650 651 if (!source->extpars) { 652 source->extpars = pmSourceExtendedParsAlloc (); 653 } 654 pmSourceExtendedPars *extpars = source->extpars; 655 656 // Assume that X_EXT Y_EXT and sigmas match the psf src so skip 657 658 // We don't have enough information to calculate the major and minor axis. Set major to 1. Should we scale this by 659 // psf size or something? 660 extpars->axes.major = 1.0; 661 extpars->axes.minor = extpars->axes.major * psMetadataLookupF32(&status, row, "F25_ARATIO"); 662 extpars->axes.theta = psMetadataLookupF32(&status, row, "F25_THETA"); 663 664 float mag = psMetadataLookupF32(&status, row, "PETRO_MAG"); 665 float magErr = psMetadataLookupF32(&status, row, "PETRO_MAG_ERR"); 666 if (isfinite(mag)) { 667 extpars->petrosianFlux = pow(10., (magOffset - mag) / 2.5); 668 if (isfinite(magErr)) { 669 extpars->petrosianFluxErr = extpars->petrosianFlux / magErr; 670 } 671 } 672 673 extpars->petrosianRadius = psMetadataLookupF32(&status, row, "PETRO_RADIUS"); 674 extpars->petrosianRadiusErr= psMetadataLookupF32(&status, row, "PETRO_RADIUS_ERR"); 675 extpars->petrosianR50 = psMetadataLookupF32(&status, row, "PETRO_RADIUS_50"); 676 extpars->petrosianR50Err = psMetadataLookupF32(&status, row, "PETRO_RADIUS_50_ERR"); 677 extpars->petrosianR90 = psMetadataLookupF32(&status, row, "PETRO_RADIUS_90"); 678 extpars->petrosianR90Err = psMetadataLookupF32(&status, row, "PETRO_RADIUS_90_ERR"); 679 extpars->petrosianFill = psMetadataLookupF32(&status, row, "PETRO_FILL"); 680 681 psVector *radSB = psMetadataLookupVector(&status, row, "PROF_SB"); 682 psVector *radFlux = psMetadataLookupVector(&status, row, "PROF_FLUX"); 683 psVector *radFill = psMetadataLookupVector(&status, row, "PROF_FILL"); 684 685 if (radSB && radSB->n > 0) { 686 extpars->radProfile = pmSourceRadialProfileAlloc(); 687 extpars->radProfile->binSB = psMemIncrRefCounter(radSB); 688 extpars->radProfile->binSum = psMemIncrRefCounter(radFlux); 689 extpars->radProfile->binFill = psMemIncrRefCounter(radFill); 690 } 691 692 psFree(row); 693 } 694 547 695 return true; 548 696 } … … 572 720 int nParamMax = 0; 573 721 for (int i = 0; i < sources->n; i++) { 574 pmSource *source = sources->data[i]; 722 // this is the source associated with this image 723 pmSource *thisSource = sources->data[i]; 724 725 // this is the "real" version of this source 726 pmSource *source = thisSource->parent ? thisSource->parent : thisSource; 727 575 728 if (source->modelFits == NULL) continue; 576 729 for (int j = 0; j < source->modelFits->n; j++) { … … 586 739 for (int i = 0; i < sources->n; i++) { 587 740 588 pmSource *source = sources->data[i]; 741 pmSource *thisSource = sources->data[i]; 742 743 // this is the "real" version of this source 744 pmSource *source = thisSource->parent ? thisSource->parent : thisSource; 589 745 590 746 // XXX if no model fits are saved, write out modelEXT? … … 641 797 642 798 // XXX these should be major and minor, not 'x' and 'y' 643 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MAJ", 0, "EXT width in x coordinate", axes.major); 644 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MIN", 0, "EXT width in y coordinate", axes.minor); 645 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_THETA", 0, "EXT orientation angle", axes.theta); 799 if (model->type == pmModelClassGetType("PS_MODEL_TRAIL")) { 800 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MAJ", 0, "EXT width (major axis), length for trail", PAR[PM_PAR_LENGTH]); 801 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MIN", 0, "EXT width (minor axis), sigma for trail", PAR[PM_PAR_SIGMA]); 802 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_THETA", 0, "EXT orientation angle", PAR[PM_PAR_THETA]); 803 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MAJ_ERR",0, "EXT width error (major axis)", dPAR[PM_PAR_LENGTH]); 804 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MIN_ERR",0, "EXT width error (minor axis)", dPAR[PM_PAR_SIGMA]); 805 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_THETA_ERR", 0, "EXT orientation angle (error)", dPAR[PM_PAR_THETA]); 806 } else { 807 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MAJ", 0, "EXT width (major axis), length for trail", axes.major); 808 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MIN", 0, "EXT width (minor axis), sigma for trail", axes.minor); 809 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_THETA", 0, "EXT orientation angle", axes.theta); 810 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MAJ_ERR",0, "EXT width error (major axis)", dPAR[PM_PAR_LENGTH]); 811 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_WIDTH_MIN_ERR",0, "EXT width error (minor axis)", dPAR[PM_PAR_SIGMA]); 812 psMetadataAddF32 (row, PS_LIST_TAIL, "EXT_THETA_ERR", 0, "EXT orientation angle (error)", dPAR[PM_PAR_THETA]); 813 } 646 814 647 815 // write out the other generic parameters … … 658 826 659 827 if (k < model->params->n) { 660 psMetadataAdd (row, PS_LIST_TAIL, name, PS_DATA_F32, "", model->params->data.F32[k]);828 psMetadataAddF32 (row, PS_LIST_TAIL, name, 0, "", model->params->data.F32[k]); 661 829 } else { 662 psMetadataAddF32 (row, PS_LIST_TAIL, name, PS_DATA_F32, "", NAN);830 psMetadataAddF32 (row, PS_LIST_TAIL, name, 0, "", NAN); 663 831 } 664 832 } 665 833 666 // XXX other parameters which may be set. 667 // XXX flag / value to define the model 668 // XXX write out the model type, fit status flags 669 834 // optionally, write out the covariance matrix values 835 // XXX do I need to pad this to match the biggest covar matrix? 836 if (model->covar) { 837 for (int iy = 0; iy < model->covar->numCols; iy++) { 838 for (int ix = iy; ix < model->covar->numCols; ix++) { 839 snprintf (name, 64, "EXT_COVAR_%02d_%02d", iy, ix); 840 psMetadataAddF32 (row, PS_LIST_TAIL, name, 0, "", model->covar->data.F32[iy][ix]); 841 842 } 843 } 844 } 670 845 psArrayAdd (table, 100, row); 671 846 psFree (row); … … 697 872 } 698 873 874 bool pmSourcesRead_CMF_@CMFMODE@_XFIT(psFits *fits, psMetadata *hduHeader, psArray *sources, long *sourceIndex) 875 { 876 PS_ASSERT_PTR_NON_NULL(fits, false); 877 PS_ASSERT_PTR_NON_NULL(sources, false); 878 879 bool status; 880 long numSources = psFitsTableSize(fits); // Number of sources in table 881 if (numSources == 0) { 882 psError(psErrorCodeLast(), false, "XFIT Table contains no entries\n"); 883 return false; 884 } 885 886 for (long i = 0; i < numSources; i++) { 887 psMetadata *row = psFitsReadTableRow(fits, i); // Table row 888 if (!row) { 889 psError(psErrorCodeLast(), false, "Unable to read row %ld of sources", i); 890 psFree(row); 891 return false; 892 } 893 // Find the source with this sequence number. 894 // XXX: I am assuming that sources is sorted in order of seq. 895 long seq = psMetadataLookupU32 (&status, row, "IPP_IDET"); 896 long j = seq < sources->n ? seq : sources->n - 1; 897 pmSource *source = NULL; 898 for (; j >= 0; j--) { 899 source = sources->data[j]; 900 if (source->seq == seq) { 901 break; 902 } 903 } 904 if (!source) { 905 psError(PS_ERR_UNKNOWN, false, "Failed to find source for row %ld sequence number %ld\n", i, seq); 906 psFree(row); 907 return false; 908 } 909 if (!source->modelFits) { 910 // XXX: where to find the number of models to expect? 911 source->modelFits = psArrayAllocEmpty(5); 912 } 913 psString modelName = psMetadataLookupStr(&status, row, "MODEL_TYPE"); 914 if (!modelName) { 915 psError(PS_ERR_UNKNOWN, true, "Failed to find model name for row %ld\n", i); 916 psFree(row); 917 return false; 918 } 919 pmModelType modelType = pmModelClassGetType(modelName); 920 if (modelType < 0) { 921 psError(PS_ERR_UNKNOWN, true, "Failed to find model type for %s\n", modelName); 922 psFree(row); 923 return false; 924 } 925 pmModel *model = pmModelAlloc(modelType); 926 927 psF32 *PAR = model->params->data.F32; 928 psF32 *dPAR = model->dparams->data.F32; 929 930 PAR[PM_PAR_XPOS] = psMetadataLookupF32(&status, row, "X_EXT"); 931 PAR[PM_PAR_YPOS] = psMetadataLookupF32(&status, row, "Y_EXT"); 932 dPAR[PM_PAR_XPOS] = psMetadataLookupF32(&status, row, "X_EXT_SIG"); 933 dPAR[PM_PAR_YPOS] = psMetadataLookupF32(&status, row, "Y_EXT_SIG"); 934 935 model->mag = psMetadataLookupF32(&status, row, "EXT_INST_MAG"); 936 model->magErr = psMetadataLookupF32(&status, row, "EXT_INST_MAG_SIG"); 937 938 psEllipseAxes axes; 939 axes.major = psMetadataLookupF32(&status, row, "EXT_WIDTH_MAJ"); 940 axes.minor = psMetadataLookupF32(&status, row, "EXT_WIDTH_MIN"); 941 axes.theta = psMetadataLookupF32(&status, row, "EXT_THETA"); 942 if (!pmPSF_AxesToModel(PAR, axes, modelType)) { 943 // Do we need to fail here or can this happen? 944 psError(PS_ERR_UNKNOWN, false, "Failed to convert psf axes to model"); 945 psFree(model); 946 psFree(row); 947 return false; 948 } 949 // XXX: clean this up 950 if (model->params->n > 7) { 951 PAR[7] = psMetadataLookupF32(&status, row, "EXT_PAR_07"); 952 } 953 // read the covariance matrix 954 int nparams = model->params->n; 955 psImage *covar = psImageAlloc(nparams, nparams, PS_TYPE_F32); 956 for (int y = 0; y < nparams; y++) { 957 for (int x = 0; x < nparams; x++) { 958 char name[64]; 959 snprintf(name, 64, "EXT_COVAR_%02d_%02d", y, x); 960 covar->data.F32[y][x] = psMetadataLookupF32(&status, row, name); 961 } 962 } 963 model->covar = covar; 964 965 psArrayAdd(source->modelFits, 1, model); 966 psFree(model); 967 968 psFree(row); 969 } 970 971 return true; 972 } 973 974 // **** write out the radial flux values for the sources for a given matched-PSF image 975 // **** how do we distinguish the matched-PSF images from the non-matched version 699 976 bool pmSourcesWrite_CMF_@CMFMODE@_XRAD(psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe) 700 977 { 978 bool status = false; 979 psArray *table; 980 psMetadata *row; 981 psF32 xPos, yPos; 982 char keyword1[80], keyword2[80]; 983 984 // perform full non-linear fits / extended source analysis? 985 if (!psMetadataLookupBool (&status, recipe, "RADIAL_APERTURES")) { 986 psLogMsg ("psphot", PS_LOG_INFO, "radial apertures were not measured, skipping\n"); 987 return true; 988 } 989 990 // create a header to hold the output data 991 psMetadata *outhead = psMetadataAlloc (); 992 993 // write the links to the image header 994 psMetadataAddStr (outhead, PS_LIST_TAIL, "EXTNAME", PS_META_REPLACE, "radial flux table extension", extname); 995 996 // we use this just to define the output vectors (which must be present for all objects) 997 psVector *radMin = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.LOWER"); 998 psVector *radMax = psMetadataLookupPtr (&status, recipe, "RADIAL.ANNULAR.BINS.UPPER"); 999 psAssert (radMax, "this must have been defined and tested earlier!"); 1000 psAssert (radMax->n, "this must have been defined and tested earlier!"); 1001 psAssert (radMin->n == radMax->n, "inconsistent annular bins"); 1002 1003 // write the radial profile apertures to header 1004 for (int i = 0; i < radMax->n; i++) { 1005 sprintf (keyword1, "RMIN_%02d", i); 1006 sprintf (keyword2, "RMAX_%02d", i); 1007 psMetadataAddF32 (outhead, PS_LIST_TAIL, keyword1, PS_META_REPLACE, "min radius for SB profile", radMin->data.F32[i]); 1008 psMetadataAddF32 (outhead, PS_LIST_TAIL, keyword2, PS_META_REPLACE, "min radius for SB profile", radMax->data.F32[i]); 1009 } 1010 1011 // the FWHM values are available if we measured a psf-matched convolved set 1012 psVector *fwhmValues = psMetadataLookupVector(&status, readout->analysis, "STACK.PSF.FWHM.VALUES"); 1013 1014 // let's write these out in S/N order 1015 sources = psArraySort (sources, pmSourceSortByFlux); 1016 1017 table = psArrayAllocEmpty (sources->n); 1018 1019 // we write out all sources, regardless of quality. the source flags tell us the state 1020 for (int i = 0; i < sources->n; i++) { 1021 1022 // this is the source associated with this image 1023 pmSource *thisSource = sources->data[i]; 1024 1025 // this is the "real" version of this source 1026 pmSource *source = thisSource->parent ? thisSource->parent : thisSource; 1027 1028 // skip sources without radial aper measurements (or insufficient) 1029 if (source->radialAper == NULL) continue; 1030 1031 // psAssert (source->radialAper->n == fwhmValues->n, "inconsistent radial aperture set"); 1032 1033 for (int entry = 0; entry < source->radialAper->n; entry++) { 1034 1035 // choose the convolved EXT model, if available, otherwise the simple one 1036 pmSourceRadialApertures *radialAper = source->radialAper->data[entry]; 1037 assert (radialAper); 1038 1039 if (pmSourcePositionUseMoments(source)) { 1040 xPos = source->moments->Mx; 1041 yPos = source->moments->My; 1042 } else { 1043 xPos = source->peak->xf; 1044 yPos = source->peak->yf; 1045 } 1046 1047 row = psMetadataAlloc (); 1048 1049 // XXX we are not writing out the mode (flags) or the type (psf, ext, etc) 1050 psMetadataAddU32 (row, PS_LIST_TAIL, "IPP_IDET", 0, "IPP detection identifier index", source->seq); 1051 psMetadataAddF32 (row, PS_LIST_TAIL, "X_APER", 0, "Center of aperture measurements", xPos); 1052 psMetadataAddF32 (row, PS_LIST_TAIL, "Y_APER", 0, "Center of aperture measurements", yPos); 1053 if (fwhmValues) { 1054 psMetadataAddF32 (row, PS_LIST_TAIL, "PSF_FWHM", 0, "FWHM of matched PSF", fwhmValues->data.F32[entry]); 1055 } else { 1056 psMetadataAddF32 (row, PS_LIST_TAIL, "PSF_FWHM", 0, "image is not FWHM-matched", NAN); 1057 } 1058 1059 // XXX if we have raw radial apertures, write them out here 1060 psVector *radFlux = psVectorAlloc(radMax->n, PS_TYPE_F32); 1061 psVector *radFluxErr = psVectorAlloc(radMax->n, PS_TYPE_F32); 1062 psVector *radFill = psVectorAlloc(radMax->n, PS_TYPE_F32); 1063 psVector *radFluxStdev = psVectorAlloc(radMax->n, PS_TYPE_F32); 1064 psVectorInit (radFlux, NAN); 1065 psVectorInit (radFluxErr, NAN); 1066 psVectorInit (radFill, NAN); 1067 if (!radialAper->flux) goto write_annuli; 1068 if (!radialAper->fill) goto write_annuli; 1069 psAssert (radialAper->flux->n <= radFlux->n, "inconsistent vector lengths"); 1070 psAssert (radialAper->fill->n <= radFlux->n, "inconsistent vector lengths"); 1071 1072 // copy the data from fluxVal (which is not guaranteed to be the full length) to radFlux 1073 for (int j = 0; j < radialAper->flux->n; j++) { 1074 radFlux->data.F32[j] = radialAper->flux->data.F32[j]; 1075 radFluxErr->data.F32[j] = radialAper->fluxErr->data.F32[j]; 1076 radFluxStdev->data.F32[j] = radialAper->fluxStdev->data.F32[j]; 1077 radFill->data.F32[j] = radialAper->fill->data.F32[j]; 1078 } 1079 1080 write_annuli: 1081 psMetadataAdd (row, PS_LIST_TAIL, "APER_FLUX", PS_DATA_VECTOR, "flux within annuli", radFlux); 1082 psMetadataAdd (row, PS_LIST_TAIL, "APER_FLUX_ERR", PS_DATA_VECTOR, "flux error in annuli", radFluxErr); 1083 psMetadataAdd (row, PS_LIST_TAIL, "APER_FLUX_STDEV", PS_DATA_VECTOR, "flux standard deviation", radFluxStdev); 1084 psMetadataAdd (row, PS_LIST_TAIL, "APER_FILL", PS_DATA_VECTOR, "fill factor of annuli", radFill); 1085 psFree (radFlux); 1086 psFree (radFluxErr); 1087 psFree (radFluxStdev); 1088 psFree (radFill); 1089 1090 psArrayAdd (table, 100, row); 1091 psFree (row); 1092 } 1093 } 1094 1095 if (table->n == 0) { 1096 if (!psFitsWriteBlank (fits, outhead, extname)) { 1097 psError(psErrorCodeLast(), false, "Unable to write empty sources file."); 1098 psFree(outhead); 1099 psFree(table); 1100 return false; 1101 } 1102 psFree (outhead); 1103 psFree (table); 1104 return true; 1105 } 1106 1107 psTrace ("pmFPAfile", 5, "writing ext data %s\n", extname); 1108 if (!psFitsWriteTable (fits, outhead, table, extname)) { 1109 psError(psErrorCodeLast(), false, "writing ext data %s\n", extname); 1110 psFree (outhead); 1111 psFree(table); 1112 return false; 1113 } 1114 psFree (outhead); 1115 psFree (table); 701 1116 return true; 702 1117 } 1118 1119 bool pmSourcesRead_CMF_@CMFMODE@_XRAD(psFits *fits, pmReadout *readout, psMetadata *hduHeader, psArray *sources, long *sourceIndex) 1120 { 1121 PS_ASSERT_PTR_NON_NULL(fits, false); 1122 PS_ASSERT_PTR_NON_NULL(sources, false); 1123 1124 bool status; 1125 long numSources = psFitsTableSize(fits); // Number of sources in table 1126 if (numSources == 0) { 1127 psError(psErrorCodeLast(), false, "XRAD Table contains no entries\n"); 1128 return false; 1129 } 1130 1131 long seq_first = -1; 1132 long seq_last = -1; 1133 psVector *fwhmValues = psVectorAllocEmpty(10, PS_TYPE_F32); 1134 long max_entries = -1; 1135 long num_entries = -1; 1136 1137 for (long i = 0; i < numSources; i++) { 1138 psMetadata *row = psFitsReadTableRow(fits, i); // Table row 1139 if (!row) { 1140 psError(psErrorCodeLast(), false, "Unable to read row %ld of sources", i); 1141 psFree(row); 1142 return false; 1143 } 1144 // Find the source with this sequence number. 1145 // XXX: I am assuming that sources is sorted in order of seq. 1146 long seq = psMetadataLookupU32 (&status, row, "IPP_IDET"); 1147 long j = seq < sources->n ? seq : sources->n - 1; 1148 pmSource *source = NULL; 1149 for (; j >= 0; j--) { 1150 source = sources->data[j]; 1151 if (source->seq == seq) { 1152 break; 1153 } 1154 } 1155 if (!source) { 1156 psError(PS_ERR_UNKNOWN, false, "Failed to find source for row %ld sequence number %ld\n", i, seq); 1157 psFree(row); 1158 return false; 1159 } 1160 if (seq_first == -1) { 1161 seq_first = seq; 1162 } 1163 if (seq == seq_first) { 1164 psF32 value = psMetadataLookupF32(&status, row, "PSF_FWHM"); 1165 psVectorAppend(fwhmValues, value); 1166 } 1167 if (seq == seq_last) { 1168 num_entries++; 1169 } else { 1170 num_entries = 1; 1171 seq_last = seq; 1172 } 1173 if (num_entries > max_entries) { 1174 max_entries = num_entries; 1175 } 1176 1177 if (!source->radialAper) { 1178 // XXX: where to find the number of models to expect? 1179 source->radialAper = psArrayAllocEmpty(5); 1180 } 1181 pmSourceRadialApertures *radialAper = pmSourceRadialAperturesAlloc(); 1182 1183 radialAper->flux = psMemIncrRefCounter(psMetadataLookupVector(&status, row, "APER_FLUX")); 1184 radialAper->fluxStdev = psMemIncrRefCounter(psMetadataLookupVector(&status, row, "APER_FLUX_STDEV")); 1185 radialAper->fluxErr = psMemIncrRefCounter(psMetadataLookupVector(&status, row, "APER_FLUX_ERR")); 1186 radialAper->fill = psMemIncrRefCounter(psMetadataLookupVector(&status, row, "APER_FILL")); 1187 1188 psArrayAdd(source->radialAper, 1, radialAper); 1189 1190 psFree(radialAper); 1191 psFree(row); 1192 } 1193 1194 // check for consistency between the length of fwhmValues and the maximum number of entries for each row 1195 if (fwhmValues->n != max_entries) { 1196 psError(PS_ERR_PROGRAMMING, true, "number of PSF_FWHM values found %ld does not match expected number: %ld\n", 1197 fwhmValues->n, max_entries); 1198 psAssert(0, "fixme"); 1199 } 1200 1201 if (!readout->analysis) { 1202 readout->analysis = psMetadataAlloc(); 1203 } 1204 1205 psMetadataAddVector(readout->analysis, PS_LIST_TAIL, "STACK.PSF.FWHM.VALUES", PS_META_REPLACE, "PSF sizes", fwhmValues); 1206 psFree(fwhmValues); 1207 1208 return true; 1209 }
Note:
See TracChangeset
for help on using the changeset viewer.
