This commit was generated by cvs2svn to compensate for changes in r2,
[u/mrichter/AliRoot.git] / GEANT321 / gcons / gplmat.F
1 *
2 * $Id$
3 *
4 * $Log$
5 * Revision 1.1.1.1  1995/10/24 10:20:15  cernlib
6 * Geant
7 *
8 *
9 #include "geant321/pilot.h"
10 *CMZ :  3.21/04 10/01/95  17.36.36  by  S.Ravndal
11 *-- Author :
12       SUBROUTINE GPLMAT(IMATES,IPART,MECAN,KDIN,TKIN,IDM)
13 C.
14 C     ******************************************************************
15 C.    *                                                                *
16 C.    *       INTERPOLATE and PLOT  the DE/DX and Cross sections       *
17 C.    *       tabulated in JMATE banks corresponding to :              *
18 C.    *       material IMATE, particle IPART, mecanism name HMECAN,    *
19 C.    *       kinetic energies TKIN                                    *
20 C.    *                                                                *
21 C.    *      The MECANism name can be :                                *
22 C.    *      'HADF'   'INEF'   'ELAF'   'FISF'   'CAPF'                *
23 C.    *      'HADG'   'INEG'   'ELAG'   'FISG'   'CAPG'                *
24 C.    *      'LOSS'   'PHOT'   'ANNI'   'COMP'   'BREM'                *
25 C.    *      'PAIR'   'DRAY'   'PFIS'   'RAYL'   'HADG'                *
26 C.    *      'MUNU'   'RANG'   'STEP'                                  *
27 C.    *                                                                *
28 C.    *       For Hadronic particles it also computes the              *
29 C.    *       hadronic cross section from FLUKA ( '***F' ) or          *
30 C.    *       GHEISHA ( '***G' ) programs:                             *
31 C.    *       HADF or HADG -- total                                    *
32 C.    *       INEF or INEG -- inelastic                                *
33 C.    *       ELAF or ELAG -- elastic                                  *
34 C.    *       FISF or FISG -- fission (0.0 for FLUKA)                  *
35 C.    *       CAPF or CAPG -- neutron capture (0.0 for FLUKA)          *
36 C.    *                                                                *
37 C.    *             Input parameters                                   *
38 C.    *  IMATE   Geant material number                                 *
39 C.    *  IPART   Geant particle number                                 *
40 C.    *  MECAN  mechanism name of the bank to be fetched               *
41 C.    *  KDIM   dimension of the arrays TKIN , VALUE                   *
42 C.    *  TKIN   array of kinetic energy of incident particle (in Gev)  *
43 C.    *  IDM    convention for histogramming mode :                    *
44 C.    *         IDM.gt.0  fill , print , keep histogram(s)             *
45 C.    *         IDM.eq.0  fill , print , delete histogram(s)           *
46 C.    *         IDM.lt.0  fill , noprint , keep histogram(s)           *
47 C.    *           The histogram IDentificator will be :                *
48 C.    *             10000*imate + 100*ipart + imeca                    *
49 C.    *          where IMECA is the link number in stucture JMATE      *
50 C.    *          (see Geant3 writeup CONS 199)                         *
51 C.    *           for 'HADG'  imeca = 17                               *
52 C.    *                                                                *
53 C.    *    ==>Called by : <USER>                                       *
54 C.    *       Authors    R.Brun, M.Maire    *********                  *
55 C.    *                                                                *
56 C.    ******************************************************************
57 C.
58 #include "geant321/gcbank.inc"
59 #include "geant321/gcnum.inc"
60 #include "geant321/gconsp.inc"
61 #include "geant321/gcunit.inc"
62 #include "geant321/gcphys.inc"
63       PARAMETER (MMX= 201,NCOL= 5)
64       CHARACTER*(*) MECAN
65       CHARACTER*4 MECA , KU(NCOL)
66       CHARACTER   NAPART*16 , NAMATE*16 , CHTITL*68
67       DIMENSION   TKIN(KDIN),VALUE(MMX),SIGT(MMX),PCUT(5)
68       DIMENSION   KI(NCOL),EK(NCOL)
69       LOGICAL     LXBARN
70 *
71 #include "geant321/gcnmec.inc"
72 *
73 *     ------------------------------------------------------------------
74 *
75       IF (KDIN.LE.0)  GO TO 999
76       KDIM = MIN(KDIN,MMX)
77       IF(IMATES.LT.0) THEN
78          LXBARN=.TRUE.
79       ELSEIF(IMATES.GT.0) THEN
80          LXBARN=.FALSE.
81       ELSE
82          GOTO 999
83       ENDIF
84       IMATE=ABS(IMATES)
85 *
86       IF (JMATE.LE.0) GO TO 999
87       IF (IMATE.GT.NMATE) GO TO 80
88       JMA = LQ(JMATE-IMATE)
89       IF  (JMA.LE.0) GO TO 80
90       CALL UHTOC(IQ(JMA+1),4,NAMATE,16)
91       IF(LXBARN) THEN
92          CMIBAR=Q(JMA+6)/(AVO*Q(JMA+8))
93       ELSE
94          CMIBAR=1.
95       ENDIF
96 *
97       IF (JPART.LE.0) GO TO 999
98       IF (IPART.LE.0) GO TO 999
99       IF (IPART.GT.NPART) GO TO 80
100       JPA = LQ(JPART-IPART)
101       IF (JPA.LE.0) GO TO 80
102       CALL UHTOC(IQ(JPA+1),4,NAPART,16)
103 *
104 * *** Print  bin meaning
105       IF (IDM.GE.0) THEN
106          CHMAIL='1'
107          CALL GMAIL(0,0)
108          CHMAIL=' '
109          CHMAIL(31:)='Kinetic energy bin meaning'
110          CALL GMAIL(0,0)
111          CHMAIL(31:)='--------------------------'
112          CALL GMAIL(0,1)
113          NROW = (KDIM-1)/NCOL + 1
114          DO 20  IR=1,NROW
115             DO 10  IC=1,NCOL
116                IKB = IR + (IC-1)*NROW
117                IF (IKB.GT.KDIM) IKB=KDIM
118                KI(IC) = IKB
119                CALL GEVKEV(TKIN(IKB),EK(IC),KU(IC))
120    10       CONTINUE
121             WRITE(CHMAIL,10200) (KI(IC),EK(IC),KU(IC),IC=1,NCOL)
122             CALL GMAIL(0,0)
123    20    CONTINUE
124       ENDIF
125 *
126       BIGINV= 1000./BIG
127       DO 30  JMX = 1, MMX
128          SIGT(JMX) = 0.
129    30 CONTINUE
130       IF(MECAN.EQ.'ALL') THEN
131          N1 = 1
132          N2 = NMECA
133       ELSE
134          N1 = 0
135          DO 40  IMECA=1,NMECA
136             IF(MECAN.EQ.CHNMEC(IMECA)) THEN
137                N1 = IMECA
138             ENDIF
139    40    CONTINUE
140          IF(N1.EQ.0) THEN
141             WRITE(CHMAIL,'('' *** GPLMAT: Mechanism '',A,
142      +      '' not implemented'')') MECAN
143             CALL GMAIL(0,0)
144             GOTO 999
145          ENDIF
146          N2 = N1
147       ENDIF
148       DO 60  IMEC = N1,N2
149 C
150          IF (MECAN.EQ.'ALL') THEN
151              IF (CHNMEC(IMEC).EQ.'RANG') GO TO 60
152              IF (CHNMEC(IMEC).EQ.'STEP') GO TO 60
153          END IF
154 C
155          IF(CHNMEC(IMEC).NE.'NULL') THEN
156             MECA = CHNMEC(IMEC)
157             CALL GFTMAT(IMATE,IPART,MECA,KDIM,TKIN,VALUE,PCUT,IXST)
158             IF(IXST.EQ.0) GO TO 60
159 *
160 * ***    Book histogram
161             ISIG = 0
162             IF (MECA.EQ.'LOSS') THEN
163                CHTITL = NAPART//' in '//NAMATE//'   dE/dx (MeV/cm)'
164             ELSEIF (MECA.EQ.'RANG') THEN
165                CHTITL = NAPART//' in '//NAMATE//'   Stopping range (cm)'
166             ELSEIF (MECA.EQ.'STEP') THEN
167                CHTITL = NAPART//' in '//NAMATE//'   continuous step '
168      +         //'(cm)'
169             ELSE
170                CHTITL = NAPART//' in '//NAMATE//'   '//MECA// ' cross '
171      +         //'section'
172                IF(LXBARN) THEN
173                   CHTITL(LNBLNK(CHTITL)+1:) = ' (barn)'
174                ELSE
175                   CHTITL(LNBLNK(CHTITL)+1:) = ' (1/cm)'
176                ENDIF
177                ISIG = 1
178             ENDIF
179 *
180             ID = 10000*IMATE + 100*IPART + IMEC
181             CALL HBOOKB(ID,CHTITL,KDIM-1,TKIN,0.)
182 *
183 * ***    Fill histogram
184 *
185             VALMI = MAX (BIGINV,VMAX(VALUE,KDIM)*1.E-8)
186             DO 50  IKB = 1,KDIM
187                IF (MECA.NE.'LOSS'.AND.MECA.NE.'RANG'
188      +             .AND.MECA.NE.'STEP')
189      +             VALUE(IKB)=VALUE(IKB)*CMIBAR
190                IF (VALUE(IKB).GE.VALMI) THEN
191                   CALL HFILL(ID,TKIN(IKB),0.,VALUE(IKB))
192                ENDIF
193                IF (ISIG.EQ.1) THEN
194                   IF(MECA(1:3).NE.'INE'.AND.MECA(1:3).NE.'ELA'.AND.
195      +            MECA(1:3).NE.'FIS'.AND.MECA(1:3).NE.'CAP'.AND.
196      +            MECA(1:3).NE.'HAD'.AND.IMEC.LT.IBLOWN) THEN
197                      SIGT(IKB) = SIGT(IKB) + VALUE(IKB)
198                   ELSE IF (MECA(1:3).EQ.'HAD') THEN
199                      IF ((MECA.EQ.'HADG'.AND.IHADR.LE.2).OR. (MECA.EQ.
200      +               'HADF'.AND.IHADR.EQ.4)) THEN
201                         SIGT(IKB) = SIGT(IKB) + VALUE(IKB)
202                      ENDIF
203                   ENDIF
204                ENDIF
205    50       CONTINUE
206             CALL HIDOPT(ID,'LOGY')
207             IF(IDM.GE.0) CALL HPHIST(ID,' ',0)
208             IF(IDM.EQ.0) CALL HDELET(ID)
209          ENDIF
210    60 CONTINUE
211 *
212 * *** plot total cross section and mean free path
213       IF (MECAN.EQ.'ALL') THEN
214          CHTITL= NAPART//' in '//NAMATE//'   total cross section'
215          IF(LXBARN) THEN
216             CHTITL(LNBLNK(CHTITL)+1:) = ' (barn)'
217          ELSE
218             CHTITL(LNBLNK(CHTITL)+1:) = ' (1/cm)'
219          ENDIF
220          ID = 10000*IMATE + 100*IPART + NMECA+1
221          CALL HBOOKB(ID,CHTITL,KDIM-1,TKIN,0.)
222 *
223          CHTITL= NAPART//' in '//NAMATE//'   total mean free path (cm)'
224          II = ID + 1
225          CALL HBOOKB(II,CHTITL,KDIM-1,TKIN,0.)
226 *
227          VALMI = MAX (BIGINV,VMAX( SIGT,KDIM)*1.E-8)
228          DO 70  IKB = 1,KDIM
229             IF (SIGT(IKB).GE.VALMI) THEN
230                CALL HFILL(ID,TKIN(IKB),0.,       SIGT(IKB))
231                CALL HFILL(II,TKIN(IKB),0.,CMIBAR/SIGT(IKB))
232             ENDIF
233    70    CONTINUE
234          CALL HIDOPT(ID,'LOGY')
235          IF(IDM.GE.0) CALL HPHIST(ID,' ',0)
236          IF(IDM.EQ.0) CALL HDELET(ID)
237 *
238          CALL HIDOPT(II,'LOGY')
239          IF(IDM.GE.0) CALL HPHIST(II,' ',0)
240          IF(IDM.EQ.0) CALL HDELET(II)
241       ENDIF
242 *
243       GO TO 999
244 *
245    80 WRITE(CHMAIL,10000) IMATE ,IPART
246       CALL GMAIL(0,0)
247 10000 FORMAT(' ***** GPLMAT error : material',I4,
248      +       '  or particle',I4,' not defined'   )
249 10100 FORMAT(6X,'BCUTE =',F6.2,A4,3X,'BCUTM =',F6.2,A4,3X,
250      +             'DCUTE =',F6.2,A4,3X,'DCUTM =',F6.2,A4,3X,
251      +            'PPCUTM =',F6.2,A4 )
252 10200 FORMAT(1X,5('   bin ',I3,' =',F7.2,A4))
253   999 END