CalorimeterLevel2/src/calol2cm.for

*****************************************************************************
      INTEGER FUNCTION CALOL2CM()
c
      IMPLICIT NONE
C
      INCLUDE 'INTEST.TXT'
C
      REAL PIANO(NPLAV), VARFIT(2)
      REAL SHIFT
      REAL TG(2), BAR(2,NPLAV)
      REAL DISTX, DISTY, Y(NPLAV), YY(NPLAV)
      REAL RIG, RMASS
      REAL CHECK, ENER, CX, CY
      REAL EINF, ESUP,ENPLANEX,ENPLANEY,TENPLAX,TENPLAY
      REAL hmemor(9000000)
c
      INTEGER gtr
      INTEGER j, m, nn, i
      INTEGER IBAR(2,NPLAV), NPFIT(2)
      INTEGER Iquest(100), IWPL(2),CHTRACK
c
      COMMON / TAGLIOEN / EINF, ESUP, ENER(2)
      SAVE / TAGLIOEN /
C
      COMMON / SHIFT / SHIFT
      SAVE / SHIFT /
C
      COMMON / ANGOLO / BAR, IBAR
      SAVE / ANGOLO /
C
      COMMON / WHERE / CX, CY, PIANO
      SAVE / WHERE /
C
      COMMON / GENERAL / RIG, RMASS
      SAVE / GENERAL /
C
      COMMON / CH / CHECK
      SAVE / CH /
C
      COMMON / CALOFIT / VARFIT, NPFIT, IWPL,CHTRACK
      SAVE / CALOFIT /
C
      COMMON / pawcd / hmemor
      save / pawcd /
C
      Common / QUESTd / Iquest
      save / questd /
C
C Begin !
C
      CALOL2CM = 0;
      RMASS = 0.938
C
C Z position of Y planes, Y1 = 0., Y2 = -8.09, Y3 = -18.18, ...
C
c      print *,'NPLA IS ',NPLA
      PIANO(1) = 0.
      DO I = 2, NPLA
         IF ( MOD(I,2).EQ.0 ) THEN 
            PIANO(I) = PIANO(I-1) - 8.09
         ELSE
            PIANO(I) = PIANO(I-1) - 10.09
         ENDIF 
      ENDDO
C
      CALL VZERO(DEXY,2*LENSEV)
      CALL VZERO(BAR,2*NPLAV)
      CALL VZERO(IBAR,2*NPLAV)
      CALL VZERO(CBAR,2*NPLAV)
      CALL VZERO(CIBAR,2*NPLAV)
      CALL VZERO(QQ,4)
      CALL VZERO(Y,NPLAV)
      CALL VZERO(YY,NPLAV)      
      VARFIT(1) = 0.
      VARFIT(2) = 0.
      NPFIT(1) = 0
      NPFIT(2) = 0
      VARCFIT(1) = 0.
      VARCFIT(2) = 0.
      NPCFIT(1) = 0
      NPCFIT(2) = 0
      ELEN = 0.
      SELEN = 0.
      NSTRIP = 0.
      QTOT = 0.
      NX22 = 0.
      QX22 = 0.
      QMAX = 0.
      GTR = 0
      CHTRACK = 0
C     
      DISTX = 0.
      DISTY = 0.
C
      PLANEMAX(1) = 0
      PLANEMAX(2) = 0
      TENPLAX = 0. 
      TENPLAY = 0. 
C
      IF (GOOD2.EQ.0.AND.TRIGTY.LT.2) goto 9696
C
C     CALCULATE QTOT, NSTRIP, QQ, NX22, QX22, QMAX AND PLANEMAX
C
      DO I = 1,NPLA
         ENPLANEY = 0.
         ENPLANEX = 0.         
         DO J = 1,96
c            if ( ESTRIP(1,I,J).gt.0)
c     &           print *,' I ',I,' J ',J,' ES1 ',ESTRIP(1,I,J)
c            if ( ESTRIP(2,I,J).gt.0)
c     &           print *,' I ',I,' J ',J,' ES2 ',ESTRIP(2,I,J)
            IF ( MOD(I,2).NE.0 ) THEN
               IF ( ESTRIP(2,I,J).GT.EMIN ) THEN
                  DEXY(2,I,J) = ESTRIP(2,I,J)
                  ENPLANEY = ENPLANEY + DEXY(2,I,J)
                  NSTRIP = NSTRIP + 1.
                  QTOT = QTOT + DEXY(2,I,J)
                  IF (I.LT.11) QQ(1) = QQ(1) + DEXY(2,I,J)
               ENDIF
               IF ( ESTRIP(1,I,J).GT.EMIN ) THEN
C               IF ( ESTRIP(1,I+1,J).GT.EMIN ) THEN
C                  
c                  print *,' i ',i,' i + 1 ',I+1
C                  DEXY(1,I,J) = ESTRIP(1,I+1,J)
                  DEXY(1,I,J) = ESTRIP(1,I,J)
C
                  ENPLANEX = ENPLANEX + DEXY(1,I,J)
                  NSTRIP = NSTRIP + 1.
                  QTOT = QTOT + DEXY(1,I,J)
                  if (i.lt.11) QQ(2) = QQ(2) + DEXY(1,I,J)
c                  IF (I.EQ.21) THEN
c                     NX22 = NX22 + 1.
c                     QX22 = QX22 + DEXY(1,I,J)
c                  ENDIF
               ENDIF
            ENDIF
            IF ( MOD(I,2).EQ.0 ) THEN
               IF (ESTRIP(2,I,J).GT.EMIN) THEN
                  DEXY(2,I,J) = ESTRIP(2,I,J)
                  ENPLANEY = ENPLANEY + DEXY(2,I,J)
                  NSTRIP = NSTRIP + 1.
                  QTOT = QTOT + DEXY(2,I,J)
                  if (i.lt.11) QQ(3) = QQ(3) + DEXY(2,I,J)
               ENDIF
               IF (ESTRIP(1,I,J).GT.EMIN) THEN
C               IF (ESTRIP(1,I-1,J).GT.EMIN) THEN
C                  
c                  print *,' i ',i,' i - 1 ',I-1
C                  DEXY(1,I,J) = ESTRIP(1,I-1,J)
                  DEXY(1,I,J) = ESTRIP(1,I,J)
C
                  ENPLANEX = ENPLANEX + DEXY(1,I,J)
                  NSTRIP = NSTRIP + 1.
                  QTOT = QTOT + DEXY(1,I,J)
                  IF (I.EQ.NPLA) THEN
                     NX22 = NX22 + 1.
                     QX22 = QX22 + DEXY(1,I,J)
                  ENDIF
                  IF (I.LT.11) QQ(4) = QQ(4) + DEXY(1,I,J)
               ENDIF
            ENDIF
C
            IF (DEXY(1,I,J).GT.QMAX) QMAX = DEXY(1,I,J)
            IF (DEXY(2,I,J).GT.QMAX) QMAX = DEXY(2,I,J)
C
         ENDDO
         IF (ENPLANEX.GT.TENPLAX) THEN
            PLANEMAX(1) = I
            TENPLAX = ENPLANEX
         ENDIF
         IF (ENPLANEY.GT.TENPLAY) THEN
            PLANEMAX(2) = I
            TENPLAY = ENPLANEY
         ENDIF
      ENDDO
C
C     FIT THE TRACK INSIDE THE CALORIMETER
C 
      if (GOOD2.EQ.1.OR.TRIGTY.GE.2) then
c         print *,' CALL CLUSTER '
         CALL CLUSTER     
c         print *,' CALL DIRECTION '
         CALL DIRECTION(TG)
         NPCFIT(1) = NPFIT(1)
         NPCFIT(2) = NPFIT(2)
C     
         DO M = 1,2
C
            IF (NPCFIT(M).GE.2) THEN
               VARCFIT(M) = VARFIT(M)
               IF (M.EQ.1) THEN
                  IMPX = (CX - XALIG)/10.
                  TANX = TG(1)
               ELSE
                  IMPY =( CY - YALIG)/10. 
                  TANY = TG(2)
               ENDIF
C               
               DO I = 1,NPLA     
C                  NN = 1
                  NN = 0
C                  IF (M.EQ.2) NN = 1
                  IF (MOD(I,2).EQ.NN) THEN
                     IF (REVERSE.EQ.0) THEN
                        SHIFT = -0.5
                     ELSE
                        SHIFT = +0.5
                     ENDIF
                  ELSE
                     IF (REVERSE.EQ.0) THEN
                        SHIFT = +0.5
                     ELSE
                        SHIFT = -0.5
                     ENDIF
                  ENDIF
C     
                  IF (M.EQ.1) THEN
                     DISTX = PIANO(I) - 5.81
c                     print *,'C X PLANE I= ',I,' Z = ',DISTX+ZALIG
                     Y(I) = DISTX * TG(1) + CX
                     BAR(M,I) = Y(I)
                     CBAR(M,I) = (Y(I) - XALIG)/10.        
C     
                  ELSE
                     DISTY = PIANO(I)                
c                     print *,'C Y PLANE I= ',I,' Z = ',DISTY+ZALIG
                     YY(I) = DISTY * TG(2) + CY
                     BAR(M,I) = YY(I)
                     CBAR(M,I) = (-YALIG + YY(I))/10.
C     
                  ENDIF
                  CALL LASTRISCIA(BAR(M,I),IBAR(M,I))
                  CIBAR(M,I) = IBAR(M,I)
                  IF (ibar(m,i).EQ.-1) THEN 
                     CHTRACK = CHTRACK + 1
                  ELSE
                     IWPL(M) = IWPL(M) + 1
                  ENDIF
               ENDDO            
            ENDIF
C
         ENDDO
      ELSE
         GOTO 9696
      ENDIF
C
C     CALCULATE ELEN AND SELEN
C
c      ELEN =  40.82 * QTOT * 1.06 * 1E-4  
C      ELEN =  260. * QTOT   
      ELEN =  QTOT / 260.   
      IF (ELEN.GT.0.) THEN
         SELEN = ELEN * (0.01183 + 0.121 / SQRT(ELEN))
      ELSE
         SELEN = 0.
      ENDIF
C
 9696 CONTINUE
C
c 45   continue

c 50   continue

      return
      END


1	mocchiut	1.1	*****************************************************************************
2			INTEGER FUNCTION CALOL2CM()
3			c
4			IMPLICIT NONE
5			C
6			INCLUDE 'INTEST.TXT'
7			C
8	mocchiut	1.8	REAL PIANO(NPLAV), VARFIT(2)
9	mocchiut	1.1	REAL SHIFT
10	mocchiut	1.8	REAL TG(2), BAR(2,NPLAV)
11			REAL DISTX, DISTY, Y(NPLAV), YY(NPLAV)
12	mocchiut	1.1	REAL RIG, RMASS
13			REAL CHECK, ENER, CX, CY
14			REAL EINF, ESUP,ENPLANEX,ENPLANEY,TENPLAX,TENPLAY
15			REAL hmemor(9000000)
16			c
17			INTEGER gtr
18			INTEGER j, m, nn, i
19	mocchiut	1.8	INTEGER IBAR(2,NPLAV), NPFIT(2)
20	mocchiut	1.5	INTEGER Iquest(100), IWPL(2),CHTRACK
21	mocchiut	1.1	c
22			COMMON / TAGLIOEN / EINF, ESUP, ENER(2)
23			SAVE / TAGLIOEN /
24			C
25			COMMON / SHIFT / SHIFT
26			SAVE / SHIFT /
27			C
28			COMMON / ANGOLO / BAR, IBAR
29			SAVE / ANGOLO /
30			C
31			COMMON / WHERE / CX, CY, PIANO
32			SAVE / WHERE /
33			C
34			COMMON / GENERAL / RIG, RMASS
35			SAVE / GENERAL /
36			C
37			COMMON / CH / CHECK
38			SAVE / CH /
39			C
40	mocchiut	1.5	COMMON / CALOFIT / VARFIT, NPFIT, IWPL,CHTRACK
41	mocchiut	1.1	SAVE / CALOFIT /
42			C
43			COMMON / pawcd / hmemor
44			save / pawcd /
45			C
46			Common / QUESTd / Iquest
47			save / questd /
48			C
49			C Begin !
50			C
51			CALOL2CM = 0;
52			RMASS = 0.938
53			C
54			C Z position of Y planes, Y1 = 0., Y2 = -8.09, Y3 = -18.18, ...
55			C
56	mocchiut	1.8	c print *,'NPLA IS ',NPLA
57	mocchiut	1.1	PIANO(1) = 0.
58	mocchiut	1.8	DO I = 2, NPLA
59	mocchiut	1.1	IF ( MOD(I,2).EQ.0 ) THEN
60			PIANO(I) = PIANO(I-1) - 8.09
61			ELSE
62			PIANO(I) = PIANO(I-1) - 10.09
63			ENDIF
64			ENDDO
65			C
66			CALL VZERO(DEXY,2*LENSEV)
67	mocchiut	1.8	CALL VZERO(BAR,2*NPLAV)
68			CALL VZERO(IBAR,2*NPLAV)
69			CALL VZERO(CBAR,2*NPLAV)
70			CALL VZERO(CIBAR,2*NPLAV)
71	mocchiut	1.1	CALL VZERO(QQ,4)
72	mocchiut	1.8	CALL VZERO(Y,NPLAV)
73			CALL VZERO(YY,NPLAV)
74	mocchiut	1.1	VARFIT(1) = 0.
75			VARFIT(2) = 0.
76	mocchiut	1.14	NPFIT(1) = 0
77			NPFIT(2) = 0
78	mocchiut	1.1	VARCFIT(1) = 0.
79			VARCFIT(2) = 0.
80	mocchiut	1.14	NPCFIT(1) = 0
81			NPCFIT(2) = 0
82	mocchiut	1.1	ELEN = 0.
83			SELEN = 0.
84			NSTRIP = 0.
85			QTOT = 0.
86			NX22 = 0.
87			QX22 = 0.
88			QMAX = 0.
89			GTR = 0
90	mocchiut	1.5	CHTRACK = 0
91	mocchiut	1.1	C
92			DISTX = 0.
93			DISTY = 0.
94			C
95			PLANEMAX(1) = 0
96			PLANEMAX(2) = 0
97			TENPLAX = 0.
98			TENPLAY = 0.
99			C
100			IF (GOOD2.EQ.0.AND.TRIGTY.LT.2) goto 9696
101			C
102			C CALCULATE QTOT, NSTRIP, QQ, NX22, QX22, QMAX AND PLANEMAX
103			C
104	mocchiut	1.8	DO I = 1,NPLA
105	mocchiut	1.1	ENPLANEY = 0.
106			ENPLANEX = 0.
107			DO J = 1,96
108	mocchiut	1.8	c if ( ESTRIP(1,I,J).gt.0)
109			c & print *,' I ',I,' J ',J,' ES1 ',ESTRIP(1,I,J)
110			c if ( ESTRIP(2,I,J).gt.0)
111			c & print *,' I ',I,' J ',J,' ES2 ',ESTRIP(2,I,J)
112	mocchiut	1.1	IF ( MOD(I,2).NE.0 ) THEN
113			IF ( ESTRIP(2,I,J).GT.EMIN ) THEN
114			DEXY(2,I,J) = ESTRIP(2,I,J)
115	mocchiut	1.7	ENPLANEY = ENPLANEY + DEXY(2,I,J)
116	mocchiut	1.1	NSTRIP = NSTRIP + 1.
117	mocchiut	1.7	QTOT = QTOT + DEXY(2,I,J)
118			IF (I.LT.11) QQ(1) = QQ(1) + DEXY(2,I,J)
119	mocchiut	1.1	ENDIF
120			IF ( ESTRIP(1,I,J).GT.EMIN ) THEN
121	mocchiut	1.7	C IF ( ESTRIP(1,I+1,J).GT.EMIN ) THEN
122			C
123			c print *,' i ',i,' i + 1 ',I+1
124			C DEXY(1,I,J) = ESTRIP(1,I+1,J)
125	mocchiut	1.1	DEXY(1,I,J) = ESTRIP(1,I,J)
126	mocchiut	1.7	C
127			ENPLANEX = ENPLANEX + DEXY(1,I,J)
128	mocchiut	1.1	NSTRIP = NSTRIP + 1.
129	mocchiut	1.7	QTOT = QTOT + DEXY(1,I,J)
130			if (i.lt.11) QQ(2) = QQ(2) + DEXY(1,I,J)
131			c IF (I.EQ.21) THEN
132			c NX22 = NX22 + 1.
133			c QX22 = QX22 + DEXY(1,I,J)
134			c ENDIF
135	mocchiut	1.1	ENDIF
136			ENDIF
137			IF ( MOD(I,2).EQ.0 ) THEN
138			IF (ESTRIP(2,I,J).GT.EMIN) THEN
139			DEXY(2,I,J) = ESTRIP(2,I,J)
140	mocchiut	1.7	ENPLANEY = ENPLANEY + DEXY(2,I,J)
141	mocchiut	1.1	NSTRIP = NSTRIP + 1.
142	mocchiut	1.7	QTOT = QTOT + DEXY(2,I,J)
143			if (i.lt.11) QQ(3) = QQ(3) + DEXY(2,I,J)
144	mocchiut	1.1	ENDIF
145			IF (ESTRIP(1,I,J).GT.EMIN) THEN
146	mocchiut	1.7	C IF (ESTRIP(1,I-1,J).GT.EMIN) THEN
147			C
148			c print *,' i ',i,' i - 1 ',I-1
149			C DEXY(1,I,J) = ESTRIP(1,I-1,J)
150	mocchiut	1.1	DEXY(1,I,J) = ESTRIP(1,I,J)
151	mocchiut	1.7	C
152			ENPLANEX = ENPLANEX + DEXY(1,I,J)
153	mocchiut	1.1	NSTRIP = NSTRIP + 1.
154	mocchiut	1.7	QTOT = QTOT + DEXY(1,I,J)
155	mocchiut	1.8	IF (I.EQ.NPLA) THEN
156	mocchiut	1.1	NX22 = NX22 + 1.
157	mocchiut	1.7	QX22 = QX22 + DEXY(1,I,J)
158	mocchiut	1.1	ENDIF
159	mocchiut	1.7	IF (I.LT.11) QQ(4) = QQ(4) + DEXY(1,I,J)
160	mocchiut	1.1	ENDIF
161			ENDIF
162			C
163			IF (DEXY(1,I,J).GT.QMAX) QMAX = DEXY(1,I,J)
164			IF (DEXY(2,I,J).GT.QMAX) QMAX = DEXY(2,I,J)
165			C
166			ENDDO
167			IF (ENPLANEX.GT.TENPLAX) THEN
168			PLANEMAX(1) = I
169			TENPLAX = ENPLANEX
170			ENDIF
171			IF (ENPLANEY.GT.TENPLAY) THEN
172			PLANEMAX(2) = I
173			TENPLAY = ENPLANEY
174			ENDIF
175			ENDDO
176			C
177			C FIT THE TRACK INSIDE THE CALORIMETER
178			C
179			if (GOOD2.EQ.1.OR.TRIGTY.GE.2) then
180	mocchiut	1.8	c print *,' CALL CLUSTER '
181	mocchiut	1.1	CALL CLUSTER
182	mocchiut	1.8	c print *,' CALL DIRECTION '
183	mocchiut	1.1	CALL DIRECTION(TG)
184			NPCFIT(1) = NPFIT(1)
185			NPCFIT(2) = NPFIT(2)
186			C
187			DO M = 1,2
188			C
189			IF (NPCFIT(M).GE.2) THEN
190			VARCFIT(M) = VARFIT(M)
191			IF (M.EQ.1) THEN
192	mocchiut	1.2	IMPX = (CX - XALIG)/10.
193	mocchiut	1.1	TANX = TG(1)
194			ELSE
195	mocchiut	1.3	IMPY =( CY - YALIG)/10.
196	mocchiut	1.1	TANY = TG(2)
197			ENDIF
198			C
199			DO I = 1,NPLA
200	mocchiut	1.10	C NN = 1
201			NN = 0
202	mocchiut	1.7	C IF (M.EQ.2) NN = 1
203	mocchiut	1.1	IF (MOD(I,2).EQ.NN) THEN
204	mocchiut	1.8	IF (REVERSE.EQ.0) THEN
205	mocchiut	1.11	SHIFT = -0.5
206			ELSE
207	mocchiut	1.8	SHIFT = +0.5
208			ENDIF
209	mocchiut	1.1	ELSE
210	mocchiut	1.8	IF (REVERSE.EQ.0) THEN
211	mocchiut	1.11	SHIFT = +0.5
212			ELSE
213	mocchiut	1.8	SHIFT = -0.5
214			ENDIF
215	mocchiut	1.1	ENDIF
216			C
217			IF (M.EQ.1) THEN
218			DISTX = PIANO(I) - 5.81
219	mocchiut	1.9	c print *,'C X PLANE I= ',I,' Z = ',DISTX+ZALIG
220	mocchiut	1.1	Y(I) = DISTX * TG(1) + CX
221			BAR(M,I) = Y(I)
222			CBAR(M,I) = (Y(I) - XALIG)/10.
223			C
224			ELSE
225			DISTY = PIANO(I)
226	mocchiut	1.9	c print *,'C Y PLANE I= ',I,' Z = ',DISTY+ZALIG
227	mocchiut	1.1	YY(I) = DISTY * TG(2) + CY
228			BAR(M,I) = YY(I)
229	mocchiut	1.3	CBAR(M,I) = (-YALIG + YY(I))/10.
230	mocchiut	1.1	C
231			ENDIF
232			CALL LASTRISCIA(BAR(M,I),IBAR(M,I))
233			CIBAR(M,I) = IBAR(M,I)
234	mocchiut	1.5	IF (ibar(m,i).EQ.-1) THEN
235			CHTRACK = CHTRACK + 1
236			ELSE
237			IWPL(M) = IWPL(M) + 1
238			ENDIF
239	mocchiut	1.1	ENDDO
240			ENDIF
241			C
242			ENDDO
243			ELSE
244			GOTO 9696
245			ENDIF
246			C
247			C CALCULATE ELEN AND SELEN
248			C
249	mocchiut	1.10	c ELEN = 40.82 * QTOT * 1.06 * 1E-4
250	mocchiut	1.12	C ELEN = 260. * QTOT
251			ELEN = QTOT / 260.
252	mocchiut	1.1	IF (ELEN.GT.0.) THEN
253			SELEN = ELEN * (0.01183 + 0.121 / SQRT(ELEN))
254			ELSE
255			SELEN = 0.
256			ENDIF
257			C
258			9696 CONTINUE
259			C
260	mocchiut	1.13	c 45 continue
261	mocchiut	1.1
262	mocchiut	1.13	c 50 continue
263	mocchiut	1.1
264			return
265			END
266
267