Bugfix in AliPoints2Memory
[u/mrichter/AliRoot.git] / GEANT321 / ggeom / gflrad.F
CommitLineData
fe4da5cc 1*
2* $Id$
3*
4* $Log$
5* Revision 1.1.1.1 1995/10/24 10:20:48 cernlib
6* Geant
7*
8*
9#include "geant321/pilot.h"
10*CMZ : 3.21/02 29/03/94 15.41.28 by S.Giani
11*-- Author :
12 SUBROUTINE GFLRAD(IAXIS,ISH,IROT,DX,PARS,CL,CH,IERR)
13C.
14C. ******************************************************************
15C. * *
16C. * ROUTINE TO COMPUTE THE LIMITS IN R FOR THE SHAPE ISH *
17C. * DISPLACED BY THE VECTOR DX AND ROTATED BY THE MATRIX IROT. *
18C. * IF IAXIS = 4 THE R IS THE XY PLANE R, IF IAXIS = 5 IT IS *
19C. * THE 3 DINEMSIONAL SPACE R. THE SHAPE HAS NPAR PARAMETERS *
20C. * IN THE ARRAY PARS. THE LOWER LIMIT IS RETURNED IN CL AND *
21C. * THE HIGHER IN CH. IF THE CALCULATION CANNOT BE PERFORMED *
22C. * IERR IS SET TO 1 OTHERWISE IT IS SET TO 0. *
23C. * *
24C. * ==>Called by : GFCLIM *
25C. * Author A.McPherson ********* *
26C. * *
27C. ******************************************************************
28C.
29#include "geant321/gcbank.inc"
30#include "geant321/gconsp.inc"
31#include "geant321/gcshno.inc"
32 DIMENSION DX(3),PARS(11),X(3),XT(3)
33C.
34C. --------------------------------------------------
35C.
36 IERR=1
37C
38C FIRST CALCULATE THE LENGTH OF THE DISPLACEMENT OF THE
39C ORIGIN.
40C
41 DXS=DX(1)*DX(1)+DX(2)*DX(2)
42 IF(IAXIS.EQ.5) DXS=DXS+DX(3)*DX(3)
43 IF(DXS.GT.0.0) DXS=SQRT(DXS)
44C
45 IF(ISH.GT.4.AND.ISH.NE.10.AND.ISH.NE.28) GO TO 40
46C
47C CUBOIDS, TRAPEZOIDS, PARALLELEPIPEDS.
48C
49 CH=0.0
50 CL=DXS
51C
52 DO 30 IP=1,8
53C
54C THIS IS A LOOP OVER THE 8 CORNERS.
55C FIRST FIND THE LOCAL COORDINATES.
56C
57 IF(ISH.EQ.28) THEN
58C
59C General twisted trapezoid.
60C
61 IL=(IP+1)/2
62 I0=IL*4+11
63 IS=(IP-IL*2)*2+1
64 X(3)=PARS(1)*IS
65 X(1)=PARS(I0)+PARS(I0+2)*X(3)
66 X(2)=PARS(I0+1)+PARS(I0+3)*X(3)
67 GO TO 20
68C
69 ENDIF
70C
71 IP3=ISH+2
72 IF(ISH.EQ.10) IP3=3
73 IF(ISH.EQ.4) IP3=1
74 X(3)=PARS(IP3)
75 IF(IP.LE.4) X(3)=-X(3)
76 IP2=3
77 IF(ISH.GT.2.AND.X(3).GT.0.0) IP2=4
78 IF(ISH.EQ.1.OR.ISH.EQ.10) IP2=2
79 IF(ISH.EQ.4) IP2=4
80 IF(ISH.EQ.4.AND.X(3).GT.0.0) IP2=8
81 X(2)=PARS(IP2)
82 IF(MOD(IP+3,4).LT.2) X(2)=-X(2)
83 IP1=1
84 IF(ISH.NE.1.AND.ISH.NE.10.AND.X(3).GT.0.0) IP1=2
85 IF(ISH.EQ.4) IP1=5
86 IF(ISH.EQ.4.AND.X(3).GT.0.0) IP1=IP1+4
87 IF(ISH.EQ.4.AND.X(2).GT.0.0) IP1=IP1+1
88 X(1)=PARS(IP1)
89 IF(MOD(IP,2).EQ.1) X(1)=-X(1)
90C
91 IF(ISH.NE.10) GO TO 10
92 X(1)=X(1)+X(2)*PARS(4)+X(3)*PARS(5)
93 X(2)=X(2)+X(3)*PARS(6)
94 10 CONTINUE
95C
96 IF(ISH.NE.4) GO TO 20
97 IP4=7
98 IF(X(3).GT.0.0) IP4=11
99 X(1)=X(1)+X(2)*PARS(IP4)+X(3)*PARS(2)
100 X(2)=X(2)+X(3)*PARS(3)
101 20 CONTINUE
102C
103C ROTATE.
104C
105 JROT=LQ(JROTM-IROT)
106 XT(1)=X(1)
107 XT(2)=X(2)
108 XT(3)=X(3)
109 IF(IROT.NE.0) CALL GINROT(X,Q(JROT+1),XT)
110C
111C NOW COMPUTE RMIN = PROJECTED R ON DX AND RMAX = R
112C AND UPDATE LIMITS IF NECESSARY.
113C
114 R2=(XT(1)+DX(1))**2+(XT(2)+DX(2))**2
115 IF(IAXIS.EQ.5) R2=R2+(XT(3)+DX(3))**2
116 R=SQRT(R2)
117 IF(R.GT.CH) CH=R
118C
119 IF(CL.LE.0.0) GO TO 30
120C
121 XPT=DX(1)*XT(1)+DX(2)*XT(2)
122 IF(IAXIS.EQ.5) XPT=XPT+DX(3)*XT(3)
123 IF(DXS.LE.1.0E-05) GO TO 30
124 RMN=DXS+XPT/DXS
125 IF(RMN.LT.CL) CL=RMN
126C
127 30 CONTINUE
128C
129 IF(CL.LE.0.0) CL=0.0
130C
131 IERR=0
132 GO TO 999
133C
134 40 CONTINUE
135 IF(ISH.GT.8.AND.ISH.NE.NSCTUB.AND.ISH.NE.13.AND.ISH.NE.14)GO TO 80
136C
137C TUBES AND CONES.
138C
139 IP3=3
140 IF(ISH.GT.6.AND.ISH.NE.NSCTUB.AND.ISH.NE.13.AND.ISH.NE.14) IP3=1
141 DZ=PARS(IP3)
142 R=PARS(2)
143 IF(ISH.EQ.NSCTUB) THEN
144 S1 = (1.0-PARS(8))*(1.0+PARS(8))
145 IF( S1 .GT. 0.0) S1 = SQRT(S1)
146 S2 = (1.0-PARS(11))*(1.0+PARS(11))
147 IF( S2 .GT. 0.0) S2 = SQRT(S2)
148 IF( S2 .GT. S1 ) S1 = S2
149 DZ = DZ+R*S1
150 ENDIF
151**
152 IF(ISH.EQ.13) THEN
153**
154** APPROXIME TO A CYLINDER WHIT RADIUS
155** EQUAL TO THE ELLIPSE MAJOR AXIS
156**
157 RMN=0.0
158 IF(PARS(1).GT.R) R=PARS(1)
159 GOTO 50
160 ENDIF
161 RMN=PARS(1)
162*
163 IF(ISH.EQ.14) THEN
164 R = SQRT(PARS(2)**2+(PARS(3)*TAN(PARS(4)*DEGRAD))**2)
165 GO TO 50
166 ENDIF
167C
168 IF(ISH.LE.6.OR.ISH.EQ.NSCTUB) GO TO 50
169C
170 R=PARS(3)
171 IF(PARS(5).GT.R) R=PARS(5)
172 RMN=PARS(2)
173 IF(PARS(4).LT.RMN) RMN=PARS(4)
174C
175 50 CONTINUE
176C
177C ROTATE THE LOCAL Z AXIS.
178C
179 X(1)=0.0
180 X(2)=0.0
181 X(3)=1.0
182 JROT=LQ(JROTM-IROT)
183 XT(1)=X(1)
184 XT(2)=X(2)
185 XT(3)=X(3)
186 IF(IROT.NE.0) CALL GINROT(X,Q(JROT+1),XT)
187C
188C COMPUTE RMIN AND RMAX ASSUMING COMPLETE TUBE HALF
189C LENGTH DZ AND RADIUS R.
190C
191 ST2=1.0
192 IF(IAXIS.EQ.4) ST2=(1+XT(3))*(1-XT(3))
193 DR=SQRT(DZ*DZ*ST2+R*R)
194 CL=DXS-DR
195 IF(CL.LT.0.0) CL=0.0
196 CH=DXS+DR
197 IF(IROT.EQ.0.AND.DXS.LT.1.0E-05) CL=RMN
198 IERR=0
199C
200 GO TO 999
201C
202 80 CONTINUE
203 IF(ISH.GT.9) GO TO 999
204C
205C SPHERE.
206C
207 CL=DXS-PARS(2)
208 IF(CL.LT.0.0) CL=0.0
209 CH=DXS+PARS(2)
210 IF(IAXIS.EQ.5.AND.DXS.LT.1.0E-05) CL=PARS(1)
211 IERR=0
212C
213 999 CONTINUE
214 END