1 |
mocchiut |
1.1 |
|
2 |
pam-fi |
1.3 |
|
3 |
mocchiut |
1.1 |
************************************************************* |
4 |
|
|
* |
5 |
pam-fi |
1.2 |
* Routines to compute the NPOINT track intersection points |
6 |
mocchiut |
1.1 |
* with planes of z-coordinate given by ZIN |
7 |
|
|
* given the track parameters |
8 |
|
|
* |
9 |
|
|
* The routine is based on GRKUTA, which computes the |
10 |
|
|
* trajectory of a charged particle in a magnetic field |
11 |
pam-fi |
1.2 |
* by solving the equations of motion with Runge-Kuta method |
12 |
mocchiut |
1.1 |
* |
13 |
|
|
* Variables that have to be assigned when the subroutine |
14 |
|
|
* is called: |
15 |
|
|
* |
16 |
|
|
* ZIN(1-NPOINT) ----> z coordinates of the planes |
17 |
|
|
* AL_P(1-5) ----> track-parameter vector |
18 |
|
|
* |
19 |
|
|
* NB !!! |
20 |
|
|
* The routine works properly only if the |
21 |
|
|
* planes are numbered in descending order |
22 |
|
|
* |
23 |
pam-fi |
1.2 |
* Routines: |
24 |
mocchiut |
1.1 |
* |
25 |
pam-fi |
1.2 |
* DOTRACK(NPOINT,ZIN,XOUT,YOUT,AL_P,IFAIL) |
26 |
|
|
* - old routine |
27 |
|
|
* |
28 |
pam-fi |
1.3 |
* DOTRACK2(NPOINT,ZIN,XOUT,YOUT,THXOUT,THYOUT,TLOUT,AL_P,IFAIL) |
29 |
pam-fi |
1.4 |
* - as the previous one, but: |
30 |
pam-fi |
1.2 |
* -- the projected angles are given as output |
31 |
|
|
* -- the track lengths are given as output: |
32 |
|
|
* ---- for planes above the reference plane, the length until |
33 |
|
|
* the lower closest plane (reference plane included) |
34 |
|
|
* ---- for planes below the reference plane, the length until |
35 |
|
|
* the higher closest plane (reference plane included) |
36 |
mocchiut |
1.1 |
* |
37 |
pam-fi |
1.3 |
* March 2008 --> optimized by Paolo |
38 |
pam-fi |
1.4 |
* |
39 |
pam-ts |
1.5 |
* DOTRACK3(NPOINT,ZIN,XOUT,YOUT,THXOUT,THYOUT,TLOUT,AL_P,ZINI,IFAIL) |
40 |
pam-fi |
1.4 |
* - as the previous one, but with the coordinate ZINI of the reference |
41 |
|
|
* plane given as input parameter |
42 |
|
|
* |
43 |
|
|
* |
44 |
mocchiut |
1.1 |
************************************************************** |
45 |
|
|
|
46 |
pam-fi |
1.4 |
|
47 |
pam-fi |
1.2 |
SUBROUTINE |
48 |
pam-fi |
1.4 |
$ DOTRACK3(NPOINT,ZIN,XOUT,YOUT,THXOUT,THYOUT,TLOUT,AL_P |
49 |
|
|
$ ,ZINI,IFAIL) |
50 |
pam-fi |
1.2 |
|
51 |
|
|
IMPLICIT DOUBLE PRECISION (A-H,O-Z) |
52 |
|
|
|
53 |
pam-fi |
1.3 |
|
54 |
|
|
DIMENSION VECT(7),VECTINI(7),VOUT(7) |
55 |
pam-fi |
1.2 |
DATA TOLL/1.d-8/ |
56 |
pam-fi |
1.3 |
|
57 |
pam-fi |
1.2 |
* tolerance in reaching the next plane during the tracking procedure |
58 |
|
|
* ----------------------------------------------- |
59 |
|
|
* I/O parameters |
60 |
pam-ts |
1.5 |
PARAMETER (NPOINT_MAX=500) |
61 |
pam-fi |
1.2 |
DIMENSION ZIN(NPOINT_MAX) |
62 |
|
|
DIMENSION XOUT(NPOINT_MAX),YOUT(NPOINT_MAX) |
63 |
|
|
DIMENSION THXOUT(NPOINT_MAX),THYOUT(NPOINT_MAX) |
64 |
|
|
DIMENSION TLOUT(NPOINT_MAX) |
65 |
|
|
DIMENSION AL_P(5) |
66 |
|
|
* ----------------------------------------------- |
67 |
pam-fi |
1.4 |
* DATA ZINI/23.5/ !z coordinate of the reference plane |
68 |
pam-fi |
1.3 |
REAL*8 L |
69 |
|
|
|
70 |
pam-fi |
1.2 |
* ================================================================== |
71 |
|
|
* divide the track in two parts: below and above the reference plane |
72 |
|
|
* ================================================================== |
73 |
|
|
IUPDOWN=0 |
74 |
|
|
DO I=1,NPOINT |
75 |
|
|
IF(ZIN(I).LT.ZINI)THEN |
76 |
|
|
if(i.ne.1)IUPDOWN=I |
77 |
|
|
GOTO 88 |
78 |
|
|
ENDIF |
79 |
|
|
IUPDOWN=NPOINT+1 |
80 |
|
|
ENDDO |
81 |
|
|
88 CONTINUE |
82 |
|
|
|
83 |
|
|
* ================================================================== |
84 |
|
|
* track from reference plane DOWN |
85 |
|
|
* ================================================================== |
86 |
|
|
* parameters for GRKUTA |
87 |
|
|
IF(AL_P(5).NE.0) CHARGE=AL_P(5)/DABS(AL_P(5)) |
88 |
|
|
IF(AL_P(5).EQ.0) CHARGE=1. |
89 |
|
|
VOUT(1)=AL_P(1) |
90 |
|
|
VOUT(2)=AL_P(2) |
91 |
|
|
VOUT(3)=ZINI |
92 |
|
|
VOUT(4)=AL_P(3)*DCOS(AL_P(4)) |
93 |
|
|
VOUT(5)=AL_P(3)*DSIN(AL_P(4)) |
94 |
|
|
VOUT(6)=-1.*DSQRT(1.-AL_P(3)**2) |
95 |
|
|
IF(AL_P(5).NE.0.) VOUT(7)=DABS(1./AL_P(5)) |
96 |
|
|
IF(AL_P(5).EQ.0.) VOUT(7)=1.E8 |
97 |
pam-fi |
1.3 |
DO I=MAX(IUPDOWN,1),NPOINT |
98 |
|
|
L = 0.0 |
99 |
|
|
10 DO J=1,7 |
100 |
pam-fi |
1.2 |
VECT(J)=VOUT(J) |
101 |
|
|
VECTINI(J)=VOUT(J) |
102 |
|
|
ENDDO |
103 |
pam-fi |
1.3 |
IF(VOUT(6).GE.0.) THEN |
104 |
|
|
IFAIL=1 |
105 |
|
|
c$$$ if(TRKVERBOSE) |
106 |
|
|
c$$$ $ PRINT *,'dofit (grkuta): WARNING ===> backward track!!' |
107 |
|
|
RETURN |
108 |
|
|
ENDIF |
109 |
|
|
step=(zin(i)-vect(3))/VOUT(6) |
110 |
pam-fi |
1.2 |
11 continue |
111 |
pam-fi |
1.3 |
CALL GRKUTA(CHARGE,STEP,VECT,VOUT) |
112 |
|
|
L = L + STEP |
113 |
pam-fi |
1.2 |
IF(VOUT(3).GT.VECT(3)) THEN |
114 |
|
|
IFAIL=1 |
115 |
pam-fi |
1.3 |
c$$$ if(TRKVERBOSE) |
116 |
|
|
c$$$ $ PRINT *,'dofit (grkuta): WARNING ===> backward track!!' |
117 |
|
|
c$$$ if(TRKVERBOSE)print*,'charge',charge |
118 |
|
|
c$$$ if(TRKVERBOSE)print*,'vect',vect |
119 |
|
|
c$$$ if(TRKVERBOSE)print*,'vout',vout |
120 |
|
|
c$$$ if(TRKVERBOSE)print*,'step',step |
121 |
|
|
c$$$ if(TRKVERBOSE)print*,'DeltaB',DELTA0,DELTA1 |
122 |
pam-fi |
1.2 |
RETURN |
123 |
|
|
ENDIF |
124 |
|
|
Z=VOUT(3) |
125 |
|
|
IF(Z.LE.ZIN(I)+TOLL.AND.Z.GE.ZIN(I)-TOLL) GOTO 100 |
126 |
|
|
IF(Z.GT.ZIN(I)+TOLL) GOTO 10 |
127 |
|
|
IF(Z.LE.ZIN(I)-TOLL) THEN |
128 |
pam-fi |
1.3 |
L = L - STEP |
129 |
pam-fi |
1.2 |
STEP=STEP*(ZIN(I)-VECT(3))/(Z-VECT(3)) |
130 |
pam-fi |
1.3 |
DO J=1,7 |
131 |
pam-fi |
1.2 |
VECT(J)=VECTINI(J) |
132 |
|
|
ENDDO |
133 |
|
|
GOTO 11 |
134 |
|
|
ENDIF |
135 |
|
|
100 XOUT(I)=VOUT(1) |
136 |
|
|
YOUT(I)=VOUT(2) |
137 |
|
|
ZIN(I)=VOUT(3) |
138 |
|
|
IF(VOUT(3).ne.0)THEN |
139 |
|
|
THXOUT(I)=DATAN(VOUT(4)/VOUT(6))*180./ACOS(-1.) |
140 |
|
|
THYOUT(I)=DATAN(VOUT(5)/VOUT(6))*180./ACOS(-1.) |
141 |
|
|
ENDIF |
142 |
pam-fi |
1.3 |
TLOUT(I) = L |
143 |
pam-fi |
1.2 |
c print*,'D ',VECT(3),VOUT(3),VOUT(8),STEP |
144 |
|
|
ENDDO |
145 |
|
|
|
146 |
|
|
* ================================================================== |
147 |
|
|
* track from reference plane UP |
148 |
|
|
* ================================================================== |
149 |
|
|
* parameters for GRKUTA: |
150 |
|
|
* -opposite charge |
151 |
|
|
* -opposite momentum direction |
152 |
|
|
IF(AL_P(5).NE.0) CHARGE=-AL_P(5)/DABS(AL_P(5)) |
153 |
|
|
IF(AL_P(5).EQ.0) CHARGE=-1. |
154 |
|
|
VOUT(1)=AL_P(1) |
155 |
|
|
VOUT(2)=AL_P(2) |
156 |
|
|
VOUT(3)=ZINI |
157 |
|
|
VOUT(4)=-AL_P(3)*DCOS(AL_P(4)) |
158 |
|
|
VOUT(5)=-AL_P(3)*DSIN(AL_P(4)) |
159 |
|
|
VOUT(6)=1.*DSQRT(1.-AL_P(3)**2) |
160 |
|
|
IF(AL_P(5).NE.0.) VOUT(7)=DABS(1./AL_P(5)) |
161 |
|
|
IF(AL_P(5).EQ.0.) VOUT(7)=1.E8 |
162 |
|
|
DO I=MIN((IUPDOWN-1),NPOINT),1,-1 |
163 |
pam-fi |
1.3 |
L = 0.0 |
164 |
|
|
|
165 |
|
|
20 DO J=1,7 |
166 |
pam-fi |
1.2 |
VECT(J)=VOUT(J) |
167 |
|
|
VECTINI(J)=VOUT(J) |
168 |
|
|
ENDDO |
169 |
pam-fi |
1.3 |
IF(VOUT(6).LE.0.) THEN |
170 |
|
|
IFAIL=1 |
171 |
|
|
c$$$ if(TRKVERBOSE) |
172 |
|
|
c$$$ $ PRINT *,'dofit (grkuta): WARNING ===> backward track!!' |
173 |
|
|
RETURN |
174 |
|
|
ENDIF |
175 |
|
|
step=(zin(i)-vect(3))/VOUT(6) |
176 |
pam-fi |
1.2 |
22 continue |
177 |
pam-fi |
1.3 |
CALL GRKUTA(CHARGE,STEP,VECT,VOUT) |
178 |
|
|
L = L + STEP |
179 |
pam-fi |
1.2 |
IF(VOUT(3).LT.VECT(3)) THEN |
180 |
|
|
IFAIL=1 |
181 |
pam-fi |
1.3 |
c$$$ if(TRKVERBOSE) |
182 |
|
|
c$$$ $ PRINT *,'dofit (grkuta): WARNING ===> backward track!!' |
183 |
|
|
c$$$ if(TRKVERBOSE)print*,'charge',charge |
184 |
|
|
c$$$ if(TRKVERBOSE)print*,'vect',vect |
185 |
|
|
c$$$ if(TRKVERBOSE)print*,'vout',vout |
186 |
|
|
c$$$ if(TRKVERBOSE)print*,'step',step |
187 |
|
|
c$$$ if(TRKVERBOSE)print*,'DeltaB',DELTA0,DELTA1 |
188 |
pam-fi |
1.2 |
RETURN |
189 |
|
|
ENDIF |
190 |
|
|
Z=VOUT(3) |
191 |
|
|
IF(Z.LE.ZIN(I)+TOLL.AND.Z.GE.ZIN(I)-TOLL) GOTO 200 |
192 |
|
|
IF(Z.LT.ZIN(I)-TOLL) GOTO 20 |
193 |
|
|
IF(Z.GE.ZIN(I)+TOLL) THEN |
194 |
pam-fi |
1.3 |
L = L - STEP |
195 |
pam-fi |
1.2 |
STEP=STEP*(ZIN(I)-VECT(3))/(Z-VECT(3)) |
196 |
pam-fi |
1.3 |
DO J=1,7 |
197 |
pam-fi |
1.2 |
VECT(J)=VECTINI(J) |
198 |
|
|
ENDDO |
199 |
|
|
GOTO 22 |
200 |
|
|
ENDIF |
201 |
|
|
200 XOUT(I)=VOUT(1) |
202 |
|
|
YOUT(I)=VOUT(2) |
203 |
|
|
ZIN(I)=VOUT(3) |
204 |
|
|
IF(VOUT(3).ne.0)THEN |
205 |
|
|
THXOUT(I)=DATAN(VOUT(4)/VOUT(6))*180./ACOS(-1.) |
206 |
|
|
THYOUT(I)=DATAN(VOUT(5)/VOUT(6))*180./ACOS(-1.) |
207 |
|
|
ENDIF |
208 |
pam-fi |
1.3 |
TLOUT(I) = L |
209 |
pam-fi |
1.2 |
ENDDO |
210 |
|
|
|
211 |
|
|
RETURN |
212 |
|
|
END |
213 |
pam-fi |
1.4 |
************************************************************************ |
214 |
|
|
* |
215 |
|
|
* |
216 |
|
|
* |
217 |
|
|
* |
218 |
|
|
************************************************************************ |
219 |
|
|
|
220 |
|
|
SUBROUTINE |
221 |
|
|
$ DOTRACK2(NPOINT,ZIN,XOUT,YOUT,THXOUT,THYOUT,TLOUT,AL_P,IFAIL) |
222 |
|
|
|
223 |
|
|
IMPLICIT DOUBLE PRECISION (A-H,O-Z) |
224 |
|
|
|
225 |
|
|
|
226 |
|
|
ccc DIMENSION VECT(7),VECTINI(7),VOUT(7) |
227 |
|
|
ccc DATA TOLL/1.d-8/ |
228 |
|
|
|
229 |
|
|
* tolerance in reaching the next plane during the tracking procedure |
230 |
|
|
* ----------------------------------------------- |
231 |
|
|
* I/O parameters |
232 |
|
|
PARAMETER (NPOINT_MAX=100) |
233 |
|
|
DIMENSION ZIN(NPOINT_MAX) |
234 |
|
|
DIMENSION XOUT(NPOINT_MAX),YOUT(NPOINT_MAX) |
235 |
|
|
DIMENSION THXOUT(NPOINT_MAX),THYOUT(NPOINT_MAX) |
236 |
|
|
DIMENSION TLOUT(NPOINT_MAX) |
237 |
|
|
DIMENSION AL_P(5) |
238 |
|
|
* ----------------------------------------------- |
239 |
|
|
DATA ZINI/23.5/ !z coordinate of the reference plane |
240 |
|
|
ccc REAL*8 L |
241 |
|
|
|
242 |
|
|
CALL DOTRACK3(NPOINT,ZIN,XOUT,YOUT,THXOUT,THYOUT,TLOUT,AL_P |
243 |
|
|
$ ,ZINI,IFAIL) |
244 |
|
|
|
245 |
|
|
ccc |
246 |
|
|
ccc OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD |
247 |
|
|
ccc |
248 |
|
|
|
249 |
|
|
c$$$* ================================================================== |
250 |
|
|
c$$$* divide the track in two parts: below and above the reference plane |
251 |
|
|
c$$$* ================================================================== |
252 |
|
|
c$$$ IUPDOWN=0 |
253 |
|
|
c$$$ DO I=1,NPOINT |
254 |
|
|
c$$$ IF(ZIN(I).LT.ZINI)THEN |
255 |
|
|
c$$$ if(i.ne.1)IUPDOWN=I |
256 |
|
|
c$$$ GOTO 88 |
257 |
|
|
c$$$ ENDIF |
258 |
|
|
c$$$ IUPDOWN=NPOINT+1 |
259 |
|
|
c$$$ ENDDO |
260 |
|
|
c$$$ 88 CONTINUE |
261 |
|
|
c$$$ |
262 |
|
|
c$$$* ================================================================== |
263 |
|
|
c$$$* track from reference plane DOWN |
264 |
|
|
c$$$* ================================================================== |
265 |
|
|
c$$$* parameters for GRKUTA |
266 |
|
|
c$$$ IF(AL_P(5).NE.0) CHARGE=AL_P(5)/DABS(AL_P(5)) |
267 |
|
|
c$$$ IF(AL_P(5).EQ.0) CHARGE=1. |
268 |
|
|
c$$$ VOUT(1)=AL_P(1) |
269 |
|
|
c$$$ VOUT(2)=AL_P(2) |
270 |
|
|
c$$$ VOUT(3)=ZINI |
271 |
|
|
c$$$ VOUT(4)=AL_P(3)*DCOS(AL_P(4)) |
272 |
|
|
c$$$ VOUT(5)=AL_P(3)*DSIN(AL_P(4)) |
273 |
|
|
c$$$ VOUT(6)=-1.*DSQRT(1.-AL_P(3)**2) |
274 |
|
|
c$$$ IF(AL_P(5).NE.0.) VOUT(7)=DABS(1./AL_P(5)) |
275 |
|
|
c$$$ IF(AL_P(5).EQ.0.) VOUT(7)=1.E8 |
276 |
|
|
c$$$ DO I=MAX(IUPDOWN,1),NPOINT |
277 |
|
|
c$$$ L = 0.0 |
278 |
|
|
c$$$ 10 DO J=1,7 |
279 |
|
|
c$$$ VECT(J)=VOUT(J) |
280 |
|
|
c$$$ VECTINI(J)=VOUT(J) |
281 |
|
|
c$$$ ENDDO |
282 |
|
|
c$$$ IF(VOUT(6).GE.0.) THEN |
283 |
|
|
c$$$ IFAIL=1 |
284 |
|
|
c$$$c$$$ if(TRKVERBOSE) |
285 |
|
|
c$$$c$$$ $ PRINT *,'dofit (grkuta): WARNING ===> backward track!!' |
286 |
|
|
c$$$ RETURN |
287 |
|
|
c$$$ ENDIF |
288 |
|
|
c$$$ step=(zin(i)-vect(3))/VOUT(6) |
289 |
|
|
c$$$ 11 continue |
290 |
|
|
c$$$ CALL GRKUTA(CHARGE,STEP,VECT,VOUT) |
291 |
|
|
c$$$ L = L + STEP |
292 |
|
|
c$$$ IF(VOUT(3).GT.VECT(3)) THEN |
293 |
|
|
c$$$ IFAIL=1 |
294 |
|
|
c$$$c$$$ if(TRKVERBOSE) |
295 |
|
|
c$$$c$$$ $ PRINT *,'dofit (grkuta): WARNING ===> backward track!!' |
296 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'charge',charge |
297 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'vect',vect |
298 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'vout',vout |
299 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'step',step |
300 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'DeltaB',DELTA0,DELTA1 |
301 |
|
|
c$$$ RETURN |
302 |
|
|
c$$$ ENDIF |
303 |
|
|
c$$$ Z=VOUT(3) |
304 |
|
|
c$$$ IF(Z.LE.ZIN(I)+TOLL.AND.Z.GE.ZIN(I)-TOLL) GOTO 100 |
305 |
|
|
c$$$ IF(Z.GT.ZIN(I)+TOLL) GOTO 10 |
306 |
|
|
c$$$ IF(Z.LE.ZIN(I)-TOLL) THEN |
307 |
|
|
c$$$ L = L - STEP |
308 |
|
|
c$$$ STEP=STEP*(ZIN(I)-VECT(3))/(Z-VECT(3)) |
309 |
|
|
c$$$ DO J=1,7 |
310 |
|
|
c$$$ VECT(J)=VECTINI(J) |
311 |
|
|
c$$$ ENDDO |
312 |
|
|
c$$$ GOTO 11 |
313 |
|
|
c$$$ ENDIF |
314 |
|
|
c$$$ 100 XOUT(I)=VOUT(1) |
315 |
|
|
c$$$ YOUT(I)=VOUT(2) |
316 |
|
|
c$$$ ZIN(I)=VOUT(3) |
317 |
|
|
c$$$ IF(VOUT(3).ne.0)THEN |
318 |
|
|
c$$$ THXOUT(I)=DATAN(VOUT(4)/VOUT(6))*180./ACOS(-1.) |
319 |
|
|
c$$$ THYOUT(I)=DATAN(VOUT(5)/VOUT(6))*180./ACOS(-1.) |
320 |
|
|
c$$$ ENDIF |
321 |
|
|
c$$$ TLOUT(I) = L |
322 |
|
|
c$$$c print*,'D ',VECT(3),VOUT(3),VOUT(8),STEP |
323 |
|
|
c$$$ ENDDO |
324 |
|
|
c$$$ |
325 |
|
|
c$$$* ================================================================== |
326 |
|
|
c$$$* track from reference plane UP |
327 |
|
|
c$$$* ================================================================== |
328 |
|
|
c$$$* parameters for GRKUTA: |
329 |
|
|
c$$$* -opposite charge |
330 |
|
|
c$$$* -opposite momentum direction |
331 |
|
|
c$$$ IF(AL_P(5).NE.0) CHARGE=-AL_P(5)/DABS(AL_P(5)) |
332 |
|
|
c$$$ IF(AL_P(5).EQ.0) CHARGE=-1. |
333 |
|
|
c$$$ VOUT(1)=AL_P(1) |
334 |
|
|
c$$$ VOUT(2)=AL_P(2) |
335 |
|
|
c$$$ VOUT(3)=ZINI |
336 |
|
|
c$$$ VOUT(4)=-AL_P(3)*DCOS(AL_P(4)) |
337 |
|
|
c$$$ VOUT(5)=-AL_P(3)*DSIN(AL_P(4)) |
338 |
|
|
c$$$ VOUT(6)=1.*DSQRT(1.-AL_P(3)**2) |
339 |
|
|
c$$$ IF(AL_P(5).NE.0.) VOUT(7)=DABS(1./AL_P(5)) |
340 |
|
|
c$$$ IF(AL_P(5).EQ.0.) VOUT(7)=1.E8 |
341 |
|
|
c$$$ DO I=MIN((IUPDOWN-1),NPOINT),1,-1 |
342 |
|
|
c$$$ L = 0.0 |
343 |
|
|
c$$$ |
344 |
|
|
c$$$ 20 DO J=1,7 |
345 |
|
|
c$$$ VECT(J)=VOUT(J) |
346 |
|
|
c$$$ VECTINI(J)=VOUT(J) |
347 |
|
|
c$$$ ENDDO |
348 |
|
|
c$$$ IF(VOUT(6).LE.0.) THEN |
349 |
|
|
c$$$ IFAIL=1 |
350 |
|
|
c$$$c$$$ if(TRKVERBOSE) |
351 |
|
|
c$$$c$$$ $ PRINT *,'dofit (grkuta): WARNING ===> backward track!!' |
352 |
|
|
c$$$ RETURN |
353 |
|
|
c$$$ ENDIF |
354 |
|
|
c$$$ step=(zin(i)-vect(3))/VOUT(6) |
355 |
|
|
c$$$ 22 continue |
356 |
|
|
c$$$ CALL GRKUTA(CHARGE,STEP,VECT,VOUT) |
357 |
|
|
c$$$ L = L + STEP |
358 |
|
|
c$$$ IF(VOUT(3).LT.VECT(3)) THEN |
359 |
|
|
c$$$ IFAIL=1 |
360 |
|
|
c$$$c$$$ if(TRKVERBOSE) |
361 |
|
|
c$$$c$$$ $ PRINT *,'dofit (grkuta): WARNING ===> backward track!!' |
362 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'charge',charge |
363 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'vect',vect |
364 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'vout',vout |
365 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'step',step |
366 |
|
|
c$$$c$$$ if(TRKVERBOSE)print*,'DeltaB',DELTA0,DELTA1 |
367 |
|
|
c$$$ RETURN |
368 |
|
|
c$$$ ENDIF |
369 |
|
|
c$$$ Z=VOUT(3) |
370 |
|
|
c$$$ IF(Z.LE.ZIN(I)+TOLL.AND.Z.GE.ZIN(I)-TOLL) GOTO 200 |
371 |
|
|
c$$$ IF(Z.LT.ZIN(I)-TOLL) GOTO 20 |
372 |
|
|
c$$$ IF(Z.GE.ZIN(I)+TOLL) THEN |
373 |
|
|
c$$$ L = L - STEP |
374 |
|
|
c$$$ STEP=STEP*(ZIN(I)-VECT(3))/(Z-VECT(3)) |
375 |
|
|
c$$$ DO J=1,7 |
376 |
|
|
c$$$ VECT(J)=VECTINI(J) |
377 |
|
|
c$$$ ENDDO |
378 |
|
|
c$$$ GOTO 22 |
379 |
|
|
c$$$ ENDIF |
380 |
|
|
c$$$ 200 XOUT(I)=VOUT(1) |
381 |
|
|
c$$$ YOUT(I)=VOUT(2) |
382 |
|
|
c$$$ ZIN(I)=VOUT(3) |
383 |
|
|
c$$$ IF(VOUT(3).ne.0)THEN |
384 |
|
|
c$$$ THXOUT(I)=DATAN(VOUT(4)/VOUT(6))*180./ACOS(-1.) |
385 |
|
|
c$$$ THYOUT(I)=DATAN(VOUT(5)/VOUT(6))*180./ACOS(-1.) |
386 |
|
|
c$$$ ENDIF |
387 |
|
|
c$$$ TLOUT(I) = L |
388 |
|
|
c$$$ ENDDO |
389 |
|
|
c$$$ |
390 |
|
|
RETURN |
391 |
|
|
END |
392 |
pam-fi |
1.2 |
|
393 |
|
|
************************************************************************ |
394 |
|
|
* |
395 |
|
|
* |
396 |
|
|
* |
397 |
|
|
* |
398 |
|
|
************************************************************************ |
399 |
mocchiut |
1.1 |
SUBROUTINE DOTRACK(NPOINT,ZIN,XOUT,YOUT,AL_P,IFAIL) |
400 |
|
|
|
401 |
|
|
IMPLICIT DOUBLE PRECISION (A-H,O-Z) |
402 |
|
|
|
403 |
|
|
* ----------------------------------------------- |
404 |
|
|
* I/O parameters |
405 |
|
|
PARAMETER (NPOINT_MAX=100) |
406 |
|
|
DIMENSION ZIN(NPOINT_MAX) |
407 |
|
|
DIMENSION XOUT(NPOINT_MAX),YOUT(NPOINT_MAX) |
408 |
pam-fi |
1.3 |
DIMENSION THXOUT(NPOINT_MAX),THYOUT(NPOINT_MAX) |
409 |
|
|
DIMENSION TLOUT(NPOINT_MAX) |
410 |
mocchiut |
1.1 |
DIMENSION AL_P(5) |
411 |
|
|
* ----------------------------------------------- |
412 |
|
|
|
413 |
pam-fi |
1.3 |
CALL DOTRACK2(NPOINT,ZIN,XOUT,YOUT,THXOUT,THYOUT,TLOUT,AL_P,IFAIL) |
414 |
mocchiut |
1.1 |
|
415 |
pam-fi |
1.3 |
ccc |
416 |
|
|
ccc OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD-OLD |
417 |
|
|
ccc |
418 |
|
|
c$$$ DIMENSION VECT(7),VECTINI(7),VOUT(7) |
419 |
|
|
c$$$ DATA TOLL/1.d-8/ |
420 |
|
|
c$$$* tolerance in reaching the next plane during the tracking procedure |
421 |
|
|
c$$$* ----------------------------------------------- |
422 |
|
|
c$$$* I/O parameters |
423 |
|
|
c$$$ PARAMETER (NPOINT_MAX=100) |
424 |
|
|
c$$$ DIMENSION ZIN(NPOINT_MAX) |
425 |
|
|
c$$$ DIMENSION XOUT(NPOINT_MAX),YOUT(NPOINT_MAX) |
426 |
|
|
c$$$ DIMENSION AL_P(5) |
427 |
|
|
c$$$* ----------------------------------------------- |
428 |
|
|
c$$$ DATA ZINI/23.5/ !z coordinate of the reference plane |
429 |
|
|
c$$$ |
430 |
|
|
c$$$* ================================================================== |
431 |
|
|
c$$$* divide the track in two parts: below and above the reference plane |
432 |
|
|
c$$$* ================================================================== |
433 |
|
|
c$$$ IUPDOWN=0 |
434 |
|
|
c$$$ DO I=1,NPOINT |
435 |
|
|
c$$$ IF(ZIN(I).LT.ZINI)THEN |
436 |
|
|
c$$$ if(i.ne.1)IUPDOWN=I |
437 |
|
|
c$$$ GOTO 88 |
438 |
|
|
c$$$ ENDIF |
439 |
|
|
c$$$ IUPDOWN=NPOINT+1 |
440 |
|
|
c$$$ ENDDO |
441 |
|
|
c$$$ 88 CONTINUE |
442 |
|
|
c$$$ |
443 |
|
|
c$$$* ================================================================== |
444 |
|
|
c$$$* track from reference plane DOWN |
445 |
|
|
c$$$* ================================================================== |
446 |
|
|
c$$$* parameters for GRKUTA |
447 |
|
|
c$$$ IF(AL_P(5).NE.0) CHARGE=AL_P(5)/DABS(AL_P(5)) |
448 |
|
|
c$$$ IF(AL_P(5).EQ.0) CHARGE=1. |
449 |
|
|
c$$$ VOUT(1)=AL_P(1) |
450 |
|
|
c$$$ VOUT(2)=AL_P(2) |
451 |
|
|
c$$$ VOUT(3)=ZINI |
452 |
|
|
c$$$ VOUT(4)=AL_P(3)*DCOS(AL_P(4)) |
453 |
|
|
c$$$ VOUT(5)=AL_P(3)*DSIN(AL_P(4)) |
454 |
|
|
c$$$ VOUT(6)=-1.*DSQRT(1.-AL_P(3)**2) |
455 |
|
|
c$$$ IF(AL_P(5).NE.0.) VOUT(7)=DABS(1./AL_P(5)) |
456 |
|
|
c$$$ IF(AL_P(5).EQ.0.) VOUT(7)=1.E8 |
457 |
|
|
c$$$ DO I=MAX(IUPDOWN,1),NPOINT |
458 |
|
|
c$$$ step=vout(3)-zin(i) |
459 |
|
|
c$$$c$$$ print*,'DOWN ',i,' - Track from ', |
460 |
|
|
c$$$c$$$ $ vout(3),' to ',zin(i),' step ',step |
461 |
|
|
c$$$ 10 DO J=1,7 |
462 |
|
|
c$$$ VECT(J)=VOUT(J) |
463 |
|
|
c$$$ VECTINI(J)=VOUT(J) |
464 |
|
|
c$$$ ENDDO |
465 |
|
|
c$$$ 11 continue |
466 |
|
|
c$$$ CALL GRKUTA(CHARGE,STEP,VECT,VOUT) |
467 |
|
|
c$$$ IF(VOUT(3).GT.VECT(3)) THEN |
468 |
|
|
c$$$ IFAIL=1 |
469 |
|
|
c$$$c$$$ PRINT *,'=== WARNING ===> tracciamento invertito (DOWN)' |
470 |
|
|
c$$$c$$$ print*,'charge',charge |
471 |
|
|
c$$$c$$$ print*,'vect',vect |
472 |
|
|
c$$$c$$$ print*,'vout',vout |
473 |
|
|
c$$$c$$$ print*,'step',step |
474 |
|
|
c$$$ RETURN |
475 |
|
|
c$$$ ENDIF |
476 |
|
|
c$$$ Z=VOUT(3) |
477 |
|
|
c$$$ IF(Z.LE.ZIN(I)+TOLL.AND.Z.GE.ZIN(I)-TOLL) GOTO 100 |
478 |
|
|
c$$$ IF(Z.GT.ZIN(I)+TOLL) GOTO 10 |
479 |
|
|
c$$$ IF(Z.LE.ZIN(I)-TOLL) THEN |
480 |
|
|
c$$$ STEP=STEP*(ZIN(I)-VECT(3))/(Z-VECT(3)) |
481 |
|
|
c$$$ DO J=1,7 |
482 |
|
|
c$$$ VECT(J)=VECTINI(J) |
483 |
|
|
c$$$ ENDDO |
484 |
|
|
c$$$ GOTO 11 |
485 |
|
|
c$$$ ENDIF |
486 |
|
|
c$$$ 100 XOUT(I)=VOUT(1) |
487 |
|
|
c$$$ YOUT(I)=VOUT(2) |
488 |
|
|
c$$$ ZIN(I)=VOUT(3) |
489 |
|
|
c$$$* THXOUT(I)= |
490 |
|
|
c$$$* THYOUT(I)= |
491 |
|
|
c$$$ ENDDO |
492 |
|
|
c$$$ |
493 |
|
|
c$$$ |
494 |
|
|
c$$$ |
495 |
|
|
c$$$* ================================================================== |
496 |
|
|
c$$$* track from reference plane UP |
497 |
|
|
c$$$* ================================================================== |
498 |
|
|
c$$$* parameters for GRKUTA: |
499 |
|
|
c$$$* -opposite charge |
500 |
|
|
c$$$* -opposite momentum direction |
501 |
|
|
c$$$ IF(AL_P(5).NE.0) CHARGE=-AL_P(5)/DABS(AL_P(5)) |
502 |
|
|
c$$$ IF(AL_P(5).EQ.0) CHARGE=-1. |
503 |
|
|
c$$$ VOUT(1)=AL_P(1) |
504 |
|
|
c$$$ VOUT(2)=AL_P(2) |
505 |
|
|
c$$$ VOUT(3)=ZINI |
506 |
|
|
c$$$ VOUT(4)=-AL_P(3)*DCOS(AL_P(4)) |
507 |
|
|
c$$$ VOUT(5)=-AL_P(3)*DSIN(AL_P(4)) |
508 |
|
|
c$$$ VOUT(6)=1.*DSQRT(1.-AL_P(3)**2) |
509 |
|
|
c$$$ IF(AL_P(5).NE.0.) VOUT(7)=DABS(1./AL_P(5)) |
510 |
|
|
c$$$ IF(AL_P(5).EQ.0.) VOUT(7)=1.E8 |
511 |
|
|
c$$$ DO I=MIN((IUPDOWN-1),NPOINT),1,-1 |
512 |
|
|
c$$$ step=vout(3)-zin(i) |
513 |
|
|
c$$$ step = -step |
514 |
|
|
c$$$c$$$ print*,'UP ',i,' - Track from ', |
515 |
|
|
c$$$c$$$ $ vout(3),' to ',zin(i),' step ',step |
516 |
|
|
c$$$ 20 DO J=1,7 |
517 |
|
|
c$$$ VECT(J)=VOUT(J) |
518 |
|
|
c$$$ VECTINI(J)=VOUT(J) |
519 |
|
|
c$$$ ENDDO |
520 |
|
|
c$$$ 22 continue |
521 |
|
|
c$$$ CALL GRKUTA(CHARGE,STEP,VECT,VOUT) |
522 |
|
|
c$$$ IF(VOUT(3).LT.VECT(3)) THEN |
523 |
|
|
c$$$ IFAIL=1 |
524 |
|
|
c$$$c$$$ PRINT *,'=== WARNING ===> tracciamento invertito (UP)' |
525 |
|
|
c$$$c$$$ print*,'charge',charge |
526 |
|
|
c$$$c$$$ print*,'vect',vect |
527 |
|
|
c$$$c$$$ print*,'vout',vout |
528 |
|
|
c$$$c$$$ print*,'step',step |
529 |
|
|
c$$$ RETURN |
530 |
|
|
c$$$ ENDIF |
531 |
|
|
c$$$ Z=VOUT(3) |
532 |
|
|
c$$$ IF(Z.LE.ZIN(I)+TOLL.AND.Z.GE.ZIN(I)-TOLL) GOTO 200 |
533 |
|
|
c$$$ IF(Z.LT.ZIN(I)-TOLL) GOTO 20 |
534 |
|
|
c$$$ IF(Z.GE.ZIN(I)+TOLL) THEN |
535 |
|
|
c$$$ STEP=STEP*(ZIN(I)-VECT(3))/(Z-VECT(3)) |
536 |
|
|
c$$$ DO J=1,7 |
537 |
|
|
c$$$ VECT(J)=VECTINI(J) |
538 |
|
|
c$$$ ENDDO |
539 |
|
|
c$$$ GOTO 22 |
540 |
|
|
c$$$ ENDIF |
541 |
|
|
c$$$ 200 XOUT(I)=VOUT(1) |
542 |
|
|
c$$$ YOUT(I)=VOUT(2) |
543 |
|
|
c$$$ ZIN(I)=VOUT(3) |
544 |
|
|
c$$$ |
545 |
|
|
c$$$ ENDDO |
546 |
mocchiut |
1.1 |
|
547 |
|
|
RETURN |
548 |
|
|
END |