[u/mrichter/AliRoot.git] / GEANT321 / gphys / glando.F

*
* $Id$
*
* $Log$
* Revision 1.1.1.1  1995/10/24 10:21:25  cernlib
* Geant
*
*
#include "geant321/pilot.h"
*CMZ :  3.21/02 29/03/94  15.41.21  by  S.Giani
*-- Author :
      SUBROUTINE GLANDO(IMODE,STEP,Z,A,RHO,P,E,XMASS,DE,IFLAG)
C.
C.    ******************************************************************
C.    *                                                                *
C.    *                                                                *
C.    *  GLANDO                                                        *
C.    *  ------                                                        *
C.    *                                                                *
C.    *  Energy straggling using Gaussian, Landau & Vavilov theories.  *
C.    *                                                                *
C.    *  Input                                                         *
C.    *  -----                                                         *
C.    *  IMODE  =  3   Landau(RANLAN)  sampling                        *
C.    *         =  4   Landau(GENLAN)  sampling                        *
C.    *         =  5   Vavilov(DINVAV) sampling                        *
C.    *         =  6   Gaussian sampling                               *
C.    *         =  2   Automatic selection of relevant distribution    *
C.    *                                                                *
C.    *  STEP   =  current step-length (cm)                            *
C.    *                                                                *
C.    *  Output                                                        *
C.    *  ------                                                        *
C.    *  DE     =  DE/DX - <DE/DX>     (GeV)                           *
C.    *                                                                *
C.    *  IFLAG  =  3   Landau(RANLAN) sampling used                    *
C.    *         =  4   Landau(GENLAN) sampling used                    *
C.    *         =  5   Vavilov(GVAVIV) sampling used                   *
C.    *         =  6   Gaussian sampling used                          *
C.    *                                                                *
C.    *  Warning                                                       *
C.    *  -------                                                       *
C.    *  Only Landau sampling should be used since this has been well  *
C.    *  tested whereas both Vavilov and Gaussian sampling are being   *
C.    *  developed.                                                    *
C.    *                                                                *
C.    *  Author      : G.N. Patrick                                    *
C.    *  Date        : 03.05.1985                                      *
C.    *  Last update : 09.09.1985                                      *
C.    *                                                                *
C.    ******************************************************************
C.
      PARAMETER (EULER=0.577215)
      PARAMETER (P1=.60715,P2=.67794,P3=.52382E-1,P4=.94753,
     +           P5=.74442,P6=1.1934)
#include "geant321/gconsp.inc"
      DIMENSION RNDM(2)
      FLAND(X) = P1+P6*X+(P2+P3*X)*EXP(P4*X+P5)
C.
C.    ------------------------------------------------------------------
C.
      IF(STEP.LT.1.E-7)THEN
         DE=0.
         IFLAG=0
         RETURN
      ENDIF
C
C     Calculate xi factor (KeV).
C
      BETA   = P/E
      GAMMA  = E/XMASS
      XI     = (153.5*Z*STEP*RHO)/(A*BETA*BETA)
C
C     Maximum energy transfer to atomic electron (KeV).
C
      ETA    = BETA*GAMMA
      ETASQ  = ETA*ETA
      RATIO  = EMASS/XMASS
      F1     = 2.*EMASS*ETASQ
      F2     = 1.+2.*RATIO*GAMMA+RATIO*RATIO
      EMAX   = F1*1.E+6/F2
C
C     Calculate Kappa significance ratio.
C
      CAPPA  = XI/EMAX
C
C     Choose correct function if IMODE set to automatic selection.
C
      IMODEI = IMODE
      IF (IMODEI.EQ.2)    THEN
         IF (CAPPA.LT.0.01)                              IMODEI = 4
*****    IF (CAPPA.GE.0.01.AND.CAPPA.LE.10.)             IMODEI = 3
         IF (CAPPA.GE.0.01.AND.CAPPA.LE.10.)             IMODEI = 5
         IF (CAPPA.GT.10.)                               IMODEI = 6
      ENDIF
C
C     +---------------------------------------------------------------+
C     I Sample lambda variable from Kolbig/Schorr Landau distribution I
C     +---------------------------------------------------------------+
C
   10 CONTINUE
      IF (IMODEI.EQ.3) THEN
         XMEAN = -BETA**2-LOG(XI/EMAX)+EULER-1.
         XLAMX = FLAND(XMEAN)
   20    CALL GRNDM(RNDM,1)
         IF( RNDM(1) .GT. 0.980 ) GO TO 20
         XLAMB = GLANDR(RNDM(1))
         IF(XLAMB.GT.XLAMX) GO TO 20
         IFLAG  = 3
C
C            +----------------------------------------------+
C            I  Sample lambda variable from James & Hancock I
C            I  Landau distribution                         I
C            +----------------------------------------------+
C
 
      ELSEIF (IMODEI.EQ.4) THEN
         XMEAN = -BETA**2-LOG(XI/EMAX)+EULER-1.
         XLAMX = FLAND(XMEAN)
  25     CALL GLANDG(XLAMB)
         IF(XLAMB.GT.XLAMX) GO TO 25
         IFLAG = 4
C
C            +---------------------------------------------------------+
C            I Sample lambda variable (Landau not Vavilov) from        I
C            I Rotondi&Montagna&Kolbig Vavilov distribution            I
C            +---------------------------------------------------------+
C
      ELSEIF (IMODEI.EQ.5) THEN
C
C            Keep within computation bounds.
C
         RKA = CAPPA
         BE2 = BETA*BETA
         CALL GRNDM(RNDM,1)
         XLAMB = GVAVIV(RKA,BE2,RNDM(1))
         IFLAG = 5
C
C     +-----------------------------------+
C     I Sample from Gaussian distribution I
C     +-----------------------------------+
C
      ELSEIF (IMODEI.EQ.6) THEN
         SIGMA  = XI*EMAX*(1.-(BETA*BETA/2.))
         SIGMA  = SQRT(SIGMA)
   30    CALL GRNDM(RNDM,2)
         IF(RNDM(1).LE.0.)GO TO 30
         F1     = -2.*LOG(RNDM(1))
         DEKEV  = SIGMA*SQRT(F1)*COS(2.*PI*RNDM(2))
         IFLAG  = 6
         GOTO 40
      ENDIF
C
C     Calculate DE/DX - <DE/DX>
C
      DEKEV  = XI*(XLAMB+BETA*BETA+LOG(XI/EMAX)-EULER+1.)
C
   40 DE=DEKEV*1.E-6
      END
Commit	Line	Data
fe4da5cc	1	*
	2	* $Id$
	3	*
	4	* $Log$
	5	* Revision 1.1.1.1 1995/10/24 10:21:25 cernlib
	6	* Geant
	7	*
	8	*
	9	#include "geant321/pilot.h"
	10	*CMZ : 3.21/02 29/03/94 15.41.21 by S.Giani
	11	*-- Author :
	12	SUBROUTINE GLANDO(IMODE,STEP,Z,A,RHO,P,E,XMASS,DE,IFLAG)
	13	C.
	14	C. ******************************************************************
	15	C. * *
	16	C. * *
	17	C. * GLANDO *
	18	C. * ------ *
	19	C. * *
	20	C. * Energy straggling using Gaussian, Landau & Vavilov theories. *
	21	C. * *
	22	C. * Input *
	23	C. * ----- *
	24	C. * IMODE = 3 Landau(RANLAN) sampling *
	25	C. * = 4 Landau(GENLAN) sampling *
	26	C. * = 5 Vavilov(DINVAV) sampling *
	27	C. * = 6 Gaussian sampling *
	28	C. * = 2 Automatic selection of relevant distribution *
	29	C. * *
	30	C. * STEP = current step-length (cm) *
	31	C. * *
	32	C. * Output *
	33	C. * ------ *
	34	C. * DE = DE/DX - <DE/DX> (GeV) *
	35	C. * *
	36	C. * IFLAG = 3 Landau(RANLAN) sampling used *
	37	C. * = 4 Landau(GENLAN) sampling used *
	38	C. * = 5 Vavilov(GVAVIV) sampling used *
	39	C. * = 6 Gaussian sampling used *
	40	C. * *
	41	C. * Warning *
	42	C. * ------- *
	43	C. * Only Landau sampling should be used since this has been well *
	44	C. * tested whereas both Vavilov and Gaussian sampling are being *
	45	C. * developed. *
	46	C. * *
	47	C. * Author : G.N. Patrick *
	48	C. * Date : 03.05.1985 *
	49	C. * Last update : 09.09.1985 *
	50	C. * *
	51	C. ******************************************************************
	52	C.
	53	PARAMETER (EULER=0.577215)
	54	PARAMETER (P1=.60715,P2=.67794,P3=.52382E-1,P4=.94753,
	55	+ P5=.74442,P6=1.1934)
	56	#include "geant321/gconsp.inc"
	57	DIMENSION RNDM(2)
	58	FLAND(X) = P1+P6X+(P2+P3X)EXP(P4X+P5)
	59	C.
	60	C. ------------------------------------------------------------------
	61	C.
	62	IF(STEP.LT.1.E-7)THEN
	63	DE=0.
	64	IFLAG=0
65	RETURN
66	ENDIF
67	C
68	C Calculate xi factor (KeV).
69	C
70	BETA = P/E
71	GAMMA = E/XMASS
72	XI = (153.5ZSTEPRHO)/(ABETA*BETA)
73	C
74	C Maximum energy transfer to atomic electron (KeV).
75	C
76	ETA = BETA*GAMMA
77	ETASQ = ETA*ETA
78	RATIO = EMASS/XMASS
79	F1 = 2.EMASSETASQ
80	F2 = 1.+2.RATIOGAMMA+RATIO*RATIO
81	EMAX = F1*1.E+6/F2
82	C
83	C Calculate Kappa significance ratio.
84	C
85	CAPPA = XI/EMAX
86	C
87	C Choose correct function if IMODE set to automatic selection.
88	C
89	IMODEI = IMODE
90	IF (IMODEI.EQ.2) THEN
91	IF (CAPPA.LT.0.01) IMODEI = 4
92	***** IF (CAPPA.GE.0.01.AND.CAPPA.LE.10.) IMODEI = 3
93	IF (CAPPA.GE.0.01.AND.CAPPA.LE.10.) IMODEI = 5
94	IF (CAPPA.GT.10.) IMODEI = 6
95	ENDIF
96	C
97	C +---------------------------------------------------------------+
98	C I Sample lambda variable from Kolbig/Schorr Landau distribution I
99	C +---------------------------------------------------------------+
100	C
101	10 CONTINUE
102	IF (IMODEI.EQ.3) THEN
103	XMEAN = -BETA**2-LOG(XI/EMAX)+EULER-1.
104	XLAMX = FLAND(XMEAN)
105	20 CALL GRNDM(RNDM,1)
106	IF( RNDM(1) .GT. 0.980 ) GO TO 20
107	XLAMB = GLANDR(RNDM(1))
108	IF(XLAMB.GT.XLAMX) GO TO 20
109	IFLAG = 3
110	C
111	C +----------------------------------------------+
112	C I Sample lambda variable from James & Hancock I
113	C I Landau distribution I
114	C +----------------------------------------------+
115	C
116
117	ELSEIF (IMODEI.EQ.4) THEN
118	XMEAN = -BETA**2-LOG(XI/EMAX)+EULER-1.
119	XLAMX = FLAND(XMEAN)
120	25 CALL GLANDG(XLAMB)
121	IF(XLAMB.GT.XLAMX) GO TO 25
122	IFLAG = 4
123	C
124	C +---------------------------------------------------------+
125	C I Sample lambda variable (Landau not Vavilov) from I
126	C I Rotondi&Montagna&Kolbig Vavilov distribution I
127	C +---------------------------------------------------------+
128	C
129	ELSEIF (IMODEI.EQ.5) THEN
130	C
131	C Keep within computation bounds.
132	C
133	RKA = CAPPA
134	BE2 = BETA*BETA
135	CALL GRNDM(RNDM,1)
136	XLAMB = GVAVIV(RKA,BE2,RNDM(1))
137	IFLAG = 5
138	C
139	C +-----------------------------------+
140	C I Sample from Gaussian distribution I
141	C +-----------------------------------+
142	C
143	ELSEIF (IMODEI.EQ.6) THEN
144	SIGMA = XIEMAX(1.-(BETA*BETA/2.))
145	SIGMA = SQRT(SIGMA)
146	30 CALL GRNDM(RNDM,2)
147	IF(RNDM(1).LE.0.)GO TO 30
148	F1 = -2.*LOG(RNDM(1))
149	DEKEV = SIGMASQRT(F1)COS(2.PIRNDM(2))
150	IFLAG = 6
151	GOTO 40
152	ENDIF
153	C
154	C Calculate DE/DX - <DE/DX>
155	C
156	DEKEV = XI(XLAMB+BETABETA+LOG(XI/EMAX)-EULER+1.)
157	C
158	40 DE=DEKEV*1.E-6
159	END