Coverage for python/lsst/pipe/tasks/computeExposureSummaryStats.py: 16%

219 statements  

« prev     ^ index     » next       coverage.py v7.4.0, created at 2024-01-04 13:31 +0000

1# This file is part of pipe_tasks. 

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__ = ["ComputeExposureSummaryStatsTask", "ComputeExposureSummaryStatsConfig"] 

23 

24import warnings 

25import numpy as np 

26from scipy.stats import median_abs_deviation as sigmaMad 

27import astropy.units as units 

28from astropy.time import Time 

29from astropy.coordinates import AltAz, SkyCoord, EarthLocation 

30from lsst.daf.base import DateTime 

31 

32import lsst.pipe.base as pipeBase 

33import lsst.pex.config as pexConfig 

34import lsst.afw.math as afwMath 

35import lsst.afw.image as afwImage 

36import lsst.geom as geom 

37from lsst.meas.algorithms import ScienceSourceSelectorTask 

38from lsst.utils.timer import timeMethod 

39 

40 

41class ComputeExposureSummaryStatsConfig(pexConfig.Config): 

42 """Config for ComputeExposureSummaryTask""" 

43 sigmaClip = pexConfig.Field( 

44 dtype=float, 

45 doc="Sigma for outlier rejection for sky noise.", 

46 default=3.0, 

47 ) 

48 clipIter = pexConfig.Field( 

49 dtype=int, 

50 doc="Number of iterations of outlier rejection for sky noise.", 

51 default=2, 

52 ) 

53 badMaskPlanes = pexConfig.ListField( 

54 dtype=str, 

55 doc="Mask planes that, if set, the associated pixel should not be included sky noise calculation.", 

56 default=("NO_DATA", "SUSPECT"), 

57 ) 

58 starSelection = pexConfig.Field( 

59 doc="Field to select full list of sources used for PSF modeling.", 

60 dtype=str, 

61 default="calib_psf_used", 

62 ) 

63 starSelector = pexConfig.ConfigurableField( 

64 target=ScienceSourceSelectorTask, 

65 doc="Selection of sources to compute PSF star statistics.", 

66 ) 

67 starShape = pexConfig.Field( 

68 doc="Base name of columns to use for the source shape in the PSF statistics computation.", 

69 dtype=str, 

70 default="slot_Shape" 

71 ) 

72 psfShape = pexConfig.Field( 

73 doc="Base name of columns to use for the PSF shape in the PSF statistics computation.", 

74 dtype=str, 

75 default="slot_PsfShape" 

76 ) 

77 psfSampling = pexConfig.Field( 

78 dtype=int, 

79 doc="Sampling rate in pixels in each dimension for the maxDistToNearestPsf metric " 

80 "caclulation grid (the tradeoff is between adequate sampling versus speed).", 

81 default=8, 

82 ) 

83 psfGridSampling = pexConfig.Field( 

84 dtype=int, 

85 doc="Sampling rate in pixels in each dimension for PSF model robustness metric " 

86 "caclulations grid (the tradeoff is between adequate sampling versus speed).", 

87 default=96, 

88 ) 

89 psfBadMaskPlanes = pexConfig.ListField( 

90 dtype=str, 

91 doc="Mask planes that, if set, the associated pixel should not be included in the PSF model " 

92 "robutsness metric calculations (namely, maxDistToNearestPsf and psfTraceRadiusDelta).", 

93 default=("BAD", "CR", "EDGE", "INTRP", "NO_DATA", "SAT", "SUSPECT"), 

94 ) 

95 

96 def setDefaults(self): 

97 super().setDefaults() 

98 

99 self.starSelector.setDefaults() 

100 self.starSelector.doFlags = True 

101 self.starSelector.doSignalToNoise = True 

102 self.starSelector.doUnresolved = False 

103 self.starSelector.doIsolated = False 

104 self.starSelector.doRequireFiniteRaDec = False 

105 self.starSelector.doRequirePrimary = False 

106 

107 self.starSelector.signalToNoise.minimum = 50.0 

108 self.starSelector.signalToNoise.maximum = 1000.0 

109 

110 self.starSelector.flags.bad = ["slot_Shape_flag", "slot_PsfFlux_flag"] 

111 # Select stars used for PSF modeling. 

112 self.starSelector.flags.good = ["calib_psf_used"] 

113 

114 self.starSelector.signalToNoise.fluxField = "slot_PsfFlux_instFlux" 

115 self.starSelector.signalToNoise.errField = "slot_PsfFlux_instFluxErr" 

116 

117 

118class ComputeExposureSummaryStatsTask(pipeBase.Task): 

119 """Task to compute exposure summary statistics. 

120 

121 This task computes various quantities suitable for DPDD and other 

122 downstream processing at the detector centers, including: 

123 - psfSigma 

124 - psfArea 

125 - psfIxx 

126 - psfIyy 

127 - psfIxy 

128 - ra 

129 - dec 

130 - zenithDistance 

131 - zeroPoint 

132 - skyBg 

133 - skyNoise 

134 - meanVar 

135 - raCorners 

136 - decCorners 

137 - astromOffsetMean 

138 - astromOffsetStd 

139 

140 These additional quantities are computed from the stars in the detector: 

141 - psfStarDeltaE1Median 

142 - psfStarDeltaE2Median 

143 - psfStarDeltaE1Scatter 

144 - psfStarDeltaE2Scatter 

145 - psfStarDeltaSizeMedian 

146 - psfStarDeltaSizeScatter 

147 - psfStarScaledDeltaSizeScatter 

148 

149 These quantities are computed based on the PSF model and image mask 

150 to assess the robustness of the PSF model across a given detector 

151 (against, e.g., extrapolation instability): 

152 - maxDistToNearestPsf 

153 - psfTraceRadiusDelta 

154 """ 

155 ConfigClass = ComputeExposureSummaryStatsConfig 

156 _DefaultName = "computeExposureSummaryStats" 

157 

158 def __init__(self, **kwargs): 

159 super().__init__(**kwargs) 

160 

161 self.makeSubtask("starSelector") 

162 

163 @timeMethod 

164 def run(self, exposure, sources, background): 

165 """Measure exposure statistics from the exposure, sources, and 

166 background. 

167 

168 Parameters 

169 ---------- 

170 exposure : `lsst.afw.image.ExposureF` 

171 sources : `lsst.afw.table.SourceCatalog` 

172 background : `lsst.afw.math.BackgroundList` 

173 

174 Returns 

175 ------- 

176 summary : `lsst.afw.image.ExposureSummary` 

177 """ 

178 self.log.info("Measuring exposure statistics") 

179 

180 summary = afwImage.ExposureSummaryStats() 

181 

182 bbox = exposure.getBBox() 

183 

184 psf = exposure.getPsf() 

185 self.update_psf_stats(summary, psf, bbox, sources, image_mask=exposure.mask) 

186 

187 wcs = exposure.getWcs() 

188 visitInfo = exposure.getInfo().getVisitInfo() 

189 self.update_wcs_stats(summary, wcs, bbox, visitInfo) 

190 

191 photoCalib = exposure.getPhotoCalib() 

192 self.update_photo_calib_stats(summary, photoCalib) 

193 

194 self.update_background_stats(summary, background) 

195 

196 self.update_masked_image_stats(summary, exposure.getMaskedImage()) 

197 

198 md = exposure.getMetadata() 

199 if 'SFM_ASTROM_OFFSET_MEAN' in md: 

200 summary.astromOffsetMean = md['SFM_ASTROM_OFFSET_MEAN'] 

201 summary.astromOffsetStd = md['SFM_ASTROM_OFFSET_STD'] 

202 

203 return summary 

204 

205 def update_psf_stats(self, summary, psf, bbox, sources=None, image_mask=None, sources_is_astropy=False): 

206 """Compute all summary-statistic fields that depend on the PSF model. 

207 

208 Parameters 

209 ---------- 

210 summary : `lsst.afw.image.ExposureSummaryStats` 

211 Summary object to update in-place. 

212 psf : `lsst.afw.detection.Psf` or `None` 

213 Point spread function model. If `None`, all fields that depend on 

214 the PSF will be reset (generally to NaN). 

215 bbox : `lsst.geom.Box2I` 

216 Bounding box of the image for which summary stats are being 

217 computed. 

218 sources : `lsst.afw.table.SourceCatalog` or `astropy.table.Table` 

219 Catalog for quantities that are computed from source table columns. 

220 If `None`, these quantities will be reset (generally to NaN). 

221 The type of this table must correspond to the 

222 ``sources_is_astropy`` argument. 

223 image_mask : `lsst.afw.image.Mask`, optional 

224 Mask image that may be used to compute distance-to-nearest-star 

225 metrics. 

226 sources_is_astropy : `bool`, optional 

227 Whether ``sources`` is an `astropy.table.Table` instance instead 

228 of an `lsst.afw.table.Catalog` instance. Default is `False` (the 

229 latter). 

230 """ 

231 nan = float("nan") 

232 summary.psfSigma = nan 

233 summary.psfIxx = nan 

234 summary.psfIyy = nan 

235 summary.psfIxy = nan 

236 summary.psfArea = nan 

237 summary.nPsfStar = 0 

238 summary.psfStarDeltaE1Median = nan 

239 summary.psfStarDeltaE2Median = nan 

240 summary.psfStarDeltaE1Scatter = nan 

241 summary.psfStarDeltaE2Scatter = nan 

242 summary.psfStarDeltaSizeMedian = nan 

243 summary.psfStarDeltaSizeScatter = nan 

244 summary.psfStarScaledDeltaSizeScatter = nan 

245 summary.maxDistToNearestPsf = nan 

246 summary.psfTraceRadiusDelta = nan 

247 

248 if psf is None: 

249 return 

250 shape = psf.computeShape(bbox.getCenter()) 

251 summary.psfSigma = shape.getDeterminantRadius() 

252 summary.psfIxx = shape.getIxx() 

253 summary.psfIyy = shape.getIyy() 

254 summary.psfIxy = shape.getIxy() 

255 im = psf.computeKernelImage(bbox.getCenter()) 

256 # The calculation of effective psf area is taken from 

257 # meas_base/src/PsfFlux.cc#L112. See 

258 # https://github.com/lsst/meas_base/blob/ 

259 # 750bffe6620e565bda731add1509507f5c40c8bb/src/PsfFlux.cc#L112 

260 summary.psfArea = float(np.sum(im.array)/np.sum(im.array**2.)) 

261 

262 if image_mask is not None: 

263 psfTraceRadiusDelta = psf_trace_radius_delta( 

264 image_mask, 

265 psf, 

266 sampling=self.config.psfGridSampling, 

267 bad_mask_bits=self.config.psfBadMaskPlanes 

268 ) 

269 summary.psfTraceRadiusDelta = float(psfTraceRadiusDelta) 

270 

271 if sources is None: 

272 # No sources are available (as in some tests and rare cases where 

273 # the selection criteria in finalizeCharacterization lead to no 

274 # good sources). 

275 return 

276 

277 # Count the total number of psf stars used (prior to stats selection). 

278 nPsfStar = sources[self.config.starSelection].sum() 

279 summary.nPsfStar = int(nPsfStar) 

280 

281 psf_mask = self.starSelector.run(sources).selected 

282 nPsfStarsUsedInStats = psf_mask.sum() 

283 

284 if nPsfStarsUsedInStats == 0: 

285 # No stars to measure statistics, so we must return the defaults 

286 # of 0 stars and NaN values. 

287 return 

288 

289 if sources_is_astropy: 

290 psf_cat = sources[psf_mask] 

291 else: 

292 psf_cat = sources[psf_mask].copy(deep=True) 

293 

294 starXX = psf_cat[self.config.starShape + '_xx'] 

295 starYY = psf_cat[self.config.starShape + '_yy'] 

296 starXY = psf_cat[self.config.starShape + '_xy'] 

297 psfXX = psf_cat[self.config.psfShape + '_xx'] 

298 psfYY = psf_cat[self.config.psfShape + '_yy'] 

299 psfXY = psf_cat[self.config.psfShape + '_xy'] 

300 

301 # Use the trace radius for the star size. 

302 starSize = np.sqrt(starXX/2. + starYY/2.) 

303 

304 starE1 = (starXX - starYY)/(starXX + starYY) 

305 starE2 = 2*starXY/(starXX + starYY) 

306 starSizeMedian = np.median(starSize) 

307 

308 # Use the trace radius for the psf size. 

309 psfSize = np.sqrt(psfXX/2. + psfYY/2.) 

310 psfE1 = (psfXX - psfYY)/(psfXX + psfYY) 

311 psfE2 = 2*psfXY/(psfXX + psfYY) 

312 

313 psfStarDeltaE1Median = np.median(starE1 - psfE1) 

314 psfStarDeltaE1Scatter = sigmaMad(starE1 - psfE1, scale='normal') 

315 psfStarDeltaE2Median = np.median(starE2 - psfE2) 

316 psfStarDeltaE2Scatter = sigmaMad(starE2 - psfE2, scale='normal') 

317 

318 psfStarDeltaSizeMedian = np.median(starSize - psfSize) 

319 psfStarDeltaSizeScatter = sigmaMad(starSize - psfSize, scale='normal') 

320 psfStarScaledDeltaSizeScatter = psfStarDeltaSizeScatter/starSizeMedian 

321 

322 summary.psfStarDeltaE1Median = float(psfStarDeltaE1Median) 

323 summary.psfStarDeltaE2Median = float(psfStarDeltaE2Median) 

324 summary.psfStarDeltaE1Scatter = float(psfStarDeltaE1Scatter) 

325 summary.psfStarDeltaE2Scatter = float(psfStarDeltaE2Scatter) 

326 summary.psfStarDeltaSizeMedian = float(psfStarDeltaSizeMedian) 

327 summary.psfStarDeltaSizeScatter = float(psfStarDeltaSizeScatter) 

328 summary.psfStarScaledDeltaSizeScatter = float(psfStarScaledDeltaSizeScatter) 

329 

330 if image_mask is not None: 

331 maxDistToNearestPsf = maximum_nearest_psf_distance( 

332 image_mask, 

333 psf_cat, 

334 sampling=self.config.psfSampling, 

335 bad_mask_bits=self.config.psfBadMaskPlanes 

336 ) 

337 summary.maxDistToNearestPsf = float(maxDistToNearestPsf) 

338 

339 def update_wcs_stats(self, summary, wcs, bbox, visitInfo): 

340 """Compute all summary-statistic fields that depend on the WCS model. 

341 

342 Parameters 

343 ---------- 

344 summary : `lsst.afw.image.ExposureSummaryStats` 

345 Summary object to update in-place. 

346 wcs : `lsst.afw.geom.SkyWcs` or `None` 

347 Astrometric calibration model. If `None`, all fields that depend 

348 on the WCS will be reset (generally to NaN). 

349 bbox : `lsst.geom.Box2I` 

350 Bounding box of the image for which summary stats are being 

351 computed. 

352 visitInfo : `lsst.afw.image.VisitInfo` 

353 Observation information used in together with ``wcs`` to compute 

354 the zenith distance. 

355 """ 

356 nan = float("nan") 

357 summary.raCorners = [nan]*4 

358 summary.decCorners = [nan]*4 

359 summary.ra = nan 

360 summary.dec = nan 

361 summary.zenithDistance = nan 

362 

363 if wcs is None: 

364 return 

365 

366 sph_pts = wcs.pixelToSky(geom.Box2D(bbox).getCorners()) 

367 summary.raCorners = [float(sph.getRa().asDegrees()) for sph in sph_pts] 

368 summary.decCorners = [float(sph.getDec().asDegrees()) for sph in sph_pts] 

369 

370 sph_pt = wcs.pixelToSky(bbox.getCenter()) 

371 summary.ra = sph_pt.getRa().asDegrees() 

372 summary.dec = sph_pt.getDec().asDegrees() 

373 

374 date = visitInfo.getDate() 

375 

376 if date.isValid(): 

377 # We compute the zenithDistance at the center of the detector 

378 # rather than use the boresight value available via the visitInfo, 

379 # because the zenithDistance may vary significantly over a large 

380 # field of view. 

381 observatory = visitInfo.getObservatory() 

382 loc = EarthLocation(lat=observatory.getLatitude().asDegrees()*units.deg, 

383 lon=observatory.getLongitude().asDegrees()*units.deg, 

384 height=observatory.getElevation()*units.m) 

385 obstime = Time(visitInfo.getDate().get(system=DateTime.MJD), 

386 location=loc, format='mjd') 

387 coord = SkyCoord( 

388 summary.ra*units.degree, 

389 summary.dec*units.degree, 

390 obstime=obstime, 

391 location=loc, 

392 ) 

393 with warnings.catch_warnings(): 

394 warnings.simplefilter('ignore') 

395 altaz = coord.transform_to(AltAz) 

396 

397 summary.zenithDistance = float(90.0 - altaz.alt.degree) 

398 

399 def update_photo_calib_stats(self, summary, photo_calib): 

400 """Compute all summary-statistic fields that depend on the photometric 

401 calibration model. 

402 

403 Parameters 

404 ---------- 

405 summary : `lsst.afw.image.ExposureSummaryStats` 

406 Summary object to update in-place. 

407 photo_calib : `lsst.afw.image.PhotoCalib` or `None` 

408 Photometric calibration model. If `None`, all fields that depend 

409 on the photometric calibration will be reset (generally to NaN). 

410 """ 

411 if photo_calib is not None: 

412 summary.zeroPoint = float(2.5*np.log10(photo_calib.getInstFluxAtZeroMagnitude())) 

413 else: 

414 summary.zeroPoint = float("nan") 

415 

416 def update_background_stats(self, summary, background): 

417 """Compute summary-statistic fields that depend only on the 

418 background model. 

419 

420 Parameters 

421 ---------- 

422 summary : `lsst.afw.image.ExposureSummaryStats` 

423 Summary object to update in-place. 

424 background : `lsst.afw.math.BackgroundList` or `None` 

425 Background model. If `None`, all fields that depend on the 

426 background will be reset (generally to NaN). 

427 

428 Notes 

429 ----- 

430 This does not include fields that depend on the background-subtracted 

431 masked image; when the background changes, it should generally be 

432 applied to the image and `update_masked_image_stats` should be called 

433 as well. 

434 """ 

435 if background is not None: 

436 bgStats = (bg[0].getStatsImage().getImage().array 

437 for bg in background) 

438 summary.skyBg = float(sum(np.median(bg[np.isfinite(bg)]) for bg in bgStats)) 

439 else: 

440 summary.skyBg = float("nan") 

441 

442 def update_masked_image_stats(self, summary, masked_image): 

443 """Compute summary-statistic fields that depend on the masked image 

444 itself. 

445 

446 Parameters 

447 ---------- 

448 summary : `lsst.afw.image.ExposureSummaryStats` 

449 Summary object to update in-place. 

450 masked_image : `lsst.afw.image.MaskedImage` or `None` 

451 Masked image. If `None`, all fields that depend 

452 on the masked image will be reset (generally to NaN). 

453 """ 

454 nan = float("nan") 

455 if masked_image is None: 

456 summary.skyNoise = nan 

457 summary.meanVar = nan 

458 return 

459 statsCtrl = afwMath.StatisticsControl() 

460 statsCtrl.setNumSigmaClip(self.config.sigmaClip) 

461 statsCtrl.setNumIter(self.config.clipIter) 

462 statsCtrl.setAndMask(afwImage.Mask.getPlaneBitMask(self.config.badMaskPlanes)) 

463 statsCtrl.setNanSafe(True) 

464 

465 statObj = afwMath.makeStatistics(masked_image, afwMath.STDEVCLIP, statsCtrl) 

466 skyNoise, _ = statObj.getResult(afwMath.STDEVCLIP) 

467 summary.skyNoise = skyNoise 

468 

469 statObj = afwMath.makeStatistics(masked_image.variance, masked_image.mask, afwMath.MEANCLIP, 

470 statsCtrl) 

471 meanVar, _ = statObj.getResult(afwMath.MEANCLIP) 

472 summary.meanVar = meanVar 

473 

474 

475def maximum_nearest_psf_distance( 

476 image_mask, 

477 psf_cat, 

478 sampling=8, 

479 bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"], 

480): 

481 """Compute the maximum distance of an unmasked pixel to its nearest PSF. 

482 

483 Parameters 

484 ---------- 

485 image_mask : `lsst.afw.image.Mask` 

486 The mask plane associated with the exposure. 

487 psf_cat : `lsst.afw.table.SourceCatalog` or `astropy.table.Table` 

488 Catalog containing only the stars used in the PSF modeling. 

489 sampling : `int` 

490 Sampling rate in each dimension to create the grid of points on which 

491 to evaluate the distance to the nearest PSF star. The tradeoff is 

492 between adequate sampling versus speed. 

493 bad_mask_bits : `list` [`str`] 

494 Mask bits required to be absent for a pixel to be considered 

495 "unmasked". 

496 

497 Returns 

498 ------- 

499 max_dist_to_nearest_psf : `float` 

500 The maximum distance (in pixels) of an unmasked pixel to its nearest 

501 PSF model star. 

502 """ 

503 mask_arr = image_mask.array[::sampling, ::sampling] 

504 bitmask = image_mask.getPlaneBitMask(bad_mask_bits) 

505 good = ((mask_arr & bitmask) == 0) 

506 

507 x = np.arange(good.shape[1]) * sampling 

508 y = np.arange(good.shape[0]) * sampling 

509 xx, yy = np.meshgrid(x, y) 

510 

511 dist_to_nearest_psf = np.full(good.shape, np.inf) 

512 for psf in psf_cat: 

513 x_psf = psf["slot_Centroid_x"] 

514 y_psf = psf["slot_Centroid_y"] 

515 dist_to_nearest_psf = np.minimum(dist_to_nearest_psf, np.hypot(xx - x_psf, yy - y_psf)) 

516 unmasked_dists = dist_to_nearest_psf * good 

517 max_dist_to_nearest_psf = np.max(unmasked_dists) 

518 

519 return max_dist_to_nearest_psf 

520 

521 

522def psf_trace_radius_delta( 

523 image_mask, 

524 image_psf, 

525 sampling=96, 

526 bad_mask_bits=["BAD", "CR", "INTRP", "SAT", "SUSPECT", "NO_DATA", "EDGE"], 

527): 

528 """Compute the delta between the maximum and minimum model PSF trace radius 

529 values evaluated on a grid of points lying in the unmasked region of the 

530 image. 

531 

532 Parameters 

533 ---------- 

534 image_mask : `lsst.afw.image.Mask` 

535 The mask plane associated with the exposure. 

536 image_psf : `lsst.afw.detection.Psf` 

537 The PSF model associated with the exposure. 

538 sampling : `int` 

539 Sampling rate in each dimension to create the grid of points at which 

540 to evaluate ``image_psf``s trace radius value. The tradeoff is between 

541 adequate sampling versus speed. 

542 bad_mask_bits : `list` [`str`] 

543 Mask bits required to be absent for a pixel to be considered 

544 "unmasked". 

545 

546 Returns 

547 ------- 

548 psf_trace_radius_delta : `float` 

549 The delta (in pixels) between the maximum and minimum model PSF trace 

550 radius values evaluated on the x,y-grid subsampled on the unmasked 

551 detector pixels by a factor of ``sampling``. If any model PSF trace 

552 radius value on the grid evaluates to NaN, then NaN is returned 

553 immediately. 

554 """ 

555 psf_trace_radius_list = [] 

556 mask_arr = image_mask.array[::sampling, ::sampling] 

557 bitmask = image_mask.getPlaneBitMask(bad_mask_bits) 

558 good = ((mask_arr & bitmask) == 0) 

559 

560 x = np.arange(good.shape[1]) * sampling 

561 y = np.arange(good.shape[0]) * sampling 

562 xx, yy = np.meshgrid(x, y) 

563 

564 for x_mesh, y_mesh, good_mesh in zip(xx, yy, good): 

565 for x_point, y_point, is_good in zip(x_mesh, y_mesh, good_mesh): 

566 if is_good: 

567 psf_trace_radius = image_psf.computeShape(geom.Point2D(x_point, y_point)).getTraceRadius() 

568 if ~np.isfinite(psf_trace_radius): 

569 return float("nan") 

570 psf_trace_radius_list.append(psf_trace_radius) 

571 

572 psf_trace_radius_delta = np.max(psf_trace_radius_list) - np.min(psf_trace_radius_list) 

573 

574 return psf_trace_radius_delta