Coverage for python / lsst / summit / utils / imageExaminer.py: 13%

333 statements  

« prev     ^ index     » next       coverage.py v7.13.5, created at 2026-05-05 19:01 +0000

1# This file is part of summit_utils. 

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__ = ["ImageExaminer"] 

23 

24from typing import cast 

25 

26import matplotlib 

27import matplotlib.patches as patches 

28import matplotlib.pyplot as plt 

29import numpy as np 

30import numpy.typing as npt 

31import scipy.ndimage as ndImage 

32from matplotlib import cm 

33from matplotlib.colors import LogNorm 

34from matplotlib.offsetbox import AnchoredText 

35from matplotlib.ticker import LinearLocator 

36from mpl_toolkits.mplot3d import Axes3D 

37from numpy.linalg import norm 

38from scipy.optimize import curve_fit 

39 

40import lsst.afw.image as afwImage 

41import lsst.geom as geom 

42import lsst.pipe.base as pipeBase 

43from lsst.pipe.tasks.quickFrameMeasurement import QuickFrameMeasurementTask, QuickFrameMeasurementTaskConfig 

44from lsst.summit.utils.utils import argMax2d, countPixels, getImageStats, quickSmooth 

45from lsst.utils.plotting.figures import make_figure 

46 

47SIGMATOFWHM = 2.0 * np.sqrt(2.0 * np.log(2.0)) 

48 

49 

50def gauss( 

51 x: float | npt.NDArray[np.float64], a: float, x0: float, sigma: float 

52) -> float | npt.NDArray[np.float64]: 

53 return a * np.exp(-((x - x0) ** 2) / (2 * sigma**2)) 

54 

55 

56class ImageExaminer: 

57 """Class for the reproducing some of the functionality of imexam. 

58 

59 For an input image create a summary plot showing: 

60 A rendering of the whole image 

61 A cutout of main source's PSF 

62 A 3d surface plot of the main star 

63 A contour plot of the main star 

64 x, y slices through the main star's centroid 

65 Radial plot of the main star 

66 Encircled energy as a function of radius 

67 A text box with assorted image statistics and measurements 

68 

69 Parameters 

70 ---------- 

71 exp : `lsst.afw.image.Exposure` 

72 The input exposure to analyze. 

73 doTweakCentroid : `bool`, optional 

74 Tweak the centroid (either the one supplied, or the one found by QFM if 

75 none supplied)? See ``tweakCentroid`` for full details of the 

76 behavior. 

77 doForceCoM : `bool`, optional 

78 Use the centre of mass inside the cutout box as the star's centroid? 

79 savePlots : `str`, optional 

80 Filename to save the plot to. Image not saved if falsey. 

81 centroid : `tuple` of `float`, optional 

82 Centroid of the star to treat as the main source. If ``None``, use 

83 ``lsst.pipe.tasks.quickFrameMeasurement.QuickFrameMeasurementTask`` to 

84 find the main source in the image. 

85 boxHalfSize : `int`, optional 

86 The half-size of the cutout to use for the star's PSF and the radius 

87 to use for the radial plots. 

88 

89 """ 

90 

91 astroMappings = { 

92 "object": "Object name", 

93 "mjd": "MJD", 

94 "expTime": "Exp Time", 

95 "filter": "Filter", 

96 "grating": "grating", 

97 "airmass": "Airmass", 

98 "rotangle": "Rotation Angle", 

99 "az": "Azimuth (deg)", 

100 "el": "Elevation (deg)", 

101 "focus": "Focus Z (mm)", 

102 } 

103 

104 imageMappings = { 

105 "centroid": "Centroid", 

106 "maxValue": "Max pixel value", 

107 "maxPixelLocation": "Max pixel location", 

108 "multipleMaxPixels": "Multiple max pixels?", 

109 "nBadPixels": "Num bad pixels", 

110 "nSatPixels": "Num saturated pixels", 

111 "percentile99": "99th percentile", 

112 "percentile9999": "99.99th percentile", 

113 "clippedMean": "Clipped mean", 

114 "clippedStddev": "Clipped stddev", 

115 } 

116 

117 cutoutMappings = { 

118 "nStatPixInBox": "nSat in cutout", 

119 "fitAmp": "Radial fitted amp", 

120 "fitGausMean": "Radial fitted position", 

121 "fitFwhm": "Radial fitted FWHM", 

122 "eeRadius50": "50% flux radius", 

123 "eeRadius80": "80% flux radius", 

124 "eeRadius90": "90% flux radius", 

125 } 

126 

127 def __init__( 

128 self, 

129 exp: afwImage.Exposure, 

130 *, 

131 doTweakCentroid: bool = True, 

132 doForceCoM: bool = False, 

133 savePlots: str | None = None, 

134 centroid: tuple[float, float] | None = None, 

135 boxHalfSize: int = 50, 

136 ): 

137 self.exp = exp 

138 self.savePlots = savePlots 

139 self.doTweakCentroid = doTweakCentroid 

140 self.doForceCoM = doForceCoM 

141 

142 self.boxHalfSize = boxHalfSize 

143 if centroid is None: 

144 qfmTaskConfig = QuickFrameMeasurementTaskConfig() 

145 qfmTask = QuickFrameMeasurementTask(config=qfmTaskConfig) 

146 result = qfmTask.run(exp) 

147 if not result.success: 

148 msg = ( 

149 "Failed to automatically find source in image. " 

150 "Either provide a centroid manually or use a new image" 

151 ) 

152 raise RuntimeError(msg) 

153 self.centroid = result.brightestObjCentroid 

154 else: 

155 self.centroid = centroid 

156 

157 self.imStats = getImageStats(self.exp) # need the background levels now 

158 

159 self.data = self.getStarBoxData() 

160 if self.doTweakCentroid: 

161 self.tweakCentroid(self.doForceCoM) 

162 self.data = self.getStarBoxData() 

163 

164 self.xx, self.yy = self.getMeshGrid(self.data) 

165 

166 self.imStats.centroid = self.centroid 

167 self.imStats.intCentroid = self.intCoords(self.centroid) 

168 self.imStats.intCentroidRounded = self.intRoundCoords(self.centroid) 

169 self.imStats.nStatPixInBox = self.nSatPixInBox 

170 

171 self.radialAverageAndFit() 

172 

173 def intCoords(self, coords: tuple[float | int, float | int]) -> npt.NDArray[np.integer]: 

174 """Get integer versions of the coordinates for dereferencing arrays. 

175 

176 Parameters are not rounded, but just cast as ints. 

177 

178 Parameters 

179 ---------- 

180 coords : `tuple` of `float` or `int` 

181 The coordinates. 

182 

183 Returns 

184 ------- 

185 intCoords : `np.array` of `int` 

186 The coordinates as integers. 

187 """ 

188 return np.asarray(coords, dtype=int) 

189 

190 def intRoundCoords(self, coords: tuple[float | int, float | int]) -> tuple[int, int]: 

191 """Get rounded integer versions of coordinates for dereferencing arrays 

192 

193 Parameters are rounded to the nearest integer value and returned. 

194 

195 Parameters 

196 ---------- 

197 coords : `tuple` of `float` or `int` 

198 The coordinates. 

199 

200 Returns 

201 ------- 

202 intCoords : Tuple[int, int] 

203 The coordinates as integers, rounded to the nearest values. 

204 """ 

205 return (int(round(coords[0])), int(round(coords[1]))) 

206 

207 def tweakCentroid(self, doForceCoM: bool) -> None: 

208 """Tweak the source centroid. Used to deal with irregular PSFs. 

209 

210 Given the star's cutout, tweak the centroid (either the one supplied 

211 manually, or the one from QFM) as follows: 

212 

213 If ``doForceCoM`` then always use the centre of mass of the cutout box 

214 as the centroid. 

215 If the star has multiple maximum values (e.g. if it is saturated and 

216 interpolated, or otherwise) then use the centre of mass of the cutout. 

217 Otherwise, use the position of the brightest pixel in the cutout. 

218 

219 Parameters 

220 ---------- 

221 doForceCoM : `bool` 

222 Forcing using the centre of mass of the cutout as the centroid? 

223 """ 

224 peak, uniquePeak, otherPeaks = argMax2d(self.data) 

225 # saturated stars don't tend to have ambiguous max pixels 

226 # due to the bunny ears left after interpolation 

227 nSatPix = self.nSatPixInBox 

228 

229 if not uniquePeak or nSatPix or doForceCoM: 

230 print("Using CoM for centroid (because was forced to, or multiple max pixels, or saturated") 

231 self.data -= self.imStats.clippedMean 

232 peak = ndImage.center_of_mass(self.data) 

233 self.data += self.imStats.clippedMean 

234 

235 offset = np.asarray(peak) - np.array((self.boxHalfSize, self.boxHalfSize)) 

236 print(f"Centroid adjusted by {offset} pixels") 

237 x = self.centroid[0] + offset[1] # yes, really, centroid is x,y offset is y,x 

238 y = self.centroid[1] + offset[0] 

239 self.centroid = (x, y) 

240 

241 def getStats(self) -> pipeBase.Struct: 

242 """Get the image stats. 

243 

244 Returns 

245 ------- 

246 stats : `dict` 

247 A dictionary of the image statistics. 

248 """ 

249 return self.imStats 

250 

251 @staticmethod 

252 def _calcMaxBoxHalfSize(centroid: tuple[float, float], chipBbox: geom.Box2I | geom.Box2D) -> int: 

253 """Calculate the maximum size the box can be without going outside the 

254 detector's bounds. 

255 

256 Returns the smallest distance between the centroid and any of the 

257 chip's edges. 

258 

259 Parameters 

260 ---------- 

261 centroid : `tuple` of `float` 

262 The centroid. 

263 chipBbox : `lsst.geom.Box` 

264 The detector's bounding box. 

265 

266 Returns 

267 ------- 

268 maxSize : `int` 

269 The maximum size for the box. 

270 """ 

271 ll = chipBbox.getBeginX() 

272 r = chipBbox.getEndX() 

273 d = chipBbox.getBeginY() 

274 u = chipBbox.getEndY() 

275 

276 x, y = np.array(centroid, dtype=int) 

277 maxSize = np.min([(x - ll), (r - x - 1), (u - y - 1), (y - d)]) # extra -1 in x because [) 

278 assert maxSize >= 0, "Box calculation went wrong" 

279 return maxSize 

280 

281 def _calcBbox(self, centroid: tuple[float, float]) -> geom.Box2I: 

282 """Get the largest valid bounding box, given the centroid and box size. 

283 

284 Parameters 

285 ---------- 

286 centroid : `tuple` of `float` 

287 The centroid 

288 

289 Returns 

290 ------- 

291 bbox : `lsst.geom.Box2I` 

292 The bounding box 

293 """ 

294 centroidPoint = geom.Point2I(centroid) 

295 extent = geom.Extent2I(1, 1) 

296 bbox = geom.Box2I(centroidPoint, extent) 

297 bbox = bbox.dilatedBy(self.boxHalfSize) 

298 bbox = bbox.clippedTo(self.exp.getBBox()) 

299 if bbox.getDimensions()[0] != bbox.getDimensions()[1]: 

300 # TODO: one day support clipped, nonsquare regions 

301 # but it's nontrivial due to all the plotting options 

302 

303 maxsize = self._calcMaxBoxHalfSize(centroid, self.exp.getBBox()) 

304 msg = ( 

305 f"With centroid at {centroid} and boxHalfSize {self.boxHalfSize} " 

306 "the selection runs off the edge of the chip. Boxsize has been " 

307 f"automatically shrunk to {maxsize} (only square selections are " 

308 "currently supported)" 

309 ) 

310 print(msg) 

311 self.boxHalfSize = maxsize 

312 return self._calcBbox(centroid) 

313 

314 return bbox 

315 

316 def getStarBoxData(self) -> npt.NDArray[np.float64]: 

317 """Get the image data for the star. 

318 

319 Calculates the maximum valid box, and uses that to return the image 

320 data, setting self.starBbox and self.nSatPixInBox as this method 

321 changes the bbox. 

322 

323 Returns 

324 ------- 

325 data : `np.array` 

326 The image data 

327 """ 

328 bbox = self._calcBbox(self.centroid) 

329 self.starBbox = bbox # needed elsewhere, so always set when calculated 

330 self.nSatPixInBox = countPixels(self.exp.maskedImage[self.starBbox], "SAT") 

331 return self.exp.image[bbox].array 

332 

333 def getMeshGrid( 

334 self, data: npt.NDArray[np.float64] 

335 ) -> tuple[npt.NDArray[np.floating], npt.NDArray[np.floating]]: 

336 """Get the meshgrid for a data array. 

337 

338 Parameters 

339 ---------- 

340 data : `np.array` 

341 The image data array. 

342 

343 Returns 

344 ------- 

345 xxyy : `tuple` of `np.array` 

346 The xx, yy as calculated by np.meshgrid 

347 """ 

348 xlen, ylen = data.shape 

349 x = np.arange(-1 * xlen / 2, xlen / 2, 1) 

350 y = np.arange(-1 * ylen / 2, ylen / 2, 1) 

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

352 return xx, yy 

353 

354 def radialAverageAndFit(self) -> None: 

355 """Calculate flux vs radius from the star's centroid and fit the width. 

356 

357 Calculate the flux vs distance from the star's centroid and fit 

358 a Gaussian to get a measurement of the width. 

359 

360 Also calculates the various encircled energy metrics. 

361 

362 Notes 

363 ----- 

364 Nothing is returned, but sets many value in the class. 

365 """ 

366 xlen, ylen = self.data.shape 

367 center = np.array([xlen / 2, ylen / 2]) 

368 # TODO: add option to move centroid to max pixel for radial (argmax 2d) 

369 

370 distances: list[float] = [] 

371 values: list[float] = [] 

372 

373 # could be much faster, but the array is tiny so its fine 

374 for i in range(xlen): 

375 for j in range(ylen): 

376 value = self.data[i, j] 

377 dist = norm((i, j) - center) 

378 fDist = float(dist) 

379 if fDist > xlen // 2: 

380 continue # clip to box size, we don't need a factor of sqrt(2) extra 

381 values.append(value) 

382 distances.append(fDist) 

383 

384 peakPos = 0 

385 amplitude = np.max(values) 

386 width = 10 

387 

388 bounds = ((0, 0, 0), (np.inf, np.inf, np.inf)) 

389 

390 try: 

391 pars, pCov = curve_fit(gauss, distances, values, [amplitude, peakPos, width], bounds=bounds) 

392 pars[0] = np.abs(pars[0]) 

393 pars[2] = np.abs(pars[2]) 

394 except RuntimeError: 

395 pars = None 

396 self.imStats.fitAmp = np.nan 

397 self.imStats.fitGausMean = np.nan 

398 self.imStats.fitFwhm = np.nan 

399 

400 if pars is not None: 

401 self.imStats.fitAmp = pars[0] 

402 self.imStats.fitGausMean = pars[1] 

403 self.imStats.fitFwhm = pars[2] * SIGMATOFWHM 

404 

405 self.radialDistances = distances 

406 self.radialValues = values 

407 

408 # calculate encircled energy metric too 

409 # sort distances and values in step by distance 

410 d = np.array([(r, v) for (r, v) in sorted(zip(self.radialDistances, self.radialValues))]) 

411 self.radii = d[:, 0] 

412 array_values = d[:, 1] 

413 self.cumFluxes = np.cumsum(array_values) 

414 self.cumFluxesNorm = self.cumFluxes / np.max(self.cumFluxes) 

415 

416 self.imStats.eeRadius50 = self.getEncircledEnergyRadius(50) 

417 self.imStats.eeRadius80 = self.getEncircledEnergyRadius(80) 

418 self.imStats.eeRadius90 = self.getEncircledEnergyRadius(90) 

419 

420 return 

421 

422 def getEncircledEnergyRadius(self, percentage: float | int) -> float: 

423 """Radius in pixels with the given percentage of encircled energy. 

424 

425 100% is at the boxHalfWidth dy definition. 

426 

427 Parameters 

428 ---------- 

429 percentage : `float` or `int` 

430 The percentage threshold to return. 

431 

432 Returns 

433 ------- 

434 radius : `float` 

435 The radius at which the ``percentage`` threshold is crossed. 

436 """ 

437 return self.radii[np.argmin(np.abs((percentage / 100) - self.cumFluxesNorm))] 

438 

439 def plotRadialAverage(self, ax: matplotlib.axes.Axes | None = None) -> None: 

440 """Make the radial average plot. 

441 

442 Parameters 

443 ---------- 

444 ax : `matplotlib.axes.Axes`, optional 

445 If ``None`` a new figure is created. Supply axes if including this 

446 as a subplot. 

447 """ 

448 plotDirect = False 

449 if not ax: 

450 ax = plt.subplot(111) 

451 plotDirect = True 

452 

453 distances = self.radialDistances 

454 values = self.radialValues 

455 pars = (self.imStats.fitAmp, self.imStats.fitGausMean, self.imStats.fitFwhm / SIGMATOFWHM) 

456 

457 fitFailed = np.isnan(pars).any() 

458 

459 ax.plot(distances, values, "x", label="Radial average") 

460 if not fitFailed: 

461 fitline = gauss(distances, *pars) # type: ignore 

462 ax.plot(distances, fitline, label="Gaussian fit") 

463 

464 ax.set_ylabel("Flux (ADU)") 

465 ax.set_xlabel("Radius (pix)") 

466 ax.set_aspect(1.0 / ax.get_data_ratio(), adjustable="box") # equal aspect for non-images 

467 ax.legend() 

468 

469 if plotDirect: 

470 plt.show() 

471 

472 def plotContours(self, ax: matplotlib.axes.Axes | None = None, nContours: int = 10) -> None: 

473 """Make the contour plot. 

474 

475 Parameters 

476 ---------- 

477 ax : `maplotlib.axes.Axes`, optional 

478 If ``None`` a new figure is created. Supply axes if including this 

479 as a subplot. 

480 nContours : `int`, optional 

481 The number of contours to use. 

482 """ 

483 plotDirect = False 

484 if not ax: 

485 fig = plt.figure(figsize=(8, 8)) # noqa F841 

486 ax = plt.subplot(111) 

487 plotDirect = True 

488 

489 vmin = np.percentile(self.data, 0.1) 

490 vmax = np.percentile(self.data, 99.9) 

491 lvls = np.linspace(vmin, vmax, nContours) 

492 intervalSize = lvls[1] - lvls[0] 

493 contourPlot = ax.contour(self.xx, self.yy, self.data, levels=lvls) # noqa F841 

494 print(f"Contoured from {vmin:,.0f} to {vmax:,.0f} using {nContours} contours of {intervalSize:.1f}") 

495 

496 ax.tick_params(which="both", direction="in", top=True, right=True, labelsize=8) 

497 ax.set_aspect("equal") 

498 

499 if plotDirect: 

500 plt.show() 

501 

502 def plotSurface(self, ax: Axes3D | None = None, useColor: bool = True) -> None: 

503 """Make the surface plot. 

504 

505 Parameters 

506 ---------- 

507 ax : `maplotlib.axes.Axes3D`, optional 

508 If ``None`` a new figure is created. Supply axes if including this 

509 as a subplot. 

510 useColor : `bool`, optional 

511 Plot at as a surface if ``True``, else plot as a wireframe. 

512 """ 

513 plotDirect = False 

514 if not ax: 

515 _, ax = plt.subplots(subplot_kw={"projection": "3d"}, figsize=(10, 10)) 

516 ax = cast(Axes3D, ax) 

517 plotDirect = True 

518 

519 if useColor: 

520 ax.plot_surface( 

521 self.xx, 

522 self.yy, 

523 self.data, 

524 cmap=cm.plasma, # type: ignore # mypy doesn't recognize dynamically created colormap attr. 

525 linewidth=1, 

526 antialiased=True, 

527 color="k", 

528 alpha=0.9, 

529 ) 

530 else: 

531 ax.plot_wireframe( 

532 self.xx, 

533 self.yy, 

534 self.data, 

535 cmap=cm.gray, # type: ignore # mypy doesn't recognize dynamically created colormap attributes 

536 linewidth=1, 

537 antialiased=True, 

538 color="k", 

539 ) 

540 

541 ax.zaxis.set_major_locator(LinearLocator(10)) 

542 ax.zaxis.set_major_formatter("{x:,.0f}") 

543 

544 if plotDirect: 

545 plt.show() 

546 

547 def plotStar(self, ax: matplotlib.axes.Axes | None = None, logScale: bool = False) -> None: 

548 """Make the PSF cutout plot. 

549 

550 Parameters 

551 ---------- 

552 ax : `maplotlib.axes`, optional 

553 If ``None`` a new figure is created. Supply axes if including this 

554 as a subplot. 

555 logScale : `bool`, optional 

556 Use a log scale? 

557 """ 

558 # TODO: display centroid in use 

559 plotDirect = False 

560 if not ax: 

561 ax = plt.subplot(111) 

562 plotDirect = True 

563 

564 interp = "none" 

565 if logScale: 

566 ax.imshow(self.data, norm=LogNorm(), origin="lower", interpolation=interp) 

567 else: 

568 ax.imshow(self.data, origin="lower", interpolation=interp) 

569 ax.tick_params(which="major", direction="in", top=True, right=True, labelsize=8) 

570 

571 xlen, ylen = self.data.shape 

572 center = np.array([xlen / 2, ylen / 2]) 

573 ax.plot(*center, "r+", markersize=10) 

574 ax.plot(*center, "rx", markersize=10) 

575 

576 if plotDirect: 

577 plt.show() 

578 

579 def plotFullExp(self, ax: matplotlib.axes.Axes | None = None) -> None: 

580 """Make the full image cutout plot. 

581 

582 Parameters 

583 ---------- 

584 ax : `maplotlib.axes`, optional 

585 If ``None`` a new figure is created. Supply axes if including this 

586 as a subplot. 

587 """ 

588 plotDirect = False 

589 if not ax: 

590 fig = plt.figure(figsize=(10, 10)) 

591 ax = fig.add_subplot(111) 

592 plotDirect = True 

593 

594 imData = quickSmooth(self.exp.image.array, 2.5) 

595 # np.percentile returns a float if imData is a 1D array, 

596 # but it can return a numpy array if imData is a 2D array 

597 # or higher. If imData is a 2D array, vmin and vmax will 

598 # be arrays, which is not what LogNorm expects. 

599 vmin = np.percentile(imData, 10) 

600 vmax = np.percentile(imData, 99.9) 

601 ax.imshow( 

602 imData, 

603 norm=LogNorm(vmin=vmin, vmax=vmax), # type: ignore 

604 origin="lower", 

605 cmap="gray_r", 

606 interpolation="bicubic", 

607 ) 

608 ax.tick_params(which="major", direction="in", top=True, right=True, labelsize=8) 

609 

610 xy0 = self.starBbox.getCorners()[0].x, self.starBbox.getCorners()[0].y 

611 width, height = self.starBbox.getWidth(), self.starBbox.getHeight() 

612 rect = patches.Rectangle(xy0, width, height, linewidth=1, edgecolor="r", facecolor="none") 

613 ax.add_patch(rect) 

614 

615 if plotDirect: 

616 plt.show() 

617 

618 def plotRowColSlices(self, ax: matplotlib.axes.Axes | None = None, logScale: bool = False) -> None: 

619 """Make the row and column slice plot. 

620 

621 Parameters 

622 ---------- 

623 ax : `maplotlib.axes`, optional 

624 If ``None`` a new figure is created. Supply axes if including this 

625 as a subplot. 

626 logScale : `bool`, optional 

627 Use a log scale? 

628 """ 

629 # TODO: display centroid in use 

630 

631 # slice through self.boxHalfSize because it's always the point being 

632 # used by definition 

633 rowSlice = self.data[self.boxHalfSize, :] 

634 colSlice = self.data[:, self.boxHalfSize] 

635 

636 plotDirect = False 

637 if not ax: 

638 ax = plt.subplot(111) 

639 plotDirect = True 

640 

641 xs = range(-1 * self.boxHalfSize, self.boxHalfSize + 1) 

642 ax.plot(xs, rowSlice, label="Row plot") 

643 ax.plot(xs, colSlice, label="Column plot") 

644 if logScale: 

645 pass 

646 # TODO: set yscale as log here also protect against negatives 

647 

648 ax.set_ylabel("Flux (ADU)") 

649 ax.set_xlabel("Radius (pix)") 

650 ax.set_aspect(1.0 / ax.get_data_ratio(), adjustable="box") # equal aspect for non-images 

651 

652 ax.legend() 

653 if plotDirect: 

654 plt.show() 

655 

656 def plotStats(self, ax: matplotlib.axes.Axes, lines: list[str]) -> None: 

657 """Make the stats box 'plot'. 

658 

659 Parameters 

660 ---------- 

661 ax : `maplotlib.axes.Axes` 

662 Axes to use. 

663 lines : `list` of `str` 

664 The data to include in the text box 

665 """ 

666 text = "\n".join([line for line in lines]) 

667 

668 stats_text = AnchoredText( 

669 text, 

670 loc="center", 

671 pad=0.5, 

672 prop=dict(size=14, ma="left", backgroundcolor="white", color="black", family="monospace"), 

673 ) 

674 ax.add_artist(stats_text) 

675 ax.axis("off") 

676 

677 def plotCurveOfGrowth(self, ax: matplotlib.axes.Axes | None = None) -> None: 

678 """Make the encircled energy plot. 

679 

680 Parameters 

681 ---------- 

682 ax : `maplotlib.axes.Axes`, optional 

683 If ``None`` a new figure is created. Supply axes if including this 

684 as a subplot. 

685 """ 

686 plotDirect = False 

687 if not ax: 

688 ax = plt.subplot(111) 

689 plotDirect = True 

690 

691 ax.plot(self.radii, self.cumFluxesNorm, markersize=10) 

692 ax.set_ylabel("Encircled flux (%)") 

693 ax.set_xlabel("Radius (pix)") 

694 

695 ax.set_aspect(1.0 / ax.get_data_ratio(), adjustable="box") # equal aspect for non-images 

696 

697 if plotDirect: 

698 plt.show() 

699 

700 def plot(self) -> matplotlib.figure.Figure: 

701 """Plot all the subplots together, including the stats box. 

702 

703 Image is saved if ``savefig`` was set. 

704 

705 Return 

706 ------ 

707 fig : `matplotlib.figure.Figure` 

708 The figure object. 

709 """ 

710 figsize = 6 

711 fig = make_figure(figsize=(figsize * 3, figsize * 2)) 

712 

713 ax1 = fig.add_subplot(331) 

714 ax2 = fig.add_subplot(332) 

715 ax3 = fig.add_subplot(333) 

716 ax4 = fig.add_subplot(334, projection="3d") 

717 ax5 = fig.add_subplot(335) 

718 ax6 = fig.add_subplot(336) 

719 ax7 = fig.add_subplot(337) 

720 ax8 = fig.add_subplot(338) 

721 ax9 = fig.add_subplot(339) 

722 

723 axExp = ax1 

724 axStar = ax2 

725 axStats1 = ax3 # noqa F841 - overwritten 

726 axSurf = ax4 

727 axCont = ax5 

728 axStats2 = ax6 # noqa F841 - overwritten 

729 axSlices = ax7 

730 axRadial = ax8 

731 axCoG = ax9 # noqa F841 - overwritten 

732 axStats1.axis("off") 

733 axStats2.axis("off") 

734 

735 self.plotFullExp(axExp) 

736 self.plotStar(axStar) 

737 self.plotSurface(axSurf) 

738 self.plotContours(axCont) 

739 self.plotRowColSlices(axSlices) 

740 self.plotRadialAverage(axRadial) 

741 

742 # overwrite three axes with this one spanning 3 rows 

743 gs = fig.add_gridspec(3, 3) 

744 axStats = fig.add_subplot(gs[0:2, 2]) # Spans rows 0-1 in column 2 

745 

746 lines = [] 

747 lines.append(" ---- Astro ----") 

748 lines.extend(self.translateStats(self.imStats, self.astroMappings)) 

749 lines.append("\n ---- Image ----") 

750 lines.extend(self.translateStats(self.imStats, self.imageMappings)) 

751 lines.append("\n ---- Cutout ----") 

752 lines.extend(self.translateStats(self.imStats, self.cutoutMappings)) 

753 self.plotStats(axStats, lines) 

754 

755 self.plotCurveOfGrowth(axCoG) 

756 

757 fig.tight_layout() 

758 if self.savePlots: 

759 print(f"Plot saved to {self.savePlots}") 

760 fig.savefig(self.savePlots) 

761 return fig 

762 

763 @staticmethod 

764 def translateStats(imStats: pipeBase.Struct, mappingDict: dict[str, str]) -> list[str]: 

765 """Create the text for the stats box from the stats themselves. 

766 

767 Parameters 

768 ---------- 

769 imStats : `lsst.pipe.base.Struct` 

770 A container with attributes containing measurements and statistics 

771 for the image. 

772 mappingDict : `dict` of `str` 

773 A mapping from attribute name to name for rendereding as text. 

774 

775 Returns 

776 ------- 

777 lines : `list` of `str` 

778 The translated lines of text. 

779 """ 

780 lines = [] 

781 for k, v in mappingDict.items(): 

782 try: 

783 value = getattr(imStats, k) 

784 except Exception: 

785 lines.append("") 

786 continue 

787 

788 # native floats are not np.floating so must check both 

789 if isinstance(value, float) or isinstance(value, np.floating): 

790 value = f"{value:,.3f}" 

791 if k == "centroid": # special case the only tuple 

792 value = f"{value[0]:.1f}, {value[1]:.1f}" 

793 lines.append(f"{v} = {value}") 

794 return lines 

795 

796 def plotAll(self) -> None: 

797 """Make each of the plots, individually. 

798 

799 Makes all the plots, full size, one by one, as opposed to plot() which 

800 creates a single image containing all the plots. 

801 """ 

802 self.plotStar() 

803 self.plotRadialAverage() 

804 self.plotContours() 

805 self.plotSurface() 

806 self.plotStar() 

807 self.plotRowColSlices()