Geant4 10.7.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4ErrorFreeTrajState.cc
Go to the documentation of this file.
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//
27// ------------------------------------------------------------
28// GEANT 4 class implementation file
29// ------------------------------------------------------------
30//
31
32#include <iomanip>
33
34#include "G4PhysicalConstants.hh"
35#include "G4SystemOfUnits.hh"
36#include "G4Field.hh"
37#include "G4FieldManager.hh"
38#include "G4TransportationManager.hh"
39#include "G4GeometryTolerance.hh"
40#include "G4Material.hh"
41#include "G4ErrorPropagatorData.hh"
42#include "G4ErrorFreeTrajState.hh"
43#include "G4ErrorFreeTrajParam.hh"
44#include "G4ErrorSurfaceTrajState.hh"
45#include "G4ErrorMatrix.hh"
46
47//------------------------------------------------------------------------
48G4ErrorFreeTrajState::G4ErrorFreeTrajState( const G4String& partName, const G4Point3D& pos, const G4Vector3D& mom, const G4ErrorTrajErr& errmat) : G4ErrorTrajState( partName, pos, mom, errmat )
49{
50 fTrajParam = G4ErrorFreeTrajParam( pos, mom );
51 Init();
52}
53
54
55//------------------------------------------------------------------------
56G4ErrorFreeTrajState::G4ErrorFreeTrajState( const G4ErrorSurfaceTrajState& tpSD ) : G4ErrorTrajState( tpSD.GetParticleType(), tpSD.GetPosition(), tpSD.GetMomentum() )
57{
58 // G4ThreeVector planeNormal = tpSD.GetPlaneNormal();
59 // G4double fPt = tpSD.GetMomentum()*planeNormal;//mom projected on normal to plane
60 // G4ErrorSurfaceTrajParam tpSDparam = tpSD.GetParameters();
61 // G4ThreeVector Psc = fPt * planeNormal + tpSDparam.GetPU()*tpSDparam.GetVectorU() + tpSD.GetPV()*tpSD.GetVectorW();
62
63 fTrajParam = G4ErrorFreeTrajParam( fPosition, fMomentum );
64 Init();
65
66 //----- Get the error matrix in SC coordinates
67 G4ErrorSurfaceTrajParam tpSDparam = tpSD.GetParameters();
68 G4double mom = fMomentum.mag();
69 G4double mom2 = fMomentum.mag2();
70 G4double TVW1 = std::sqrt( mom2 / ( mom2 + tpSDparam.GetPV()*tpSDparam.GetPV() + tpSDparam.GetPW()*tpSDparam.GetPW()) );
71 G4ThreeVector vTVW( TVW1, tpSDparam.GetPV()/mom * TVW1, tpSDparam.GetPW()/mom * TVW1 );
72 G4Vector3D vectorU = tpSDparam.GetVectorV().cross( tpSDparam.GetVectorW() );
73 G4Vector3D vTN = vTVW.x()*vectorU + vTVW.y()*tpSDparam.GetVectorV() + vTVW.z()*tpSDparam.GetVectorW();
74
75#ifdef G4EVERBOSE
76 if( iverbose >= 5){
77 G4double pc2 = std::asin( vTN.z() );
78 G4double pc3 = std::atan (vTN.y()/vTN.x());
79
80 G4cout << " CHECK: pc2 " << pc2 << " = " << GetParameters().GetLambda() << " diff " << pc2-GetParameters().GetLambda() << G4endl;
81 G4cout << " CHECK: pc3 " << pc3 << " = " << GetParameters().GetPhi() << " diff " << pc3-GetParameters().GetPhi() << G4endl;
82 }
83#endif
84
85 //--- Get the unit vectors perp to P
86 G4double cosl = std::cos( GetParameters().GetLambda() );
87 if (cosl < 1.E-30) cosl = 1.E-30;
88 G4double cosl1 = 1./cosl;
89 G4Vector3D vUN(-vTN.y()*cosl1, vTN.x()*cosl1, 0. );
90 G4Vector3D vVN(-vTN.z()*vUN.y(), vTN.z()*vUN.x(), cosl );
91
92 G4Vector3D vUperp = G4Vector3D( -fMomentum.y(), fMomentum.x(), 0.);
93 G4Vector3D vVperp = vUperp.cross( fMomentum );
94 vUperp *= 1./vUperp.mag();
95 vVperp *= 1./vVperp.mag();
96
97#ifdef G4EVERBOSE
98 if( iverbose >= 5 ){
99 G4cout << " CHECK: vUN " << vUN << " = " << vUperp << " diff " << (vUN-vUperp).mag() << G4endl;
100 G4cout << " CHECK: vVN " << vVN << " = " << vVperp << " diff " << (vVN-vVperp).mag() << G4endl;
101 }
102#endif
103
104 //get the dot products of vectors perpendicular to direction and vector defining SD plane
105 G4double dUU = vUperp * tpSD.GetVectorV();
106 G4double dUV = vUperp * tpSD.GetVectorW();
107 G4double dVU = vVperp * tpSD.GetVectorV();
108 G4double dVV = vVperp * tpSD.GetVectorW();
109
110 //--- Get transformation first
111 G4ErrorMatrix transfM(5, 5, 1 );
112 //--- Get magnetic field
113 const G4Field* field = G4TransportationManager::GetTransportationManager()->GetFieldManager()->GetDetectorField();
114 G4ThreeVector dir = fTrajParam.GetDirection();
115 G4double invCosTheta = 1./std::cos( dir.theta() );
116 G4cout << " dir="<<dir<<" invCosTheta "<<invCosTheta << G4endl;
117
118 if( fCharge != 0 && field ) {
119 G4double pos1[3]; pos1[0] = fPosition.x()*cm; pos1[1] = fPosition.y()*cm; pos1[2] = fPosition.z()*cm;
120 G4double h1[3];
121 field->GetFieldValue( pos1, h1 );
122 G4ThreeVector HPre = G4ThreeVector( h1[0], h1[1], h1[2] ) / tesla *10.;
123 G4double magHPre = HPre.mag();
124 G4double invP = 1./fMomentum.mag();
125 G4double magHPreM = magHPre * invP;
126 if( magHPre != 0. ) {
127 G4double magHPreM2 = fCharge / magHPre;
128
129 G4double Q = -magHPreM * c_light;
130 G4double sinz = -HPre*vUperp * magHPreM2;
131 G4double cosz = HPre*vVperp * magHPreM2;
132
133 transfM[1][3] = -Q*dir.y()*sinz;
134 transfM[1][4] = -Q*dir.z()*sinz;
135 transfM[2][3] = -Q*dir.y()*cosz*invCosTheta;
136 transfM[2][4] = -Q*dir.z()*cosz*invCosTheta;
137 }
138 }
139
140 transfM[0][0] = 1.;
141 transfM[1][1] = dir.x()*dVU;
142 transfM[1][2] = dir.x()*dVV;
143 transfM[2][1] = dir.x()*dUU*invCosTheta;
144 transfM[2][2] = dir.x()*dUV*invCosTheta;
145 transfM[3][3] = dUU;
146 transfM[3][4] = dUV;
147 transfM[4][3] = dVU;
148 transfM[4][4] = dVV;
149
150 fError = G4ErrorTrajErr( tpSD.GetError().similarity( transfM ) );
151
152#ifdef G4EVERBOSE
153 if( iverbose >= 1) G4cout << "error matrix SD2SC " << fError << G4endl;
154 if( iverbose >= 4) G4cout << "G4ErrorFreeTrajState from SD " << *this << G4endl;
155#endif
156}
157
158
159//------------------------------------------------------------------------
160void G4ErrorFreeTrajState::Init()
161{
162 theTSType = G4eTS_FREE;
163 BuildCharge();
164 theTransfMat = G4ErrorMatrix(5,5,0);
165 theFirstStep = true;
166}
167
168//------------------------------------------------------------------------
169void G4ErrorFreeTrajState::Dump( std::ostream& out ) const
170{
171 out << *this;
172}
173
174//------------------------------------------------------------------------
175G4int G4ErrorFreeTrajState::Update( const G4Track* aTrack )
176{
177 G4int ierr = 0;
178 fTrajParam.Update( aTrack );
179 UpdatePosMom( aTrack->GetPosition(), aTrack->GetMomentum() );
180 return ierr;
181}
182
183
184//------------------------------------------------------------------------
185std::ostream& operator<<(std::ostream& out, const G4ErrorFreeTrajState& ts)
186{
187 std::ios::fmtflags orig_flags = out.flags();
188
189 out.setf(std::ios::fixed,std::ios::floatfield);
190
191 ts.DumpPosMomError( out );
192
193 out << " G4ErrorFreeTrajState: Params: " << ts.fTrajParam << G4endl;
194
195 out.flags(orig_flags);
196
197 return out;
198}
199
200
201//------------------------------------------------------------------------
202G4int G4ErrorFreeTrajState::PropagateError( const G4Track* aTrack )
203{
204 G4double stepLengthCm = aTrack->GetStep()->GetStepLength()/cm;
205 if( G4ErrorPropagatorData::GetErrorPropagatorData()->GetStage() == G4ErrorStage_Deflation ) stepLengthCm *= -1.;
206
207 G4double kCarTolerance = G4GeometryTolerance::GetInstance()->GetSurfaceTolerance();
208
209 if( std::fabs(stepLengthCm) <= kCarTolerance/cm ) return 0;
210
211#ifdef G4EVERBOSE
212 if( iverbose >= 2 )G4cout << " G4ErrorFreeTrajState::PropagateError " << G4endl;
213 G4cout << "G4EP: iverbose="<< iverbose << G4endl;
214#endif
215
216 // * *** ERROR PROPAGATION ON A HELIX ASSUMING SC VARIABLES
217 G4Point3D vposPost = aTrack->GetPosition()/cm;
218 G4Vector3D vpPost = aTrack->GetMomentum()/GeV;
219 // G4Point3D vposPre = fPosition/cm;
220 // G4Vector3D vpPre = fMomentum/GeV;
221 G4Point3D vposPre = aTrack->GetStep()->GetPreStepPoint()->GetPosition()/cm;
222 G4Vector3D vpPre = aTrack->GetStep()->GetPreStepPoint()->GetMomentum()/GeV;
223 //correct to avoid propagation along Z
224 if( vpPre.mag() == vpPre.z() ) vpPre.setX( 1.E-6*MeV );
225 if( vpPost.mag() == vpPost.z() ) vpPost.setX( 1.E-6*MeV );
226
227 G4double pPre = vpPre.mag();
228 G4double pPost = vpPost.mag();
229#ifdef G4EVERBOSE
230 if( iverbose >= 2 ) {
231 G4cout << "G4EP: vposPre " << vposPre << G4endl
232 << "G4EP: vposPost " << vposPost << G4endl;
233 G4cout << "G4EP: vpPre " << vpPre << G4endl
234 << "G4EP: vpPost " << vpPost << G4endl;
235 G4cout << " err start step " << fError << G4endl;
236 G4cout << "G4EP: stepLengthCm " << stepLengthCm << G4endl;
237 }
238#endif
239
240 if( pPre == 0. || pPost == 0 ) return 2;
241 G4double pInvPre = 1./pPre;
242 G4double pInvPost = 1./pPost;
243 G4double deltaPInv = pInvPost - pInvPre;
244 if( iverbose >= 2 ) G4cout << "G4EP: pInvPre" << pInvPre<< " pInvPost:" << pInvPost<<" deltaPInv:" << deltaPInv<< G4endl;
245
246
247 G4Vector3D vpPreNorm = vpPre * pInvPre;
248 G4Vector3D vpPostNorm = vpPost * pInvPost;
249 if( iverbose >= 2 ) G4cout << "G4EP: vpPreNorm " << vpPreNorm << " vpPostNorm " << vpPostNorm << G4endl;
250 //return if propagation along Z??
251 if( 1. - std::fabs(vpPreNorm.z()) < kCarTolerance ) return 4;
252 if( 1. - std::fabs(vpPostNorm.z()) < kCarTolerance ) return 4;
253 G4double sinpPre = std::sin( vpPreNorm.theta() ); //cosine perpendicular to pPre = sine pPre
254 G4double sinpPost = std::sin( vpPostNorm.theta() ); //cosine perpendicular to pPost = sine pPost
255 G4double sinpPostInv = 1./std::sin( vpPostNorm.theta() );
256
257#ifdef G4EVERBOSE
258 if( iverbose >= 2 ) G4cout << "G4EP: cosl " << sinpPre << " cosl0 " << sinpPost << G4endl;
259#endif
260 //* *** DEFINE TRANSFORMATION MATRIX BETWEEN X1 AND X2 FOR
261 //* *** NEUTRAL PARTICLE OR FIELDFREE REGION
262 G4ErrorMatrix transf(5, 5, 0 );
263
264 transf[3][2] = stepLengthCm * sinpPost;
265 transf[4][1] = stepLengthCm;
266 for( size_t ii=0;ii < 5; ii++ ){
267 transf[ii][ii] = 1.;
268 }
269#ifdef G4EVERBOSE
270 if( iverbose >= 2 ) {
271 G4cout << "G4EP: transf matrix neutral " << transf;
272 }
273#endif
274
275 // charge X propagation direction
276 G4double charge = aTrack->GetDynamicParticle()->GetCharge();
277 if( G4ErrorPropagatorData::GetErrorPropagatorData()->GetMode() == G4ErrorMode_PropBackwards ) {
278 charge *= -1.;
279 }
280 // G4cout << " charge " << charge << G4endl;
281 //t check if particle has charge
282 //t if( charge == 0 ) goto 45;
283 // check if the magnetic field is = 0.
284
285 //position is from geant4, it is assumed to be in mm (for debugging, eventually it will not be transformed)
286 //it is assumed vposPre[] is in cm and pos1[] is in mm.
287 G4double pos1[3]; pos1[0] = vposPre.x()*cm; pos1[1] = vposPre.y()*cm; pos1[2] = vposPre.z()*cm;
288 G4double pos2[3]; pos2[0] = vposPost.x()*cm; pos2[1] = vposPost.y()*cm; pos2[2] = vposPost.z()*cm;
289 G4double h1[3], h2[3];
290
291 const G4Field* field = G4TransportationManager::GetTransportationManager()->GetFieldManager()->GetDetectorField();
292 if( !field ) return 0; //goto 45
293
294
295
296 // calculate transformation except it NEUTRAL PARTICLE OR FIELDFREE REGION
297 if( charge != 0. && field ) {
298
299 field->GetFieldValue( pos1, h1 ); //here pos1[], pos2[] are in mm, not changed
300 field->GetFieldValue( pos2, h2 );
301 G4ThreeVector HPre = G4ThreeVector( h1[0], h1[1], h1[2] ) / tesla *10.; //10. is to get same dimensions as GEANT3 (kilogauss)
302 G4ThreeVector HPost= G4ThreeVector( h2[0], h2[1], h2[2] ) / tesla *10.;
303 G4double magHPre = HPre.mag();
304 G4double magHPost = HPost.mag();
305#ifdef G4EVERBOSE
306 if( iverbose >= 2 ) {
307 G4cout << "G4EP: h1 = "
308 << h1[0] << ", " << h1[1] << ", " << h1[2] << G4endl;
309 G4cout << "G4EP: pos1/mm = "
310 << pos1[0] << ", " << pos1[1] << ", " << pos1[2] << G4endl;
311 G4cout << "G4EP: pos2/mm = "
312 << pos2[0] << ", " << pos2[1] << ", " << pos2[2] << G4endl;
313 G4cout << "G4EP: B-filed in KGauss HPre " << HPre << G4endl
314 << "G4EP: in KGauss HPost " << HPost << G4endl;
315 }
316#endif
317
318 if( magHPre + magHPost != 0. ) {
319
320 //* *** CHECK WHETHER H*ALFA/P IS TOO DIFFERENT AT X1 AND X2
321 G4double gam;
322 if( magHPost != 0. ){
323 gam = HPost * vpPostNorm / magHPost;
324 }else {
325 gam = HPre * vpPreNorm / magHPre;
326 }
327
328 // G4eMagneticLimitsProcess will limit the step, but based on an straight line trajectory
329 G4double alphaSqr = 1. - gam * gam;
330 G4double diffHSqr = ( HPre * pInvPre - HPost * pInvPost ).mag2();
331 G4double delhp6Sqr = 300.*300.;
332#ifdef G4EVERBOSE
333 if( iverbose >= 2 ) {
334 G4cout << " G4EP: gam " << gam << " alphaSqr " << alphaSqr
335 << " diffHSqr " << diffHSqr << G4endl;
336 G4cout << " alpha= " << std::sqrt(alphaSqr) << G4endl;
337 }
338#endif
339 if( diffHSqr * alphaSqr > delhp6Sqr ) return 3;
340
341
342 //* *** DEFINE AVERAGE MAGNETIC FIELD AND GRADIENT
343 G4double pInvAver = 1./(pInvPre + pInvPost );
344 G4double CFACT8 = 2.997925E-4;
345 //G4double HAver
346 G4ThreeVector vHAverNorm( (HPre*pInvPre + HPost*pInvPost ) * pInvAver * charge * CFACT8 );
347 G4double HAver = vHAverNorm.mag();
348 G4double invHAver = 1./HAver;
349 vHAverNorm *= invHAver;
350#ifdef G4EVERBOSE
351 if( iverbose >= 2 ) G4cout << " G4EP: HaverNorm " << vHAverNorm << " magHAver " << HAver << " charge " << charge<< G4endl;
352#endif
353
354 G4double pAver = (pPre+pPost)*0.5;
355 G4double QAver = -HAver/pAver;
356 G4double thetaAver = QAver * stepLengthCm;
357 G4double sinThetaAver = std::sin(thetaAver);
358 G4double cosThetaAver = std::cos(thetaAver);
359 G4double gamma = vHAverNorm * vpPostNorm;
360 G4ThreeVector AN2 = vHAverNorm.cross( vpPostNorm );
361
362#ifdef G4EVERBOSE
363 if( iverbose >= 2 ) G4cout << " G4EP: AN2 " << AN2 << " gamma:"<<gamma<< " theta="<< thetaAver<<G4endl;
364#endif
365 G4double AU = 1./vpPreNorm.perp();
366 //t G4ThreeVector vU( vpPreNorm.cross( G4ThreeVector(0.,0.,1.) ) * AU );
367 G4ThreeVector vUPre( -AU*vpPreNorm.y(),
368 AU*vpPreNorm.x(),
369 0. );
370 G4ThreeVector vVPre( -vpPreNorm.z()*vUPre.y(),
371 vpPreNorm.z()*vUPre.x(),
372 vpPreNorm.x()*vUPre.y() - vpPreNorm.y()*vUPre.x() );
373
374 //
375 AU = 1./vpPostNorm.perp();
376 //t G4ThreeVector vU( vpPostNorm.cross( G4ThreeVector(0.,0.,1.) ) * AU );
377 G4ThreeVector vUPost( -AU*vpPostNorm.y(),
378 AU*vpPostNorm.x(),
379 0. );
380 G4ThreeVector vVPost( -vpPostNorm.z()*vUPost.y(),
381 vpPostNorm.z()*vUPost.x(),
382 vpPostNorm.x()*vUPost.y() - vpPostNorm.y()*vUPost.x() );
383#ifdef G4EVERBOSE
384 G4cout << " vpPostNorm " << vpPostNorm << G4endl;
385 if( iverbose >= 2 ) G4cout << " G4EP: AU " << AU << " vUPre " << vUPre << " vVPre " << vVPre << " vUPost " << vUPost << " vVPost " << vVPost << G4endl;
386#endif
387 G4Point3D deltaPos( vposPre - vposPost );
388
389 // * *** COMPLETE TRANSFORMATION MATRIX BETWEEN ERRORS AT X1 AND X2
390 // * *** FIELD GRADIENT PERPENDICULAR TO TRACK IS PRESENTLY NOT
391 // * *** TAKEN INTO ACCOUNT
392
393 G4double QP = QAver * pAver; // = -HAver
394#ifdef G4EVERBOSE
395 if( iverbose >= 2) G4cout << " G4EP: QP " << QP << " QAver " << QAver << " pAver " << pAver << G4endl;
396#endif
397 G4double ANV = -( vHAverNorm.x()*vUPost.x() + vHAverNorm.y()*vUPost.y() );
398 G4double ANU = ( vHAverNorm.x()*vVPost.x() + vHAverNorm.y()*vVPost.y() + vHAverNorm.z()*vVPost.z() );
399 G4double OMcosThetaAver = 1. - cosThetaAver;
400#ifdef G4EVERBOSE
401 if( iverbose >= 2) G4cout << "G4EP: OMcosThetaAver " << OMcosThetaAver << " cosThetaAver " << cosThetaAver << " thetaAver " << thetaAver << " QAver " << QAver << " stepLengthCm " << stepLengthCm << G4endl;
402#endif
403 G4double TMSINT = thetaAver - sinThetaAver;
404#ifdef G4EVERBOSE
405 if( iverbose >= 2 ) G4cout << " G4EP: ANV " << ANV << " ANU " << ANU << G4endl;
406#endif
407
408 G4ThreeVector vHUPre( -vHAverNorm.z() * vUPre.y(),
409 vHAverNorm.z() * vUPre.x(),
410 vHAverNorm.x() * vUPre.y() - vHAverNorm.y() * vUPre.x() );
411#ifdef G4EVERBOSE
412 // if( iverbose >= 2 ) G4cout << "G4EP: HUPre(1) " << vHUPre.x() << " " << vHAverNorm.z() << " " << vUPre.y() << G4endl;
413#endif
414 G4ThreeVector vHVPre( vHAverNorm.y() * vVPre.z() - vHAverNorm.z() * vVPre.y(),
415 vHAverNorm.z() * vVPre.x() - vHAverNorm.x() * vVPre.z(),
416 vHAverNorm.x() * vVPre.y() - vHAverNorm.y() * vVPre.x() );
417#ifdef G4EVERBOSE
418 if( iverbose >= 2 ) G4cout << " G4EP: HUPre " << vHUPre << " HVPre " << vHVPre << G4endl;
419#endif
420
421 //------------------- COMPUTE MATRIX
422 //---------- 1/P
423
424 transf[0][0] = 1.-deltaPInv*pAver*(1.+(vpPostNorm.x()*deltaPos.x()+vpPostNorm.y()*deltaPos.y()+vpPostNorm.z()*deltaPos.z())/stepLengthCm)
425 +2.*deltaPInv*pAver;
426
427 transf[0][1] = -deltaPInv/thetaAver*
428 ( TMSINT*gamma*(vHAverNorm.x()*vVPre.x()+vHAverNorm.y()*vVPre.y()+vHAverNorm.z()*vVPre.z()) +
429 sinThetaAver*(vVPre.x()*vpPostNorm.x()+vVPre.y()*vpPostNorm.y()+vVPre.z()*vpPostNorm.z()) +
430 OMcosThetaAver*(vHVPre.x()*vpPostNorm.x()+vHVPre.y()*vpPostNorm.y()+vHVPre.z()*vpPostNorm.z()) );
431
432 transf[0][2] = -sinpPre*deltaPInv/thetaAver*
433 ( TMSINT*gamma*(vHAverNorm.x()*vUPre.x()+vHAverNorm.y()*vUPre.y() ) +
434 sinThetaAver*(vUPre.x()*vpPostNorm.x()+vUPre.y()*vpPostNorm.y() ) +
435 OMcosThetaAver*(vHUPre.x()*vpPostNorm.x()+vHUPre.y()*vpPostNorm.y()+vHUPre.z()*vpPostNorm.z()) );
436
437 transf[0][3] = -deltaPInv/stepLengthCm*(vUPre.x()*vpPostNorm.x()+vUPre.y()*vpPostNorm.y() );
438
439 transf[0][4] = -deltaPInv/stepLengthCm*(vVPre.x()*vpPostNorm.x()+vVPre.y()*vpPostNorm.y()+vVPre.z()*vpPostNorm.z());
440
441 // *** Lambda
442 transf[1][0] = -QP*ANV*(vpPostNorm.x()*deltaPos.x()+vpPostNorm.y()*deltaPos.y()+vpPostNorm.z()*deltaPos.z())
443 *(1.+deltaPInv*pAver);
444#ifdef G4EVERBOSE
445 if(iverbose >= 3) G4cout << "ctransf10= " << transf[1][0] << " " << -QP<< " " << ANV<< " " << vpPostNorm.x()<< " " << deltaPos.x()<< " " << vpPostNorm.y()<< " " << deltaPos.y()<< " " << vpPostNorm.z()<< " " << deltaPos.z()
446 << " " << deltaPInv<< " " << pAver << G4endl;
447#endif
448
449 transf[1][1] = cosThetaAver*(vVPre.x()*vVPost.x()+vVPre.y()*vVPost.y()+vVPre.z()*vVPost.z()) +
450 sinThetaAver*(vHVPre.x()*vVPost.x()+vHVPre.y()*vVPost.y()+vHVPre.z()*vVPost.z()) +
451 OMcosThetaAver*(vHAverNorm.x()*vVPre.x()+vHAverNorm.y()*vVPre.y()+vHAverNorm.z()*vVPre.z())*
452 (vHAverNorm.x()*vVPost.x()+vHAverNorm.y()*vVPost.y()+vHAverNorm.z()*vVPost.z()) +
453 ANV*( -sinThetaAver*(vVPre.x()*vpPostNorm.x()+vVPre.y()*vpPostNorm.y()+vVPre.z()*vpPostNorm.z()) +
454 OMcosThetaAver*(vVPre.x()*AN2.x()+vVPre.y()*AN2.y()+vVPre.z()*AN2.z()) -
455 TMSINT*gamma*(vHAverNorm.x()*vVPre.x()+vHAverNorm.y()*vVPre.y()+vHAverNorm.z()*vVPre.z()) );
456
457 transf[1][2] = cosThetaAver*(vUPre.x()*vVPost.x()+vUPre.y()*vVPost.y() ) +
458 sinThetaAver*(vHUPre.x()*vVPost.x()+vHUPre.y()*vVPost.y()+vHUPre.z()*vVPost.z()) +
459 OMcosThetaAver*(vHAverNorm.x()*vUPre.x()+vHAverNorm.y()*vUPre.y() )*
460 (vHAverNorm.x()*vVPost.x()+vHAverNorm.y()*vVPost.y()+vHAverNorm.z()*vVPost.z()) +
461 ANV*( -sinThetaAver*(vUPre.x()*vpPostNorm.x()+vUPre.y()*vpPostNorm.y() ) +
462 OMcosThetaAver*(vUPre.x()*AN2.x()+vUPre.y()*AN2.y() ) -
463 TMSINT*gamma*(vHAverNorm.x()*vUPre.x()+vHAverNorm.y()*vUPre.y() ) );
464 transf[1][2] = sinpPre*transf[1][2];
465
466 transf[1][3] = -QAver*ANV*(vUPre.x()*vpPostNorm.x()+vUPre.y()*vpPostNorm.y() );
467
468 transf[1][4] = -QAver*ANV*(vVPre.x()*vpPostNorm.x()+vVPre.y()*vpPostNorm.y()+vVPre.z()*vpPostNorm.z());
469
470 // *** Phi
471
472 transf[2][0] = -QP*ANU*(vpPostNorm.x()*deltaPos.x()+vpPostNorm.y()*deltaPos.y()+vpPostNorm.z()*deltaPos.z())*sinpPostInv
473 *(1.+deltaPInv*pAver);
474#ifdef G4EVERBOSE
475 if(iverbose >= 3)G4cout <<"ctransf20= " << transf[2][0] <<" "<< -QP<<" "<<ANU<<" "<<vpPostNorm.x()<<" "<<deltaPos.x()<<" "<<vpPostNorm.y()<<" "<<deltaPos.y()<<" "<<vpPostNorm.z()<<" "<<deltaPos.z()<<" "<<sinpPostInv
476 <<" "<<deltaPInv<<" "<<pAver<< G4endl;
477#endif
478 transf[2][1] = cosThetaAver*(vVPre.x()*vUPost.x()+vVPre.y()*vUPost.y() ) +
479 sinThetaAver*(vHVPre.x()*vUPost.x()+vHVPre.y()*vUPost.y() ) +
480 OMcosThetaAver*(vHAverNorm.x()*vVPre.x()+vHAverNorm.y()*vVPre.y()+vHAverNorm.z()*vVPre.z())*
481 (vHAverNorm.x()*vUPost.x()+vHAverNorm.y()*vUPost.y() ) +
482 ANU*( -sinThetaAver*(vVPre.x()*vpPostNorm.x()+vVPre.y()*vpPostNorm.y()+vVPre.z()*vpPostNorm.z()) +
483 OMcosThetaAver*(vVPre.x()*AN2.x()+vVPre.y()*AN2.y()+vVPre.z()*AN2.z()) -
484 TMSINT*gamma*(vHAverNorm.x()*vVPre.x()+vHAverNorm.y()*vVPre.y()+vHAverNorm.z()*vVPre.z()) );
485 transf[2][1] = sinpPostInv*transf[2][1];
486
487 transf[2][2] = cosThetaAver*(vUPre.x()*vUPost.x()+vUPre.y()*vUPost.y() ) +
488 sinThetaAver*(vHUPre.x()*vUPost.x()+vHUPre.y()*vUPost.y() ) +
489 OMcosThetaAver*(vHAverNorm.x()*vUPre.x()+vHAverNorm.y()*vUPre.y() )*
490 (vHAverNorm.x()*vUPost.x()+vHAverNorm.y()*vUPost.y() ) +
491 ANU*( -sinThetaAver*(vUPre.x()*vpPostNorm.x()+vUPre.y()*vpPostNorm.y() ) +
492 OMcosThetaAver*(vUPre.x()*AN2.x()+vUPre.y()*AN2.y() ) -
493 TMSINT*gamma*(vHAverNorm.x()*vUPre.x()+vHAverNorm.y()*vUPre.y() ) );
494 transf[2][2] = sinpPostInv*sinpPre*transf[2][2];
495
496 transf[2][3] = -QAver*ANU*(vUPre.x()*vpPostNorm.x()+vUPre.y()*vpPostNorm.y() )*sinpPostInv;
497#ifdef G4EVERBOSE
498 if(iverbose >= 3)G4cout <<"ctransf23= " << transf[2][3] <<" "<< -QAver<<" "<<ANU<<" "<<vUPre.x()<<" "<<vpPostNorm.x()<<" "<< vUPre.y()<<" "<<vpPostNorm.y()<<" "<<sinpPostInv<<G4endl;
499#endif
500
501 transf[2][4] = -QAver*ANU*(vVPre.x()*vpPostNorm.x()+vVPre.y()*vpPostNorm.y()+vVPre.z()*vpPostNorm.z())*sinpPostInv;
502
503 // *** Yt
504
505 transf[3][0] = pAver*(vUPost.x()*deltaPos.x()+vUPost.y()*deltaPos.y() )
506 *(1.+deltaPInv*pAver);
507#ifdef G4EVERBOSE
508 if(iverbose >= 3) G4cout <<"ctransf30= " << transf[3][0] <<" "<< pAver<<" "<<vUPost.x()<<" "<<deltaPos.x()<<" "<<vUPost.y()<<" "<<deltaPos.y()
509 <<" "<<deltaPInv<<" "<<pAver<<G4endl;
510#endif
511
512 transf[3][1] = ( sinThetaAver*(vVPre.x()*vUPost.x()+vVPre.y()*vUPost.y() ) +
513 OMcosThetaAver*(vHVPre.x()*vUPost.x()+vHVPre.y()*vUPost.y() ) +
514 TMSINT*(vHAverNorm.x()*vUPost.x()+vHAverNorm.y()*vUPost.y() )*
515 (vHAverNorm.x()*vVPre.x()+vHAverNorm.y()*vVPre.y()+vHAverNorm.z()*vVPre.z()) )/QAver;
516
517 transf[3][2] = ( sinThetaAver*(vUPre.x()*vUPost.x()+vUPre.y()*vUPost.y() ) +
518 OMcosThetaAver*(vHUPre.x()*vUPost.x()+vHUPre.y()*vUPost.y() ) +
519 TMSINT*(vHAverNorm.x()*vUPost.x()+vHAverNorm.y()*vUPost.y() )*
520 (vHAverNorm.x()*vUPre.x()+vHAverNorm.y()*vUPre.y() ) )*sinpPre/QAver;
521#ifdef G4EVERBOSE
522 if(iverbose >= 3) G4cout <<"ctransf32= " << transf[3][2] <<" "<< sinThetaAver<<" "<<vUPre.x()<<" "<<vUPost.x()<<" "<<vUPre.y()<<" "<<vUPost.y() <<" "<<
523 OMcosThetaAver<<" "<<vHUPre.x()<<" "<<vUPost.x()<<" "<<vHUPre.y()<<" "<<vUPost.y() <<" "<<
524 TMSINT<<" "<<vHAverNorm.x()<<" "<<vUPost.x()<<" "<<vHAverNorm.y()<<" "<<vUPost.y() <<" "<<
525 vHAverNorm.x()<<" "<<vUPre.x()<<" "<<vHAverNorm.y()<<" "<<vUPre.y() <<" "<<sinpPre<<" "<<QAver<<G4endl;
526#endif
527
528 transf[3][3] = (vUPre.x()*vUPost.x()+vUPre.y()*vUPost.y() );
529
530 transf[3][4] = (vVPre.x()*vUPost.x()+vVPre.y()*vUPost.y() );
531
532 // *** Zt
533 transf[4][0] = pAver*(vVPost.x()*deltaPos.x()+vVPost.y()*deltaPos.y()+vVPost.z()*deltaPos.z())
534 *(1.+deltaPInv*pAver);
535
536 transf[4][1] = ( sinThetaAver*(vVPre.x()*vVPost.x()+vVPre.y()*vVPost.y()+vVPre.z()*vVPost.z()) +
537 OMcosThetaAver*(vHVPre.x()*vVPost.x()+vHVPre.y()*vVPost.y()+vHVPre.z()*vVPost.z()) +
538 TMSINT*(vHAverNorm.x()*vVPost.x()+vHAverNorm.y()*vVPost.y()+vHAverNorm.z()*vVPost.z())*
539 (vHAverNorm.x()*vVPre.x()+vHAverNorm.y()*vVPre.y()+vHAverNorm.z()*vVPre.z()) )/QAver;
540#ifdef G4EVERBOSE
541 if(iverbose >= 3)G4cout <<"ctransf41= " << transf[4][1] <<" "<< sinThetaAver<<" "<< OMcosThetaAver <<" "<<TMSINT<<" "<< vVPre <<" "<<vVPost <<" "<<vHVPre<<" "<<vHAverNorm <<" "<< QAver<<G4endl;
542#endif
543
544 transf[4][2] = ( sinThetaAver*(vUPre.x()*vVPost.x()+vUPre.y()*vVPost.y() ) +
545 OMcosThetaAver*(vHUPre.x()*vVPost.x()+vHUPre.y()*vVPost.y()+vHUPre.z()*vVPost.z()) +
546 TMSINT*(vHAverNorm.x()*vVPost.x()+vHAverNorm.y()*vVPost.y()+vHAverNorm.z()*vVPost.z())*
547 (vHAverNorm.x()*vUPre.x()+vHAverNorm.y()*vUPre.y() ) )*sinpPre/QAver;
548
549 transf[4][3] = (vUPre.x()*vVPost.x()+vUPre.y()*vVPost.y() );
550
551 transf[4][4] = (vVPre.x()*vVPost.x()+vVPre.y()*vVPost.y()+vVPre.z()*vVPost.z());
552 // if(iverbose >= 3) G4cout <<"ctransf44= " << transf[4][4] <<" "<< vVPre.x() <<" "<<vVPost.x() <<" "<< vVPre.y() <<" "<< vVPost.y() <<" "<< vVPre.z() <<" "<< vVPost.z() << G4endl;
553
554
555#ifdef G4EVERBOSE
556 if( iverbose >= 1 ) G4cout << "G4EP: transf matrix computed " << transf << G4endl;
557#endif
558 /* for( G4int ii=0;ii<5;ii++){
559 for( G4int jj=0;jj<5;jj++){
560 G4cout << transf[ii][jj] << " ";
561 }
562 G4cout << G4endl;
563 } */
564 }
565 }
566 // end of calculate transformation except it NEUTRAL PARTICLE OR FIELDFREE REGION
567 /* if( iverbose >= 1 ) G4cout << "G4EP: transf not updated but initialized " << theFirstStep << G4endl;
568 if( theFirstStep ) {
569 theTransfMat = transf;
570 theFirstStep = false;
571 }else{
572 theTransfMat = theTransfMat * transf;
573 if( iverbose >= 1 ) G4cout << "G4EP: transf matrix accumulated" << theTransfMat << G4endl;
574 }
575 */
576 theTransfMat = transf;
577#ifdef G4EVERBOSE
578 if( iverbose >= 1 ) G4cout << "G4EP: error matrix before transformation " << fError << G4endl;
579 if( iverbose >= 2 ) G4cout << " tf * err " << theTransfMat * fError << G4endl
580 << " transf matrix " << theTransfMat.T() << G4endl;
581#endif
582
583 fError = fError.similarity(theTransfMat).T();
584 //- fError = transf * fError * transf.T();
585#ifdef G4EVERBOSE
586 if( iverbose >= 1 ) G4cout << "G4EP: error matrix propagated " << fError << G4endl;
587#endif
588
589 //? S = B*S*BT S.similarity(B)
590 //? R = S
591 // not needed * *** TRANSFORM ERROR MATRIX FROM INTERNAL TO EXTERNAL VARIABLES;
592
593 PropagateErrorMSC( aTrack );
594
595 PropagateErrorIoni( aTrack );
596
597 return 0;
598}
599
600
601//------------------------------------------------------------------------
602G4int G4ErrorFreeTrajState::PropagateErrorMSC( const G4Track* aTrack )
603{
604 G4ThreeVector vpPre = aTrack->GetMomentum()/GeV;
605 G4double pPre = vpPre.mag();
606 G4double pBeta = pPre*pPre / (aTrack->GetTotalEnergy()/GeV);
607 G4double stepLengthCm = aTrack->GetStep()->GetStepLength()/cm;
608
609 G4Material* mate = aTrack->GetVolume()->GetLogicalVolume()->GetMaterial();
610 G4double effZ, effA;
611 CalculateEffectiveZandA( mate, effZ, effA );
612
613#ifdef G4EVERBOSE
614 if( iverbose >= 4 ) G4cout << "material " << mate->GetName()
615 //<< " " << mate->GetZ() << " " << mate->GetA()
616 << " effZ:" << effZ << " effA:" << effA
617 << " dens(g/mole):" << mate->GetDensity()/g*mole << " Radlen/cm:" << mate->GetRadlen()/cm << " nuclLen/cm" << mate->GetNuclearInterLength()/cm << G4endl;
618#endif
619
620 G4double RI = stepLengthCm / (mate->GetRadlen()/cm);
621#ifdef G4EVERBOSE
622 if( iverbose >= 4 ) G4cout << std::setprecision(6) << std::setw(6) << "G4EP:MSC: RI=X/X0 " << RI << " stepLengthCm " << stepLengthCm << " radlen/cm " << (mate->GetRadlen()/cm) << " RI*1.e10:" << RI*1.e10 << G4endl;
623#endif
624 G4double charge = aTrack->GetDynamicParticle()->GetCharge();
625 G4double DD = 1.8496E-4*RI*(charge/pBeta * charge/pBeta );
626#ifdef G4EVERBOSE
627 if( iverbose >= 3 ) G4cout << "G4EP:MSC: D*1E6= " << DD*1.E6 <<" pBeta " << pBeta << G4endl;
628#endif
629 G4double S1 = DD*stepLengthCm*stepLengthCm/3.;
630 G4double S2 = DD;
631 G4double S3 = DD*stepLengthCm/2.;
632
633 G4double CLA = std::sqrt( vpPre.x() * vpPre.x() + vpPre.y() * vpPre.y() )/pPre;
634#ifdef G4EVERBOSE
635 if( iverbose >= 2 ) G4cout << std::setw(6) << "G4EP:MSC: RI " << RI << " S1 " << S1 << " S2 " << S2 << " S3 " << S3 << " CLA " << CLA << G4endl;
636#endif
637 fError[1][1] += S2;
638 fError[1][4] -= S3;
639 fError[2][2] += S2/CLA/CLA;
640 fError[2][3] += S3/CLA;
641 fError[3][3] += S1;
642 fError[4][4] += S1;
643
644#ifdef G4EVERBOSE
645 if( iverbose >= 2 ) G4cout << "G4EP:MSC: error matrix propagated msc " << fError << G4endl;
646#endif
647
648 return 0;
649}
650
651
652//------------------------------------------------------------------------
653void G4ErrorFreeTrajState::CalculateEffectiveZandA( const G4Material* mate, G4double& effZ, G4double& effA )
654{
655 effZ = 0.;
656 effA = 0.;
657 G4int ii, nelem = mate->GetNumberOfElements();
658 const G4double* fracVec = mate->GetFractionVector();
659 for(ii=0; ii < nelem; ii++ ) {
660 effZ += mate->GetElement( ii )->GetZ() * fracVec[ii];
661 effA += mate->GetElement( ii )->GetA() * fracVec[ii] /g*mole;
662 }
663
664}
665
666
667//------------------------------------------------------------------------
668G4int G4ErrorFreeTrajState::PropagateErrorIoni( const G4Track* aTrack )
669{
670 G4double stepLengthCm = aTrack->GetStep()->GetStepLength()/cm;
671#ifdef G4EVERBOSE
672 G4double DEDX2;
673 if( stepLengthCm < 1.E-7 ) {
674 DEDX2=0.;
675 }
676#endif
677 // * Calculate xi factor (KeV).
678 G4Material* mate = aTrack->GetVolume()->GetLogicalVolume()->GetMaterial();
679 G4double effZ, effA;
680 CalculateEffectiveZandA( mate, effZ, effA );
681
682 G4double Etot = aTrack->GetTotalEnergy()/GeV;
683 G4double beta = aTrack->GetMomentum().mag()/GeV / Etot;
684 G4double mass = aTrack->GetDynamicParticle()->GetMass() / GeV;
685 G4double gamma = Etot / mass;
686
687 // * Calculate xi factor (keV).
688 G4double XI = 153.5*effZ*stepLengthCm*(mate->GetDensity()/mg*mole) /
689 (effA*beta*beta);
690
691#ifdef G4EVERBOSE
692 if( iverbose >= 2 ){
693 G4cout << "G4EP:IONI: XI/keV " << XI << " beta " << beta << " gamma " << gamma << G4endl;
694 G4cout << " density " << (mate->GetDensity()/mg*mole) << " effA " << effA << " step " << stepLengthCm << G4endl;
695 }
696#endif
697 // * Maximum energy transfer to atomic electron (KeV).
698 G4double eta = beta*gamma;
699 G4double etasq = eta*eta;
700 G4double eMass = 0.51099906/GeV;
701 G4double massRatio = eMass / mass;
702 G4double F1 = 2*eMass*etasq;
703 G4double F2 = 1. + 2. * massRatio * gamma + massRatio * massRatio;
704 G4double Emax = 1.E+6*F1/F2; // now in keV
705
706 // * *** and now sigma**2 in GeV
707 G4double dedxSq = XI*Emax*(1.-(beta*beta/2.))*1.E-12; // now in GeV^2
708 /*The above formula for var(1/p) good for dens scatterers. However, for MIPS passing
709 through a gas it leads to overestimation. Further more for incident electrons
710 the Emax is almost equal to incident energy.
711 This leads to k=Xi/Emax as small as e-6 and gradually the cov matrix explodes.
712
713 http://www2.pv.infn.it/~rotondi/kalman_1.pdf
714
715 Since I do not have enough info at the moment to implement Landau & sub-Landau models for k=Xi/Emax <0.01 I'll saturate k at this value for now
716 */
717
718 if (XI/Emax<0.01) dedxSq *=XI/Emax*100 ;// Quench for low Elos, see above: newVar=odVar *k/0.01
719
720#ifdef G4EVERBOSE
721 if( iverbose >= 2 ) G4cout << "G4EP:IONI: DEDX^2(GeV^2) " << dedxSq << " emass/GeV: " << eMass << " Emax/keV: " << Emax
722 <<" k=Xi/Emax="<< XI/Emax<< G4endl;
723
724#endif
725
726 G4double pPre6 = (aTrack->GetStep()->GetPreStepPoint()->GetMomentum()/GeV).mag();
727 pPre6 = std::pow(pPre6, 6 );
728 //Apply it to error
729 fError[0][0] += Etot*Etot*dedxSq / pPre6;
730#ifdef G4EVERBOSE
731 if( iverbose >= 2 ) G4cout << "G4:IONI Etot/GeV: " << Etot << " err_dedx^2/GeV^2: " << dedxSq << " p^6: " << pPre6 << G4endl;
732 if( iverbose >= 2 ) G4cout << "G4EP:IONI: error2_from_ionisation " << (Etot*Etot*dedxSq) / pPre6 << G4endl;
733#endif
734
735 return 0;
736}
kCarTolerance
const G4double kCarTolerance
Definition: G4BoundingEnvelope.cc:37
operator<<
std::ostream & operator<<(std::ostream &out, const G4ErrorFreeTrajState &ts)
Definition: G4ErrorFreeTrajState.cc:185
G4ErrorMode_PropBackwards
@ G4ErrorMode_PropBackwards
Definition: G4ErrorPropagatorData.hh:44
G4ErrorStage_Deflation
@ G4ErrorStage_Deflation
Definition: G4ErrorPropagatorData.hh:60
G4ErrorTrajErr
G4ErrorSymMatrix G4ErrorTrajErr
Definition: G4ErrorTrajErr.hh:43
G4eTS_FREE
@ G4eTS_FREE
Definition: G4ErrorTrajState.hh:46
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:36
G4double
double G4double
Definition: G4Types.hh:83
G4int
int G4int
Definition: G4Types.hh:85
G4Vector3D
HepGeom::Vector3D< G4double > G4Vector3D
Definition: G4Vector3D.hh:34
G4endl
#define G4endl
Definition: G4ios.hh:57
G4cout
G4GLOB_DLL std::ostream G4cout
CLHEP::Hep3Vector::z
double z() const
CLHEP::Hep3Vector::theta
double theta() const
CLHEP::Hep3Vector::x
double x() const
CLHEP::Hep3Vector::y
double y() const
CLHEP::Hep3Vector::cross
Hep3Vector cross(const Hep3Vector &) const
CLHEP::Hep3Vector::mag
double mag() const
G4DynamicParticle::GetMass
G4double GetMass() const
G4DynamicParticle::GetCharge
G4double GetCharge() const
G4Element::GetZ
G4double GetZ() const
Definition: G4Element.hh:130
G4Element::GetA
G4double GetA() const
Definition: G4Element.hh:138
G4ErrorFreeTrajParam::Update
void Update(const G4Track *aTrack)
Definition: G4ErrorFreeTrajParam.cc:78
G4ErrorFreeTrajParam::GetLambda
G4double GetLambda() const
Definition: G4ErrorFreeTrajParam.hh:80
G4ErrorFreeTrajParam::GetDirection
G4Vector3D GetDirection() const
Definition: G4ErrorFreeTrajParam.hh:77
G4ErrorFreeTrajState::Update
virtual G4int Update(const G4Track *aTrack)
Definition: G4ErrorFreeTrajState.cc:175
G4ErrorFreeTrajState::PropagateError
virtual G4int PropagateError(const G4Track *aTrack)
Definition: G4ErrorFreeTrajState.cc:202
G4ErrorFreeTrajState::GetParameters
G4ErrorFreeTrajParam GetParameters() const
Definition: G4ErrorFreeTrajState.hh:107
G4ErrorFreeTrajState::Dump
virtual void Dump(std::ostream &out=G4cout) const
Definition: G4ErrorFreeTrajState.cc:169
G4ErrorMatrix::T
G4ErrorMatrix T() const
Definition: G4ErrorMatrix.cc:373
G4ErrorPropagatorData::GetErrorPropagatorData
static G4ErrorPropagatorData * GetErrorPropagatorData()
Definition: G4ErrorPropagatorData.cc:54
G4ErrorSurfaceTrajState::GetParameters
G4ErrorSurfaceTrajParam GetParameters() const
Definition: G4ErrorSurfaceTrajState.hh:92
G4ErrorSymMatrix::similarity
G4ErrorSymMatrix similarity(const G4ErrorMatrix &m1) const
Definition: G4ErrorSymMatrix.cc:588
G4ErrorSymMatrix::T
G4ErrorSymMatrix T() const
G4ErrorTrajState::DumpPosMomError
void DumpPosMomError(std::ostream &out=G4cout) const
Definition: G4ErrorTrajState.cc:163
G4ErrorTrajState::fError
G4ErrorTrajErr fError
Definition: G4ErrorTrajState.hh:135
G4ErrorTrajState::UpdatePosMom
void UpdatePosMom(const G4Point3D &pos, const G4Vector3D &mom)
Definition: G4ErrorTrajState.cc:124
G4ErrorTrajState::GetError
G4ErrorTrajErr GetError() const
Definition: G4ErrorTrajState.hh:111
G4FieldManager::GetDetectorField
const G4Field * GetDetectorField() const
G4Field::GetFieldValue
virtual void GetFieldValue(const G4double Point[4], G4double *fieldArr) const =0
G4GeometryTolerance::GetSurfaceTolerance
G4double GetSurfaceTolerance() const
Definition: G4GeometryTolerance.cc:78
G4GeometryTolerance::GetInstance
static G4GeometryTolerance * GetInstance()
Definition: G4GeometryTolerance.cc:64
G4LogicalVolume::GetMaterial
G4Material * GetMaterial() const
Definition: G4LogicalVolume.cc:429
G4Material::GetDensity
G4double GetDensity() const
Definition: G4Material.hh:178
G4Material::GetElement
const G4Element * GetElement(G4int iel) const
Definition: G4Material.hh:200
G4Material::GetFractionVector
const G4double * GetFractionVector() const
Definition: G4Material.hh:192
G4Material::GetNumberOfElements
size_t GetNumberOfElements() const
Definition: G4Material.hh:184
G4Material::GetRadlen
G4double GetRadlen() const
Definition: G4Material.hh:218
G4Material::GetName
const G4String & GetName() const
Definition: G4Material.hh:175
G4Material::GetNuclearInterLength
G4double GetNuclearInterLength() const
Definition: G4Material.hh:221
G4StepPoint::GetMomentum
G4ThreeVector GetMomentum() const
G4StepPoint::GetPosition
const G4ThreeVector & GetPosition() const
G4Step::GetPreStepPoint
G4StepPoint * GetPreStepPoint() const
G4Step::GetStepLength
G4double GetStepLength() const
G4Track::GetVolume
G4VPhysicalVolume * GetVolume() const
G4Track::GetPosition
const G4ThreeVector & GetPosition() const
G4Track::GetMomentum
G4ThreeVector GetMomentum() const
G4Track::GetDynamicParticle
const G4DynamicParticle * GetDynamicParticle() const
G4Track::GetTotalEnergy
G4double GetTotalEnergy() const
G4Track::GetStep
const G4Step * GetStep() const
G4TransportationManager::GetTransportationManager
static G4TransportationManager * GetTransportationManager()
Definition: G4TransportationManager.cc:104
G4TransportationManager::GetFieldManager
G4FieldManager * GetFieldManager() const
G4VPhysicalVolume::GetLogicalVolume
G4LogicalVolume * GetLogicalVolume() const
HepGeom::BasicVector3D::cross
BasicVector3D< T > cross(const BasicVector3D< T > &v) const
Definition: BasicVector3D.h:281