I am
fascinated by space exploration. The other day, inspired by the
recent anniversary of the Apollo 16 mission, I was watching a short documentary about Moon exploration.
The
documentary had an interesting note about the Lunar ReconnaissanceOrbiter, a space probe orbiting the Moon since 2009, and performing
different measurement to prepare for future Moon expeditions. One of
the six on-board instruments of the orbiter is called Lunar Orbiter
Laser Altimeter, or LOLA. It measures the precise height of the lunar
surface under the probe by centimeter accuracy, and it is used to
build a topographic map of the Moon.
That
has turned on my brain and I was carried away. What could be more
cool than play around with the real 3D topographical map of the Moon!
Then I figured, that NASA is a public institute, they publish
practically everything they measure or record during their missions.
So maybe, if I dig deep enough, I can find LOLA data out in the
Internet wilderness.
And
here we go! The orbiter itself has a great and informative website.
They provide links not only to interesting images and multimedia, but
also to actual real measurements and datasets via NASA's PlanetaryData System archives (sounds mind blowing!).
The
actual measured raw data can be downloaded from the Lunar OrbitalData Explorer.
For example, a single set of measurements can be imported into R like
this:
raw.data<-read .table="" a="" href="http://pds-geosciences.wustl.edu/lro/lro-l-lola-3-rdr-v1/lrolol_1xxx/data/lola_radr/laser1/lro_no_01/lolaradr_092582345.tab">http://pds-geosciences.wustl.edu/lro/lro-l-lola-3-rdr-v1/lrolol_1xxx/data/lola_radr/laser1/lro_no_01/lolaradr_092582345.tab-read>
")
head(raw.data)
This is
a table file which can be nicely imported into an R data frame. Its
structure is explained here.
The
main point is that the first two columns contain longitude/latitude
coordinates, while column 10 is the distance of the surface measured
from the orbiter. It is a simplification, because the measuring laser
is not always pointing “down”, but for our purposes that will do.
Because after investigating a few of the files (and there are many of
them), it is clear that they contain data segmented by collection
time, meaning that in a file there are reads from a segment of a
single orbit. Since I am more interested in all the measurements from
a given lunar location, I was digging further on the website.
Fortunately,
they offer also a tool to download all the elevation reads from a
location, and it can be selected using simply the coordinates of the
desired region: http://ode.rsl.wustl.edu/moon/indextools.aspx
So I
have fired up Google Moon and picked up one of my favorite regions: the Posidonius Crater. There are multiple crater rings in this complex, and I have randomly
picked up crater “J” and its surrounding and downloaded data
between 30 to 32 E; and 33 to 35 N coordinates.
The
data is kindly generated as a CSV file, which is very easy to handle
with R:
dat.file<- font="" opofull_csv_table.csv="">->
my.data<-read .csv="" dat.file="" font="">-read>
head(my.data)
summary(my.data$topography)
The
important information is stured in the first three columns of the
created data frame: the longitudes, the latitudes, and the
“topography”, which can be simply interpreted as elevation
compared to a hypothetical “sea level”.
After a
brief investigation of this topography data, it can be seen that it
is actually “below” that level, so all the elevations are
negative:
summary(my.data$topography)
It is
easy to use the plot3d() function to draw a 3D scatterplot of
the measured data. I wanted to use a color code which would resemble
map colors, with high values as brown and low areas as green. The
terrain.colors() function offers a handy solution for that:
zlim<-range font="" my.data="" topography="">-range>
zlen<-zlim font="" zlim="">-zlim>
colorlut<-terrain .colors="" font="" zlen="">-terrain>
col<-colorlut font="" my.data="" topography-zlim="">-colorlut>
With
the colors prepared, plotting is a matter of providing the proper
parameters of plot3d():
plot3d(my.data$Pt_Longitude,
my.data$Pt_Latitude,
my.data$topography,
col=col,
zlim=c(-4000,11000),
xlab="Longitude",
ylab="Latitude",
zlab="Elevation",
main="Posidonius J crater on Moon",
sub="Lunar Reconnaissance Orbiter LOLA data")
This is
excellent, the passings of orbiter are clearly visible, and with a
bit of turning the plot this and that ways, we can clearly see the
shape of the terrain. Still, it is not what I wanted as a nice 3D
topography. For that, the measurements taken under the orbital
passings should e converted to a regular grid format, where the
elevation is recorded in a grid like matrix. It is not a trivial task
to convert our existing data to grid format as it involves heavy
interpolation. It is a computationally expensive calculation, but
nothing that a modern computer could not handle. R offers the akima
package, which implements the interpolation of H. Akima
(http://www.iue.tuwien.ac.at/phd/rottinger/node60.html
).
library(akima)
int.dat <- font="" interp="" my.data="" t_longitude="">->
my.data$Pt_Latitude,
my.data$topography,
xo=seq(30,32,length=100),
yo=seq(33,35,length=100))
It
takes some 10-15 minutes to run the interpolation, but the results
are suitable for
visualization with the rgl
package. First, I wanted to convert the elevation values so that they
scale close to 1:10 with the horizontal dimensions of the plot.
elevation<-int .dat="" font="" z="">-int>
elevation<-elevation 1:10="" distortion="" font="" of="" sizes="" vertical="">-elevation>
Again, the colors are recalculated so the usual topographical colors
are scaled to this data:
zlim<-range elevation="" na.rm="TRUE)</font">-range>
zlen<-zlim font="" zlim="">-zlim>
colorlut<-terrain .colors="" font="" zlen="">-terrain>
col<-colorlut elevation-zlim="" font="">-colorlut>
And finally, we are ready to see the surface of the Moon!
library(rgl)
rgl.open()
rgl.surface(100:1,1:100,elevation,color=col,back="lines")
Incredible!
A piece of Moon
on my computer screen!