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