Building Surface Topography Files

From Planets
Revision as of 16:04, 26 October 2023 by Martin.turbet (talk | contribs) (Convert the data into a netcdf surface topography map)

Jump to: navigation, search

Here we present a simple tutorial to build a new topography map that can be used in the Generic PCM.

Surface topography maps are netcdf files stored in your /datagcm/surface_data/ repository. Such files contain information on terrain elevations, as well as albedo and thermal inertia.

As a reminder, our database of surface topography files is available here: https://web.lmd.jussieu.fr/~lmdz/planets/LMDZ.GENERIC/datagcm/surface_data/

Getting topography and surface data

A few options to do that:

You can also derive data files by digitalizing color maps. See below an example:

HERE WE ADD AN EXAMPLE

Convert the data into a netcdf surface topography map

We have developed Python and fortran routines to easily convert surface data files (elevation, albedo, thermal inertia) into GCM-compatible surface topography files.

Please find below examples:

from    numpy import *
import  numpy                 as        np
import  matplotlib.pyplot     as        mpl
import  math
import xarray as xr

# 1. WE GET THE SURFACE TOPOGRAPHY DATA
nc = xr.open_dataset('surface_mars.nc',decode_times=False) # can be any netcdf file (e.g. start/startfi.nc files)

# 2. WE READ THE VARIABLES
lat=nc['latitude']
lon=nc['longitude']
albedo=nc['albedo']
thermal_inertia=nc['thermal']
elevation=nc['zMOL']

# BELOW THE VARIABLES WE WANT TO UPDATE
new_elevation = np.empty((len(lat),len(lon)))
new_thermal_inertia = np.empty((len(lat),len(lon)))
new_albedo = np.empty((len(lat),len(lon)))

# LOOP TO MODIFY THE VARIABLES
for i in range(0,len(lat),1):
    for j in range(0,len(lon),1):
        new_elevation[i,j]=max(elevation[i,j],-3.9) # here you put whatever you want ; in this example, we take the (MOLA) topography of present-day Mars and modify it to fill the topographic depressions with an ocean at an altitude of -3.9km
        new_albedo[i,j]=albedo[i,j] # here you put whatever you want
        new_thermal_inertia[i,j]=thermal_inertia[i,j] # here you put whatever you want

# SANITY CHECK PLOTS

fig1 = mpl.figure(1)
mpl.contourf(lon,lat,elevation)
mpl.xlabel('longitude (deg)')
mpl.ylabel('latitude (deg)')
mpl.show()

fig2 = mpl.figure(2)
mpl.contourf(lon,lat,new_elevation)
mpl.xlabel('longitude (deg)')
mpl.ylabel('latitude (deg)')
mpl.show()

# CREATE THE NEW NETCDF TOPOGRAPHY FILE

nc['albedo'].values = new_albedo
nc['zMOL'].values = new_elevation
nc['thermal'].values = new_thermal_inertia
nc.to_netcdf('new_topography.nc')

Using a new topography in the model

follow the tutorial here: https://lmdz-forge.lmd.jussieu.fr/mediawiki/Planets/index.php/Advanced_Use_of_the_GCM