]>
Commit | Line | Data |
---|---|---|
fe4da5cc | 1 | * |
2 | * $Id$ | |
3 | * | |
4 | * $Log$ | |
5 | * Revision 1.1.1.1 1995/10/24 10:20:54 cernlib | |
6 | * Geant | |
7 | * | |
8 | * | |
9 | #include "geant321/pilot.h" | |
10 | *CMZ : 3.21/02 29/03/94 15.41.30 by S.Giani | |
11 | *-- Author : | |
12 | SUBROUTINE GNTRP (X, PAR, IACT, SNEXT, SNXT, SAFE) | |
13 | C. | |
14 | C. ****************************************************************** | |
15 | C. * * | |
16 | C. * COMPUTE DISTANCE UP TO INTERSECTION WITH 'TRAP' VOLUME, * | |
17 | C. * FROM INSIDE POINT X(1-3) ALONG DIRECTION X(4-6) * | |
18 | C. * * | |
19 | C. * PAR (input) : volume parameters * | |
20 | C. * IACT (input) : action flag * | |
21 | C. * = 0 Compute SAFE only * | |
22 | C. * = 1 Compute SAFE, and SNXT only if SNEXT .GT.new SAFE * | |
23 | C. * = 2 Compute both SAFE and SNXT * | |
24 | C. * = 3 Compute SNXT only * | |
25 | C. * SNEXT (input) : see IACT = 1 * | |
26 | C. * SNXT (output) : distance to volume boundary * | |
27 | C. * SAFE (output) : shortest distance to any boundary * | |
28 | C. * * | |
29 | C. * ==>Called by : GNEXT, GTNEXT * | |
30 | C. * Author A.McPherson ********* * | |
31 | C. * * | |
32 | C. ****************************************************************** | |
33 | C. | |
34 | #include "geant321/gconsp.inc" | |
35 | DIMENSION X(6), PAR(11) | |
36 | #if defined(CERNLIB_SINGLE) | |
37 | DIMENSION XT(6), P(11) | |
38 | #endif | |
39 | #if !defined(CERNLIB_SINGLE) | |
40 | DOUBLE PRECISION XT(6), P(11) | |
41 | DOUBLE PRECISION A,B,C,DSN,TOP1,TOP2,BOT1,BOT2 | |
42 | DOUBLE PRECISION DX1,DX2,X1,X2,X3,DCX1,DCX2,DX,DY | |
43 | DOUBLE PRECISION SP1,SP2,DCX,DCY,TANG1,TANG2,A1,B1,A2,B2 | |
44 | DOUBLE PRECISION XL1,XL2,XL3,XL4,CY,SY,SX,CX,C1,C2 | |
45 | DOUBLE PRECISION TX,H0,DHDZ,CH | |
46 | DOUBLE PRECISION CHTAN,SHTAN | |
47 | DOUBLE PRECISION SN1,SN2,SN3,SN4,SN5,SN6,SN7,SN8,P15 | |
48 | #endif | |
49 | C. | |
50 | C. -------------------------------------- | |
51 | C. | |
52 | SNXT = BIG | |
53 | SAFE = 0.0 | |
54 | DO 1 I=1,6 | |
55 | 1 XT(I)=X(I) | |
56 | DO 2 I=1,11 | |
57 | 2 P(I)=PAR(I) | |
58 | ||
59 | P15=0.5/P(1) | |
60 | DHDZ=(P(4)-P(8))*P15 | |
61 | ||
62 | IF (IACT .LT. 3) THEN | |
63 | ||
64 | ||
65 | C ------------------------------------------------- | |
66 | C | Compute safety-distance 'SAFE' (McPherson) | | |
67 | C ------------------------------------------------- | |
68 | ||
69 | ||
70 | C CALCULATE RECTANGLE ON FACE AT Z=-DZ. | |
71 | ||
72 | DX1=P(4)*P(7) | |
73 | X1=-DX1-P(5) | |
74 | X2=DX1-P(6) | |
75 | IF(X2.GT.X1) X1=X2 | |
76 | X2=-DX1+P(5) | |
77 | X3=DX1+P(6) | |
78 | IF(X3.LT.X2) X2=X3 | |
79 | DCX1=(X1+X2)*0.5 | |
80 | DX1=(X2-X1)*0.5 | |
81 | C | |
82 | C CALCULATE RECTANGLE ON FACE AT Z=DZ. | |
83 | C | |
84 | DX2=P(8)*P(11) | |
85 | X1=-DX2-P(9) | |
86 | X2=DX2-P(10) | |
87 | IF(X2.GT.X1) X1=X2 | |
88 | X2=-DX2+P(9) | |
89 | X3=DX2+P(10) | |
90 | IF(X3.LT.X2) X2=X3 | |
91 | DCX2=(X1+X2)*0.5 | |
92 | DX2=(X2-X1)*0.5 | |
93 | ||
94 | C | |
95 | TX=P(2)+(DCX2-DCX1)*P15 | |
96 | C | |
97 | C CALCULATE LOCAL RECTANGLE. | |
98 | C | |
99 | SP1=(P(1)-XT(3))*P15 | |
100 | SP2=(P(1)+XT(3))*P15 | |
101 | C | |
102 | DY=P(4)*SP1+P(8)*SP2 | |
103 | DX=DX1*SP1+DX2*SP2 | |
104 | DCX=(DCX1+DCX2)*0.5+TX*XT(3) | |
105 | DCY=P(3)*XT(3) | |
106 | C | |
107 | C CHECK POINT IS INSIDE RECTANGLE. | |
108 | C | |
109 | IF(ABS(XT(1)-DCX).GT.DX) GO TO 10 | |
110 | C | |
111 | C CALCULATE ANGLE OF YZ PLANES. | |
112 | C | |
113 | TANG1=TX+(DX1-DX2)*P15 | |
114 | TANG2=TX-(DX1-DX2)*P15 | |
115 | C | |
116 | C CALCULATE SAFETY FROM YZ PLANES. | |
117 | C | |
118 | SAF1=(XT(1)-DCX+DX)/SQRT(1.0+TANG1*TANG1) | |
119 | SAF2=(DX-XT(1)+DCX)/SQRT(1.0+TANG2*TANG2) | |
120 | C | |
121 | C CALCULATE ANGLE OF XZ PLANES. | |
122 | C | |
123 | TANG1=P(3)+DHDZ | |
124 | TANG2=P(3)-DHDZ | |
125 | C | |
126 | C CALCULATE SAFETY FROM XZ PLANES. | |
127 | C | |
128 | SAF3=(XT(2)-DCY+DY)/SQRT(1.0+TANG1*TANG1) | |
129 | SAF4=(DY-XT(2)+DCY)/SQRT(1.0+TANG2*TANG2) | |
130 | C | |
131 | C CALCULATE SAFETY FROM XY PLANES. | |
132 | C | |
133 | SAF5=P(1)-ABS(XT(3)) | |
134 | C | |
135 | C OVERALL SAFETY. | |
136 | C | |
137 | SAFE = MIN(SAF1,SAF2,SAF3,SAF4,SAF5) | |
138 | ||
139 | ||
140 | IF (IACT .EQ. 0) GO TO 999 | |
141 | IF (IACT .EQ. 1) THEN | |
142 | IF (SNEXT .LT. SAFE) GO TO 999 | |
143 | ENDIF | |
144 | ENDIF | |
145 | ||
146 | ||
147 | 10 CONTINUE | |
148 | ||
149 | C ------------------------------------------------ | |
150 | C | Compute vector-distance 'SNXT' (McPherson) | | |
151 | C ------------------------------------------------ | |
152 | ||
153 | ||
154 | ||
155 | C FIRST FIND S TO Z LIMITS. | |
156 | ||
157 | SN1=BIG | |
158 | SN2=BIG | |
159 | SN3=BIG | |
160 | SN4=BIG | |
161 | C | |
162 | IF(XT(6).EQ.0.0) GOTO 15 | |
163 | SN1=(P(1)-XT(3))/XT(6) | |
164 | SN2=-(P(1)+XT(3))/XT(6) | |
165 | IF(SN1.GE.0.) GOTO 15 | |
166 | SNX=SN1 | |
167 | SN1=SN2 | |
168 | SN2=SNX | |
169 | 15 CONTINUE | |
170 | C | |
171 | C NOW Y LIMITS. | |
172 | C | |
173 | H0=(P(4)+P(8))*0.5 | |
174 | C | |
175 | TOP1=H0-DHDZ*XT(3) | |
176 | TOP2=XT(2)-XT(3)*P(3) | |
177 | BOT1=XT(5)-XT(6)*P(3) | |
178 | BOT2=DHDZ*XT(6) | |
179 | C | |
180 | IF(BOT1+BOT2.NE.0.0) SN3=(TOP1-TOP2)/(BOT1+BOT2) | |
181 | IF(BOT2-BOT1.NE.0.0) SN4=(TOP1+TOP2)/(BOT2-BOT1) | |
182 | C | |
183 | C NOW X LIMTS. | |
184 | C | |
185 | XL1=(P(5)+P(6)+P(9)+P(10))*0.25 | |
186 | XL2=(P(5)-P(6)+P(9)-P(10))*0.25 | |
187 | XL3=(P(5)+P(6)-P(9)-P(10))*0.25 | |
188 | XL4=(P(5)-P(6)-P(9)+P(10))*0.25 | |
189 | C | |
190 | CY=XT(2)-XT(3)*P(3) | |
191 | SY=XT(5)-XT(6)*P(3) | |
192 | CH=H0-DHDZ*XT(3) | |
193 | C | |
194 | A1=XL4*XT(6)*SY+XL3*XT(6)*XT(6)*DHDZ | |
195 | B1=XL4*(XT(6)*CY+XT(3)*SY)-XL3*XT(6)*(CH-XT(3)*DHDZ)-XL2*P(1)*SY- | |
196 | +XL1*XT(6)*DHDZ*P(1) | |
197 | C1=XL4*CY*XT(3)-XL3*XT(3)*CH-XL2*CY*P(1)+XL1*CH*P(1) | |
198 | C | |
199 | CHTAN=(P(7)*P(4)+P(11)*P(8))*0.5 | |
200 | SHTAN=-(P(7)*P(4)-P(11)*P(8))*P15 | |
201 | CHTAN=CHTAN+SHTAN*XT(3) | |
202 | SHTAN=SHTAN*XT(6) | |
203 | CX=XT(1)-XT(3)*P(2) | |
204 | SX=XT(4)-XT(6)*P(2) | |
205 | C | |
206 | A2=-P(1)*(DHDZ*XT(6)*SX+SY*SHTAN) | |
207 | B2=P(1)*(CH*SX-DHDZ*XT(6)*CX-CY*SHTAN-SY*CHTAN) | |
208 | C2=P(1)*(CH*CX-CY*CHTAN) | |
209 | C | |
210 | ISN56=0 | |
211 | A=A1+A2 | |
212 | B=B1+B2 | |
213 | C=C1+C2 | |
214 | IF(B*B-A*C*4.0.LT.0.0) GO TO 20 | |
215 | IF(ABS(A).LT.1.0E-7) GO TO 19 | |
216 | DSN=SQRT(B*B-A*C*4.0) | |
217 | SN5=(-B+DSN)*0.5/A | |
218 | SN6=(-B-DSN)*0.5/A | |
219 | ISN56=1 | |
220 | GO TO 20 | |
221 | 19 CONTINUE | |
222 | IF(ABS(B).LT.1.0E-5) GO TO 20 | |
223 | ISN56=1 | |
224 | SN5=-C/B | |
225 | SN6=-C/B | |
226 | 20 CONTINUE | |
227 | C | |
228 | ISN78=0 | |
229 | A=A1-A2 | |
230 | B=B1-B2 | |
231 | C=C1-C2 | |
232 | IF(B*B-A*C*4.0.LT.0.0) GO TO 30 | |
233 | DSN=SQRT(B*B-A*C*4.0) | |
234 | IF(ABS(A).LT.1.0E-7) GO TO 25 | |
235 | SN7=(-B+DSN)*0.5/A | |
236 | SN8=(-B-DSN)*0.5/A | |
237 | ISN78=1 | |
238 | GO TO 30 | |
239 | 25 CONTINUE | |
240 | IF(ABS(B).LT.1.0E-5) GO TO 30 | |
241 | ISN78=1 | |
242 | SN7=-C/B | |
243 | SN8=-C/B | |
244 | 30 CONTINUE | |
245 | C | |
246 | IF(SN2.GT.0.0.AND.SN2.LT.SN1) SN1=SN2 | |
247 | IF(SN3.GT.0.0.AND.SN3.LT.SN1) SN1=SN3 | |
248 | IF(SN4.GT.0.0.AND.SN4.LT.SN1) SN1=SN4 | |
249 | IF(ISN56.EQ.0) GO TO 40 | |
250 | IF(SN5.GT.0.0.AND.SN5.LT.SN1) SN1=SN5 | |
251 | IF(SN6.GT.0.0.AND.SN6.LT.SN1) SN1=SN6 | |
252 | 40 CONTINUE | |
253 | IF(ISN78.EQ.0) GO TO 50 | |
254 | IF(SN7.GT.0.0.AND.SN7.LT.SN1) SN1=SN7 | |
255 | IF(SN8.GT.0.0.AND.SN8.LT.SN1) SN1=SN8 | |
256 | 50 CONTINUE | |
257 | C | |
258 | SNXT=SN1 | |
259 | C | |
260 | 999 CONTINUE | |
261 | END |