      program basin_volume
      implicit none
      integer natoms,attract_label
      integer, dimension (:), allocatable :: nat
      double precision, dimension (:), allocatable :: xat,yat,zat
      double precision xmin,xmax,ymin,ymax,zmin,zmax,dv
      double precision xminold,xmaxold,yminold,ymaxold,zminold,zmaxold
      double precision trans(3,3),dx,dy,dz
      double precision, dimension (:), allocatable :: volume 
      double precision, dimension (:,:), allocatable :: vol_bas 
      double precision xc, yc, zc, dist
      double precision v(64,3)
      character*16, dimension (:), allocatable :: basin_name
      character*4, dimension (:), allocatable :: atom_name
      character*40 file_in
      character*1 molecule(40),blank
      character*3 file_type
      integer, dimension (:,:), allocatable :: index_basin,aux_label
      integer indexi,indexk
      integer npx,npy,npz,nop,sym(3,3,3)
      integer npxold,npyold,npzold
      integer i, j, k, l, m, n, maxv, iop, np,npts
      integer ix1,ix2,iy1,iy2,iz1,iz2,ix,iy,iz,ijk
      logical if_file
      integer*2, dimension (:), allocatable :: iaux
      real, dimension (:), allocatable :: aux
      real r(8,3)
      double precision, dimension (:), allocatable :: x,y,z
      character*3 name_op(3)
      common /mat/xc,yc,zc,trans,sym,name_op
      equivalence (molecule(1),file_in)
      data blank/' '/
      write(*,*) ' input  ebas_sbf file '
      read(*,'(a)') file_in
      inquire(file=file_in,exist=if_file)
      if(.not.if_file) then
        write(*,'(2x,"requested file does not exist:",a40)') file_in
        stop
      endif
      open(unit=1,file=file_in,status='old',form='unformatted')
      read(1) npx,npy,npz
      read(1) xmin,xmax,ymin,ymax,zmin,zmax
      npts=npx*npy*npz
      dv=(xmax-xmin)*(ymax-ymin)*(zmax-zmin)/npts
      read(1) natoms,attract_label
      allocate(xat(natoms),yat(natoms),zat(natoms),nat(natoms))
      allocate(atom_name(natoms))
      allocate(volume(attract_label),basin_name(attract_label))
      allocate(index_basin(attract_label,5),aux_label(attract_label,64))
      allocate(x(attract_label),y(attract_label),z(attract_label))
      allocate(vol_bas(attract_label,4))
      do i=1,natoms
         read(1) nat(i),atom_name(i),xat(i),yat(i),zat(i)
      enddo
      do i=1,attract_label
         read(1) index_basin(i,1),basin_name(i),volume(i),
     .           x(i),y(i),z(i)
      enddo
      allocate(iaux(npts))
      read(1)(iaux(i),i=1,npts)
      close(unit=1)
      write(*,*) ' input  rho_sbf file '
      read(*,'(a)') file_in
      inquire(file=file_in,exist=if_file)
      if(.not.if_file) then
        write(*,'(2x,"requested file does not exist:",a40)') file_in
        stop
      endif
      open(unit=1,file=file_in,status='old',form='unformatted')
      read(1) npxold,npyold,npzold
      read(1) xminold,xmaxold,yminold,ymaxold,zminold,zmaxold
      allocate(aux(npts))
      read(1)(aux(i),i=1,npts)
      close(unit=1)
      do i=1,attract_label
         do j=1,4
            vol_bas(i,j)=0.0d0
         enddo
      enddo
      do i=1,npts
         k=iaux(i)
         if(k.gt.0) then
         if(aux(i).ge.1.0d-2) then
           vol_bas(k,1)=vol_bas(k,1)+dv
         endif
         if(aux(i).ge.1.0d-3) then
           vol_bas(k,2)=vol_bas(k,2)+dv
         endif
         if(aux(i).ge.1.0d-4) then
           vol_bas(k,3)=vol_bas(k,3)+dv
         endif
         if(aux(i).ge.1.0d-5) then
           vol_bas(k,4)=vol_bas(k,4)+dv
         endif
         endif
      enddo
      do i=1,attract_label
         write(*,'(2x,a16,4f12.3)') basin_name(i),(vol_bas(i,j),j=1,4)
      enddo
      end
