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 LoadIndexedReferenceObjectsTask, 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 refObjLoader = pexConfig.ConfigurableField(
51 target=LoadIndexedReferenceObjectsTask,
52 doc=
"Reference object loader for photometry",
53 deprecated=
"refObjLoader is deprecated, and will be removed after v24.",
55 filterMap = pexConfig.DictField(
56 doc=
"Mapping from physicalFilter label to reference filter name.",
61 applyColorTerms = pexConfig.Field(
62 doc=(
"Apply photometric color terms to reference stars?"
63 "Requires that colorterms be set to a ColorTermLibrary"),
67 colorterms = pexConfig.ConfigField(
68 doc=
"Library of photometric reference catalog name to color term dict.",
69 dtype=ColortermLibrary,
71 referenceSelector = pexConfig.ConfigurableField(
72 target=ReferenceSourceSelectorTask,
73 doc=
"Selection of reference sources",
79 msg =
'Must set filterMap'
80 raise pexConfig.FieldValidationError(FgcmLoadReferenceCatalogConfig.filterMap, self, msg)
82 msg =
"applyColorTerms=True requires the `colorterms` field be set to a ColortermLibrary."
83 raise pexConfig.FieldValidationError(FgcmLoadReferenceCatalogConfig.colorterms, self, msg)
88 Load multi-band reference objects from a reference catalog.
92 refObjLoader : `lsst.meas.algorithms.ReferenceObjectLoader`
93 Reference object loader.
95 Name of reference catalog (
for color term lookups).
97 ConfigClass = FgcmLoadReferenceCatalogConfig
98 _DefaultName = 'fgcmLoadReferenceCatalog'
100 def __init__(self, refObjLoader=None, refCatName=None, **kwargs):
101 """Construct an FgcmLoadReferenceCatalogTask
103 pipeBase.Task.__init__(self, **kwargs)
107 if refObjLoader
is None or refCatName
is None:
108 raise RuntimeError(
"FgcmLoadReferenceCatalogTask requires a refObjLoader and refCatName.")
110 self.makeSubtask(
'referenceSelector')
117 Get a reference catalog that overlaps a healpix pixel, using multiple
118 filters. In addition, apply colorterms if available.
120 Return format
is a numpy recarray
for use
with fgcm,
with the format:
122 dtype = ([(
'ra', `np.float64`),
123 (
'dec', `np.float64`),
124 (
'refMag', `np.float32`, len(filterList)),
125 (
'refMagErr', `np.float32`, len(filterList)])
127 Reference magnitudes (AB) will be 99
for non-detections.
132 Healpix nside of pixel to load
134 Healpix pixel of pixel to load
136 list of `str` of camera filter names.
137 nest: `bool`, optional
138 Is the pixel
in nest format? Default
is False.
142 fgcmRefCat: `np.recarray`
146 theta, phi = hp.pix2ang(nside, pixel, nest=nest)
147 center = lsst.geom.SpherePoint(phi * lsst.geom.radians, (np.pi/2. - theta) * lsst.geom.radians)
149 corners = hp.boundaries(nside, pixel, step=1, nest=nest)
150 theta_phi = hp.vec2ang(np.transpose(corners))
152 radius = 0.0 * lsst.geom.radians
153 for ctheta, cphi
in zip(*theta_phi):
154 rad = center.separation(lsst.geom.SpherePoint(cphi * lsst.geom.radians,
155 (np.pi/2. - ctheta) * lsst.geom.radians))
161 center.getDec().asDegrees(),
164 catPix = hp.ang2pix(nside, np.radians(90.0 - fgcmRefCat[
'dec']),
165 np.radians(fgcmRefCat[
'ra']), nest=nest)
167 inPix, = np.where(catPix == pixel)
169 return fgcmRefCat[inPix]
173 Get a reference catalog that overlaps a circular sky region, using
174 multiple filters. In addition, apply colorterms if available.
176 Return format
is a numpy recarray
for use
with fgcm.
178 dtype = ([(
'ra', `np.float64`),
179 (
'dec', `np.float64`),
180 (
'refMag', `np.float32`, len(filterList)),
181 (
'refMagErr', `np.float32`, len(filterList)])
183 Reference magnitudes (AB) will be 99
for non-detections.
188 ICRS right ascension, degrees.
190 ICRS declination, degrees.
192 Radius to search, degrees.
194 list of `str` of camera filter names.
198 fgcmRefCat: `np.recarray`
201 center = lsst.geom.SpherePoint(ra * lsst.geom.degrees, dec * lsst.geom.degrees)
208 radius * lsst.geom.degrees,
211 if not skyCircle.refCat.isContiguous():
212 refCat = skyCircle.refCat.copy(deep=
True)
214 refCat = skyCircle.refCat
217 goodSources = self.referenceSelector.selectSources(refCat)
218 selected = goodSources.selected
220 fgcmRefCat = np.zeros(np.sum(selected), dtype=[(
'ra',
'f8'),
222 (
'refMag',
'f4', len(filterList)),
223 (
'refMagErr',
'f4', len(filterList))])
224 if fgcmRefCat.size == 0:
234 conv = refCat[0][
'coord_ra'].asDegrees() / float(refCat[0][
'coord_ra'])
235 fgcmRefCat[
'ra'] = refCat[
'coord_ra'][selected] * conv
236 fgcmRefCat[
'dec'] = refCat[
'coord_dec'][selected] * conv
239 fgcmRefCat[
'refMag'][:, :] = 99.0
240 fgcmRefCat[
'refMagErr'][:, :] = 99.0
242 if self.config.applyColorTerms:
244 if fluxField
is None:
247 self.log.debug(
"Applying color terms for filtername=%r" % (filterName))
249 colorterm = self.config.colorterms.getColorterm(filterName, self.
refCatName, doRaise=
True)
251 refMag, refMagErr = colorterm.getCorrectedMagnitudes(refCat)
257 good, = np.where((np.nan_to_num(refMag[selected]) < 90.0)
258 & (np.nan_to_num(refMagErr[selected]) < 90.0)
259 & (np.nan_to_num(refMagErr[selected]) > 0.0))
261 fgcmRefCat[
'refMag'][good, i] = refMag[selected][good]
262 fgcmRefCat[
'refMagErr'][good, i] = refMagErr[selected][good]
270 good, = np.where((np.nan_to_num(refCat[fluxField][selected]) > 0.0)
271 & (np.nan_to_num(refCat[fluxField+
'Err'][selected]) > 0.0))
272 refMag = (refCat[fluxField][selected][good] * units.nJy).to_value(units.ABmag)
273 refMagErr = abMagErrFromFluxErr(refCat[fluxField+
'Err'][selected][good],
274 refCat[fluxField][selected][good])
275 fgcmRefCat[
'refMag'][good, i] = refMag
276 fgcmRefCat[
'refMagErr'][good, i] = refMagErr
280 def _determine_flux_fields(self, center, filterList):
282 Determine the flux field names for a reference catalog.
288 center: `lsst.geom.SpherePoint`
289 The center around which to load test sources.
291 list of `str` of camera filter names.
299 foundReferenceFilter =
False
300 for filterName
in filterList:
301 refFilterName = self.config.filterMap.get(filterName)
302 if refFilterName
is None:
307 0.05 * lsst.geom.degrees,
309 foundReferenceFilter =
True
317 if not foundReferenceFilter:
318 raise RuntimeError(
"Could not find any valid flux field(s) %s" %
319 (
", ".join(filterList)))
323 for filterName
in filterList:
326 refFilterName = self.config.filterMap.get(filterName)
328 if refFilterName
is not None:
330 fluxField = getRefFluxField(results.refCat.schema, filterName=refFilterName)
335 if fluxField
is None:
336 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)