function [ret] = geodist ( LAT1, LONG1, LAT2, LONG2 )
  LAT1r=LAT1*pi/180;
  LONG1r=LONG1*pi/180;
  LAT2r=LAT2*pi/180;
  LONG2r=LONG2*pi/180;
  x=acos( sin(LAT1r)*sin(LAT2r)+cos(LAT1r)*cos(LAT2r)*cos(LONG1r-LONG2r))*3959; 
  ret=x;
endfunction
