Air mass: Difference between revisions

=={{header|BASIC}}==

==={{header|BASIC256}}===

{{trans|FreeBASIC}}

<lang freebasic>global RE, dd, LIM

RE = 6371000 #Earth radius in meters

dd = 0.001 #integrate in this fraction of the distance already covered

LIM = 10000000 #integrate only to a height of 10000km, effectively inLIMity

print "Angle 0 m 13700 m"

print "------------------------------------"

for z = 0 to 90 step 5

print rjust(z,2); " "; ljust(airmass(0, z),13,"0"); " "; ljust(airmass(13700, z),13,"0")

next z

end

function max(a, b)

if a > b then return a else return b

end function

function rho(a)

#the density of air as a function of height above sea level

return exp(-a/8500.0)

end function

function height(a, z, d)

#a = altitude of observer

#z = zenith angle (in degrees)

#d = distance along line of sight

AA = RE + a

HH = sqr(AA^2 + d^2 - 2*d*AA*cos(radians(180-z)))

return HH - RE

end function

function column_density(a, z)

#integrates density along the line of sight

sum = 0.0

d = 0.0

while d < LIM

delta = max(dd, (dd)*d) #adaptive step size to avoid it taking forever:

sum += rho(height(a, z, d+0.5*delta)) * delta

d += delta

end while

return sum

end function

function airmass(a, z)

return column_density(a, z) / column_density(a, 0)

end function</lang>

==={{header|FreeBASIC}}===

<lang freebasic>

#define DEG 0.017453292519943295769236907684886127134 'degrees to radians

#define RE 6371000 'Earth radius in meters

#define dd 0.001 'integrate in this fraction of the distance already covered

#define FIN 10000000 'integrate only to a height of 10000km, effectively infinity

#define max(a, b) iif(a>b,a,b)

function rho(a as double) as double

'the density of air as a function of height above sea level

return exp(-a/8500.0)

end function

function height( a as double, z as double, d as double ) as double

'a = altitude of observer

'z = zenith angle (in degrees)

'd = distance along line of sight

dim as double AA = RE + a, HH

HH = sqr( AA^2 + d^2 - 2*d*AA*cos((180-z)*DEG) )

return HH - RE

end function

function column_density( a as double, z as double ) as double

'integrates density along the line of sight

dim as double sum = 0.0, d = 0.0, delta

while d<FIN

delta = max(dd, (dd)*d) 'adaptive step size to avoid it taking forever:

sum += rho(height(a, z, d+0.5*delta))*delta

d += delta

wend

return sum

end function

function airmass( a as double, z as double ) as double

return column_density( a, z ) / column_density( a, 0 )

end function

print "Angle 0 m 13700 m"

print "------------------------------------"

for z as double = 0 to 90 step 5.0

print using "## ##.######## ##.########";z;airmass(0, z);airmass(13700, z)

next z

</lang>

{{out}}

<pre>

Angle 0 m 13700 m

------------------------------------

0 1.00000000 1.00000000

5 1.00380963 1.00380965

10 1.01538466 1.01538475

15 1.03517744 1.03517765

20 1.06399053 1.06399093

25 1.10305937 1.10306005

30 1.15418974 1.15419083

35 1.21998076 1.21998246

40 1.30418931 1.30419190

45 1.41234169 1.41234567

50 1.55280404 1.55281025

55 1.73875921 1.73876915

60 1.99212000 1.99213665

65 2.35199740 2.35202722

70 2.89531368 2.89537287

75 3.79582352 3.79596149

80 5.53885809 5.53928113

85 10.07896219 10.08115981

90 34.32981136 34.36666557

</pre>

==={{header|True BASIC}}===

{{trans|FreeBASIC}}

<lang qbasic>FUNCTION max(a, b)

IF a > b then LET max = a else LET max = b

END FUNCTION

FUNCTION rho(a)

!the density of air as a function of height above sea level

LET rho = exp(-a/8500)

END FUNCTION

FUNCTION height(a, z, d)

!a = altitude of observer

!z = zenith angle (in degrees)

!d = distance along line of sight

LET aa = re+a

LET hh = sqr(aa^2+d^2-2*d*aa*cos((180-z)*deg))

LET height = hh-re

END FUNCTION

FUNCTION columndensity(a, z)

!integrates density along the line of sight

LET sum = 0

LET d = 0

DO while d < lim

LET delta = max(dd, (dd)*d) !adaptive step size to avoid it taking forever:

LET sum = sum+rho(height(a, z, d+.5*delta))*delta

LET d = d+delta

LOOP

LET columndensity = sum

END FUNCTION

FUNCTION airmass(a, z)

LET airmass = columndensity(a, z)/columndensity(a, 0)

END FUNCTION

LET deg = .0174532925199433 !degrees to radians

LET re = 6371000 !Earth radius in meters

LET dd = .001 !integrate in this fraction of the distance already covered

LET lim = 10000000 !integrate only to a height of 10000km, effectively infinity

PRINT "Angle 0 m 13700 m"

PRINT "------------------------------------"

FOR z = 0 to 90 step 5

PRINT using "## ##.######## ##.########": z, airmass(0, z), airmass(13700, z)

NEXT z

END</lang>