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