GEANT321/ggeom/gvdcar.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:56  cernlib
   9 * Geant
  10 *
  11 *
  12 #include "geant321/pilot.h"
  13 *CMZ :  3.21/03 10/10/94  20.01.58  by  S.Giani
  14 *-- Author :
  15       SUBROUTINE GVDCAR(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 : GVDLIM                                      *
  27 C.    *         Author  S.Giani  ********                             *
  28 C.    *                                                               *
  29 C.    *****************************************************************
  30 C.
  31 #include "geant321/gcbank.inc"
  32 #include "geant321/gconsp.inc"
  33 #include "geant321/gcshno.inc"
  34       DIMENSION PARS(100),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       JROT=LQ(JROTM-IROT)
 102       XT(1)=X(1)
 103       XT(2)=X(2)
 104       XT(3)=X(3)
 105       IF(IROT.NE.0) CALL GINROT(X,Q(JROT+1),XT)
 106 C
 107 C          UPDATE LIMITS IF NECESSARY.
 108 C
 109       IF(XT(IAXIS).LT.CL) CL=XT(IAXIS)
 110       IF(XT(IAXIS).GT.CH) CH=XT(IAXIS)
 111 C
 112    30 CONTINUE
 113 C
 114       GO TO 999
 115 C
 116    40 CONTINUE
 117       IF(ISH.EQ.9) GO TO 90
 118 C
 119 C              TUBES , CONES, POLYGONS, POLYCONES.
 120 C              AND CUT TUBES.
 121 C
 122       MYFLAG=0
 123       IF((ISH.EQ.11.OR.ISH.EQ.12).AND.(IAXIS.LT.3))THEN
 124         MYFLAG=1
 125       ENDIF
 126       X(1)=0.0
 127       X(2)=0.0
 128       X(3)=1.0
 129       JROT=LQ(JROTM-IROT)
 130       XT(1)=X(1)
 131       XT(2)=X(2)
 132       XT(3)=X(3)
 133       IF(IROT.NE.0) CALL GINROT(X,Q(JROT+1),XT)
 134 C
 135 C          XT IS Z AXIS ROTATED.
 136 C
 137       IF(MYFLAG.EQ.0)THEN
 138        IF(ABS(XT(IAXIS)).LT.0.99) GO TO 50
 139       ELSE
 140        IF(ABS(XT(3)).LT.0.99) GO TO 50
 141       ENDIF
 142       IF(ISH.EQ.11)GO TO 45
 143       IF(ISH.EQ.12)GO TO 46
 144 C
 145 C           PARALLEL.
 146 C
 147       IP=3
 148       IF(ISH.GT.6.AND.ISH.NE.NSCTUB.AND.ISH.NE.13.AND.ISH.NE.14) IP=1
 149       CL=-PARS(IP)
 150       CH=PARS(IP)
 151       IERR=0
 152 C
 153       GO TO 999
 154   45  IF(MYFLAG.EQ.0)THEN
 155        NZLAST=PARS(4)
 156        IZLAST=2+3*NZLAST
 157        CL=PARS(5)
 158        GO TO 49
 159       ELSEIF(MYFLAG.EQ.1)THEN
 160        NZLAST=PARS(4)
 161        IZLAST=2+3*NZLAST
 162        TMPRAD=0.
 163        DO 145 I=7,IZLAST+2,3
 164          IF(PARS(I).GT.TMPRAD)TMPRAD=PARS(I)
 165  145   CONTINUE
 166        PHIMIN=PARS(1)
 167        PHIMAX=PHIMIN+PARS(2)
 168        AANG=ABS(PHIMAX-PHIMIN)
 169        NANG=PARS(3)
 170        AATMAX=NANG*360./AANG
 171        LATMAX=AATMAX
 172        ALA=AATMAX-LATMAX
 173        IF(ALA.GT..5)LATMAX=LATMAX+1
 174        AFINV=1./COS(PI/LATMAX)
 175        FINV=ABS(AFINV)
 176        R=TMPRAD*FINV
 177        CL=-R
 178        CH= R
 179        IERR=0
 180        GOTO 999
 181       ENDIF
 182 C
 183   46  IF(MYFLAG.EQ.0)THEN
 184        NZLAST=PARS(3)
 185        IZLAST=1+3*NZLAST
 186        CL=PARS(4)
 187       ELSEIF(MYFLAG.EQ.1)THEN
 188        NZLAST=PARS(3)
 189        IZLAST=1+3*NZLAST
 190        TMPRAD=0.
 191        DO 146 I=6,IZLAST+2,3
 192          IF(PARS(I).GT.TMPRAD)TMPRAD=PARS(I)
 193  146   CONTINUE
 194        CL=-TMPRAD
 195        CH= TMPRAD
 196        IERR=0
 197        GOTO 999
 198       ENDIF
 199 C
 200   49  CH=PARS(IZLAST)
 201       IF ( ABS(XT(IAXIS)-X(IAXIS)) .GT.1.) THEN
 202          TEMP = CL
 203          CL = -CH
 204          CH = -TEMP
 205       ENDIF
 206       IERR=0
 207       GO TO 999
 208 C
 209    50 CONTINUE
 210 **
 211       IF(ISH.EQ.13) THEN
 212          CL=-PARS(IAXIS)
 213          CH=PARS(IAXIS)
 214          IERR=0
 215          GOTO 999
 216       ENDIF
 217 **
 218       IF(ISH.EQ.14) THEN
 219 C     for hyperboloid, use escribed cylinder
 220          CH = SQRT(PARS(2)**2+(PARS(3)*TAN(PARS(4)*DEGRAD))**2)
 221          CL = -CH
 222          IERR=0
 223          GOTO 999
 224       ENDIF
 225 **
 226       IF(ISH.GT.10.AND.ISH.NE.NSCTUB)GO TO 999
 227       IF(ABS(XT(IAXIS)).GT.0.01) GO TO 70
 228 C
 229 C         Z AXIS PERPENDICULAR TO IAXIS. ASSUME COMPLETE TUBE OR
 230 C         CONE (I.E. IGNORE PHI SEGMENTATION).
 231 C
 232       IF(ISH.GT.6.AND.ISH.NE.NSCTUB) GO TO 60
 233 C
 234       CL=-PARS(2)
 235       CH=PARS(2)
 236       IERR=0
 237       IF(ISH.EQ.6)THEN
 238         RMIN=PARS(1)
 239         RMAX=PARS(2)
 240         IF(IROT.NE.0)THEN
 241          IF(Q(JROT+15).EQ.0.)THEN
 242            PHI1=(PARS(4)+Q(JROT+12))*DEGRAD
 243            PHI2=(PARS(5)+Q(JROT+12))*DEGRAD
 244          ELSEIF(Q(JROT+15).EQ.180.)THEN
 245            PHI1=(PARS(4)+Q(JROT+12)-(PARS(5)-PARS(4)))*DEGRAD
 246            PHI2=(PARS(5)+Q(JROT+12)-(PARS(5)-PARS(4)))*DEGRAD
 247          ELSE
 248            GOTO 999
 249          ENDIF
 250         ELSE
 251          PHI1=PARS(4)*DEGRAD
 252          PHI2=PARS(5)*DEGRAD
 253         ENDIF
 254         IF(IAXIS.EQ.1)THEN
 255          IF(PHI1.GE.0..AND.PHI2.LE.PI)THEN
 256           XMIN1=RMIN*COS(PHI2)
 257           XMIN2=RMAX*COS(PHI2)
 258           CL=MIN(XMIN1,XMIN2)
 259           XMAX1=RMIN*COS(PHI1)
 260           XMAX2=RMAX*COS(PHI1)
 261           CH=MAX(XMAX1,XMAX2)
 262          ELSEIF(PHI1.GE.PI.AND.PHI2.LE.TWOPI.OR.
 263      +          PHI1.GE.-PI.AND.PHI2.LE.0.)THEN
 264           XMIN1=RMIN*COS(PHI1)
 265           XMIN2=RMAX*COS(PHI1)
 266           CL=MIN(XMIN1,XMIN2)
 267           XMAX1=RMIN*COS(PHI2)
 268           XMAX2=RMAX*COS(PHI2)
 269           CH=MAX(XMAX1,XMAX2)
 270          ELSEIF(PHI1.LT.0..AND.PHI2.GT.0..AND.
 271      +         (PHI2-PHI1).LE.PI)THEN
 272           XMIN1=RMIN*COS(PHI2)
 273           XMIN2=RMIN*COS(PHI1)
 274           CL1=MIN(XMIN1,XMIN2)
 275           XMIN3=RMAX*COS(PHI2)
 276           XMIN4=RMAX*COS(PHI1)
 277           CL2=MIN(XMIN3,XMIN4)
 278           CL=MIN(CL1,CL2)
 279           CH=RMAX
 280          ELSEIF(PHI1.LT.PI.AND.PHI2.GT.PI.AND.
 281      +         (PHI2-PHI1).LE.PI)THEN
 282           CL=-RMAX
 283           XMAX1=RMIN*COS(PHI2)
 284           XMAX2=RMIN*COS(PHI1)
 285           CH1=MAX(XMAX1,XMAX2)
 286           XMAX3=RMAX*COS(PHI2)
 287           XMAX4=RMAX*COS(PHI1)
 288           CH2=MAX(XMAX3,XMAX4)
 289           CH=MAX(CH1,CH2)
 290          ENDIF
 291         ELSEIF(IAXIS.EQ.2)THEN
 292          IF(PHI1.GE.(-PI*.5).AND.PHI2.LE.(PI*.5))THEN
 293           YMIN1=RMIN*SIN(PHI1)
 294           YMIN2=RMAX*SIN(PHI1)
 295           CL=MIN(YMIN1,YMIN2)
 296           YMAX1=RMIN*SIN(PHI2)
 297           YMAX2=RMAX*SIN(PHI2)
 298           CH=MAX(YMAX1,YMAX2)
 299          ELSEIF(PHI1.GE.(PI*.5).AND.PHI2.LE.(PI*3*.5))THEN
 300           YMIN1=RMIN*SIN(PHI2)
 301           YMIN2=RMAX*SIN(PHI2)
 302           CL=MIN(YMIN1,YMIN2)
 303           YMAX1=RMIN*SIN(PHI1)
 304           YMAX2=RMAX*SIN(PHI1)
 305           CH=MAX(YMAX1,YMAX2)
 306          ELSEIF(PHI1.LT.(PI*.5).AND.PHI2.GT.(PI*.5).AND.
 307      +         (PHI2-PHI1).LE.PI)THEN
 308           YMIN1=RMIN*SIN(PHI2)
 309           YMIN2=RMIN*SIN(PHI1)
 310           CL1=MIN(YMIN1,YMIN2)
 311           YMIN3=RMAX*SIN(PHI2)
 312           YMIN4=RMAX*SIN(PHI1)
 313           CL2=MIN(YMIN3,YMIN4)
 314           CL=MIN(CL1,CL2)
 315           CH=RMAX
 316          ELSEIF(((PHI1.LT.(PI*3*.5).AND.PHI2.GT.(PI*3*.5)).OR.
 317      +          (PHI1.LT.-(PI*.5).AND.PHI2.GT.-(PI*.5)))
 318      +          .AND.(PHI2-PHI1).LE.PI)THEN
 319           CL=-RMAX
 320           YMAX1=RMIN*SIN(PHI2)
 321           YMAX2=RMIN*SIN(PHI1)
 322           CH1=MAX(YMAX1,YMAX2)
 323           YMAX3=RMAX*SIN(PHI2)
 324           YMAX4=RMAX*SIN(PHI1)
 325           CH2=MAX(YMAX3,YMAX4)
 326           CH=MAX(CH1,CH2)
 327          ENDIF
 328         ENDIF
 329       ENDIF
 330 C
 331       GO TO 999
 332 C
 333    60 CONTINUE
 334 C
 335       RM=PARS(3)
 336       IF(PARS(5).GT.PARS(3)) RM=PARS(5)
 337 C
 338       CL=-RM
 339       CH=RM
 340       IERR=0
 341 C
 342       GO TO 999
 343 C
 344    70 CONTINUE
 345 C
 346 C           ARBITRARY ROTATION.
 347 C
 348       DZ=PARS(3)
 349       RM=PARS(2)
 350       IF(ISH.EQ.13) THEN
 351 **
 352 **       approxime to a cylinder whit radius
 353 **       equal to the ellipse major axis
 354 **
 355          IF(PARS(1).GT.RM) RM=PARS(1)
 356          GOTO 80
 357       ENDIF
 358 **
 359       IF(ISH.EQ.14) THEN
 360         RM = SQRT(PARS(2)**2+(PARS(3)*TAN(PARS(4)*DEGRAD))**2)
 361         GO TO 80
 362       ENDIF
 363 *
 364       IF(ISH.EQ.NSCTUB) THEN
 365         S1 = (1.0-PARS(8))*(1.0+PARS(8))
 366         IF( S1 .GT. 0.0) S1 = SQRT(S1)
 367         S2 = (1.0-PARS(11))*(1.0+PARS(11))
 368         IF( S2 .GT. 0.0) S2 = SQRT(S2)
 369         IF( S2 .GT. S1 ) S1 = S2
 370         DZ = DZ+RM*S1
 371       ENDIF
 372       IF(ISH.LE.6) GO TO 80
 373 C
 374       DZ=PARS(1)
 375       RM=PARS(3)
 376       IF(PARS(5).GT.RM) RM=PARS(5)
 377 C
 378    80 CONTINUE
 379 C
 380       COST=ABS(XT(IAXIS))
 381       SINT=(1+COST)*(1-COST)
 382       IF(SINT.GT.0.0) SINT=SQRT(SINT)
 383 C
 384       CH=COST*DZ+SINT*RM
 385       CL=-CH
 386       IERR=0
 387 C
 388       GO TO 999
 389    90 CONTINUE
 390 C
 391 C           SPHERE - ASSUME COMPLETE SPHERE, TAKE OUTER RADIUS.
 392 C
 393       IERR=0
 394       CL=-PARS(2)
 395       CH=PARS(2)
 396 C
 397   999 CONTINUE
 398       END