| Geant4 Cross Reference |
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // G4FieldTrack implementation
27 //
28 // Author: John Apostolakis, CERN - First version, 14.10.1996
29 // -------------------------------------------------------------------
30
31 #include "G4FieldTrack.hh"
32
33 std::ostream& operator<<( std::ostream& os, const G4FieldTrack& SixVec)
34 {
35 const G4double* SixV = SixVec.SixVector;
36 const G4int precPos= 9; // For position
37 const G4int precEp= 9; // For Energy / momentum
38 const G4int precLen= 12; // For Length along track
39 const G4int precSpin= 9; // For polarisation
40 const G4int precTime= 6; // For time of flight
41 const G4long oldpr= os.precision(precPos);
42 os << " ( ";
43 os << " X= " << SixV[0] << " " << SixV[1] << " "
44 << SixV[2] << " "; // Position
45 os.precision(precEp);
46 os << " P= " << SixV[3] << " " << SixV[4] << " "
47 << SixV[5] << " "; // Momentum
48 os << " Pmag= "
49 << G4ThreeVector(SixV[3], SixV[4], SixV[5]).mag(); // mom magnitude
50 os << " Ekin= " << SixVec.fKineticEnergy ;
51 os.precision(precLen);
52 os << " l= " << SixVec.GetCurveLength();
53 os.precision(6);
54 os << " m0= " << SixVec.fRestMass_c2;
55 os << " (Pdir-1)= " << SixVec.fMomentumDir.mag()-1.0;
56 if( SixVec.fLabTimeOfFlight > 0.0 )
57 {
58 os.precision(precTime);
59 }
60 else
61 {
62 os.precision(3);
63 }
64 os << " t_lab= " << SixVec.fLabTimeOfFlight;
65 os << " t_proper= " << SixVec.fProperTimeOfFlight ;
66 G4ThreeVector pol= SixVec.GetPolarization();
67 if( pol.mag2() > 0.0 )
68 {
69 os.precision(precSpin);
70 os << " PolV= " << pol; // SixVec.GetPolarization();
71 }
72 else
73 {
74 os << " PolV= (0,0,0) ";
75 }
76 os << " ) ";
77 os.precision(oldpr);
78 return os;
79 }
80
81 G4FieldTrack::G4FieldTrack( const G4ThreeVector& pPosition,
82 G4double LaboratoryTimeOfFlight,
83 const G4ThreeVector& pMomentumDirection,
84 G4double kineticEnergy,
85 G4double restMass_c2,
86 G4double charge,
87 const G4ThreeVector& vecPolarization,
88 G4double magnetic_dipole_moment,
89 G4double curve_length,
90 G4double pdgSpin )
91 : fDistanceAlongCurve(curve_length),
92 fKineticEnergy(kineticEnergy),
93 fRestMass_c2(restMass_c2),
94 fLabTimeOfFlight(LaboratoryTimeOfFlight),
95 fProperTimeOfFlight(0.),
96 // fMomentumDir(pMomentumDirection),
97 fChargeState( charge, magnetic_dipole_moment, pdgSpin )
98 // fChargeState( charge, magnetic_dipole_moment ) ,
99 // fPDGSpin( pdgSpin )
100 {
101 UpdateFourMomentum( kineticEnergy, pMomentumDirection );
102 // Sets momentum direction as well.
103
104 SetPosition( pPosition );
105 SetPolarization( vecPolarization );
106 }
107
108 G4FieldTrack::G4FieldTrack( const G4ThreeVector& pPosition,
109 const G4ThreeVector& pMomentumDirection,
110 G4double curve_length,
111 G4double kineticEnergy,
112 const G4double restMass_c2,
113 G4double, // velocity
114 G4double pLaboratoryTimeOfFlight,
115 G4double pProperTimeOfFlight,
116 const G4ThreeVector* pPolarization,
117 G4double pdgSpin )
118 : fDistanceAlongCurve(curve_length),
119 fKineticEnergy(kineticEnergy),
120 fRestMass_c2(restMass_c2),
121 fLabTimeOfFlight(pLaboratoryTimeOfFlight),
122 fProperTimeOfFlight(pProperTimeOfFlight),
123 fChargeState( DBL_MAX, DBL_MAX, -1.0 ) // charge not set
124 {
125 UpdateFourMomentum( kineticEnergy, pMomentumDirection );
126 // Sets momentum direction as well.
127
128 SetPosition( pPosition );
129 fChargeState.SetPDGSpin( pdgSpin );
130
131 G4ThreeVector PolarVec(0.0, 0.0, 0.0);
132 if( pPolarization != nullptr ) { PolarVec= *pPolarization; }
133 SetPolarization( PolarVec );
134 }
135
136 G4FieldTrack::G4FieldTrack( char ) // Nothing is set !!
137 : fKineticEnergy(0.), fRestMass_c2(0.), fLabTimeOfFlight(0.),
138 fProperTimeOfFlight(0.), fChargeState( DBL_MAX , DBL_MAX, -1 )
139 {
140 G4ThreeVector Zero(0.0, 0.0, 0.0);
141 SetCurvePnt( Zero, Zero, 0.0 );
142 SetPolarization( Zero );
143 // fInitialMomentumMag = 0.00; // Invalid
144 // fLastMomentumMag = 0.0;
145 }
146
147 void G4FieldTrack::
148 SetChargeAndMoments(G4double charge,
149 G4double magnetic_dipole_moment, // default = DBL_MAX
150 G4double electric_dipole_moment, // ditto
151 G4double magnetic_charge ) // ditto
152 {
153 fChargeState.SetChargesAndMoments( charge,
154 magnetic_dipole_moment,
155 electric_dipole_moment,
156 magnetic_charge );
157
158 // NOTE: Leaves Spin unchanged !
159 //
160 // G4double pdgSpin= fChargeState.GetSpin();
161 // New Property of ChargeState (not well documented! )
162
163 // IDEA: Improve the implementation using handles
164 // -- and handle to the old one (which can be shared by other copies) and
165 // must not be left to hang loose
166 //
167 // fpChargeState= new G4ChargeState( charge, magnetic_dipole_moment,
168 // electric_dipole_moment, magnetic_charge );
169 }
170
171 // Load values from array
172 //
173 // Note that momentum direction must-be/is normalised
174 //
175 void G4FieldTrack::LoadFromArray(const G4double valArrIn[ncompSVEC],
176 G4int noVarsIntegrated)
177 {
178 // Fill the variables not integrated with zero -- so it's clear !!
179 //
180 G4double valArr[ncompSVEC];
181 for(G4int i=0; i<noVarsIntegrated; ++i)
182 {
183 valArr[i] = valArrIn[i];
184 }
185 for(G4int i=noVarsIntegrated; i<ncompSVEC; ++i)
186 {
187 valArr[i] = 0.0;
188 }
189
190 SixVector[0] = valArr[0];
191 SixVector[1] = valArr[1];
192 SixVector[2] = valArr[2];
193 SixVector[3] = valArr[3];
194 SixVector[4] = valArr[4];
195 SixVector[5] = valArr[5];
196
197 G4ThreeVector Momentum(valArr[3],valArr[4],valArr[5]);
198
199 G4double momentum_square= Momentum.mag2();
200 fMomentumDir= Momentum.unit();
201
202 fKineticEnergy = momentum_square
203 / (std::sqrt(momentum_square+fRestMass_c2*fRestMass_c2)
204 + fRestMass_c2 );
205 // The above equation is stable for small and large momenta
206
207 // The following components may or may not be
208 // integrated over -- integration is optional
209 // fKineticEnergy = valArr[6];
210
211 fLabTimeOfFlight = valArr[7];
212 fProperTimeOfFlight = valArr[8];
213 G4ThreeVector vecPolarization= G4ThreeVector(valArr[9],valArr[10],valArr[11]);
214 SetPolarization( vecPolarization );
215
216 // fMomentumDir=G4ThreeVector(valArr[13],valArr[14],valArr[15]);
217 // fDistanceAlongCurve= valArr[];
218 }
219