Mars GIS Data

This notebook uses data obtained from the USGS via the following URL.

http://astrogeology.usgs.gov/Teams/Geomatics/photogrammetry/topography_mars.html

SetDirectory[
        ToFileName[("FileName" /. 
                NotebookInformation[
                  EvaluationNotebook[]])[[1]]]];

Some header data was manually stripped off of the data prior to importing and is then reversed so that it is orineted correctly.

data = 
        Reverse[Import["mtm-05282e_dem.asc", "TABLE"]];

Dimensions[data]

Convert the data from an array of heights into a list of {x,y,z} triples using integer values for x and y.

newdata = 
        Flatten[Table[{i, j, data[[j, i]]}, {j, 1, Length[data]}, {i, 1, 
              Length[data[[1]]]}], 1];

Now the data is prepared for interpolation.

datafun = Interpolation[newdata]

Next, we want to apply some nice color transitions to the surface.  First, we need to define a function that can create a smooth blend of specified colors to be used as a color function.

Define some possible color heads.

ccolor = (_RGBColor | _Hue | _GrayLevel | _CMYKColor);

Define the blending function which accepts sequences of colors of the types described above.

ablendsequence[c : {ccolor, ccolor .. }, u_?NumericQ] /; 
        And[SameQ[Head[#] & /@ c], 
          SameQ[Length[#] & /@ c]] := Head[
              First[c]] @@ (({1 - #, #} &[
                    FractionalPart[#]]) . (List @@@ 
                    Take[RotateLeft[c, Floor[#]], 2])) &[
        u*(-1 + Length[c])]

Create a surface plot of the interpolated data using flat colors initially via ColorFunction.  Colors range from blue for the lowest areas to white for the highest areas.  The PlotPoints value is set rather high to create a very detailed surface (better seen in the next step).

gr = 
      Graphics3D[
        Plot3D[datafun[x, y], {x, 50, 1490}, {y, 50, 1490}, 
          Mesh -> False, BoxRatios -> {1, 1,  .01}, 
          ImageSize -> 1280, 
          ColorFunction -> (ablendsequence[{RGBColor[0, 0,  .5], 
                    RGBColor[ .3,  .5,  .3], Brown, RGBColor[1, 1, 1], 
                    RGBColor[1, 1, 1]}, #] &), 
          PlotPoints -> 600]]

Now, we can transform the flat colors into SurfaceColors which respond to the specified lighting conditions. We also zoom in a bit to show more detail.

newgr = 
      Show[Graphics3D[
          gr[[1]] /. 
            RGBColor[r_, g_, b_] -> SurfaceColor[RGBColor[r, g, b]]], 
        BoxRatios -> {1, 1,  .1},  
        LightSources -> {{{1, 1, 1}, White}}, ImageSize -> 1280, 
        Background -> Black, Boxed -> False, 
        PlotRegion -> {{- .3, 1 +  .3}, {- .3, 
              1 +  .1}}]