*-- Author : Stephen Burke 07/05/92
SUBROUTINE FVZWM(INFTKR,NFTKR,ZNOM,LFIRST,FVVEC,IERR)
*-----------------------------------------Updates 26/07/93-------
**: FVZWM 30907 RP. Farm changes.
*-----------------------------------------Updates 30/10/92-------
**: FVZWM 30907 SB. New debug histogram numbers.
*-----------------------------------------Updates 06/05/92-------
**: FVZWM 30907 SB. New deck to take a weighted mean of z values.
*-----------------------------------------Updates----------------
**********************************************************************
* *
* Calculate z-vertex by weighted mean *
* *
* INPUT; *
* INFTKR - FTKR bank index *
* NFTKR - The number of FTKR rows (= size of work bank) *
* ZNOM - the nominal z-vertex position *
* LFIRST - TRUE for the first call (primary vertex) *
* *
* OUTPUT; *
* FVVEC - the four words of the FTGR bank *
* IERR - non-zero if the weighted mean fails *
* *
**********************************************************************
DIMENSION FVVEC(4)
LOGICAL LFIRST,LRJCT,LPRIM
*KEEP,FVSTEE.
LOGICAL LTRUTH,LCUT,LRESID
COMMON /FVSTEE/ IDIAG,LUN,LUNHB,LTRUTH,LCUT,LRESID
*KEEP,FVPAR.
DOUBLE PRECISION ZWALL1,ZWALL2,RADLEN
COMMON /FVPAR/ ZWALL1,ZWALL2,RADLEN,MINHTP,MINHTR,ZSQMAX
&, PMIN,DCAMAX,Z0MAX,CHIMAX
*KEEP,FVSCAL.
* Various counters
PARAMETER (NSCAL=16)
COMMON /FVSCAL/ NNEVNT,NNVTX,NNFTKR,NNXTR,NNFIT,NNOUT,NNSIN
&, NNFTKP,NNXTRP,NNFITP,NNOUTP,NNSINP
&, NNVTXC,NNSINC,NNFVNC,NNFSNC
*KEEP,FVWBI.
* Work bank indices
PARAMETER (NFVWBI=2)
COMMON /FVWBI/ INFTPR,INFVWK
*KEEP,BCS.
INTEGER NHROW,NHCOL,NHLEN
PARAMETER (NHROW = 2, NHCOL = 1, NHLEN=2)
INTEGER NBOSIW
PARAMETER (NBOSIW=1000000)
INTEGER IW(NBOSIW)
REAL RW(NBOSIW)
COMMON /BCS/ IW
EQUIVALENCE (RW(1),IW(1))
SAVE /BCS/
*KEEP,STFUNCT.
* index of element before row number IROW
INDR(IND,IROW)=IND+2+IW(IND+1)*(IROW-1)
* index of L'th element of row number IROW
INDCR(IND,L,IROW)=INDR(IND,IROW) + L
* L'th integer element of the IROW'th row of bank with index IND
IBTAB(IND,L,IROW)=IW(INDCR(IND,L,IROW))
* L'th real element of the IROW'th row of bank with index IND
RBTAB(IND,L,IROW)=RW(INDCR(IND,L,IROW))
*KEND.
**********************************************************************
* Default is failure
IERR = 1
*
* Take a weighted mean of the z0 values
*
SIGZW = 0.
SIGW = 0.
NDF = -1
DO 200 JFT=1,NFTKR-1,2
WZ0 = RW(INFVWK+JFT+1)
* Rejected tracks have weight zero
IF (WZ0.LE.0.) GOTO 200
Z0 = RW(INFVWK+JFT)
IF (SIGW.GT.0.) THEN
ZMEAN = SIGZW/SIGW
CHI = (Z0 - ZMEAN)**2*WZ0
ELSE
CHI = 0.
ENDIF
IF (CHI.LE.10*CHIMAX) THEN
* Accept new value if reasonably consistent with current estimate
SIGZW = SIGZW + Z0*WZ0
SIGW = SIGW + WZ0
NDF = NDF + 1
* Mark value used by setting weight negative
RW(INFVWK+JFT+1) = -WZ0
ELSEIF (ABS(Z0-ZNOM).LT.ABS(ZMEAN-ZNOM)) THEN
* New value is closer to nominal z, so discard previous estimate
SIGZW = Z0*WZ0
SIGW = WZ0
NDF = 0
* Mark used by setting weight negative
RW(INFVWK+JFT+1) = -WZ0
* Reset used flag for previous tracks
DO 100 JJFT=1,JFT-2,2
WZ0 = RW(INFVWK+JJFT+1)
IF (WZ0.LT.0.) RW(INFVWK+JJFT+1) = -WZ0
100 CONTINUE
ENDIF
200 CONTINUE
* First guess at z0 and error
IF (SIGW.LE.0.) RETURN
CZMEAN = 1./SIGW
ZMEAN = SIGZW*CZMEAN
300 CONTINUE
*
* Now work out a chi-squared, and throw away any tracks
* which are too far from the mean
*
CHISQ = 0.
NOUTP = 0
LRJCT = .FALSE.
DO 400 JFT=1,NFTKR-1,2
WZ0 = RW(INFVWK+JFT+1)
* Used tracks now have negative weight
IF (WZ0.GE.0.) GOTO 400
WZ0 = -WZ0
Z0 = RW(INFVWK+JFT)
* Primary flag for diagnostics
LPRIM = IBTAB(INFTKR,7,JFT).EQ.0
IF (LPRIM) NOUTP = NOUTP + 1
IF (NDF.EQ.0) THEN
* If there's only one track, accept it
IF (LFIRST) THEN
NNSIN = NNSIN + 1
IF (LPRIM) NNSINP = NNSINP + 1
ENDIF
CHISQ = 0.
ELSE
CHI = (Z0 - ZMEAN)**2*WZ0
IF (CHI.LE.CHIMAX) THEN
* Accept if close enough to the mean
CHISQ = CHISQ + CHI
ELSE
* Remove this track (but don't yet update ZMEAN)
SIGZW = SIGZW - Z0*WZ0
SIGW = SIGW - WZ0
NDF = NDF - 1
* Reset used flag
RW(INFVWK+JFT+1) = WZ0
LRJCT = .TRUE.
ENDIF
IF (LCUT .AND. LPRIM) THEN
CALL HFILL(213,CHI,0.,1.)
ELSEIF (LCUT) THEN
CALL HFILL(214,CHI,0.,1.)
ENDIF
ENDIF
400 CONTINUE
* Iterate if necessary
IF (LRJCT) THEN
IF (SIGW.LE.0.) RETURN
CZMEAN = 1./SIGW
ZMEAN = SIGZW*CZMEAN
IF (NDF.GT.0) GOTO 300
ENDIF
IF (LFIRST) THEN
NNOUTP = NNOUTP + NOUTP
* Create FTGX bank
INFTGX = NBANK('FTGX',0,2+NFTKR/2)
IF (INFTGX.LE.0) THEN
CALL ERRLOG(512,'S:FVZWM: Unable to create FTGX')
RETURN
ENDIF
IW(INFTGX+1) = 1
IW(INFTGX+2) = 0
ENDIF
* Remove all tracks used for this vertex
DO 500 JFT=1,NFTKR-1,2
WZ0 = RW(INFVWK+JFT+1)
IF (WZ0.LT.0.) THEN
RW(INFVWK+JFT+1) = 0.
IW(INFTGX+2) = IW(INFTGX+2) + 1
IW(INDCR(INFTGX,1,IW(INFTGX+2))) = JFT
ENDIF
500 CONTINUE
* Fill the output vector
FVVEC(1) = ZMEAN
FVVEC(2) = SQRT(CZMEAN)
FVVEC(3) = CHISQ
CALL UCOPY(NDF,FVVEC(4),1)
IERR = 0
IF (LCUT) THEN
PRAT = FLOAT(NOUTP)/FLOAT(NDF+1)
IF (PRAT.GT.0.999) PRAT = 0.999
CALL HFILL(300,PRAT,0.,1.)
ENDIF
RETURN
END
*