      subroutine potts_zln(nlink,namin,beta0,b,ha,hasum,iopt,
     &                     zln,aln,a2ln) ! bernd berg, may 14 2002.
c potts model reweighting with jackknife errors.
c ln of the partition function z and related variables.
      include '../../libs/fortran/implicit.sta'
      include '../../libs/fortran/constants.par'
      dimension b(0:nlink),ha(0:nlink),hasum(0:nlink)
      data a/zero/,zln1_max/zero/,iact_old/0/ ! avoids compiler warning.
      ltest=.true.
      ltest=.false.
      if(ltest) print*,"potts_zln: beta0 =",beta0
c
      iloop_ha=0 ! counts runs through the ha(iact).gt.half loop.
      if(namin.le.0) stop "potts_zln: namin false."
      do iact=0,nlink
        ha(iact)=hasum(iact)-ha(iact) ! jackknife histogram.
        if(ha(iact).gt.half) then
          iloop_ha=iloop_ha+1
          if(iloop_ha.eq.1) then
            a=-(beta0-b(iact))*namin ! ferdinand-fisher normalization.
            zln=log(ha(iact))+two*((beta0-b(iact))*iact+a)
            zln1_max=zln
            aln=zln+log(iact*one)
            a2ln=zln+two*log(iact*one)
            if(ltest) print*,"iact,a,aln,zln =",iact,a,aln,zln
          else
            a=a+(b(iact)-b(iact_old))*iact_old
            zln1=log(ha(iact))+two*((beta0-b(iact))*iact+a)
            zln1_max=max(zln1_max,zln1)
            aln1=zln1+log(iact*one)
            a2ln1=zln1+two*log(iact*one)
            zln=addln(zln,zln1)
            aln=addln(aln,aln1)
            a2ln=addln(a2ln,a2ln1)
          end if
          iact_old=iact
        end if
      end do
c
      if(iopt.ne.1) return
c
c iact probability density at beta0: (zln1_max used only here.)
      iloop_ha=0
      hsum=zero
      do iact=0,nlink
        if(ha(iact).gt.half) then
          iloop_ha=iloop_ha+1
          if(iloop_ha.eq.1) then
            a=-(beta0-b(iact))*namin
          else
            a=a+(b(iact)-b(iact_old))*iact_old
          end if
          zln1=log(ha(iact))+two*((beta0-b(iact))*iact+a)-zln1_max
          if((zln1/zln1_max).gt.one+em12) stop "potts_zln: zln1." 
          ha(iact)=exp(zln1) ! maximum value at this point one.
          hsum=hsum+ha(iact)
          iact_old=iact
        else
          ha(iact)=zero
        end if
      end do
      do iact=0,nlink ! normalization:
        ha(iact)=ha(iact)/hsum
      end do
c
      return
      end
