| Geant4 Cross Reference |
1 // -*- C++ -*-
2 // ---------------------------------------------------------------------------
3 //
4 // This file is a part of the CLHEP - a Class Library for High Energy Physics.
5 //
6 // History:
7 // 12.06.01 E.Chernyaev - CLHEP-1.7: initial version
8 // 14.03.03 E.Chernyaev - CLHEP-1.9: template version
9 //
10
11 #ifndef BASIC_VECTOR3D_H
12 #define BASIC_VECTOR3D_H
13
14 #include <iosfwd>
15 #include <type_traits>
16 #include "CLHEP/Vector/ThreeVector.h"
17
18 namespace HepGeom {
19 /**
20 * Base class for Point3D<T>, Vector3D<T> and Normal3D<T>.
21 * It defines only common functionality for those classes and
22 * should not be used as separate class.
23 *
24 * @author Evgeni Chernyaev <Evgueni.Tcherniaev@cern.ch>
25 * @ingroup geometry
26 */
27 template<class T> class BasicVector3D {
28 protected:
29 T v_[3];
30
31 /**
32 * Default constructor.
33 * It is protected - this class should not be instantiated directly.
34 */
35 BasicVector3D() { v_[0] = 0; v_[1] = 0; v_[2] = 0; }
36
37 public:
38 /**
39 * Safe indexing of the coordinates when using with matrices, arrays, etc.
40 */
41 enum {
42 X = 0, /**< index for x-component */
43 Y = 1, /**< index for y-component */
44 Z = 2, /**< index for z-component */
45 NUM_COORDINATES = 3, /**< number of components */
46 SIZE = NUM_COORDINATES /**< number of components */
47 };
48
49 /**
50 * Constructor from three numbers. */
51 BasicVector3D(T x1, T y1, T z1) { v_[0] = x1; v_[1] = y1; v_[2] = z1; }
52
53 /**
54 * Copy constructor. */
55 BasicVector3D(const BasicVector3D<T> &) = default;
56
57 /**
58 * Constructor for BasicVector3D<double> from BasicVector3D<float>. */
59 template<typename U = T,
60 typename = typename std::enable_if<!std::is_same<U,float>::value >::type>
61 BasicVector3D(const BasicVector3D<float> & v) {
62 v_[0] = v.x(); v_[1] = v.y(); v_[2] = v.z();
63 }
64
65 /**
66 * Move constructor. */
67 BasicVector3D(BasicVector3D<T> &&) = default;
68
69 /**
70 * Destructor. */
71 virtual ~BasicVector3D() = default;
72
73 // -------------------------
74 // Interface to "good old C"
75 // -------------------------
76
77 /**
78 * Conversion (cast) to ordinary array. */
79 operator T * () { return v_; }
80
81 /**
82 * Conversion (cast) to ordinary const array. */
83 operator const T * () const { return v_; }
84
85 /**
86 * Conversion (cast) to CLHEP::Hep3Vector.
87 * This operator is needed only for backward compatibility and
88 * in principle should not exit.
89 */
90 operator CLHEP::Hep3Vector () const { return CLHEP::Hep3Vector(x(),y(),z()); }
91
92 // -----------------------------
93 // General arithmetic operations
94 // -----------------------------
95
96 /**
97 * Assignment. */
98 BasicVector3D<T> & operator= (const BasicVector3D<T> &) = default;
99 /**
100 * Move assignment. */
101 BasicVector3D<T> & operator= (BasicVector3D<T> &&) = default;
102 /**
103 * Addition. */
104 BasicVector3D<T> & operator+=(const BasicVector3D<T> & v) {
105 v_[0] += v.v_[0]; v_[1] += v.v_[1]; v_[2] += v.v_[2]; return *this;
106 }
107 /**
108 * Subtraction. */
109 BasicVector3D<T> & operator-=(const BasicVector3D<T> & v) {
110 v_[0] -= v.v_[0]; v_[1] -= v.v_[1]; v_[2] -= v.v_[2]; return *this;
111 }
112 /**
113 * Multiplication by scalar. */
114 BasicVector3D<T> & operator*=(double a) {
115 v_[0] *= a; v_[1] *= a; v_[2] *= a; return *this;
116 }
117 /**
118 * Division by scalar. */
119 BasicVector3D<T> & operator/=(double a) {
120 v_[0] /= a; v_[1] /= a; v_[2] /= a; return *this;
121 }
122
123 // ------------
124 // Subscripting
125 // ------------
126
127 /**
128 * Gets components by index. */
129 T operator()(int i) const { return v_[i]; }
130 /**
131 * Gets components by index. */
132 T operator[](int i) const { return v_[i]; }
133
134 /**
135 * Sets components by index. */
136 T & operator()(int i) { return v_[i]; }
137 /**
138 * Sets components by index. */
139 T & operator[](int i) { return v_[i]; }
140
141 // ------------------------------------
142 // Cartesian coordinate system: x, y, z
143 // ------------------------------------
144
145 /**
146 * Gets x-component in cartesian coordinate system. */
147 T x() const { return v_[0]; }
148 /**
149 * Gets y-component in cartesian coordinate system. */
150 T y() const { return v_[1]; }
151 /**
152 * Gets z-component in cartesian coordinate system. */
153 T z() const { return v_[2]; }
154
155 /**
156 * Sets x-component in cartesian coordinate system. */
157 void setX(T a) { v_[0] = a; }
158 /**
159 * Sets y-component in cartesian coordinate system. */
160 void setY(T a) { v_[1] = a; }
161 /**
162 * Sets z-component in cartesian coordinate system. */
163 void setZ(T a) { v_[2] = a; }
164
165 /**
166 * Sets components in cartesian coordinate system. */
167 void set(T x1, T y1, T z1) { v_[0] = x1; v_[1] = y1; v_[2] = z1; }
168
169 // ------------------------------------------
170 // Cylindrical coordinate system: rho, phi, z
171 // ------------------------------------------
172
173 /**
174 * Gets transverse component squared. */
175 T perp2() const { return x()*x()+y()*y(); }
176 /**
177 * Gets transverse component. */
178 T perp() const { return std::sqrt(perp2()); }
179 /**
180 * Gets rho-component in cylindrical coordinate system */
181 T rho() const { return perp(); }
182
183 /**
184 * Sets transverse component keeping phi and z constant. */
185 void setPerp(T rh) {
186 T factor = perp();
187 if (factor > 0) {
188 factor = rh/factor; v_[0] *= factor; v_[1] *= factor;
189 }
190 }
191
192 // ------------------------------------------
193 // Spherical coordinate system: r, phi, theta
194 // ------------------------------------------
195
196 /**
197 * Gets magnitude squared of the vector. */
198 T mag2() const { return x()*x()+y()*y()+z()*z(); }
199 /**
200 * Gets magnitude of the vector. */
201 T mag() const { return std::sqrt(mag2()); }
202 /**
203 * Gets r-component in spherical coordinate system */
204 T r() const { return mag(); }
205 /**
206 * Gets azimuth angle. */
207 T phi() const {
208 return x() == 0 && y() == 0 ? 0 : std::atan2(y(),x());
209 }
210 /**
211 * Gets polar angle. */
212 T theta() const {
213 return x() == 0 && y() == 0 && z() == 0 ? 0 : std::atan2(perp(),z());
214 }
215 /**
216 * Gets cosine of polar angle. */
217 T cosTheta() const { T ma = mag(); return ma == 0 ? 1 : z()/ma; }
218
219 /**
220 * Gets r-component in spherical coordinate system */
221 T getR() const { return r(); }
222 /**
223 * Gets phi-component in spherical coordinate system */
224 T getPhi() const { return phi(); }
225 /**
226 * Gets theta-component in spherical coordinate system */
227 T getTheta() const { return theta(); }
228
229 /**
230 * Sets magnitude. */
231 void setMag(T ma) {
232 T factor = mag();
233 if (factor > 0) {
234 factor = ma/factor; v_[0] *= factor; v_[1] *= factor; v_[2] *= factor;
235 }
236 }
237 /**
238 * Sets r-component in spherical coordinate system. */
239 void setR(T ma) { setMag(ma); }
240 /**
241 * Sets phi-component in spherical coordinate system. */
242 void setPhi(T ph) { T xy = perp(); setX(xy*std::cos(ph)); setY(xy*std::sin(ph)); }
243 /**
244 * Sets theta-component in spherical coordinate system. */
245 void setTheta(T th) {
246 T ma = mag();
247 T ph = phi();
248 set(ma*std::sin(th)*std::cos(ph), ma*std::sin(th)*std::sin(ph), ma*std::cos(th));
249 }
250
251 // ---------------
252 // Pseudo rapidity
253 // ---------------
254
255 /**
256 * Gets pseudo-rapidity: -ln(tan(theta/2)) */
257 T pseudoRapidity() const;
258 /**
259 * Gets pseudo-rapidity. */
260 T eta() const { return pseudoRapidity(); }
261 /**
262 * Gets pseudo-rapidity. */
263 T getEta() const { return pseudoRapidity(); }
264
265 /**
266 * Sets pseudo-rapidity, keeping magnitude and phi fixed. */
267 void setEta(T a);
268
269 // -------------------
270 // Combine two vectors
271 // -------------------
272
273 /**
274 * Scalar product. */
275 T dot(const BasicVector3D<T> & v) const {
276 return x()*v.x()+y()*v.y()+z()*v.z();
277 }
278
279 /**
280 * Vector product. */
281 BasicVector3D<T> cross(const BasicVector3D<T> & v) const {
282 return BasicVector3D<T>(y()*v.z()-v.y()*z(),
283 z()*v.x()-v.z()*x(),
284 x()*v.y()-v.x()*y());
285 }
286
287 /**
288 * Returns transverse component w.r.t. given axis squared. */
289 T perp2(const BasicVector3D<T> & v) const {
290 T tot = v.mag2(), s = dot(v);
291 return tot > 0 ? mag2()-s*s/tot : mag2();
292 }
293
294 /**
295 * Returns transverse component w.r.t. given axis. */
296 T perp(const BasicVector3D<T> & v) const {
297 return std::sqrt(perp2(v));
298 }
299
300 /**
301 * Returns angle w.r.t. another vector. */
302 T angle(const BasicVector3D<T> & v) const;
303
304 // ---------------
305 // Related vectors
306 // ---------------
307
308 /**
309 * Returns unit vector parallel to this. */
310 BasicVector3D<T> unit() const {
311 T len = mag();
312 return (len > 0) ?
313 BasicVector3D<T>(x()/len, y()/len, z()/len) : BasicVector3D<T>();
314 }
315
316 /**
317 * Returns orthogonal vector. */
318 BasicVector3D<T> orthogonal() const {
319 T dx = x() < 0 ? -x() : x();
320 T dy = y() < 0 ? -y() : y();
321 T dz = z() < 0 ? -z() : z();
322 if (dx < dy) {
323 return dx < dz ?
324 BasicVector3D<T>(0,z(),-y()) : BasicVector3D<T>(y(),-x(),0);
325 }else{
326 return dy < dz ?
327 BasicVector3D<T>(-z(),0,x()) : BasicVector3D<T>(y(),-x(),0);
328 }
329 }
330
331 // ---------
332 // Rotations
333 // ---------
334
335 /**
336 * Rotates around x-axis. */
337 BasicVector3D<T> & rotateX(T a);
338 /**
339 * Rotates around y-axis. */
340 BasicVector3D<T> & rotateY(T a);
341 /**
342 * Rotates around z-axis. */
343 BasicVector3D<T> & rotateZ(T a);
344 /**
345 * Rotates around the axis specified by another vector. */
346 BasicVector3D<T> & rotate(T a, const BasicVector3D<T> & v);
347 };
348
349 /*************************************************************************
350 * *
351 * Non-member functions for BasicVector3D<float> *
352 * *
353 *************************************************************************/
354
355 /**
356 * Output to stream.
357 * @relates BasicVector3D
358 */
359 std::ostream &
360 operator<<(std::ostream &, const BasicVector3D<float> &);
361
362 /**
363 * Input from stream.
364 * @relates BasicVector3D
365 */
366 std::istream &
367 operator>>(std::istream &, BasicVector3D<float> &);
368
369 /**
370 * Unary plus.
371 * @relates BasicVector3D
372 */
373 inline BasicVector3D<float>
374 operator+(const BasicVector3D<float> & v) { return v; }
375
376 /**
377 * Addition of two vectors.
378 * @relates BasicVector3D
379 */
380 inline BasicVector3D<float>
381 operator+(const BasicVector3D<float> & a, const BasicVector3D<float> & b) {
382 return BasicVector3D<float>(a.x()+b.x(), a.y()+b.y(), a.z()+b.z());
383 }
384
385 /**
386 * Unary minus.
387 * @relates BasicVector3D
388 */
389 inline BasicVector3D<float>
390 operator-(const BasicVector3D<float> & v) {
391 return BasicVector3D<float>(-v.x(), -v.y(), -v.z());
392 }
393
394 /**
395 * Subtraction of two vectors.
396 * @relates BasicVector3D
397 */
398 inline BasicVector3D<float>
399 operator-(const BasicVector3D<float> & a, const BasicVector3D<float> & b) {
400 return BasicVector3D<float>(a.x()-b.x(), a.y()-b.y(), a.z()-b.z());
401 }
402
403 /**
404 * Multiplication vector by scalar.
405 * @relates BasicVector3D
406 */
407 inline BasicVector3D<float>
408 operator*(const BasicVector3D<float> & v, double a) {
409 return BasicVector3D<float>(v.x()*static_cast<float>(a), v.y()*static_cast<float>(a), v.z()*static_cast<float>(a));
410 }
411
412 /**
413 * Scalar product of two vectors.
414 * @relates BasicVector3D
415 */
416 inline float
417 operator*(const BasicVector3D<float> & a, const BasicVector3D<float> & b) {
418 return a.dot(b);
419 }
420
421 /**
422 * Multiplication scalar by vector.
423 * @relates BasicVector3D
424 */
425 inline BasicVector3D<float>
426 operator*(double a, const BasicVector3D<float> & v) {
427 return BasicVector3D<float>(static_cast<float>(a)*v.x(), static_cast<float>(a)*v.y(), static_cast<float>(a)*v.z());
428 }
429
430 /**
431 * Division vector by scalar.
432 * @relates BasicVector3D
433 */
434 inline BasicVector3D<float>
435 operator/(const BasicVector3D<float> & v, double a) {
436 return BasicVector3D<float>(v.x()/static_cast<float>(a), v.y()/static_cast<float>(a), v.z()/static_cast<float>(a));
437 }
438
439 /**
440 * Comparison of two vectors for equality.
441 * @relates BasicVector3D
442 */
443 inline bool
444 operator==(const BasicVector3D<float> & a, const BasicVector3D<float> & b) {
445 return (a.x()==b.x() && a.y()==b.y() && a.z()==b.z());
446 }
447
448 /**
449 * Comparison of two vectors for inequality.
450 * @relates BasicVector3D
451 */
452 inline bool
453 operator!=(const BasicVector3D<float> & a, const BasicVector3D<float> & b) {
454 return (a.x()!=b.x() || a.y()!=b.y() || a.z()!=b.z());
455 }
456
457 /*************************************************************************
458 * *
459 * Non-member functions for BasicVector3D<double> *
460 * *
461 *************************************************************************/
462
463 /**
464 * Output to stream.
465 * @relates BasicVector3D
466 */
467 std::ostream &
468 operator<<(std::ostream &, const BasicVector3D<double> &);
469
470 /**
471 * Input from stream.
472 * @relates BasicVector3D
473 */
474 std::istream &
475 operator>>(std::istream &, BasicVector3D<double> &);
476
477 /**
478 * Unary plus.
479 * @relates BasicVector3D
480 */
481 inline BasicVector3D<double>
482 operator+(const BasicVector3D<double> & v) { return v; }
483
484 /**
485 * Addition of two vectors.
486 * @relates BasicVector3D
487 */
488 inline BasicVector3D<double>
489 operator+(const BasicVector3D<double> & a,const BasicVector3D<double> & b) {
490 return BasicVector3D<double>(a.x()+b.x(), a.y()+b.y(), a.z()+b.z());
491 }
492
493 /**
494 * Unary minus.
495 * @relates BasicVector3D
496 */
497 inline BasicVector3D<double>
498 operator-(const BasicVector3D<double> & v) {
499 return BasicVector3D<double>(-v.x(), -v.y(), -v.z());
500 }
501
502 /**
503 * Subtraction of two vectors.
504 * @relates BasicVector3D
505 */
506 inline BasicVector3D<double>
507 operator-(const BasicVector3D<double> & a,const BasicVector3D<double> & b) {
508 return BasicVector3D<double>(a.x()-b.x(), a.y()-b.y(), a.z()-b.z());
509 }
510
511 /**
512 * Multiplication vector by scalar.
513 * @relates BasicVector3D
514 */
515 inline BasicVector3D<double>
516 operator*(const BasicVector3D<double> & v, double a) {
517 return BasicVector3D<double>(v.x()*a, v.y()*a, v.z()*a);
518 }
519
520 /**
521 * Scalar product of two vectors.
522 * @relates BasicVector3D
523 */
524 inline double
525 operator*(const BasicVector3D<double> & a,const BasicVector3D<double> & b) {
526 return a.dot(b);
527 }
528
529 /**
530 * Multiplication scalar by vector.
531 * @relates BasicVector3D
532 */
533 inline BasicVector3D<double>
534 operator*(double a, const BasicVector3D<double> & v) {
535 return BasicVector3D<double>(a*v.x(), a*v.y(), a*v.z());
536 }
537
538 /**
539 * Division vector by scalar.
540 * @relates BasicVector3D
541 */
542 inline BasicVector3D<double>
543 operator/(const BasicVector3D<double> & v, double a) {
544 return BasicVector3D<double>(v.x()/a, v.y()/a, v.z()/a);
545 }
546
547 /**
548 * Comparison of two vectors for equality.
549 * @relates BasicVector3D
550 */
551 inline bool
552 operator==(const BasicVector3D<double> & a, const BasicVector3D<double> & b)
553 {
554 return (a.x()==b.x() && a.y()==b.y() && a.z()==b.z());
555 }
556
557 /**
558 * Comparison of two vectors for inequality.
559 * @relates BasicVector3D
560 */
561 inline bool
562 operator!=(const BasicVector3D<double> & a, const BasicVector3D<double> & b)
563 {
564 return (a.x()!=b.x() || a.y()!=b.y() || a.z()!=b.z());
565 }
566 } /* namespace HepGeom */
567
568 #endif /* BASIC_VECTOR3D_H */
569