23"""Load reference catalog objects for input to FGCM.
25This task will load multi-band reference objects and apply color terms (if
26configured). This wrapper around LoadReferenceObjects task also allows loading
27by healpix pixel (the native pixelization of fgcm), and is self-contained so
28the task can be called by third-party code.
33from astropy
import units
35import lsst.pex.config
as pexConfig
36import lsst.pipe.base
as pipeBase
37from lsst.meas.algorithms
import ReferenceSourceSelectorTask
38from lsst.meas.algorithms
import getRefFluxField
39from lsst.pipe.tasks.colorterms
import ColortermLibrary
40from lsst.afw.image
import abMagErrFromFluxErr
44__all__ = [
'FgcmLoadReferenceCatalogConfig',
'FgcmLoadReferenceCatalogTask']
48 """Config for FgcmLoadReferenceCatalogTask"""
50 filterMap = pexConfig.DictField(
51 doc=
"Mapping from physicalFilter label to reference filter name.",
56 applyColorTerms = pexConfig.Field(
57 doc=(
"Apply photometric color terms to reference stars?"
58 "Requires that colorterms be set to a ColorTermLibrary"),
62 colorterms = pexConfig.ConfigField(
63 doc=
"Library of photometric reference catalog name to color term dict.",
64 dtype=ColortermLibrary,
66 referenceSelector = pexConfig.ConfigurableField(
67 target=ReferenceSourceSelectorTask,
68 doc=
"Selection of reference sources",
74 msg =
'Must set filterMap'
75 raise pexConfig.FieldValidationError(FgcmLoadReferenceCatalogConfig.filterMap, self, msg)
77 msg =
"applyColorTerms=True requires the `colorterms` field be set to a ColortermLibrary."
78 raise pexConfig.FieldValidationError(FgcmLoadReferenceCatalogConfig.colorterms, self, msg)
83 Load multi-band reference objects from a reference catalog.
87 refObjLoader : `lsst.meas.algorithms.ReferenceObjectLoader`
88 Reference object loader.
90 Name of reference catalog (
for color term lookups).
92 ConfigClass = FgcmLoadReferenceCatalogConfig
93 _DefaultName = 'fgcmLoadReferenceCatalog'
95 def __init__(self, refObjLoader=None, refCatName=None, **kwargs):
96 """Construct an FgcmLoadReferenceCatalogTask
98 pipeBase.Task.__init__(self, **kwargs)
102 if refObjLoader
is None or refCatName
is None:
103 raise RuntimeError(
"FgcmLoadReferenceCatalogTask requires a refObjLoader and refCatName.")
105 self.makeSubtask(
'referenceSelector')
112 Get a reference catalog that overlaps a healpix pixel, using multiple
113 filters. In addition, apply colorterms if available.
115 Return format
is a numpy recarray
for use
with fgcm,
with the format:
117 dtype = ([(
'ra', `np.float64`),
118 (
'dec', `np.float64`),
119 (
'refMag', `np.float32`, len(filterList)),
120 (
'refMagErr', `np.float32`, len(filterList)])
122 Reference magnitudes (AB) will be 99
for non-detections.
127 Healpix nside of pixel to load
129 Healpix pixel of pixel to load
131 list of `str` of camera filter names.
132 nest: `bool`, optional
133 Is the pixel
in nest format? Default
is False.
137 fgcmRefCat: `np.recarray`
141 lon, lat = hpg.pixel_to_angle(nside, pixel, nest=nest, degrees=
False)
142 center = lsst.geom.SpherePoint(lon * lsst.geom.degrees, lat * lsst.geom.radians)
144 theta_phi = hpg.boundaries(nside, pixel, step=1, nest=nest, lonlat=
False)
146 radius = 0.0 * lsst.geom.radians
147 for ctheta, cphi
in zip(*theta_phi):
148 rad = center.separation(lsst.geom.SpherePoint(cphi * lsst.geom.radians,
149 (np.pi/2. - ctheta) * lsst.geom.radians))
155 center.getDec().asDegrees(),
158 catPix = hpg.angle_to_pixel(nside, fgcmRefCat[
'ra'], fgcmRefCat[
'dec'], nest=nest)
160 inPix, = np.where(catPix == pixel)
162 return fgcmRefCat[inPix]
166 Get a reference catalog that overlaps a circular sky region, using
167 multiple filters. In addition, apply colorterms if available.
169 Return format
is a numpy recarray
for use
with fgcm.
171 dtype = ([(
'ra', `np.float64`),
172 (
'dec', `np.float64`),
173 (
'refMag', `np.float32`, len(filterList)),
174 (
'refMagErr', `np.float32`, len(filterList)])
176 Reference magnitudes (AB) will be 99
for non-detections.
181 ICRS right ascension, degrees.
183 ICRS declination, degrees.
185 Radius to search, degrees.
187 list of `str` of camera filter names.
191 fgcmRefCat: `np.recarray`
194 center = lsst.geom.SpherePoint(ra * lsst.geom.degrees, dec * lsst.geom.degrees)
201 radius * lsst.geom.degrees,
204 if not skyCircle.refCat.isContiguous():
205 refCat = skyCircle.refCat.copy(deep=
True)
207 refCat = skyCircle.refCat
210 goodSources = self.referenceSelector.selectSources(refCat)
211 selected = goodSources.selected
213 fgcmRefCat = np.zeros(np.sum(selected), dtype=[(
'ra',
'f8'),
215 (
'refMag',
'f4', len(filterList)),
216 (
'refMagErr',
'f4', len(filterList))])
217 if fgcmRefCat.size == 0:
227 conv = refCat[0][
'coord_ra'].asDegrees() / float(refCat[0][
'coord_ra'])
228 fgcmRefCat[
'ra'] = refCat[
'coord_ra'][selected] * conv
229 fgcmRefCat[
'dec'] = refCat[
'coord_dec'][selected] * conv
232 fgcmRefCat[
'refMag'][:, :] = 99.0
233 fgcmRefCat[
'refMagErr'][:, :] = 99.0
235 if self.config.applyColorTerms:
237 if fluxField
is None:
240 self.log.debug(
"Applying color terms for filtername=%r" % (filterName))
242 colorterm = self.config.colorterms.getColorterm(filterName, self.
refCatName, doRaise=
True)
244 refMag, refMagErr = colorterm.getCorrectedMagnitudes(refCat)
250 good, = np.where((np.nan_to_num(refMag[selected]) < 90.0)
251 & (np.nan_to_num(refMagErr[selected]) < 90.0)
252 & (np.nan_to_num(refMagErr[selected]) > 0.0))
254 fgcmRefCat[
'refMag'][good, i] = refMag[selected][good]
255 fgcmRefCat[
'refMagErr'][good, i] = refMagErr[selected][good]
263 good, = np.where((np.nan_to_num(refCat[fluxField][selected]) > 0.0)
264 & (np.nan_to_num(refCat[fluxField+
'Err'][selected]) > 0.0))
265 refMag = (refCat[fluxField][selected][good] * units.nJy).to_value(units.ABmag)
266 refMagErr = abMagErrFromFluxErr(refCat[fluxField+
'Err'][selected][good],
267 refCat[fluxField][selected][good])
268 fgcmRefCat[
'refMag'][good, i] = refMag
269 fgcmRefCat[
'refMagErr'][good, i] = refMagErr
273 def _determine_flux_fields(self, center, filterList):
275 Determine the flux field names for a reference catalog.
281 center: `lsst.geom.SpherePoint`
282 The center around which to load test sources.
284 list of `str` of camera filter names.
292 foundReferenceFilter =
False
293 for filterName
in filterList:
294 refFilterName = self.config.filterMap.get(filterName)
295 if refFilterName
is None:
300 0.05 * lsst.geom.degrees,
302 foundReferenceFilter =
True
310 if not foundReferenceFilter:
311 raise RuntimeError(
"Could not find any valid flux field(s) %s" %
312 (
", ".join(filterList)))
316 for filterName
in filterList:
319 refFilterName = self.config.filterMap.get(filterName)
321 if refFilterName
is not None:
323 fluxField = getRefFluxField(results.refCat.schema, filterName=refFilterName)
328 if fluxField
is None:
329 self.log.warning(f
'No reference flux field for camera filter {filterName}')
def __init__(self, refObjLoader=None, refCatName=None, **kwargs)
def getFgcmReferenceStarsSkyCircle(self, ra, dec, radius, filterList)
def _determine_flux_fields(self, center, filterList)
def getFgcmReferenceStarsHealpix(self, nside, pixel, filterList, nest=False)