Added protection. In case IROT=0 the address Q(LQ(JROTM-IROT)) should not
[u/mrichter/AliRoot.git] / GEANT321 / ggeom / gflcar.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1999/05/18 15:55:17  fca
6 * AliRoot sources
7 *
8 * Revision 1.1.1.1  1995/10/24 10:20:48  cernlib
9 * Geant
10 *
11 *
12 #include "geant321/pilot.h"
13 *CMZ :  3.21/02 29/03/94  15.41.28  by  S.Giani
14 *-- Author :
15       SUBROUTINE GFLCAR(IAXIS,ISH,IROT,PARS,CL,CH,IERR)
16 C.
17 C.    *****************************************************************
18 C.    *                                                               *
19 C.    *    ROUTINE TO FIND THE LIMITS ALONG AXIS IAXIS IN CARTESIAN   *
20 C.    *    COORDINATES FOR VOLUME OF SHAPE ISH ROTATED BY THE         *
21 C.    *    ROTATION MATRIX IROT. THE SHAPE HAS NPAR PARAMETERS IN     *
22 C.    *    THE ARRAY PARS. THE LOWER LIMIT IS RETURNED IN CL, THE     *
23 C.    *    HIGHER IN CH. IF THE CALCULATION CANNOT BE MADE IERR IS    *
24 C.    *    SET TO 1 OTHERWISE IT IS SET TO 0.                         *
25 C.    *                                                               *
26 C.    *    ==>Called by : GFCLIM                                      *
27 C.    *         Author  A.McPherson  *********                        *
28 C.    *                                                               *
29 C.    *****************************************************************
30 C.
31 #include "geant321/gcbank.inc"
32 #include "geant321/gconsp.inc"
33 #include "geant321/gcshno.inc"
34       DIMENSION PARS(11),X(3),XT(3)
35 C.
36 C.          ---------------------------------------------------
37 C.
38       IERR=1
39       IF (ISH.GT.4.AND.ISH.NE.10.AND.ISH.NE.28) GO TO 40
40 C
41 C           CUBOIDS, TRAPEZOIDS, PARALLELEPIPEDS.
42 C
43 C
44       IERR=0
45       CL=0
46       CH=0
47 C
48       DO 30 IP=1,8
49 C
50 C           THIS IS A LOOP OVER THE 8 CORNERS.
51 C           FIRST FIND THE LOCAL COORDINATES.
52 C
53       IF(ISH.EQ.28) THEN
54 C
55 C            General twisted trapezoid.
56 C
57          IL=(IP+1)/2
58          I0=IL*4+11
59          IS=(IP-IL*2)*2+1
60          X(3)=PARS(1)*IS
61          X(1)=PARS(I0)+PARS(I0+2)*X(3)
62          X(2)=PARS(I0+1)+PARS(I0+3)*X(3)
63          GO TO 20
64 C
65       ENDIF
66 C
67       IP3=ISH+2
68       IF(ISH.EQ.10) IP3=3
69       IF(ISH.EQ.4) IP3=1
70       X(3)=PARS(IP3)
71       IF(IP.LE.4) X(3)=-X(3)
72       IP2=3
73       IF(ISH.GT.2.AND.X(3).GT.0.0) IP2=4
74       IF(ISH.EQ.1.OR.ISH.EQ.10) IP2=2
75       IF(ISH.EQ.4) IP2=4
76       IF(ISH.EQ.4.AND.X(3).GT.0.0) IP2=8
77       X(2)=PARS(IP2)
78       IF(MOD(IP+3,4).LT.2) X(2)=-X(2)
79       IP1=1
80       IF(ISH.NE.1.AND.ISH.NE.10.AND.X(3).GT.0.0) IP1=2
81       IF(ISH.EQ.4) IP1=5
82       IF(ISH.EQ.4.AND.X(3).GT.0.0) IP1=IP1+4
83       IF(ISH.EQ.4.AND.X(2).GT.0.0) IP1=IP1+1
84       X(1)=PARS(IP1)
85       IF(MOD(IP,2).EQ.1) X(1)=-X(1)
86 C
87       IF(ISH.NE.10) GO TO 10
88       X(1)=X(1)+X(2)*PARS(4)+X(3)*PARS(5)
89       X(2)=X(2)+X(3)*PARS(6)
90    10 CONTINUE
91 C
92       IF(ISH.NE.4) GO TO 20
93       IP4=7
94       IF(X(3).GT.0.0) IP4=11
95       X(1)=X(1)+X(2)*PARS(IP4)+X(3)*PARS(2)
96       X(2)=X(2)+X(3)*PARS(3)
97    20 CONTINUE
98 C
99 C          ROTATE.
100 C
101       IF(IROT.NE.0) THEN
102          JROT=LQ(JROTM-IROT)
103          CALL GINROT(X,Q(JROT+1),XT)
104       ELSE
105          XT(1)=X(1)
106          XT(2)=X(2)
107          XT(3)=X(3)
108       ENDIF
109 C
110 C          UPDATE LIMITS IF NECESSARY.
111 C
112       IF(XT(IAXIS).LT.CL) CL=XT(IAXIS)
113       IF(XT(IAXIS).GT.CH) CH=XT(IAXIS)
114 C
115    30 CONTINUE
116 C
117       GO TO 999
118 C
119    40 CONTINUE
120       IF(ISH.EQ.9) GO TO 90
121 C
122 C              TUBES , CONES, POLYGONS, POLYCONES.
123 C              AND CUT TUBES.
124 C
125       X(1)=0.0
126       X(2)=0.0
127       X(3)=1.0
128       IF(IROT.NE.0) THEN
129          JROT=LQ(JROTM-IROT)
130          CALL GINROT(X,Q(JROT+1),XT)
131       ELSE
132          XT(1)=X(1)
133          XT(2)=X(2)
134          XT(3)=X(3)
135       ENDIF
136 C
137 C          XT IS Z AXIS ROTATED.
138 C
139       IF(ABS(XT(IAXIS)).LT.0.99) GO TO 50
140       IF(ISH.EQ.11)GO TO 45
141       IF(ISH.EQ.12)GO TO 46
142 C
143 C           PARALLEL.
144 C
145       IP=3
146       IF(ISH.GT.6.AND.ISH.NE.NSCTUB.AND.ISH.NE.13.AND.ISH.NE.14) IP=1
147       CL=-PARS(IP)
148       CH=PARS(IP)
149       IERR=0
150 C
151       GO TO 999
152   45  NZLAST=PARS(4)
153       IZLAST=2+3*NZLAST
154       CL=PARS(5)
155       GO TO 49
156 C
157   46  NZLAST=PARS(3)
158       IZLAST=1+3*NZLAST
159       CL=PARS(4)
160 C
161   49  CH=PARS(IZLAST)
162       IF ( ABS(XT(IAXIS)-X(IAXIS)) .GT.1.) THEN
163          TEMP = CL
164          CL = -CH
165          CH = -TEMP
166       ENDIF
167       IERR=0
168       GO TO 999
169 C
170    50 CONTINUE
171 **
172       IF(ISH.EQ.13) THEN
173          CL=-PARS(IAXIS)
174          CH=PARS(IAXIS)
175          IERR=0
176          GOTO 999
177       ENDIF
178 **
179       IF(ISH.EQ.14) THEN
180 C     for hyperboloid, use escribed cylinder
181          CH = SQRT(PARS(2)**2+(PARS(3)*TAN(PARS(4)*DEGRAD))**2)
182          CL = -CH
183          IERR=0
184          GOTO 999
185       ENDIF
186 **
187       IF(ISH.GT.10.AND.ISH.NE.NSCTUB)GO TO 999
188       IF(ABS(XT(IAXIS)).GT.0.01) GO TO 70
189 C
190 C         Z AXIS PERPENDICULAR TO IAXIS. ASSUME COMPLETE TUBE OR
191 C         CONE (I.E. IGNORE PHI SEGMENTATION).
192 C
193       IF(ISH.GT.6.AND.ISH.NE.NSCTUB) GO TO 60
194 C
195       CL=-PARS(2)
196       CH=PARS(2)
197       IERR=0
198 C
199       GO TO 999
200 C
201    60 CONTINUE
202 C
203       RM=PARS(3)
204       IF(PARS(5).GT.PARS(3)) RM=PARS(5)
205 C
206       CL=-RM
207       CH=RM
208       IERR=0
209 C
210       GO TO 999
211 C
212    70 CONTINUE
213 C
214 C           ARBITRARY ROTATION.
215 C
216       DZ=PARS(3)
217       RM=PARS(2)
218       IF(ISH.EQ.13) THEN
219 **
220 **       approxime to a cylinder whit radius
221 **       equal to the ellipse major axis
222 **
223          IF(PARS(1).GT.RM) RM=PARS(1)
224          GOTO 80
225       ENDIF
226 **
227       IF(ISH.EQ.14) THEN
228         RM = SQRT(PARS(2)**2+(PARS(3)*TAN(PARS(4)*DEGRAD))**2)
229         GO TO 80
230       ENDIF
231 *
232       IF(ISH.EQ.NSCTUB) THEN
233         S1 = (1.0-PARS(8))*(1.0+PARS(8))
234         IF( S1 .GT. 0.0) S1 = SQRT(S1)
235         S2 = (1.0-PARS(11))*(1.0+PARS(11))
236         IF( S2 .GT. 0.0) S2 = SQRT(S2)
237         IF( S2 .GT. S1 ) S1 = S2
238         DZ = DZ+RM*S1
239       ENDIF
240       IF(ISH.LE.6) GO TO 80
241 C
242       DZ=PARS(1)
243       RM=PARS(3)
244       IF(PARS(5).GT.RM) RM=PARS(5)
245 C
246    80 CONTINUE
247 C
248       COST=ABS(XT(IAXIS))
249       SINT=(1+COST)*(1-COST)
250       IF(SINT.GT.0.0) SINT=SQRT(SINT)
251 C
252       CH=COST*DZ+SINT*RM
253       CL=-CH
254       IERR=0
255 C
256       GO TO 999
257    90 CONTINUE
258 C
259 C           SPHERE - ASSUME COMPLETE SPHERE, TAKE OUTER RADIUS.
260 C
261       IERR=0
262       CL=-PARS(2)
263       CH=PARS(2)
264 C
265   999 CONTINUE
266       END