set dummy x, y
set format x "% h"
set format y "% h"
set grid nopolar
set grid xtics nomxtics ytics nomytics noztics nomztics nortics nomrtics \
 nox2tics nomx2tics noy2tics nomy2tics nocbtics nomcbtics
set grid layerdefault   lt 0 linecolor 0 linewidth 0.500,  lt 0 linecolor 0 linewidth 0.500
unset raxis
set key notitle
set key fixed left top vertical Left reverse enhanced autotitle box lt black linewidth 1.000 dashtype solid
set key noinvert samplen 2 spacing 1 width 2 height 0
set style textbox  transparent margins  1.0,  1.0 border  lt -1 linewidth  1.0
set offsets 0, 0, 0, 0
set pointsize 1
set term svg
set pointintervalbox 1
set encoding default
set cntrlabel  format '%8.3g' font '' start 5 interval 20
set title "Viewing angle and air masses"
set xlabel "Zenit Angle (grad)"
set xrange [ 80.0000 : 90.0000 ] noreverse writeback
set xtick 1
set ylabel "Air masses (AM)"
set yrange [ 0.00000 : 40.0000 ] noreverse writeback
set zero 1e-08
set lmargin  -1
set bmargin  -1
set rmargin  -1
set tmargin  -1
sec(x) = 1/cos(x/180*pi)
yAi(x) = sec(x) * (1 - 0.0012*(sec(x)*sec(x) - 1))
hardie(x) = sec(x) - 0.0018167 * (sec(x) - 1) - 0.002875 * (sec(x) - 1)*(sec(x) - 1) - 0.0008083*(sec(x) - 1)*(sec(x) - 1)*(sec(x) - 1)
rozenberg(x) = 1/(cos(x/180*pi) + 0.025*exp(-11*cos(x/180*pi)))
pow(a, b) = exp(b*log(a))
kay(x) = 1/(cos(x/180*pi) + 0.50571*pow(6.07995 + 90 - x, -1.6364))
z(x) = x/180*pi
young(x) = (1.002432*cos(z(x))*cos(z(x)) + 0.148386*cos(z(x)) + 0.0096467)/(cos(z(x))*cos(z(x))*cos(z(x)) + 0.149864*cos(z(x))*cos(z(x)) + 0.0102963*cos(z(x)) + 0.000303978)
pickering(x) = 1/(sin(z(90 - x + 244/(165 + 47*pow(90 - x, 1.1)))))
spherical(x) = 6371.0/10.0*sqrt(cos(z(x))*cos(z(x)) + 2*10.0/6371.0 + 10.0*10.0/6371.0/6371.0) - 6371.0/10.0*cos(z(x))
plot sec(x) title "Plane Parallel model" lc "black", rozenberg(x) title "Rozenberg (1966)" lc "red", kay(x) title "Kasten and Young (1989)" lc "blue", "sph.dat" title "Buger (1729, spherical)" with lines lc "orange", young(x) title "Young (1994, without refraction)" lc "green"