function d = gc_dist(lon, lat)
%Returns the cumulative great circle distance from a series of ([lon] [lat])
%points. Units will be in units of R (currently km)

R = 6378.1;

if(size(lon) ~= size(lat))
    error('lon and lat arrays must be the same size');
end

lon = lon*pi/180;
lat = lat*pi/180;

d = zeros(size(lon));

for(i=2:length(d))
    dlon = lon(i) - lon(i-1);
    dlat = lat(i) - lat(i-1);
    a = (sin(dlat/2)).^2 + cos(lat(i-1)) * cos(lat(i)) * (sin(dlon/2))^2;
    c = 2 * asin(min(1,sqrt(a)));
    d(i) = d(i-1) + (R * c);
end