Download or view MarsAtmospherePlot.frink in plain text format
/** This plots atmospheric data for Mars based on the Mars-GRAM data.
It is used to create the much better MarsAtmosphere.frink model
See the Mars-GRAM User's Guide at:
https://ntrs.nasa.gov/api/citations/20210023957/downloads/Mars-GRAM%20User%20Guide.pdf
(jump to page 62 of the PDF file which is labeled page 49 on the document itself)
It currently is modeled between the latitudes between -30 deg and 30 deg
where most landings occur.
For a map of landing spots, see:
https://www.planetary.org/space-images/mars_landing_site_map_lakdawalla
*/
use Grid.frink
use MarsAtmosphereOld.frink
use curveFit.frink
use functionUtils.frink
use LeastSquares.frink
url = filenameToURL["/home/eliasen/builds/Mars-GRAM 2010/Release1.0_Nov10/txtFiles/tpdloy11.txt"]
maxlat = 30 deg
//minheight = -5 km
//maxheight = 200 km
minheight = -5 km
maxheight = 170 km
g = new graphics
dtemp = new dict
ddensity = new dict
dpressure = new dict
dummyx =!= dummyx
AITemp[h is height] :=
{
h1 = h / km
return K (242 - 0.087 h1 + 1.15 × 10^-4 h^2)
}
for line = rest[lines[url]] // Skips first line of headers
{
line = trim[line]
fields = eval[split[%r/\s+/, line]]
long = fields@0 degrees
height = fields@1 km
lat = fields@2 deg
temp = fields@3 K
density = fields@8 kg/m^3
pressure = fields@13 Pa
// Change these to focus on a certain latitude.
// See
// https://www.planetary.org/space-images/mars_landing_site_map_lakdawalla
if abs[lat] <= maxlat and height <= maxheight and height >= minheight
{
dtemp.addToList[height, temp]
ddensity.addToList[height, density]
dpressure.addToList[height, pressure]
}
}
// Medium altitude file
url = filenameToURL["/home/eliasen/builds/Mars-GRAM 2010/Release1.0_Nov10/txtFiles/tpdms101.txt"]
for line = rest[lines[url]] // Skips first line of headers
{
line = trim[line]
fields = eval[split[%r/\s+/, line]]
long = fields@0 degrees
height = fields@1 km
lat = fields@2 deg
temp = fields@3 K
pressure = fields@8 Pa
density = fields@13 kg/m^3
// Change these to focus on a certain latitude.
// See
// https://www.planetary.org/space-images/mars_landing_site_map_lakdawalla
if abs[lat] <= maxlat and height <= maxheight and height >= minheight
{
dtemp.addToList[height, temp]
ddensity.addToList[height, density]
dpressure.addToList[height, pressure]
}
}
println["\nTemperature"]
meantemp = new dict
ptemp = new polyline
ptemp1 = new polyline
for alt = sort[dtemp.keys[]]
{
meantemp@alt = meanAndSD[dtemp@alt, false]@0
println["$alt\t" + meantemp@alt]
ptemp.addPoint[meantemp@alt/K, -alt/km]
[t,p,d] = MarsAtmosphereOld.getTPD[alt]
ptemp1.addPoint[t/K, -alt/km]
}
// Fit the values to a piecewise polynomial of the specified degree
[tf1, rtf1] = LeastSquares.fitDegreePiecewise[meantemp, 5]
println[inputForm[tf1]]
println["rValue: $rtf1"]
[tempFunc, rtemp] = quadraticFitFunctionPiecewise3[meantemp]
//tempFunc = quadraticFitFunction[meantemp]
//rtemp = rValue[tempFunc, toArray[meantemp]]
//println[inputForm[functionBody[tempFunc]]]
println[inputForm[tempFunc]]
println["rvalue = $rtemp"]
ptempfit = new polyline
ptempfit1 = new polyline
for alt = sort[meantemp.keys[]]
{
ptempfit.addPoint[tempFunc[alt]/K, -alt/km]
ptempfit1.addPoint[tf1[alt]/K, -alt/km]
}
gtemp = new graphics
gtemp.font["Monospaced", .05]
gtemp.add[ptemp]
gtemp.color[0,0,1]
gtemp.add[ptemp1]
gtemp.color[1,0,0]
gtemp.add[ptempfit]
gtemp.color[0,1,0]
gtemp.add[ptempfit1]
grid = new Grid
grid.setUnits[1, -1]
grid.auto[gtemp]
gtemp.add[grid.getGrid[]]
gtemp.caption["Temperature, K", "top"]
gtemp.caption["altitude (km)", "right", -90 deg]
gtemp.show[]
gtemp.write["MarsAtmosphereTemperature.svg", 1000, 800]
gtemp.write["MarsAtmosphereTemperature.png", 1000, 800]
println["\nPressure"]
meanpressure = new dict
ppressure = new polyline
ppressure1 = new polyline
for alt = sort[dpressure.keys[]]
{
meanpressure@alt = meanAndSD[dpressure@alt, false]@0
println["$alt\t" + meanpressure@alt]
ppressure.addPoint[log[meanpressure@alt / Pa] dummyx, -alt/km]
[t,p,d] = MarsAtmosphereOld.getTPD[alt]
ppressure1.addPoint[log[p / Pa] dummyx, -alt/km]
}
[pressureFunc, rpressure] = exponentialFitFunctionPiecewise[meanpressure]
//rpressure = rValue[pressureFunc, toArray[meanpressure]]
//println[inputForm[functionBody[pressureFunc]]]
println[inputForm[pressureFunc]]
println["rvalue = $rpressure"]
ppressurefit = new polyline
for alt = sort[meanpressure.keys[]]
ppressurefit.addPoint[log[pressureFunc[alt] / Pa] dummyx, -alt/km]
gpressure = new graphics
gpressure.add[ppressure]
gpressure.color[0,0,1]
gpressure.add[ppressure1]
gpressure.color[1,0,0]
gpressure.add[ppressurefit]
grid = new Grid
grid.setUnits[dummyx, -1]
grid.auto[gpressure]
gpressure.add[grid.getGrid[]]
gpressure.caption["Pressure", "top"]
gpressure.caption["log10[pressure/pascal]"]
gpressure.show[]
gpressure.write["MarsAtmospherePressure.svg", 1000, 800]
gpressure.write["MarsAtmospherePressure.png", 1000, 800]
println["\nDensity"]
meandensity = new dict
pdensity = new polyline
pdensity1 = new polyline
for alt = sort[ddensity.keys[]]
{
meandensity@alt = meanAndSD[ddensity@alt, false]@0
println["$alt\t" + meandensity@alt]
pdensity.addPoint[log[meandensity@alt / (kg/m^3)] dummyx, -alt/km]
[t,p,d] = MarsAtmosphereOld.getTPD[alt]
if d != 0 kg/m^3
pdensity1.addPoint[log[d / (kg/m^3)] dummyx, -alt/km]
}
[densityFunc, rdensity] = exponentialFitFunctionPiecewise[meandensity]
println[inputForm[densityFunc]]
println["rvalue = $rdensity"]
pdensityfit = new polyline
for alt = sort[meandensity.keys[]]
pdensityfit.addPoint[log[densityFunc[alt] / (kg/m^3)] dummyx, -alt/km]
gdensity = new graphics
gdensity.add[pdensity]
gdensity.color[0,0,1]
gdensity.add[pdensity1]
gdensity.color[1,0,0]
gdensity.add[pdensityfit]
grid = new Grid
grid.setUnits[dummyx, -1]
grid.auto[gdensity]
gdensity.add[grid.getGrid[]]
gdensity.caption["Density", "top"]
gdensity.caption["log10[density/(kg/m^3)]"]
gdensity.show[]
gdensity.write["MarsAtmosphereDensity.svg", 1000, 800]
gdensity.write["MarsAtmosphereDensity.png", 1000, 800]
Download or view MarsAtmospherePlot.frink in plain text format
This is a program written in the programming language Frink.
For more information, view the Frink
Documentation or see More Sample Frink Programs.
Alan Eliasen was born 20217 days, 23 hours, 47 minutes ago.