This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / GEANT321 / ggeom / gintco.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1995/10/24 10:20:51  cernlib
6 * Geant
7 *
8 *
9 #include "geant321/pilot.h"
10 *CMZ :  3.21/02 29/03/94  15.41.29  by  S.Giani
11 *-- Author :
12       SUBROUTINE GINTCO (X, RLEFT, RRIGHT, DZ, TAU, TAUL)
13  
14 C        ********************************************
15 C        *  This subroutine finds the intersection  *
16 C        *  of a given ray (described by array X)   *
17 C        *  with a given cone (described by radii   *
18 C        *  RLEFT and RRIGHT and half-distance DZ). *
19 C        *  Output parameter is TAU, and inter-     *
20 C        *  section point is X = XP + TAU*XD, etc.  *
21 C        *                                          *
22 C        *  Called by GNOCON                        *
23 C        *  Programmed by:  Patrick Weidhaas        *
24 C        *          CERN,   March 1988              *
25 C        ********************************************
26  
27 #include "geant321/gconsp.inc"
28       DIMENSION X(6)
29 #if !defined(CERNLIB_SINGLE)
30       DOUBLE PRECISION XP,YP,ZP,XD,YD,ZD,S,T,U,V,W
31       DOUBLE PRECISION DISCR,SQDISC
32 #endif
33 C----------------------------------------------------
34  
35 C...... Point of origin of ray:
36  
37       XP = X(1)
38       YP = X(2)
39       ZP = X(3)
40  
41 C...... Direction cosines:
42  
43       XD = X(4)
44       YD = X(5)
45       ZD = X(6)
46  
47       TAU  = BIG
48       TAUL = BIG
49  
50       S = 0.5 * (RLEFT + RRIGHT)
51       T = (RLEFT - RRIGHT) / DZ
52  
53 C......  Cone equation is:    x**2 + y**2 - Az**2 + Bz + C = 0
54  
55       A = 0.25 * T*T
56       B = S * T
57       C = -S*S
58  
59 C......  To obtain "TAU", we must solve the quadratic equation
60 C......  Ut**2 + Vt + W = 0 .
61  
62       U = XD**2 + YD**2 - A*ZD**2
63       V = 2.0 * (XP*XD + YP*YD - A*ZP*ZD) + B*ZD
64       W = XP**2 + YP**2 - A*ZP**2 + B*ZP + C
65  
66       DISCR = V*V - 4.0*U*W
67       IF (DISCR .LE. 0.0) GO TO 999
68       IF(U.EQ.0.)GO TO 999
69       SQDISC = SQRT (DISCR)
70       TAU1 = (-V + SQDISC) / (2.0*U)
71       TAU2 = (-V - SQDISC) / (2.0*U)
72  
73  
74 C......  Set TAU to the smallest positive root;
75 C......  otherwise let TAU = BIG .
76 C
77 C......  If both roots are positive, set TAUL to
78 C......  the larger one: it may be needed in the
79 C......  case of a PHI-segmented cone.
80  
81       IF (TAU1 .LT. 0.0) THEN
82         IF (TAU2 .LT. 0.0) GO TO 999
83         TAU = TAU2
84       ELSE
85         TAU = TAU1
86         IF (TAU2 .GT. 0.0) THEN
87           TAUL = TAU2
88           IF (TAU2.LT.TAU1) THEN
89             TAU = TAU2
90             TAUL = TAU1
91           ENDIF
92         ENDIF
93       ENDIF
94  
95   999 CONTINUE
96       END
97