LSST Applications g05325017a3+be04ab5bc8,g2079a07aa2+3e9fd84d81,g20cdd03214+600de29eb2,g2305ad1205+5e1b0c67c7,g2bbee38e9b+69e041a28d,g337abbeb29+69e041a28d,g3a166c0a6a+69e041a28d,g4322eb9e3a+9670871960,g4e26c3aaff+9125fe5495,g50ff169b8f+96c6868917,g52b1c1532d+585e252eca,g591dd9f2cf+608c162c81,g5c3423f6d4+600de29eb2,g607f77f49a+600de29eb2,g69c4c69b53+69e041a28d,g858d7b2824+600de29eb2,g99cad8db69+a7f8c65a2d,g9ddcbc5298+9a081db1e4,ga1e77700b3+2cbb763275,ga8d6e2f9c9+7a336fa28a,gae0086650b+585e252eca,gb0e22166c9+0e73c8378f,gb7c702c636+aec2f3b930,gbb8dafda3b+7d1409958d,gc120e1dc64+acd995d0d8,gc28159a63d+69e041a28d,gcdd4ae20e8+5ef4c697da,gcf0d15dbbd+5ef4c697da,gd1535ee943+1ef0950dc9,gdaeeff99f8+f9a426f77a,gddc38dedce+585e252eca,ge79ae78c31+69e041a28d,gf905fbdb4b+471d9f2c03,gf907d663e3+01e3c6626b,gfbcc870c63+b310236976,w.2024.25
LSST Data Management Base Package
Loading...
Searching...
No Matches
GaussianPsf.cc
Go to the documentation of this file.
1// -*- LSST-C++ -*-
2/*
3 * LSST Data Management System
4 * Copyright 2008-2014 LSST Corporation.
5 *
6 * This product includes software developed by the
7 * LSST Project (http://www.lsst.org/).
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 LSST License Statement and
20 * the GNU General Public License along with this program. If not,
21 * see <http://www.lsstcorp.org/LegalNotices/>.
22 */
23
24#include <memory>
25
26#include "ndarray/eigen.h"
27
33
34namespace lsst {
35namespace afw {
36
37template std::shared_ptr<detection::GaussianPsf> table::io::PersistableFacade<
38 detection::GaussianPsf>::dynamicCast(std::shared_ptr<table::io::Persistable> const&);
39
40namespace detection {
41
42namespace {
43
44// Read-only singleton struct containing the schema and keys that a GaussianPsf is mapped
45// to in record persistence.
46struct GaussianPsfPersistenceHelper {
47 afw::table::Schema schema;
48 afw::table::PointKey<int> dimensions;
49 afw::table::Key<double> sigma;
50
51 static GaussianPsfPersistenceHelper const& get() {
52 static GaussianPsfPersistenceHelper instance;
53 return instance;
54 }
55
56 // No copying
57 GaussianPsfPersistenceHelper(const GaussianPsfPersistenceHelper&) = delete;
58 GaussianPsfPersistenceHelper& operator=(const GaussianPsfPersistenceHelper&) = delete;
59
60 // No moving
61 GaussianPsfPersistenceHelper(GaussianPsfPersistenceHelper&&) = delete;
62 GaussianPsfPersistenceHelper& operator=(GaussianPsfPersistenceHelper&&) = delete;
63
64private:
65 GaussianPsfPersistenceHelper()
66 : schema(),
67 dimensions(afw::table::PointKey<int>::addFields(schema, "dimensions",
68 "width/height of realization of Psf", "pixel")),
69 sigma(schema.addField<double>("sigma", "radius of Gaussian", "pixel")) {}
70};
71
72class GaussianPsfFactory : public afw::table::io::PersistableFactory {
73public:
74 std::shared_ptr<afw::table::io::Persistable> read(InputArchive const& archive,
75 CatalogVector const& catalogs) const override {
76 static GaussianPsfPersistenceHelper const& keys = GaussianPsfPersistenceHelper::get();
77 LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
78 LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
79 afw::table::BaseRecord const& record = catalogs.front().front();
80 LSST_ARCHIVE_ASSERT(record.getSchema() == keys.schema);
81 return std::make_shared<GaussianPsf>(record.get(keys.dimensions.getX()),
82 record.get(keys.dimensions.getY()), record.get(keys.sigma));
83 }
84
85 GaussianPsfFactory(std::string const& name) : afw::table::io::PersistableFactory(name) {}
86};
87
88GaussianPsfFactory registration("GaussianPsf");
89
90void checkDimensions(lsst::geom::Extent2I const& dimensions) {
91 if (dimensions.getX() % 2 == 0 || dimensions.getY() % 2 == 2) {
92 throw LSST_EXCEPT(pex::exceptions::InvalidParameterError, "GaussianPsf dimensions must be odd");
93 }
94}
95
96} // namespace
97
98GaussianPsf::GaussianPsf(int width, int height, double sigma)
99 : Psf(true), _dimensions(width, height), _sigma(sigma) {
100 checkDimensions(_dimensions);
101}
102
104 : Psf(true), _dimensions(dimensions), _sigma(sigma) {
105 checkDimensions(_dimensions);
106}
107
108GaussianPsf::GaussianPsf(GaussianPsf const&) = default;
110GaussianPsf::~GaussianPsf() = default;
111
113 return std::make_shared<GaussianPsf>(_dimensions, _sigma);
114}
115
117 return std::make_shared<GaussianPsf>(lsst::geom::Extent2I(width, height), _sigma);
118}
119
120std::string GaussianPsf::getPersistenceName() const { return "GaussianPsf"; }
121
122std::string GaussianPsf::getPythonModule() const { return "lsst.afw.detection"; }
123
125 static GaussianPsfPersistenceHelper const& keys = GaussianPsfPersistenceHelper::get();
126 afw::table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
127 std::shared_ptr<afw::table::BaseRecord> record = catalog.addNew();
128 (*record)[keys.dimensions.getX()] = _dimensions.getX();
129 (*record)[keys.dimensions.getY()] = _dimensions.getY();
130 (*record)[keys.sigma] = getSigma();
131 handle.saveCatalog(catalog);
132}
133
135 image::Color const&) const {
137 Image::Array array = r->getArray();
138 double sum = 0.0;
139 for (int yIndex = 0, y = r->getY0(); yIndex < _dimensions.getY(); ++yIndex, ++y) {
140 Image::Array::Reference row = array[yIndex];
141 for (int xIndex = 0, x = r->getX0(); xIndex < _dimensions.getX(); ++xIndex, ++x) {
142 sum += row[xIndex] = std::exp(-0.5 * (x * x + y * y) / (_sigma * _sigma));
143 }
144 }
145 ndarray::asEigenMatrix(array) /= sum;
146 return r;
147}
148
149double GaussianPsf::doComputeApertureFlux(double radius, lsst::geom::Point2D const& position,
150 image::Color const& color) const {
151 return 1.0 - std::exp(-0.5 * radius * radius / (_sigma * _sigma));
152}
153
155 image::Color const& color) const {
156 return geom::ellipses::Quadrupole(_sigma * _sigma, _sigma * _sigma, 0.0);
157}
158
160 image::Color const& color) const {
161 return lsst::geom::Box2I(lsst::geom::Point2I(-_dimensions / 2),
162 _dimensions); // integer truncation intentional
163}
164} // namespace detection
165} // namespace afw
166} // namespace lsst
#define LSST_EXCEPT(type,...)
Create an exception with a given type.
Definition Exception.h:48
afw::table::Key< double > sigma
int y
Definition SpanSet.cc:48
#define LSST_ARCHIVE_ASSERT(EXPR)
An assertion macro used to validate the structure of an InputArchive.
Definition Persistable.h:48
A circularly symmetric Gaussian Psf class with no spatial variation, intended mostly for testing purp...
Definition GaussianPsf.h:42
std::shared_ptr< afw::detection::Psf > resized(int width, int height) const override
Return a clone with specified kernel dimensions.
std::string getPersistenceName() const override
Return the unique name used to persist this object and look up its factory.
geom::ellipses::Quadrupole doComputeShape(lsst::geom::Point2D const &position, image::Color const &color) const override
double getSigma() const
Return the radius of the Gaussian.
Definition GaussianPsf.h:77
void write(OutputArchiveHandle &handle) const override
Write the object to one or more catalogs.
std::string getPythonModule() const override
Return the fully-qualified Python module that should be imported to guarantee that its factory is reg...
double doComputeApertureFlux(double radius, lsst::geom::Point2D const &position, image::Color const &color) const override
std::shared_ptr< Image > doComputeKernelImage(lsst::geom::Point2D const &position, image::Color const &color) const override
These virtual member functions are private, not protected, because we only want derived classes to im...
lsst::geom::Box2I doComputeBBox(lsst::geom::Point2D const &position, image::Color const &color) const override
std::shared_ptr< afw::detection::Psf > clone() const override
Polymorphic deep copy; should usually be unnecessary because Psfs are immutable.
GaussianPsf(int width, int height, double sigma)
Constructor for a GaussianPsf.
A polymorphic base class for representing an image's Point Spread Function.
Definition Psf.h:82
lsst::geom::Box2I computeBBox(lsst::geom::Point2D position, image::Color color=image::Color()) const
Return the bounding box of the image returned by computeKernelImage()
Definition Psf.cc:127
virtual lsst::geom::Point2D getAveragePosition() const
Return the average position of the stars used to construct the Psf.
Definition Psf.cc:189
image::Image< Pixel > Image
Image type returned by computeImage.
Definition Psf.h:89
An ellipse core with quadrupole moments as parameters.
Definition Quadrupole.h:47
Describe the colour of a source.
Definition Color.h:25
typename ndarray::Array< PixelT, 2, 1 > Array
A mutable ndarray representation of the image.
Definition ImageBase.h:149
An object passed to Persistable::write to allow it to persist itself.
void saveCatalog(BaseCatalog const &catalog)
Save a catalog in the archive.
BaseCatalog makeCatalog(Schema const &schema)
Return a new, empty catalog with the given schema.
An integer coordinate rectangle.
Definition Box.h:55
T exp(T... args)
int row
Definition CR.cc:145
std::shared_ptr< table::io::Persistable > read(table::io::InputArchive const &archive, table::io::CatalogVector const &catalogs) const override