24__all__ = [
"getRefFluxField",
"getRefFluxKeys",
"LoadReferenceObjectsTask",
"LoadReferenceObjectsConfig",
25 "ReferenceObjectLoader",
"ReferenceObjectLoaderBase"]
30from deprecated.sphinx
import deprecated
40from lsst
import sphgeom
46 """Return True if this name/units combination corresponds to an
47 "old-style" reference catalog flux field.
49 unitsCheck = units != 'nJy'
50 isFlux = name.endswith(
'_flux')
51 isFluxSigma = name.endswith(
'_fluxSigma')
52 isFluxErr = name.endswith(
'_fluxErr')
53 return (isFlux
or isFluxSigma
or isFluxErr)
and unitsCheck
58 """Return True if the units of all flux and fluxErr are correct (nJy).
67 """"Return the format version stored in a reference catalog header.
72 Reference catalog to inspect.
77 Format verison integer. Returns `0` if the catalog has no metadata
78 or the metadata does
not include a
"REFCAT_FORMAT_VERSION" key.
82 deprecation_msg =
"Support for version 0 refcats (pre-nJy fluxes) will be removed after v25."
83 md = refCat.getMetadata()
85 warnings.warn(deprecation_msg)
88 return md.getScalar(
"REFCAT_FORMAT_VERSION")
90 warnings.warn(deprecation_msg)
95@deprecated(reason="Support for version 0 refcats (pre-nJy fluxes) will be removed after v25.
",
96 version="v24.0", category=FutureWarning)
98 """Convert fluxes in a catalog from jansky to nanojansky.
103 The catalog to convert.
105 Log to send messages to.
106 doConvert : `bool`, optional
107 Return a converted catalog,
or just identify the fields that need to be converted?
108 This supports the
"write=False" mode of `bin/convert_to_nJy.py`.
113 The converted catalog,
or None if ``doConvert``
is False.
117 Support
for old units
in reference catalogs will be removed after the
118 release of late calendar year 2019.
119 Use `meas_algorithms/bin/convert_to_nJy.py` to update your reference catalog.
123 mapper = afwTable.SchemaMapper(catalog.schema, shareAliasMap=
False)
124 mapper.addMinimalSchema(afwTable.SimpleTable.makeMinimalSchema())
127 for field
in catalog.schema:
128 oldName = field.field.getName()
129 oldUnits = field.field.getUnits()
133 if oldName.endswith(
'_fluxSigma'):
134 name = oldName.replace(
'_fluxSigma',
'_fluxErr')
137 newField = afwTable.Field[field.dtype](name, field.field.getDoc(), units)
138 mapper.addMapping(field.getKey(), newField)
139 input_fields.append(field.field)
140 output_fields.append(newField)
142 mapper.addMapping(field.getKey())
144 fluxFieldsStr =
'; '.join(
"(%s, '%s')" % (field.getName(), field.getUnits())
for field
in input_fields)
147 newSchema = mapper.getOutputSchema()
148 output = afwTable.SimpleCatalog(newSchema)
149 output.reserve(len(catalog))
150 output.extend(catalog, mapper=mapper)
151 for field
in output_fields:
152 output[field.getName()] *= 1e9
153 log.info(
"Converted refcat flux fields to nJy (name, units): %s", fluxFieldsStr)
156 log.info(
"Found old-style refcat flux fields (name, units): %s", fluxFieldsStr)
161 """This is a private helper class which filters catalogs by
162 row based on the row being inside the region used to initialize
168 The spatial region which all objects should lie within
174 """This call method on an instance of this class takes in a reference
175 catalog, and the region
from which the catalog was generated.
177 If the catalog region
is entirely contained within the region used to
178 initialize this
class, then all the entries
in the catalog must be
179 within the region
and so the whole catalog
is returned.
181 If the catalog region
is not entirely contained, then the location
for
182 each record
is tested against the region used to initialize the
class.
183 Records which fall inside this region are added to a new catalog,
and
184 this catalog
is then returned.
189 SourceCatalog to be filtered.
191 Region
in which the catalog was created
193 if catRegion.isWithin(self.
region):
197 filteredRefCat = type(refCat)(refCat.table)
198 for record
in refCat:
199 if self.
region.contains(record.getCoord().getVector()):
200 filteredRefCat.append(record)
201 return filteredRefCat
205 pixelMargin = pexConfig.RangeField(
206 doc=
"Padding to add to 4 all edges of the bounding box (pixels)",
211 anyFilterMapsToThis = pexConfig.Field(
212 doc=(
"Always use this reference catalog filter, no matter whether or what filter name is "
213 "supplied to the loader. Effectively a trivial filterMap: map all filter names to this filter."
214 " This can be set for purely-astrometric catalogs (e.g. Gaia DR2) where there is only one "
215 "reasonable choice for every camera filter->refcat mapping, but not for refcats used for "
216 "photometry, which need a filterMap and/or colorterms/transmission corrections."),
221 filterMap = pexConfig.DictField(
222 doc=(
"Mapping of camera filter name: reference catalog filter name; "
223 "each reference filter must exist in the refcat."
224 " Note that this does not perform any bandpass corrections: it is just a lookup."),
229 requireProperMotion = pexConfig.Field(
230 doc=
"Require that the fields needed to correct proper motion "
231 "(epoch, pm_ra and pm_dec) are present?",
235 ref_dataset_name = pexConfig.Field(
236 doc=
"Deprecated; do not use. Added for easier transition from LoadIndexedReferenceObjectsConfig to "
237 "LoadReferenceObjectsConfig",
240 deprecated=
'This field is not used. It will be removed after v25.',
246 msg =
"`filterMap` and `anyFilterMapsToThis` are mutually exclusive"
247 raise pexConfig.FieldValidationError(LoadReferenceObjectsConfig.anyFilterMapsToThis,
252 """This class facilitates loading reference catalogs.
254 The QuantumGraph generation will create a list of datasets that may
255 possibly overlap a given region. These datasets are then used to construct
256 an instance of this class. The
class instance should then be passed into
257 a task which needs reference catalogs. These tasks should then determine
258 the exact region of the sky reference catalogs will be loaded
for,
and
259 call a corresponding method to load the reference objects.
263 dataIds : iterable of `lsst.daf.butler.DataCoordinate`
264 An iterable object of data IDs that point to reference catalogs.
265 refCats : iterable of `lsst.daf.butler.DeferredDatasetHandle`
266 Handles to load refCats on demand.
267 name : `str`, optional
268 The name of the refcat that this object will load. This name
is used
269 for applying colorterms,
for example.
270 config : `LoadReferenceObjectsConfig`
271 Configuration of this reference loader.
272 log : `
lsst.log.Log`, `logging.Logger`
or `
None`, optional
273 Logger object used to write out messages. If `
None` a default
276 ConfigClass = LoadReferenceObjectsConfig
278 def __init__(self, dataIds, refCats, name=None, log=None, config=None, **kwargs):
280 warnings.warn(
"Instantiating ReferenceObjectLoader with additional kwargs is deprecated "
281 "and will be removed after v25.0", FutureWarning, stacklevel=2)
289 self.
log = log
or logging.getLogger(__name__).getChild(
"ReferenceObjectLoader")
292 """Apply proper motion correction to a reference catalog.
294 Adjust position and position error
in the ``catalog``
295 for proper motion to the specified ``epoch``,
296 modifying the catalog
in place.
301 Catalog of positions, containing at least these fields:
303 - Coordinates, retrieved by the table
's coordinate key.
304 - ``coord_raErr`` : Error in Right Ascension (rad).
305 - ``coord_decErr`` : Error
in Declination (rad).
306 - ``pm_ra`` : Proper motion
in Right Ascension (rad/yr,
308 - ``pm_raErr`` : Error
in ``pm_ra`` (rad/yr), optional.
309 - ``pm_dec`` : Proper motion
in Declination (rad/yr,
311 - ``pm_decErr`` : Error
in ``pm_dec`` (rad/yr), optional.
312 - ``epoch`` : Mean epoch of object (an astropy.time.Time)
313 epoch : `astropy.time.Time`
314 Epoch to which to correct proper motion.
315 If
None, do
not apply PM corrections
or raise if
316 ``config.requireProperMotion``
is True.
321 Raised
if ``config.requireProperMotion``
is set but we cannot
322 apply the proper motion correction
for some reason.
325 if self.
config.requireProperMotion:
326 raise RuntimeError(
"requireProperMotion=True but epoch not provided to loader.")
328 self.
log.debug(
"No epoch provided: not applying proper motion corrections to refcat.")
332 if (
"pm_ra" in catalog.schema
333 and not isinstance(catalog.schema[
"pm_ra"].asKey(), afwTable.KeyAngle)):
334 if self.
config.requireProperMotion:
335 raise RuntimeError(
"requireProperMotion=True but refcat pm_ra field is not an Angle.")
337 self.
log.warning(
"Reference catalog pm_ra field is not an Angle; cannot apply proper motion.")
340 if (
"epoch" not in catalog.schema
or "pm_ra" not in catalog.schema):
341 if self.
config.requireProperMotion:
342 raise RuntimeError(
"requireProperMotion=True but PM data not available from catalog.")
344 self.
log.warning(
"Proper motion correction not available for this reference catalog.")
351 filterMap=None, centroids=False):
352 """This function takes in a reference catalog and returns a new catalog
353 with additional columns defined
from the remaining function arguments.
358 Reference catalog to map to new catalog
359 anyFilterMapsToThis : `str`, optional
360 Always use this reference catalog filter.
361 Mutually exclusive
with `filterMap`
362 filterMap : `dict` [`str`,`str`], optional
363 Mapping of camera filter name: reference catalog filter name.
364 centroids : `bool`, optional
365 Add centroid fields to the loaded Schema. ``loadPixelBox`` expects
366 these fields to exist.
371 Deep copy of input reference catalog
with additional columns added
373 if anyFilterMapsToThis
or filterMap:
374 ReferenceObjectLoader._addFluxAliases(refCat.schema, anyFilterMapsToThis, filterMap)
376 mapper = afwTable.SchemaMapper(refCat.schema,
True)
377 mapper.addMinimalSchema(refCat.schema,
True)
378 mapper.editOutputSchema().disconnectAliases()
385 mapper.editOutputSchema().addField(
"centroid_x", type=float, doReplace=
True)
386 mapper.editOutputSchema().addField(
"centroid_y", type=float, doReplace=
True)
387 mapper.editOutputSchema().addField(
"hasCentroid", type=
"Flag", doReplace=
True)
388 mapper.editOutputSchema().getAliasMap().set(
"slot_Centroid",
"centroid")
390 expandedCat = afwTable.SimpleCatalog(mapper.getOutputSchema())
391 expandedCat.setMetadata(refCat.getMetadata())
392 expandedCat.extend(refCat, mapper=mapper)
398 """Add aliases for camera filter fluxes to the schema.
400 For each camFilter: refFilter in filterMap, adds these aliases:
401 <camFilter>_camFlux: <refFilter>_flux
402 <camFilter>_camFluxErr: <refFilter>_fluxErr,
if the latter exists
403 or sets `anyFilterMapsToThis`
in the schema.
408 Schema
for reference catalog.
409 anyFilterMapsToThis : `str`, optional
410 Always use this reference catalog filter.
411 Mutually exclusive
with `filterMap`.
412 filterMap : `dict` [`str`,`str`], optional
413 Mapping of camera filter name: reference catalog filter name.
414 Mutually exclusive
with `anyFilterMapsToThis`.
419 Raised
if any required reference flux field
is missing
from the
423 if anyFilterMapsToThis
and filterMap:
424 raise ValueError(
"anyFilterMapsToThis and filterMap are mutually exclusive!")
426 aliasMap = schema.getAliasMap()
428 if anyFilterMapsToThis
is not None:
429 refFluxName = anyFilterMapsToThis +
"_flux"
430 if refFluxName
not in schema:
431 msg = f
"Unknown reference filter for anyFilterMapsToThis='{refFluxName}'"
432 raise RuntimeError(msg)
433 aliasMap.set(
"anyFilterMapsToThis", refFluxName)
436 def addAliasesForOneFilter(filterName, refFilterName):
437 """Add aliases for a single filter
441 filterName : `str` (optional)
442 Camera filter name. The resulting alias name is
444 refFilterName : `str`
445 Reference catalog filter name; the field
446 <refFilterName>_flux must exist.
448 camFluxName = filterName + "_camFlux"
449 refFluxName = refFilterName +
"_flux"
450 if refFluxName
not in schema:
451 raise RuntimeError(
"Unknown reference filter %s" % (refFluxName,))
452 aliasMap.set(camFluxName, refFluxName)
453 refFluxErrName = refFluxName +
"Err"
454 if refFluxErrName
in schema:
455 camFluxErrName = camFluxName +
"Err"
456 aliasMap.set(camFluxErrName, refFluxErrName)
458 if filterMap
is not None:
459 for filterName, refFilterName
in filterMap.items():
460 addAliasesForOneFilter(filterName, refFilterName)
476 outerLocalBBox.grow(BBoxPadding)
477 innerLocalBBox.grow(-1*BBoxPadding)
489 innerBoxCorners = innerLocalBBox.getCorners()
490 innerSphCorners = [wcs.pixelToSky(corner).getVector()
for corner
in innerBoxCorners]
493 outerBoxCorners = outerLocalBBox.getCorners()
494 outerSphCorners = [wcs.pixelToSky(corner).getVector()
for corner
in outerBoxCorners]
497 return innerSkyRegion, outerSkyRegion, innerSphCorners, outerSphCorners
501 """Compute on-sky center and radius of search region.
508 WCS; used to convert pixel positions to sky coordinates.
510 Padding to add to 4 all edges of the bounding box (pixels).
514 results : `lsst.pipe.base.Struct`
518 ICRS center of the search region.
520 Radius of the search region.
522 Bounding box used to compute the circle.
525 bbox.grow(pixelMargin)
526 coord = wcs.pixelToSky(bbox.getCenter())
527 radius = max(coord.separation(wcs.pixelToSky(pp))
for pp
in bbox.getCorners())
528 return pipeBase.Struct(coord=coord, radius=radius, bbox=bbox)
532 """Return metadata about the loaded reference catalog, in an on-sky
535 This metadata is used
for reloading the catalog (e.g.
for
536 reconstituting a normalized match list).
541 ICRS center of the search region.
543 Radius of the search region.
545 Name of the camera filter.
546 epoch : `astropy.time.Time`
or `
None`, optional
547 Epoch that proper motion
and parallax were corrected to,
or `
None`
548 if no such corrections were applied.
553 Metadata about the catalog.
556 md.add('RA', coord.getRa().asDegrees(),
'field center in degrees')
557 md.add(
'DEC', coord.getDec().asDegrees(),
'field center in degrees')
558 md.add(
'RADIUS', radius.asDegrees(),
'field radius in degrees, minimum')
561 md.add(
'SMATCHV', 2,
'SourceMatchVector version number')
562 md.add(
'FILTER', filterName,
'camera filter name for photometric data')
563 md.add(
'TIMESYS',
"TAI",
"time scale of time keywords")
564 md.add(
'JEPOCH',
None if epoch
is None else epoch.tai.jyear,
565 'Julian epoch (TAI Julian Epoch year) for catalog')
569 bboxToSpherePadding=100):
570 """Return metadata about the loaded reference catalog, in an
573 This metadata is used
for reloading the catalog (e.g.,
for
574 reconstituting a normalised match list).
579 Bounding box
for the pixels.
581 The WCS object associated
with ``bbox``.
583 Name of the camera filter.
584 epoch : `astropy.time.Time`
or `
None`, optional
585 Epoch that proper motion
and parallax were corrected to,
or `
None`
586 if no such corrections were applied.
587 bboxToSpherePadding : `int`, optional
588 Padding
in pixels to account
for translating a set of corners into
589 a spherical (convex) boundary that
is certain to encompass the
590 enitre area covered by the box.
595 The metadata detailing the search parameters used
for this
599 md = self.getMetadataCircle(circle.coord, circle.radius, filterName, epoch=epoch)
601 paddedBbox = circle.bbox
602 _, _, innerCorners, outerCorners = self._makeBoxRegion(paddedBbox, wcs, bboxToSpherePadding)
603 for box, corners
in zip((
"INNER",
"OUTER"), (innerCorners, outerCorners)):
604 for (name, corner)
in zip((
"UPPER_LEFT",
"UPPER_RIGHT",
"LOWER_LEFT",
"LOWER_RIGHT"),
606 md.add(f
"{box}_{name}_RA",
geom.SpherePoint(corner).getRa().asDegrees(), f
"{box}_corner")
607 md.add(f
"{box}_{name}_DEC",
geom.SpherePoint(corner).getDec().asDegrees(), f
"{box}_corner")
611 bboxToSpherePadding=100):
612 """Load reference objects that are within a pixel-based rectangular
615 This algorithm works by creating a spherical box whose corners
616 correspond to the WCS converted corners of the input bounding box
617 (possibly padded). It then defines a filtering function which looks at
618 the pixel position of the reference objects and accepts only those that
619 lie within the specified bounding box.
621 The spherical box region
and filtering function are passed to the
622 generic loadRegion method which loads
and filters the reference objects
623 from the datastore
and returns a single catalog containing the filtered
624 set of reference objects.
629 Box which bounds a region
in pixel space.
631 Wcs object defining the pixel to sky (
and inverse) transform
for
632 the supplied ``bbox``.
634 Name of camera filter.
635 epoch : `astropy.time.Time`
or `
None`, optional
636 Epoch to which to correct proper motion
and parallax,
or `
None`
637 to
not apply such corrections.
638 bboxToSpherePadding : `int`, optional
639 Padding to account
for translating a set of corners into a
640 spherical (convex) boundary that
is certain to encompase the
641 enitre area covered by the box.
645 output : `lsst.pipe.base.Struct`
646 Results struct
with attributes:
649 Catalog containing reference objects inside the specified
650 bounding box (padded by self.
config.pixelMargin).
652 Name of the field containing the flux associated
with
658 Raised
if no reference catalogs could be found
for the specified
661 Raised
if the loaded reference catalogs do
not have matching
665 paddedBbox.grow(self.config.pixelMargin)
666 innerSkyRegion, outerSkyRegion, _, _ = self._makeBoxRegion(paddedBbox, wcs, bboxToSpherePadding)
668 def _filterFunction(refCat, region):
679 refCat = preFiltFunc(refCat, region)
685 afwTable.updateRefCentroids(wcs, refCat)
688 if innerSkyRegion.contains(region):
692 filteredRefCat = type(refCat)(refCat.table)
693 centroidKey = afwTable.Point2DKey(refCat.schema[
'centroid'])
694 for record
in refCat:
695 pixCoords = record[centroidKey]
697 filteredRefCat.append(record)
698 return filteredRefCat
699 return self.
loadRegion(outerSkyRegion, filterName, filtFunc=_filterFunction, epoch=epoch)
701 def loadRegion(self, region, filterName, filtFunc=None, epoch=None):
702 """Load reference objects within a specified region.
704 This function loads the DataIds used to construct an instance of this
705 class which intersect or are contained within the specified region. The
706 reference catalogs which intersect but are
not fully contained within
707 the input region are further filtered by the specified filter function.
708 This function returns a single source catalog containing all reference
709 objects inside the specified region.
715 should define the spatial region
for which reference objects are to
717 filtFunc : callable
or `
None`, optional
718 This optional parameter should be a callable object that takes a
719 reference catalog
and its corresponding region
as parameters,
720 filters the catalog by some criteria
and returns the filtered
721 reference catalog. If `
None`, an internal filter function
is used
722 which filters according to
if a reference object falls within the
725 Name of camera filter.
726 epoch : `astropy.time.Time`
or `
None`, optional
727 Epoch to which to correct proper motion
and parallax,
or `
None` to
728 not apply such corrections.
732 output : `lsst.pipe.base.Struct`
733 Results struct
with attributes:
736 Catalog containing reference objects which intersect the
737 input region, filtered by the specified filter function.
739 Name of the field containing the flux associated
with
745 Raised
if no reference catalogs could be found
for the specified
748 Raised
if the loaded reference catalogs do
not have matching
751 regionLat = region.getBoundingBox().getLat()
752 regionLon = region.getBoundingBox().getLon()
753 self.log.info("Loading reference objects from %s in region bounded by "
754 "[%.8f, %.8f], [%.8f, %.8f] RA Dec",
756 regionLon.getA().asDegrees(), regionLon.getB().asDegrees(),
757 regionLat.getA().asDegrees(), regionLat.getB().asDegrees())
766 intersects = dataId.region.intersects(region)
768 intersects = region.intersects(dataId.region)
771 overlapList.append((dataId, refCat))
773 if len(overlapList) == 0:
774 raise RuntimeError(
"No reference tables could be found for input region")
776 firstCat = overlapList[0][1].get()
777 refCat = filtFunc(firstCat, overlapList[0][0].region)
778 trimmedAmount = len(firstCat) - len(refCat)
781 for dataId, inputRefCat
in overlapList[1:]:
782 tmpCat = inputRefCat.get()
784 if tmpCat.schema != firstCat.schema:
785 raise TypeError(
"Reference catalogs have mismatching schemas")
787 filteredCat = filtFunc(tmpCat, dataId.region)
788 refCat.extend(filteredCat)
789 trimmedAmount += len(tmpCat) - len(filteredCat)
791 self.
log.debug(
"Trimmed %d refCat objects lying outside padded region, leaving %d",
792 trimmedAmount, len(refCat))
793 self.
log.info(
"Loaded %d reference objects", len(refCat))
796 if not refCat.isContiguous():
797 refCat = refCat.copy(deep=
True)
805 self.
log.warning(
"Found version 0 reference catalog with old style units in schema.")
806 self.
log.warning(
"run `meas_algorithms/bin/convert_refcat_to_nJy.py` to convert fluxes to nJy.")
807 self.
log.warning(
"See RFC-575 for more details.")
811 anyFilterMapsToThis=self.
config.anyFilterMapsToThis,
812 filterMap=self.
config.filterMap)
815 if not expandedCat.isContiguous():
816 expandedCat = expandedCat.copy(deep=
True)
819 return pipeBase.Struct(refCat=expandedCat, fluxField=fluxField)
822 """Load reference objects that lie within a circular region on the sky.
824 This method constructs a circular region from an input center
and
825 angular radius, loads reference catalogs which are contained
in or
826 intersect the circle,
and filters reference catalogs which intersect
827 down to objects which lie within the defined circle.
832 Point defining the center of the circular region.
834 Defines the angular radius of the circular region.
836 Name of camera filter.
837 epoch : `astropy.time.Time`
or `
None`, optional
838 Epoch to which to correct proper motion
and parallax,
or `
None` to
839 not apply such corrections.
843 output : `lsst.pipe.base.Struct`
844 Results struct
with attributes:
847 Catalog containing reference objects inside the specified
850 Name of the field containing the flux associated
with
853 centerVector = ctrCoord.getVector()
856 return self.
loadRegion(circularRegion, filterName, epoch=epoch)
860 """Get the name of a flux field from a schema.
865 Reference catalog schema.
867 Name of camera filter.
871 fluxFieldName : `str`
876 Return the alias of ``anyFilterMapsToThis``, if present
877 else,
return ``*filterName*_camFlux``
if present,
878 else,
return ``*filterName*_flux``
if present (camera filter name
879 matches reference filter name),
else raise an exception.
884 Raised
if an appropriate field
is not found.
886 if not isinstance(schema, afwTable.Schema):
887 raise RuntimeError(
"schema=%s is not a schema" % (schema,))
889 return schema.getAliasMap().get(
"anyFilterMapsToThis")
893 fluxFieldList = [filterName +
"_camFlux", filterName +
"_flux"]
894 for fluxField
in fluxFieldList:
895 if fluxField
in schema:
898 raise RuntimeError(
"Could not find flux field(s) %s" % (
", ".join(fluxFieldList)))
902 """Return keys for flux and flux error.
907 Reference catalog schema.
909 Name of camera filter.
917 - flux error key, if present,
else None
922 If flux field
not found.
925 fluxErrField = fluxField + "Err"
926 fluxKey = schema[fluxField].asKey()
928 fluxErrKey = schema[fluxErrField].asKey()
931 return (fluxKey, fluxErrKey)
934@deprecated(reason=(
"This task is used in gen2 only; it will be removed after v25. "
935 "See DM-35671 for details on updating code to avoid this warning."),
936 version=
"v25.0", category=FutureWarning)
938 """Abstract gen2 base class to load objects from reference catalogs.
940 _DefaultName = "LoadReferenceObjects"
943 def loadSkyCircle(self, ctrCoord, radius, filterName, epoch=None, centroids=False):
944 """Load reference objects that overlap a circular sky region.
949 ICRS center of search region.
951 Radius of search region.
953 Name of filter. This can be used for flux limit comparisons.
954 epoch : `astropy.time.Time`
or `
None`, optional
955 Epoch to which to correct proper motion
and parallax,
or `
None` to
956 not apply such corrections.
957 centroids : `bool`, optional
958 Add centroid fields to the loaded Schema. ``loadPixelBox`` expects
959 these fields to exist.
963 results : `lsst.pipe.base.Struct`
964 A `~lsst.pipe.base.Struct` containing the following fields:
967 A catalog of reference objects
with the standard
968 schema,
as documented
in the main doc string
for
969 `LoadReferenceObjects`.
970 The catalog
is guaranteed to be contiguous.
971 (`lsst.afw.catalog.SimpleCatalog`)
973 Name of flux field
for specified `filterName`. (`str`)
977 Note that subclasses are responsible
for performing the proper motion
978 correction, since this
is the lowest-level interface
for retrieving
984@deprecated(reason=
"Base class only used for gen2 interface, and will be removed after v25.0. "
985 "Please use ReferenceObjectLoader directly.",
986 version=
"v25.0", category=FutureWarning)
988 """Stub of a deprecated class.
992 config : `lsst.pex.config.Config`
993 Configuration for the loader.
1000 """Apply proper motion correction to a reference catalog.
1002 Adjust position and position error
in the ``catalog``
1003 for proper motion to the specified ``epoch``,
1004 modifying the catalog
in place.
1009 Log object to write to.
1011 Catalog of positions, containing:
1013 - Coordinates, retrieved by the table
's coordinate key.
1014 - ``coord_raErr`` : Error in Right Ascension (rad).
1015 - ``coord_decErr`` : Error
in Declination (rad).
1016 - ``pm_ra`` : Proper motion
in Right Ascension (rad/yr,
1018 - ``pm_raErr`` : Error
in ``pm_ra`` (rad/yr), optional.
1019 - ``pm_dec`` : Proper motion
in Declination (rad/yr,
1021 - ``pm_decErr`` : Error
in ``pm_dec`` (rad/yr), optional.
1022 - ``epoch`` : Mean epoch of object (an astropy.time.Time)
1023 epoch : `astropy.time.Time`
1024 Epoch to which to correct proper motion.
1026 if "epoch" not in catalog.schema
or "pm_ra" not in catalog.schema
or "pm_dec" not in catalog.schema:
1027 log.warning(
"Proper motion correction not available from catalog")
1029 if not catalog.isContiguous():
1030 raise RuntimeError(
"Catalog must be contiguous")
1031 catEpoch = astropy.time.Time(catalog[
"epoch"], scale=
"tai", format=
"mjd")
1032 log.info(
"Correcting reference catalog for proper motion to %r", epoch)
1034 timeDiffsYears = (epoch.tai - catEpoch).to(astropy.units.yr).value
1035 coordKey = catalog.table.getCoordKey()
1038 pmRaRad = catalog[
"pm_ra"]
1039 pmDecRad = catalog[
"pm_dec"]
1040 offsetsRaRad = pmRaRad*timeDiffsYears
1041 offsetsDecRad = pmDecRad*timeDiffsYears
1049 offsetBearingsRad = numpy.arctan2(pmDecRad*1e6, pmRaRad*1e6)
1050 offsetAmountsRad = numpy.hypot(offsetsRaRad, offsetsDecRad)
1051 for record, bearingRad, amountRad
in zip(catalog, offsetBearingsRad, offsetAmountsRad):
1052 record.set(coordKey,
1053 record.get(coordKey).offset(bearing=bearingRad*geom.radians,
1054 amount=amountRad*geom.radians))
1057 if "coord_raErr" in catalog.schema:
1058 catalog[
"coord_raErr"] = numpy.hypot(catalog[
"coord_raErr"],
1059 catalog[
"pm_raErr"]*timeDiffsYears)
1060 if "coord_decErr" in catalog.schema:
1061 catalog[
"coord_decErr"] = numpy.hypot(catalog[
"coord_decErr"],
1062 catalog[
"pm_decErr"]*timeDiffsYears)
def __call__(self, refCat, catRegion)
def __init__(self, region)
pexConfig anyFilterMapsToThis
def loadSkyCircle(self, ctrCoord, radius, filterName, epoch=None, centroids=False)
def __init__(self, config=None, *args, **kwargs)
def _calculateCircle(bbox, wcs, pixelMargin)
LoadReferenceObjectsConfig ConfigClass
def __init__(self, dataIds, refCats, name=None, log=None, config=None, **kwargs)
def _makeBoxRegion(BBox, wcs, BBoxPadding)
def applyProperMotions(self, catalog, epoch)
def loadSkyCircle(self, ctrCoord, radius, filterName, epoch=None)
def _remapReferenceCatalogSchema(refCat, *anyFilterMapsToThis=None, filterMap=None, centroids=False)
def loadRegion(self, region, filterName, filtFunc=None, epoch=None)
def getMetadataBox(self, bbox, wcs, filterName, epoch=None, bboxToSpherePadding=100)
def _addFluxAliases(schema, anyFilterMapsToThis=None, filterMap=None)
def loadPixelBox(self, bbox, wcs, filterName, epoch=None, bboxToSpherePadding=100)
def getMetadataCircle(coord, radius, filterName, epoch=None)
def hasNanojanskyFluxUnits(schema)
def getRefFluxKeys(schema, filterName)
def convertToNanojansky(catalog, log, doConvert=True)
def isOldFluxField(name, units)
def getRefFluxField(schema, filterName)
def getFormatVersionFromRefCat(refCat)
def applyProperMotionsImpl(log, catalog, epoch)