Coverage for python / lsst / analysis / tools / actions / plot / propertyMapPlot.py: 9%

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/>. 

21from __future__ import annotations 

22 

23__all__ = ( 

24 "PerTractPropertyMapPlot", 

25 "SurveyWidePropertyMapPlot", 

26) 

27 

28import logging 

29from collections.abc import Mapping 

30 

31import matplotlib.patheffects as mpl_path_effects 

32import numpy as np 

33import skyproj 

34from healsparse.healSparseMap import HealSparseMap 

35from matplotlib import cm, rc_context 

36from matplotlib.figure import Figure 

37from matplotlib.legend_handler import HandlerTuple 

38from mpl_toolkits.axes_grid1 import make_axes_locatable 

39from mpl_toolkits.axes_grid1.axes_size import AxesX, AxesY, Fraction 

40 

41import lsst.pex.config as pexConfig 

42import lsst.sphgeom as sphgeom 

43from lsst.analysis.tools.tasks.propertyMapAnalysis import ( 

44 PerTractPropertyMapAnalysisConfig, 

45 SurveyWidePropertyMapAnalysisConfig, 

46) 

47from lsst.skymap.tractInfo import ExplicitTractInfo 

48from lsst.utils.plotting import make_figure, set_rubin_plotstyle 

49 

50from ...interfaces import KeyedData, PlotAction 

51 

52_LOG = logging.getLogger(__name__) 

53 

54# Holds unit renames to match style guidelines: {"old_unit": "new_unit"}. 

55unitRenameDict = {"mag(AB)": r"$\rm mag_{AB}$"} 

56 

57 

58def getZoomedExtent(box, n): 

59 """Get a zoomed plot extent in degrees for a bounding box. 

60 

61 Parameters 

62 ---------- 

63 box : `lsst.sphgeom.NormalizedAngleInterval` 

64 The bounding box to get an extent for. 

65 n : `float` 

66 Zoom factor; for instance, n=2 means zooming in 2 times at the 

67 center. Must be positive. 

68 

69 Returns 

70 ------- 

71 `tuple` [`float`] 

72 New extent as (new_lon_min, new_lon_max, new_lat_min, new_lat_max). 

73 """ 

74 lon = box.getLon() 

75 lat = box.getLat() 

76 if n > 1: 

77 factor = (1.0 - 1.0 / n) / 2.0 

78 lon.erodeBy(sphgeom.NormalizedAngle(factor * box.getWidth())) 

79 lat.erodeBy(sphgeom.NormalizedAngle(factor * box.getHeight())) 

80 elif n < 1: 

81 factor = (1.0 / n - 1) / 2.0 

82 lon.dilateBy(sphgeom.NormalizedAngle(factor * box.getWidth())) 

83 lat.dilateBy(sphgeom.NormalizedAngle(factor * box.getHeight())) 

84 

85 extent = ( 

86 lon.getA().asDegrees(), 

87 lon.getB().asDegrees(), 

88 lat.getA().asDegrees(), 

89 lat.getB().asDegrees(), 

90 ) 

91 return extent 

92 

93 

94def getLongestSuffixMatch(s, options): 

95 """Find the longest suffix in the provided list that matches the end of 

96 the given string. 

97 

98 Parameters 

99 ---------- 

100 s : `str` 

101 The target string for which we want to find a matching suffix. 

102 options : `list` [`str`] 

103 A list of potential suffix strings to match against the target 

104 string `s`. 

105 

106 Returns 

107 ------- 

108 `str` 

109 The longest matching suffix from the `options` list. If no match is 

110 found, returns `None`. 

111 """ 

112 return next((opt for opt in sorted(options, key=len, reverse=True) if s.endswith(opt)), None) 

113 

114 

115def addTextToColorbar( 

116 cb, text, orientation="vertical", color="black", fontsize=14, fontweight="bold", alpha=0.8 

117): 

118 """Helper method to add text inside the horizontal colorbar. 

119 

120 Parameters 

121 ---------- 

122 cb : `~matplotlib.colorbar.Colorbar` 

123 The colorbar object. 

124 text : `str` 

125 The text to add. 

126 orientation : `str`, optional 

127 The orientation of the colorbar. Can be either "vertical" or 

128 "horizontal". 

129 fontsize : `int`, optional 

130 The fontsize of the text. 

131 fontweight : `str`, optional 

132 The fontweight of the text. 

133 alpha : `float`, optional 

134 The alpha value of the text. 

135 

136 Returns 

137 ------- 

138 `None` 

139 The text is added to the colorbar in place. 

140 """ 

141 if color is None: 

142 color = "black" 

143 vmid = (cb.vmin + cb.vmax) / 2 

144 positions = {"vertical": (0.5, vmid), "horizontal": (vmid, 0.5)} 

145 cbtext = cb.ax.text( 

146 *positions[orientation], 

147 text, 

148 color=color, 

149 va="center", 

150 ha="center", 

151 fontsize=fontsize, 

152 fontweight=fontweight, 

153 rotation=orientation, 

154 alpha=alpha, 

155 ) 

156 # Add a distinct outline around the text for better visibility in 

157 # various backgrounds. 

158 cbtext.set_path_effects( 

159 [mpl_path_effects.Stroke(linewidth=4, foreground="white", alpha=0.8), mpl_path_effects.Normal()] 

160 ) 

161 

162 

163class CustomHandler(HandlerTuple): 

164 """Custom legend handler to overlay multiple patches for a single 

165 legend entry. 

166 

167 This handler class inherits from `HandlerTuple` and is designed to 

168 handle cases where multiple artists (e.g., patches) need to be overlaid 

169 on top of each other in a single legend entry, as opposed to 

170 side-by-side which is the default behavior of `HandlerTuple`. 

171 

172 Methods 

173 ------- 

174 create_artists: 

175 Override the `create_artists` method of `HandlerTuple` to modify 

176 the positioning of the artists so that they overlay directly on top 

177 of one another in the legend. 

178 

179 Example 

180 ------- 

181 # Plot some data. 

182 line, = ax.plot(x, y, label="Sample Line") 

183 

184 # Use CustomHandler for overlaid patches and also include the regular 

185 # line legend if desired. 

186 handles = [(patch1, patch2), line] 

187 labels = ['Overlaid Patches', line.get_label()] 

188 leg = ax.legend( 

189 handles, labels, handler_map={tuple: CustomHandler()}, loc="best" 

190 ) 

191 """ 

192 

193 def create_artists(self, legend, orig_handle, xdescent, ydescent, width, height, fontsize, trans): 

194 artists = HandlerTuple.create_artists( 

195 self, legend, orig_handle, xdescent, ydescent, width, height, fontsize, trans 

196 ) 

197 # Overlay the two patches. 

198 for a in artists: 

199 a.set_transform(trans) 

200 return artists 

201 

202 

203class PerTractPropertyMapPlot(PlotAction): 

204 draw_patch_bounds = pexConfig.Field[bool]( 

205 doc="Whether to draw patch inner boundaries or not", 

206 default=False, 

207 ) 

208 label_patches = pexConfig.ChoiceField[str]( 

209 doc="Which patches to label by ID", 

210 default="none", 

211 allowed={ 

212 "all": "All patches", 

213 "edge": "Edge patches only", 

214 "none": "No labels", 

215 }, 

216 ) 

217 plotName = pexConfig.Field[str](doc="The name for the plotting task.", optional=True) 

218 

219 def __call__( 

220 self, 

221 data: KeyedData, 

222 tractInfo: ExplicitTractInfo, 

223 plotConfig: PerTractPropertyMapAnalysisConfig, 

224 plotInfo: Mapping[str, Mapping[str, str] | str | int], 

225 **kwargs, 

226 ) -> Mapping[str, Figure]: 

227 self._validateInput(data, tractInfo, plotConfig, plotInfo) 

228 return self.makePlot(data, tractInfo, plotConfig, plotInfo) 

229 

230 def _validateInput( 

231 self, 

232 data: KeyedData, 

233 tractInfo: ExplicitTractInfo, 

234 plotConfig: PerTractPropertyMapAnalysisConfig, 

235 plotInfo: Mapping[str, Mapping[str, str] | str | int], 

236 ) -> None: 

237 """Validate the input data.""" 

238 

239 if not isinstance(tractInfo, ExplicitTractInfo): 

240 raise TypeError(f"Input `tractInfo` type must be {ExplicitTractInfo} not {type(tractInfo)}.") 

241 

242 if not isinstance(plotConfig, PerTractPropertyMapAnalysisConfig): 

243 raise TypeError( 

244 "`plotConfig` must be a " 

245 "`lsst.analysis.tools.tasks.propertyMapAnalysis.PerTractPropertyMapAnalysisConfig`. " 

246 f"Got {type(plotConfig)}." 

247 ) 

248 

249 if not isinstance(plotInfo, dict): 

250 raise TypeError("`plotConfig` must be a dictionary.") 

251 

252 if not plotConfig.publicationStyle: 

253 zoomFactors = plotConfig.zoomFactors 

254 isListOfFloats = isinstance(zoomFactors, pexConfig.listField.List) and all( 

255 isinstance(zf, float) for zf in zoomFactors 

256 ) 

257 if not (isListOfFloats and len(zoomFactors) == 2) or any(zf <= 1 for zf in zoomFactors): 

258 raise TypeError( 

259 "`zoomFactors` must be a two-element `lsst.pex.config.listField.List` of floats > 1." 

260 ) 

261 

262 for atool in plotConfig.atools: 

263 if not isinstance(atool.nBinsHist, int) or atool.nBinsHist <= 0: 

264 raise ValueError( 

265 f"`nBinsHist` for property `{atool.process.buildActions.data.mapKey}` must be a positive " 

266 f"integer. Got {atool.nBinsHist}." 

267 ) 

268 

269 # Identify any invalid entries in `data`. 

270 invalidEntries = { 

271 key: pytype 

272 for key, pytype in {k: v.ref.datasetType.storageClass.pytype for k, v in data.items()}.items() 

273 if pytype != HealSparseMap 

274 } 

275 

276 # If any invalid entries are found, raise a TypeError with details. 

277 if invalidEntries: 

278 errorMessage = "; ".join( 

279 f"`{key}` should be {HealSparseMap}, got {type_}" for key, type_ in invalidEntries.items() 

280 ) 

281 raise TypeError(f"Invalid input types found in `data`: {errorMessage}") 

282 

283 def addPlotInfo( 

284 self, 

285 fig: Figure, 

286 plotInfo: Mapping[str, Mapping[str, str] | str | int], 

287 toolName: str, 

288 ) -> Figure: 

289 """Add useful information to the plot. 

290 

291 Parameters 

292 ---------- 

293 fig : `matplotlib.figure.Figure` 

294 The figure to add the information to. 

295 plotInfo : `dict` 

296 A dictionary of the plot information. 

297 toolName : `str` 

298 The name of the tool used to generate the plot. 

299 

300 Returns 

301 ------- 

302 fig : `matplotlib.figure.Figure` 

303 The figure with the information added. 

304 """ 

305 

306 run = plotInfo["run"] 

307 tableType = f"\nTable: {plotInfo['tableNames'][toolName]}" 

308 

309 dataIdText = f"Tract: {plotInfo['tract']}, Band: {plotInfo['band']}" 

310 propertyDescription = plotInfo["description"] 

311 # Lowercase the first letter unless the string starts with more than 

312 # one uppercase letter in which case we keep it as is, e.g. PSF, DCR. 

313 if propertyDescription[0].isupper() and propertyDescription[1].islower(): 

314 propertyDescription = propertyDescription[0].lower() + propertyDescription[1:] 

315 mapText = ( 

316 f", Property: {propertyDescription}, " 

317 f"Unit: {plotInfo['unit']}, " 

318 f"Operation: {plotInfo['operation']}, " 

319 f"Coadd: {plotInfo['coaddName']}" 

320 ) 

321 geomText = ( 

322 f", Valid area: {plotInfo['valid_area']:.2f} sq. deg., " 

323 + f"NSIDE: {plotInfo['nside']}, projection: {plotInfo['projection']}" 

324 ) 

325 infoText = f"\n{dataIdText}{mapText}" 

326 

327 fig.text( 

328 0.04, 

329 0.965, 

330 f'{plotInfo["plotName"]} of {plotInfo["property"]}', 

331 fontsize=19, 

332 transform=fig.transFigure, 

333 ha="left", 

334 va="top", 

335 ) 

336 t = fig.text( 

337 0.04, 

338 0.942, 

339 f"{run}{tableType}{geomText}{infoText}", 

340 fontsize=15, 

341 transform=fig.transFigure, 

342 alpha=0.6, 

343 ha="left", 

344 va="top", 

345 ) 

346 t.set_linespacing(1.4) 

347 

348 return fig 

349 

350 def makePlot( 

351 self, 

352 data: KeyedData, 

353 tractInfo: ExplicitTractInfo, 

354 plotConfig: PerTractPropertyMapAnalysisConfig, 

355 plotInfo: Mapping[str, Mapping[str, str] | str | int], 

356 ) -> Mapping[str, Figure]: 

357 """Make the survey property map plot. 

358 

359 Parameters 

360 ---------- 

361 data : `KeyedData` 

362 The HealSparseMap to plot the points from. 

363 tractInfo: `~lsst.skymap.tractInfo.ExplicitTractInfo` 

364 The tract info object. 

365 plotConfig : 

366 `~lsst.analysis.tools.tasks.perTractPropertyMapAnalysis. 

367 PerTractPropertyMapAnalysisConfig` 

368 The configuration for the plot. 

369 plotInfo : `dict` 

370 A dictionary of information about the data being plotted. 

371 

372 Returns 

373 ------- 

374 fig : `~matplotlib.figure.Figure` 

375 The resulting figure. 

376 """ 

377 

378 set_rubin_plotstyle() 

379 

380 # 'plotName' by default is constructed from the attribute specified in 

381 # 'atools.<attribute>' in the pipeline YAML. If it is explicitly 

382 # set in `~lsst.analysis.tools.atools.propertyMap.PropertyMapTool`, 

383 # it will override this default. 

384 if self.plotName: 

385 # Set the plot name using 'produce.plot.plotName' from 

386 # PropertyMapTool's instance. 

387 plotInfo["plotName"] = self.plotName 

388 

389 # Plotting customization. 

390 colorbarLocation = "right" 

391 colorbarTickLabelSize = 14 

392 colorBarAspect = 16 

393 colorBarFraction = 0.07 # w.r.t. the axes size. 

394 colorBarPad = 0 # w.r.t. the axes size. 

395 rcparams = { 

396 "axes.labelsize": 22 if plotConfig.publicationStyle else 19, 

397 "axes.linewidth": 1.8, 

398 "xtick.labelsize": 15 if plotConfig.publicationStyle else 14, 

399 "ytick.labelsize": 15 if plotConfig.publicationStyle else 14, 

400 } 

401 

402 colorbarKwargs = dict(plotConfig.colorbarKwargs) 

403 projectionKwargs = dict(plotConfig.projectionKwargs) 

404 

405 if plotConfig.publicationStyle: 

406 zoomFactors = [] 

407 _LOG.info( 

408 "Zoom factors are not used in publication-style plots. " 

409 "Only the full-tract map will be plotted." 

410 ) 

411 if "cmap" in colorbarKwargs and colorbarKwargs["cmap"] != "viridis": 

412 _LOG.warning( 

413 "Color map set to 'viridis' for publication style plots. " 

414 f"The color map '{colorbarKwargs['cmap']}' set in the config will be ignored." 

415 ) 

416 colorbarKwargs["cmap"] = "viridis" 

417 cmap = cm.get_cmap(colorbarKwargs["cmap"]) 

418 # Colorbar text color only, not used for any histogram. 

419 histColors = [cmap(0.1)] 

420 labelpad = 14 

421 else: 

422 zoomFactors = plotConfig.zoomFactors 

423 colorbarKwargs["cmap"] = colorbarKwargs.get("cmap", "viridis") 

424 # Muted green for the full map, and muted red and blue for the two 

425 # zoomed-in maps. Used for boxes, colorbar texts and histograms. 

426 histColors = ["#265D40", "#8B0000", "#00008B"] 

427 labelpad = 11 

428 

429 toolName = data["data"].ref.datasetType.name 

430 mapName = toolName.replace("_map_", "_") 

431 mapData = data["data"].get() 

432 

433 with rc_context(rcparams): 

434 if plotConfig.publicationStyle: 

435 # Make a single plot for the full tract. 

436 fig = make_figure(figsize=(8, 8)) 

437 ax1 = fig.add_subplot(111) 

438 fig.subplots_adjust(left=0.133, right=0.888, bottom=0.1, top=0.855) 

439 else: 

440 # Make a 2x2 grid of subplots for the full tract, two zoomed-in 

441 # views, and a histogram of values. 

442 fig = make_figure(figsize=(16, 16)) 

443 ax1 = fig.add_subplot(221) 

444 ax2 = fig.add_subplot(222) 

445 ax3 = fig.add_subplot(223) 

446 ax4 = fig.add_subplot(224) 

447 

448 # Reduce whitespace but leave some room at the top for info. 

449 fig.subplots_adjust(left=0.064, right=0.96, top=0.855, bottom=0.07, wspace=0.275, hspace=0.24) 

450 

451 # Get the values for the valid pixels of the full tract. 

452 values = mapData[mapData.valid_pixels] 

453 goodValues = np.isfinite(values) 

454 values = values[goodValues] # As a precaution. 

455 

456 # Make a concise human-readable label for the plot. 

457 plotInfo["unit"] = "N/A" # Unless overridden. 

458 if hasattr(mapData, "metadata") and all( 

459 key in mapData.metadata for key in ["DESCRIPTION", "OPERATION", "UNIT"] 

460 ): 

461 hasMetadata = True 

462 metadata = mapData.metadata 

463 plotInfo["description"] = metadata["DESCRIPTION"] 

464 plotInfo["operation"] = metadata["OPERATION"] 

465 if metadata["UNIT"]: 

466 plotInfo["unit"] = metadata["UNIT"] 

467 elif metadata["UNIT"] == "": 

468 plotInfo["unit"] = "dimensionless" 

469 else: 

470 hasMetadata = False 

471 plotInfo["operation"] = getLongestSuffixMatch( 

472 mapName, ["min", "max", "mean", "weighted_mean", "sum"] 

473 ).replace("_", " ") 

474 plotInfo["coaddName"] = mapName.split("Coadd_")[0] 

475 plotInfo["operation"] = plotInfo["operation"].replace("minimum", "min").replace("maximum", "max") 

476 propertyName = mapName[len(f"{plotInfo['coaddName']}Coadd_") : -len(plotInfo["operation"])].strip( 

477 "_" 

478 ) 

479 if not hasMetadata: 

480 # Infer the property description from the map name (all 

481 # lower case), and properly handle formatting. 

482 plotInfo["description"] = ( 

483 propertyName.replace("_", " ") 

484 .replace("psf", "PSF") 

485 .replace("dcr", "DCR") 

486 .replace("dra", r"$\Delta$RA") 

487 .replace("ddec", r"$\Delta$Dec") 

488 ) 

489 plotInfo["property"] = ( 

490 propertyName.replace("_", " ") 

491 .title() # Capitalize and handle edge cases below. 

492 .replace("Psf", "PSF") 

493 .replace("Dcr", "DCR") 

494 .replace("Dra", r"$\Delta$RA") 

495 .replace("Ddec", r"$\Delta$Dec") 

496 .replace("E1", "e1") 

497 .replace("E2", "e2") 

498 ) 

499 

500 # Handle unit renaming. 

501 if plotInfo["unit"] in unitRenameDict: 

502 plotInfo["unit"] = unitRenameDict[plotInfo["unit"]] 

503 

504 atool = getattr(plotConfig.atools, toolName) 

505 nBinsHist = atool.nBinsHist 

506 

507 # Use the full tract bounding box to set the default extent. 

508 box = tractInfo.getOuterSkyPolygon().getBoundingBox() 

509 ctr_lonlat = box.getCenter() 

510 ctr_lon = ctr_lonlat.getLon().asDegrees() 

511 ctr_lat = ctr_lonlat.getLat().asDegrees() 

512 

513 # Prepare plot elements for the full tract and optional zoomed-in 

514 # views. 

515 if plotConfig.publicationStyle: 

516 skyprojAxes = [ax1] 

517 else: 

518 skyprojAxes = [ax1, ax3, ax4] 

519 

520 zoomIdx = [] 

521 for ax, zoomFactor, histColor in zip(skyprojAxes, [1.0, *zoomFactors], histColors): 

522 extent = getZoomedExtent(box, zoomFactor) 

523 

524 sp = skyproj.GnomonicSkyproj( 

525 ax=ax, 

526 lon_0=ctr_lon, 

527 lat_0=ctr_lat, 

528 extent=extent, 

529 **projectionKwargs, 

530 ) 

531 

532 sp.draw_hspmap(mapData, zoom=False, cmap=colorbarKwargs["cmap"]) 

533 

534 sp.ax.set_xlabel("R.A.", labelpad=labelpad, fontsize=rcparams["axes.labelsize"]) 

535 sp.ax.set_ylabel("Dec.", labelpad=labelpad, fontsize=rcparams["axes.labelsize"]) 

536 

537 if self.draw_patch_bounds: 

538 label_all = self.label_patches == "all" 

539 label_edge = self.label_patches == "edge" 

540 if label_all or label_edge: 

541 patchids_max = tuple(x - 1 for x in tractInfo.getNumPatches()) 

542 lon_min, lon_max = min(extent[:2]), max(extent[:2]) 

543 lat_min, lat_max = min(extent[2:]), max(extent[2:]) 

544 

545 for patchInfo in tractInfo: 

546 vertices = patchInfo.getInnerSkyPolygon().getVertices() 

547 clipped = tractInfo.inner_sky_region.clipTo( 

548 sphgeom.Box(sphgeom.LonLat(vertices[0]), sphgeom.LonLat(vertices[2])) 

549 ) 

550 lonlats_patch = np.array( 

551 [ 

552 [x.asDegrees() for x in (lonlat.getA(), lonlat.getB())] 

553 for lonlat in (clipped.getLon(), clipped.getLat()) 

554 ] 

555 ) 

556 lons_patch = np.concat((lonlats_patch[0, :], lonlats_patch[0, ::-1])) 

557 lats_patch = np.repeat(lonlats_patch[1, :], 2) 

558 

559 sp.draw_polygon(lons_patch, lats_patch, edgecolor="gray") 

560 label_id = label_all or ( 

561 label_edge 

562 and ( 

563 (patchInfo.index[0] == 0) 

564 or (patchInfo.index[1] == 0) 

565 or (patchInfo.index[0] == patchids_max[0]) 

566 or (patchInfo.index[1] == patchids_max[1]) 

567 ) 

568 ) 

569 if label_id: 

570 lon_label = clipped.getCenter().getLon().asDegrees() 

571 lat_label = clipped.getCenter().getLat().asDegrees() 

572 if (lon_min < lon_label < lon_max) and (lat_min < lat_label < lat_max): 

573 sp.ax.text( 

574 lon_label, 

575 lat_label, 

576 patchInfo.getSequentialIndex(), 

577 ha="center", 

578 va="center", 

579 size=7, 

580 c=[0, 0, 0, 0.5], 

581 ) 

582 

583 # Specify the size and padding of the colorbar axes with 

584 # respect to the main axes. 

585 refAx = AxesX(ax) if colorbarLocation in ("left", "right") else AxesY(ax) 

586 cbsize = Fraction(colorBarFraction, refAx) 

587 cbpad = Fraction(colorBarPad, refAx) 

588 

589 # Make divider and a colorbar axes to be attached to the main 

590 # axes. 

591 divider = make_axes_locatable(ax) 

592 cax = divider.append_axes(colorbarLocation, size=cbsize, pad=cbpad) 

593 

594 cbar = sp.draw_colorbar( 

595 **{ 

596 "location": colorbarLocation, 

597 "aspect": colorBarAspect, 

598 "fraction": colorBarFraction, 

599 "pad": colorBarPad, 

600 "cax": cax, 

601 **colorbarKwargs, 

602 } 

603 ) 

604 cbar.ax.tick_params(labelsize=colorbarTickLabelSize) 

605 if plotConfig.publicationStyle: 

606 if plotInfo["unit"] not in ["dimensionless", "N/A"]: 

607 cbarText = f"{plotInfo['property']} ({plotInfo['unit']})" 

608 else: 

609 cbarText = plotInfo["property"] 

610 else: 

611 cbarText = ( 

612 "Full Tract" if zoomFactor == 1.0 else f"{self.prettyPrintFloat(zoomFactor)}x Zoom" 

613 ) 

614 addTextToColorbar(cbar, cbarText, color=histColor) 

615 if zoomFactor == 1.0: 

616 # Store the "full tract" map so that we can overplot 

617 # the zoom rectangles. 

618 spf = sp 

619 else: 

620 # Create a rectangle for the zoomed-in region. 

621 x0, x1, y0, y1 = extent 

622 for c, ls, lw in zip(["white", histColor], ["solid", "dashed"], [3.5, 1.5]): 

623 spf.draw_polygon( 

624 [x0, x0, x1, x1], 

625 [y0, y1, y1, y0], 

626 facecolor="none", 

627 edgecolor=c, 

628 linestyle=ls, 

629 linewidth=lw, 

630 alpha=0.8, 

631 ) 

632 zoomText = spf.ax.text( 

633 (x0 + x1) / 2, 

634 y0, 

635 f"{self.prettyPrintFloat(zoomFactor)}x", 

636 color=histColor, 

637 fontsize=14, 

638 fontweight="bold", 

639 alpha=0.8, 

640 ha="center", 

641 va="bottom", 

642 ) 

643 # Add a distinct outline around the text for better 

644 # visibility in various backgrounds. 

645 zoomText.set_path_effects( 

646 [ 

647 mpl_path_effects.Stroke(linewidth=4, foreground="white", alpha=0.8), 

648 mpl_path_effects.Normal(), 

649 ] 

650 ) 

651 # Get the indices of pixels in the zoomed-in region. 

652 pos = mapData.valid_pixels_pos() 

653 # Reversed axes consideration. 

654 xmin, xmax = sorted([x0, x1]) 

655 idx = (pos[0] > xmin) & (pos[0] < xmax) & (pos[1] > y0) & (pos[1] < y1) 

656 zoomIdx.append(idx[goodValues]) 

657 

658 # Calculate weights for each bin to ensure that the peak of the 

659 # histogram reaches 1. 

660 if len(values) == 0 or np.ptp(values) < 1e-12: 

661 # np.histogram cannot use 100 bins if no variance in values 

662 weights = np.ones_like(values) / values.size 

663 nBinsHist = 1 

664 _LOG.info(f"No variance in {toolName}; set histogram to a single bin.") 

665 else: 

666 weights = np.ones_like(values) / np.histogram(values, bins=nBinsHist)[0].max() 

667 

668 if not plotConfig.publicationStyle: 

669 # Compute full-tract histogram and get its bins. 

670 # NOTE: `exposure_time` histograms are quantized and look more 

671 # bar-like, so they look better with fewer bins. 

672 bins = ax2.hist( 

673 values, 

674 bins=nBinsHist, 

675 label="Full Tract", 

676 color=histColors[0], 

677 weights=weights, 

678 alpha=0.7, 

679 )[1] 

680 

681 # Align the histogram (top right panel) with the skyproj plots. 

682 pos1 = spf.ax.get_position() # Top left. 

683 pos4 = sp.ax.get_position() # Bottom right. 

684 cbarWidth = cbar.ax.get_position().height / colorBarAspect 

685 # NOTE: cbarWidth != cbar.ax.get_position().width 

686 ax2.set_position([pos4.x0, pos1.y0, pos4.width + cbarWidth, pos1.height]) 

687 

688 # Overplot the histograms for the zoomed-in plots. 

689 for zoomFactor, zidx, color, linestyle, hatch in zip( 

690 zoomFactors, zoomIdx, histColors[1:], ["solid", "dotted"], ["//", "xxxx"] 

691 ): 

692 weights = np.ones_like(values[zidx]) / np.histogram(values[zidx], bins=bins)[0].max() 

693 histLabel = f"{self.prettyPrintFloat(zoomFactor)}x Zoom" 

694 histValues = ax2.hist( 

695 values[zidx], 

696 bins=bins, 

697 label=histLabel, 

698 color=color, 

699 weights=weights, 

700 histtype="step", 

701 linewidth=2, 

702 linestyle=linestyle, 

703 alpha=0.6, 

704 )[0] 

705 # Fill the area under the step. 

706 ax2.fill_between( 

707 bins[:-1], 

708 histValues, 

709 step="post", 

710 color=color, 

711 alpha=0.2, 

712 hatch=hatch, 

713 label="hidden", 

714 ) 

715 

716 # Set labels and legend. 

717 xlabel = plotInfo["property"] 

718 if plotInfo["unit"] not in ["dimensionless", "N/A"]: 

719 xlabel += f" ({plotInfo['unit']})" 

720 xtext = ax2.set_xlabel(xlabel, labelpad=labelpad) 

721 ytext = ax2.set_ylabel("Normalized Count", labelpad=labelpad) 

722 xtext.set_fontsize(rcparams["axes.labelsize"]) 

723 ytext.set_fontsize(rcparams["axes.labelsize"]) 

724 

725 # Get handles and labels from the axis. 

726 handles, labels = ax2.get_legend_handles_labels() 

727 

728 # Add a legend with custom handler that combines the handle 

729 # pairs for the zoomed-in cases. 

730 handles = [handles[0], (handles[1], handles[2]), (handles[3], handles[4])] 

731 while "hidden" in labels: 

732 labels.remove("hidden") 

733 legend = ax2.legend( 

734 handles, 

735 labels, 

736 handler_map={tuple: CustomHandler()}, 

737 loc="best", 

738 frameon=False, 

739 fontsize=15, 

740 ) 

741 

742 for line, text in zip(handles, legend.get_texts()): 

743 if isinstance(line, tuple): 

744 # Use the first handle to get the color. 

745 line = line[0] 

746 color = line.get_edgecolor() if line.get_facecolor()[-1] == 0 else line.get_facecolor() 

747 text.set_color(color) 

748 

749 # Add extra info to plotInfo. 

750 plotInfo["projection"] = "Gnomonic" 

751 plotInfo["nside"] = mapData.nside_sparse 

752 plotInfo["valid_area"] = mapData.get_valid_area() 

753 

754 # Add useful information to the plot. 

755 self.addPlotInfo(fig, plotInfo, toolName) 

756 style = "" 

757 else: 

758 style = "publication-style " 

759 fig.suptitle( 

760 f"{plotInfo['description']} {plotInfo['operation']} map", 

761 fontsize=18.5, 

762 ha="center", 

763 va="top", 

764 y=0.985, 

765 ) 

766 fig.text( 

767 0.5, 

768 0.925, 

769 f"Tract: {plotInfo['tract']}, Band: {plotInfo['band']}, Coadd: {plotInfo['coaddName']}", 

770 ha="center", 

771 fontsize=15.5, 

772 ) 

773 

774 _LOG.info( 

775 f"Made {style}per-tract property map plot for dataset type '{toolName}', " 

776 f"tract: {plotInfo['tract']}, band: '{plotInfo['band']}'." 

777 ) 

778 

779 return fig 

780 

781 @staticmethod 

782 def prettyPrintFloat(n): 

783 if n.is_integer(): 

784 return str(int(n)) 

785 return str(n) 

786 

787 

788class SurveyWidePropertyMapPlot(PlotAction): 

789 plotName = pexConfig.Field[str](doc="The name for the plotting task.", optional=True) 

790 

791 def __call__( 

792 self, 

793 data: KeyedData, 

794 plotConfig: SurveyWidePropertyMapAnalysisConfig, 

795 plotInfo: Mapping[str, Mapping[str, str] | str | int], 

796 **kwargs, 

797 ) -> Mapping[str, Figure]: