LSSTApplications  10.0-2-g4f67435,11.0.rc2+1,11.0.rc2+12,11.0.rc2+3,11.0.rc2+4,11.0.rc2+5,11.0.rc2+6,11.0.rc2+7,11.0.rc2+8
LSSTDataManagementBasePackage
applyApCorr.py
Go to the documentation of this file.
1 from __future__ import absolute_import, division
2 #
3 # LSST Data Management System
4 # Copyright 2015 LSST Corporation.
5 #
6 # This product includes software developed by the
7 # LSST Project (http://www.lsst.org/).
8 #
9 # This program is free software: you can redistribute it and/or modify
10 # it under the terms of the GNU General Public License as published by
11 # the Free Software Foundation, either version 3 of the License, or
12 # (at your option) any later version.
13 #
14 # This program is distributed in the hope that it will be useful,
15 # but WITHOUT ANY WARRANTY; without even the implied warranty of
16 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
17 # GNU General Public License for more details.
18 #
19 # You should have received a copy of the LSST License Statement and
20 # the GNU General Public License along with this program. If not,
21 # see <http://www.lsstcorp.org/LegalNotices/>.
22 #
23 import math
24 
25 import numpy
26 
27 import lsst.pex.config
29 import lsst.afw.image
30 import lsst.pipe.base
31 from .apCorrRegistry import getApCorrNameSet
32 
33 # If True then scale flux sigma by apCorr; if False then use a more complex computation
34 # that over-estimates flux error (often grossly so) because it double-counts photon noise.
35 # This flag is intended to be temporary until we figure out a better way to compute
36 # the effects of aperture correction on flux uncertainty
37 UseNaiveFluxSigma = True
38 
39 __all__ = ("ApplyApCorrConfig", "ApplyApCorrTask")
40 
41 class ApCorrInfo(object):
42  """!Catalog field names and keys needed to aperture correct a particular flux
43  """
44  def __init__(self, schema, name):
45  """!Construct an ApCorrInfo and add fields to the schema
46 
47  @param[in,out] schema source catalog schema;
48  three fields are used to generate keys:
49  - {name}_flux
50  - {name}_fluxSigma
51  - {name}_fluxFlag
52  three fields are added:
53  - {name}_apCorr
54  - {name}_apCorrSigma
55  - {name}_flag_apCorr
56  @param[in] name field name prefix for flux needing aperture correction, e.g. "base_PsfFlux"
57 
58  ApCorrInfo has the following attributes:
59  - name: field name prefix for flux needing aperture correction
60  - fluxName: name of flux field
61  - fluxSigmaName: name of flux sigma field
62  - fluxKey: key to flux field
63  - fluxSigmaKey: key to flux sigma field
64  - fluxFlagKey: key to flux flag field
65  - apCorrKey: key to new aperture correction field
66  - apCorrSigmaKey: key to new aperture correction sigma field
67  - apCorrFlagKey: key to new aperture correction flag field
68  """
69  self.name = name
70  self.fluxName = name + "_flux"
71  self.fluxSigmaName = name + "_fluxSigma"
72  self.fluxKey = schema.find(self.fluxName).key
73  self.fluxSigmaKey = schema.find(self.fluxSigmaName).key
74  self.fluxFlagKey = schema.find(name + "_flag").key
75  self.apCorrKey = schema.addField(
76  name + "_apCorr",
77  doc = "aperture correction applied to %s" % (name,),
78  type = float,
79  )
80  self.apCorrSigmaKey = schema.addField(
81  name + "_apCorrSigma",
82  doc = "aperture correction applied to %s" % (name,),
83  type = float,
84  )
85  self.apCorrFlagKey = schema.addField(
86  name + "_flag_apCorr",
87  doc = "set if unable to aperture correct %s" % (name,),
88  type = "Flag",
89  )
90 
91 class ApplyApCorrConfig(lsst.pex.config.Config):
92  ignoreList = lsst.pex.config.ListField(
93  doc = "flux measurement algorithms in getApCorrNameSet() to ignore;" +
94  " if a name is listed that does not appear in getApCorrNameSet() then a warning is logged",
95  dtype = str,
96  optional = False,
97  default = (),
98  )
99  doFlagApCorrFailures = lsst.pex.config.Field(
100  doc = "set the general failure flag for a flux when it cannot be aperture-corrected?",
101  dtype = bool,
102  default = True,
103  )
104 
105 
106 class ApplyApCorrTask(lsst.pipe.base.Task):
107  """!Apply aperture corrections
108  """
109  ConfigClass = ApplyApCorrConfig
110  _DefaultName = "applyApCorr"
111 
112  def __init__(self, schema, **kwds):
113  """Construct an instance of this task
114  """
115  lsst.pipe.base.Task.__init__(self, **kwds)
116 
117  self.apCorrInfoDict = dict()
118  apCorrNameSet = getApCorrNameSet()
119  ignoreSet = set(self.config.ignoreList)
120  missingNameSet = ignoreSet - set(apCorrNameSet)
121  if missingNameSet:
122  self.log.warn("Fields in ignoreList that are not in fluxCorrectList: %s" %
123  (sorted(list(missingNameSet)),))
124  for name in apCorrNameSet - ignoreSet:
125  if name + "_flux" in schema:
126  self.apCorrInfoDict[name] = ApCorrInfo(schema=schema, name=name)
127  else:
128  # if a field in the registry is missing from the schema, silently ignore it.
129  pass
130 
131  def run(self, catalog, apCorrMap):
132  """Apply aperture corrections to a catalog of sources
133 
134  @param[in,out] catalog catalog of sources
135  @param[in] apCorrMap aperture correction map (an lsst.afw.image.ApCorrMap)
136 
137  If you show debug-level log messages then you will see statistics for the effects of
138  aperture correction.
139  """
140  self.log.info("Applying aperture corrections to %d flux fields" % (len(self.apCorrInfoDict),))
141  if UseNaiveFluxSigma:
142  self.log.info("Use naive flux sigma computation")
143  else:
144  self.log.info("Use complex flux sigma computation that double-counts photon noise "
145  " and thus over-estimates flux uncertainty")
146  for apCorrInfo in self.apCorrInfoDict.itervalues():
147  apCorrModel = apCorrMap.get(apCorrInfo.fluxName)
148  apCorrSigmaModel = apCorrMap.get(apCorrInfo.fluxSigmaName)
149  if None in (apCorrModel, apCorrSigmaModel):
150  missingNames = [(apCorrInfo.fluxName, apCorrInfo.fluxSigmaName)[i]
151  for i, model in enumerate((apCorrModel, apCorrSigmaModel)) if model is None]
152  self.log.warn("Could not find %s in apCorrMap" % (" or ".join(missingNames),))
153  for source in catalog:
154  source.set(apCorrInfo.apCorrFlagKey, True)
155  continue
156 
157  for source in catalog:
158  center = source.getCentroid()
159  # say we've failed when we start; we'll unset these flags when we succeed
160  source.set(apCorrInfo.apCorrFlagKey, True)
161  oldFluxFlagState = False
162  if self.config.doFlagApCorrFailures:
163  oldFluxFlagState = source.get(apCorrInfo.fluxFlagKey)
164  source.set(apCorrInfo.fluxFlagKey, True)
165 
166  apCorr = 1.0
167  apCorrSigma = 0.0
168  try:
169  apCorr = apCorrModel.evaluate(center)
170  if not UseNaiveFluxSigma:
171  apCorrSigma = apCorrSigmaModel.evaluate(center)
172  except lsst.pex.exceptions.DomainError:
173  continue
174 
175  source.set(apCorrInfo.apCorrKey, apCorr)
176  source.set(apCorrInfo.apCorrSigmaKey, apCorrSigma)
177  if apCorr <= 0.0 or apCorrSigma < 0.0:
178  continue
179 
180  flux = source.get(apCorrInfo.fluxKey)
181  fluxSigma = source.get(apCorrInfo.fluxSigmaKey)
182  source.set(apCorrInfo.fluxKey, flux*apCorr)
183  if UseNaiveFluxSigma:
184  source.set(apCorrInfo.fluxSigmaKey, fluxSigma*apCorr)
185  else:
186  a = fluxSigma/flux
187  b = apCorrSigma/apCorr
188  source.set(apCorrInfo.fluxSigmaKey, abs(flux*apCorr)*math.sqrt(a*a + b*b))
189  source.set(apCorrInfo.apCorrFlagKey, False)
190  if self.config.doFlagApCorrFailures:
191  source.set(apCorrInfo.fluxFlagKey, oldFluxFlagState)
192 
193  if self.log.getThreshold() <= self.log.DEBUG:
194  # log statistics on the effects of aperture correction
195  apCorrArr = numpy.array([s.get(apCorrInfo.apCorrKey) for s in catalog])
196  apCorrSigmaArr = numpy.array([s.get(apCorrInfo.apCorrSigmaKey) for s in catalog])
197  self.log.logdebug("For flux field %r: mean apCorr=%s, stdDev apCorr=%s,"
198  " mean apCorrSigma=%s, stdDev apCorrSigma=%s for %s sources" %
199  (apCorrInfo.name, apCorrArr.mean(), apCorrArr.std(),
200  apCorrSigmaArr.mean(), apCorrSigmaArr.std(), len(catalog)))
201 
202 
def __init__
Construct an ApCorrInfo and add fields to the schema.
Definition: applyApCorr.py:44
def getApCorrNameSet
Return a copy of the set of field name prefixes for fluxes that should be aperture corrected...
Catalog field names and keys needed to aperture correct a particular flux.
Definition: applyApCorr.py:41