Create a gist now

Instantly share code, notes, and snippets.

What would you like to do?
set terminal pngcairo font "arial,8" size 640,350
set output 'simple_distance.png'
set multiplot layout 1, 2
r(d) = d*pi/180
dist_hav(x,y,x0,y0) = 2*asin(sqrt(sin(r(y-y0)/2)**2 + cos(r(y0))*cos(r(y))*sin(r(x-x0)/2)**2))*6371
dist_sim(x,y,x0,y0) = sqrt(r(y-y0)**2 + (cos(r(y0))*r(x-x0))**2)*6371
error(x,y,x0,y0) = 100*abs(dist_hav(x,y,x0,y0)-dist_sim(x,y,x0,y0))/dist_hav(x,y,x0,y0)
set view map
set palette cubehelix negative
set contour
set isosamples 200, 200
set cntrparam level discrete 5
set title "Relative error at 89°N"
set xrange [-.05/cos(r(89)): .05/cos(r(89))]
set yrange [89-.05:89+.05]
splot error(x,y,0,89) notitle with pm3d, \
dist_hav(x,y,0,89) nosurface notitle
set title "Relative error at 65°N"
set xrange [-.05/cos(r(65)):.05/cos(r(65))]
set yrange [64.95:65.05]
splot error(x,y,0,65) with pm3d notitle, \
dist_hav(x,y,0,65) nosurface notitle
unset multiplot
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment