Merge branch 'tickets/DM-30829'

COPYRIGHT

@@ -1,7 +1,7 @@
-Copyright 2014-2019 University of Washington
-Copyright 2015-2019 Association of Universities for Research in Astronomy
-Copyright 2014-2019 The Trustees of Princeton University
-Copyright 2018 The Board of Trustees of the Leland Stanford Junior University, through SLAC National Accelerator Laboratory
+Copyright 2014-2021 University of Washington
+Copyright 2015-2021 Association of Universities for Research in Astronomy, Inc. (AURA)
+Copyright 2014-2021 The Trustees of Princeton University
+Copyright 2018,2020 The Board of Trustees of the Leland Stanford Junior University, through SLAC National Accelerator Laboratory
 Copyright 2014-2017 The Regents of the University of California
 Copyright 2016 University of Illinois Board of Trustees
 Copyright 2008-2014 LSST Corporation

include/lsst/meas/algorithms/WarpedPsf.h

@@ -1,11 +1,13 @@
 // -*- lsst-c++ -*-
 
 /*
- * LSST Data Management System
- * Copyright 2008, 2009, 2010 LSST Corporation.
+ * This file is part of meas_algorithms.
  *
- * This product includes software developed by the
- * LSST Project (http://www.lsst.org/).
+ * Developed for the LSST Data Management System.
+ * This product includes software developed by the LSST Project
+ * (https://www.lsst.org).
+ * See the COPYRIGHT file at the top-level directory of this distribution
+ * for details of code ownership.
  *
  * This program is free software: you can redistribute it and/or modify
  * it under the terms of the GNU General Public License as published by
@@ -17,25 +19,25 @@
  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  * GNU General Public License for more details.
  *
- * You should have received a copy of the LSST License Statement and
- * the GNU General Public License along with this program.  If not,
- * see <http://www.lsstcorp.org/LegalNotices/>.
+ * You should have received a copy of the GNU General Public License
+ * along with this program.  If not, see <https://www.gnu.org/licenses/>.
  */
 
+#ifndef LSST_MEAS_ALGORITHMS_WARPEDPSF_H
+#define LSST_MEAS_ALGORITHMS_WARPEDPSF_H
+
 #include "lsst/geom/Box.h"
 #include "lsst/afw/geom/Transform.h"
 #include "lsst/afw/math/warpExposure.h"
+#include "lsst/afw/table/io/Persistable.h"
 #include "lsst/meas/algorithms/ImagePsf.h"
 
-#ifndef LSST_AFW_DETECTION_WARPEDPSF_H
-#define LSST_AFW_DETECTION_WARPEDPSF_H
-
 namespace lsst {
 namespace meas {
 namespace algorithms {
 
 /**
- * @brief A Psf class that maps an arbitrary Psf through a coordinate transformation
+ * A Psf class that maps an arbitrary Psf through a coordinate transformation
  *
  * If K_0(x,x') is the unwarped PSF, and f is the coordinate transform, then the
  * warped PSF is defined by
@@ -47,47 +49,55 @@ namespace algorithms {
  * transformation, since the afw convention is that PSF's are normalized to
  * have integral 1 anyway.
  */
-class WarpedPsf : public ImagePsf {
+class WarpedPsf final : public afw::table::io::PersistableFacade<WarpedPsf>, public ImagePsf {
 public:
     /**
-     * @brief Construct WarpedPsf from unwarped psf and distortion.
+     * Construct WarpedPsf from unwarped psf and distortion.
      *
      * If p is the nominal pixel position, and p' is the true position on the sky, then our
      * convention for the transform is that p' = distortion.applyForward(p)
      */
-    WarpedPsf(CONST_PTR(afw::detection::Psf) undistortedPsf,
-              CONST_PTR(afw::geom::TransformPoint2ToPoint2) distortion,
-              CONST_PTR(afw::math::WarpingControl) control);
-    WarpedPsf(CONST_PTR(afw::detection::Psf) undistortedPsf,
-              CONST_PTR(afw::geom::TransformPoint2ToPoint2) distortion,
+    WarpedPsf(std::shared_ptr<afw::detection::Psf const> undistortedPsf,
+              std::shared_ptr<afw::geom::TransformPoint2ToPoint2 const> distortion,
+              std::shared_ptr<afw::math::WarpingControl const> control);
+    WarpedPsf(std::shared_ptr<afw::detection::Psf const> undistortedPsf,
+              std::shared_ptr<afw::geom::TransformPoint2ToPoint2 const> distortion,
               std::string const& kernelName = "lanczos3", unsigned int cache = 10000);
 
     /**
-     *  @brief Return the average of the positions of the stars that went into this Psf.
+     *  Return the average of the positions of the stars that went into this Psf.
      *
      *  For WarpedPsf, this is just the transform of the undistorted Psf's average position.
      */
-    virtual geom::Point2D getAveragePosition() const;
+    geom::Point2D getAveragePosition() const override;
 
     /// Polymorphic deep copy.  Usually unnecessary, as Psfs are immutable.
-    virtual PTR(afw::detection::Psf) clone() const;
+    std::shared_ptr<afw::detection::Psf> clone() const override;
 
     /// Return a clone with specified kernel dimensions
-    virtual PTR(afw::detection::Psf) resized(int width, int height) const;
+    std::shared_ptr<afw::detection::Psf> resized(int width, int height) const override;
 
-protected:
-    virtual PTR(afw::detection::Psf::Image)
-            doComputeKernelImage(geom::Point2D const& position, afw::image::Color const& color) const;
+    bool isPersistable() const noexcept override {
+        return _undistortedPsf->isPersistable() && _distortion->isPersistable() &&
+               _warpingControl->isPersistable();
+    }
 
 protected:
-    PTR(afw::detection::Psf const) _undistortedPsf;
-    PTR(afw::geom::TransformPoint2ToPoint2 const) _distortion;
+    std::shared_ptr<afw::detection::Psf::Image> doComputeKernelImage(
+            geom::Point2D const& position, afw::image::Color const& color) const override;
+
+    std::string getPersistenceName() const override;
+    std::string getPythonModule() const override;
+    void write(OutputArchiveHandle& handle) const override;
+
+    std::shared_ptr<afw::detection::Psf const> _undistortedPsf;
+    std::shared_ptr<afw::geom::TransformPoint2ToPoint2 const> _distortion;
 
 private:
     void _init();
-    CONST_PTR(afw::math::WarpingControl) _warpingControl;
+    std::shared_ptr<afw::math::WarpingControl const> _warpingControl;
 
-    virtual geom::Box2I doComputeBBox(geom::Point2D const& position, afw::image::Color const& color) const;
+    geom::Box2I doComputeBBox(geom::Point2D const& position, afw::image::Color const& color) const override;
 };
 
 }  // namespace algorithms

python/lsst/meas/algorithms/warpedPsf.cc

@@ -21,18 +21,22 @@
  */
 #include "pybind11/pybind11.h"
 
+#include "lsst/utils/python/PySharedPtr.h"
+
+#include "lsst/meas/algorithms/WarpedPsf.h"
+
 namespace py = pybind11;
 using namespace pybind11::literals;
 
-#include "lsst/meas/algorithms/WarpedPsf.h"
+using lsst::utils::python::PySharedPtr;
 
 namespace lsst {
 namespace meas {
 namespace algorithms {
 namespace {
 
 PYBIND11_MODULE(warpedPsf, mod) {
-    py::class_<WarpedPsf, std::shared_ptr<WarpedPsf>, ImagePsf> clsWarpedPsf(mod, "WarpedPsf");
+    py::class_<WarpedPsf, PySharedPtr<WarpedPsf>, ImagePsf> clsWarpedPsf(mod, "WarpedPsf", py::is_final());
 
     /* Constructors */
     clsWarpedPsf.def(py::init<std::shared_ptr<afw::detection::Psf const>,

src/WarpedPsf.cc

@@ -28,8 +28,23 @@
 #include "lsst/meas/algorithms/WarpedPsf.h"
 #include "lsst/afw/math/warpExposure.h"
 #include "lsst/afw/image/Image.h"
+#include "lsst/afw/table/io/InputArchive.h"
+#include "lsst/afw/table/io/OutputArchive.h"
+#include "lsst/afw/table/io/CatalogVector.h"
+#include "lsst/afw/table/io/Persistable.cc"
 
 namespace lsst {
+namespace afw {
+namespace table {
+namespace io {
+
+template std::shared_ptr<meas::algorithms::WarpedPsf>
+PersistableFacade<meas::algorithms::WarpedPsf>::dynamicCast(std::shared_ptr<Persistable> const &);
+
+}  // namespace io
+}  // namespace table
+}  // namespace afw
+
 namespace meas {
 namespace algorithms {
 
@@ -45,12 +60,12 @@ inline double max4(double a, double b, double c, double d) {
 
 // TODO: make this routine externally callable and more generic using templates
 //  (also useful in e.g. math/offsetImage.cc)
-PTR(afw::detection::Psf::Image) zeroPadImage(afw::detection::Psf::Image const &im, int xPad, int yPad) {
+std::shared_ptr<afw::detection::Psf::Image> zeroPadImage(afw::detection::Psf::Image const &im, int xPad,
+                                                         int yPad) {
     int nx = im.getWidth();
     int ny = im.getHeight();
 
-    PTR(afw::detection::Psf::Image)
-    out = std::make_shared<afw::detection::Psf::Image>(nx + 2 * xPad, ny + 2 * yPad);
+    auto out = std::make_shared<afw::detection::Psf::Image>(nx + 2 * xPad, ny + 2 * yPad);
     out->setXY0(im.getX0() - xPad, im.getY0() - yPad);
 
     geom::Box2I box(geom::Point2I(xPad, yPad), geom::Extent2I(nx, ny));
@@ -110,9 +125,9 @@ geom::Box2I computeBBoxFromTransform(geom::Box2I const bbox, geom::AffineTransfo
  *
  * The input image is assumed zero-padded.
  */
-PTR(afw::detection::Psf::Image)
-warpAffine(afw::detection::Psf::Image const &im, geom::AffineTransform const &srcToDest,
-           afw::math::WarpingControl const &wc) {
+std::shared_ptr<afw::detection::Psf::Image> warpAffine(afw::detection::Psf::Image const &im,
+                                                       geom::AffineTransform const &srcToDest,
+                                                       afw::math::WarpingControl const &wc) {
     std::shared_ptr<afw::geom::TransformPoint2ToPoint2> srcToDestTransform =
             afw::geom::makeTransform(srcToDest);
 
@@ -123,10 +138,10 @@ warpAffine(afw::detection::Psf::Image const &im, geom::AffineTransform const &sr
 
     // allocate output image
     geom::Box2I bbox = computeBBoxFromTransform(im.getBBox(), srcToDest);
-    PTR(afw::detection::Psf::Image) ret = std::make_shared<afw::detection::Psf::Image>(bbox);
+    auto ret = std::make_shared<afw::detection::Psf::Image>(bbox);
 
     // zero-pad input image
-    PTR(afw::detection::Psf::Image) im_padded = zeroPadImage(im, xPad, yPad);
+    std::shared_ptr<afw::detection::Psf::Image> im_padded = zeroPadImage(im, xPad, yPad);
 
     // warp it!
     afw::math::warpImage(*ret, *im_padded, *srcToDestTransform, wc, 0.0);
@@ -135,19 +150,19 @@ warpAffine(afw::detection::Psf::Image const &im, geom::AffineTransform const &sr
 
 }  // namespace
 
-WarpedPsf::WarpedPsf(PTR(afw::detection::Psf const) undistortedPsf,
-                     PTR(afw::geom::TransformPoint2ToPoint2 const) distortion,
-                     CONST_PTR(afw::math::WarpingControl) control)
+WarpedPsf::WarpedPsf(std::shared_ptr<afw::detection::Psf const> undistortedPsf,
+                     std::shared_ptr<afw::geom::TransformPoint2ToPoint2 const> distortion,
+                     std::shared_ptr<afw::math::WarpingControl const> control)
         : ImagePsf(false),
           _undistortedPsf(undistortedPsf),
           _distortion(distortion),
           _warpingControl(control) {
     _init();
 }
 
-WarpedPsf::WarpedPsf(PTR(afw::detection::Psf const) undistortedPsf,
-                     PTR(afw::geom::TransformPoint2ToPoint2 const) distortion, std::string const &kernelName,
-                     unsigned int cache)
+WarpedPsf::WarpedPsf(std::shared_ptr<afw::detection::Psf const> undistortedPsf,
+                     std::shared_ptr<afw::geom::TransformPoint2ToPoint2 const> distortion,
+                     std::string const &kernelName, unsigned int cache)
         : ImagePsf(false),
           _undistortedPsf(undistortedPsf),
           _distortion(distortion),
@@ -173,27 +188,27 @@ geom::Point2D WarpedPsf::getAveragePosition() const {
     return _distortion->applyForward(_undistortedPsf->getAveragePosition());
 }
 
-PTR(afw::detection::Psf) WarpedPsf::clone() const {
+std::shared_ptr<afw::detection::Psf> WarpedPsf::clone() const {
     return std::make_shared<WarpedPsf>(_undistortedPsf->clone(), _distortion, _warpingControl);
 }
 
-PTR(afw::detection::Psf) WarpedPsf::resized(int width, int height) const {
+std::shared_ptr<afw::detection::Psf> WarpedPsf::resized(int width, int height) const {
     // For a given set of requested dimensions and distortion, it is not guaranteed that a
     // _undistortedPsf would exist to manifest those dimensions after distortion
     // Not possible to implement with member data currently in WarpedPsf
     throw LSST_EXCEPT(pex::exceptions::LogicError, "Not Implemented");
 }
 
-PTR(afw::detection::Psf::Image)
-WarpedPsf::doComputeKernelImage(geom::Point2D const &position, afw::image::Color const &color) const {
+std::shared_ptr<afw::detection::Psf::Image> WarpedPsf::doComputeKernelImage(
+        geom::Point2D const &position, afw::image::Color const &color) const {
     geom::AffineTransform t = afw::geom::linearizeTransform(*_distortion->inverted(), position);
     geom::Point2D tp = t(position);
 
-    PTR(Image) im = _undistortedPsf->computeKernelImage(tp, color);
+    std::shared_ptr<Image> im = _undistortedPsf->computeKernelImage(tp, color);
 
     // Go to the warped coordinate system with 'p' at the origin
     auto srcToDest = geom::AffineTransform(t.inverted().getLinear());
-    PTR(afw::detection::Psf::Psf::Image) ret = warpAffine(*im, srcToDest, *_warpingControl);
+    std::shared_ptr<afw::detection::Psf::Psf::Image> ret = warpAffine(*im, srcToDest, *_warpingControl);
 
     double normFactor = 1.0;
     //
@@ -223,6 +238,66 @@ geom::Box2I WarpedPsf::doComputeBBox(geom::Point2D const &position, afw::image::
     return ret;
 }
 
+namespace {
+
+struct PersistenceHelper {
+    afw::table::Schema schema;
+    afw::table::Key<int> psfIndex;
+    afw::table::Key<int> transformIndex;
+    afw::table::Key<int> controlIndex;
+
+    static PersistenceHelper const &get() {
+        static PersistenceHelper const instance;
+        return instance;
+    }
+
+private:
+    PersistenceHelper()
+            : schema(),
+              psfIndex(schema.addField<int>("psfIndex", "archive ID of nested Psf object")),
+              transformIndex(schema.addField<int>("transformIndex", "archive ID of nested Transform object")),
+              controlIndex(
+                      schema.addField<int>("controlIndex", "archive ID of nested WarpingControl object")) {}
+};
+
+std::string _getPersistenceName() { return "WarpedPsf"; }
+
+class : public afw::table::io::PersistableFactory {
+public:
+    virtual std::shared_ptr<afw::table::io::Persistable> read(
+            afw::table::io::InputArchive const &archive,
+            afw::table::io::CatalogVector const &catalogs) const {
+        static PersistenceHelper const &keys = PersistenceHelper::get();
+        LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
+        LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
+        afw::table::BaseRecord const &record = catalogs.front().front();
+        LSST_ARCHIVE_ASSERT(record.getSchema() == keys.schema);
+        return std::make_shared<WarpedPsf>(
+                archive.get<afw::detection::Psf>(record.get(keys.psfIndex)),
+                archive.get<afw::geom::TransformPoint2ToPoint2>(record.get(keys.transformIndex)),
+                archive.get<afw::math::WarpingControl>(record.get(keys.controlIndex)));
+    }
+
+    using afw::table::io::PersistableFactory::PersistableFactory;
+} warpedPsfFactory(_getPersistenceName());
+
+}  // namespace
+
+std::string WarpedPsf::getPersistenceName() const { return _getPersistenceName(); }
+std::string WarpedPsf::getPythonModule() const { return "lsst.meas.algorithms"; }
+
+void WarpedPsf::write(OutputArchiveHandle &handle) const {
+    static PersistenceHelper const &keys = PersistenceHelper::get();
+    afw::table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
+    PTR(afw::table::BaseRecord) record = catalog.addNew();
+
+    record->set(keys.psfIndex, handle.put(_undistortedPsf));
+    record->set(keys.transformIndex, handle.put(_distortion));
+    record->set(keys.controlIndex, handle.put(_warpingControl));
+
+    handle.saveCatalog(catalog);
+}
+
 }  // namespace algorithms
 }  // namespace meas
 }  // namespace lsst