| Summary: | computed error of |P| too large | ||
|---|---|---|---|
| Product: | Examples/Extended | Reporter: | janstar1122 |
| Component: | errorpropagation | Assignee: | Pedro.Arce |
| Status: | RESOLVED FIXED | ||
| Severity: | trivial | ||
| Priority: | P5 | ||
| Version: | 9.6 | ||
| Hardware: | All | ||
| OS: | All | ||
Follow up....
G4int G4ErrorFreeTrajState::PropagateErrorIoni( const G4Track* aTrack ) {
uses formula for var(1/p) good for dens scatterers. However, for MIPS passing through a gas it leads to overestimation. Further more for incident electrons the Emax is almost equal to incident energy.
This leads to k=Xi/Emax as small as e-6 and gradually the cov matrix explodes
http://www2.pv.infn.it/~rotondi/kalman_1.pdf
I do not have enough info at the moment to implement Landau & sub-Landau models for k=Xi/Emax <0.01
|
Hi, I'm questioning how error of 1/P component of the covariance matrix is computed. My starting point is your example use case: geant4.9.6.p01/examples/extended/errorpropagation I copied it, and make only one change: moved the target plane closer to the start, to 10cm and later to 20 cm. Set G4E to propagate forward: export G4ERROR_TARGET=PLANE_SURFACE export G4ERROR_MODE=FORWARDS export G4ERROR_PROP=UNTIL_TARGET and printed position, momentum and momentum error (for the free state) at start and at this target location. This is what I see (I run my code twice for 10 & 20 cm X-target location), still Air if the only material traversed. ---- At START ------ Position (mm): (0,0,0) phi/deg=0 Rxy/cm=0 Momentum (MeV): (20000,0,0) phi/deg=0 PT/MeV=20000 freeErr: |P| (MeV) 20000.00 +/- 0.00 Looks good --- at Target X=10 cm------ Position (mm): (100,-0.07494821925585556,0) phi/deg=-0.04294215841327705 Rxy/cm=10.00000280861739 Momentum (MeV): (19999.94417956696,-29.97923330181213,0) phi/deg=-0.08588435244242193 PT/MeV=19999.96664847778 freeErr: |P| (MeV) 19999.97 +/- 2440.97 --- at Target X=20 cm------ Position (mm): (200,-0.299793548487774,0) phi/deg=-0.08588446094331748 Rxy/cm=20.00002246903031 Momentum (MeV): (19999.84342105256,-59.95842910881996,0) phi/deg=-0.1717690768833353 PT/MeV=19999.93329688478 freeErr: |P| (MeV) 19999.93 +/- 3452.04 My conclusions: the mean energy=momentum for initial 20 GeV muon are reduced by: - 0.06 MeV after 10 cm in air. - 0.16 MeV after 20 cm in air. what makes sense. But the error of momentum , computed as shown below, is humongous and irrational: err(P) =2.4 GeV after 10 cm in air. err(P)=3.5 GeV after 20 cm in air. It can't be the Landau energy loss width is of 3.5 GeV if the mean of it is just 0.16 MeV. Sth is realy wrong . Could you please have a look at it? Thanks Jan --- code used ------ void printFreeTraj(G4ErrorFreeTrajState *g4ErrorTrajState, G4String txt){ G4cout << "\n $$$---------------- FreeTraj " << txt<<G4endl; // extract current state info G4Point3D posEnd = g4ErrorTrajState->GetPosition(); G4Normal3D momEnd = g4ErrorTrajState->GetMomentum(); G4ErrorTrajErr errorEnd = g4ErrorTrajState->GetError(); G4cout << " Position (mm): " << posEnd <<" phi/deg="<<posEnd.phi()/deg<<" Rxy/cm="<<posEnd.perp()/cm<< G4endl << " Momentum (MeV): " << momEnd <<" phi/deg="<<momEnd.phi()/deg<<" PT/MeV="<<momEnd.perp()/MeV<< G4endl; double rsigP= sqrt(errorEnd(1,1))* (momEnd.mag()); printf("freeErr: |P| (MeV) %.2f +/- %.2f rsigP=%.3f\n",momEnd.mag(),rsigP*momEnd.mag() ,rsigP); }