anchor.diff:         end if
   c      maximum radius
          if(rmax.lt.0) then
  -         rmax = 0.5*MINVAL(latv)
  +         message = "please provide rmax"
  +         goto 911
          end if
   c      number of bins
          if(nb.lt.0) then
  
positive.code:      >             real(k)*dr,gr(k)
         end do
         stop
  
  911    continue       
  c      if you reach this you are in trouble
         write(*,*) "Emergency STOP"
         write(*,*) message
  
negative.code: c      -----------------
  30     FORMAT(180A1)
  20     FORMAT(A180)
         end
  
         subroutine compute_default_rmax(rmax,lat)
         implicit none
         double precision rmax,lat(3,3)
  c       
         double precision d(3)
  
         d(1) = sqrt(SUM(lat(:,1)*lat(:,1)))
         d(2) = sqrt(SUM(lat(:,2)*lat(:,2)))
         d(3) = sqrt(SUM(lat(:,3)*lat(:,3)))
         rmax = 0.5*MINVAL(d)
         return
         end
  
         subroutine xyz_read_natoms(n,fn_in)
  