| 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 // This is the implementation of the Hep2Vector class.
7 //
8 //-------------------------------------------------------------
9
10 #include "CLHEP/Vector/TwoVector.h"
11 #include "CLHEP/Vector/ThreeVector.h"
12
13 #include <cmath>
14 #include <iostream>
15
16 namespace CLHEP {
17
18 double Hep2Vector::tolerance = Hep2Vector::ZMpvToleranceTicks * 2.22045e-16;
19
20 double Hep2Vector::setTolerance (double tol) {
21 // Set the tolerance for Hep2Vectors to be considered near one another
22 double oldTolerance (tolerance);
23 tolerance = tol;
24 return oldTolerance;
25 }
26
27 double Hep2Vector::operator () (int i) const {
28 if (i == 0) {
29 return x();
30 }else if (i == 1) {
31 return y();
32 }else{
33 // std::cerr << "Hep2Vector::operator () - "
34 // << "Hep2Vector::operator(): bad index" << std::endl;
35 return 0.0;
36 }
37 }
38
39 double & Hep2Vector::operator () (int i) {
40 static double dummy;
41 switch(i) {
42 case X:
43 return dx;
44 case Y:
45 return dy;
46 default:
47 // std::cerr << "Hep2Vector::operator () - "
48 // << "Hep2Vector::operator() : bad index" << std::endl;
49 return dummy;
50 }
51 }
52
53 void Hep2Vector::rotate(double aangle) {
54 double ss = std::sin(aangle);
55 double cc = std::cos(aangle);
56 double xx = dx;
57 dx = cc*xx - ss*dy;
58 dy = ss*xx + cc*dy;
59 }
60
61 Hep2Vector operator/ (const Hep2Vector & p, double a) {
62 // if (a==0) {
63 // std::cerr << "Hep2Vector operator/ () - "
64 // << "Division of Hep2Vector by zero" << std::endl;
65 // }
66 return Hep2Vector(p.x()/a, p.y()/a);
67 }
68
69 std::ostream & operator << (std::ostream & os, const Hep2Vector & q) {
70 os << "(" << q.x() << ", " << q.y() << ")";
71 return os;
72 }
73
74 void ZMinput2doubles ( std::istream & is, const char * type,
75 double & x, double & y );
76
77 std::istream & operator>>(std::istream & is, Hep2Vector & p) {
78 double x, y;
79 ZMinput2doubles ( is, "Hep2Vector", x, y );
80 p.set(x, y);
81 return is;
82 } // operator>>()
83
84 Hep2Vector::operator Hep3Vector () const {
85 return Hep3Vector ( dx, dy, 0.0 );
86 }
87
88 int Hep2Vector::compare (const Hep2Vector & v) const {
89 if ( dy > v.dy ) {
90 return 1;
91 } else if ( dy < v.dy ) {
92 return -1;
93 } else if ( dx > v.dx ) {
94 return 1;
95 } else if ( dx < v.dx ) {
96 return -1;
97 } else {
98 return 0;
99 }
100 } /* Compare */
101
102
103 bool Hep2Vector::operator > (const Hep2Vector & v) const {
104 return (compare(v) > 0);
105 }
106 bool Hep2Vector::operator < (const Hep2Vector & v) const {
107 return (compare(v) < 0);
108 }
109 bool Hep2Vector::operator>= (const Hep2Vector & v) const {
110 return (compare(v) >= 0);
111 }
112 bool Hep2Vector::operator<= (const Hep2Vector & v) const {
113 return (compare(v) <= 0);
114 }
115
116 bool Hep2Vector::isNear(const Hep2Vector & p, double epsilon) const {
117 double limit = dot(p)*epsilon*epsilon;
118 return ( (*this - p).mag2() <= limit );
119 } /* isNear() */
120
121 double Hep2Vector::howNear(const Hep2Vector & p ) const {
122 double d = (*this - p).mag2();
123 double pdp = dot(p);
124 if ( (pdp > 0) && (d < pdp) ) {
125 return std::sqrt (d/pdp);
126 } else if ( (pdp == 0) && (d == 0) ) {
127 return 0;
128 } else {
129 return 1;
130 }
131 } /* howNear */
132
133 double Hep2Vector::howParallel (const Hep2Vector & v) const {
134 // | V1 x V2 | / | V1 dot V2 |
135 // Of course, the "cross product" is fictitious but the math is valid
136 double v1v2 = std::fabs(dot(v));
137 if ( v1v2 == 0 ) {
138 // Zero is parallel to no other vector except for zero.
139 return ( (mag2() == 0) && (v.mag2() == 0) ) ? 0 : 1;
140 }
141 double abscross = std::fabs ( dx * v.y() - dy - v.x() );
142 if ( abscross >= v1v2 ) {
143 return 1;
144 } else {
145 return abscross/v1v2;
146 }
147 } /* howParallel() */
148
149 bool Hep2Vector::isParallel (const Hep2Vector & v,
150 double epsilon) const {
151 // | V1 x V2 | <= epsilon * | V1 dot V2 |
152 // Of course, the "cross product" is fictitious but the math is valid
153 double v1v2 = std::fabs(dot(v));
154 if ( v1v2 == 0 ) {
155 // Zero is parallel to no other vector except for zero.
156 return ( (mag2() == 0) && (v.mag2() == 0) );
157 }
158 double abscross = std::fabs ( dx * v.y() - dy - v.x() );
159 return ( abscross <= epsilon * v1v2 );
160 } /* isParallel() */
161
162 double Hep2Vector::howOrthogonal (const Hep2Vector & v) const {
163 // | V1 dot V2 | / | V1 x V2 |
164 // Of course, the "cross product" is fictitious but the math is valid
165 double v1v2 = std::fabs(dot(v));
166 if ( v1v2 == 0 ) {
167 return 0; // Even if one or both are 0, they are considered orthogonal
168 }
169 double abscross = std::fabs ( dx * v.y() - dy - v.x() );
170 if ( v1v2 >= abscross ) {
171 return 1;
172 } else {
173 return v1v2/abscross;
174 }
175 } /* howOrthogonal() */
176
177 bool Hep2Vector::isOrthogonal (const Hep2Vector & v,
178 double epsilon) const {
179 // | V1 dot V2 | <= epsilon * | V1 x V2 |
180 // Of course, the "cross product" is fictitious but the math is valid
181 double v1v2 = std::fabs(dot(v));
182 double abscross = std::fabs ( dx * v.y() - dy - v.x() );
183 return ( v1v2 <= epsilon * abscross );
184 } /* isOrthogonal() */
185
186 } // namespace CLHEP
187