00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011 #ifdef GNUPRAGMA
00012 #pragma implementation
00013 #endif
00014
00015 #include "CLHEP/Vector/LorentzVector.h"
00016
00017 #include <cmath>
00018
00019 namespace CLHEP {
00020
00021
00022
00023
00024
00025 ZMpvMetric_t HepLorentzVector::setMetric( ZMpvMetric_t met ) {
00026 ZMpvMetric_t oldMetric = (metric > 0) ? TimePositive : TimeNegative;
00027 if ( met == TimeNegative ) {
00028 metric = -1.0;
00029 } else {
00030 metric = 1.0;
00031 }
00032 return oldMetric;
00033 }
00034
00035 ZMpvMetric_t HepLorentzVector::getMetric() {
00036 return ( (metric > 0) ? TimePositive : TimeNegative );
00037 }
00038
00039
00040
00041
00042
00043
00044 double HepLorentzVector::plus (const Hep3Vector & ref) const {
00045 double r = ref.mag();
00046 if (r == 0) {
00047 std::cerr << "HepLorentzVector::plus() - "
00048 << "A zero vector used as reference to LorentzVector plus-part"
00049 << std::endl;
00050 return ee;
00051 }
00052 return ee + pp.dot(ref)/r;
00053 }
00054
00055 double HepLorentzVector::minus (const Hep3Vector & ref) const {
00056 double r = ref.mag();
00057 if (r == 0) {
00058 std::cerr << "HepLorentzVector::minus() - "
00059 << "A zero vector used as reference to LorentzVector minus-part"
00060 << std::endl;
00061 return ee;
00062 }
00063 return ee - pp.dot(ref)/r;
00064 }
00065
00066 HepLorentzVector HepLorentzVector::rest4Vector() const {
00067 return HepLorentzVector (0, 0, 0, (t() < 0.0 ? -m() : m()));
00068 }
00069
00070
00071
00072
00073
00074
00075 double HepLorentzVector::beta() const {
00076 if (ee == 0) {
00077 if (pp.mag2() == 0) {
00078 return 0;
00079 } else {
00080 std::cerr << "HepLorentzVector::beta() - "
00081 << "beta computed for HepLorentzVector with t=0 -- infinite result"
00082 << std::endl;
00083 return 1./ee;
00084 }
00085 }
00086
00087
00088
00089
00090
00091 return std::sqrt (pp.mag2() / (ee*ee)) ;
00092 }
00093
00094 double HepLorentzVector::gamma() const {
00095 double v2 = pp.mag2();
00096 double t2 = ee*ee;
00097 if (ee == 0) {
00098 if (pp.mag2() == 0) {
00099 return 1;
00100 } else {
00101 std::cerr << "HepLorentzVector::gamma() - "
00102 << "gamma computed for HepLorentzVector with t=0 -- zero result"
00103 << std::endl;
00104 return 0;
00105 }
00106 }
00107 if (t2 < v2) {
00108 std::cerr << "HepLorentzVector::gamma() - "
00109 << "gamma computed for a spacelike HepLorentzVector -- imaginary result"
00110 << std::endl;
00111
00112 return 0;
00113
00114
00115
00116
00117 }
00118 return 1./std::sqrt(1. - v2/t2 );
00119 }
00120
00121
00122
00123
00124
00125
00126
00127
00128 double HepLorentzVector::rapidity() const {
00129 register double z1 = pp.getZ();
00130
00131
00132
00133
00134
00135 if (std::fabs(ee) < std::fabs(z1)) {
00136 std::cerr << "HepLorentzVector::rapidity() - "
00137 << "rapidity for spacelike 4-vector with |E| < |Pz| -- undefined"
00138 << std::endl;
00139 return 0;
00140 }
00141 double q = (ee + z1) / (ee - z1);
00142
00143
00144 return .5 * std::log(q);
00145 }
00146
00147 double HepLorentzVector::rapidity(const Hep3Vector & ref) const {
00148 register double r = ref.mag2();
00149 if (r == 0) {
00150 std::cerr << "HepLorentzVector::rapidity() - "
00151 << "A zero vector used as reference to LorentzVector rapidity"
00152 << std::endl;
00153 return 0;
00154 }
00155 register double vdotu = pp.dot(ref)/std::sqrt(r);
00156
00157
00158
00159
00160
00161 if (std::fabs(ee) < std::fabs(vdotu)) {
00162 std::cerr << "HepLorentzVector::rapidity() - "
00163 << "rapidity for spacelike 4-vector with |E| < |P*ref| -- undefined "
00164 << std::endl;
00165 return 0;
00166 }
00167 double q = (ee + vdotu) / (ee - vdotu);
00168 return .5 * std::log(q);
00169 }
00170
00171 double HepLorentzVector::coLinearRapidity() const {
00172 register double v1 = pp.mag();
00173
00174
00175
00176
00177
00178 if (std::fabs(ee) < std::fabs(v1)) {
00179 std::cerr << "HepLorentzVector::coLinearRapidity() - "
00180 << "co-linear rapidity for spacelike 4-vector -- undefined"
00181 << std::endl;
00182 return 0;
00183 }
00184 double q = (ee + v1) / (ee - v1);
00185 return .5 * std::log(q);
00186 }
00187
00188
00189
00190
00191
00192 double HepLorentzVector::invariantMass(const HepLorentzVector & w) const {
00193 double m1 = invariantMass2(w);
00194 if (m1 < 0) {
00195
00196 if ( ee * w.ee < 0 ) {
00197 std::cerr << "HepLorentzVector::invariantMass() - "
00198 << "invariant mass meaningless: \n"
00199 << "a negative-mass input led to spacelike 4-vector sum" << std::endl;
00200 return 0;
00201 } else if ( (isSpacelike() && !isLightlike()) ||
00202 (w.isSpacelike() && !w.isLightlike()) ) {
00203 std::cerr << "HepLorentzVector::invariantMass() - "
00204 << "invariant mass meaningless because of spacelike input"
00205 << std::endl;
00206 return 0;
00207 } else {
00208
00209
00210
00211
00212
00213 return 0;
00214 }
00215 }
00216 return (ee+w.ee >=0 ) ? std::sqrt(m1) : - std::sqrt(m1);
00217 }
00218
00219
00220
00221
00222
00223 Hep3Vector HepLorentzVector::findBoostToCM() const {
00224 return -boostVector();
00225 }
00226
00227 Hep3Vector HepLorentzVector::findBoostToCM (const HepLorentzVector & w) const {
00228 double t1 = ee + w.ee;
00229 Hep3Vector v1 = pp + w.pp;
00230 if (t1 == 0) {
00231 if (v1.mag2() == 0) {
00232 return Hep3Vector(0,0,0);
00233 } else {
00234 std::cerr << "HepLorentzVector::findBoostToCM() - "
00235 << "boostToCM computed for two 4-vectors with combined t=0 -- "
00236 << "infinite result" << std::endl;
00237 return Hep3Vector(v1*(1./t1));
00238 }
00239 }
00240
00241
00242
00243
00244
00245
00246 return Hep3Vector(v1 * (-1./t1));
00247 }
00248
00249 }
00250