LSSTApplications  18.0.0+106,18.0.0+50,19.0.0,19.0.0+1,19.0.0+10,19.0.0+11,19.0.0+13,19.0.0+17,19.0.0+2,19.0.0-1-g20d9b18+6,19.0.0-1-g425ff20,19.0.0-1-g5549ca4,19.0.0-1-g580fafe+6,19.0.0-1-g6fe20d0+1,19.0.0-1-g7011481+9,19.0.0-1-g8c57eb9+6,19.0.0-1-gb5175dc+11,19.0.0-1-gdc0e4a7+9,19.0.0-1-ge272bc4+6,19.0.0-1-ge3aa853,19.0.0-10-g448f008b,19.0.0-12-g6990b2c,19.0.0-2-g0d9f9cd+11,19.0.0-2-g3d9e4fb2+11,19.0.0-2-g5037de4,19.0.0-2-gb96a1c4+3,19.0.0-2-gd955cfd+15,19.0.0-3-g2d13df8,19.0.0-3-g6f3c7dc,19.0.0-4-g725f80e+11,19.0.0-4-ga671dab3b+1,19.0.0-4-gad373c5+3,19.0.0-5-ga2acb9c+2,19.0.0-5-gfe96e6c+2,w.2020.01
LSSTDataManagementBasePackage
FunctionLibrary.cc
Go to the documentation of this file.
1 // -*- lsst-c++ -*-
2 
3 #include <memory>
4 
10 
11 namespace lsst {
12 namespace afw {
13 namespace math {
14 
15 namespace {
16 
17 // Singleton persistence schema for 2-d Gaussians
18 struct GaussianFunction2PersistenceHelper {
19  table::Schema schema;
20  table::Key<double> sigma1;
21  table::Key<double> sigma2;
22  table::Key<double> angle;
23 
24  static GaussianFunction2PersistenceHelper const& get() {
25  static GaussianFunction2PersistenceHelper instance;
26  return instance;
27  }
28 
29  // No copying
30  GaussianFunction2PersistenceHelper(const GaussianFunction2PersistenceHelper&) = delete;
31  GaussianFunction2PersistenceHelper& operator=(const GaussianFunction2PersistenceHelper&) = delete;
32 
33  // No moving
34  GaussianFunction2PersistenceHelper(GaussianFunction2PersistenceHelper&&) = delete;
35  GaussianFunction2PersistenceHelper& operator=(GaussianFunction2PersistenceHelper&&) = delete;
36 
37 private:
38  GaussianFunction2PersistenceHelper()
39  : schema(),
40  sigma1(schema.addField<double>("sigma1", "sigma along axis 1")),
41  sigma2(schema.addField<double>("sigma2", "sigma along axis 2")),
42  angle(schema.addField<double>("angle", "angle of axis 1 in rad (along x=0, y=pi/2)")) {}
43 };
44 
45 // Singleton persistence schema for 2-d circular DoubleGaussians
46 struct DoubleGaussianFunction2PersistenceHelper {
47  table::Schema schema;
48  table::Key<double> sigma1;
49  table::Key<double> sigma2;
50  table::Key<double> ampl2;
51 
52  static DoubleGaussianFunction2PersistenceHelper const& get() {
53  static DoubleGaussianFunction2PersistenceHelper instance;
54  return instance;
55  }
56 
57  // No copying
58  DoubleGaussianFunction2PersistenceHelper(const DoubleGaussianFunction2PersistenceHelper&) = delete;
59  DoubleGaussianFunction2PersistenceHelper& operator=(const DoubleGaussianFunction2PersistenceHelper&) =
60  delete;
61 
62  // No moving
63  DoubleGaussianFunction2PersistenceHelper(DoubleGaussianFunction2PersistenceHelper&&) = delete;
64  DoubleGaussianFunction2PersistenceHelper& operator=(DoubleGaussianFunction2PersistenceHelper&&) = delete;
65 
66 private:
67  DoubleGaussianFunction2PersistenceHelper()
68  : schema(),
69  sigma1(schema.addField<double>("sigma1", "sigma of first Gaussian")),
70  sigma2(schema.addField<double>("sigma2", "sigma of second Gaussian")),
71  ampl2(schema.addField<double>("ampl2", "peak of second Gaussian relative to peak of first")) {}
72 };
73 
74 // Persistence schema for 2-d polynomials; not a singleton because it depends on the order.
75 struct PolynomialFunction2PersistenceHelper {
76  table::Schema schema;
77  table::Key<table::Array<double> > coefficients;
78 
79  explicit PolynomialFunction2PersistenceHelper(int nCoefficients)
80  : schema(),
81  coefficients(schema.addField<table::Array<double> >(
82  "coefficients",
83  "polynomial coefficients, ordered (x,y) [0,0; 1,0; 0,1; 2,0; 1,1; 0,2; ...]",
84  nCoefficients)) {}
85 
86  explicit PolynomialFunction2PersistenceHelper(table::Schema const& schema_)
87  : schema(schema_), coefficients(schema["coefficients"]) {}
88 };
89 
90 // Persistance schema for 2-d Chebyshevs; not a singleton because it depends on the order.
91 struct Chebyshev1Function2PersistenceHelper : public PolynomialFunction2PersistenceHelper {
92  table::PointKey<double> min;
93  table::PointKey<double> max;
94 
95  explicit Chebyshev1Function2PersistenceHelper(int nCoefficients)
96  : PolynomialFunction2PersistenceHelper(nCoefficients),
97  min(table::PointKey<double>::addFields(schema, "min", "minimum point for function's bbox",
98  "pixel")),
99  max(table::PointKey<double>::addFields(schema, "max", "maximum point for function's bbox",
100  "pixel")) {}
101 
102  explicit Chebyshev1Function2PersistenceHelper(table::Schema const& schema_)
103  : PolynomialFunction2PersistenceHelper(schema_), min(schema["min"]), max(schema["max"]) {}
104 };
105 
106 template <typename ReturnT>
107 class GaussianFunction2Factory : public table::io::PersistableFactory {
108 public:
109  std::shared_ptr<table::io::Persistable> read(InputArchive const& archive,
110  CatalogVector const& catalogs) const override {
111  LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
112  LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
113  GaussianFunction2PersistenceHelper const& keys = GaussianFunction2PersistenceHelper::get();
114  LSST_ARCHIVE_ASSERT(catalogs.front().getSchema().contains(keys.schema));
115  table::BaseRecord const& record = catalogs.front().front();
116  return std::make_shared<GaussianFunction2<ReturnT> >(record.get(keys.sigma1), record.get(keys.sigma2),
117  record.get(keys.angle));
118  }
119 
120  GaussianFunction2Factory(std::string const& name) : table::io::PersistableFactory(name) {}
121 };
122 
123 template <typename ReturnT>
124 class DoubleGaussianFunction2Factory : public table::io::PersistableFactory {
125 public:
126  std::shared_ptr<table::io::Persistable> read(InputArchive const& archive,
127  CatalogVector const& catalogs) const override {
128  LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
129  LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
130  DoubleGaussianFunction2PersistenceHelper const& keys =
131  DoubleGaussianFunction2PersistenceHelper::get();
132  LSST_ARCHIVE_ASSERT(catalogs.front().getSchema().contains(keys.schema));
133  table::BaseRecord const& record = catalogs.front().front();
134  return std::make_shared<DoubleGaussianFunction2<ReturnT> >(
135  record.get(keys.sigma1), record.get(keys.sigma2), record.get(keys.ampl2));
136  }
137 
138  DoubleGaussianFunction2Factory(std::string const& name) : table::io::PersistableFactory(name) {}
139 };
140 
141 template <typename ReturnT>
142 class PolynomialFunction2Factory : public table::io::PersistableFactory {
143 public:
144  std::shared_ptr<table::io::Persistable> read(InputArchive const& archive,
145  CatalogVector const& catalogs) const override {
146  LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
147  LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
148  PolynomialFunction2PersistenceHelper const keys(catalogs.front().getSchema());
149  return std::make_shared<PolynomialFunction2<ReturnT> >(
150  keys.coefficients.extractVector(catalogs.front().front()));
151  }
152 
153  PolynomialFunction2Factory(std::string const& name) : table::io::PersistableFactory(name) {}
154 };
155 
156 template <typename ReturnT>
157 class Chebyshev1Function2Factory : public table::io::PersistableFactory {
158 public:
159  std::shared_ptr<table::io::Persistable> read(InputArchive const& archive,
160  CatalogVector const& catalogs) const override {
161  LSST_ARCHIVE_ASSERT(catalogs.size() == 1u);
162  LSST_ARCHIVE_ASSERT(catalogs.front().size() == 1u);
163  Chebyshev1Function2PersistenceHelper keys(catalogs.front().getSchema());
164  table::BaseRecord const& record = catalogs.front().front();
165  lsst::geom::Box2D bbox(record.get(keys.min), record.get(keys.max));
166  return std::make_shared<Chebyshev1Function2<ReturnT> >(keys.coefficients.extractVector(record), bbox);
167  }
168 
169  Chebyshev1Function2Factory(std::string const& name) : table::io::PersistableFactory(name) {}
170 };
171 
172 GaussianFunction2Factory<float> registrationGaussian2F("GaussianFunction2F");
173 GaussianFunction2Factory<double> registrationGaussian2D("GaussianFunction2D");
174 
175 DoubleGaussianFunction2Factory<float> registrationDoubleGaussian2F("DoubleGaussianFunction2F");
176 DoubleGaussianFunction2Factory<double> registrationDoubleGaussian2D("DoubleGaussianFunction2D");
177 
178 PolynomialFunction2Factory<float> registrationPoly2F("PolynomialFunction2F");
179 PolynomialFunction2Factory<double> registrationPoly2D("PolynomialFunction2D");
180 
181 Chebyshev1Function2Factory<float> registrationCheb2F("Chebyshev1Function2F");
182 Chebyshev1Function2Factory<double> registrationCheb2D("Chebyshev1Function2D");
183 
184 template <typename T>
185 struct Suffix;
186 template <>
187 struct Suffix<float> {
188  static std::string get() { return "F"; }
189 };
190 template <>
191 struct Suffix<double> {
192  static std::string get() { return "D"; }
193 };
194 
195 } // namespace
196 
197 template <typename ReturnT>
199  return "GaussianFunction2" + Suffix<ReturnT>::get();
200 }
201 
202 template <typename ReturnT>
204  return "DoubleGaussianFunction2" + Suffix<ReturnT>::get();
205 }
206 
207 template <typename ReturnT>
209  return "PolynomialFunction2" + Suffix<ReturnT>::get();
210 }
211 
212 template <typename ReturnT>
214  return "Chebyshev1Function2" + Suffix<ReturnT>::get();
215 }
216 
217 template <typename ReturnT>
219  GaussianFunction2PersistenceHelper const& keys = GaussianFunction2PersistenceHelper::get();
220  table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
221  std::shared_ptr<table::BaseRecord> record = catalog.addNew();
222  record->set(keys.sigma1, this->getParameters()[0]);
223  record->set(keys.sigma2, this->getParameters()[1]);
224  record->set(keys.angle, this->getParameters()[2]);
225  handle.saveCatalog(catalog);
226 }
227 
228 template <typename ReturnT>
230  DoubleGaussianFunction2PersistenceHelper const& keys = DoubleGaussianFunction2PersistenceHelper::get();
231  table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
232  std::shared_ptr<table::BaseRecord> record = catalog.addNew();
233  record->set(keys.sigma1, this->getParameters()[0]);
234  record->set(keys.sigma2, this->getParameters()[1]);
235  record->set(keys.ampl2, this->getParameters()[2]);
236  handle.saveCatalog(catalog);
237 }
238 
239 template <typename ReturnT>
241  PolynomialFunction2PersistenceHelper const keys(this->getNParameters());
242  table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
243  keys.coefficients.assignVector(*catalog.addNew(), this->getParameters());
244  handle.saveCatalog(catalog);
245 }
246 
247 template <typename ReturnT>
249  Chebyshev1Function2PersistenceHelper const keys(this->getNParameters());
250  table::BaseCatalog catalog = handle.makeCatalog(keys.schema);
251  std::shared_ptr<table::BaseRecord> record = catalog.addNew();
252  keys.coefficients.assignVector(*record, this->getParameters());
253  lsst::geom::Box2D bbox = getXYRange();
254  record->set(keys.min, bbox.getMin());
255  record->set(keys.max, bbox.getMax());
256  handle.saveCatalog(catalog);
257 }
258 
259 // Explicit instantiation
260 #define INSTANTIATE(TYPE) \
261  template class IntegerDeltaFunction1<TYPE>; \
262  template class IntegerDeltaFunction2<TYPE>; \
263  template class GaussianFunction1<TYPE>; \
264  template class GaussianFunction2<TYPE>; \
265  template class DoubleGaussianFunction2<TYPE>; \
266  template class PolynomialFunction1<TYPE>; \
267  template class PolynomialFunction2<TYPE>; \
268  template class Chebyshev1Function1<TYPE>; \
269  template class Chebyshev1Function2<TYPE>; \
270  template class LanczosFunction1<TYPE>; \
271  template class LanczosFunction2<TYPE>;
272 
273 INSTANTIATE(float);
274 INSTANTIATE(double);
275 } // namespace math
276 } // namespace afw
277 } // namespace lsst
AmpInfoBoxKey bbox
Definition: Amplifier.cc:117
void write(afw::table::io::OutputArchiveHandle &handle) const override
Write the object to one or more catalogs.
A floating-point coordinate rectangle geometry.
Definition: Box.h:413
An object passed to Persistable::write to allow it to persist itself.
#define LSST_ARCHIVE_ASSERT(EXPR)
An assertion macro used to validate the structure of an InputArchive.
Definition: Persistable.h:48
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.
table::Key< double > angle
std::string getPersistenceName() const override
Return the unique name used to persist this object and look up its factory.
Point2D const getMin() const noexcept
Definition: Box.h:513
table::Key< double > sigma2
STL class.
void write(afw::table::io::OutputArchiveHandle &handle) const override
Write the object to one or more catalogs.
#define INSTANTIATE(TYPE)
A base class for image defects.
table::PointKey< double > min
table::PointKey< double > max
std::string getPersistenceName() const override
Return the unique name used to persist this object and look up its factory.
BaseCatalog makeCatalog(Schema const &schema)
Return a new, empty catalog with the given schema.
Point2D const getMax() const noexcept
Definition: Box.h:517
table::Key< table::Array< double > > coefficients
table::Key< double > ampl2
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 saveCatalog(BaseCatalog const &catalog)
Save a catalog in the archive.
table::Schema schema
table::Key< double > sigma1
std::shared_ptr< RecordT > addNew()
Create a new record, add it to the end of the catalog, and return a pointer to it.
Definition: Catalog.h:485