      subroutine bmha_tabfig(nbm1,nbm2,tabbm,area,ymin,ymax,
     &  iud1,iud2,nareas) ! a. bazavov and b.a. berg, may 26 2005. 
c prepare bmha areas figure.
c nbm1 -- number of steps on x (abscissa).
c nbm2 -- number of steps on y (ordinate).
c tabbm -- 2d table.
c area -- array with calculated probability areas (not sorted).
c nareas -- number of areas.
      include '../../libs/fortran/implicit.sta'
      include '../../libs/fortran/constants.par'
      dimension tabbm(nbm2,nbm1),area(nbm1*nbm2)

      open(iud1,file="bmha01.d",form='formatted',status='unknown')
      open(iud2,file="bmha02.d",form='formatted',status='unknown')
      write(iud1,'(2f12.5)') +half,ymin
      write(iud1,'(2f12.5)') +half,ymax
      pm_max=ymax
      pm_min=ymin
      xa=+half
      nareas=0
      do ibm1=1,nbm1 ! loop on abscissa.
c output vertical lines equally spaced (into first file):
        write(iud1,'(2f12.5)') ibm1+half,pm_max
        write(iud1,'(2f12.5)') ibm1+half,pm_min
c swap up and down:
        pm_one=pm_max
	pm_max=pm_min
	pm_min=pm_one
c increment x
        xb=ibm1+half
        ya=ymin
c loop on y
        do ibm2=1,nbm2-1
c read current y value from the table
          yb=tabbm(ibm2,ibm1)
          write(iud2,'(2f12.5)') xa,yb
          write(iud2,'(2f12.5)') xb,yb
c calculate current area and add to the array
          nareas=nareas+1
          area(nareas)=abs(xb-xa)*(yb-ya)
c swap left and right, up and down
          ya=yb
          vsave=xa
          xa=xb
          xb=vsave
        end do
c calculate current area and add to the array
        nareas=nareas+1
        area(nareas)=abs(xb-xa)*(ymax-ya)
        xa=ibm1+half
      end do
      close(iud2)
      close(iud1)

      return
      end
