source: trunk/Documentation/geant4/UserDocumentation/UsersGuides/ForApplicationDeveloper/html/ch05s08.html @ 901

Last change on this file since 901 was 901, checked in by garnier, 16 years ago

Add Geant4 Documentation at 8.12.2008
File size: 21.8 KB
Line 
1<html><head><meta http-equiv="Content-Type" content="text/html; charset=ISO-8859-1"><title>5.8.  Track Error Propagation</title><link rel="stylesheet" href="../xml/XSLCustomizationLayer/G4HTMLStylesheet.css" type="text/css"><meta name="generator" content="DocBook XSL Stylesheets V1.71.1"><link rel="start" href="index.html" title="Geant4 User's Guide for Application Developers"><link rel="up" href="ch05.html" title="Chapter 5.  Tracking and Physics"><link rel="prev" href="ch05s07.html" title="5.7.  User Limits"><link rel="next" href="ch06.html" title="Chapter 6.  User Actions"><script language="JavaScript">
2function remote_win(fName)
3{
4   var url = "AllResources/Detector/geometry.src/" + fName;
5   RemoteWin=window.open(url,"","resizable=no,toolbar=0,location=0,directories=0,status=0,menubar=0,scrollbars=0,copyhistory=0,width=520,height=520")
6   RemoteWin.creator=self
7}
8</script></head><body bgcolor="white" text="black" link="#0000FF" vlink="#840084" alink="#0000FF"><div class="navheader"><table width="100%" summary="Navigation header"><tr><th colspan="3" align="center">5.8. 
9Track Error Propagation
10</th></tr><tr><td width="20%" align="left"><a accesskey="p" href="ch05s07.html"><img src="AllResources/IconsGIF/prev.gif" alt="Prev"></a> </td><th width="60%" align="center">Chapter 5. 
11Tracking and Physics
12</th><td width="20%" align="right"> <a accesskey="n" href="ch06.html"><img src="AllResources/IconsGIF/next.gif" alt="Next"></a></td></tr></table><hr></div><div class="sect1" lang="en"><div class="titlepage"><div><div><h2 class="title" style="clear: both"><a name="sect.G4Error"></a>5.8. 
13Track Error Propagation
14</h2></div></div></div><p>
15The error propagation package serves to propagate one particle together with its error from a given trajectory state until a user-defined target is reached (a surface, a volume, a given track length,...).
16</p><div class="sect2" lang="en"><div class="titlepage"><div><div><h3 class="title"><a name="sect.G4Error.Physics"></a>5.8.1. 
17Physics
18</h3></div></div></div><p>
19The error propagator package computes the average trajectory that a particle would follow. 
20This means that the physics list must have the following characteristics:
21</p><div class="itemizedlist"><ul type="disc"><li> No multiple scattering </li><li> No random fluctuations for energy loss </li><li> No creation of secondary tracks </li><li> No hadronic processes </li></ul></div><p>
22It has also to be taken into account that when the propagation is done backwards (in the direction opposed to the one the original track traveled) the energy loss has to be changed into an energy gain.
23</p><p>
24All this is done in the <code class="literal">G4ErrorPhysicsList</code> class, that is automatically set by <code class="literal">G4ErrorPropagatorManager</code> as the GEANT4 physics list.
25It sets <code class="literal">G4ErrorEnergyLoss</code> as unique electromagnetic process. This process uses the GEANT4 class <code class="literal">G4EnergyLossForExtrapolator</code> to compute the average energy loss for forwards or backwards propagation.
26To avoid getting too different energy loss calculation when the propagation is done forwards (when the energy at the beginning of the step is used) or backwards (when the energy at the end of the step is used, always smaller than at the beginning) <code class="literal">G4ErrorEnergyLoss</code> computes once the energy loss and then replaces the original energy loss by subtracting/adding half of this value (what is approximately the same as computing the energy loss with the energy at the middle of the step).
27In this way, a better calculation of the energy loss is obtained with a minimal impact on the total CPU time.
28</p><p>
29The user may use his/her own physics list instead of <code class="literal">G4ErrorPhysicsList</code>. As it is not needed to define a physics list when running this package, the user may have not realized that somewhere else in his/her application it has been defined; therefore a warning will be sent to advert the user that he is using a physics list different to <code class="literal">G4ErrorPhysicsList</code>. If a new physics list is used, it should also initialize the <code class="literal">G4ErrorMessenger</code> with the classes that serve to limit the step:
30
31</p><div class="informalexample"><pre class="programlisting">
32  G4ErrorEnergyLoss* eLossProcess = new G4ErrorEnergyLoss;
33  G4ErrorStepLengthLimitProcess* stepLengthLimitProcess = new G4ErrorStepLengthLimitProcess;
34  G4ErrorMagFieldLimitProcess* magFieldLimitProcess = new G4ErrorMagFieldLimitProcess;
35  new G4ErrorMessenger( stepLengthLimitProcess, magFieldLimitProcess, eLossProcess );
36</pre></div><p>
37
38</p><p>
39To ease the use of this package in the reconstruction code, the physics list, whether <code class="literal">G4ErrorPhysicsList</code> or the user's one, will be automatically initialized before starting the track propagation if it has not been done by the user.
40</p></div><div class="sect2" lang="en"><div class="titlepage"><div><div><h3 class="title"><a name="sect.G4Error.TrajState"></a>5.8.2. 
41Trajectory state
42</h3></div></div></div><p>
43The user has to provide the particle trajectory state at the initial point.
44To do this it has to create an object of one of the children classes of <code class="literal">G4ErrorTrajState</code>, providing:
45</p><div class="itemizedlist"><ul type="disc"><li> Particle type </li><li> Position </li><li> Momentum </li><li> Trajectory error matrix </li></ul></div><div class="informalexample"><pre class="programlisting">
46 G4ErrorTrajState( const G4String&amp; partType,
47                   const G4Point3D&amp; pos,
48                   const G4Vector3D&amp; mom,
49                   const G4ErrorTrajErr&amp; errmat = G4ErrorTrajErr(5,0) );
50</pre></div><p>
51A particle trajectory is characterized by five independent variables as a function of one parameter (e.g. the path length).
52Among the five variables, one is related to the curvature (to the absolute value of the momentum), two are related to the direction of the particle and the other two are related to the spatial location.
53</p><p>
54There are two possible representations of these five parameters in the error propagator package: as a free trajectory state, class <code class="literal">G4ErrorTrajStateFree</code>, or as a trajectory state on a surface, class <code class="literal">G4ErrorTrajStateonSurface</code>.
55</p><div class="sect3" lang="en"><div class="titlepage"><div><div><h4 class="title"><a name="sect.G4Error.TrajState.FreTrajState"></a>5.8.2.1. 
56Free trajectory state
57</h4></div></div></div>
58
59In the free trajectory state representation the five trajectory parameters are
60
61<div class="itemizedlist"><ul type="disc"><li> G4double fInvP </li><li> G4double fLambda </li><li> G4double fPhi </li><li> G4double fYPerp </li><li> G4double fZPerp </li></ul></div><p>
62where <code class="literal">fInvP</code> is the inverse of the momentum.
63<code class="literal">fLambda</code> and <code class="literal">fPhi</code> are the dip and azimuthal angles related to the momentum components in the following way:
64</p><code class="literal">
65            p_x = p cos(lambda) cos(phi)
66            p_y = p cos(lambda) sin(phi)   
67            p_z = p sin(lambda)
68</code><p>
69that is, <code class="literal">lambda = 90 - theta</code>, where <code class="literal">theta</code> is the usual angle with respect to the Z axis.
70</p><p>
71<code class="literal">fYperp</code> and <code class="literal">fZperp</code> are the coordinates of the trajectory in a local orthonormal reference frame with the X axis along the particle direction, the Y axis being parallel to the X-Y plane (obtained by the vectorial product of the global Z axis and the momentum).
72</p></div><div class="sect3" lang="en"><div class="titlepage"><div><div><h4 class="title"><a name="sect.G4Error.TrajState.SurfaceTrajState"></a>5.8.2.2. 
73Trajectory state on a surface
74</h4></div></div></div><p>
75In the trajectory state on a surface representation the five trajectory parameters are
76</p><div class="itemizedlist"><ul type="disc"><li> G4double fInvP </li><li> G4double fPV </li><li> G4double fPW </li><li> G4double fV </li><li> G4double fW </li></ul></div><p>
77where <code class="literal">fInvP</code> is the inverse of the momentum;
78<code class="literal">fPV</code> and <code class="literal">fPW</code> are the momentum components in an orthonormal coordinate system with axis U, V and W;
79<code class="literal">fV</code> and <code class="literal">fW</code> are the position components on this coordinate system.
80</p>
81
82For this representation the user has to provide the plane where the parameters are calculated.
83This can be done by providing two vectors, V and W, contained in the plane:
84
85<div class="informalexample"><pre class="programlisting">
86  G4ErrorSurfaceTrajState( const G4String&amp; partType,
87                           const G4Point3D&amp; pos,
88                           const G4Vector3D&amp; mom,
89                           const G4Vector3D&amp; vecV,
90                           const G4Vector3D&amp; vecW,
91                           const G4ErrorTrajErr&amp; errmat = G4ErrorTrajErr(5,0) );
92</pre></div>
93
94 or by providing a plane
95
96<div class="informalexample"><pre class="programlisting">
97  G4ErrorSurfaceTrajState( const G4String&amp; partType,
98                           const G4Point3D&amp; pos,
99                           const G4Vector3D&amp; mom,
100                           const G4Plane3D&amp; plane,
101                           const G4ErrorTrajErr&amp; errmat = G4ErrorTrajErr(5,0) );
102</pre></div>
103
104In this second case the vector V is calculated as the vector in the plane perpendicular to the global vector X (if the plane normal is equal to X, Z is used instead) and W is calculated as the vector in the plane perpendicular to V.
105
106</div></div><div class="sect2" lang="en"><div class="titlepage"><div><div><h3 class="title"><a name="sect.G4Error.Err"></a>5.8.3. 
107Trajectory state error
108</h3></div></div></div><p>
109The 5X5 error matrix should also be provided at the creation of the trajectory state as a <code class="literal">G4ErrorTrajErr</code> object.
110If it is not provided a default object will be created filled with null values.
111</p><p>
112Currently the <code class="literal">G4ErrorTrajErr</code> is a <code class="literal">G4ErrorSymMatrix</code>, a simplified version of <code class="literal">CLHEP HepSymMatrix</code>.
113</p><p>
114The error matrix is given in units of GeV and cm. Therefore you should do the conversion if your code is using other units.
115</p></div><div class="sect2" lang="en"><div class="titlepage"><div><div><h3 class="title"><a name="sect.G4Error.Target"></a>5.8.4. 
116Targets
117</h3></div></div></div><p>
118The user has to define up to where the propagation must be done: the target.
119The target can be a surface <code class="literal">G4ErrorSurfaceTarget</code>, which is not part of the GEANT4 geometry.
120It can also be the surface of a GEANT4 volume  <code class="literal">G4ErrorGeomVolumeTarget</code>, so that the particle will be stopped when it enters this volume.
121Or it can be that the particle is stopped when a certain track length is reached, by implementing a <code class="literal">G4ErrorTrackLengthTarget</code>.
122</p><div class="sect3" lang="en"><div class="titlepage"><div><div><h4 class="title"><a name="sect.G4Error.Target.SurfaceTarget"></a>5.8.4.1. 
123Surface target
124</h4></div></div></div>
125
126When the user chooses a <code class="literal">G4ErrorSurfaceTarget</code> as target, the track is propagated until the surface is reached.
127This surface is not part of GEANT4 geometry, but usually traverses many GEANT4 volumes.
128The class <code class="literal">G4ErrorNavigator</code> takes care of the double navigation: for each step the step length is calculated as the minimum of the step length in the full geometry (up to a GEANT4 volume surface) and the distance to the user-defined surface.
129To do it, <code class="literal">G4ErrorNavigator</code> inherits from <code class="literal">G4Navigator</code> and overwrites the methods <code class="literal">ComputeStep()</code> and <code class="literal">ComputeSafety()</code>.
130
131Two types of surface are currently supported (more types could be easily implemented at user request): plane and cylindrical.
132
133<div class="sect4" lang="en"><div class="titlepage"><div><div><h5 class="title"><a name="sect.G4Error.Target.SurfaceTarget.PlaneSurfaceTarget"></a>5.8.4.1.1. 
134Plane surface target
135</h5></div></div></div><code class="literal">G4ErrorPlaneSurfaceTarget</code> implements an infinite plane surface.
136The surface can be given as the four coefficients of the plane equation <code class="literal">ax+by+cz+d = 0</code>:
137
138<div class="informalexample"><pre class="programlisting">
139    G4ErrorPlaneSurfaceTarget(G4double a=0,
140                              G4double b=0,
141                              G4double c=0,
142                              G4double d=0);
143</pre></div>
144
145or as the normal to the plane and a point contained in it:
146
147<div class="informalexample"><pre class="programlisting">
148    G4ErrorPlaneSurfaceTarget(const G4Normal3D &amp;n,
149                              const G4Point3D &amp;p);
150</pre></div>
151
152or as three points contained in it:
153
154<div class="informalexample"><pre class="programlisting">
155    G4ErrorPlaneSurfaceTarget(const G4Point3D &amp;p1,
156                              const G4Point3D &amp;p2,
157                              const G4Point3D &amp;p3);
158</pre></div></div><div class="sect4" lang="en"><div class="titlepage"><div><div><h5 class="title"><a name="sect.G4Error.Target.SurfaceTarget.CylSurfaceTarget"></a>5.8.4.1.2. 
159Cylindrical surface target
160</h5></div></div></div><code class="literal">G4ErrorCylSurfaceTarget</code> implements an infinite-length cylindrical surface (a cylinder without end-caps).
161The surface can be given as the radius, the translation and the rotation
162
163<div class="informalexample"><pre class="programlisting">
164    G4ErrorCylSurfaceTarget( const G4double&amp; radius,
165                             const G4ThreeVector&amp; trans=G4ThreeVector(),
166                             const G4RotationMatrix&amp; rotm=G4RotationMatrix() );
167</pre></div>
168
169or as the radius and the affine transformation
170
171<div class="informalexample"><pre class="programlisting">
172    G4ErrorCylSurfaceTarget( const G4double&amp; radius,
173                             const G4AffineTransform&amp; trans );
174</pre></div></div></div><div class="sect3" lang="en"><div class="titlepage"><div><div><h4 class="title"><a name="sect.G4Error.Target.VolumeTarget"></a>5.8.4.2. 
175Geometry volume target
176</h4></div></div></div><p>
177When the user chooses a <code class="literal">G4ErrorGeomVolumeTarget</code> as target, the track is propagated until the surface of a GEANT4 volume is reached. User can choose if the track will be stopped only when the track enters the volume, only when the track exits the volume or in both cases.
178</p>
179
180The object has to be instantiated giving the name of a logical volume existing in the geometry:
181
182<div class="informalexample"><pre class="programlisting">
183  G4ErrorGeomVolumeTarget( const G4String&amp; name );
184</pre></div></div><div class="sect3" lang="en"><div class="titlepage"><div><div><h4 class="title"><a name="sect.G4Error.Target.TrkLenTarget"></a>5.8.4.3. 
185Track Length target
186</h4></div></div></div><p>
187When the user chooses a <code class="literal">G4ErrorTrackLengthTarget</code> as target, the track is propagated until the given track length is reached.
188</p>
189
190The object has to be instantiated giving the value of the track length:
191<div class="informalexample"><pre class="programlisting">
192  G4ErrorTrackLengthTarget(const G4double maxTrkLength );
193</pre></div><p>
194It is implemented as a <code class="literal">G4VDiscreteProcess</code> and it limits the step in <code class="literal">PostStepGetPhysicalInteractionLength</code>.
195To ease its use, the process is registered to all particles in the constructor.
196</p></div></div><div class="sect2" lang="en"><div class="titlepage"><div><div><h3 class="title"><a name="sect.G4Error.Propagation"></a>5.8.5. 
197Managing the track propagation
198</h3></div></div></div><p>
199The user needs to propagate just one track, so there is no need of run and events. neither of <code class="literal">G4VPrimaryGeneratorAction</code>.
200
201<code class="literal">G4ErrorPropagator</code> creates a track from the information given in the <code class="literal">G4ErrorTrajState</code> and manages the step propagation.
202The propagation is done by the standard GEANT4 methods, invoking <code class="literal">G4SteppingManager::Stepping() </code> to propagate each step.
203</p><p>
204After one step is propagated, <code class="literal">G4ErrorPropagator</code> takes cares of propagating the track errors for this step, what is done by <code class="literal">G4ErrorTrajStateFree::PropagateError()</code>.
205The equations of error propagation are only implemented in the representation of  <code class="literal">G4ErrorTrajStateFree</code>.
206Therefore if the user has provided instead a <code class="literal">G4ErrorTrajStateOnSurface</code> object, it will be transformed into a <code class="literal">G4ErrorTrajStateFree</code> at the beginning of tracking, and at the end it is converted back into <code class="literal">G4ErrorTrajStateOnSurface</code> on the target surface (on the normal plane to the surface at the final point).
207</p><p>
208The user <code class="literal">G4VUserTrackingAction::PreUserTrackingAction( const G4Track* )</code> and  <code class="literal">G4VUserTrackingAction::PreUserTrackingAction( const G4Track* )</code> are also invoked at the beginning and at the end of the track propagation.
209</p><p>
210<code class="literal">G4ErrorPropagator</code> stops the tracking when one of the three conditions is true:
211
212</p><div class="itemizedlist"><ul type="disc"><li> Energy is exhausted </li><li> World boundary is reached </li><li> User-defined target is reached  </li></ul></div><p>
213
214In case the defined target is not reached, <code class="literal">G4ErrorPropagator::Propagate()</code> returns a negative value.
215</p><p>
216The propagation of a trajectory state until a user defined target can be done by invoking the method of <code class="literal">G4ErrorPropagatorManager</code>
217</p><div class="informalexample"><pre class="programlisting">
218  G4int Propagate( G4ErrorTrajState* currentTS, const G4ErrorTarget* target,
219                   G4ErrorMode mode = G4ErrorMode_PropForwards );
220</pre></div><p>
221</p>
222
223You can get the pointer to the only instance of <code class="literal">G4ErrorPropagatorManager</code> with
224
225<div class="informalexample"><pre class="programlisting">
226  G4ErrorPropagatorManager* g4emgr = G4ErrorPropagatorManager::GetErrorPropagatorManager();
227</pre></div>
228
229Another possibility is to invoke the propagation step by step, returning control to the user after each step. This can be done with the method
230
231<div class="informalexample"><pre class="programlisting">
232  G4int PropagateOneStep( G4ErrorTrajState* currentTS,
233                          G4ErrorMode mode = G4ErrorMode_PropForwards );
234</pre></div>
235
236In this case you should register the target first with the command
237
238<div class="informalexample"><pre class="programlisting">
239  G4ErrorPropagatorData::GetG4ErrorPropagatorData()-&gt;SetTarget( theG4eTarget );
240</pre></div><div class="sect3" lang="en"><div class="titlepage"><div><div><h4 class="title"><a name="sect.G4Error.Propagation.ErrorPropagation"></a>5.8.5.1. 
241Error propagation
242</h4></div></div></div>
243
244As in the GEANT3-based GEANE package, the error propagation is based on the equations of the European Muon Collaboration, that take into account:
245
246<div class="itemizedlist"><ul type="disc"><li> Error from curved trajectory in magnetic field </li><li> Error from multiple scattering </li><li> Error from ionization </li></ul></div><p>
247The formulas assume propagation along an helix.
248This means that it is necessary to make steps small enough to assure magnetic field constantness and not too big energy loss.
249</p></div></div><div class="sect2" lang="en"><div class="titlepage"><div><div><h3 class="title"><a name="sect.G4Error.StepLimits"></a>5.8.6. 
250Limiting the step
251</h3></div></div></div>
252
253There are three ways to limit the step.
254The first one is by using a fixed length value. This can be set by invoking the user command :
255
256<div class="informalexample"><pre class="programlisting">
257G4UImanager::GetUIpointer()-&gt;ApplyCommand("/geant4e/limits/stepLength MY_VALUE MY_UNIT");
258</pre></div>
259
260The second one is by setting the maximum percentage of energy loss in the step (or energy gain is propagation is backwards).
261This can be set by invoking the user command :
262
263<div class="informalexample"><pre class="programlisting">
264G4UImanager::GetUIpointer()-&gt;ApplyCommand("/geant4e/limits/energyLoss MY_VALUE");
265</pre></div>
266
267The last one is by setting the maximum difference between the value of the magnetic field at the beginning and at the end of the step. Indeed what is limited is the curvature, or exactly the value of the magnetic field divided by the value of the momentum transversal to the field.
268This can be set by invoking the user command :
269
270<div class="informalexample"><pre class="programlisting">
271G4UImanager::GetUIpointer()-&gt;ApplyCommand("/geant4e/limits/magField MY_VALUE");
272</pre></div>
273
274The classes that limit the step are implemented as GEANT4 processes. Therefore, the invocation of the above-mentioned commands should only be done after the initialization (for example after <code class="literal">G4ErrorPropagatorManager::InitGeant4e()</code>.
275
276</div></div><div class="navfooter"><hr><table width="100%" summary="Navigation footer"><tr><td width="40%" align="left"><a accesskey="p" href="ch05s07.html"><img src="AllResources/IconsGIF/prev.gif" alt="Prev"></a> </td><td width="20%" align="center"><a accesskey="u" href="ch05.html"><img src="AllResources/IconsGIF/up.gif" alt="Up"></a></td><td width="40%" align="right"> <a accesskey="n" href="ch06.html"><img src="AllResources/IconsGIF/next.gif" alt="Next"></a></td></tr><tr><td width="40%" align="left" valign="top">5.7. 
277User Limits
278 </td><td width="20%" align="center"><a accesskey="h" href="index.html"><img src="AllResources/IconsGIF/home.gif" alt="Home"></a></td><td width="40%" align="right" valign="top"> Chapter 6. 
279User Actions
280</td></tr></table></div></body></html>
Note: See TracBrowser for help on using the repository browser.