Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/externals/clhep/src/ThreeVector.cc

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

Diff markup

Differences between /externals/clhep/src/ThreeVector.cc (Version 11.3.0) and /externals/clhep/src/ThreeVector.cc (Version 9.5.p1)


  1 // -*- C++ -*-                                      1 // -*- C++ -*-
                                                   >>   2 // $Id:$
  2 // -------------------------------------------      3 // ---------------------------------------------------------------------------
  3 //                                                  4 //
  4 // This file is a part of the CLHEP - a Class       5 // This file is a part of the CLHEP - a Class Library for High Energy Physics.
  5 //                                                  6 //
  6 // This is the implementation of the Hep3Vecto      7 // This is the implementation of the Hep3Vector class.
  7 //                                                  8 //
  8 // See also ThreeVectorR.cc for implementation      9 // See also ThreeVectorR.cc for implementation of Hep3Vector methods which 
  9 // would couple in all the HepRotation methods     10 // would couple in all the HepRotation methods.
 10 //                                                 11 //
 11                                                    12 
                                                   >>  13 #ifdef GNUPRAGMA
                                                   >>  14 #pragma implementation
                                                   >>  15 #endif
                                                   >>  16 
 12 #include "CLHEP/Vector/ThreeVector.h"              17 #include "CLHEP/Vector/ThreeVector.h"
 13 #include "CLHEP/Units/PhysicalConstants.h"         18 #include "CLHEP/Units/PhysicalConstants.h"
 14                                                    19 
 15 #include <cmath>                                   20 #include <cmath>
 16 #include <iostream>                                21 #include <iostream>
 17                                                    22 
 18 namespace CLHEP  {                                 23 namespace CLHEP  {
 19                                                    24 
 20 void Hep3Vector::setMag(double ma) {               25 void Hep3Vector::setMag(double ma) {
 21   double factor = mag();                           26   double factor = mag();
 22   if (factor == 0) {                               27   if (factor == 0) {
 23     std::cerr << "Hep3Vector::setMag() - "         28     std::cerr << "Hep3Vector::setMag() - "
 24               << "zero vector can't be stretch     29               << "zero vector can't be stretched" << std::endl;
 25   }else{                                           30   }else{
 26     factor = ma/factor;                            31     factor = ma/factor;
 27     setX(x()*factor);                              32     setX(x()*factor);
 28     setY(y()*factor);                              33     setY(y()*factor);
 29     setZ(z()*factor);                              34     setZ(z()*factor);
 30   }                                                35   }
 31 }                                                  36 }
 32                                                    37 
                                                   >>  38 double Hep3Vector::operator () (int i) const {
                                                   >>  39   switch(i) {
                                                   >>  40   case X:
                                                   >>  41     return x();
                                                   >>  42   case Y:
                                                   >>  43     return y();
                                                   >>  44   case Z:
                                                   >>  45     return z();
                                                   >>  46   default:
                                                   >>  47     std::cerr << "Hep3Vector::operator () - "
                                                   >>  48               << "Hep3Vector subscripting: bad index (" << i << ")"
                                                   >>  49               << std::endl;
                                                   >>  50   }
                                                   >>  51   return 0.;
                                                   >>  52 }
                                                   >>  53 
                                                   >>  54 double & Hep3Vector::operator () (int i) {
                                                   >>  55   static double dummy;
                                                   >>  56   switch(i) {
                                                   >>  57   case X:
                                                   >>  58     return dx;
                                                   >>  59   case Y:
                                                   >>  60     return dy;
                                                   >>  61   case Z:
                                                   >>  62     return dz;
                                                   >>  63   default:
                                                   >>  64     std::cerr
                                                   >>  65       << "Hep3Vector::operator () - "
                                                   >>  66       << "Hep3Vector subscripting: bad index (" << i << ")"
                                                   >>  67       << std::endl;
                                                   >>  68     return dummy;
                                                   >>  69   }
                                                   >>  70 }
                                                   >>  71 
 33 Hep3Vector & Hep3Vector::rotateUz(const Hep3Ve     72 Hep3Vector & Hep3Vector::rotateUz(const Hep3Vector& NewUzVector) {
 34   // NewUzVector must be normalized !              73   // NewUzVector must be normalized !
 35                                                    74 
 36   double u1 = NewUzVector.x();                     75   double u1 = NewUzVector.x();
 37   double u2 = NewUzVector.y();                     76   double u2 = NewUzVector.y();
 38   double u3 = NewUzVector.z();                     77   double u3 = NewUzVector.z();
 39   double up = u1*u1 + u2*u2;                       78   double up = u1*u1 + u2*u2;
 40                                                    79 
 41   if (up > 0) {                                <<  80   if (up>0) {
 42     up = std::sqrt(up);                        <<  81       up = std::sqrt(up);
 43     double px = (u1 * u3 * x() - u2 * y()) / u <<  82       double px = dx,  py = dy,  pz = dz;
 44     double py = (u2 * u3 * x() + u1 * y()) / u <<  83       dx = (u1*u3*px - u2*py)/up + u1*pz;
 45     double pz = -up * x() + u3 * z();          <<  84       dy = (u2*u3*px + u1*py)/up + u2*pz;
 46     set(px, py, pz);                           <<  85       dz =    -up*px +             u3*pz;
 47   } else if (u3 < 0.) {                        <<  86     }
 48     setX(-x());                                <<  87   else if (u3 < 0.) { dx = -dx; dz = -dz; }      // phi=0  teta=pi
 49     setZ(-z());                                <<  88   else {};
 50   } // phi=0  teta=pi                          << 
 51                                                << 
 52   return *this;                                    89   return *this;
 53 }                                                  90 }
 54                                                    91 
 55 double Hep3Vector::pseudoRapidity() const {        92 double Hep3Vector::pseudoRapidity() const {
 56   double m1 = mag();                           <<  93   double m = mag();
 57   if ( m1==  0   ) return  0.0;                <<  94   if ( m==  0   ) return  0.0;   
 58   if ( m1==  z() ) return  1.0E72;             <<  95   if ( m==  z() ) return  1.0E72;
 59   if ( m1== -z() ) return -1.0E72;             <<  96   if ( m== -z() ) return -1.0E72;
 60   return 0.5*std::log( (m1+z())/(m1-z()) );    <<  97   return 0.5*std::log( (m+z())/(m-z()) );
 61 }                                                  98 }
 62                                                    99 
 63 std::ostream & operator<< (std::ostream & os,     100 std::ostream & operator<< (std::ostream & os, const Hep3Vector & v) {
 64   return os << "(" << v.x() << "," << v.y() <<    101   return os << "(" << v.x() << "," << v.y() << "," << v.z() << ")";
 65 }                                                 102 }
 66                                                   103 
 67 void ZMinput3doubles ( std::istream & is, cons    104 void ZMinput3doubles ( std::istream & is, const char * type,
 68                        double & x, double & y,    105                        double & x, double & y, double & z );
 69                                                   106 
 70 std::istream & operator>>(std::istream & is, H    107 std::istream & operator>>(std::istream & is, Hep3Vector & v) {
 71   double x, y, z;                                 108   double x, y, z;
 72   ZMinput3doubles ( is, "Hep3Vector", x, y, z     109   ZMinput3doubles ( is, "Hep3Vector", x, y, z );
 73   v.set(x, y, z);                                 110   v.set(x, y, z);
 74   return  is;                                     111   return  is;
 75 }  // operator>>()                                112 }  // operator>>()
 76                                                   113 
 77 const Hep3Vector HepXHat(1.0, 0.0, 0.0);          114 const Hep3Vector HepXHat(1.0, 0.0, 0.0);
 78 const Hep3Vector HepYHat(0.0, 1.0, 0.0);          115 const Hep3Vector HepYHat(0.0, 1.0, 0.0);
 79 const Hep3Vector HepZHat(0.0, 0.0, 1.0);          116 const Hep3Vector HepZHat(0.0, 0.0, 1.0);
 80                                                   117 
 81 //-------------------                             118 //-------------------
 82 //                                                119 //
 83 // New methods introduced when ZOOM PhysicsVec    120 // New methods introduced when ZOOM PhysicsVectors was merged in:
 84 //                                                121 //
 85 //-------------------                             122 //-------------------
 86                                                   123 
 87 Hep3Vector & Hep3Vector::rotateX (double phi1) << 124 Hep3Vector & Hep3Vector::rotateX (double phi) {
 88   double sinphi = std::sin(phi1);              << 125   double sinphi = std::sin(phi);
 89   double cosphi = std::cos(phi1);              << 126   double cosphi = std::cos(phi);
 90   double ty = y() * cosphi - z() * sinphi;     << 127   double ty;
 91   double tz = z() * cosphi + y() * sinphi;     << 128   ty = dy * cosphi - dz * sinphi;
 92   setY(ty);                                    << 129   dz = dz * cosphi + dy * sinphi;
 93   setZ(tz);                                    << 130   dy = ty;
 94   return *this;                                   131   return *this;
 95 } /* rotateX */                                   132 } /* rotateX */
 96                                                   133 
 97 Hep3Vector & Hep3Vector::rotateY (double phi1) << 134 Hep3Vector & Hep3Vector::rotateY (double phi) {
 98   double sinphi = std::sin(phi1);              << 135   double sinphi = std::sin(phi);
 99   double cosphi = std::cos(phi1);              << 136   double cosphi = std::cos(phi);
100   double tx = x() * cosphi + z() * sinphi;     << 137   double tz;
101   double tz = z() * cosphi - x() * sinphi;     << 138   tz = dz * cosphi - dx * sinphi;
102   setX(tx);                                    << 139   dx = dx * cosphi + dz * sinphi;
103   setZ(tz);                                    << 140   dz = tz;
104   return *this;                                   141   return *this;
105 } /* rotateY */                                   142 } /* rotateY */
106                                                   143 
107 Hep3Vector & Hep3Vector::rotateZ (double phi1) << 144 Hep3Vector & Hep3Vector::rotateZ (double phi) {
108   double sinphi = std::sin(phi1);              << 145   double sinphi = std::sin(phi);
109   double cosphi = std::cos(phi1);              << 146   double cosphi = std::cos(phi);
110   double tx = x() * cosphi - y() * sinphi;     << 147   double tx;
111   double ty = y() * cosphi + x() * sinphi;     << 148   tx = dx * cosphi - dy * sinphi;
112   setX(tx);                                    << 149   dy = dy * cosphi + dx * sinphi;
113   setY(ty);                                    << 150   dx = tx;
114   return *this;                                   151   return *this;
115 } /* rotateZ */                                   152 } /* rotateZ */
116                                                   153 
117 bool Hep3Vector::isNear(const Hep3Vector & v,     154 bool Hep3Vector::isNear(const Hep3Vector & v, double epsilon) const {
118   double limit = dot(v)*epsilon*epsilon;          155   double limit = dot(v)*epsilon*epsilon;
119   return ( (*this - v).mag2() <= limit );         156   return ( (*this - v).mag2() <= limit );
120 } /* isNear() */                                  157 } /* isNear() */
121                                                   158 
122 double Hep3Vector::howNear(const Hep3Vector &     159 double Hep3Vector::howNear(const Hep3Vector & v ) const {
123   // | V1 - V2 | **2  / V1 dot V2, up to 1        160   // | V1 - V2 | **2  / V1 dot V2, up to 1
124   double d   = (*this - v).mag2();                161   double d   = (*this - v).mag2();
125   double vdv = dot(v);                            162   double vdv = dot(v);
126   if ( (vdv > 0) && (d < vdv)  ) {                163   if ( (vdv > 0) && (d < vdv)  ) {
127     return std::sqrt (d/vdv);                     164     return std::sqrt (d/vdv);
128   } else if ( (vdv == 0) && (d == 0) ) {          165   } else if ( (vdv == 0) && (d == 0) ) {
129     return 0;                                     166     return 0;
130   } else {                                        167   } else {
131     return 1;                                     168     return 1;
132   }                                               169   }
133 } /* howNear */                                   170 } /* howNear */
134                                                   171 
135 double Hep3Vector::deltaPhi  (const Hep3Vector    172 double Hep3Vector::deltaPhi  (const Hep3Vector & v2) const {
136   double dphi = v2.getPhi() - getPhi();           173   double dphi = v2.getPhi() - getPhi();
137   if ( dphi > CLHEP::pi ) {                       174   if ( dphi > CLHEP::pi ) {
138     dphi -= CLHEP::twopi;                         175     dphi -= CLHEP::twopi;
139   } else if ( dphi <= -CLHEP::pi ) {              176   } else if ( dphi <= -CLHEP::pi ) {
140     dphi += CLHEP::twopi;                         177     dphi += CLHEP::twopi;
141   }                                               178   }
142   return dphi;                                    179   return dphi;
143 } /* deltaPhi */                                  180 } /* deltaPhi */
144                                                   181 
145 double Hep3Vector::deltaR ( const Hep3Vector &    182 double Hep3Vector::deltaR ( const Hep3Vector & v ) const {
146   double a = eta() - v.eta();                     183   double a = eta() - v.eta();
147   double b = deltaPhi(v);                         184   double b = deltaPhi(v); 
148   return std::sqrt ( a*a + b*b );                 185   return std::sqrt ( a*a + b*b );
149 } /* deltaR */                                    186 } /* deltaR */
150                                                   187 
151 double Hep3Vector::cosTheta(const Hep3Vector &    188 double Hep3Vector::cosTheta(const Hep3Vector & q) const {
152   double arg;                                     189   double arg;
153   double ptot2 = mag2()*q.mag2();                 190   double ptot2 = mag2()*q.mag2();
154   if(ptot2 <= 0) {                                191   if(ptot2 <= 0) {
155     arg = 0.0;                                    192     arg = 0.0;
156   }else{                                          193   }else{
157     arg = dot(q)/std::sqrt(ptot2);                194     arg = dot(q)/std::sqrt(ptot2);
158     if(arg >  1.0) arg =  1.0;                    195     if(arg >  1.0) arg =  1.0;
159     if(arg < -1.0) arg = -1.0;                    196     if(arg < -1.0) arg = -1.0;
160   }                                               197   }
161   return arg;                                     198   return arg;
162 }                                                 199 }
163                                                   200 
164 double Hep3Vector::cos2Theta(const Hep3Vector     201 double Hep3Vector::cos2Theta(const Hep3Vector & q) const {
165   double arg;                                     202   double arg;
166   double ptot2 = mag2();                          203   double ptot2 = mag2();
167   double qtot2 = q.mag2();                        204   double qtot2 = q.mag2();
168   if ( ptot2 == 0 || qtot2 == 0 )  {              205   if ( ptot2 == 0 || qtot2 == 0 )  {
169     arg = 1.0;                                    206     arg = 1.0;
170   }else{                                          207   }else{
171     double pdq = dot(q);                          208     double pdq = dot(q);
172     arg = (pdq/ptot2) * (pdq/qtot2);              209     arg = (pdq/ptot2) * (pdq/qtot2);
173         // More naive methods overflow on vect    210         // More naive methods overflow on vectors which can be squared
174         // but can't be raised to the 4th powe    211         // but can't be raised to the 4th power.
175     if(arg >  1.0) arg =  1.0;                    212     if(arg >  1.0) arg =  1.0;
176  }                                                213  }
177  return arg;                                      214  return arg;
178 }                                                 215 }
179                                                   216 
180 void Hep3Vector::setEta (double eta1) {        << 217 void Hep3Vector::setEta (double eta) {
181   double phi1 = 0;                             << 218   double phi = 0;
182   double r1;                                   << 219   double r;
183   if ( (x() == 0) && (y() == 0) ) {            << 220   if ( (dx == 0) && (dy == 0) ) {
184     if (z() == 0) {                            << 221     if (dz == 0) {
185       std::cerr << "Hep3Vector::setEta() - "      222       std::cerr << "Hep3Vector::setEta() - "
186                 << "Attempt to set eta of zero    223                 << "Attempt to set eta of zero vector -- vector is unchanged"
187                 << std::endl;                     224                 << std::endl;
188       return;                                     225       return;
189     }                                             226     }
190   std::cerr << "Hep3Vector::setEta() - "          227   std::cerr << "Hep3Vector::setEta() - "
191             << "Attempt to set eta of vector a    228             << "Attempt to set eta of vector along Z axis -- will use phi = 0"
192             << std::endl;                         229             << std::endl;
193     r1 = std::fabs(z());                       << 230     r = std::fabs(dz);
194   } else {                                        231   } else {
195     r1 = getR();                               << 232     r = getR();
196     phi1 = getPhi();                           << 233     phi = getPhi();
197   }                                               234   }
198   double tanHalfTheta = std::exp ( -eta1 );    << 235   double tanHalfTheta = std::exp ( -eta );
199   double cosTheta1 =                           << 236   double cosTheta =
200         (1 - tanHalfTheta*tanHalfTheta) / (1 +    237         (1 - tanHalfTheta*tanHalfTheta) / (1 + tanHalfTheta*tanHalfTheta);
201   double rho1 = r1*std::sqrt(1 - cosTheta1*cos << 238   dz = r * cosTheta;
202   setZ(r1 * cosTheta1);                        << 239   double rho = r*std::sqrt(1 - cosTheta*cosTheta);
203   setY(rho1 * std::sin (phi1));                << 240   dy = rho * std::sin (phi);
204   setX(rho1 * std::cos (phi1));                << 241   dx = rho * std::cos (phi);
205   return;                                         242   return;
206 }                                                 243 }
207                                                   244 
208 void Hep3Vector::setCylTheta (double theta1) { << 245 void Hep3Vector::setCylTheta (double theta) {
209                                                   246 
210   // In cylindrical coords, set theta while ke    247   // In cylindrical coords, set theta while keeping rho and phi fixed
211                                                   248 
212   if ( (x() == 0) && (y() == 0) ) {            << 249   if ( (dx == 0) && (dy == 0) ) {
213     if (z() == 0) {                            << 250     if (dz == 0) {
214       std::cerr << "Hep3Vector::setCylTheta()     251       std::cerr << "Hep3Vector::setCylTheta() - "
215                 << "Attempt to set cylTheta of    252                 << "Attempt to set cylTheta of zero vector -- vector is unchanged"
216                 << std::endl;                     253                 << std::endl;
217       return;                                     254       return;
218     }                                             255     }
219     if (theta1 == 0) {                         << 256     if (theta == 0) {
220       setZ(std::fabs(z()));                    << 257       dz = std::fabs(dz);
221       return;                                     258       return;
222     }                                             259     }
223     if (theta1 == CLHEP::pi) {                 << 260     if (theta == CLHEP::pi) {
224       setZ(-std::fabs(z()));                   << 261       dz = -std::fabs(dz);
225       return;                                     262       return;
226     }                                             263     }
227     std::cerr << "Hep3Vector::setCylTheta() -     264     std::cerr << "Hep3Vector::setCylTheta() - "
228       << "Attempt set cylindrical theta of vec    265       << "Attempt set cylindrical theta of vector along Z axis "
229       << "to a non-trivial value, while keepin    266       << "to a non-trivial value, while keeping rho fixed -- "
230       << "will return zero vector" << std::end    267       << "will return zero vector" << std::endl;
231     setZ(0.0);                                 << 268     dz = 0;
232     return;                                       269     return;
233   }                                               270   }
234   if ( (theta1 < 0) || (theta1 > CLHEP::pi) )  << 271   if ( (theta < 0) || (theta > CLHEP::pi) ) {
235     std::cerr << "Hep3Vector::setCylTheta() -     272     std::cerr << "Hep3Vector::setCylTheta() - "
236       << "Setting Cyl theta of a vector based     273       << "Setting Cyl theta of a vector based on a value not in [0, PI]"
237       << std::endl;                               274       << std::endl;
238         // No special return needed if warning    275         // No special return needed if warning is ignored.
239   }                                               276   }
240   double phi1 (getPhi());                      << 277   double phi (getPhi());
241   double rho1 = getRho();                      << 278   double rho = getRho();
242   if ( (theta1 == 0) || (theta1 == CLHEP::pi)  << 279   if ( (theta == 0) || (theta == CLHEP::pi) ) {
243     std::cerr << "Hep3Vector::setCylTheta() -     280     std::cerr << "Hep3Vector::setCylTheta() - "
244       << "Attempt to set cylindrical theta to     281       << "Attempt to set cylindrical theta to 0 or PI "
245       << "while keeping rho fixed -- infinite     282       << "while keeping rho fixed -- infinite Z will be computed"
246       << std::endl;                               283       << std::endl;
247       setZ((theta1==0) ? 1.0E72 : -1.0E72);    << 284       dz = (theta==0) ? 1.0E72 : -1.0E72;
248     return;                                       285     return;
249   }                                               286   }
250   setZ(rho1 / std::tan (theta1));              << 287   dz = rho / std::tan (theta);
251   setY(rho1 * std::sin (phi1));                << 288   dy = rho * std::sin (phi);
252   setX(rho1 * std::cos (phi1));                << 289   dx = rho * std::cos (phi);
253                                                   290 
254 } /* setCylTheta */                               291 } /* setCylTheta */
255                                                   292 
256 void Hep3Vector::setCylEta (double eta1) {     << 293 void Hep3Vector::setCylEta (double eta) {
257                                                   294 
258   // In cylindrical coords, set eta while keep    295   // In cylindrical coords, set eta while keeping rho and phi fixed
259                                                   296 
260   double theta1 = 2 * std::atan ( std::exp (-e << 297   double theta = 2 * std::atan ( std::exp (-eta) );
261                                                   298 
262         //-| The remaining code is similar to     299         //-| The remaining code is similar to setCylTheta,  The reason for
263         //-| using a copy is so as to be able     300         //-| using a copy is so as to be able to change the messages in the
264         //-| ZMthrows to say eta rather than t    301         //-| ZMthrows to say eta rather than theta.  Besides, we assumedly
265         //-| need not check for theta of 0 or     302         //-| need not check for theta of 0 or PI.
266                                                   303 
267   if ( (x() == 0) && (y() == 0) ) {            << 304   if ( (dx == 0) && (dy == 0) ) {
268     if (z() == 0) {                            << 305     if (dz == 0) {
269       std::cerr << "Hep3Vector::setCylEta() -     306       std::cerr << "Hep3Vector::setCylEta() - "
270         << "Attempt to set cylEta of zero vect    307         << "Attempt to set cylEta of zero vector -- vector is unchanged"
271         << std::endl;                             308         << std::endl;
272       return;                                     309       return;
273     }                                             310     }
274     if (theta1 == 0) {                         << 311     if (theta == 0) {
275       setZ(std::fabs(z()));                    << 312       dz = std::fabs(dz);
276       return;                                     313       return;
277     }                                             314     }
278     if (theta1 == CLHEP::pi) {                 << 315     if (theta == CLHEP::pi) {
279       setZ(-std::fabs(z()));                   << 316       dz = -std::fabs(dz);
280       return;                                     317       return;
281     }                                             318     }
282     std::cerr << "Hep3Vector::setCylEta() - "     319     std::cerr << "Hep3Vector::setCylEta() - "
283       << "Attempt set cylindrical eta of vecto    320       << "Attempt set cylindrical eta of vector along Z axis "
284       << "to a non-trivial value, while keepin    321       << "to a non-trivial value, while keeping rho fixed -- "
285       << "will return zero vector" << std::end    322       << "will return zero vector" << std::endl;
286     setZ(0.0);                                 << 323     dz = 0;
287     return;                                       324     return;
288   }                                               325   }
289   double phi1 (getPhi());                      << 326   double phi (getPhi());
290   double rho1 = getRho();                      << 327   double rho = getRho();
291   setZ(rho1 / std::tan (theta1));              << 328   dz = rho / std::tan (theta);
292   setY(rho1 * std::sin (phi1));                << 329   dy = rho * std::sin (phi);
293   setX(rho1 * std::cos (phi1));                << 330   dx = rho * std::cos (phi);
294                                                   331 
295 } /* setCylEta */                                 332 } /* setCylEta */
296                                                   333 
297                                                   334 
298 Hep3Vector operator/ ( const Hep3Vector & v1,  << 335 Hep3Vector operator/  ( const Hep3Vector & v1, double c ) {
299 //  if (c == 0) {                                 336 //  if (c == 0) {
300 //    std::cerr << "Hep3Vector::operator/ () -    337 //    std::cerr << "Hep3Vector::operator/ () - "
301 //      << "Attempt to divide vector by 0 -- "    338 //      << "Attempt to divide vector by 0 -- "
302 //      << "will produce infinities and/or NAN    339 //      << "will produce infinities and/or NANs" << std::endl;
303 //  }                                             340 //  } 
304   return v1 * (1.0/c);                         << 341   double   oneOverC = 1.0/c;
                                                   >> 342   return Hep3Vector  (  v1.x() * oneOverC,
                                                   >> 343                         v1.y() * oneOverC,
                                                   >> 344                         v1.z() * oneOverC );
305 } /* v / c */                                     345 } /* v / c */
306                                                   346 
307 Hep3Vector & Hep3Vector::operator/= (double c)    347 Hep3Vector & Hep3Vector::operator/= (double c) {
308 //  if (c == 0) {                                 348 //  if (c == 0) {
309 //    std::cerr << "Hep3Vector::operator/ () -    349 //    std::cerr << "Hep3Vector::operator/ () - "
310 //      << "Attempt to do vector /= 0 -- "        350 //      << "Attempt to do vector /= 0 -- "
311 //      << "division by zero would produce inf    351 //      << "division by zero would produce infinite or NAN components"
312 //      << std::endl;                             352 //      << std::endl;
313 //  }                                             353 //  }
314   *this *= 1.0/c;                              << 354   double oneOverC = 1.0/c;
                                                   >> 355   dx *= oneOverC;
                                                   >> 356   dy *= oneOverC;
                                                   >> 357   dz *= oneOverC;
315   return *this;                                   358   return *this;
316 }                                                 359 }
317                                                   360 
318 double Hep3Vector::tolerance = Hep3Vector::Tol    361 double Hep3Vector::tolerance = Hep3Vector::ToleranceTicks * 2.22045e-16;
319                                                   362 
320 }  // namespace CLHEP                             363 }  // namespace CLHEP
321                                                   364