lsst.jointcal  16.0-3-g21cc1d5+6
Gtransfo.h
Go to the documentation of this file.
1 // -*- LSST-C++ -*-
2 #ifndef LSST_JOINTCAL_GTRANSFO_H
3 #define LSST_JOINTCAL_GTRANSFO_H
4 
5 #include <iostream>
6 #include <memory>
7 #include <string>
8 #include <sstream>
9 #include <vector>
10 
11 #include "Eigen/Core"
12 
13 #include "lsst/pex/exceptions.h"
14 #include "lsst/afw/geom/SkyWcs.h"
15 #include "lsst/jointcal/FatPoint.h"
16 #include "lsst/jointcal/Frame.h"
17 
19 
20 namespace lsst {
21 namespace jointcal {
22 
23 class StarMatchList;
24 class Frame;
25 class GtransfoLin;
26 
28 
42 class Gtransfo {
43 public:
45  virtual void apply(const double xIn, const double yIn, double &xOut, double &yOut) const = 0;
46 
48  void apply(Point const &in, Point &out) const { apply(in.x, in.y, out.x, out.y); }
49 
52  Point apply(Point const &in) const {
53  double xout, yout;
54  apply(in.x, in.y, xout, yout);
55  return Point(xout, yout);
56  }
57 
66  Frame apply(Frame const &inputframe, bool inscribed) const;
67 
69  virtual void dump(std::ostream &stream = std::cout) const = 0;
70 
73  dump(s);
74  return s.str();
75  }
76 
78 
82  virtual double fit(StarMatchList const &starMatchList) = 0;
83 
85  void transformStar(FatPoint &in) const { transformPosAndErrors(in, in); }
86 
88  virtual double getJacobian(Point const &point) const { return getJacobian(point.x, point.y); }
89 
91  virtual std::unique_ptr<Gtransfo> clone() const = 0;
92 
110 
112  virtual double getJacobian(const double x, const double y) const;
113 
119  virtual void computeDerivative(Point const &where, GtransfoLin &derivative,
120  const double step = 0.01) const;
121 
123  virtual GtransfoLin linearApproximation(Point const &where, const double step = 0.01) const;
124 
125  virtual void transformPosAndErrors(const FatPoint &in, FatPoint &out) const;
126 
128  virtual void transformErrors(Point const &where, const double *vIn, double *vOut) const;
129 
131 
133  virtual std::unique_ptr<Gtransfo> inverseTransfo(const double precision, const Frame &region) const;
134 
136  void getParams(double *params) const;
137 
139  void offsetParams(Eigen::VectorXd const &delta);
140 
142  virtual double paramRef(const int i) const;
143 
145  virtual double &paramRef(const int i);
146 
149  virtual void paramDerivatives(Point const &where, double *dx, double *dy) const;
150 
152 
154  virtual std::unique_ptr<Gtransfo> roughInverse(const Frame &region) const;
155 
157  virtual int getNpar() const { return 0; }
158 
167  throw std::logic_error("toAstMap is not implemented for this class.");
168  }
169 
170  void write(const std::string &fileName) const;
171 
172  virtual void write(std::ostream &stream) const;
173 
174  virtual ~Gtransfo(){};
175 };
176 
179 
192 
193 /*=============================================================*/
195 
196 class GtransfoIdentity : public Gtransfo {
197 public:
200 
202  void apply(const double xIn, const double yIn, double &xOut, double &yOut) const override {
203  xOut = xIn;
204  yOut = yIn;
205  } // to speed up
206 
207  double fit(StarMatchList const &starMatchList) override {
208  throw pexExcept::TypeError(
209  "GtransfoIdentity is the identity transformation: it cannot be fit to anything.");
210  }
211 
213  std::unique_ptr<Gtransfo> composeAndReduce(Gtransfo const &right) const override { return right.clone(); }
214 
215  void dump(std::ostream &stream = std::cout) const override { stream << "x' = x\ny' = y" << std::endl; }
216 
217  int getNpar() const override { return 0; }
218 
219  std::unique_ptr<Gtransfo> clone() const override {
221  }
222 
223  void computeDerivative(Point const &where, GtransfoLin &derivative,
224  const double step = 0.01) const override;
225 
227  virtual GtransfoLin linearApproximation(Point const &where, const double step = 0.01) const override;
228 
230  std::shared_ptr<ast::Mapping> toAstMap(jointcal::Frame const &domain) const override;
231 
232  void write(std::ostream &s) const override;
233 
234  void read(std::istream &s);
235 
236  // ClassDef(GtransfoIdentity,1)
237 };
238 
246 
248 bool isIntegerShift(const Gtransfo *gtransfo);
249 
250 /*==================== GtransfoPoly =======================*/
251 
253 class GtransfoPoly : public Gtransfo {
254 public:
260  GtransfoPoly(const unsigned order = 1);
261 
263  GtransfoPoly(const Gtransfo *gtransfo, const Frame &frame, unsigned order, unsigned nPoint = 1000);
264 
274  unsigned const order, unsigned const nSteps = 50);
275 
277  void setOrder(const unsigned order);
279  unsigned getOrder() const { return _order; }
280 
281  using Gtransfo::apply; // to unhide Gtransfo::apply(Point const &)
282 
283  void apply(const double xIn, const double yIn, double &xOut, double &yOut) const override;
284 
286  void computeDerivative(Point const &where, GtransfoLin &derivative,
287  const double step = 0.01) const override;
288 
290  virtual void transformPosAndErrors(const FatPoint &in, FatPoint &out) const override;
291 
293  int getNpar() const override { return 2 * _nterms; }
294 
296  void dump(std::ostream &stream = std::cout) const override;
297 
299  double fit(StarMatchList const &starMatchList) override;
300 
302  GtransfoPoly operator*(GtransfoPoly const &right) const;
303 
305  GtransfoPoly operator+(GtransfoPoly const &right) const;
306 
308  GtransfoPoly operator-(GtransfoPoly const &right) const;
309 
310  using Gtransfo::composeAndReduce; // to unhide Gtransfo::composeAndReduce(Gtransfo const &)
313 
314  std::unique_ptr<Gtransfo> clone() const override {
315  return std::unique_ptr<Gtransfo>(new GtransfoPoly(*this));
316  }
317 
319  double coeff(const unsigned powX, const unsigned powY, const unsigned whichCoord) const;
320 
322  double &coeff(const unsigned powX, const unsigned powY, const unsigned whichCoord);
323 
325  double coeffOrZero(const unsigned powX, const unsigned powY, const unsigned whichCoord) const;
326 
327  double determinant() const;
328 
330  double paramRef(const int i) const override;
331 
333  double &paramRef(const int i) override;
334 
337  void paramDerivatives(Point const &where, double *dx, double *dy) const override;
338 
340  std::shared_ptr<ast::Mapping> toAstMap(jointcal::Frame const &domain) const override;
341 
342  void write(std::ostream &s) const override;
343  void read(std::istream &s);
344 
345 private:
346  double computeFit(StarMatchList const &starMatchList, Gtransfo const &InTransfo, const bool UseErrors);
347 
348  unsigned _order; // The highest sum of exponents of the largest monomial.
349  unsigned _nterms; // number of parameters per coordinate
350  std::vector<double> _coeffs; // the actual coefficients
351  // both polynomials in a single vector to speed up allocation and copies
352 
353  /* use std::vector rather than double * to avoid
354  writing copy constructor and "operator =".
355  Vect would work as well but introduces a dependence
356  that can be avoided */
357 
358  /* This routine take a double * for the vector because the array can
359  then be allocated on the execution stack, which speeds thing
360  up. However this uses Variable Length Array (VLA) which is not
361  part of C++, but gcc implements it. */
362  void computeMonomials(double xIn, double yIn, double *monomial) const;
363 
369  ndarray::Array<double, 2, 2> toAstPolyMapCoefficients() const;
370 };
371 
384  double const precision, int const maxOrder = 9,
385  unsigned const nSteps = 50);
386 
388 
389 /*=============================================================*/
391 class GtransfoLin : public GtransfoPoly {
392 public:
393  using GtransfoPoly::apply; // to unhide Gtransfo::apply(Point const &)
394 
397 
399  explicit GtransfoLin(GtransfoPoly const &gtransfoPoly);
400 
402  GtransfoLin operator*(GtransfoLin const &right) const;
403 
405  GtransfoLin invert() const;
406 
407  // useful? double jacobian(const double x, const double y) const { return determinant();}
408 
410  void computeDerivative(Point const &where, GtransfoLin &derivative, const double step = 0.01) const;
412  GtransfoLin linearApproximation(Point const &where, const double step = 0.01) const;
413 
414  // void dump(std::ostream &stream = std::cout) const;
415 
416  // double fit(StarMatchList const &starMatchList);
417 
419  GtransfoLin(const double ox, const double oy, const double aa11, const double aa12, const double aa21,
420  const double aa22);
421 
424 
426 
427  std::unique_ptr<Gtransfo> inverseTransfo(const double precision, const Frame &region) const;
428 
429  double A11() const { return coeff(1, 0, 0); }
430  double A12() const { return coeff(0, 1, 0); }
431  double A21() const { return coeff(1, 0, 1); }
432  double A22() const { return coeff(0, 1, 1); }
433  double Dx() const { return coeff(0, 0, 0); }
434  double Dy() const { return coeff(0, 0, 1); }
435 
436 protected:
437  double &a11() { return coeff(1, 0, 0); }
438  double &a12() { return coeff(0, 1, 0); }
439  double &a21() { return coeff(1, 0, 1); }
440  double &a22() { return coeff(0, 1, 1); }
441  double &dx() { return coeff(0, 0, 0); }
442  double &dy() { return coeff(0, 0, 1); }
443 
444  friend class Gtransfo;
445  friend class GtransfoIdentity; // for Gtransfo::Derivative
446  friend class GtransfoPoly; // // for Gtransfo::Derivative
447 
448 private:
449  void setOrder(const unsigned order); // to hide GtransfoPoly::setOrder
450 };
451 
452 /*=============================================================*/
453 
456 public:
457  using Gtransfo::apply; // to unhide Gtransfo::apply(Point const &)
459  GtransfoLinShift(double ox = 0., double oy = 0.) : GtransfoLin(ox, oy, 1., 0., 0., 1.) {}
460  GtransfoLinShift(Point const &point) : GtransfoLin(point.x, point.y, 1., 0., 0., 1.){};
461  double fit(StarMatchList const &starMatchList);
462 
463  int getNpar() const { return 2; }
464 };
465 
466 /*=============================================================*/
468 class GtransfoLinRot : public GtransfoLin {
469 public:
470  using Gtransfo::apply; // to unhide apply(const Point&)
471 
473  GtransfoLinRot(const double angleRad, const Point *center = nullptr, const double scaleFactor = 1.0);
474  double fit(StarMatchList const &starMatchList);
475 
476  int getNpar() const { return 4; }
477 };
478 
479 /*=============================================================*/
480 
483 public:
484  using Gtransfo::apply; // to unhide apply(const Point&)
486  GtransfoLinScale(const double scale = 1) : GtransfoLin(0.0, 0.0, scale, 0., 0., scale){};
488  GtransfoLinScale(const double scaleX, const double scaleY)
489  : GtransfoLin(0.0, 0.0, scaleX, 0., 0., scaleY){};
490 
491  int getNpar() const { return 2; }
492 };
493 
504 class GtransfoSkyWcs : public Gtransfo {
505 public:
507 
508  using Gtransfo::apply;
509 
510  // Input is x, y pixels; output is ICRS RA, Dec in degrees
511  void apply(const double xIn, const double yIn, double &xOut, double &yOut) const override;
512 
513  void dump(std::ostream &stream = std::cout) const override;
514 
516  double fit(const StarMatchList &starMatchList) override;
517 
518  std::unique_ptr<Gtransfo> clone() const override;
519 
520  std::shared_ptr<afw::geom::SkyWcs> getSkyWcs() const { return _skyWcs; }
521 
522 private:
524 };
525 
526 /*==================WCS's transfo's =====================================*/
527 
528 class BaseTanWcs : public Gtransfo {
529 public:
530  using Gtransfo::apply; // to unhide apply(const Point&)
531 
532  BaseTanWcs(GtransfoLin const &pix2Tan, Point const &tangentPoint,
533  const GtransfoPoly *corrections = nullptr);
534 
535  BaseTanWcs(const BaseTanWcs &original);
536 
537  void operator=(const BaseTanWcs &original);
538 
540  void apply(const double xIn, const double yIn, double &xOut, double &yOut) const;
541 
543  Point getTangentPoint() const;
544 
546  GtransfoLin getLinPart() const;
547 
549  const GtransfoPoly *getCorr() const { return corr.get(); }
550 
552  void setCorrections(std::unique_ptr<GtransfoPoly> corrections);
553 
555  Point getCrPix() const;
556 
559  virtual GtransfoPoly getPix2TangentPlane() const = 0;
560 
562  virtual void pix2TP(double xPixel, double yPixel, double &xTangentPlane, double &yTangentPlane) const = 0;
563 
564  ~BaseTanWcs();
565 
566 protected:
567  GtransfoLin linPix2Tan; // transform from pixels to tangent plane (degrees)
568  // a linear approximation centered at the pixel and sky origins
570  double ra0, dec0; // sky origin (radians)
571  double cos0, sin0; // cos(dec0), sin(dec0)
572 };
573 
574 class TanRaDec2Pix; // the inverse of TanPix2RaDec.
575 
577 class TanPix2RaDec : public BaseTanWcs {
578 public:
579  using Gtransfo::apply; // to unhide apply(const Point&)
582  TanPix2RaDec(GtransfoLin const &pix2Tan, Point const &tangentPoint,
583  const GtransfoPoly *corrections = nullptr);
584 
586  GtransfoPoly getPix2TangentPlane() const;
587 
589  virtual void pix2TP(double xPixel, double yPixel, double &xTangentPlane, double &yTangentPlane) const;
590 
591  TanPix2RaDec();
592 
594  TanPix2RaDec operator*(GtransfoLin const &right) const;
595 
596  using Gtransfo::composeAndReduce; // to unhide Gtransfo::composeAndReduce(Gtransfo const &)
599 
601  TanRaDec2Pix invert() const;
602 
604  std::unique_ptr<Gtransfo> roughInverse(const Frame &region) const;
605 
608  std::unique_ptr<Gtransfo> inverseTransfo(const double precision, const Frame &region) const;
609 
611 
612  void dump(std::ostream &stream) const;
613 
615  double fit(StarMatchList const &starMatchList);
616 };
617 
619 class TanSipPix2RaDec : public BaseTanWcs {
620 public:
623  TanSipPix2RaDec(GtransfoLin const &pix2Tan, Point const &tangentPoint,
624  const GtransfoPoly *corrections = nullptr);
625 
627  GtransfoPoly getPix2TangentPlane() const;
628 
630  virtual void pix2TP(double xPixel, double yPixel, double &xTangentPlane, double &yTangentPlane) const;
631 
632  TanSipPix2RaDec();
633 
636  std::unique_ptr<Gtransfo> inverseTransfo(const double precision, const Frame &region) const;
637 
639 
640  void dump(std::ostream &stream) const;
641 
643  double fit(StarMatchList const &starMatchList);
644 };
645 
647 
653 class TanRaDec2Pix : public Gtransfo {
654 public:
655  using Gtransfo::apply; // to unhide apply(const Point&)
656 
658  TanRaDec2Pix(GtransfoLin const &tan2Pix, Point const &tangentPoint);
659 
661  TanRaDec2Pix();
662 
664  GtransfoLin getLinPart() const;
665 
667  void setTangentPoint(Point const &tangentPoint);
668 
670  Point getTangentPoint() const;
671 
673  void apply(const double xIn, const double yIn, double &xOut, double &yOut) const;
674 
676  void transformPosAndErrors(const FatPoint &in, FatPoint &out) const;
677 
679  TanPix2RaDec invert() const;
680 
682  std::unique_ptr<Gtransfo> roughInverse(const Frame &region) const;
683 
685  std::unique_ptr<Gtransfo> inverseTransfo(const double precision, const Frame &region) const;
686 
687  void dump(std::ostream &stream) const;
688 
690 
691  double fit(StarMatchList const &starMatchList);
692 
693 private:
694  double ra0, dec0; // tangent point (radians)
695  double cos0, sin0;
696  GtransfoLin linTan2Pix; // tangent plane (probably degrees) to pixels
697 };
698 
700 typedef void(GtransfoFun)(const double, const double, double &, double &, const void *);
701 
703 class UserTransfo : public Gtransfo {
704 public:
705  using Gtransfo::apply; // to unhide apply(const Point&)
706 
708  UserTransfo(GtransfoFun &fun, const void *userData);
709 
710  void apply(const double xIn, const double yIn, double &xOut, double &yOut) const;
711 
712  void dump(std::ostream &stream = std::cout) const;
713 
714  double fit(StarMatchList const &starMatchList);
715 
717 
718 private:
719  GtransfoFun *_userFun;
720  const void *_userData;
721 };
722 
727 } // namespace jointcal
728 } // namespace lsst
729 
730 #endif // LSST_JOINTCAL_GTRANSFO_H
virtual std::unique_ptr< Gtransfo > roughInverse(const Frame &region) const
Rough inverse.
Definition: Gtransfo.cc:166
Implements the (forward) SIP distorsion scheme.
Definition: Gtransfo.h:619
implements the linear transformations (6 real coefficients).
Definition: Gtransfo.h:391
A Gtransfo that holds a SkyWcs.
Definition: Gtransfo.h:504
void getParams(double *params) const
params should be at least Npar() long
Definition: Gtransfo.cc:186
GtransfoLin(GtransfoIdentity const &)
Handy converter:
Definition: Gtransfo.h:423
virtual void computeDerivative(Point const &where, GtransfoLin &derivative, const double step=0.01) const
Computes the local Derivative of a transfo, w.r.t.
Definition: Gtransfo.cc:92
std::ostream & operator<<(std::ostream &out, CcdImageKey const &key)
Definition: CcdImage.cc:28
void transformStar(FatPoint &in) const
allows to write MyTransfo(MyStar)
Definition: Gtransfo.h:85
virtual double getJacobian(Point const &point) const
returns the local jacobian.
Definition: Gtransfo.h:88
A point in a plane.
Definition: Point.h:13
void apply(Point const &in, Point &out) const
applies the tranfo to in and writes into out. Is indeed virtual.
Definition: Gtransfo.h:48
Relationship invert(Relationship r)
void offsetParams(Eigen::VectorXd const &delta)
Definition: Gtransfo.cc:191
void() GtransfoFun(const double, const double, double &, double &, const void *)
signature of the user-provided routine that actually does the coordinate transfo for UserTransfo...
Definition: Gtransfo.h:700
the transformation that handles pix to sideral transfos (Gnomonic, possibly with polynomial distortio...
Definition: Gtransfo.h:577
GtransfoLinShift(double ox=0., double oy=0.)
Add ox and oy.
Definition: Gtransfo.h:459
std::unique_ptr< Gtransfo > clone() const override
returns a copy (allocated by new) of the transformation.
Definition: Gtransfo.h:219
virtual int getNpar() const
returns the number of parameters (to compute chi2&#39;s)
Definition: Gtransfo.h:157
GtransfoLinScale(const double scale=1)
Definition: Gtransfo.h:486
void dump(std::ostream &stream=std::cout) const override
dumps the transfo coefficients to stream.
Definition: Gtransfo.h:215
def scale(algorithm, min, max=None, frame=None)
T endl(T... args)
virtual void transformErrors(Point const &where, const double *vIn, double *vOut) const
transform errors (represented as double[3] in order V(xx),V(yy),Cov(xy))
Definition: Gtransfo.cc:133
T right(T... args)
constexpr Angle operator+(Angle a, Angle d) noexcept
Point apply(Point const &in) const
All these apply(..) shadow the virtual one in derived classes, unless one writes "using Gtransfo::app...
Definition: Gtransfo.h:52
void write(const std::string &fileName) const
Definition: Gtransfo.cc:215
virtual double fit(StarMatchList const &starMatchList)=0
fits a transfo to a std::list of Point pairs (p1,p2, the Point fields in StarMatch).
Polynomial transformation class.
Definition: Gtransfo.h:253
std::shared_ptr< GtransfoPoly > inversePolyTransfo(Gtransfo const &forward, Frame const &domain, double const precision, int const maxOrder=9, unsigned const nSteps=50)
Approximate the inverse by a polynomial, to some precision.
Definition: Gtransfo.cc:1072
GtransfoLin normalizeCoordinatesTransfo(const Frame &frame)
Returns the transformation that maps the input frame along both axes to [-1,1].
Definition: Gtransfo.cc:764
STL class.
A Point with uncertainties.
Definition: FatPoint.h:11
int const step
STL class.
double x
coordinate
Definition: Point.h:18
GtransfoLin()
the default constructor constructs the do-nothing transformation.
Definition: Gtransfo.h:396
bool isIntegerShift(const Gtransfo *gtransfo)
Shorthand test to tell if a transfo is a simple integer shift.
Definition: Gtransfo.cc:31
std::unique_ptr< Gtransfo > clone() const
returns a copy (allocated by new) of the transformation.
Definition: Gtransfo.h:425
just here to provide a specialized constructor, and fit.
Definition: Gtransfo.h:468
constexpr Angle operator-(Angle a, Angle d) noexcept
GtransfoLinScale(const double scaleX, const double scaleY)
Definition: Gtransfo.h:488
rectangle with sides parallel to axes.
Definition: Frame.h:15
Class for a simple mapping implementing a generic Gtransfo.
int getNpar() const
total number of parameters
Definition: Gtransfo.h:463
int getNpar() const override
returns the number of parameters (to compute chi2&#39;s)
Definition: Gtransfo.h:217
GtransfoIdentity()
constructor.
Definition: Gtransfo.h:199
constexpr Angle operator*(Angle a, Angle d) noexcept
std::shared_ptr< afw::geom::SkyWcs > getSkyWcs() const
Definition: Gtransfo.h:520
T str(T... args)
std::string __str__()
Definition: Gtransfo.h:71
GtransfoLinShift(Point const &point)
Definition: Gtransfo.h:460
virtual std::unique_ptr< Gtransfo > composeAndReduce(Gtransfo const &right) const
Return a reduced composition of newTransfo = this(right()), or nullptr if it cannot be reduced...
Definition: Gtransfo.cc:68
virtual double paramRef(const int i) const
Definition: Gtransfo.cc:196
virtual void transformPosAndErrors(const FatPoint &in, FatPoint &out) const
Definition: Gtransfo.cc:116
virtual std::shared_ptr< ast::Mapping > toAstMap(jointcal::Frame const &domain) const
Create an equivalent AST mapping for this transformation, including an analytic inverse if possible...
Definition: Gtransfo.h:166
const GtransfoPoly * getCorr() const
Get a non-owning pointer to the correction transform polynomial.
Definition: Gtransfo.h:549
unsigned getOrder() const
Returns the polynomial order.
Definition: Gtransfo.h:279
A do-nothing transformation. It anyway has dummy routines to mimick a Gtransfo.
Definition: Gtransfo.h:196
just here to provide a specialized constructor, and fit.
Definition: Gtransfo.h:455
void apply(const double xIn, const double yIn, double &xOut, double &yOut) const override
xOut = xIn; yOut = yIn !
Definition: Gtransfo.h:202
STL class.
This one is the Tangent Plane (called gnomonic) projection (from celestial sphere to tangent plane) ...
Definition: Gtransfo.h:653
STL class.
just here to provide specialized constructors. GtransfoLin fit routine.
Definition: Gtransfo.h:482
a virtual (interface) class for geometric transformations.
Definition: Gtransfo.h:42
std::unique_ptr< Gtransfo > composeAndReduce(Gtransfo const &right) const override
Return a reduced composition of newTransfo = this(right()), or nullptr if it cannot be reduced...
Definition: Gtransfo.h:213
double x
std::unique_ptr< GtransfoPoly > corr
Definition: Gtransfo.h:569
virtual void paramDerivatives(Point const &where, double *dx, double *dy) const
Derivative w.r.t parameters.
Definition: Gtransfo.cc:205
std::unique_ptr< Gtransfo > clone() const override
returns a copy (allocated by new) of the transformation.
Definition: Gtransfo.h:314
double fit(StarMatchList const &starMatchList) override
fits a transfo to a std::list of Point pairs (p1,p2, the Point fields in StarMatch).
Definition: Gtransfo.h:207
a run-time transfo that allows users to define a Gtransfo with minimal coding (just the transfo routi...
Definition: Gtransfo.h:703
std::unique_ptr< Gtransfo > gtransfoRead(const std::string &fileName)
The virtual constructor from a file.
Definition: Gtransfo.cc:1694
virtual std::unique_ptr< Gtransfo > clone() const =0
returns a copy (allocated by new) of the transformation.
virtual void dump(std::ostream &stream=std::cout) const =0
dumps the transfo coefficients to stream.
table::Key< int > transform
int getNpar() const
total number of parameters
Definition: Gtransfo.h:476
void apply(const double xIn, const double yIn, double &xOut, double &yOut) const override
Definition: Gtransfo.cc:522
STL class.
int getNpar() const
total number of parameters
Definition: Gtransfo.h:491
T forward(T... args)
virtual GtransfoLin linearApproximation(Point const &where, const double step=0.01) const
linear (local) approximation.
Definition: Gtransfo.cc:109
int y
int getNpar() const override
total number of parameters
Definition: Gtransfo.h:293
virtual std::unique_ptr< Gtransfo > inverseTransfo(const double precision, const Frame &region) const
returns an inverse transfo. Numerical if not overloaded.
Definition: Gtransfo.cc:268
virtual void apply(const double xIn, const double yIn, double &xOut, double &yOut) const =0
std::unique_ptr< Gtransfo > gtransfoCompose(Gtransfo const &left, Gtransfo const &right)
Returns a pointer to a composition of gtransfos, representing left(right()).
Definition: Gtransfo.cc:384