casacore
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Gaussian3D.h
Go to the documentation of this file.
1 //# Gaussian3D.h: A three-dimensional Gaussian class
2 //# Copyright (C) 1995,1996,1997,2001,2002,2005
3 //# Associated Universities, Inc. Washington DC, USA.
4 //#
5 //# This library is free software; you can redistribute it and/or modify it
6 //# under the terms of the GNU Library General Public License as published by
7 //# the Free Software Foundation; either version 2 of the License, or (at your
8 //# option) any later version.
9 //#
10 //# This library is distributed in the hope that it will be useful, but WITHOUT
11 //# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
12 //# FITNESS FOR A PARTICULAR PURPOSE. See the GNU Library General Public
13 //# License for more details.
14 //#
15 //# You should have received a copy of the GNU Library General Public License
16 //# along with this library; if not, write to the Free Software Foundation,
17 //# Inc., 675 Massachusetts Ave, Cambridge, MA 02139, USA.
18 //#
19 //# Correspondence concerning AIPS++ should be addressed as follows:
20 //# Internet email: aips2-request@nrao.edu.
21 //# Postal address: AIPS++ Project Office
22 //# National Radio Astronomy Observatory
23 //# 520 Edgemont Road
24 //# Charlottesville, VA 22903-2475 USA
25 //#
26 //# $Id$
27 
28 #ifndef SCIMATH_GAUSSIAN3D_H
29 #define SCIMATH_GAUSSIAN3D_H
30 
31 #include <casacore/casa/aips.h>
37 
38 namespace casacore { //# NAMESPACE CASACORE - BEGIN
39 
40 // <summary> A three dimensional Gaussian class.</summary>
41 
42 // <use visibility=export>
43 
44 // <reviewed reviewer="" date="" tests="tGaussian3D"
45 // demos="">
46 // </reviewed>
47 
48 // <prerequisite>
49 // <li> <linkto class="Gaussian3DParam">Gaussian3DParam</linkto>
50 // <li> <linkto class="Function">Function</linkto>
51 // </prerequisite>
52 
53 // <etymology>
54 // A Gaussian3D functional is designed exclusively for calculating a
55 // Gaussian (or Normal) distribution in three dimensions. Other classes exist
56 // for calculating these functions in one
57 // (<linkto class=Gaussian1D>Gaussian1D</linkto>), two
58 // (<linkto class=Gaussian2D>Gaussian2D</linkto>), and N
59 // (<linkto class=GaussianND>GaussianND</linkto>) dimensions.
60 // </etymology>
61 
62 // <synopsis>
63 // A <src>Gaussian3D</src> is described by a height, center, and width,
64 // and position angles. Its fundamental operation is evaluating itself
65 // at some <src>(x,y,z)</src>
66 // coordinate. Its parameters (height, center and width, position angles) may
67 // be changed at run time.
68 
69 // The width of the Gaussian is now specified in terms of the full width
70 // at half maximum (FWHM), like the 2D and 1D Gaussian functional classes.
71 
72 // The three axis values refer to the x, y, and z axes, and unlike with the
73 // 2D Gaussian any of the three axes may be the longest; instead, the position
74 // angles are restricted. The first position angle, theta, is the longitudinal
75 // angle, referring to the rotation (counterclockwise) around the z-axis. The
76 // second, phi, is the latidudinal angle, referring to the rotation around
77 // the theta-rotated y axis. The domain of both angles is -pi/4 < A < pi/4,
78 // although the angles are not constrained when fitting and can be set outside
79 // the domain by setting the parameters directly using Functional operator[].
80 // (Note that the use of theta and phi corresponds to the mathematics
81 // convention for these angles, not the physics convention.)
82 
83 // The parameter interface (see
84 // <linkto class="Gaussian3DParam">Gaussian3DParam</linkto> class),
85 // is used to provide an interface to the
86 // <linkto module="Fitting">Fitting</linkto> classes.
87 //
88 // There are 9 parameters that are used to describe the Gaussian:
89 // <ol>
90 // <li> The height of the Gaussian. This is identical to the value
91 // returned using the <src> height </src> member function.
92 // <li> The center of the Gaussian in the x direction. This is identical to
93 // the value returned using the <src> xCenter </src> member function.
94 // <li> The center of the Gaussian in the y direction. This is identical to
95 // the value returned using the <src> yCenter </src> member function.
96 // <li> The center of the Gaussian in the z direction. This is identical to
97 // the value returned using the <src> zCenter </src> member function.
98 // <li> The width of the Gaussian along the x-axis.
99 // <li> The width of the Gaussian along the y-axis.
100 // <li> The width of the Gaussian along the z-axis.
101 // <li> The longitudinal position angle, theta (in radians)
102 // <li> The latitudinal position angle, phi (also in radians).
103 // </ol>
104 
105 // An enumeration for the parameter index is provided, enabling the setting
106 // and reading of parameters with the <src>[]</src> operator. The
107 // <src>mask()</src> methods can be used to check and set the parameter masks.
108 //
109 // </synopsis>
110 
111 // <example>
112 // <srcblock>
113 // Gaussian3D<Double> g(9.0, 0.0, 0.0, 0.0, 1.0, 1.0, 1.0, 0.0, 0.0);
114 // Vector<Double> x(3);
115 // x(0) = 1.0; x(1) = 0.5; x(2) = 0.0
116 // cout << "g(" << x(0) << "," << x(1) << "," << x(2) << ")=" << g(x) << endl;
117 // </srcblock>
118 // </example>
119 
120 // <motivation>
121 // The GaussianND class does not contain explicit derivatives
122 // and was insufficient for fitting 3D Gaussians to data.
123 // </motivation>
124 
125 // <templating arg=T>
126 // <li> T should have standard numerical operators and exp() function. Current
127 // implementation only tested for real types.
128 // <li> To obtain derivatives, the derivatives should be defined.
129 // </templating>
130 
131 // <thrown>
132 // <li> Assertion in debug mode if attempt is made to set a negative width
133 // <li> AipsError if incorrect parameter number specified.
134 // <li> Assertion in debug mode if operator(Vector<>) with empty Vector
135 // </thrown>
136 
137 // <todo asof="2002/07/22">
138 // <li> Optimize derivative calculations for faster fitting?
139 // </todo>
140 
141 
142 
143 template<class T> class Gaussian3D : public Gaussian3DParam<T>
144 {
145 public:
146  // A functional for a rotated, 3D Gaussian. Similar to Gaussian2D, but
147  // the xWidth, yWidth, and zWidth parameters are not adjusted for FWHM;
148  // they are identical to the parameters used in the function.
149 
150  // Constructs the three-dimensional Gaussians. Defaults:
151  // height = 1, center = {0,0,0}, width = {1,1,1}, theta = phi = 0.
152  // The center and width vectors must have three elements.
153  // <group>
154  Gaussian3D();
155  Gaussian3D(T height, const Vector<T>& center,
156  const Vector<T>& width, T theta, T phi);
157  Gaussian3D(T &height, T &xCenter, T &yCenter, T &zCenter,
158  T &xWidth, T &yWidth, T &zWidth, T &theta, T &phi);
159  // </group>
160 
161  // Copy constructor
162  // <group>
163  Gaussian3D(const Gaussian3D<T> &other);
164  template <class W>
165  Gaussian3D(const Gaussian3D<W> &other) : Gaussian3DParam<T>(other) {}
166  // </group>
167 
168  // Destructor
169  virtual ~Gaussian3D();
170 
171  // Assignment operator
172  Gaussian3D<T> &operator=(const Gaussian3D<T> &other);
173 
174  // Evaluate the Gaussian at <src>x</src>.
175  virtual T eval(typename Function<T>::FunctionArg x) const;
176 
177  // Return a copy of this object from the heap. The caller is responsible
178  // for deleting this pointer.
179  // <group>
180  virtual Function<T> *clone() const;
185  // </group>
186 
187 private:
188  // AutoDiff does not have a square() function, so one is provided here.
189  T sq(T v) const;
190 
191  //# Make members of parent classes known.
192 protected:
204 public:
205  using Gaussian3DParam<T>::H;
216 };
217 
218 
219 // AUTODIFF SPECIALIZATION
220 
221 #define Gaussian3D_PS Gaussian3D
222 
223 // <summary> Partial specialization of Gaussian3D for <src>AutoDiff</src>
224 // </summary>
225 
226 // <synopsis>
227 // <note role=warning> The name <src>Gaussian3D_PS</src> is only for cxx2html
228 // documentation problems. Use <src>Gaussian3D</src> in your code.</note>
229 // </synopsis>
230 
232 template <class T> class Gaussian3D_PS<AutoDiff<T> > : public Gaussian3DParam<AutoDiff<T> >
233 {
234 public:
235  Gaussian3D_PS();
236  Gaussian3D_PS(const AutoDiff<T> &height,
237  const Vector<AutoDiff<T> >& center,
238  const Vector<AutoDiff<T> >& width,
239  const AutoDiff<T>& theta,
240  const AutoDiff<T>& phi);
241  Gaussian3D_PS(AutoDiff<T>& height, AutoDiff<T>& xCenter,
242  AutoDiff<T>& yCenter, AutoDiff<T>& zCenter,
243  AutoDiff<T>& xWidth, AutoDiff<T>& yWidth,
244  AutoDiff<T>& zWidth, AutoDiff<T>& theta,
245  AutoDiff<T>& phi);
246  Gaussian3D_PS(const Gaussian3D_PS<AutoDiff<T> > &other);
247  template <class W>
248  Gaussian3D_PS(const Gaussian3D_PS<W> &other) :
249  Gaussian3DParam<AutoDiff<T> >(other) {}
250  virtual ~Gaussian3D_PS();
251 //
252  Gaussian3D_PS<AutoDiff<T> > &operator=(const Gaussian3D_PS<AutoDiff<T> > &other);
253 //
254  virtual AutoDiff<T> eval(typename Function<AutoDiff<T> >::FunctionArg x) const;
255  virtual Function<AutoDiff<T> > *clone() const;
257  *cloneAD() const {
258  return new Gaussian3D<typename FunctionTraits<AutoDiff<T> >::DiffType>
259  (*this); }
261  *cloneNonAD() const {
262  return new Gaussian3D<typename FunctionTraits<AutoDiff<T> >::BaseType>
263  (*this); }
264 
265 private:
266  T sq(T v) const;
267 
268  //# Make members of parent classes known.
269 protected:
270  using Gaussian3DParam<AutoDiff<T> >::param_p;
271  using Gaussian3DParam<AutoDiff<T> >::stoT_p;
272  using Gaussian3DParam<AutoDiff<T> >::stoP_p;
273  using Gaussian3DParam<AutoDiff<T> >::cosT_p;
274  using Gaussian3DParam<AutoDiff<T> >::cosP_p;
275  using Gaussian3DParam<AutoDiff<T> >::sinT_p;
276  using Gaussian3DParam<AutoDiff<T> >::sinP_p;
277  using Gaussian3DParam<AutoDiff<T> >::cosTcosP_p;
278  using Gaussian3DParam<AutoDiff<T> >::cosTsinP_p;
279  using Gaussian3DParam<AutoDiff<T> >::sinTcosP_p;
280  using Gaussian3DParam<AutoDiff<T> >::sinTsinP_p;
281 public:
282  using Gaussian3DParam<AutoDiff<T> >::H;
283  using Gaussian3DParam<AutoDiff<T> >::CX;
284  using Gaussian3DParam<AutoDiff<T> >::CY;
285  using Gaussian3DParam<AutoDiff<T> >::CZ;
286  using Gaussian3DParam<AutoDiff<T> >::AX;
287  using Gaussian3DParam<AutoDiff<T> >::AY;
288  using Gaussian3DParam<AutoDiff<T> >::AZ;
289  using Gaussian3DParam<AutoDiff<T> >::THETA;
290  using Gaussian3DParam<AutoDiff<T> >::PHI;
291  using Gaussian3DParam<AutoDiff<T> >::fwhm2int;
292  using Gaussian3DParam<AutoDiff<T> >::settrigvals;
293 };
294 
295 #undef Gaussian3D_PS
296 
297 
298 } //# NAMESPACE CASACORE - END
299 
300 #ifndef CASACORE_NO_AUTO_TEMPLATES
301 #include <casacore/scimath/Functionals/Gaussian3D.tcc>
302 #include <casacore/scimath/Functionals/Gaussian3D2.tcc>
303 #endif //# CASACORE_NO_AUTO_TEMPLATES
304 #endif
305 
306 
307 
308 
309 
T height() const
Get or set the peak height of the Gaussian.
virtual T eval(typename Function< T >::FunctionArg x) const
Evaluate the Gaussian at x.
T theta() const
Get or set the rotation angles of the Gaussian.
virtual Function< T > * clone() const
Return a copy of this object from the heap.
PtrHolder< T > & operator=(const PtrHolder< T > &other)
Gaussian3D(const Gaussian3D< W > &other)
Definition: Gaussian3D.h:165
virtual Function< typename FunctionTraits< T >::DiffType > * cloneAD() const
Definition: Gaussian3D.h:181
Vector< T > center() const
Get or cet the center coordinates of the Gaussian.
Parameter handling for 3 dimensional Gaussian class.
virtual Function< typename FunctionTraits< T >::BaseType > * cloneNonAD() const
Definition: Gaussian3D.h:183
Gaussian3D()
A functional for a rotated, 3D Gaussian.
Numerical functional interface class.
Definition: GenericL2Fit.h:46
Vector< T > width() const
Get or set the sigma-width of the Gaussian.
virtual ~Gaussian3D()
Destructor.
Class that computes partial derivatives by automatic differentiation.
Definition: AutoDiff.h:257
A three dimensional Gaussian class.
Definition: Gaussian3D.h:143
Gaussian3D< T > & operator=(const Gaussian3D< T > &other)
Assignment operator.
T sq(T v) const
AutoDiff does not have a square() function, so one is provided here.
#define Gaussian3D_PS
AUTODIFF SPECIALIZATION.
Definition: Gaussian3D.h:221