Logo ROOT   6.08/07
Reference Guide
Functions.h
Go to the documentation of this file.
1 // @(#)root/smatrix:$Id$
2 // Authors: T. Glebe, L. Moneta 2005
3 
4 #ifndef ROOT_Math_Functions
5 #define ROOT_Math_Functions
6 // ********************************************************************
7 //
8 // source:
9 //
10 // type: source code
11 //
12 // created: 16. Mar 2001
13 //
14 // author: Thorsten Glebe
15 // HERA-B Collaboration
16 // Max-Planck-Institut fuer Kernphysik
17 // Saupfercheckweg 1
18 // 69117 Heidelberg
19 // Germany
20 // E-mail: T.Glebe@mpi-hd.mpg.de
21 //
22 // Description: Functions which are not applied like a unary operator
23 //
24 // changes:
25 // 16 Mar 2001 (TG) creation
26 // 03 Apr 2001 (TG) minimum added, doc++ comments added
27 // 07 Apr 2001 (TG) Lmag2, Lmag added
28 // 19 Apr 2001 (TG) added #include <cmath>
29 // 24 Apr 2001 (TG) added sign()
30 // 26 Jul 2001 (TG) round() added
31 // 27 Sep 2001 (TG) added Expr declaration
32 //
33 // ********************************************************************
34 #include <cmath>
35 
36 #ifndef ROOT_Math_Expression
37 #include "Math/Expression.h"
38 #endif
39 
40 /**
41  @defgroup TempFunction Generic Template Functions
42  @ingroup SMatrixGroup
43 
44  These functions apply for any type T, such as a scalar, a vector or a matrix.
45  */
46 /**
47  @defgroup VectFunction Vector Template Functions
48  @ingroup SMatrixGroup
49 
50  These functions apply to SVector types (and also to Vector expressions) and can
51  return a vector expression or
52  a scalar, like in the Dot product, or a matrix, like in the Tensor product
53  */
54 
55 
56 namespace ROOT {
57 
58  namespace Math {
59 
60 
61 
62 template <class T, unsigned int D> class SVector;
63 
64 
65 /** square
66  Template function to compute \f$x\cdot x \f$, for any type T returning a type T
67 
68  @ingroup TempFunction
69  @author T. Glebe
70 */
71 //==============================================================================
72 // square: x*x
73 //==============================================================================
74 template <class T>
75 inline const T Square(const T& x) { return x*x; }
76 
77 /** maximum.
78  Template to find max(a,b) where a,b are of type T
79 
80  @ingroup TempFunction
81  @author T. Glebe
82 */
83 //==============================================================================
84 // maximum
85 //==============================================================================
86 template <class T>
87 inline const T Maximum(const T& lhs, const T& rhs) {
88  return (lhs > rhs) ? lhs : rhs;
89 }
90 
91 /** minimum.
92  Template to find min(a,b) where a,b are of type T
93 
94  @ingroup TempFunction
95  @author T. Glebe
96 */
97 //==============================================================================
98 // minimum
99 //==============================================================================
100 template <class T>
101 inline const T Minimum(const T& lhs, const T& rhs) {
102  return (lhs < rhs) ? lhs : rhs;
103 }
104 
105 /** round.
106  Template to compute nearest integer value for any type T
107  @ingroup TempFunction
108  @author T. Glebe
109 */
110 //==============================================================================
111 // round
112 //==============================================================================
113 template <class T>
114 inline int Round(const T& x) {
115  return (x-static_cast<int>(x) < 0.5) ? static_cast<int>(x) : static_cast<int>(x+1);
116 }
117 
118 
119 /** sign.
120  Template to compute the sign of a number
121 
122  @ingroup TempFunction
123  @author T. Glebe
124 */
125 //==============================================================================
126 // sign
127 //==============================================================================
128 template <class T>
129 inline int Sign(const T& x) { return (x==0)? 0 : (x<0)? -1 : 1; }
130 
131 //==============================================================================
132 // meta_dot
133 //==============================================================================
134 template <unsigned int I>
135 struct meta_dot {
136  template <class A, class B, class T>
137  static inline T f(const A& lhs, const B& rhs, const T& x) {
138  return lhs.apply(I) * rhs.apply(I) + meta_dot<I-1>::f(lhs,rhs,x);
139  }
140 };
141 
142 
143 //==============================================================================
144 // meta_dot<0>
145 //==============================================================================
146 template <>
147 struct meta_dot<0> {
148  template <class A, class B, class T>
149  static inline T f(const A& lhs, const B& rhs, const T& /*x */) {
150  return lhs.apply(0) * rhs.apply(0);
151  }
152 };
153 
154 
155 /**
156  Vector dot product.
157  Template to compute \f$\vec{a}\cdot\vec{b} = \sum_i a_i\cdot b_i \f$.
158 
159  @ingroup VectFunction
160  @author T. Glebe
161 */
162 //==============================================================================
163 // dot
164 //==============================================================================
165 template <class T, unsigned int D>
166 inline T Dot(const SVector<T,D>& lhs, const SVector<T,D>& rhs) {
167  return meta_dot<D-1>::f(lhs,rhs, T());
168 }
169 
170 //==============================================================================
171 // dot
172 //==============================================================================
173 template <class A, class T, unsigned int D>
174 inline T Dot(const SVector<T,D>& lhs, const VecExpr<A,T,D>& rhs) {
175  return meta_dot<D-1>::f(lhs,rhs, T());
176 }
177 
178 //==============================================================================
179 // dot
180 //==============================================================================
181 template <class A, class T, unsigned int D>
182 inline T Dot(const VecExpr<A,T,D>& lhs, const SVector<T,D>& rhs) {
183  return meta_dot<D-1>::f(lhs,rhs, T());
184 }
185 
186 
187 //==============================================================================
188 // dot
189 //==============================================================================
190 template <class A, class B, class T, unsigned int D>
191 inline T Dot(const VecExpr<A,T,D>& lhs, const VecExpr<B,T,D>& rhs) {
192  return meta_dot<D-1>::f(rhs,lhs, T());
193 }
194 
195 
196 //==============================================================================
197 // meta_mag
198 //==============================================================================
199 template <unsigned int I>
200 struct meta_mag {
201  template <class A, class T>
202  static inline T f(const A& rhs, const T& x) {
203  return Square(rhs.apply(I)) + meta_mag<I-1>::f(rhs, x);
204  }
205 };
206 
207 
208 //==============================================================================
209 // meta_mag<0>
210 //==============================================================================
211 template <>
212 struct meta_mag<0> {
213  template <class A, class T>
214  static inline T f(const A& rhs, const T& ) {
215  return Square(rhs.apply(0));
216  }
217 };
218 
219 
220 /**
221  Vector magnitude square
222  Template to compute \f$|\vec{v}|^2 = \sum_iv_i^2 \f$.
223 
224  @ingroup VectFunction
225  @author T. Glebe
226 */
227 //==============================================================================
228 // mag2
229 //==============================================================================
230 template <class T, unsigned int D>
231 inline T Mag2(const SVector<T,D>& rhs) {
232  return meta_mag<D-1>::f(rhs, T());
233 }
234 
235 //==============================================================================
236 // mag2
237 //==============================================================================
238 template <class A, class T, unsigned int D>
239 inline T Mag2(const VecExpr<A,T,D>& rhs) {
240  return meta_mag<D-1>::f(rhs, T());
241 }
242 
243 /**
244  Vector magnitude (Euclidian norm)
245  Compute : \f$ |\vec{v}| = \sqrt{\sum_iv_i^2} \f$.
246 
247  @ingroup VectFunction
248  @author T. Glebe
249 */
250 //==============================================================================
251 // mag
252 //==============================================================================
253 template <class T, unsigned int D>
254 inline T Mag(const SVector<T,D>& rhs) {
255  return std::sqrt(Mag2(rhs));
256 }
257 
258 //==============================================================================
259 // mag
260 //==============================================================================
261 template <class A, class T, unsigned int D>
262 inline T Mag(const VecExpr<A,T,D>& rhs) {
263  return std::sqrt(Mag2(rhs));
264 }
265 
266 
267 /** Lmag2: Square of Minkowski Lorentz-Vector norm (only for 4D Vectors)
268  Template to compute \f$ |\vec{v}|^2 = v_0^2 - v_1^2 - v_2^2 -v_3^2 \f$.
269 
270  @ingroup VectFunction
271  @author T. Glebe
272 */
273 //==============================================================================
274 // Lmag2
275 //==============================================================================
276 template <class T>
277 inline T Lmag2(const SVector<T,4>& rhs) {
278  return Square(rhs[0]) - Square(rhs[1]) - Square(rhs[2]) - Square(rhs[3]);
279 }
280 
281 //==============================================================================
282 // Lmag2
283 //==============================================================================
284 template <class A, class T>
285 inline T Lmag2(const VecExpr<A,T,4>& rhs) {
286  return Square(rhs.apply(0))
287  - Square(rhs.apply(1)) - Square(rhs.apply(2)) - Square(rhs.apply(3));
288 }
289 
290 /** Lmag: Minkowski Lorentz-Vector norm (only for 4-dim vectors)
291  Length of a vector Lorentz-Vector:
292  \f$ |\vec{v}| = \sqrt{v_0^2 - v_1^2 - v_2^2 -v_3^2} \f$.
293 
294  @ingroup VectFunction
295  @author T. Glebe
296 */
297 //==============================================================================
298 // Lmag
299 //==============================================================================
300 template <class T>
301 inline T Lmag(const SVector<T,4>& rhs) {
302  return std::sqrt(Lmag2(rhs));
303 }
304 
305 //==============================================================================
306 // Lmag
307 //==============================================================================
308 template <class A, class T>
309 inline T Lmag(const VecExpr<A,T,4>& rhs) {
310  return std::sqrt(Lmag2(rhs));
311 }
312 
313 
314 /** Vector Cross Product (only for 3-dim vectors)
315  \f$ \vec{c} = \vec{a}\times\vec{b} \f$.
316 
317  @ingroup VectFunction
318  @author T. Glebe
319 */
320 //==============================================================================
321 // cross product
322 //==============================================================================
323 template <class T>
324 inline SVector<T,3> Cross(const SVector<T,3>& lhs, const SVector<T,3>& rhs) {
325  return SVector<T,3>(lhs.apply(1)*rhs.apply(2) -
326  lhs.apply(2)*rhs.apply(1),
327  lhs.apply(2)*rhs.apply(0) -
328  lhs.apply(0)*rhs.apply(2),
329  lhs.apply(0)*rhs.apply(1) -
330  lhs.apply(1)*rhs.apply(0));
331 }
332 
333 //==============================================================================
334 // cross product
335 //==============================================================================
336 template <class A, class T>
337 inline SVector<T,3> Cross(const VecExpr<A,T,3>& lhs, const SVector<T,3>& rhs) {
338  return SVector<T,3>(lhs.apply(1)*rhs.apply(2) -
339  lhs.apply(2)*rhs.apply(1),
340  lhs.apply(2)*rhs.apply(0) -
341  lhs.apply(0)*rhs.apply(2),
342  lhs.apply(0)*rhs.apply(1) -
343  lhs.apply(1)*rhs.apply(0));
344 }
345 
346 //==============================================================================
347 // cross product
348 //==============================================================================
349 template <class T, class A>
350 inline SVector<T,3> Cross(const SVector<T,3>& lhs, const VecExpr<A,T,3>& rhs) {
351  return SVector<T,3>(lhs.apply(1)*rhs.apply(2) -
352  lhs.apply(2)*rhs.apply(1),
353  lhs.apply(2)*rhs.apply(0) -
354  lhs.apply(0)*rhs.apply(2),
355  lhs.apply(0)*rhs.apply(1) -
356  lhs.apply(1)*rhs.apply(0));
357 }
358 
359 //==============================================================================
360 // cross product
361 //==============================================================================
362 template <class A, class B, class T>
363 inline SVector<T,3> Cross(const VecExpr<A,T,3>& lhs, const VecExpr<B,T,3>& rhs) {
364  return SVector<T,3>(lhs.apply(1)*rhs.apply(2) -
365  lhs.apply(2)*rhs.apply(1),
366  lhs.apply(2)*rhs.apply(0) -
367  lhs.apply(0)*rhs.apply(2),
368  lhs.apply(0)*rhs.apply(1) -
369  lhs.apply(1)*rhs.apply(0));
370 }
371 
372 
373 /** Unit.
374  Return a vector of unit length: \f$ \vec{e}_v = \vec{v}/|\vec{v}| \f$.
375 
376  @ingroup VectFunction
377  @author T. Glebe
378 */
379 //==============================================================================
380 // unit: returns a unit vector
381 //==============================================================================
382 template <class T, unsigned int D>
383 inline SVector<T,D> Unit(const SVector<T,D>& rhs) {
384  return SVector<T,D>(rhs).Unit();
385 }
386 
387 //==============================================================================
388 // unit: returns a unit vector
389 //==============================================================================
390 template <class A, class T, unsigned int D>
391 inline SVector<T,D> Unit(const VecExpr<A,T,D>& rhs) {
392  return SVector<T,D>(rhs).Unit();
393 }
394 
395 #ifdef XXX
396 //==============================================================================
397 // unit: returns a unit vector (worse performance)
398 //==============================================================================
399 template <class T, unsigned int D>
401  unit(const SVector<T,D>& lhs) {
402  typedef BinaryOp<DivOp<T>, SVector<T,D>, Constant<T>, T> DivOpBinOp;
403  return VecExpr<DivOpBinOp,T,D>(DivOpBinOp(DivOp<T>(),lhs,Constant<T>(mag(lhs))));
404 }
405 
406 //==============================================================================
407 // unit: returns a unit vector (worse performance)
408 //==============================================================================
409 template <class A, class T, unsigned int D>
410 inline VecExpr<BinaryOp<DivOp<T>, VecExpr<A,T,D>, Constant<T>, T>, T, D>
411  unit(const VecExpr<A,T,D>& lhs) {
412  typedef BinaryOp<DivOp<T>, VecExpr<A,T,D>, Constant<T>, T> DivOpBinOp;
413  return VecExpr<DivOpBinOp,T,D>(DivOpBinOp(DivOp<T>(),lhs,Constant<T>(mag(lhs))));
414 }
415 #endif
416 
417 
418  } // namespace Math
419 
420 } // namespace ROOT
421 
422 
423 
424 #endif /* ROOT_Math_Functions */
T Dot(const SVector< T, D > &lhs, const SVector< T, D > &rhs)
Vector dot product.
Definition: Functions.h:166
static double B[]
Constant expression class A class representing constant expressions (literals) in the parse tree...
Definition: Expression.h:400
T Lmag(const SVector< T, 4 > &rhs)
Lmag: Minkowski Lorentz-Vector norm (only for 4-dim vectors) Length of a vector Lorentz-Vector: ...
Definition: Functions.h:301
This namespace contains pre-defined functions to be used in conjuction with TExecutor::Map and TExecu...
Definition: StringConv.hxx:21
Division (element-wise) Operation Class.
double T(double x)
Definition: ChebyshevPol.h:34
T Mag(const SVector< T, D > &rhs)
Vector magnitude (Euclidian norm) Compute : .
Definition: Functions.h:254
static T f(const A &lhs, const B &rhs, const T &x)
Definition: Functions.h:137
static T f(const A &rhs, const T &x)
Definition: Functions.h:202
int Round(const T &x)
round.
Definition: Functions.h:114
const T Minimum(const T &lhs, const T &rhs)
minimum.
Definition: Functions.h:101
const T Maximum(const T &lhs, const T &rhs)
maximum.
Definition: Functions.h:87
static double A[]
int Sign(const T &x)
sign.
Definition: Functions.h:129
double sqrt(double)
Double_t x[n]
Definition: legend1.C:17
static T f(const A &lhs, const B &rhs, const T &)
Definition: Functions.h:149
BinaryOperation class A class representing binary operators in the parse tree.
Definition: Expression.h:234
T apply(unsigned int i) const
access the parse tree. Index starts from zero
Definition: SVector.icc:530
T Mag2(const SVector< T, D > &rhs)
Vector magnitude square Template to compute .
Definition: Functions.h:231
SVector< T, 3 > Cross(const SVector< T, 3 > &lhs, const SVector< T, 3 > &rhs)
Vector Cross Product (only for 3-dim vectors) .
Definition: Functions.h:324
T Lmag2(const SVector< T, 4 > &rhs)
Lmag2: Square of Minkowski Lorentz-Vector norm (only for 4D Vectors) Template to compute ...
Definition: Functions.h:277
Namespace for new Math classes and functions.
Expression wrapper class for Vector objects.
Definition: Expression.h:64
static T f(const A &rhs, const T &)
Definition: Functions.h:214
const T Square(const T &x)
square Template function to compute , for any type T returning a type T
Definition: Functions.h:75
SVector< T, D > Unit(const SVector< T, D > &rhs)
Unit.
Definition: Functions.h:383
#define I(x, y, z)
SVector: a generic fixed size Vector class.
T apply(unsigned int i) const
Definition: Expression.h:77