Skip to content

Commit c28450e

Browse files
authored
Merge pull request #479 from lsst/tickets/DM-55220
DM-55220: add accessors for template provenance to CoaddPsf
2 parents fbb632c + d5a9214 commit c28450e

3 files changed

Lines changed: 84 additions & 33 deletions

File tree

include/lsst/meas/algorithms/CoaddPsf.h

Lines changed: 28 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -72,7 +72,8 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
7272
* @param[in] warpingKernelName Name of warping kernel
7373
* @param[in] cacheSize Warping kernel cache size
7474
*/
75-
explicit CoaddPsf(afw::table::ExposureCatalog const& catalog, afw::geom::SkyWcs const& coaddWcs,
75+
explicit CoaddPsf(afw::table::ExposureCatalog const& catalog,
76+
std::shared_ptr<afw::geom::SkyWcs const> coaddWcs,
7677
std::string const& weightFieldName = "weight",
7778
std::string const& warpingKernelName = "lanczos3", int cacheSize = 10000);
7879

@@ -90,12 +91,12 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
9091
* @param[in] weightFieldName Field name that contains the weight of the exposure in the coadd;
9192
* defaults to "weight".
9293
*/
93-
CoaddPsf(afw::table::ExposureCatalog const& catalog, afw::geom::SkyWcs const& coaddWcs,
94+
CoaddPsf(afw::table::ExposureCatalog const& catalog, std::shared_ptr<afw::geom::SkyWcs const> coaddWcs,
9495
CoaddPsfControl const& ctrl, std::string const& weightFieldName = "weight")
9596
: CoaddPsf(catalog, coaddWcs, weightFieldName, ctrl.warpingKernelName, ctrl.cacheSize) {}
9697

9798
explicit CoaddPsf(afw::table::ExposureCatalog const& catalog, ///< Unpersisted catalog
98-
afw::geom::SkyWcs const& coaddWcs, ///< WCS for the coadd
99+
std::shared_ptr<afw::geom::SkyWcs const> coaddWcs, ///< WCS for the coadd
99100
geom::Point2D const& averagePosition, ///< Default position for accessors
100101
std::string const& warpingKernelName = "lanczos3", ///< Warping kernel name
101102
int cacheSize = 10000 ///< Kernel cache size
@@ -115,7 +116,7 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
115116
geom::Point2D getAveragePosition() const override { return _averagePosition; }
116117

117118
/// Return the Wcs of the coadd (defines the coordinate system of the Psf).
118-
afw::geom::SkyWcs getCoaddWcs() { return _coaddWcs; }
119+
std::shared_ptr<afw::geom::SkyWcs const> getCoaddWcs() const { return _coaddWcs; }
119120

120121
/// Return the number of component Psfs in this CoaddPsf
121122
int getComponentCount() const;
@@ -127,7 +128,7 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
127128
* @returns Corresponding Psf.
128129
* @throws RangeError Index of component is out of range.
129130
*/
130-
std::shared_ptr<afw::detection::Psf const> getPsf(int index);
131+
std::shared_ptr<afw::detection::Psf const> getPsf(int index) const;
131132

132133
/**
133134
* Get the Wcs of the component image at index.
@@ -136,7 +137,7 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
136137
* @returns Corresponding Wcs.
137138
* @throws RangeError Index of component is out of range.
138139
*/
139-
afw::geom::SkyWcs getWcs(int index);
140+
std::shared_ptr<afw::geom::SkyWcs const> getWcs(int index) const;
140141

141142
/**
142143
* Get the weight of the component image at index.
@@ -145,7 +146,7 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
145146
* @returns Corresponding weight.
146147
* @throws RangeError Index of component is out of range.
147148
*/
148-
double getWeight(int index);
149+
double getWeight(int index) const;
149150

150151
/**
151152
* Get the exposure ID of the component image at index.
@@ -154,7 +155,7 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
154155
* @returns Corresponding exposure ID.
155156
* @throws RangeError Index of component is out of range.
156157
*/
157-
afw::table::RecordId getId(int index);
158+
afw::table::RecordId getId(int index) const;
158159

159160
/**
160161
* Get the bounding box (in component image Pixel coordinates) of the component image at index.
@@ -163,7 +164,7 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
163164
* @returns Corresponding bounding box.
164165
* @throws RangeError Index of component is out of range.
165166
*/
166-
geom::Box2I getBBox(int index);
167+
geom::Box2I getBBox(int index) const;
167168

168169
/**
169170
* Get the validPolygon (in component image Pixel coordinates) of the component image at index.
@@ -172,7 +173,21 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
172173
* @returns Corresponding validPolygon.
173174
* @throws RangeError Index of component is out of range.
174175
*/
175-
std::shared_ptr<afw::geom::polygon::Polygon const> getValidPolygon(int index);
176+
std::shared_ptr<afw::geom::polygon::Polygon const> getValidPolygon(int index) const;
177+
178+
/**
179+
* Get the tract ID of the component image at the given index.
180+
*
181+
* @throws lsst::pex::exceptions::NotFoundError There is no tract column.
182+
*/
183+
int getTract(int index) const;
184+
185+
/**
186+
* Get the patch ID of the component image at the given index.
187+
*
188+
* @throws lsst::pex::exceptions::NotFoundError There is no patch column.
189+
*/
190+
int getPatch(int index) const;
176191

177192
/**
178193
* @brief Return true if the CoaddPsf persistable (always true).
@@ -206,8 +221,10 @@ class CoaddPsf : public afw::table::io::PersistableFacade<CoaddPsf>, public Imag
206221

207222
private:
208223
afw::table::ExposureCatalog _catalog;
209-
afw::geom::SkyWcs _coaddWcs;
224+
std::shared_ptr<afw::geom::SkyWcs const> _coaddWcs;
210225
afw::table::Key<double> _weightKey;
226+
afw::table::Key<int> _tractKey;
227+
afw::table::Key<int> _patchKey;
211228
geom::Point2D _averagePosition;
212229
std::string _warpingKernelName; // could be removed if we could get this from _warpingControl (#2949)
213230
std::shared_ptr<afw::math::WarpingControl const> _warpingControl;

python/lsst/meas/algorithms/coaddPsf/coaddPsf.cc

Lines changed: 5 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -49,15 +49,15 @@ void declareCoaddPsf(lsst::cpputils::python::WrapperCollection &wrappers) {
4949

5050
auto clsCoaddPsf = wrappers.wrapType(PyCoaddPsf(wrappers.module, "CoaddPsf"), [](auto &mod, auto &cls) {
5151
/* Constructors */
52-
cls.def(py::init<afw::table::ExposureCatalog const &, afw::geom::SkyWcs const &,
52+
cls.def(py::init<afw::table::ExposureCatalog const &, std::shared_ptr<afw::geom::SkyWcs>,
5353
std::string const &, std::string const &, int>(),
5454
"catalog"_a, "coaddWcs"_a, "weightFieldName"_a = "weight",
5555
"warpingKernelName"_a = "lanczos3", "cacheSize"_a = 10000);
56-
cls.def(py::init<afw::table::ExposureCatalog const &, afw::geom::SkyWcs const &,
56+
cls.def(py::init<afw::table::ExposureCatalog const &, std::shared_ptr<afw::geom::SkyWcs>,
5757
geom::Point2D const &, std::string const &, int>(),
5858
"catalog"_a, "coaddWcs"_a, "averagePosition"_a,
5959
"warpingKernelName"_a = "lanczos3", "cacheSize"_a = 10000);
60-
cls.def(py::init<afw::table::ExposureCatalog const &, afw::geom::SkyWcs const &,
60+
cls.def(py::init<afw::table::ExposureCatalog const &, std::shared_ptr<afw::geom::SkyWcs>,
6161
CoaddPsfControl const &, std::string const &>(),
6262
"catalog"_a, "coaddWcs"_a, "ctrl"_a, "weightFieldName"_a = "weight");
6363

@@ -72,6 +72,8 @@ void declareCoaddPsf(lsst::cpputils::python::WrapperCollection &wrappers) {
7272
cls.def("getId", &CoaddPsf::getId);
7373
cls.def("getBBox", &CoaddPsf::getBBox);
7474
cls.def("getValidPolygon", &CoaddPsf::getValidPolygon);
75+
cls.def("getTract", &CoaddPsf::getTract);
76+
cls.def("getPatch", &CoaddPsf::getPatch);
7577
cls.def("isPersistable", &CoaddPsf::isPersistable);
7678
});
7779

src/CoaddPsf.cc

Lines changed: 51 additions & 19 deletions
Original file line numberDiff line numberDiff line change
@@ -139,7 +139,8 @@ geom::Point2D computeAveragePosition(afw::table::ExposureCatalog const &catalog,
139139

140140
} // namespace
141141

142-
CoaddPsf::CoaddPsf(afw::table::ExposureCatalog const &catalog, afw::geom::SkyWcs const &coaddWcs,
142+
CoaddPsf::CoaddPsf(afw::table::ExposureCatalog const &catalog,
143+
std::shared_ptr<afw::geom::SkyWcs const> coaddWcs,
143144
std::string const &weightFieldName, std::string const &warpingKernelName, int cacheSize)
144145
: _coaddWcs(coaddWcs),
145146
_warpingKernelName(warpingKernelName),
@@ -165,17 +166,29 @@ CoaddPsf::CoaddPsf(afw::table::ExposureCatalog const &catalog, afw::geom::SkyWcs
165166
record->assign(*i, mapper);
166167
_catalog.push_back(record);
167168
}
168-
_averagePosition = computeAveragePosition(_catalog, _coaddWcs, _weightKey);
169+
_averagePosition = computeAveragePosition(_catalog, *_coaddWcs, _weightKey);
169170
}
170171

171-
CoaddPsf::CoaddPsf(afw::table::ExposureCatalog const &catalog, afw::geom::SkyWcs const &coaddWcs,
172+
CoaddPsf::CoaddPsf(afw::table::ExposureCatalog const &catalog,
173+
std::shared_ptr<afw::geom::SkyWcs const> coaddWcs,
172174
geom::Point2D const &averagePosition, std::string const &warpingKernelName, int cacheSize)
173175
: _catalog(catalog),
174176
_coaddWcs(coaddWcs),
175177
_weightKey(_catalog.getSchema()["weight"]),
178+
_tractKey(),
179+
_patchKey(),
176180
_averagePosition(averagePosition),
177181
_warpingKernelName(warpingKernelName),
178-
_warpingControl(new afw::math::WarpingControl(warpingKernelName, "", cacheSize)) {}
182+
_warpingControl(new afw::math::WarpingControl(warpingKernelName, "", cacheSize)) {
183+
try {
184+
_tractKey = _catalog.getSchema()["tract"];
185+
} catch (pex::exceptions::NotFoundError &) {
186+
}
187+
try {
188+
_patchKey = _catalog.getSchema()["patch"];
189+
} catch (pex::exceptions::NotFoundError &) {
190+
}
191+
}
179192

180193
std::shared_ptr<afw::detection::Psf> CoaddPsf::clone() const { return std::make_shared<CoaddPsf>(*this); }
181194

@@ -227,7 +240,7 @@ void addToImage(std::shared_ptr<afw::image::Image<double>> image,
227240
} // anonymous
228241

229242
geom::Box2I CoaddPsf::doComputeBBox(geom::Point2D const &ccdXY, afw::image::Color const &color) const {
230-
afw::table::ExposureCatalog subcat = _catalog.subsetContaining(ccdXY, _coaddWcs, true);
243+
afw::table::ExposureCatalog subcat = _catalog.subsetContaining(ccdXY, *_coaddWcs, true);
231244
if (subcat.empty()) {
232245
throw LSST_EXCEPT(
233246
lsst::afw::detection::InvalidPsfError,
@@ -238,7 +251,7 @@ geom::Box2I CoaddPsf::doComputeBBox(geom::Point2D const &ccdXY, afw::image::Colo
238251
geom::Box2I ret;
239252
for (auto const &exposureRecord : subcat) {
240253
// compute transform from exposure pixels to coadd pixels
241-
auto exposureToCoadd = afw::geom::makeWcsPairTransform(*exposureRecord.getWcs(), _coaddWcs);
254+
auto exposureToCoadd = afw::geom::makeWcsPairTransform(*exposureRecord.getWcs(), *_coaddWcs);
242255
WarpedPsf warpedPsf = WarpedPsf(exposureRecord.getPsf(), exposureToCoadd, _warpingControl);
243256
geom::Box2I componentBBox = warpedPsf.computeBBox(ccdXY, color);
244257
ret.include(componentBBox);
@@ -250,7 +263,7 @@ geom::Box2I CoaddPsf::doComputeBBox(geom::Point2D const &ccdXY, afw::image::Colo
250263
std::shared_ptr<afw::detection::Psf::Image>
251264
CoaddPsf::doComputeKernelImage(geom::Point2D const &ccdXY, afw::image::Color const &color) const {
252265
// Get the subset of exposures which contain our coordinate within their validPolygons.
253-
afw::table::ExposureCatalog subcat = _catalog.subsetContaining(ccdXY, _coaddWcs, true);
266+
afw::table::ExposureCatalog subcat = _catalog.subsetContaining(ccdXY, *_coaddWcs, true);
254267
if (subcat.empty()) {
255268
throw LSST_EXCEPT(
256269
lsst::afw::detection::InvalidPsfError,
@@ -267,7 +280,7 @@ CoaddPsf::doComputeKernelImage(geom::Point2D const &ccdXY, afw::image::Color con
267280

268281
for (auto const &exposureRecord : subcat) {
269282
// compute transform from exposure pixels to coadd pixels
270-
auto exposureToCoadd = afw::geom::makeWcsPairTransform(*exposureRecord.getWcs(), _coaddWcs);
283+
auto exposureToCoadd = afw::geom::makeWcsPairTransform(*exposureRecord.getWcs(), *_coaddWcs);
271284
std::shared_ptr<afw::image::Image<double>> componentImg;
272285
try {
273286
WarpedPsf warpedPsf = WarpedPsf(exposureRecord.getPsf(), exposureToCoadd, _warpingControl);
@@ -295,48 +308,68 @@ CoaddPsf::doComputeKernelImage(geom::Point2D const &ccdXY, afw::image::Color con
295308

296309
int CoaddPsf::getComponentCount() const { return _catalog.size(); }
297310

298-
std::shared_ptr<afw::detection::Psf const> CoaddPsf::getPsf(int index) {
311+
std::shared_ptr<afw::detection::Psf const> CoaddPsf::getPsf(int index) const {
299312
if (index < 0 || index >= getComponentCount()) {
300313
throw LSST_EXCEPT(pex::exceptions::RangeError, "index of CoaddPsf component out of range");
301314
}
302315
return _catalog[index].getPsf();
303316
}
304317

305-
afw::geom::SkyWcs CoaddPsf::getWcs(int index) {
318+
std::shared_ptr<afw::geom::SkyWcs const> CoaddPsf::getWcs(int index) const {
306319
if (index < 0 || index >= getComponentCount()) {
307320
throw LSST_EXCEPT(pex::exceptions::RangeError, "index of CoaddPsf component out of range");
308321
}
309-
return *_catalog[index].getWcs();
322+
return _catalog[index].getWcs();
310323
}
311324

312-
std::shared_ptr<afw::geom::polygon::Polygon const> CoaddPsf::getValidPolygon(int index) {
325+
std::shared_ptr<afw::geom::polygon::Polygon const> CoaddPsf::getValidPolygon(int index) const {
313326
if (index < 0 || index >= getComponentCount()) {
314327
throw LSST_EXCEPT(pex::exceptions::RangeError, "index of CoaddPsf component out of range");
315328
}
316329
return _catalog[index].getValidPolygon();
317330
}
318331

319-
double CoaddPsf::getWeight(int index) {
332+
double CoaddPsf::getWeight(int index) const {
320333
if (index < 0 || index >= getComponentCount()) {
321334
throw LSST_EXCEPT(pex::exceptions::RangeError, "index of CoaddPsf component out of range");
322335
}
323336
return _catalog[index].get(_weightKey);
324337
}
325338

326-
afw::table::RecordId CoaddPsf::getId(int index) {
339+
afw::table::RecordId CoaddPsf::getId(int index) const {
327340
if (index < 0 || index >= getComponentCount()) {
328341
throw LSST_EXCEPT(pex::exceptions::RangeError, "index of CoaddPsf component out of range");
329342
}
330343
return _catalog[index].getId();
331344
}
332345

333-
geom::Box2I CoaddPsf::getBBox(int index) {
346+
geom::Box2I CoaddPsf::getBBox(int index) const {
334347
if (index < 0 || index >= getComponentCount()) {
335348
throw LSST_EXCEPT(pex::exceptions::RangeError, "index of CoaddPsf component out of range");
336349
}
337350
return _catalog[index].getBBox();
338351
}
339352

353+
int CoaddPsf::getTract(int index) const {
354+
if (index < 0 || index >= getComponentCount()) {
355+
throw LSST_EXCEPT(pex::exceptions::RangeError, "index of CoaddPsf component out of range");
356+
}
357+
if (!_tractKey.isValid()) {
358+
throw LSST_EXCEPT(pex::exceptions::NotFoundError, "CoaddPsf has no tract column");
359+
}
360+
return _catalog[index][_tractKey];
361+
}
362+
363+
int CoaddPsf::getPatch(int index) const {
364+
if (index < 0 || index >= getComponentCount()) {
365+
throw LSST_EXCEPT(pex::exceptions::RangeError, "index of CoaddPsf component out of range");
366+
}
367+
if (!_patchKey.isValid()) {
368+
throw LSST_EXCEPT(pex::exceptions::NotFoundError, "CoaddPsf has no patch column");
369+
}
370+
return _catalog[index][_patchKey];
371+
}
372+
340373
// ---------- Persistence -----------------------------------------------------------------------------------
341374

342375
// For persistence of CoaddPsf, we have two catalogs: the first has just one record, and contains
@@ -388,7 +421,7 @@ class CoaddPsf::Factory : public afw::table::io::PersistableFactory {
388421
afw::table::BaseRecord const &record1 = catalogs.front().front();
389422
return std::shared_ptr<CoaddPsf>(
390423
new CoaddPsf(afw::table::ExposureCatalog::readFromArchive(archive, catalogs.back()),
391-
*archive.get<afw::geom::SkyWcs>(record1.get(keys1.coaddWcs)),
424+
archive.get<afw::geom::SkyWcs>(record1.get(keys1.coaddWcs)),
392425
record1.get(keys1.averagePosition), record1.get(keys1.warpingKernelName),
393426
record1.get(keys1.cacheSize)));
394427
}
@@ -413,7 +446,7 @@ class CoaddPsf::Factory : public afw::table::io::PersistableFactory {
413446
} catch (pex::exceptions::NotFoundError &) {
414447
}
415448
auto averagePos = computeAveragePosition(internalCat, *coaddWcs, weightKey);
416-
return std::shared_ptr<CoaddPsf>(new CoaddPsf(internalCat, *coaddWcs, averagePos));
449+
return std::shared_ptr<CoaddPsf>(new CoaddPsf(internalCat, coaddWcs, averagePos));
417450
}
418451

419452
Factory(std::string const &name) : afw::table::io::PersistableFactory(name) {}
@@ -435,8 +468,7 @@ void CoaddPsf::write(OutputArchiveHandle &handle) const {
435468
CoaddPsfPersistenceHelper const &keys1 = CoaddPsfPersistenceHelper::get();
436469
afw::table::BaseCatalog cat1 = handle.makeCatalog(keys1.schema);
437470
std::shared_ptr<afw::table::BaseRecord> record1 = cat1.addNew();
438-
auto coaddWcsPtr = std::make_shared<afw::geom::SkyWcs>(_coaddWcs);
439-
record1->set(keys1.coaddWcs, handle.put(coaddWcsPtr));
471+
record1->set(keys1.coaddWcs, handle.put(_coaddWcs));
440472
record1->set(keys1.cacheSize, _warpingControl->getCacheSize());
441473
record1->set(keys1.averagePosition, _averagePosition);
442474
record1->set(keys1.warpingKernelName, _warpingKernelName);

0 commit comments

Comments
 (0)