LSST Applications 24.1.6,g063fba187b+56b85ce14a,g0f08755f38+df8a265115,g12f32b3c4e+891a09f10d,g1524ad2192+7a5d7b3fbd,g1653933729+a8ce1bb630,g168dd56ebc+a8ce1bb630,g28da252d5a+07cb1400be,g2bbee38e9b+ae03bbfc84,g2bc492864f+ae03bbfc84,g3156d2b45e+6e55a43351,g347aa1857d+ae03bbfc84,g35bb328faa+a8ce1bb630,g3a166c0a6a+ae03bbfc84,g3e281a1b8c+c5dd892a6c,g414038480c+6b9177ef31,g41af890bb2+8f257c4c0b,g781aacb6e4+a8ce1bb630,g7af13505b9+7137b3b17d,g80478fca09+6df6903293,g82479be7b0+091ce1d07f,g858d7b2824+df8a265115,g89c8672015+f4add4ffd5,g9125e01d80+a8ce1bb630,g9726552aa6+414189b318,ga5288a1d22+4a2bca08d7,gacef1a1666+c9a8ff65f4,gb58c049af0+d64f4d3760,gbcfae0f0a0+de1d42d831,gc28159a63d+ae03bbfc84,gcf0d15dbbd+72117bf34e,gda6a2b7d83+72117bf34e,gdaeeff99f8+1711a396fd,ge500cccec5+c8c9c9af63,ge79ae78c31+ae03bbfc84,gf0baf85859+c1f95f4921,gfa517265be+df8a265115,gfa999e8aa5+17cd334064,gfb92a5be7c+df8a265115
LSST Data Management Base Package
Loading...
Searching...
No Matches
FunctionLibrary.cc
Go to the documentation of this file.
1// -*- lsst-c++ -*-
2
3#include <memory>
4
9
10namespace lsst {
11namespace afw {
12namespace math {
13
14namespace {
15
16// Singleton persistence schema for 2-d Gaussians
17struct GaussianFunction2PersistenceHelper {
18 table::Schema schema;
19 table::Key<double> sigma1;
20 table::Key<double> sigma2;
21 table::Key<double> angle;
22
23 static GaussianFunction2PersistenceHelper const& get() {
24 static GaussianFunction2PersistenceHelper instance;
25 return instance;
26 }
27
28 // No copying
29 GaussianFunction2PersistenceHelper(const GaussianFunction2PersistenceHelper&) = delete;
30 GaussianFunction2PersistenceHelper& operator=(const GaussianFunction2PersistenceHelper&) = delete;
31
32 // No moving
33 GaussianFunction2PersistenceHelper(GaussianFunction2PersistenceHelper&&) = delete;
34 GaussianFunction2PersistenceHelper& operator=(GaussianFunction2PersistenceHelper&&) = delete;
35
36private:
37 GaussianFunction2PersistenceHelper()
38 : schema(),
39 sigma1(schema.addField<double>("sigma1", "sigma along axis 1")),
40 sigma2(schema.addField<double>("sigma2", "sigma along axis 2")),
41 angle(schema.addField<double>("angle", "angle of axis 1 in rad (along x=0, y=pi/2)")) {}
42};
43
44// Singleton persistence schema for 2-d circular DoubleGaussians
45struct DoubleGaussianFunction2PersistenceHelper {
46 table::Schema schema;
47 table::Key<double> sigma1;
48 table::Key<double> sigma2;
49 table::Key<double> ampl2;
50
51 static DoubleGaussianFunction2PersistenceHelper const& get() {
52 static DoubleGaussianFunction2PersistenceHelper instance;
53 return instance;
54 }
55
56 // No copying
57 DoubleGaussianFunction2PersistenceHelper(const DoubleGaussianFunction2PersistenceHelper&) = delete;
58 DoubleGaussianFunction2PersistenceHelper& operator=(const DoubleGaussianFunction2PersistenceHelper&) =
59 delete;
60
61 // No moving
62 DoubleGaussianFunction2PersistenceHelper(DoubleGaussianFunction2PersistenceHelper&&) = delete;
63 DoubleGaussianFunction2PersistenceHelper& operator=(DoubleGaussianFunction2PersistenceHelper&&) = delete;
64
65private:
66 DoubleGaussianFunction2PersistenceHelper()
67 : schema(),
68 sigma1(schema.addField<double>("sigma1", "sigma of first Gaussian")),
69 sigma2(schema.addField<double>("sigma2", "sigma of second Gaussian")),
70 ampl2(schema.addField<double>("ampl2", "peak of second Gaussian relative to peak of first")) {}
71};
72
73// Persistence schema for 2-d polynomials; not a singleton because it depends on the order.
74struct PolynomialFunction2PersistenceHelper {
75 table::Schema schema;
76 table::Key<table::Array<double> > coefficients;
77
78 explicit PolynomialFunction2PersistenceHelper(int nCoefficients)
79 : schema(),
80 coefficients(schema.addField<table::Array<double> >(
81 "coefficients",
82 "polynomial coefficients, ordered (x,y) [0,0; 1,0; 0,1; 2,0; 1,1; 0,2; ...]",
83 nCoefficients)) {}
84
85 explicit PolynomialFunction2PersistenceHelper(table::Schema const& schema_)
86 : schema(schema_), coefficients(schema["coefficients"]) {}
87};
88
89// Persistance schema for 2-d Chebyshevs; not a singleton because it depends on the order.
90struct Chebyshev1Function2PersistenceHelper : public PolynomialFunction2PersistenceHelper {
91 table::PointKey<double> min;
92 table::PointKey<double> max;
93
94 explicit Chebyshev1Function2PersistenceHelper(int nCoefficients)
95 : PolynomialFunction2PersistenceHelper(nCoefficients),
96 min(table::PointKey<double>::addFields(schema, "min", "minimum point for function's bbox",
97 "pixel")),
98 max(table::PointKey<double>::addFields(schema, "max", "maximum point for function's bbox",
99 "pixel")) {}
100
101 explicit Chebyshev1Function2PersistenceHelper(table::Schema const& schema_)
102 : PolynomialFunction2PersistenceHelper(schema_), min(schema["min"]), max(schema["max"]) {}
103};
104
105template <typename ReturnT>
106class GaussianFunction2Factory : public table::io::PersistableFactory {
107public:
108 std::shared_ptr<table::io::Persistable> read(InputArchive const& archive,
109 CatalogVector const& catalogs) const override {
110 LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
111 LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
112 GaussianFunction2PersistenceHelper const& keys = GaussianFunction2PersistenceHelper::get();
113 LSST_ARCHIVE_ASSERT(catalogs.front().getSchema().contains(keys.schema));
114 table::BaseRecord const& record = catalogs.front().front();
115 return std::make_shared<GaussianFunction2<ReturnT> >(record.get(keys.sigma1), record.get(keys.sigma2),
116 record.get(keys.angle));
117 }
118
119 GaussianFunction2Factory(std::string const& name) : table::io::PersistableFactory(name) {}
120};
121
122template <typename ReturnT>
123class DoubleGaussianFunction2Factory : public table::io::PersistableFactory {
124public:
125 std::shared_ptr<table::io::Persistable> read(InputArchive const& archive,
126 CatalogVector const& catalogs) const override {
127 LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
128 LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
129 DoubleGaussianFunction2PersistenceHelper const& keys =
130 DoubleGaussianFunction2PersistenceHelper::get();
131 LSST_ARCHIVE_ASSERT(catalogs.front().getSchema().contains(keys.schema));
132 table::BaseRecord const& record = catalogs.front().front();
133 return std::make_shared<DoubleGaussianFunction2<ReturnT> >(
134 record.get(keys.sigma1), record.get(keys.sigma2), record.get(keys.ampl2));
135 }
136
137 DoubleGaussianFunction2Factory(std::string const& name) : table::io::PersistableFactory(name) {}
138};
139
140template <typename ReturnT>
141class PolynomialFunction2Factory : public table::io::PersistableFactory {
142public:
143 std::shared_ptr<table::io::Persistable> read(InputArchive const& archive,
144 CatalogVector const& catalogs) const override {
145 LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
146 LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
147 PolynomialFunction2PersistenceHelper const keys(catalogs.front().getSchema());
148 return std::make_shared<PolynomialFunction2<ReturnT> >(
149 keys.coefficients.extractVector(catalogs.front().front()));
150 }
151
152 PolynomialFunction2Factory(std::string const& name) : table::io::PersistableFactory(name) {}
153};
154
155template <typename ReturnT>
156class Chebyshev1Function2Factory : public table::io::PersistableFactory {
157public:
158 std::shared_ptr<table::io::Persistable> read(InputArchive const& archive,
159 CatalogVector const& catalogs) const override {
160 LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
161 LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
162 Chebyshev1Function2PersistenceHelper keys(catalogs.front().getSchema());
163 table::BaseRecord const& record = catalogs.front().front();
164 lsst::geom::Box2D bbox(record.get(keys.min), record.get(keys.max));
165 return std::make_shared<Chebyshev1Function2<ReturnT> >(keys.coefficients.extractVector(record), bbox);
166 }
167
168 Chebyshev1Function2Factory(std::string const& name) : table::io::PersistableFactory(name) {}
169};
170
171GaussianFunction2Factory<float> registrationGaussian2F("GaussianFunction2F");
172GaussianFunction2Factory<double> registrationGaussian2D("GaussianFunction2D");
173
174DoubleGaussianFunction2Factory<float> registrationDoubleGaussian2F("DoubleGaussianFunction2F");
175DoubleGaussianFunction2Factory<double> registrationDoubleGaussian2D("DoubleGaussianFunction2D");
176
177PolynomialFunction2Factory<float> registrationPoly2F("PolynomialFunction2F");
178PolynomialFunction2Factory<double> registrationPoly2D("PolynomialFunction2D");
179
180Chebyshev1Function2Factory<float> registrationCheb2F("Chebyshev1Function2F");
181Chebyshev1Function2Factory<double> registrationCheb2D("Chebyshev1Function2D");
182
183template <typename T>
184struct Suffix;
185template <>
186struct Suffix<float> {
187 static std::string get() { return "F"; }
188};
189template <>
190struct Suffix<double> {
191 static std::string get() { return "D"; }
192};
193
194} // namespace
195
196template <typename ReturnT>
198 return "GaussianFunction2" + Suffix<ReturnT>::get();
199}
200
201template <typename ReturnT>
203 return "DoubleGaussianFunction2" + Suffix<ReturnT>::get();
204}
205
206template <typename ReturnT>
208 return "PolynomialFunction2" + Suffix<ReturnT>::get();
209}
210
211template <typename ReturnT>
213 return "Chebyshev1Function2" + Suffix<ReturnT>::get();
214}
215
216template <typename ReturnT>
218 GaussianFunction2PersistenceHelper const& keys = GaussianFunction2PersistenceHelper::get();
219 table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
220 std::shared_ptr<table::BaseRecord> record = catalog.addNew();
221 record->set(keys.sigma1, this->getParameters()[0]);
222 record->set(keys.sigma2, this->getParameters()[1]);
223 record->set(keys.angle, this->getParameters()[2]);
224 handle.saveCatalog(catalog);
225}
226
227template <typename ReturnT>
229 DoubleGaussianFunction2PersistenceHelper const& keys = DoubleGaussianFunction2PersistenceHelper::get();
230 table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
231 std::shared_ptr<table::BaseRecord> record = catalog.addNew();
232 record->set(keys.sigma1, this->getParameters()[0]);
233 record->set(keys.sigma2, this->getParameters()[1]);
234 record->set(keys.ampl2, this->getParameters()[2]);
235 handle.saveCatalog(catalog);
236}
237
238template <typename ReturnT>
240 PolynomialFunction2PersistenceHelper const keys(this->getNParameters());
241 table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
242 keys.coefficients.assignVector(*catalog.addNew(), this->getParameters());
243 handle.saveCatalog(catalog);
244}
245
246template <typename ReturnT>
248 Chebyshev1Function2PersistenceHelper const keys(this->getNParameters());
249 table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
250 std::shared_ptr<table::BaseRecord> record = catalog.addNew();
251 keys.coefficients.assignVector(*record, this->getParameters());
252 lsst::geom::Box2D bbox = getXYRange();
253 record->set(keys.min, bbox.getMin());
254 record->set(keys.max, bbox.getMax());
255 handle.saveCatalog(catalog);
256}
257
258// Explicit instantiation
259#define INSTANTIATE(TYPE) \
260 template class IntegerDeltaFunction1<TYPE>; \
261 template class IntegerDeltaFunction2<TYPE>; \
262 template class GaussianFunction1<TYPE>; \
263 template class GaussianFunction2<TYPE>; \
264 template class DoubleGaussianFunction2<TYPE>; \
265 template class PolynomialFunction1<TYPE>; \
266 template class PolynomialFunction2<TYPE>; \
267 template class Chebyshev1Function1<TYPE>; \
268 template class Chebyshev1Function2<TYPE>; \
269 template class LanczosFunction1<TYPE>; \
270 template class LanczosFunction2<TYPE>;
271
272INSTANTIATE(float);
273INSTANTIATE(double);
274} // namespace math
275} // namespace afw
276} // namespace lsst
AmpInfoBoxKey bbox
Definition Amplifier.cc:117
int min
int max
ndarray::Array< double const, 2, 2 > coefficients
#define INSTANTIATE(FROMSYS, TOSYS)
Definition Detector.cc:509
table::Key< double > angle
table::Key< double > ampl2
table::Key< double > sigma2
table::Key< double > sigma1
#define LSST_ARCHIVE_ASSERT(EXPR)
An assertion macro used to validate the structure of an InputArchive.
Definition Persistable.h:48
table::Schema schema
Definition python.h:134
std::string getPersistenceName() const override
Return the unique name used to persist this object and look up its factory.
void write(afw::table::io::OutputArchiveHandle &handle) const override
Write the object to one or more catalogs.
void write(afw::table::io::OutputArchiveHandle &handle) const override
Write the object to one or more catalogs.
std::string getPersistenceName() const override
Return the unique name used to persist this object and look up its factory.
void write(afw::table::io::OutputArchiveHandle &handle) const override
Write the object to one or more catalogs.
std::string getPersistenceName() const override
Return the unique name used to persist this object and look up its factory.
std::string getPersistenceName() const override
Return the unique name used to persist this object and look up its factory.
void write(afw::table::io::OutputArchiveHandle &handle) const override
Write the object to one or more catalogs.
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.
A floating-point coordinate rectangle geometry.
Definition Box.h:413
std::shared_ptr< table::io::Persistable > read(table::io::InputArchive const &archive, table::io::CatalogVector const &catalogs) const override