1============================================================================ 2 3 4Minor Issues arising when merging ZOOM PhysicsVectors into CLHEP/Vector 5 (and their resolution) 6 7 ---------- 8 Hep3Vector 9 ---------- 10 11 12 136 - Spherical coordinate setting: 14 15Even though the direct spherical coordinate constructors are no longer 16in the main Hep3Vector class, we want to make it possible for the user 17to use the set() methods in spherical coordinates. But (see item 4.1) 18the keywords DERGREES RADIANS ETA are not available in Hep3Vector, so 19instead, we provide a few new methods, for example: 20setRThetaPhi(double r, double theta, double phi); 21In these, alll angles are always treated as measured in RADIANS. 22 23The full set of original signatures for set are still implemented via 24the SpaceVector.h header for compatibility with ZOOM use. 25 267 - Use of UnitVector 27 28Some code at Fermilab uses UnitVector, and I have also found some 29situations where the efficiency can make a significant difference. 30I provide UnitVector as a class for ZOOM users, in the zmpv namespace, 31via a file UnitVector.h. This no longer inherits from SpaceVector; that 32avoids the design flaw of "specialization by non-virtual inheritance." 33 34The entire implementation is inline (by casting to Hep3Vector to get its 35implementations) so no additional code need go into the CLHEP library. 36 37Since UnitVector is not in CLHEP, of necessity no CLHEP classes (in particular 38not Hep3Vector and not HepRotation) will be aware of the existance of 39UnitVector. 40 415.1 X_HAT Y_HAT Z_HAT: 42 43These were in the original physicsVectors package as UnitVectors. 44They remain so. Pure CLHEP users do have the Hep3Vectors 45HepXHat, HepYHat, HepZHat if they need such vectors. 46Note, by the way, that X_HAT, Y_HAT, Z_HAT are static to each code unit 47that includes UnitVector.h. This is safer from the viewpoint of subtle 48order-of-initialization issues. Coders should not use HepXHat in 49initialization code that will execute prior to main() because it is not 50guaranteed that they will be initialized before they are used. (This 51"should not" applies equally to original CLHEP.) 52 5311 - Treatment of ambiguous angles 54 55CLHEP often, when faced with an ambiguous angle such as the cosine of the 56angle between vector v and a zero vector, simply returns 1. ZOOM had 57issued a ZMthrow, and returned 1. This was OK because a user could set 58up to ignore the ZMthorw. I retain the original CLHEP behavior on those 59routines already in CLHEP. For some that are not originally in CLHEP, 60such as cos2Theta, we leave the ZOOM behavior intact. 61 6212 - The names angle and theta 63 64These are synonmous. So why not drop one? Well, even if I said let's 65use the CLHEP form, CLHEP uses "theta" for the angle with respect to the 66Z axis, and "angle" for angle wrt an arbitrary direction. Since in ZOOM 67these concepts are somewhat unified, we we simply stick with angle and 68theta both being accepted. 69 7013 - R, rho nomenclature confusion 71 72In ThreeVector, CLHEP, --IN COMMENTS ONLY-- refers to the radius in spherical 73coordinates as "rho", and the distance from the Z axis in cylindrical 74coordinates as R. ZOOM uses the oppopsite names, and in fact uses these 75names in method names. 76 77I did not change the ZOOM method names. I have altered the CLHEP 78comments to make the nomenclature consistent; the meanings are still 79clear from the comments. For example: 80 inline double perp() const; 81 // The transverse component (rho in cylindrical coordinate system). 82 8314 - Signatures of Hep3Vector::rotate 84 85For equivalent rotate() methods, ZOOM takes (axis, delta) where 86CLHEP takes (delta, axis). There is no harm in leaving this alone: 87both signatures will be accepted. 88 89It is amazing that there were not more cases of identical routiens be with 90different choices for order ofthe arguments. 91 9215 - An inefficiency in Hep3Vector::rotate(delta, axis) 93 94CLHEP has implemented a rotate(delta, axis) method on a 3-vector in a very 95inefficent way, first forming an identity Rotation, then rotating that 96by axis and angle (in itself quite a task) then douing vector *= rotation. 97 98For now, I leave the CLHEP code alone -- people are of course free to use 99the ZOOM-originated method with signature rotate(axis,delta), which I 100believe will be faster. 101 10216 - Return types for rotateX, rotateY, rotateZ 103 104CLHEP and PhysicsVectors each have these three methods, and they are 105identical except that the ZOOM version returns a reference to *this, 106while in CLHEP they return void. 107 108Having methods that alter an object return a reference to that object 109is convenient for certain chained operations, and costs nothing. 110 111I don't wish to potentially break ZOOM user code for no good reason, so I 112have made these CLHEP method conform to this convention. There are a couple 113of other CLHEP methods, rotate and rotateUz, which use the void return type, 114but since these methods or signatures don't appear in the original ZOOM 115classes, this can't break any code, so I leave the void return type alone 116for those. 117 118After discussion with A.P. and others, I have modified the return types 119of other CLHEP methods which return void and would by the same reasoning 120be better returning *this. 121 122These include rotate and boost methods in LorentzVector.h. 123 124 ---------------- 125 HepLorentzVector 126 ---------------- 127 12820 - "explicit" keyword: 129 130ZOOM allows constcuting a LorentzVector from just p (with e = 0) 131and from just e (with p=0). ZOOM uses the "explicit" keyword to 132protect against possible mishaps with such construcutions. Will 133that give some compilers trouble? No: CLHEP already uses the 134"explicit" keyword (in the Matrix package). So we will retain it here. 135 13622 - boostVector 137 138In CLHEP there is no checking on whether it makes sense. It does not make 139sense if this vector is non-timelike, and leads to division by zero if this 140vector has ee=0. The ZOOM routine does check, which takes very little time. 141I think the zoom implementation is therefore better. This, by the way, puts 142boostVector into the .cc file rather than the .icc. 143 14423 - Rotation methods: 145 146We should separate methods that force the load of the Rotation class. 147For practical purposes, that implies that methods like rotate() are no 148longer inline, and that as in the ThreeVector class we separate the .cc 149files. Also, again we have the rotation methods returning 150HepLorentzVector& rather than void, so things can be chained. 151 15224 - Boost methods: 153 154We have the boost methods returning HepLorentzVector& rather than void, 155so things can be chained. Also, we feel the boost methods along an axis, 156boostZ in particular, really ought to be in the main part of the header. 157ZOOM does several checks to see that the boost vector is not tachyonic. 158However, we will forego these and use the CLHEP implementations. 159 16025 - Methods acting on containers of LorentzVectors: 161 162These are present in ZOOM, for example, 163 164 list<LorentzVector> s; 165 double m = s.invariantMass(); 166 167At least for now, we will omit this, so as not to introduce template 168complications into CLHEP. If necessary, we can put this back in the 169LorentzVector.h file in the ZOOM area that typedefs LorentzVector from 170HepLorentzVector. And if there is a tremendous desire for this from the 171CLHEP crowd, we can provide this capability be suppying a header file which 172users can optionally include. (Since these are templated global methods, 173all information can naturally go into a header file which need not be included 174if these methods are not called.) 175 17626 - Unit 4-Vectors along each axis: 177 178We provide X_HAT4, Y_HAT4, Z_HAT4, T_HAT4 as statics, in the header. 179We use the inline constructor in original CLHEP to avoid the need to pull 180in any additional code. 181 18226.5 - Et 183 184We provide a new methode, et() which is the transverse energy. 185This is defined as E sin theta, or E Pperp/P. For completeness, 186we also have et(v) the transverse energy with respect to a given 187direction: p.et(v) = p.t() * p.V().perp(v) / p.V().mag(). 188 189 ---------- 190 Hep2Vector 191 ---------- 192 19327 - New methods relative to PlaneVector: 194 195Mostly, we have taken the PlaneVector class from ZOOM, intact, as the 196Hep2Vector class in CLHEP. As usual, Hep2Vector is in the global namespace, 197while PlaneVector is (if namespaces are enabled) placed in the zmpv space. 198 199For the sake of similarity with the CLHEP Hep3Vector and HepLorentzVector 200classes, we have added several simple methods to the PlaneVector class: 201 Operator[] 202 L-value operator() and [] 203 isParallel, isOrthogonal, howParallel, howOrthogonal 204 compare(), operator < > <= >= 205 setPolar(r, phi) 206 20728 - Physical Coupling of Hep2Vector to Hep3Vector: 208 209A couple of methods of PlaneVector take SpaceVector arguments. 210 211We forward declare Hep3Vector, and any methods involving Hep3Vector are 212implemented in the .cc file, so a user using Hep2Vector need not pull in 213Hep3Vector.h. 214 21529 - Use of CMATH_NAMESPACE:: 216 217In PlaneVector.h there were many instances of code like 218 y = CMATH_NAMESPACE::sqrt(x); 219This takes advantage of ISOcxx portability, which handles all the of places 220where things like sqrt are found. CLHEP has its own congfigure mechanism, 221for this, and till ISOcxx is in place, it is best to just remove the 222CMATH_NAMESPACE accomodation. 223 224(Later we may want to put it back in when ISOcxx is in CLHEP.) 225 226 ---------- 227 UnitVector 228 ---------- 229 23031 - C-style or modern casts: 231 232We have gone with old C-style casts for the purposes of using the Hep3Vector 233implementations, where needed, in the UnitVector class. In particular we 234use syntax like 235 double phi() const {return ((Hep3Vector*)this)->phi();} 236 237This maximizes the number of compilers which will handle the code, and 238minimizes the chance of some compilers using the conversion operator and 239doing an unnecessary temporary copy. 240 24134 - rotationOf global methods with UnitVector: 242 243Since there is a conversion operator from UnitVector to Hep3Vector, 244we could have omitted the global rotationOf(UnitVector, ...) methods. 245But because when such a method is applied the answer is known to be a 246UnitVector, we chose to provide methods with those signatures, returning 247HepUnit3Vector instead of Hep3Vector. These simply treat the UnitVector 248as a Hep3Vector, perform the rotation, and return the result as a unit 249vector. The only advantage is that you end up with a known UnitVector 250rather than a Hep3Vector which might later be normalized. 251 252 ------------ 253 Hep3Rotation 254 ------------ 255 25635 - Rotation Classes: 257 258Along with Rotation, ZOOM has the concepts of RotationX, RotationY, RotationZ. 259Along with LorentzTransformation, ZOOM has LorentzBoost, LorentzBoostX, 260LorentzBoostY, LorentzBoostZ. 261 262These are useful (particularly the pure LorentzBoost concept) and will become 263part of CLHEP. To keep LorentzRotation.h small, the boost classes will be in 264their own header files. 265 26637 - HEP_SHORT_NAMES: 267 268It looks like if HEP_SHORT_NAMES is defined, then Rotation will collide 269with the ZOOM Rotation class. (That is, the Rotaion.h wrapper in the 270PhysicsVectors package normally says 271 class Rotation : public HepRotation { 272but with HEP_SHORT_NAMES defined, this becomes 273 class Rotation : public Rotation { 274which won't compile. 275 276We will (for now) just say that you can't have both. 277 278We could do something fancier like have the ZOOM header check for and undef 279HEP_SHORT_NAMES when including the CLHEP Rotation.h but I will leave that for 280future consideration. 281 28239 - Header file completeness: 283 284It is necessary (so that CLHEP users can continue to rely on the header for 285finding out what methods are available) that all the methods available for 286HepRotation appear in that class' definition. No problem, since most 287are pure virtual in the base class anyway (and thus had to be present anyway.) 288 28940 - Virtual destructor in Rotation classes 290 291The Rotation and LorentzRotations classes will have virtual destructors. 292 293Since the classes have to inherit from a virtual interface (for ZOOM 294compatibility) they have Vtables anyway. For various technical reasons, 295the virtual mechanism is useful. By the way, UNLIKE the case for 3-vectors 296and 4-vectors, the overhead cost is small, and the number of distinct 297Rotation objects in our key jobs is not in the millions, so the virtual 298mechanism is not at all costly. 299 30042 - The tolerance member in Rotation and LorentzRotation: 301 302The variable tolerance is used as a default for isNear. It is shared by 303all Rotation and LorentzRotation classes. 304 305It is right to use the same tolerance for RotationInterface that it gets by 306inheritance from LTI. Providing distinct instances of this class static in 307a base and a derived class would lead to thorny dynamic-vs-static behavior 308differences. 309 31043 - LTI::decompose() method: 311 312There is a method in ZOOM to decompose a LT into a pure rotation times 313a pure boost; this needs to go into CLHEP. The decompose() method does 314not really belong in the LTI interface, even though it is a matter of getting 315information about a LT, and even though they can be defined for any sort 316of specialized LT. The reason is that the interface class should not need 317to know about two of its concrete descendants! Instead, decompose(R,B) will 318be a method of HepLorentzRotation. 319 320This is a subtle change from the ZOOM structure but I think it VERY unlikely 321that anybody was doing decompose() on an anonymous LTI. 322 32344 - Dictionary ordering: compare() for Rotations and LorentzRotations 324 325The compare() method belongs in LTI and RI. However, compare() really needs 326decompose! And if I take compare() out of the interface, then all the 327dictionary-ordering comparisons will need to be repeated in each derived 328class. 329 330Solution: There is an ADDITIONAL method 331 decompose (Hep3Vector& boost, HepAxisAngle& rotation). 332The LTI and RI interfaces do know about Hep3Vector and HepAxisAngle, so 333that is fine. And comparison is of course just as valid using these 334(in fact, that is how it was being done anyway in ZOOM, deep down). 335 336Now compare() becomes inline in the base class, and decompose() is the 337virtual method it utilizes. 338 33945 - Operator() and [][]: 340 341ZOOM did not have these. CLHEP does, so they must be present in HepRotation. 342Since these are not easy to implement with the same semantics as in CLHEP, 343for an arbitrary Hep3RotationInterface, and since no compatibility or use 344case agruments tell us they must be in RI, we will keep these only in 345HepRotation. Similarly, they will appear in HepLorentzRotation but not 346in HepBoost. 347 34848 - PhiX, thetaX, etc for RotationX, RotationY, RotationZ: 349 350These routines were added for Geant4, and that code won't be using the 351specialized rotations. Still, for completeness we should provide these. 352 353In these special cases, we could fairly easily avoid taking an atan2 or acos. 354For example, RotationX::thetaZ can just do return d. In some of the other 355cases you have to branch according to the quadrant of delta, but in no case 356do you need an inverse trig function. However, for maximal safety, I simply 357use the same code as was in CLHEP originally. Why introduce a chance for 358error to creep in? 359 360That begs the question of why not put these up into Hep3RotationInterface? 361I think eventually we may want to use the much simpler equivalent algorithms. 362By putting the methods into each class, I have been able to add comments 363containing the improved code. Plus, these methods are not so fundamental 364that I would like to see them in the interface. I'm not so sure my decision 365here is best, but I don't think it is worththe time to re-think it since 366it is such a minor point. 367 36849 - CLHEP/ZOOM "equivalent" methods for Rotations: 369 370In principle you don't want to repeat code. So for example, the rotateX 371method (in CLHEP) could be implemented as *this= RotationX(delta)* *this. 372However, except for absolutely trivial cases, I feel it unacceptable to 373take any risk of "breaking" existing CLHEP methods (in the sense of altering 374them such that they yield different results). Similarly for ZOOM methods - 375unless they have an exact CLHEP match, I prefer to use the existing ZOOM 376code. 377 378So for instance I leave the implementation of 379 getAngleAxis(double &, Hep3Vector &) 380even though it could be written as simply as 381 { HepAxisAngle z = axisAngle(); delta=a.delta; axis=a.axis; } 382 383Careful note will be made of all such "dual-implementation" situations, and 384we should take advantage by testing for equivalence in our test suites. 385 38650 - Operator *+ with a 3- or 4-vector and a 3- or 4-rotation: 387 388For example, v *= R where v is a Hep3Vector and R a HepRotation. 389In C++ semantics, z *= b usually means a = a * b, so this leads you to expect 390the behaviour v = v*R. But physicists know that they would rarely apply 391a rotation on the right like that, and thus CLHEP defined v *= R as v = R*v. 392 393This is a case where ZOOM decided there was no misleading way to treat 394the operator and thus omitted it, while CLHEP has just decided which way 395to treat it. We must retain this meaning for compatibility. We will 396point out the unnatural syntax in both the header file and the documentation, 397to warn users that they will get the useful physics operation rather than 398the expected meaning of operator *=. 399 400-------------------- 401LorentzTransformation 402-------------------- 403 40452 - set from (bx,by,bz): 405 406The merged classes have to contain a constructor of LorentzRotation from 407three boost components because CLHEP does. They must contain a set() 408method taking a boost vector because ZOOM does. They need not, a 409contain LorentzRotation::set(double,double,Hepdouble), but we include 410that because we have the constructor from three doubles. 411 41253 - transform (Rotation r): 413 414While CLHEP has a separate method for transform taking a HepRotation, this 415is not strictly needed since HepRotation is a Hep3RotationInterface which in 416turn is a Hep4RotationInterface. If you know the argument is 3x3 you could 417do 3x3 times 4x4 (in principle faster) and for transform, the original CLHEP 418would do this. On the other hand, for multiplication this "efficiency" was 419not done; there is no operator* (HepRotation). In reallity, the time needed 420for 3x3 times 4x4 is not significantly less, certainly not enough to warrant 421the additional complication. So we simply provide the transform method taking 422a 4 rotation interface, which will break no code. 423 424But for pure rotations around coordinate axes and boosts along axes, we do 425provide such methods as rotateZ and boostZ. Here the efficiency difference 426is marked. Besides, CLHEP has these explicitly. 427 428We also must keep the technically superfluous rotate(delta,axis) and 429boost(3 doubles or 3-vector) methods. 430 43154 - rotateX, Y and Z: 432 433Although these are inline in CLHEP, we make them normal now. The reason is 434that we can thus take advantage of rather big savings in algorithm (by 435utilizing the fact that the rotation only involves an X rotation or whatever). 436 43755 - ZMxpv in LorentzRotation.cc 438 439As a reminder, the ZMthrow macro is defined as outputting its argument and 440aborting. (ZMthrowC does not abort). That said, should we use ZMthrow 441(defined in ZMxpv.h) in the basic LorentzRotation class? I think we have to: 442the method set (double, double, double) -- which was present in CLHEP -- had 443a flaw in that it would take sqrt of a negative number if the arguments 444supplied had mag2 > 1. We really should behave more sensibly here ahd ZMthrow 445is the best we can do. 446 44756 - Physical coupling: 448 449BoostX needs LorentzRotation.h and so forth, because the product of 450a boost times an Hep4RotationInterface is a LorentzRotation. We could 451probably have lessened this sort of coupling but maybe not without sacrificing 452speed. At any rate, I short-circuited design though on this issue in favor 453of getting the merge done! 454 45557 - Efficiency of distance2 for boosts: 456 457We could, by providing distance2(HepLorentzRotation) and the other concrete 458classes, significantly speed up this computation, because we don't have to go 459back and forth through HepAxisAngle for the rotation decompositin and 460Hep3Vector for the boost. However, this would add another 15 methods to each 461header (assuming we would want these gains in isNear and HowNear also); 462I decided it would not be worth it. 463 464