--- DarthVader/CalorimeterLevel2/src/direction.for 2006/05/19 13:15:50 1.1.1.1 +++ DarthVader/CalorimeterLevel2/src/direction.for 2007/01/22 09:17:01 1.2 @@ -19,18 +19,21 @@ INCLUDE 'INTEST.TXT' + REAL FLIMIT INTEGER NVAR,NPAR PARAMETER (NVAR=1,NPAR=2) - + PARAMETER (FLIMIT=10000) REAL X(NPLA),Y(NPLA),W(NPLA) REAL DEVIA(2),TG(2) REAL BAR(2,NPLA) REAL A, B, VAR + REAL SA, SB, SVAR + REAL SSA, SSB, SSVAR REAL PIANO(22) REAL VARFIT(2) - INTEGER ICONT, M, I - INTEGER NPFIT(2) + INTEGER ICONT, M, I, SNPFIT,SSNPFIT + INTEGER NPFIT(2), CHTRACK,IWPL(2) INTEGER IBAR(2,NPLA) COMMON/ANGOLO/BAR,IBAR @@ -40,7 +43,7 @@ COMMON/WHERE/CX,CY,PIANO SAVE /WHERE/ - COMMON/CALOFIT/VARFIT,NPFIT + COMMON/CALOFIT/VARFIT,NPFIT,IWPL,CHTRACK SAVE/CALOFIT/ C ICONT = 0 @@ -59,6 +62,11 @@ C HIGHEST DETECTED ENERGY IS USED . C DO M = 1,2 + SA = 0 + SB = 0 + SVAR = 0 + FMODE(M) = 0 + 10 CONTINUE ICONT = 0 CALL VZERO(X,NPLA) CALL VZERO(Y,NPLA) @@ -70,27 +78,118 @@ ELSE X(I) = PIANO(I) ENDIF -C X(I) = -(I-1) * PIANO Y(I) = 0. - IF (NCL(M,I).GE.1) THEN - ICONT = ICONT + 1 - Y(I) = CLUS(M,I,1) - W(I) = ((CLUS(M,I,1+NCHA/2)**0.79))**2. - IF (W(I).GT.2000.) W(I) = 2000. +C +C +C + IF (FMODE(M).EQ.0) THEN + IF (NCL(M,I).GE.1) THEN + ICONT = ICONT + 1 + Y(I) = CLUS(M,I,1) + W(I) = ((CLUS(M,I,1+NCHA/2)**0.79))**2. + IF (W(I).GT.2000.) THEN + W(I) = 2000. +C Y(I) = 0. +C ICONT = ICONT -1 + ENDIF + ENDIF + ENDIF +C +C More importance to weigth and less cluster per plane +C + IF(FMODE(M).EQ.1) THEN + IF (NCL(M,I).GE.1.AND.NCL(M,I).LE.2) THEN + ICONT = ICONT + 1 + Y(I) = CLUS(M,I,1) + W(I) = ((CLUS(M,I,1+NCHA/2)**0.29))**2. + ENDIF + IF (W(I).GT.2000.) THEN + W(I) = 2000. +c Y(I) = 0. +c ICONT = ICONT -1 + ENDIF + ENDIF +c +c One cluster per plane +c + IF(FMODE(M).EQ.2) THEN + IF (NCL(M,I).EQ.1) THEN + ICONT = ICONT + 1 + Y(I) = CLUS(M,I,1) + W(I) = ((CLUS(M,I,1+NCHA/2)**0.79))**2. + ENDIF ENDIF - if (Y(I).eq.0.or.W(I).eq.2000.) NPFIT(M) = NPFIT(M) - 1 +c if (Y(I).eq.0.or.W(I).eq.2000.) NPFIT(M) = NPFIT(M) - 1 + if (Y(I).eq.0.) NPFIT(M) = NPFIT(M) - 1 ENDDO C A = 0. B = 0. - IF (ICONT.GE.2) CALL LFITW(X,Y,W,NPLA,0,A,B,VAR) + VAR = 0. +c print *,' icont = ',icont + IF (ICONT.GE.2) THEN + CALL LFITW(X,Y,W,NPLA,0,A,B,VAR) + ELSE + NPFIT(M) = 0 + ENDIF +C +C TRY A BETTER FIT IF VAR IS BIG ON THE FIRST TRY +C + IF ((ABS(VAR).GT.FLIMIT.OR.VAR.EQ.0.).AND.FMODE(M).EQ.0) THEN + SVAR = VAR + SB = B + SA = A + SNPFIT = NPFIT(M) + NPFIT(M) = NPLA + FMODE(M) = 1 + GOTO 10 + ENDIF +C +C TRY AGAIN IF STILL BAD AND WE HAVE AT LEAST 5 POINT TO MAKE THE FIT +C + IF ((ABS(VAR).GT.FLIMIT.OR.VAR.EQ.0.).AND.FMODE(M).EQ.1 + & .AND.NPFIT(M).GT.4) THEN + SSVAR = VAR + SSB = B + SSA = A + SSNPFIT = NPFIT(M) + NPFIT(M) = NPLA + FMODE(M) = 2 + GOTO 10 + ENDIF +C +C IF WE ARE AT THE THIRD TRY AND THE FIT IS WORST THAN THE SECOND TRY TAKE THE FIRST TRY +C + IF (FMODE(M).EQ.2) THEN + IF (ABS(VAR).GT.ABS(SSVAR).OR.ABS(VAR).EQ.0. + & .OR.NPFIT(M).LE.4) THEN + FMODE(M) = 1 + VAR=SSVAR + B = SSB + A = SSA + NPFIT(M) = SSNPFIT + ENDIF + ENDIF +C +C IF THE SECOND TRY AND THE FIT IS WORST THAN THE FIRST TRY TAKE THE FIRST TRY +C + IF (FMODE(M).EQ.1) THEN + IF (ABS(VAR).GT.ABS(SVAR).OR.ABS(VAR).EQ.0.) THEN + FMODE(M) = 0 + VAR=SVAR + B = SB + A = SA + NPFIT(M) = SNPFIT + ENDIF + ENDIF +C VARFIT(M) = VAR IF (M.EQ.1) CX = B IF (M.EQ.2) CY = B -C DEVIA(M) = A TG(M) = A +C ENDDO C