*********************************************************************** * * Subroutine READVME * Routine to read a single event: use C subroutine readevnt. * Use conditional compilation with cpp preprocessor, so this file must * end in .F and cannot be included with INCLUDE directive (use #include) * Input: IEV event # to read * Output: fill ADC and DSP variables * * Rewritten from a previous routine from an HP machine * S. Piccardi Nov. 1998 * * $Id: readvme.F,v 1.1 1999/12/15 10:05:56 piccardi Exp $ * ***************************************************************************** PROGRAM NNPLOT_NT5 c IMPLICIT NONE C INCLUDE 'INTEST.TXT' C CHARACTER*6 file_name, file_name2 CHARACTER*40 file_name3 C integer icount,icount2 c integer RUNERROR, PEDERROR C C Normal variables definition C INTEGER FFD C INTEGER l, j, kk, m, mm, ll, il, ii, nn, ival, iev INTEGER i, istat, icycle, ierr, icheck, icheck2, icheck3 INTEGER LVMIN INTEGER*2 length integer IP INTEGER ipre, ipre2, jc INTEGER ic, k, jj INTEGER status INTEGER inf, isup, idet INTEGER LEVEL REAL ADCERR C$$$ INTEGER*4 lenght(4) INTEGER lundata, iosop INTEGER*2 st, st2 CHARACTER*2 nr c integer*4 evnum, systime, len integer nev C integer nset,nword(10) C INTEGER*2 STRINGA(200000) C REAL auto(11) integer lun(4) REAL SUM(6), SUM2(6), COUNT(6), & BAS(4,11,6), DIS(6), DDIS(6), BASE(2) REAL SSUM, SSUM2, CCOUNT, MME(2), SSIG REAL VAL REAL PED(4,11,96), GOOD(4,11,96) COMMON / PEDE / PED, GOOD REAL MIP1(11,96), MIP2(11,96), MIP3(11,96), MIP4(11,96) REAL SPOS1(3,11), SPOS2(3,11), SPOS3(3,11), SPOS4(3,11), DET REAL D1(11,96), D2(11,96), D3(11,96), D4(11,96) REAL DI, DS REAL ME(6), SIG(6) REAL TG(2) REAL SHIFT COMMON /SHIFT/ SHIFT REAL BAR(2,NPLA) INTEGER IBAR(2,NPLA) COMMON/ANGOLO/BAR,IBAR REAL DISTX, DISTY, Y(NPLA), YY(NPLA) REAL CX, CY COMMON/WHERE/CX,CY REAL RIG, PLANEMAX, RMASS COMMON/GENERAL/RIG,RMASS INTEGER IPLANE, NNX, NNY, INFX, INFY, ISUPX, ISUPY INTEGER NLK REAL RNSS, QTOTT, RQT real dedx1(11,96),dedx2(11,96),dedx3(11,96),dedx4(11,96) real qpl(96), ene(96), qst(96) INTEGER nin REAL CHECK COMMON / CH / CHECK REAL QMAXX, QMAXY REAL VMAX real nstrip, qtot, ncore, qcore, impx, impy, tanx, tany, nint REAL ncyl, qcyl, qtrack, qmax, cher, lead, nx22, qx22, qq(4) INTEGER event COMMON / Planes / nstrip, qtot, ncore, qcore, impx, impy, & tanx, tany, nint, ncyl, qcyl, qtrack, qmax, cher, lead, nx22, & qx22, qq, event INTEGER lrec PARAMETER (lrec=6144) REAL hmemor(3000000) COMMON /pawc/hmemor CALL HLIMIT(3000000) C C Begin ! C RUNERROR = 0 ! error variable PEDERROR = 0 PRINT *,'Data file to read (only data and number included)? ' READ(*,905)file_name CALL PEDESTAL(PEDERROR) IF (PEDERROR.EQ.1) THEN PRINT *,'Problems with pedestal file' GOTO 50 ENDIF C 905 FORMAT(A6) lundata=33 OPEN(UNIT=lundata,FILE='/wizard3/pamela/tb-0903/data/output_0309' + //file_name//'.dat',STATUS='OLD', + FORM='UNFORMATTED',ERR=50,IOSTAT=IOSOP) C C C PRINT *,'File to save? ' READ(*,904)file_name3 904 FORMAT(A40) C C Histos creation C CALL HROPEN(59,'Event','/wizard3/pamela/tb-0903/paw/'// + file_name3,'n',lrec,istat) CALL HBNT(1,'Pamela Calo',' ') CALL HBSET('BSIZE',lrec,ierr) *** /* Book ntuple variables */ CALL HBNAME(1,'calo',nstrip,'nstrip:R,qtot:R, & ncore:R,qcore:R,impx:R,impy:R,tanx:R,tany:R, & nint:R,ncyl:R,qcyl:R,qtrack:R,qmax:R, & cher:R,lead:R,nx22:R,qx22:R,qq(4):R,event:I') C 906 FORMAT(A2) C OPEN(71,FILE='/wizard3/pamela/tb-0903/source/mip1.dat', & STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP) OPEN(72,FILE='/wizard3/pamela/tb-0903/source/mip2.dat', & STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP) OPEN(73,FILE='/wizard3/pamela/tb-0903/source/nmip3.dat', & STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP) OPEN(74,FILE='/wizard3/pamela/tb-0903/source/mip4.dat', & STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP) C OPEN(81,FILE='/wizard3/pamela/tb-0903/source/spos1.dat', & STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP) OPEN(82,FILE='/wizard3/pamela/tb-0903/source/spos4.dat', & STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP) OPEN(83,FILE='/wizard3/pamela/tb-0903/source/spos3.dat', & STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP) OPEN(84,FILE='/wizard3/pamela/tb-0903/source/spos2.dat', & STATUS='OLD',FORM='FORMATTED',ERR=50,IOSTAT=IOSOP) C DO I = 1,11 READ (71,102)(mip1(I,J),J=1,96) READ (72,102)(mip2(I,J),J=1,96) READ (73,102)(mip3(I,J),J=1,96) READ (74,102)(mip4(I,J),J=1,96) ENDDO DO I = 1,3 READ (81,103)(SPOS1(I,J),J=1,11) READ (82,103)(SPOS2(I,J),J=1,11) READ (83,103)(SPOS3(I,J),J=1,11) READ (84,103)(SPOS4(I,J),J=1,11) print *,'I :',I,' J : 2 spos1 :',SPOS1(I,2) ENDDO CLOSE(21) CLOSE(71) CLOSE(72) CLOSE(73) CLOSE(74) CLOSE(81) CLOSE(82) CLOSE(83) CLOSE(84) 102 FORMAT (4(2X,F9.3)) 103 FORMAT (11(2X,F7.5)) C C Take the file descriptor number (unix) C PRINT *,'Momentum:' READ *,RIG C PRINT *,'Mass of the particle:' READ *,RMASS C PRINT *,'How many events to read:' READ *,nev C ffd=FNum(lundata) ! take the file descriptor number C C Read event: use C readevent routine C icount = 0 icount2 = 0 event = 0 iev = 0. CALL VZERO(BAS,66*4) c 10 continue C event = event + 1 ival = 0 DO I = 1,11 DO J = 1,96 DEDX1(I,J) = 0. DEDX2(I,J) = 0. DEDX3(I,J) = 0. DEDX4(I,J) = 0. ENDDO ENDDO C call readev(evnum,systime,len,stringa,RUNERROR,ffd) ! XXX if (runerror.eq.-1.or.runerror.eq.1) goto 50 C c write (*,31) evnum C write (*,32) systime C write (*,33) len C print *,'event number:',evnum C print *,'systime:',systime C print *,'event length:',len 31 format(1x,'event number :',1x,Z8) 32 format(1x,'systime :',1x,Z8) 33 format(1x,'event length :',1x,Z8) C$$$ write(6,*)nset,nword C$$ print *,'lenght 1 :',lenght(1) C$$ print *,'lenght 2 :',lenght(2) C$$ print *,'lenght 3 :',lenght(3) if (runerror.eq.-1.or.runerror.eq.1) goto 50 ic = 0 ic = ic + 4 length = ic do kk = 1,4 C if (kk.eq.3) goto 45 ic = length + 1 C$$$ do k = ic-5,ic+5 k = ic C write (*,41) stringa(k) st2 = IAND(stringa(k),'FF00'x) status = ISHFT(st2,-8) C write (*,41) status IF (KK.eq.1.and.status.ne.234) then print *,'Problems with view:',KK GOTO 45 endif IF (KK.eq.2.and.status.ne.246) then print *,'Problems with view:',KK GOTO 45 endif IF (KK.eq.3.and.status.ne.241) then print *,'Problems with view:',KK GOTO 45 endif IF (KK.eq.4.and.status.ne.237) then print *,'Problems with view:',KK GOTO 45 endif C$$$ enddo if (status.ne.-30460) then c print *,'problems with DSP, event:',event endif C 41 format(1x,'status :',1x,Z4) C ic = ic + 1 lun(kk) = stringa(ic) length = length + (stringa(ic) + 2) C print *,'Lenght1 :',length c write (*,41) stringa(ic) C print *,'Lenght2 :',stringa(ic) C$$$ if (stringa(ic).ne.1064) then C$$$ print *,'problems with DSP' C$$$ endif C do i = 1,11 ic = ic + 1 auto(i) = stringa(ic) enddo C do i = 1,11 do j = 1,96 ic = ic + 1 if (kk.eq.1) then d1(i,j) = stringa(ic) - & ped(kk,i,j) if (i.eq.4.and.j.gt.64) then c print *,'d1 :',d1(i,j),' ped :',ped(kk,i,j), c & ' stringa :',stringa(ic),' j :',j endif endif if (kk.eq.2) then d2(i,j) = stringa(ic) - & ped(kk,i,j) endif if (kk.eq.3) then d3(i,j) = stringa(ic) - & ped(kk,i,j) endif if (kk.eq.4) then d4(i,j) = stringa(ic) - & ped(kk,i,j) C$ if (i.eq.2.and.j.gt.32.and.j.lt.49.and.event.eq.42) C$ & then C$ print *,'d4 :',d4(i,j),' ped :',ped(kk,i,j), C$ & ' stringa :',stringa(ic),' j :',j C$ endif endif enddo enddo C DO I = 1,11 C DO M = 1,6 INF = (M - 1) * 16 + 1 ISUP = M * 16 CHECK = 0. C$ IF (M.EQ.3.AND.I.EQ.2.AND.KK.EQ.4.and.event.eq.42) C$ & CHECK = 1. IF (KK.EQ.1) CALL CALCBASE(KK,INF,ISUP,I,D1,BAS(KK,I,M)) c IF (KK.EQ.1.and.i.eq.4) PRINT *,'BASE 1:',BAS(KK,I,M), c & ' m :',m IF (KK.EQ.2) CALL CALCBASE(KK,INF,ISUP,I,D2,BAS(KK,I,M)) IF (KK.EQ.3) CALL CALCBASE(KK,INF,ISUP,I,D3,BAS(KK,I,M)) c IF (M.EQ.4.AND.I.EQ.4.AND.KK.EQ.4) c & CHECK = 1. IF (KK.EQ.4) CALL CALCBASE(KK,INF,ISUP,I,D4,BAS(KK,I,M)) ENDDO C C do j = 1,96 IP = INT((J - 1) / 16) + 1 IF (KK.EQ.1) THEN IF (BAS(1,I,IP).GT.0..AND.MIP1(I,J).GT.0.) THEN DEDX1(I,J) = (D1(I,J) - BAS(1,I,IP)) / MIP1(I,J) C if (i.eq.1) print *,'dedx :',dedx1(I,j), C & ' base :',bas(1,i,ip),' mip :',mip1(i,j) if (i.eq.4.and.j.gt.64) then c print *,'dedx1 :',dedx1(i,j),' i :',i,' j :',j endif ELSE DEDX1(I,J) = 0. ENDIF ENDIF IF (KK.EQ.2) THEN IF (BAS(2,I,IP).GT.0..AND.MIP2(I,J).GT.0.) THEN DEDX2(I,J) = (D2(I,J) - BAS(2,I,IP)) / MIP2(I,J) ELSE DEDX2(I,J) = 0. ENDIF ENDIF c IF (KK.EQ.3) THEN c IF (BAS(3,I,IP).GT.0.) THEN c DEDX3(I,J) = (D3(I,J) - BAS(3,I,IP)) / MIP3(I,J) c ELSE c DEDX3(I,J) = 0. c ENDIF c ENDIF IF (KK.EQ.4) THEN IF (BAS(4,I,IP).GT.0..AND.MIP4(I,J).GT.0.) THEN DEDX4(I,J) = (D4(I,J) - BAS(4,I,IP)) / MIP4(I,J) c IF (I.EQ.4.and.j.gt.48.and.j.lt.65) c PRINT *,'DEDX4 :',DEDX4(I,J),' BASE :', c & BAS(4,I,IP),' D4 :',D4(I,J),' J :',J ELSE DEDX4(I,J) = 0. ENDIF ENDIF enddo c enddo C ic = ic + 1 C 45 continue enddo C DO M = 1,4 DO K = 1,11 C DO I = 1,96 IF (M.EQ.1.AND.MIP1(K,I).NE.0.) QPL(I) = DEDX1(K,I) IF (M.EQ.2.AND.MIP2(K,I).NE.0.) QPL(I) = DEDX2(K,I) IF (M.EQ.3.AND.MIP3(K,I).NE.0.) QPL(I) = DEDX3(K,I) IF (M.EQ.4.AND.MIP4(K,I).NE.0.) QPL(I) = DEDX4(K,I) ENE(I) = QPL(I) QST(I) = QPL(I) ENDDO C DO I = 1,96 IPRE = INT((I-1)/16) + 1 JC = MOD(I,32) IF (QPL(I).GT.EMIN) THEN C INF = (IPRE - 1) * 16 + 1 ISUP = IPRE * 16 DO J = INF,ISUP C$$$ ENE(J) = ENE(J) + QPL(I) * 0.00963 ENE(J) = ENE(J) + QPL(I) * 0.00478 ENDDO C IF (MOD(IPRE,2).EQ.0) IPRE2 = IPRE - 1 IF (MOD(IPRE,2).EQ.1) IPRE2 = IPRE + 1 INF = (IPRE2 - 1) * 16 + 1 ISUP = IPRE2 * 16 C$$$ DO J = INF,ISUP C$$$ ENE(J) = ENE(J) + QPL(I) * 0.00485 C$$$ ENDDO C IF (JC.NE.1) ENE(I-1) = ENE(I-1) - QPL(I)*0.01581 IF (JC.NE.0) ENE(I+1) = ENE(I+1) - QPL(I)*0.01581 C ENDIF ENDDO C DO I = 1,96 IPRE = INT((I-1)/16) + 1 JC = MOD(I,32) IF (ENE(I).GT.EMIN) THEN C INF = (IPRE - 1) * 16 + 1 ISUP = IPRE * 16 IDET = INT((I - 1) / 32) + 1 IF (M.EQ.1) DET = SPOS1(IDET,K) IF (M.EQ.2) DET = SPOS2(IDET,K) IF (M.EQ.3) DET = SPOS3(IDET,K) IF (M.EQ.4) DET = SPOS4(IDET,K) DO J = INF,ISUP C$$$ QST(J) = QST(J) + ENE(I) * 0.00963 QST(J) = QST(J) + ENE(I) * (0.00478 + & DET) ENDDO C IF (MOD(IPRE,2).EQ.0) IPRE2 = IPRE - 1 IF (MOD(IPRE,2).EQ.1) IPRE2 = IPRE + 1 INF = (IPRE2 - 1) * 16 + 1 ISUP = IPRE2 * 16 DO J = INF,ISUP QST(J) = QST(J) + ENE(I) * DET ENDDO C IF (JC.NE.1) QST(I-1) = QST(I-1) - ENE(I)*0.01581 IF (JC.NE.0) QST(I+1) = QST(I+1) - ENE(I)*0.01581 C ENDIF ENDDO C DO I = 1,96 C IF (M.EQ.1.AND.MIP1(K,I).NE.0.) DEDX1(K,I) = QST(I) c IF (M.EQ.2.AND.MIP2(K,I).NE.0.) DEDX2(K,I) = QST(I) c IF (M.EQ.3.AND.MIP3(K,I).NE.0.) DEDX3(K,I) = QST(I) c IF (M.EQ.4.AND.MIP4(K,I).NE.0.) DEDX4(K,I) = QST(I) IF (M.EQ.4.AND.K.EQ.5) DEDX4(K,I) = 0. ENDDO C ENDDO C ENDDO C CALL VZERO(DEXY,2*LENSEV) CALL VZERO(QQ,4) NSTRIP = 0. QTOT = 0. NX22 = 0. QX22 = 0. DO I = 1,11 DO J = 1,96 ICHECK3 = 1 IF (DEDX1(I,J).GT.EMIN) & THEN DEXY(2,2*I-1,97-J) = DEDX1(I,J) NSTRIP = NSTRIP + 1. QTOT = QTOT + DEDX1(I,J) IF (I.NE.1) QQ(1) = QQ(1) + DEDX1(I,J) c if (i.eq.1) print *,'qtot1 :',qtot,' dedx :',dedx1(i,j) ENDIF IF (DEDX2(I,J).GT.EMIN) THEN DEXY(1,2*I-1,J) = DEDX2(I,J) NSTRIP = NSTRIP + 1 QTOT = QTOT + DEDX2(I,J) IF (I.EQ.11) THEN NX22 = NX22 + 1. QX22 = QX22 + DEDX2(I,J) ENDIF if (i.lt.11) QQ(2) = QQ(2) + DEDX2(I,J) c print *,'qtot2 :',qtot ENDIF IF (DEDX3(I,J).GT.EMIN) THEN C DEXY(2,2*I,J) = DEDX3(I,J) C NSTRIP = NSTRIP + 1 C QTOT = QTOT + DEDX3(I,J) C QQ(3) = QQ(3) + DEDX3(I,J) c print *,'qtot3 :',qtot ENDIF IF (DEDX4(I,J).GT.EMIN) THEN DEXY(1,2*I,J) = DEDX4(I,J) NSTRIP = NSTRIP + 1 QTOT = QTOT + DEDX4(I,J) IF (I.LT.11) QQ(4) = QQ(4) + DEDX4(I,J) c if (i.eq.5) PRINT *,'DEDX4 :',DEDX4(I,J),' qtot :', c & qtot,' J :',J,' I :',i c print *,'qtot4 :',qtot ENDIF ENDDO ENDDO C CALL CLUSTER C CALL DIRECTION(TG) IMPX = CX - 241. / 2. IMPY = CY - 241. / 2. SHIFT = -0.5 CALL LASTRISCIA(CX,II) IMPX = II SHIFT = +0.5 CALL LASTRISCIA(CY,II) IMPY = II TANX = TG(1) TANY = TG(2) C DO M = 1,2 DO I = 1,NPLA C NN = 0 IF (M.EQ.2) NN = 1 IF (MOD(I,2).EQ.NN) THEN SHIFT = +0.5 ELSE SHIFT = -0.5 ENDIF C DISTY = -PIANO * (I - 1.) DISTX = -PIANO * (I - 1.) C IF (M.EQ.1) THEN Y(I) = DISTX * TG(1) + CX BAR(M,I) = Y(I) C ELSE YY(I) = DISTY * TG(2) + CY BAR(M,I) = YY(I) C ENDIF CALL LASTRISCIA(BAR(M,I),IBAR(M,I)) ENDDO ENDDO C NCORE = 0. QCORE = 0. PLANEMAX = 1.01*(LOG(ABS(RIG)/0.0081)-1.) IPLANE = INT(ANINT(PLANEMAX)) + 5 IF (IPLANE.GT.NPLA) IPLANE=NPLA DO J = 1,IPLANE C NNX = IBAR(1,J) NNY = IBAR(2,J) IF (NNX.LT.9) NNX = 9 IF (NNY.LT.9) NNY = 9 IF (NNX.GT.88) NNX = 88 IF (NNY.GT.88) NNY = 88 INFX = NNX - 8 INFY = NNY - 8 C C 8 STRIPS ARE 2.88 cm , A MOLIERE RADIUS IS ABOUT 0.7 cm . C ISUPX = NNX + 8 ISUPY = NNY + 8 RNSS = 0. QTOTT = 0. DO I = INFX,ISUPX IF (DEXY(1,J,I).GE.EMIN) THEN RNSS = RNSS + 1 QTOTT = QTOTT + DEXY(1,J,I) ENDIF ENDDO DO I = INFY,ISUPY IF (DEXY(2,J,I).GE.EMIN) THEN RNSS = RNSS + 1 QTOTT = QTOTT + DEXY(1,J,I) ENDIF ENDDO NCORE = RNSS * FLOAT(J) + NCORE QCORE = QTOTT * FLOAT(J) + QCORE ENDDO C NINT = 0. CALL NOINT(NIN) ! if NINT=1 not interacting particle NINT = FLOAT(NIN) C C QCYL = 0. NCYL = 0. C C QCYL = DETECTED ENERGY AND NCYL = NUMBER OF HIT STRIPS IN A CYLINDER oF C RADIUS 8.5 STRIPS WITH AXIS DEFINED BY THE DIRECTION OF THE INCOMING C PARTICLE . C DO J = 1,NPLA C NNX = IBAR(1,J) NNY = IBAR(2,J) IF (NNX.LT.9) NNX = 9 IF (NNY.LT.9) NNY = 9 IF (NNX.GT.88) NNX = 88 IF (NNY.GT.88) NNY = 88 INFX = NNX - 8 INFY = NNY - 8 C C 8 STRIPS ARE 2.88 cm , A MOLIERE RADIUS IS ABOUT 0.7 cm . C ISUPX = NNX + 8 ISUPY = NNY + 8 DO I = INFX,ISUPX IF (DEXY(1,J,I).LT.EMIN) GO TO 710 NCYL = NCYL + 1 QCYL = QCYL + DEXY(1,J,I) 710 ENDDO DO I=INFY,ISUPY IF (DEXY(2,J,I).LT.EMIN) GO TO 810 NCYL = NCYL + 1 QCYL = QCYL + DEXY(1,J,I) 810 ENDDO ENDDO C QTRACK = 0. CALL LATERALE(QTRACK,RQT) C QMAX = 0. DO M = 1,2 DO I = 1,NPLA DO J = 1,NCHA IF (DEXY(M,I,J).GT.QMAX) QMAX = DEXY(M,I,J) ENDDO ENDDO ENDDO C if (ival.eq.1) iev = iev + 1 icount2 = icount2 + 1 if (icount2.ge.nev) goto 50 c print *,'qtot :',qtot,' nstrip :',nstrip,' ncore :',ncore, c & ' qcore :',qcore call hfnt(1) C goto 10 50 continue print *,'Processed events :',icount2 print *,'Events with no calculated baseline :',iev *** /* Save histograms */ CALL HROUT(0,icycle,' ') CALL HREND('Event') CLOSE (59) STOP END C C------------------------------------------------------------------------- SUBROUTINE CALCBASE(K,IA,IB,IPL,DX,BASE) IMPLICIT NONE INTEGER IA, IB, IPL, K REAL DX(11,96) REAL S1, S2 REAL SIG REAL C1 REAL V(16) REAL RINF, RSUP REAL BASE, BAS REAL CHECK INTEGER J, I REAL PED(4,11,96),GOOD(4,11,96) COMMON / PEDE / PED, GOOD COMMON / CH / CHECK C$$$ print *,'ia :',ia,' ib :',ib,' di :',di,' ds :',ds S1 = 0. S2 = 0. C1 = 0. DO J = IA,IB IF (GOOD(K,IPL,J).EQ.0) THEN S1 = S1 + DX(IPL,J) S2 = S2 + DX(IPL,J)*DX(IPL,J) C1 = C1 + 1. ENDIF ENDDO C BAS = 0. SIG = 0. IF (C1.GT.1) THEN BAS = S1 / C1 SIG = SQRT(ABS((S2 - C1 * BAS*BAS) /(C1 - 1.))) ENDIF C RINF = BAS - 3. * SIG RSUP = BAS + 3. * SIG C S1 = 0. S2 = 0. C1 = 0. DO J = IA,IB IF (DX(IPL,J).GT.RINF.AND.DX(IPL,J).LT.RSUP & .AND.GOOD(K,IPL,J).EQ.0) THEN S1 = S1 + DX(IPL,J) S2 = S2 + DX(IPL,J)*DX(IPL,J) C1 = C1 + 1. ENDIF ENDDO C BAS = 0. SIG = 0. IF (C1.GT.1) THEN BAS = S1 / C1 SIG = SQRT(ABS((S2 - C1 * BAS*BAS) /(C1 - 1.))) ENDIF C IF (CHECK.EQ.1) PRINT *,'BAS :',BAS,' SIG :',SIG IF (SIG.LT.30) THEN BASE = BAS RETURN ENDIF C RINF = 1.E6 DO J = IA,IB IF (DX(IPL,J).GT.0.AND.DX(IPL,J).LT.RINF.AND. & GOOD(K,IPL,J).EQ.0) THEN RINF = DX(IPL,J) I = J IF (CHECK.EQ.1) PRINT *,'RINF :',RINF,' I :',I ENDIF ENDDO C RINF = 1.E6 DO J = IA,IB IF (DX(IPL,J).GT.0.AND.DX(IPL,J).LT.RINF.AND. & GOOD(K,IPL,J).EQ.0.AND.J.NE.I) THEN RINF = DX(IPL,J) if (check.eq.1) print *,'rinf2 :',rinf ENDIF ENDDO C S1 = 0. C1 = 0. DO J = IA,IB IF (DX(IPL,J).LT.(RINF+40.).AND.DX(IPL,J).GT.(RINF-40.) & .AND.GOOD(K,IPL,J).EQ.0) THEN S1 = S1 + DX(IPL,J) C1 = C1 + 1. if (check.eq.1) print *,'s1 :',s1,' c1 :',c1 ENDIF ENDDO IF (C1.GT.1) BASE = S1 / C1 if (check.eq.1) print *,'BASE :',BASE C PRINT *,'BASE :',BASE,' C1 :',C1 RETURN END C C--------------------------------------------------------------------- SUBROUTINE LATERALE(RQT1,RQT2) C--------------------------------------------------------------------- C RQT1 (IT WILL BE CALLED QTRACK IN THE N-TUPLE) IS THE SUM OF THE DETECTED C ENERGY IN THE STRIP ALONG THE TRACK AND THE TWO CLOSEST STRIPS . FOR ALL THE C LAYERS . RQT2 (IS NOT USED IN THE N-TUPLA) IS THE TOTAL ENERGY MINUS RQT1 . C INCLUDE 'INTEST.TXT' REAL RQT1 INTEGER A,B REAL BAR(2,NPLA) REAL Q(0:NPLA) INTEGER IBAR(2,NPLA) COMMON/ANGOLO/BAR,IBAR RQT2=0. INPIA = 1 C QQQ=0 MAX=0 Q(MAX)=0 C DO I = INPIA,NPLA A = IBAR(1,I) B = IBAR(2,I) IF (A.LE.2) A = 3 IF (B.LE.2) B = 3 IF (A.GE.(NCHA-1)) A = NCHA - 2 IF (B.GE.(NCHA-1)) B = NCHA - 2 DO J = A-1,A+1 IF (DEXY(1,I,J).GE.EMIN) RQT1 = RQT1 + DEXY(1,I,J) 600 ENDDO C DO J = B-1,B+1 IF (DEXY(2,I,J).GE.EMIN) RQT1 = RQT1 + DEXY(2,I,J) ENDDO C DO J=1,A-2 PXY = DEXY(1,I,J) IF (PXY.GE.EMIN) RQT2 = RQT2 + PXY 650 ENDDO C DO J=A+2,NCHA PXY = DEXY(1,I,J) IF (PXY.GE.EMIN) RQT2 = RQT2 + PXY 700 ENDDO C DO J=1,B-2 PXY = DEXY(2,I,J) IF (PXY.GE.EMIN) RQT2 = RQT2 + PXY 750 ENDDO C DO J=B+2,NCHA PXY = DEXY(2,I,J) IF (PXY.GE.EMIN) RQT2 = RQT2 + PXY 800 ENDDO C ENDDO C C 400 RETURN END C------------------------------------------------ SUBROUTINE PEDESTAL(ERROR) C------------------------------------------------ IMPLICIT NONE C CHARACTER*6 file_name C integer icount,icount2 c integer*4 evnum, systime, len integer nev integer RUNERROR, ERROR C C Normal variables definition C INTEGER FFD C integer event C INTEGER l, j, kk, m, mm, ll, ival, iev INTEGER i, istat, icycle, ierr, icheck, icheck2, icheck3 INTEGER LVMIN integer nset,nword(10) C INTEGER*2 STRINGA(20000) C INTEGER ic, k, jj INTEGER status INTEGER inf, isup, idet INTEGER*2 length INTEGER lundata, iosop INTEGER*2 st, st2 CHARACTER*2 nr CHARACTER*1 alfa REAL auto(7) REAL RP,RI REAL SUM(6), SUM2(6), COUNT(6), & BAS(4,11), DIS(6), DDIS(6) REAL SSUM, SSUM2, CCOUNT, MME, SSIG REAL VAL REAL DI, DS REAL ME(6), SIG(6) REAL PED(4,11,96), GOOD(4,11,96) COMMON / PEDE / PED, GOOD C C Begin ! C RUNERROR=0 ! error variable PRINT *,'Pedestal file to read (only day and number included)? ' READ(*,905)file_name 905 FORMAT(A6) lundata=44 OPEN(UNIT=lundata,FILE='/wizard3/pamela/tb-0903/data/output_0309' + //file_name//'.dat',STATUS='OLD', + FORM='UNFORMATTED',ERR=50,IOSTAT=IOSOP) C ffd=FNum(lundata) ! take the file descriptor number C C Read event: use C readevent routine C 10 continue C ival = 0 C call readev(evnum,systime,len,stringa,RUNERROR,ffd) if (runerror.eq.-1.or.runerror.eq.1) goto 50 ic = 0 ic = ic + 1 length = ic do kk = 1,4 ic = length + 1 k = ic write (*,41) stringa(k) st2 = IAND(stringa(k),'FF00'x) status = ISHFT(st2,-8) write (*,41) status status = ISHFT(st2,-8) IF (KK.eq.1.and.status.ne.170) then print *,'Problems with view:',KK ERROR = 1 GOTO 50 endif IF (KK.eq.2.and.status.ne.182) then print *,'Problems with view:',KK ERROR = 1 GOTO 50 endif IF (KK.eq.3.and.status.ne.177) then print *,'Problems with view:',KK ERROR = 1 GOTO 50 endif IF (KK.eq.4.and.status.ne.173) then print *,'Problems with view:',KK ERROR = 1 GOTO 50 endif C 41 format(1x,'status :',1x,Z4) C ic = ic + 1 length = length + (stringa(ic) + 2) print *,'Lenght1 :',length print *,'Lenght2 :',stringa(ic) if (stringa(ic).ne.4629) then print *,'problems with view',KK ERROR = 1 GOTO 50 endif C do i = 1,11 do j = 1,96 ic = ic + 1 ped(kk,i,j) = stringa(ic) good(kk,i,j) = stringa(ic+1) if (kk.eq.1.and.i.eq.4.and.j.eq.77) & print *,'ped :',stringa(ic),' g/b :', % stringa(ic+1) ic = ic + 1 enddo enddo C enddo C 50 continue CLOSE (44) RETURN END