anchor.diff:         allocate(c1(n,3))
          allocate(c2(n,3))
          allocate(gr(nb))
  +       allocate(gr0(nb))
          write(*,*) "number of atoms", n
   c       call xyz_read_coords_byname(n,c,atom1_tag,fn_in)
  +c      -----------------------------
  +c      start looping over the frames
  +c      -----------------------------
  +       nf = 0
  +       gr = 0
   30     continue
  +       do
          n1=n
          n2=n
          call xyz_read_coords_byname2(n1,c1,atom1_tag,
        +                              n2,c2,atom2_tag,
        +                              fn_in)
  -       if(n1.eq.0.or.n2.eq.0) goto 31
  +       if(n1.eq.0.or.n2.eq.0) then
  +          write(*,*) "exiting because there are no coordinates",nf
  +          exit
  +       end if
          do i=1,n1
            write(34,*) (c1(i,k),k=1,3)
          end do
          do i=1,n2
            write(35,*) (c2(i,k),k=1,3)
          end do
  -c
  +
  +       gr0 = 0
          call rdf_compute(n, n1,c1,atom1_tag,
        +                  n2,c2,atom2_tag,
        +                  nb,rmax,lat,
  -     +                  gr,
  +     +                  gr0,
        +                  fn_in)
  -       goto 30
  +       gr = gr + gr0
  +       nf = nf+1
  +       if(nfrm.gt.0.and.nf.gt.nfrm) then
  +          write(*,*) "exiting because exceeded number of frames",nf,nfrm
  +          exit
  +       end if
  +       end do
  +c       goto 30
  +c      ---------------------------
  +c      end looping over the frames
  +c      ---------------------------
   31     continue
  -       write(*,*) "came to the end of the file"
  +c      average RDF over the frames
  +       gr = gr/real(nf)
  +       write(*,*) "came to the end of the file",nf
  +       dr = rmax/nb
  +       do k=1,nb
  +         write(fn_out,'(2F12.6)') 
  +     >             real(k)*dr,gr(k)
  +       end do
          stop
   
          call pair_compute(rmax,lat,i1,atag,nb,infile)
  
positive.code:        character*30 pname 
  c
         pname = "xyz_read_coords_byname"
         write(*,*) "in "//pname
         j1 = 0
         j2 = 0
  
negative.code: c       ------------------------------------
  c       read coordinates(skipping empty lines)
  c       ------------------------------------
         write(*,*) "reading coordinates"
         do i=1,n
          read(fn_in,20) bigbuf
          l=l+1
  