How to: GDAL with Martian Data
I’ve seen a lot of tutorial on working with GDAL with Earth data, but not a lot for planetary data, Mars included. As I use GDAL a lot (and in need of reminder myself), why not write a how-to with Martian data.
MOLA, but in different projections (https://core2.gsfc.nasa.gov/PGDA/images/MarsTopography_GMM3.jpeg, https://astrogeology.usgs.gov/search/map/Mars/GlobalSurveyor/MOLA/Mars_MGS_MOLA_ClrShade_merge_global_463m)
There is a reason why there is not a lot of discussion about GDAL in planetary community. It is because a lot of data processing, conversion, projection, are done in ISIS3. However, ISIS3 doesn’t work that well with HRSC data (my main Martian data), VICAR-formatted data, and GeoTIFF you produce yourself. Because of that, I suggest to have gdal installed as well.
Preparation: install GDAL. You can download GDAL here. Installing GDAL on Mac here (don’t forget to read the readme). Installing GDAL on Python on Mac here. To work with VICAR data, gdal version has to be higher than 2.0.0. Sadly if you run RHEL, you can’t go further than 1.11.0. Also don’t forget to add gdal on your PATH.
Checking gdal version (also useful to check whether gdal is properly installed and added to PATH)
Preparation: set projection information.
Compared to Earth where proj4 for projection are easily available as well as definable with just a few characters like WGS84 and EPSG:4926, you need to define yours for Mars. Instead of typing a new CRS/proj4 every time you convert a strip, why not save them in some files
Generally, we can have MOLA Sphere (a,b,c=3396000) or areoid (a,b=3396190 c=3376200). For the projection, we can have latlon, equirectangular, sinusoidal, or stereographic.
Create a new textfile (in Notepad, gedit, vi, or Atom), write the proj4, save it under name you easily remember (remember, no space!) such as MarsSouthStereo, MarsEQC, MarsMOLAEQC, etc.
Equirectangular areoid projection, centered in Lat=0 and Lon=0
+proj=eqc +lat_ts=0 +lat_0=0 +lon_0=0 +x_0=0 +y_0=0 +a=3396190 +b=3396190 +units=m +no_defs
LatLon projection, MOLA sphere
+proj=longlat +a=3396000 +b=3396000 +no_defs
Is useful to convert IMG to GeoTIFF, set no data value, change bit depth
Converting an HRSC vicar file to GeoTIFF with world file (TFW), world file is defined by projection information
gdal_translate H2288_0000_DT4.vic -co TFW=YES h2288_000_dt4.tif
Changing bit depth to 8bit and defining no data value as -32768
gdal_translate H2288_0000_DT4.tif –a_nodata –32768 –ot Int H2288_0000_8bit_DT4.tif
Removing georeference information
gdal_translate H2288_0000_DT4.tif –co PROFILE=BASELINE H2288_0000_DT4_base.tif
gdal_translate -srcwin 0 0 256 256 H2288_0000_DT4.tif H2288_0000_DT4_crop.tif
Crop according to Lat/Lon (projwin_srs available from gdal 2.0.0 and above)
gdal_translate -projwin -89.2 1.2 -89.1 0.9 -projwin_srs “+proj=longlat +a=3396000 +b=3396000 +no_defs” H2288_0000_DT4.tif H2288_0000_DT4_crop.tif
or if you have defined your proj4 as a file named marsMOLAlatlon
gdal_translate -projwin -89.2 1.2 -89.1 0.9 -projwin_srs MarsMOLAlatlon H2288_0000_DT4.tif H2288_0000_DT4_crop.tif
Obtain projection information, size, bit depth, no data
The same with gdalinfo, with minimum and maximum values
gdalinfo -stats H2288_0000_DT4.tif
Useful to change projection and fitting a strip to another strip
Example (with MArsMOLAEQC and MarsSouthStereo as projection files in the same folder with the image strip)
gdalwarp -s_srs MarsMOLAEQC -t_srs MarsSouthStereo H2288_0000_DT4.tif H2288_0000_D4_south_DT4.tif
(Will add more when I’m working with gdal)