Coverage for python / lsst / obs / base / defineVisits.py: 25%

430 statements  

« prev     ^ index     » next       coverage.py v7.13.2, created at 2026-02-04 01:12 -0800

1# This file is part of obs_base. 

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 <http://www.gnu.org/licenses/>. 

21 

22from __future__ import annotations 

23 

24__all__ = [ 

25 "DefineVisitsConfig", 

26 "DefineVisitsTask", 

27 "GroupExposuresConfig", 

28 "GroupExposuresTask", 

29 "VisitDefinitionData", 

30 "VisitSystem", 

31] 

32 

33import cmath 

34import dataclasses 

35import enum 

36import math 

37import operator 

38from abc import ABCMeta, abstractmethod 

39from collections import defaultdict 

40from collections.abc import Callable, Iterable, Sequence 

41from typing import Any, ClassVar, TypeVar, cast 

42 

43import lsst.geom 

44from lsst.afw.cameraGeom import FOCAL_PLANE, PIXELS 

45from lsst.daf.butler import Butler, DataId, DimensionRecord, Progress, Timespan 

46from lsst.geom import Box2D 

47from lsst.pex.config import Config, Field, makeRegistry, registerConfigurable 

48from lsst.pipe.base import Struct, Task 

49from lsst.sphgeom import ConvexPolygon, Region, UnitVector3d 

50from lsst.utils.introspection import get_full_type_name 

51 

52from ._instrument import Instrument, loadCamera 

53 

54 

55class VisitSystem(enum.Enum): 

56 """Enumeration used to label different visit systems.""" 

57 

58 ONE_TO_ONE = 0 

59 """Each exposure is assigned to its own visit.""" 

60 

61 BY_GROUP_METADATA = 1 

62 """Visit membership is defined by the value of the group dimension or, for 

63 older dimension universes, exposure.group_id.""" 

64 

65 BY_SEQ_START_END = 2 

66 """Visit membership is defined by the values of the ``exposure.day_obs``, 

67 ``exposure.seq_start``, and ``exposure.seq_end`` values. 

68 """ 

69 

70 @classmethod 

71 def all(cls) -> frozenset[VisitSystem]: 

72 """Return a `frozenset` containing all members.""" 

73 return frozenset(cls.__members__.values()) 

74 

75 @classmethod 

76 def from_name(cls, external_name: str) -> VisitSystem: 

77 """Construct the enumeration from given name.""" 

78 name = external_name.upper() 

79 name = name.replace("-", "_") 

80 try: 

81 return cls.__members__[name] 

82 except KeyError: 

83 raise KeyError(f"Visit system named '{external_name}' not known.") from None 

84 

85 @classmethod 

86 def from_names(cls, names: Iterable[str] | None) -> frozenset[VisitSystem]: 

87 """Return a `frozenset` of all the visit systems matching the supplied 

88 names. 

89 

90 Parameters 

91 ---------- 

92 names : iterable of `str`, or `None` 

93 Names of visit systems. Case insensitive. If `None` or empty, all 

94 the visit systems are returned. 

95 

96 Returns 

97 ------- 

98 systems : `frozenset` of `VisitSystem` 

99 The matching visit systems. 

100 """ 

101 if not names: 

102 return cls.all() 

103 

104 return frozenset({cls.from_name(name) for name in names}) 

105 

106 def __str__(self) -> str: 

107 name = self.name.lower() 

108 name = name.replace("_", "-") 

109 return name 

110 

111 

112@dataclasses.dataclass 

113class VisitDefinitionData: 

114 """Struct representing a group of exposures that will be used to define a 

115 visit. 

116 """ 

117 

118 instrument: str 

119 """Name of the instrument this visit will be associated with. 

120 """ 

121 

122 id: int 

123 """Integer ID of the visit. 

124 

125 This must be unique across all visit systems for the instrument. 

126 """ 

127 

128 name: str 

129 """String name for the visit. 

130 

131 This must be unique across all visit systems for the instrument. 

132 """ 

133 

134 visit_systems: set[VisitSystem] 

135 """All the visit systems associated with this visit.""" 

136 

137 exposures: list[DimensionRecord] = dataclasses.field(default_factory=list) 

138 """Dimension records for the exposures that are part of this visit. 

139 """ 

140 

141 

142@dataclasses.dataclass 

143class _VisitRecords: 

144 """Struct containing the dimension records associated with a visit.""" 

145 

146 visit: DimensionRecord 

147 """Record for the 'visit' dimension itself. 

148 """ 

149 

150 visit_definition: list[DimensionRecord] 

151 """Records for 'visit_definition', which relates 'visit' to 'exposure'. 

152 """ 

153 

154 visit_detector_region: list[DimensionRecord] 

155 """Records for 'visit_detector_region', which associates the combination 

156 of a 'visit' and a 'detector' with a region on the sky. 

157 """ 

158 

159 visit_system_membership: list[DimensionRecord] 

160 """Records relating visits to an associated visit system.""" 

161 

162 

163class GroupExposuresConfig(Config): 

164 """Configure exposure grouping.""" 

165 

166 

167class GroupExposuresTask(Task, metaclass=ABCMeta): 

168 """Abstract base class for the subtask of `DefineVisitsTask` that is 

169 responsible for grouping exposures into visits. 

170 

171 Subclasses should be registered with `GroupExposuresTask.registry` to 

172 enable use by `DefineVisitsTask`, and should generally correspond to a 

173 particular 'visit_system' dimension value. They are also responsible for 

174 defining visit IDs and names that are unique across all visit systems in 

175 use by an instrument. 

176 

177 Parameters 

178 ---------- 

179 config : `GroupExposuresConfig` 

180 Configuration information. 

181 **kwargs 

182 Additional keyword arguments forwarded to the `lsst.pipe.baseTask` 

183 constructor. 

184 """ 

185 

186 def __init__(self, config: GroupExposuresConfig, **kwargs: Any): 

187 Task.__init__(self, config=config, **kwargs) 

188 

189 ConfigClass = GroupExposuresConfig 

190 

191 _DefaultName = "groupExposures" 

192 

193 registry = makeRegistry( 

194 doc="Registry of algorithms for grouping exposures into visits.", 

195 configBaseType=GroupExposuresConfig, 

196 ) 

197 

198 @abstractmethod 

199 def find_missing( 

200 self, exposures: list[DimensionRecord], registry: lsst.daf.butler.Registry 

201 ) -> list[DimensionRecord]: 

202 """Determine, if possible, which exposures might be missing. 

203 

204 Parameters 

205 ---------- 

206 exposures : `list` of `lsst.daf.butler.DimensionRecord` 

207 The exposure records to analyze. 

208 registry : `lsst.daf.butler.Registry` 

209 A butler registry that contains these exposure records. 

210 

211 Returns 

212 ------- 

213 missing : `list` of `lsst.daf.butler.DimensionRecord` 

214 Any exposure records present in registry that were related to 

215 the given exposures but were missing from that list and deemed 

216 to be relevant. 

217 

218 Notes 

219 ----- 

220 Only some grouping schemes are able to find missing exposures. It 

221 is acceptable to return an empty list. 

222 """ 

223 raise NotImplementedError() 

224 

225 @abstractmethod 

226 def group_exposures(self, exposures: list[DimensionRecord]) -> dict[Any, list[DimensionRecord]]: 

227 """Group the exposures in a way most natural for this visit definition. 

228 

229 Parameters 

230 ---------- 

231 exposures : `list` of `lsst.daf.butler.DimensionRecord` 

232 The exposure records to group. 

233 

234 Returns 

235 ------- 

236 groups : `dict` [Any, `list` [ `lsst.daf.butler.DimensionRecord` ] ] 

237 Groupings of exposure records. The key type is relevant to the 

238 specific visit definition and could be a string or a tuple. 

239 """ 

240 raise NotImplementedError() 

241 

242 @abstractmethod 

243 def group( 

244 self, exposures: list[DimensionRecord], instrument: Instrument 

245 ) -> Iterable[VisitDefinitionData]: 

246 """Group the given exposures into visits. 

247 

248 Parameters 

249 ---------- 

250 exposures : `list` [ `lsst.daf.butler.DimensionRecord` ] 

251 DimensionRecords (for the 'exposure' dimension) describing the 

252 exposures to group. 

253 instrument : `~lsst.obs.base.Instrument` 

254 Instrument specification that can be used to optionally support 

255 some visit ID definitions. 

256 

257 Returns 

258 ------- 

259 visits : `Iterable` [ `VisitDefinitionData` ] 

260 Structs identifying the visits and the exposures associated with 

261 them. This may be an iterator or a container. 

262 """ 

263 raise NotImplementedError() 

264 

265 def getVisitSystems(self) -> set[VisitSystem]: 

266 """Return identifiers for the 'visit_system' dimension this 

267 algorithm implements. 

268 

269 Returns 

270 ------- 

271 visit_systems : `Set` [`VisitSystem`] 

272 The visit systems used by this algorithm. 

273 """ 

274 raise NotImplementedError() 

275 

276 

277class ComputeVisitRegionsConfig(Config): 

278 """Configure visit region calculations.""" 

279 

280 padding: Field[int] = Field( 

281 dtype=int, 

282 default=250, 

283 doc=( 

284 "Pad raw image bounding boxes with specified number of pixels " 

285 "when calculating their (conservatively large) region on the " 

286 "sky. Note that the config value for pixelMargin of the " 

287 "reference object loaders in meas_algorithms should be <= " 

288 "the value set here." 

289 ), 

290 ) 

291 

292 

293class ComputeVisitRegionsTask(Task, metaclass=ABCMeta): 

294 """Abstract base class for the subtask of `DefineVisitsTask` that is 

295 responsible for extracting spatial regions for visits and visit+detector 

296 combinations. 

297 

298 Subclasses should be registered with `ComputeVisitRegionsTask.registry` to 

299 enable use by `DefineVisitsTask`. 

300 

301 Parameters 

302 ---------- 

303 config : `ComputeVisitRegionsConfig` 

304 Configuration information. 

305 butler : `lsst.daf.butler.Butler` 

306 The butler to use. 

307 **kwargs 

308 Additional keyword arguments forwarded to the `Task` constructor. 

309 """ 

310 

311 def __init__(self, config: ComputeVisitRegionsConfig, *, butler: Butler, **kwargs: Any): 

312 Task.__init__(self, config=config, **kwargs) 

313 self.butler = butler 

314 self.instrumentMap: dict[str, Instrument] = {} 

315 

316 ConfigClass = ComputeVisitRegionsConfig 

317 

318 _DefaultName = "computeVisitRegions" 

319 

320 registry = makeRegistry( 

321 doc="Registry of algorithms for computing on-sky regions for visits and visit+detector combinations.", 

322 configBaseType=ComputeVisitRegionsConfig, 

323 ) 

324 

325 def getInstrument(self, instrumentName: str) -> Instrument: 

326 """Retrieve an `~lsst.obs.base.Instrument` associated with this 

327 instrument name. 

328 

329 Parameters 

330 ---------- 

331 instrumentName : `str` 

332 The name of the instrument. 

333 

334 Returns 

335 ------- 

336 instrument : `~lsst.obs.base.Instrument` 

337 The associated instrument object. 

338 

339 Notes 

340 ----- 

341 The result is cached. 

342 """ 

343 instrument = self.instrumentMap.get(instrumentName) 

344 if instrument is None: 

345 instrument = Instrument.fromName(instrumentName, self.butler.registry) 

346 self.instrumentMap[instrumentName] = instrument 

347 return instrument 

348 

349 @abstractmethod 

350 def compute( 

351 self, 

352 visit: VisitDefinitionData, 

353 *, 

354 collections: Sequence[str] | str | None = None, 

355 ) -> tuple[Region, dict[int, Region]]: 

356 """Compute regions for the given visit and all detectors in that visit. 

357 

358 Parameters 

359 ---------- 

360 visit : `VisitDefinitionData` 

361 Struct describing the visit and the exposures associated with it. 

362 collections : `Sequence` [ `str` ] or `str` or `None` 

363 Collections to be searched for camera geometry, overriding 

364 ``self.butler.collections.defaults``. Can be any of the types 

365 supported by the ``collections`` argument to butler construction. 

366 

367 Returns 

368 ------- 

369 visitRegion : `lsst.sphgeom.Region` 

370 Region for the full visit. 

371 visitDetectorRegions : `dict` [ `int`, `lsst.sphgeom.Region` ] 

372 Dictionary mapping detector ID to the region for that detector. 

373 Should include all detectors in the visit. 

374 """ 

375 raise NotImplementedError() 

376 

377 

378class DefineVisitsConfig(Config): 

379 """Configure visit definition.""" 

380 

381 groupExposures = GroupExposuresTask.registry.makeField( 

382 doc="Algorithm for grouping exposures into visits.", 

383 default="one-to-one-and-by-counter", 

384 ) 

385 computeVisitRegions = ComputeVisitRegionsTask.registry.makeField( 

386 doc="Algorithm from computing visit and visit+detector regions.", 

387 default="single-raw-wcs", 

388 ) 

389 ignoreNonScienceExposures: Field[bool] = Field( 

390 doc=( 

391 "If True, silently ignore input exposures that do not have " 

392 "observation_type=SCIENCE. If False, raise an exception if one " 

393 "encountered." 

394 ), 

395 dtype=bool, 

396 optional=False, 

397 default=True, 

398 ) 

399 updateObsCoreTable: Field[bool] = Field( 

400 doc=( 

401 "If True, update exposure regions in obscore table after visits " 

402 "are defined. If False, do not update obscore table." 

403 ), 

404 dtype=bool, 

405 default=True, 

406 ) 

407 

408 

409class DefineVisitsTask(Task): 

410 """Driver Task for defining visits (and their spatial regions) in Gen3 

411 Butler repositories. 

412 

413 Parameters 

414 ---------- 

415 config : `DefineVisitsConfig` 

416 Configuration for the task. 

417 butler : `~lsst.daf.butler.Butler` 

418 Writeable butler instance. Will be used to read ``camera`` datasets 

419 and insert/sync dimension data. 

420 **kwargs 

421 Additional keyword arguments are forwarded to the `lsst.pipe.base.Task` 

422 constructor. 

423 

424 Notes 

425 ----- 

426 Each instance of `DefineVisitsTask` reads from / writes to the same Butler. 

427 Each invocation of `DefineVisitsTask.run` processes an independent group of 

428 exposures into one or more new visits, all belonging to the same visit 

429 system and instrument. 

430 

431 The actual work of grouping exposures and computing regions is delegated to 

432 pluggable subtasks (`GroupExposuresTask` and `ComputeVisitRegionsTask`), 

433 respectively. The defaults are to create one visit for every exposure, and 

434 to use exactly one (arbitrary) detector-level raw dataset's WCS along with 

435 camera geometry to compute regions for all detectors, but the raw WCS is 

436 recomputed from the ``exposure`` dimension record's rotation angle and 

437 boresight rather than by loading the ``raw.wcs`` dataset directly. Other 

438 implementations can be created and configured for instruments for which 

439 these choices are unsuitable (e.g. because visits and exposures are not 

440 one-to-one, or because ``raw.wcs`` datasets for different detectors may not 

441 be consistent with camera geometry). 

442 

443 Defining the same visit the same way multiple times (e.g. via multiple 

444 invocations of this task on the same exposures, with the same 

445 configuration) is safe, but it may be inefficient, as most of the work must 

446 be done before new visits can be compared to existing visits. 

447 """ 

448 

449 def __init__(self, config: DefineVisitsConfig, *, butler: Butler, **kwargs: Any): 

450 config.validate() # Not a CmdlineTask nor PipelineTask, so have to validate the config here. 

451 super().__init__(config, **kwargs) 

452 self.butler = butler 

453 self.universe = self.butler.dimensions 

454 self.progress = Progress("obs.base.DefineVisitsTask") 

455 self.makeSubtask("groupExposures") 

456 self.makeSubtask("computeVisitRegions", butler=self.butler) 

457 

458 def _reduce_kwargs(self) -> dict: 

459 # Add extra parameters to pickle 

460 return dict(**super()._reduce_kwargs(), butler=self.butler) 

461 

462 ConfigClass: ClassVar[type[Config]] = DefineVisitsConfig 

463 

464 _DefaultName: ClassVar[str] = "defineVisits" 

465 

466 config: DefineVisitsConfig 

467 groupExposures: GroupExposuresTask 

468 computeVisitRegions: ComputeVisitRegionsTask 

469 

470 def _buildVisitRecords( 

471 self, definition: VisitDefinitionData, *, collections: Sequence[str] | str | None = None 

472 ) -> _VisitRecords: 

473 """Build the DimensionRecords associated with a visit. 

474 

475 Parameters 

476 ---------- 

477 definition : `VisitDefinitionData` 

478 Struct with identifiers for the visit and records for its 

479 constituent exposures. 

480 collections : `Sequence` [ `str` ] or `str` or `None` 

481 Collections to be searched for camera geometry, overriding 

482 ``self.butler.collections.defaults``. Can be any of the types 

483 supported by the ``collections`` argument to butler construction. 

484 

485 Returns 

486 ------- 

487 records : `_VisitRecords` 

488 Struct containing DimensionRecords for the visit, including 

489 associated dimension elements. 

490 """ 

491 dimension = self.universe["visit"] 

492 

493 # Some registries support additional items. 

494 supported = {meta.name for meta in dimension.metadata} 

495 

496 # Compute all regions. 

497 visitRegion, visitDetectorRegions = self.computeVisitRegions.compute( 

498 definition, collections=collections 

499 ) 

500 # Aggregate other exposure quantities. 

501 timespan = Timespan( 

502 begin=_reduceOrNone(min, (e.timespan.begin for e in definition.exposures)), 

503 end=_reduceOrNone(max, (e.timespan.end for e in definition.exposures)), 

504 ) 

505 exposure_time = _reduceOrNone(operator.add, (e.exposure_time for e in definition.exposures)) 

506 physical_filter = _reduceOrNone(_value_if_equal, (e.physical_filter for e in definition.exposures)) 

507 target_name = _reduceOrNone(_value_if_equal, (e.target_name for e in definition.exposures)) 

508 science_program = _reduceOrNone(_value_if_equal, (e.science_program for e in definition.exposures)) 

509 

510 # observing day for a visit is defined by the earliest observation 

511 # of the visit 

512 observing_day = _reduceOrNone(min, (e.day_obs for e in definition.exposures)) 

513 observation_reason = _reduceOrNone( 

514 _value_if_equal, (e.observation_reason for e in definition.exposures) 

515 ) 

516 if observation_reason is None: 

517 # Be explicit about there being multiple reasons 

518 observation_reason = "various" 

519 

520 # Use the mean zenith angle as an approximation 

521 zenith_angle = _reduceOrNone(operator.add, (e.zenith_angle for e in definition.exposures)) 

522 if zenith_angle is not None: 

523 zenith_angle /= len(definition.exposures) 

524 

525 # New records that may not be supported. 

526 extras: dict[str, Any] = {} 

527 if "seq_num" in supported: 

528 extras["seq_num"] = _reduceOrNone(min, (e.seq_num for e in definition.exposures)) 

529 if "azimuth" in supported: 

530 # Must take into account 0/360 problem. 

531 extras["azimuth"] = _calc_mean_angle([e.azimuth for e in definition.exposures]) 

532 

533 # visit_system handling changed. This is the logic for visit/exposure 

534 # that has support for seq_start/seq_end. 

535 if "seq_num" in supported: 

536 # Map visit to exposure. 

537 visit_definition = [ 

538 self.universe["visit_definition"].RecordClass( 

539 instrument=definition.instrument, 

540 visit=definition.id, 

541 exposure=exposure.id, 

542 ) 

543 for exposure in definition.exposures 

544 ] 

545 

546 # Map visit to visit system. 

547 visit_system_membership = [] 

548 for visit_system in self.groupExposures.getVisitSystems(): 

549 if visit_system in definition.visit_systems: 

550 record = self.universe["visit_system_membership"].RecordClass( 

551 instrument=definition.instrument, 

552 visit=definition.id, 

553 visit_system=visit_system.value, 

554 ) 

555 visit_system_membership.append(record) 

556 

557 else: 

558 # The old approach can only handle one visit system at a time. 

559 # If we have been configured with multiple options, prefer the 

560 # one-to-one. 

561 visit_systems = self.groupExposures.getVisitSystems() 

562 if len(visit_systems) > 1: 

563 one_to_one = VisitSystem.from_name("one-to-one") 

564 if one_to_one not in visit_systems: 

565 raise ValueError( 

566 f"Multiple visit systems specified ({visit_systems}) for use with old" 

567 " dimension universe but unable to find one-to-one." 

568 ) 

569 visit_system = one_to_one 

570 else: 

571 visit_system = visit_systems.pop() 

572 

573 extras["visit_system"] = visit_system.value 

574 

575 # The old visit_definition included visit system. 

576 visit_definition = [ 

577 self.universe["visit_definition"].RecordClass( 

578 instrument=definition.instrument, 

579 visit=definition.id, 

580 exposure=exposure.id, 

581 visit_system=visit_system.value, 

582 ) 

583 for exposure in definition.exposures 

584 ] 

585 

586 # This concept does not exist in old schema. 

587 visit_system_membership = [] 

588 

589 # Construct the actual DimensionRecords. 

590 return _VisitRecords( 

591 visit=dimension.RecordClass( 

592 instrument=definition.instrument, 

593 id=definition.id, 

594 name=definition.name, 

595 physical_filter=physical_filter, 

596 target_name=target_name, 

597 science_program=science_program, 

598 observation_reason=observation_reason, 

599 day_obs=observing_day, 

600 zenith_angle=zenith_angle, 

601 exposure_time=exposure_time, 

602 timespan=timespan, 

603 region=visitRegion, 

604 # TODO: no seeing value in exposure dimension records, so we 

605 # can't set that here. But there are many other columns that 

606 # both dimensions should probably have as well. 

607 **extras, 

608 ), 

609 visit_definition=visit_definition, 

610 visit_system_membership=visit_system_membership, 

611 visit_detector_region=[ 

612 self.universe["visit_detector_region"].RecordClass( 

613 instrument=definition.instrument, 

614 visit=definition.id, 

615 detector=detectorId, 

616 region=detectorRegion, 

617 ) 

618 for detectorId, detectorRegion in visitDetectorRegions.items() 

619 ], 

620 ) 

621 

622 def run( 

623 self, 

624 dataIds_or_records: Iterable[DataId | DimensionRecord], 

625 *, 

626 collections: Sequence[str] | str | None = None, 

627 update_records: bool = False, 

628 incremental: bool = False, 

629 ) -> Struct: 

630 """Add visit definitions to the registry for the given exposures. 

631 

632 Parameters 

633 ---------- 

634 dataIds_or_records : `Iterable` [ `dict` or \ 

635 `~lsst.daf.butler.DataCoordinate` or \ 

636 `~lsst.daf.butler.DimensionRecord` ] 

637 Exposure-level data IDs or explicit exposure records. These must 

638 all correspond to the same instrument, and are expected to be 

639 on-sky science exposures. 

640 collections : `Sequence` [ `str` ] or `str` or `None` 

641 Collections to be searched for camera geometry, overriding 

642 ``self.butler.collections.defaults``. Can be any of the types 

643 supported by the ``collections`` argument to butler construction. 

644 update_records : `bool`, optional 

645 If `True` (`False` is default), update existing ``visit`` records 

646 and ``visit_detector_region`` records. THIS IS AN ADVANCED OPTION 

647 THAT SHOULD ONLY BE USED TO FIX REGIONS AND/OR METADATA THAT ARE 

648 KNOWN TO BE BAD, AND IT CANNOT BE USED TO REMOVE EXPOSURES OR 

649 DETECTORS FROM A VISIT. 

650 incremental : `bool`, optional 

651 If `True` indicate that exposures are being ingested incrementally 

652 and visit definition will be run on partial visits. This will 

653 allow the ``visit`` record to be updated if it already exists, but 

654 (unlike ``update_records=True``) it will only update the 

655 ``visit_detector_region`` records if the ``visit`` record's region 

656 changes. If there is any risk that files are being ingested 

657 incrementally it is critical that this parameter is set to `True` 

658 and not to rely on ``update_records``. 

659 

660 Returns 

661 ------- 

662 result : `lsst.pipe.base.Struct` 

663 Structure with the following attributes (all `int`): 

664 

665 - n_visits: total number of visits defined 

666 - n_skipped: number of visits that were already present left alone 

667 - n_new: number of new visits inserted 

668 - n_fully_updated: number of existing visits fully updated 

669 - n_partially_updated: number of visits with non-geometry updates. 

670 

671 Raises 

672 ------ 

673 lsst.daf.butler.registry.ConflictingDefinitionError 

674 Raised if a visit ID conflict is detected and the existing visit 

675 differs from the new one. 

676 """ 

677 # Normalize, expand, and deduplicate data IDs. 

678 self.log.info("Preprocessing data IDs.") 

679 dimensions = self.universe.conform(["exposure"]) 

680 

681 exposures: list[DimensionRecord] = [] 

682 instruments: set[str] = set() 

683 instrument_cls_name: str | None = None 

684 instrument_record: DimensionRecord | None = None 

685 

686 # Go through the supplied dataset extracting records. 

687 # Check that only a single instrument is being used. 

688 for external in dataIds_or_records: 

689 if isinstance(external, DimensionRecord): 

690 record = external 

691 if str(record.definition) != "exposure": 

692 raise ValueError(f"Can only define visits from exposure records, not {record}.") 

693 else: 

694 data_id = self.butler.registry.expandDataId(external, dimensions=dimensions) 

695 exp_record = data_id.records["exposure"] 

696 assert exp_record is not None, "Guaranteed by expandDataIds call earlier." 

697 record = exp_record 

698 instrument_record = data_id.records["instrument"] 

699 

700 # LSSTCam data can assign ra/dec to flats, and dome-closed 

701 # engineering tests. Do not assign a visit if we know that 

702 # can_see_sky is False. Treat None as True for this test. 

703 can_see_sky = getattr(record, "can_see_sky", True) 

704 if ( 

705 record.tracking_ra is None 

706 or record.tracking_dec is None 

707 or record.sky_angle is None 

708 or can_see_sky is False 

709 ): 

710 if self.config.ignoreNonScienceExposures: 

711 continue 

712 else: 

713 raise RuntimeError( 

714 f"Input exposure {external} has observation_type " 

715 f"{record.observation_type}, but is not on sky." 

716 ) 

717 instrument_name = record.instrument 

718 instruments.add(instrument_name) 

719 exposures.append(record) 

720 if not exposures: 

721 self.log.info("No on-sky exposures found after filtering.") 

722 return Struct(n_visits=0, n_skipped=0, n_new=0, n_partially_updated=0, n_fully_updated=0) 

723 if len(instruments) > 1: 

724 raise RuntimeError( 

725 "All data IDs passed to DefineVisitsTask.run must be " 

726 f"from the same instrument; got {instruments}." 

727 ) 

728 (instrument,) = instruments 

729 

730 # Might need the instrument class for later depending on universe 

731 # and grouping scheme. 

732 if instrument_cls_name is None: 

733 if instrument_record is None: 

734 # We were given a DimensionRecord instead of a DataCoordinate. 

735 

736 instrument_records = self.butler.query_dimension_records( 

737 "instrument", instrument=instrument, limit=1 

738 ) 

739 if len(instrument_records) != 1: 

740 raise RuntimeError( 

741 f"Instrument {instrument} found in dimension record but unknown to butler." 

742 ) 

743 instrument_record = instrument_records[0] 

744 instrument_cls_name = instrument_record.class_name 

745 assert instrument_cls_name is not None, "Instrument must be defined by this point" 

746 instrument_helper = Instrument.from_string(instrument_cls_name) 

747 

748 # Ensure the visit_system our grouping algorithm uses is in the 

749 # registry, if it wasn't already. 

750 visitSystems = self.groupExposures.getVisitSystems() 

751 for visitSystem in visitSystems: 

752 self.log.info("Registering visit_system %d: %s.", visitSystem.value, visitSystem) 

753 self.butler.registry.syncDimensionData( 

754 "visit_system", 

755 {"instrument": instrument, "id": visitSystem.value, "name": str(visitSystem)}, 

756 ) 

757 

758 # In true incremental we will be given the second snap on its 

759 # own on the assumption that the previous snap was already handled. 

760 # For correct grouping we need access to the other exposures in the 

761 # visit. 

762 if incremental: 

763 exposures.extend(self.groupExposures.find_missing(exposures, self.butler.registry)) 

764 

765 # Group exposures into visits, delegating to subtask. 

766 self.log.info("Grouping %d exposure(s) into visits.", len(exposures)) 

767 definitions = list(self.groupExposures.group(exposures, instrument_helper)) 

768 # Iterate over visits, compute regions, and insert dimension data, one 

769 # transaction per visit. If a visit already exists, we skip all other 

770 # inserts. 

771 self.log.info("Computing regions and other metadata for %d visit(s).", len(definitions)) 

772 n_skipped: int = 0 

773 n_new: int = 0 

774 n_fully_updated: int = 0 

775 n_partially_updated: int = 0 

776 for visitDefinition in self.progress.wrap( 

777 definitions, total=len(definitions), desc="Computing regions and inserting visits" 

778 ): 

779 visitRecords = self._buildVisitRecords(visitDefinition, collections=collections) 

780 with self.butler.registry.transaction(): 

781 inserted_or_updated = self.butler.registry.syncDimensionData( 

782 "visit", 

783 visitRecords.visit, 

784 update=(update_records or incremental), 

785 ) 

786 if inserted_or_updated or update_records: 

787 if inserted_or_updated is True: 

788 # This is a new visit, not an update to an existing 

789 # one, so insert visit definition. 

790 # We don't allow visit definitions to change even when 

791 # asked to update, because we'd have to delete the old 

792 # visit_definitions first and also worry about what 

793 # this does to datasets that already use the visit. 

794 self.butler.registry.insertDimensionData( 

795 "visit_definition", *visitRecords.visit_definition 

796 ) 

797 if visitRecords.visit_system_membership: 

798 self.butler.registry.insertDimensionData( 

799 "visit_system_membership", *visitRecords.visit_system_membership 

800 ) 

801 elif incremental and len(visitRecords.visit_definition) > 1: 

802 # The visit record was modified. This could happen 

803 # if a multi-snap visit was redefined with an 

804 # additional snap so play it safe and allow for the 

805 # visit definition to be updated. We use update=False 

806 # here since there should not be any rows updated, 

807 # just additional rows added. update=True does not work 

808 # correctly with multiple records. In incremental mode 

809 # we assume that the caller wants the visit definition 

810 # to be updated and has no worries about provenance 

811 # with the previous definition. 

812 for definition in visitRecords.visit_definition: 

813 self.butler.registry.syncDimensionData("visit_definition", definition) 

814 if inserted_or_updated is True: 

815 # Insert visit-detector regions if the visit is new. 

816 self.butler.registry.insertDimensionData( 

817 "visit_detector_region", 

818 *visitRecords.visit_detector_region, 

819 replace=False, 

820 ) 

821 self.log.verbose( 

822 "Inserted %s visit_detector_region records for new visit %s.", 

823 len(visitRecords.visit_detector_region), 

824 visitRecords.visit.id, 

825 ) 

826 n_new += 1 

827 # Cast below is because MyPy can't determine that 

828 # inserted_or_updated can only be False if update_records 

829 # is True. 

830 elif update_records or "region" in cast(dict, inserted_or_updated): 

831 # Replace visit-detector regions if we were told to 

832 # update records explicitly, or if the visit region 

833 # changed in an incremental=True update. 

834 self.butler.registry.insertDimensionData( 

835 "visit_detector_region", 

836 *visitRecords.visit_detector_region, 

837 replace=True, 

838 ) 

839 self.log.verbose( 

840 "Re-inserted %s visit_detector_region records for updated visit %s.", 

841 len(visitRecords.visit_detector_region), 

842 visitRecords.visit.id, 

843 ) 

844 n_fully_updated += 1 

845 else: 

846 self.log.verbose( 

847 "Updated visit %s without modifying visit_detector_region records.", 

848 visitRecords.visit.id, 

849 ) 

850 n_partially_updated += 1 

851 

852 # Update obscore exposure records with region information 

853 # from corresponding visits. 

854 if self.config.updateObsCoreTable: 

855 if obscore_manager := self.butler.registry.obsCoreTableManager: 

856 obscore_updates: list[tuple[int, int, Region]] = [] 

857 exposure_ids = [rec.exposure for rec in visitRecords.visit_definition] 

858 for record in visitRecords.visit_detector_region: 

859 obscore_updates += [ 

860 (exposure, record.detector, record.region) for exposure in exposure_ids 

861 ] 

862 if obscore_updates: 

863 obscore_manager.update_exposure_regions(instrument, obscore_updates) 

864 else: 

865 self.log.verbose("Skipped already-existing visit %s.", visitRecords.visit.id) 

866 n_skipped += 1 

867 self.log.info( 

868 "Finished writing database records for %d visit(s): %s left unchanged, %s new, " 

869 "%s updated with new detector regions, %s updated without new detector regions.", 

870 len(definitions), 

871 n_skipped, 

872 n_new, 

873 n_fully_updated, 

874 n_partially_updated, 

875 ) 

876 return Struct( 

877 n_visits=len(definitions), 

878 n_skipped=n_skipped, 

879 n_new=n_new, 

880 n_fully_updated=n_fully_updated, 

881 n_partially_updated=n_partially_updated, 

882 ) 

883 

884 

885_T = TypeVar("_T") 

886 

887 

888def _reduceOrNone(func: Callable[[_T, _T], _T | None], iterable: Iterable[_T | None]) -> _T | None: 

889 """Apply a binary function to pairs of elements in an iterable until a 

890 single value is returned, but return `None` if any element is `None` or 

891 there are no elements. 

892 """ 

893 r: _T | None = None 

894 for v in iterable: 

895 if v is None: 

896 return None 

897 if r is None: 

898 r = v 

899 else: 

900 r = func(r, v) 

901 return r 

902 

903 

904def _value_if_equal(a: _T, b: _T) -> _T | None: 

905 """Return either argument if they are equal, or `None` if they are not.""" 

906 return a if a == b else None 

907 

908 

909def _calc_mean_angle(angles: list[float]) -> float: 

910 """Calculate the mean angle, taking into account 0/360 wrapping. 

911 

912 Parameters 

913 ---------- 

914 angles : `list` [`float`] 

915 Angles to average together, in degrees. 

916 

917 Returns 

918 ------- 

919 average : `float` 

920 Average angle in degrees. 

921 """ 

922 # Save on all the math if we only have one value. 

923 if len(angles) == 1: 

924 return angles[0] 

925 

926 # Convert polar coordinates of unit circle to complex values. 

927 # Average the complex values. 

928 # Convert back to a phase angle. 

929 return math.degrees(cmath.phase(sum(cmath.rect(1.0, math.radians(d)) for d in angles) / len(angles))) 

930 

931 

932class _GroupExposuresOneToOneConfig(GroupExposuresConfig): 

933 visitSystemId: Field[int] = Field( 

934 doc="Integer ID of the visit_system implemented by this grouping algorithm.", 

935 dtype=int, 

936 default=0, 

937 deprecated="No longer used. Replaced by enum.", 

938 ) 

939 visitSystemName: Field[str] = Field( 

940 doc="String name of the visit_system implemented by this grouping algorithm.", 

941 dtype=str, 

942 default="one-to-one", 

943 deprecated="No longer used. Replaced by enum.", 

944 ) 

945 

946 

947@registerConfigurable("one-to-one", GroupExposuresTask.registry) 

948class _GroupExposuresOneToOneTask(GroupExposuresTask, metaclass=ABCMeta): 

949 """An exposure grouping algorithm that simply defines one visit for each 

950 exposure, reusing the exposures identifiers for the visit. 

951 """ 

952 

953 ConfigClass = _GroupExposuresOneToOneConfig 

954 

955 def find_missing( 

956 self, exposures: list[DimensionRecord], registry: lsst.daf.butler.Registry 

957 ) -> list[DimensionRecord]: 

958 # By definition no exposures can be missing. 

959 return [] 

960 

961 def group_exposures(self, exposures: list[DimensionRecord]) -> dict[Any, list[DimensionRecord]]: 

962 # No grouping. 

963 return {exposure.id: [exposure] for exposure in exposures} 

964 

965 def group( 

966 self, exposures: list[DimensionRecord], instrument: Instrument 

967 ) -> Iterable[VisitDefinitionData]: 

968 # Docstring inherited from GroupExposuresTask. 

969 visit_systems = {VisitSystem.from_name("one-to-one")} 

970 for exposure in exposures: 

971 yield VisitDefinitionData( 

972 instrument=exposure.instrument, 

973 id=exposure.id, 

974 name=exposure.obs_id, 

975 exposures=[exposure], 

976 visit_systems=visit_systems, 

977 ) 

978 

979 def getVisitSystems(self) -> set[VisitSystem]: 

980 # Docstring inherited from GroupExposuresTask. 

981 return set(VisitSystem.from_names(["one-to-one"])) 

982 

983 

984class _GroupExposuresByGroupMetadataConfig(GroupExposuresConfig): 

985 visitSystemId: Field[int] = Field( 

986 doc="Integer ID of the visit_system implemented by this grouping algorithm.", 

987 dtype=int, 

988 default=1, 

989 deprecated="No longer used. Replaced by enum.", 

990 ) 

991 visitSystemName: Field[str] = Field( 

992 doc="String name of the visit_system implemented by this grouping algorithm.", 

993 dtype=str, 

994 default="by-group-metadata", 

995 deprecated="No longer used. Replaced by enum.", 

996 ) 

997 

998 

999@registerConfigurable("by-group-metadata", GroupExposuresTask.registry) 

1000class _GroupExposuresByGroupMetadataTask(GroupExposuresTask, metaclass=ABCMeta): 

1001 """An exposure grouping algorithm that uses the exposure group. 

1002 

1003 This algorithm uses the ``group`` dimension for modern universes and the 

1004 ``exposure.group_id`` for older universes. 

1005 

1006 This algorithm *assumes* group ID (generally populated from 

1007 `astro_metadata_translator.ObservationInfo.visit_id`) is not just unique, 

1008 but disjoint from all `ObservationInfo.exposure_id` values - if it isn't, 

1009 it will be impossible to ever use both this grouping algorithm and the 

1010 one-to-one algorithm for a particular camera in the same data repository. 

1011 """ 

1012 

1013 ConfigClass = _GroupExposuresByGroupMetadataConfig 

1014 

1015 def find_missing( 

1016 self, exposures: list[DimensionRecord], registry: lsst.daf.butler.Registry 

1017 ) -> list[DimensionRecord]: 

1018 groups = self.group_exposures(exposures) 

1019 # Determine which group implementation we are using. 

1020 if "group" in registry.dimensions["exposure"].implied: 

1021 group_key = "group" 

1022 else: 

1023 group_key = "group_name" 

1024 missing_exposures: list[DimensionRecord] = [] 

1025 for exposures_in_group in groups.values(): 

1026 # We can not tell how many exposures are expected to be in the 

1027 # visit so we have to query every time. 

1028 first = exposures_in_group[0] 

1029 records = set( 

1030 registry.queryDimensionRecords( 

1031 "exposure", 

1032 where=f"exposure.{group_key} = groupnam", 

1033 bind={"groupnam": getattr(first, group_key)}, 

1034 instrument=first.instrument, 

1035 ) 

1036 ) 

1037 records.difference_update(set(exposures_in_group)) 

1038 missing_exposures.extend(list(records)) 

1039 return missing_exposures 

1040 

1041 def group_exposures(self, exposures: list[DimensionRecord]) -> dict[Any, list[DimensionRecord]]: 

1042 groups = defaultdict(list) 

1043 group_key = "group" 

1044 if exposures and hasattr(exposures[0], "group_name"): 

1045 group_key = "group_name" 

1046 for exposure in exposures: 

1047 groups[getattr(exposure, group_key)].append(exposure) 

1048 return groups 

1049 

1050 def group( 

1051 self, exposures: list[DimensionRecord], instrument: Instrument 

1052 ) -> Iterable[VisitDefinitionData]: 

1053 # Docstring inherited from GroupExposuresTask. 

1054 visit_systems = {VisitSystem.from_name("by-group-metadata")} 

1055 groups = self.group_exposures(exposures) 

1056 has_group_dimension: bool | None = None 

1057 for visitName, exposuresInGroup in groups.items(): 

1058 instrument_name = exposuresInGroup[0].instrument 

1059 assert instrument_name == instrument.getName(), "Inconsistency in instrument name" 

1060 visit_ids: set[int] = set() 

1061 if has_group_dimension is None: 

1062 has_group_dimension = hasattr(exposuresInGroup[0], "group") 

1063 if has_group_dimension: 

1064 visit_ids = {instrument.group_name_to_group_id(e.group) for e in exposuresInGroup} 

1065 else: 

1066 visit_ids = {e.group_id for e in exposuresInGroup} 

1067 assert len(visit_ids) == 1, "Grouping by exposure group does not yield consistent group IDs" 

1068 yield VisitDefinitionData( 

1069 instrument=instrument_name, 

1070 id=visit_ids.pop(), 

1071 name=visitName, 

1072 exposures=exposuresInGroup, 

1073 visit_systems=visit_systems, 

1074 ) 

1075 

1076 def getVisitSystems(self) -> set[VisitSystem]: 

1077 # Docstring inherited from GroupExposuresTask. 

1078 return set(VisitSystem.from_names(["by-group-metadata"])) 

1079 

1080 

1081class _GroupExposuresByCounterAndExposuresConfig(GroupExposuresConfig): 

1082 visitSystemId: Field[int] = Field( 

1083 doc="Integer ID of the visit_system implemented by this grouping algorithm.", 

1084 dtype=int, 

1085 default=2, 

1086 deprecated="No longer used. Replaced by enum.", 

1087 ) 

1088 visitSystemName: Field[str] = Field( 

1089 doc="String name of the visit_system implemented by this grouping algorithm.", 

1090 dtype=str, 

1091 default="by-counter-and-exposures", 

1092 deprecated="No longer used. Replaced by enum.", 

1093 ) 

1094 

1095 

1096@registerConfigurable("one-to-one-and-by-counter", GroupExposuresTask.registry) 

1097class _GroupExposuresByCounterAndExposuresTask(GroupExposuresTask, metaclass=ABCMeta): 

1098 """An exposure grouping algorithm that uses the sequence start and 

1099 sequence end metadata to create multi-exposure visits, but also 

1100 creates one-to-one visits. 

1101 

1102 This algorithm uses the exposure.seq_start and 

1103 exposure.seq_end fields to collect related snaps. 

1104 It also groups single exposures. 

1105 """ 

1106 

1107 ConfigClass = _GroupExposuresByCounterAndExposuresConfig 

1108 

1109 def find_missing( 

1110 self, exposures: list[DimensionRecord], registry: lsst.daf.butler.Registry 

1111 ) -> list[DimensionRecord]: 

1112 """Analyze the exposures and return relevant exposures known to 

1113 registry. 

1114 """ 

1115 groups = self.group_exposures(exposures) 

1116 missing_exposures: list[DimensionRecord] = [] 

1117 for exposures_in_group in groups.values(): 

1118 sorted_exposures = sorted(exposures_in_group, key=lambda e: e.seq_num) 

1119 first = sorted_exposures[0] 

1120 

1121 # Only need to look for the seq_nums that we don't already have. 

1122 seq_nums = set(range(first.seq_start, first.seq_end + 1)) 

1123 seq_nums.difference_update({exp.seq_num for exp in sorted_exposures}) 

1124 

1125 if seq_nums: 

1126 # Missing something. Check registry. 

1127 records = list( 

1128 registry.queryDimensionRecords( 

1129 "exposure", 

1130 where="exposure.seq_start = seq_start AND exposure.seq_end = seq_end AND " 

1131 "exposure.seq_num IN (seq_nums)", 

1132 bind={"seq_start": first.seq_start, "seq_end": first.seq_end, "seq_nums": seq_nums}, 

1133 instrument=first.instrument, 

1134 ) 

1135 ) 

1136 missing_exposures.extend(records) 

1137 

1138 return missing_exposures 

1139 

1140 def group_exposures(self, exposures: list[DimensionRecord]) -> dict[Any, list[DimensionRecord]]: 

1141 groups = defaultdict(list) 

1142 for exposure in exposures: 

1143 groups[exposure.day_obs, exposure.seq_start, exposure.seq_end].append(exposure) 

1144 return groups 

1145 

1146 def group( 

1147 self, exposures: list[DimensionRecord], instrument: Instrument 

1148 ) -> Iterable[VisitDefinitionData]: 

1149 # Docstring inherited from GroupExposuresTask. 

1150 system_one_to_one = VisitSystem.from_name("one-to-one") 

1151 system_seq_start_end = VisitSystem.from_name("by-seq-start-end") 

1152 

1153 groups = self.group_exposures(exposures) 

1154 for visit_key, exposures_in_group in groups.items(): 

1155 instrument_name = exposures_in_group[0].instrument 

1156 

1157 # It is possible that the first exposure in a visit has not 

1158 # been ingested. This can be determined and if that is the case 

1159 # we can not reliably define the multi-exposure visit. 

1160 skip_multi = False 

1161 sorted_exposures = sorted(exposures_in_group, key=lambda e: e.seq_num) 

1162 first = sorted_exposures.pop(0) 

1163 if first.seq_num != first.seq_start: 

1164 # Special case seq_num == 0 since that implies that the 

1165 # instrument has no counters and therefore no multi-exposure 

1166 # visits. 

1167 if first.seq_num != 0: 

1168 self.log.warning( 

1169 "First exposure for visit %s is not present. Skipping the multi-snap definition.", 

1170 visit_key, 

1171 ) 

1172 skip_multi = True 

1173 

1174 multi_exposure = False 

1175 if first.seq_start != first.seq_end: 

1176 # This is a multi-exposure visit regardless of the number 

1177 # of exposures present. 

1178 multi_exposure = True 

1179 

1180 # Define the one-to-one visits. 

1181 for exposure in exposures_in_group: 

1182 # Default is to use the exposure ID and name unless 

1183 # this is the first exposure in a multi-exposure visit. 

1184 visit_name = exposure.obs_id 

1185 visit_id = exposure.id 

1186 visit_systems = {system_one_to_one} 

1187 

1188 if not multi_exposure: 

1189 # This is also a by-counter visit. 

1190 # It will use the same visit_name and visit_id. 

1191 visit_systems.add(system_seq_start_end) 

1192 

1193 elif not skip_multi and exposure == first: 

1194 # This is the first legitimate exposure in a multi-exposure 

1195 # visit. It therefore needs a modified visit name and ID 

1196 # so it does not clash with the multi-exposure visit 

1197 # definition. 

1198 visit_name = f"{visit_name}_first" 

1199 visit_id = int(f"9{visit_id}") 

1200 

1201 yield VisitDefinitionData( 

1202 instrument=instrument_name, 

1203 id=visit_id, 

1204 name=visit_name, 

1205 exposures=[exposure], 

1206 visit_systems=visit_systems, 

1207 ) 

1208 

1209 # Multi-exposure visit. 

1210 if not skip_multi and multi_exposure: 

1211 # Define the visit using the first exposure 

1212 visit_name = first.obs_id 

1213 visit_id = first.id 

1214 

1215 yield VisitDefinitionData( 

1216 instrument=instrument_name, 

1217 id=visit_id, 

1218 name=visit_name, 

1219 exposures=exposures_in_group, 

1220 visit_systems={system_seq_start_end}, 

1221 ) 

1222 

1223 def getVisitSystems(self) -> set[VisitSystem]: 

1224 # Docstring inherited from GroupExposuresTask. 

1225 # Using a Config for this is difficult because what this grouping 

1226 # algorithm is doing is using two visit systems. 

1227 # One is using metadata (but not by-group) and the other is the 

1228 # one-to-one. For now hard-code in class. 

1229 return set(VisitSystem.from_names(["one-to-one", "by-seq-start-end"])) 

1230 

1231 

1232class _ComputeVisitRegionsFromSingleRawWcsConfig(ComputeVisitRegionsConfig): 

1233 mergeExposures: Field[bool] = Field( 

1234 doc=( 

1235 "If True, merge per-detector regions over all exposures in a " 

1236 "visit (via convex hull) instead of using the first exposure and " 

1237 "assuming its regions are valid for all others." 

1238 ), 

1239 dtype=bool, 

1240 default=False, 

1241 ) 

1242 detectorId: Field[int | None] = Field( 

1243 doc=( 

1244 "Load the WCS for the detector with this ID. If None, use an " 

1245 "arbitrary detector (the first found in a query of the data " 

1246 "repository for each exposure (or all exposures, if " 

1247 "mergeExposures is True)." 

1248 ), 

1249 dtype=int, 

1250 optional=True, 

1251 default=None, 

1252 ) 

1253 requireVersionedCamera: Field[bool] = Field( 

1254 doc=( 

1255 "If True, raise LookupError if version camera geometry cannot be " 

1256 "loaded for an exposure. If False, use the nominal camera from " 

1257 "the Instrument class instead." 

1258 ), 

1259 dtype=bool, 

1260 optional=False, 

1261 default=False, 

1262 ) 

1263 

1264 

1265@registerConfigurable("single-raw-wcs", ComputeVisitRegionsTask.registry) 

1266class _ComputeVisitRegionsFromSingleRawWcsTask(ComputeVisitRegionsTask): 

1267 """A visit region calculator that uses a single raw WCS (recomputed from 

1268 the ``exposure`` dimension record) and a camera to project the bounding 

1269 boxes of all detectors onto the sky, relating different detectors by their 

1270 positions in focal plane coordinates. 

1271 

1272 Notes 

1273 ----- 

1274 Most instruments should have their raw WCSs determined from a combination 

1275 of boresight angle, rotator angle, and camera geometry, and hence this 

1276 algorithm should produce stable results regardless of which detector the 

1277 raw corresponds to. If this is not the case (e.g. because a per-file FITS 

1278 WCS is used instead), either the ID of the detector should be fixed (see 

1279 the ``detectorId`` config parameter) or a different algorithm used. 

1280 """ 

1281 

1282 ConfigClass = _ComputeVisitRegionsFromSingleRawWcsConfig 

1283 config: _ComputeVisitRegionsFromSingleRawWcsConfig 

1284 

1285 def computeExposureBounds( 

1286 self, exposure: DimensionRecord, *, collections: Any = None 

1287 ) -> dict[int, list[UnitVector3d]]: 

1288 """Compute the lists of unit vectors on the sphere that correspond to 

1289 the sky positions of detector corners. 

1290 

1291 Parameters 

1292 ---------- 

1293 exposure : `DimensionRecord` 

1294 Dimension record for the exposure. 

1295 collections : Any, optional 

1296 Collections to be searched for raws and camera geometry, overriding 

1297 ``self.butler.collections.defaults``. 

1298 Can be any of the types supported by the ``collections`` argument 

1299 to butler construction. 

1300 

1301 Returns 

1302 ------- 

1303 bounds : `dict` 

1304 Dictionary mapping detector ID to a list of unit vectors on the 

1305 sphere representing that detector's corners projected onto the sky. 

1306 """ 

1307 if collections is None: 

1308 collections = list(self.butler.collections.defaults) 

1309 camera, versioned = loadCamera(self.butler, exposure.dataId, collections=collections) 

1310 if not versioned and self.config.requireVersionedCamera: 

1311 raise LookupError(f"No versioned camera found for exposure {exposure.dataId}.") 

1312 

1313 orientation = lsst.geom.Angle(exposure.sky_angle, lsst.geom.degrees) 

1314 radec = lsst.geom.SpherePoint( 

1315 lsst.geom.Angle(exposure.tracking_ra, lsst.geom.degrees), 

1316 lsst.geom.Angle(exposure.tracking_dec, lsst.geom.degrees), 

1317 ) 

1318 

1319 if self.config.detectorId is None: 

1320 detectorId = next(camera.getIdIter()) 

1321 else: 

1322 detectorId = self.config.detectorId 

1323 wcsDetector = camera[detectorId] 

1324 

1325 # Ask the raw formatter to create the relevant WCS 

1326 # This allows flips to be taken into account 

1327 instrument = self.getInstrument(exposure.instrument) 

1328 rawFormatter = instrument.getRawFormatter({"detector": detectorId}) 

1329 

1330 try: 

1331 wcs = rawFormatter.makeRawSkyWcsFromBoresight(radec, orientation, wcsDetector) # type: ignore 

1332 except AttributeError: 

1333 raise TypeError( 

1334 f"Raw formatter is {get_full_type_name(rawFormatter)} but visit" 

1335 " definition requires it to support 'makeRawSkyWcsFromBoresight'" 

1336 ) from None 

1337 

1338 fpToSky = wcsDetector.getTransform(FOCAL_PLANE, PIXELS).then(wcs.getTransform()) 

1339 bounds = {} 

1340 for detector in camera: 

1341 pixelsToSky = detector.getTransform(PIXELS, FOCAL_PLANE).then(fpToSky) 

1342 pixCorners = Box2D(detector.getBBox().dilatedBy(self.config.padding)).getCorners() 

1343 bounds[detector.getId()] = [ 

1344 skyCorner.getVector() for skyCorner in pixelsToSky.applyForward(pixCorners) 

1345 ] 

1346 return bounds 

1347 

1348 def compute( 

1349 self, visit: VisitDefinitionData, *, collections: Any = None 

1350 ) -> tuple[Region, dict[int, Region]]: 

1351 # Docstring inherited from ComputeVisitRegionsTask. 

1352 if self.config.mergeExposures: 

1353 detectorBounds: dict[int, list[UnitVector3d]] = defaultdict(list) 

1354 for exposure in visit.exposures: 

1355 exposureDetectorBounds = self.computeExposureBounds(exposure, collections=collections) 

1356 for detectorId, bounds in exposureDetectorBounds.items(): 

1357 detectorBounds[detectorId].extend(bounds) 

1358 else: 

1359 detectorBounds = self.computeExposureBounds(visit.exposures[0], collections=collections) 

1360 visitBounds = [] 

1361 detectorRegions = {} 

1362 for detectorId, bounds in detectorBounds.items(): 

1363 detectorRegions[detectorId] = ConvexPolygon.convexHull(bounds) 

1364 visitBounds.extend(bounds) 

1365 return ConvexPolygon.convexHull(visitBounds), detectorRegions