00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039 #include "G4LorentzConvertor.hh"
00040 #include "G4ThreeVector.hh"
00041 #include "G4HadronicException.hh"
00042 #include "G4InuclParticle.hh"
00043
00044
00045 const G4double G4LorentzConvertor::small = 1.0e-10;
00046
00047 G4LorentzConvertor::G4LorentzConvertor()
00048 : verboseLevel(0), v2(0.), ecm_tot(0.), valong(0.), degenerated(false) {}
00049
00050 G4LorentzConvertor::
00051 G4LorentzConvertor(const G4LorentzVector& bmom, G4double bmass,
00052 const G4LorentzVector& tmom, G4double tmass)
00053 : verboseLevel(0), v2(0.), ecm_tot(0.), valong(0.), degenerated(false) {
00054 setBullet(bmom, bmass);
00055 setTarget(tmom, tmass);
00056 }
00057
00058 G4LorentzConvertor::
00059 G4LorentzConvertor(const G4InuclParticle* bullet,
00060 const G4InuclParticle* target)
00061 : verboseLevel(0), v2(0.), ecm_tot(0.), valong(0.), degenerated(false) {
00062 setBullet(bullet);
00063 setTarget(target);
00064 }
00065
00066
00067 void G4LorentzConvertor::setBullet(const G4InuclParticle* bullet) {
00068 setBullet(bullet->getMomentum());
00069 }
00070
00071 void G4LorentzConvertor::setTarget(const G4InuclParticle* target) {
00072 setTarget(target->getMomentum());
00073 }
00074
00075
00076
00077
00078 void G4LorentzConvertor::toTheCenterOfMass() {
00079 if (verboseLevel > 2)
00080 G4cout << " >>> G4LorentzConvertor::toTheCenterOfMass" << G4endl;
00081
00082 velocity = (target_mom+bullet_mom).boostVector();
00083 if (verboseLevel > 3) G4cout << " boost " << velocity << G4endl;
00084
00085
00086 scm_momentum = target_mom;
00087 scm_momentum.boost(-velocity);
00088 scm_momentum.setVect(-scm_momentum.vect());
00089
00090 if (verboseLevel > 3) G4cout << " pscm " << scm_momentum.vect() << G4endl;
00091
00092 fillKinematics();
00093 }
00094
00095 void G4LorentzConvertor::toTheTargetRestFrame() {
00096 if (verboseLevel > 2)
00097 G4cout << " >>> G4LorentzConvertor::toTheTargetRestFrame" << G4endl;
00098
00099 velocity = target_mom.boostVector();
00100 if (verboseLevel > 3) G4cout << " boost " << velocity << G4endl;
00101
00102
00103 scm_momentum = bullet_mom;
00104 scm_momentum.boost(-velocity);
00105
00106 if (verboseLevel > 3) G4cout << " pseudo-pscm " << scm_momentum.vect() << G4endl;
00107
00108 fillKinematics();
00109 }
00110
00111
00112
00113 void G4LorentzConvertor::fillKinematics() {
00114 ecm_tot = (target_mom+bullet_mom).m();
00115
00116 scm_direction = scm_momentum.vect().unit();
00117 valong = velocity.dot(scm_direction);
00118
00119 v2 = velocity.mag2();
00120
00121 G4double pvsq = v2 - valong*valong;
00122 if (verboseLevel > 3) G4cout << " pvsq " << pvsq << G4endl;
00123
00124 degenerated = (pvsq < small);
00125 if (degenerated && verboseLevel > 2)
00126 G4cout << " degenerated case (already along Z) " << G4endl;
00127
00128 if (verboseLevel > 3) {
00129 G4cout << " v2 " << v2 << " valong " << valong
00130 << " valong*valong " << valong*valong << G4endl;
00131 }
00132 }
00133
00134 G4LorentzVector
00135 G4LorentzConvertor::backToTheLab(const G4LorentzVector& mom) const {
00136 if (verboseLevel > 2)
00137 G4cout << " >>> G4LorentzConvertor::backToTheLab" << G4endl;
00138
00139 if (verboseLevel > 3)
00140 G4cout << " at rest: px " << mom.x() << " py " << mom.y() << " pz "
00141 << mom.z() << " e " << mom.e() << G4endl
00142 << " v2 " << v2 << G4endl;
00143
00144 G4LorentzVector mom1 = mom;
00145 if (v2 > small) mom1.boost(velocity);
00146
00147 if (verboseLevel > 3)
00148 G4cout << " at lab: px " << mom1.x() << " py " << mom1.y() << " pz "
00149 << mom1.z() << G4endl;
00150
00151 return mom1;
00152 }
00153
00154
00155
00156
00157 G4double G4LorentzConvertor::getKinEnergyInTheTRS() const {
00158 if (verboseLevel > 2)
00159 G4cout << " >>> G4LorentzConvertor::getKinEnergyInTheTRS" << G4endl;
00160
00161 G4LorentzVector bmom = bullet_mom;
00162 bmom.boost(-target_mom.boostVector());
00163 return bmom.e()-bmom.m();
00164 }
00165
00166 G4double G4LorentzConvertor::getTRSMomentum() const {
00167 if (verboseLevel > 2)
00168 G4cout << " >>> G4LorentzConvertor::getTRSMomentum" << G4endl;
00169
00170 G4LorentzVector bmom = bullet_mom;
00171 bmom.boost(-target_mom.boostVector());
00172 return bmom.rho();
00173 }
00174
00175 G4LorentzVector G4LorentzConvertor::rotate(const G4LorentzVector& mom) const {
00176 if (verboseLevel > 2)
00177 G4cout << " >>> G4LorentzConvertor::rotate(G4LorentzVector)" << G4endl;
00178
00179 if (verboseLevel > 3) {
00180 G4cout << " valong " << valong << " degenerated " << degenerated << G4endl
00181 << " before rotation: px " << mom.x() << " py " << mom.y()
00182 << " pz " << mom.z() << G4endl;
00183 }
00184
00185 G4LorentzVector mom_rot = mom;
00186 if (!degenerated) {
00187 if (verboseLevel > 2)
00188 G4cout << " rotating to align with reference z axis " << G4endl;
00189
00190 G4ThreeVector vscm = velocity - valong*scm_direction;
00191 G4ThreeVector vxcm = scm_direction.cross(velocity);
00192
00193 if (vscm.mag() > small && vxcm.mag() > small) {
00194 if (verboseLevel > 3) {
00195 G4cout << " reference z axis " << scm_direction
00196 << " vscm " << vscm << " vxcm " << vxcm << G4endl;
00197 }
00198
00199 mom_rot.setVect(mom.x()*vscm.unit() + mom.y()*vxcm.unit() +
00200 mom.z()*scm_direction);
00201 } else {
00202 if (verboseLevel)
00203 G4cerr << ">>> G4LorentzVector::rotate zero with !degenerated" << G4endl;
00204 }
00205 }
00206
00207 if (verboseLevel > 3) {
00208 G4cout << " after rotation: px " << mom_rot.x() << " py " << mom_rot.y()
00209 << " pz " << mom_rot.z() << G4endl;
00210 }
00211
00212 return mom_rot;
00213 }
00214
00215 G4LorentzVector G4LorentzConvertor::rotate(const G4LorentzVector& mom1,
00216 const G4LorentzVector& mom) const {
00217 if (verboseLevel > 2)
00218 G4cout << " >>> G4LorentzConvertor::rotate(G4LorentzVector,G4LorentzVector)"
00219 << G4endl;
00220
00221 if (verboseLevel > 3) {
00222 G4cout << " before rotation: px " << mom.x() << " py " << mom.y()
00223 << " pz " << mom.z() << G4endl;
00224 }
00225
00226 G4ThreeVector mom1_dir = mom1.vect().unit();
00227 G4double pv = velocity.dot(mom1_dir);
00228
00229 G4double vperp = v2 - pv*pv;
00230 if (verboseLevel > 3) {
00231 G4cout << " vperp " << vperp << " small? " << (vperp <= small) << G4endl;
00232 }
00233
00234 G4LorentzVector mom_rot = mom;
00235
00236 if (vperp > small) {
00237 if (verboseLevel > 2)
00238 G4cout << " rotating to align with first z axis " << G4endl;
00239
00240 G4ThreeVector vmom1 = velocity - pv*mom1_dir;
00241 G4ThreeVector vxm1 = mom1_dir.cross(velocity);
00242
00243 if (vmom1.mag() > small && vxm1.mag() > small) {
00244 if (verboseLevel > 3) {
00245 G4cout << " first z axis " << mom1_dir << G4endl
00246 << " vmom1 " << vmom1 << " vxm1 " << vxm1 << G4endl;
00247 }
00248
00249 mom_rot.setVect(mom.x()*vmom1.unit() + mom.y()*vxm1.unit() +
00250 mom.z()*mom1_dir );
00251 } else {
00252 if (verboseLevel)
00253 G4cerr << ">>> G4LorentzVector::rotate zero with !degenerated" << G4endl;
00254 }
00255 }
00256
00257 if (verboseLevel > 3) {
00258 G4cout << " after rotation: px " << mom_rot.x() << " py " << mom_rot.y()
00259 << " pz " << mom_rot.z() << G4endl;
00260 }
00261
00262 return mom_rot;
00263 }
00264
00265 G4bool G4LorentzConvertor::reflectionNeeded() const {
00266 if (verboseLevel > 2)
00267 G4cout << " >>> G4LorentzConvertor::reflectionNeeded (query)" << G4endl;
00268
00269 if (verboseLevel > 3) {
00270 G4cout << " v2 = " << v2 << " SCM z = " << scm_momentum.z()
00271 << " degenerated? " << degenerated << G4endl;
00272 }
00273
00274 if (v2 < small && !degenerated)
00275 throw G4HadronicException(__FILE__, __LINE__, "G4LorentzConvertor::reflectionNeeded - return value undefined");
00276
00277 if (verboseLevel > 2) {
00278 G4cout << " reflection across XY is"
00279 << ((v2>=small && (!degenerated || scm_momentum.z()<0.0))?"":" NOT")
00280 << " needed" << G4endl;
00281 }
00282
00283 return (v2>=small && (!degenerated || scm_momentum.z()<0.0));
00284 }
00285
00286
00287
00288
00289 void G4LorentzConvertor::printBullet() const {
00290 G4cout << " G4LC bullet: px " << bullet_mom.px() << " py " << bullet_mom.py()
00291 << " pz " << bullet_mom.pz() << " e " << bullet_mom.e()
00292 << " mass " << bullet_mom.m() << G4endl;
00293 }
00294
00295 void G4LorentzConvertor::printTarget() const {
00296 G4cout << " G4LC target: px " << target_mom.px() << " py " << target_mom.py()
00297 << " pz " << target_mom.pz() << " e " << target_mom.e()
00298 << " mass " << target_mom.m() << G4endl;
00299 }
00300