Coverage for python / lsst / analysis / tools / tasks / astrometricCatalogMatch.py: 31%

1# This file is part of analysis_tools. 

2# 

3# Developed for the LSST Data Management System. 

4# This product includes software developed by the LSST Project 

5# (https://www.lsst.org). 

6# See the COPYRIGHT file at the top-level directory of this distribution 

7# for details of code ownership. 

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 GNU General Public License 

20# along with this program. If not, see <https://www.gnu.org/licenses/>. 

21 

22__all__ = ( 

23 "AstrometricCatalogMatchConfig", 

24 "AstrometricCatalogMatchTask", 

25 "AstrometricCatalogMatchVisitConfig", 

26 "AstrometricCatalogMatchVisitTask", 

27) 

28 

29import lsst.geom 

30import lsst.pex.config as pexConfig 

31import lsst.pipe.base as pipeBase 

32import numpy as np 

33from astropy.table import Table 

34from lsst.pex.config.configurableActions import ConfigurableActionStructField 

35from lsst.pipe.tasks.loadReferenceCatalog import LoadReferenceCatalogTask 

36 

37from ..actions.vector import StarSelector, VisitPlotFlagSelector 

38from ..tasks.catalogMatch import CatalogMatchConfig, CatalogMatchConnections, CatalogMatchTask 

39 

40 

41class AstrometricCatalogMatchConnections( 

42 CatalogMatchConnections, 

43 dimensions=("tract", "skymap"), 

44 defaultTemplates={"targetCatalog": "objectTable_tract", "refCatalog": "gaia_dr3_20230707"}, 

45): 

46 matchedCatalog = pipeBase.connectionTypes.Output( 

47 doc="Catalog with matched target and reference objects with separations", 

48 name="{targetCatalog}_{refCatalog}_match_astrom", 

49 storageClass="ArrowAstropy", 

50 dimensions=("tract", "skymap"), 

51 ) 

52 

53 

54class AstrometricCatalogMatchConfig( 

55 CatalogMatchConfig, pipelineConnections=AstrometricCatalogMatchConnections 

56): 

57 bands = pexConfig.ListField[str]( 

58 doc="The bands to persist downstream", 

59 default=["u", "g", "r", "i", "z", "y"], 

60 ) 

61 

62 def setDefaults(self): 

63 super().setDefaults() 

64 self.matchesRefCat = True 

65 self.referenceCatalogLoader.doApplyColorTerms = False 

66 self.referenceCatalogLoader.refObjLoader.requireProperMotion = True 

67 self.referenceCatalogLoader.refObjLoader.anyFilterMapsToThis = "phot_g_mean" 

68 self.connections.refCatalog = "gaia_dr3_20230707" 

69 

70 

71class AstrometricCatalogMatchTask(CatalogMatchTask): 

72 """Match a tract-level catalog to a reference catalog""" 

73 

74 ConfigClass = AstrometricCatalogMatchConfig 

75 _DefaultName = "analysisToolsAstrometricCatalogMatch" 

76 

77 def runQuantum(self, butlerQC, inputRefs, outputRefs): 

78 # Docs inherited from base class 

79 

80 inputs = butlerQC.get(inputRefs) 

81 columns = self.prepColumns(self.config.bands) 

82 table = inputs["catalog"].get(parameters={"columns": columns}) 

83 

84 tract = butlerQC.quantum.dataId["tract"] 

85 

86 loaderTask = LoadReferenceCatalogTask( 

87 config=self.config.referenceCatalogLoader, 

88 dataIds=[ref.dataId for ref in inputRefs.refCat], 

89 name=inputs["refCat"][0].ref.datasetType.name, 

90 refCats=inputs["refCat"], 

91 ) 

92 

93 skymap = inputs.pop("skymap") 

94 loadedRefCat = self._loadRefCat(loaderTask, skymap[tract]) 

95 outputs = self.run(targetCatalog=table, refCatalog=loadedRefCat, bands=self.config.bands) 

96 

97 butlerQC.put(outputs, outputRefs) 

98 

99 

100class AstrometricCatalogMatchVisitConnections( 

101 pipeBase.PipelineTaskConnections, 

102 dimensions=("visit",), 

103 defaultTemplates={"targetCatalog": "sourceTable_visit", "refCatalog": "gaia_dr3_20230707"}, 

104): 

105 catalog = pipeBase.connectionTypes.Input( 

106 doc="The visit-wide catalog to make plots from.", 

107 storageClass="ArrowAstropy", 

108 name="{targetCatalog}", 

109 dimensions=("visit",), 

110 deferLoad=True, 

111 ) 

112 

113 refCat = pipeBase.connectionTypes.PrerequisiteInput( 

114 doc="The astrometry reference catalog to match to loaded input catalog sources.", 

115 name="{refCatalog}", 

116 storageClass="SimpleCatalog", 

117 dimensions=("skypix",), 

118 deferLoad=True, 

119 multiple=True, 

120 ) 

121 

122 visitSummaryTable = pipeBase.connectionTypes.Input( 

123 doc="A summary table of the ccds in the visit", 

124 storageClass="ExposureCatalog", 

125 name="finalVisitSummary", 

126 dimensions=("visit",), 

127 ) 

128 

129 matchedCatalog = pipeBase.connectionTypes.Output( 

130 doc="Catalog with matched target and reference objects with separations", 

131 name="{targetCatalog}_{refCatalog}_match_astrom", 

132 storageClass="ArrowAstropy", 

133 dimensions=("visit",), 

134 ) 

135 

136 

137class AstrometricCatalogMatchVisitConfig( 

138 AstrometricCatalogMatchConfig, pipelineConnections=AstrometricCatalogMatchVisitConnections 

139): 

140 selectorActions = ConfigurableActionStructField( 

141 doc="Which selectors to use to narrow down the data for QA plotting.", 

142 default={"flagSelector": VisitPlotFlagSelector}, 

143 ) 

144 

145 extraColumns = pexConfig.ListField[str]( 

146 doc="Other catalog columns to persist to downstream tasks", 

147 default=["psfFlux", "psfFluxErr"], 

148 ) 

149 

150 bands = pexConfig.ListField[str]( 

151 doc="The bands to persist downstream", 

152 default=[], 

153 ) 

154 

155 def setDefaults(self): 

156 self.matchesRefCat = True 

157 self.idColumn = "sourceId" 

158 # sourceSelectorActions.sourceSelector is StarSelector 

159 self.sourceSelectorActions.sourceSelector = StarSelector() 

160 self.sourceSelectorActions.sourceSelector.vectorKey = "extendedness" 

161 # extraColumnSelectors.selector1 is SnSelector 

162 self.extraColumnSelectors.selector1.fluxType = "psfFlux" 

163 # extraColumnSelectors.selector2 is GalaxySelector 

164 self.extraColumnSelectors.selector2.vectorKey = "extendedness" 

165 self.extraColumnSelectors.selector3.vectorKey = "extendedness" 

166 self.extraColumnSelectors.selector4 = VisitPlotFlagSelector 

167 self.referenceCatalogLoader.doApplyColorTerms = False 

168 self.referenceCatalogLoader.refObjLoader.requireProperMotion = False 

169 self.referenceCatalogLoader.refObjLoader.anyFilterMapsToThis = "phot_g_mean" 

170 

171 

172class AstrometricCatalogMatchVisitTask(AstrometricCatalogMatchTask): 

173 """Match a visit-level catalog to a reference catalog""" 

174 

175 ConfigClass = AstrometricCatalogMatchVisitConfig 

176 _DefaultName = "analysisToolsAstrometricCatalogMatchVisit" 

177 

178 def runQuantum(self, butlerQC, inputRefs, outputRefs): 

179 # Docs inherited from base class 

180 

181 inputs = butlerQC.get(inputRefs) 

182 

183 columns = [self.config.idColumn, "coord_ra", "coord_dec", "detector"] 

184 columns.extend(self.config.extraColumns) 

185 for selectorAction in [ 

186 self.config.selectorActions, 

187 self.config.sourceSelectorActions, 

188 self.config.extraColumnSelectors, 

189 ]: 

190 for selector in selectorAction: 

191 selectorSchema = selector.getFormattedInputSchema() 

192 columns += [s[0] for s in selectorSchema] 

193 

194 table = inputs["catalog"].get(parameters={"columns": columns}) 

195 

196 loaderTask = LoadReferenceCatalogTask( 

197 config=self.config.referenceCatalogLoader, 

198 dataIds=[ref.dataId for ref in inputRefs.refCat], 

199 name=inputs["refCat"][0].ref.datasetType.name, 

200 refCats=inputs["refCat"], 

201 ) 

202 

203 visitSummaryTable = inputs.pop("visitSummaryTable") 

204 loadedRefCat, badDetectors = self._loadRefCat(loaderTask, visitSummaryTable) 

205 for badDetector in badDetectors: 

206 mask = table["detector"] != badDetector 

207 self.log.info( 

208 "Trimming %d sources from detector %s with incomplete or missing astrometry.", 

209 len(table) - sum(mask), 

210 badDetector, 

211 ) 

212 table = table[mask] 

213 outputs = self.run(targetCatalog=table, refCatalog=loadedRefCat, bands=self.config.bands) 

214 

215 butlerQC.put(outputs, outputRefs) 

216 

217 def _loadRefCat(self, loaderTask, visitSummaryTable): 

218 """Make a reference catalog with coordinates in degrees 

219 

220 Parameters 

221 ---------- 

222 visitSummaryTable : `lsst.afw.table.ExposureCatalog` 

223 The table of visit information 

224 

225 Returns 

226 ------- 

227 refCat : `astropy.table.Table` 

228 Reference catalog. 

229 badDetectors : `list` [ `int` ] 

230 IDs of detectors with invalid visit summaries (should be skipped). 

231 """ 

232 # Get convex hull around the detectors, then get its center and radius 

233 corners = [] 

234 badDetectors = [] 

235 for visSum in visitSummaryTable: 

236 for ra, dec in zip(visSum["raCorners"], visSum["decCorners"]): 

237 # If the coordinates are nan skip this detector. 

238 if not np.isfinite(ra) or not np.isfinite(dec): 

239 badDetectors.append(visSum.getId()) 

240 break 

241 corners.append(lsst.geom.SpherePoint(ra, dec, units=lsst.geom.degrees).getVector()) 

242 if not corners: 

243 raise pipeBase.NoWorkFound("No visit summary corners were finite.") 

244 visitBoundingCircle = lsst.sphgeom.ConvexPolygon.convexHull(corners).getBoundingCircle() 

245 center = lsst.geom.SpherePoint(visitBoundingCircle.getCenter()) 

246 radius = visitBoundingCircle.getOpeningAngle() 

247 

248 # Get the observation date of the visit 

249 epoch = visSum.getVisitInfo().getDate().toAstropy() 

250 

251 # Load the reference catalog in the skyCircle of the detectors, then 

252 # convert the coordinates to degrees and convert the catalog to a 

253 # dataframe 

254 

255 filterName = self.config.referenceCatalogLoader.refObjLoader.anyFilterMapsToThis 

256 loadedRefCat = loaderTask.getSkyCircleCatalog(center, radius, filterName, epoch=epoch) 

257 

258 return Table(loadedRefCat), badDetectors