Usage

The functionalities of gcm_toolkit is likely to expand in the future, since gcm_toolkit is in active development. However, gcm_toolkit can already be used to perform a variety of tasks. These docs explain how gcm_toolkit can be used to load data and do plots with the data.

GCMT class

The GCMT class is the main point of user interaction. It is used to load data into xarray datasets, which can then be used for analysis.

from gcm_toolkit import GCMT

# Load data
tools = GCMT()
tools.read_raw(..., tag='model_name')
# tools.read_raw(..., tag='model_name_2')  # if you have multiple models
ds = tools['model_name']

The GCMT class can also be used for quick plotting of data, since it wraps a few plotting functions, which are outlined in Plotting.

Note

You can also load netcdf data:

tools.read_reduced(..., tag='model_name')

class gcm_toolkit.GCMT(p_unit='bar', time_unit='day', write='on')

The main gcm_toolkit class with which the user can interact.

Attributes:

modelsdict

Shortcut to return all GCMs in memory.

Methods

read_raw(self, gcm, data_path, tag=None):	Read in the raw data output from a GCM.
read_reduced(self, data_path, tag=None):	Read in the previously reduced GCM

__init__(p_unit='bar', time_unit='day', write='on')

Constructor for the gcm_toolkit class.

Parameters:

p_unit: str, optional

Set the unit that is used internally for pressure related things

time_unit: str, optional

Set the unit that is used internally for time related things

write: str, optional

Set the level of logging information. Check core.writer for more infos

property models

Shortcut to return all GCMs in memory.

Returns:

selected_modelsGCMDatasetCollection or xarray Dataset

All models in self._models

get(tag, default=None)

Pythonic get, will return default if tag is not available.

Parameters:

tag: str

Name of the model that should be returned

default: Any

Default value in case that the tag is not available in the dataset

Returns:

dsi or default, depending on whether dsi is available

get_models(tag=None, always_dict=False)

Function return all GCMs in memory. If a tag is given, only return this one.

Parameters:

tagstr

Name of the model that should be returned.

always_dict: bool

Force result to be a dictionary (if tag is None)

Returns:

selected_modelsGCMDatasetCollection or xarray Dataset

All models in self._models, or only the one with the right tag. Will definetly be GCMDatasetCollection if always_dict=True

read_raw(gcm, data_path, iters='last', load_existing=False, tag=None, **kwargs)

General read in function for GCM data

Parameters:

gcmstr

Type of GCM, must be ‘MITgcm’.

data_pathstr

Folder path to the standard output of the GCM.

iterslist, str

The iteration (time step) of the input files to be read. If None, no data will be read. If ‘last’ (default), only the last iteration will be read. If ‘all’, all iterations will be read.

load_existing: bool

Set to false if you want to overwrite already loaded data Set to true if you want to increment already loaded data

tagstr

Tag to reference the simulation in the collection of models.

kwargs: dict

Additional options passed down to read functions

read_reduced(data_path, tag=None, time_unit_in=None, p_unit_in=None)

Read in function for GCM data that has been reduced a nd saved according to the gcm_toolkit GCMDataset format.

Parameters:

data_pathstr

Folder path to the reduced (gcm_toolkit) data.

time_unit_in: str, None

units of time dimension in input dataset. If None, try to read from nc file (ds.attr.time_unit)

p_unit_in: str, None

units of pressure dimensions in input dataset, If None, try to read from nc file (ds.attr.p_unit)

tagstr

Tag to reference the simulation in the collection of models.

save(direct, method='nc', update_along_time=False, tag=None)

Save function to store current member variables.

Parameters:

directstr

directory at which the gcm_toolkit datasets should be stored.

methodstr, optional

Datasets can be stored as ‘.zarr’ or ‘.nc’. Decide which type you prefer. Defaults to ‘.nc’.

update_along_time: bool, optional

Decide if you want to update already saved datasets along the timedimension. This only works with method=’zarr’.

tag: str, optional

tag of the model that should be loaded. Will save all available models by default.

Returns:

NoneType

None

load(direct, method='nc', tag=None)

Load function to load stored member variables.

Parameters:

directstr

directory at which the gcm_toolkit datasets are stored

methodstr, optional

Should be the same method with which you stored the data

tag: str, optional

tag of the model that should be loaded. Will load all available models by default.

Postprocessing

The GCMT class has quite a few methods that can calculate diagnostics such as total energy, location of the RCB and so on. You find their documentation here.

Example for the total energy:

from gcm_toolkit import GCMT
import gcm_toolkit.gcm_plotting as gcmp

# Load data
tools = GCMT()
tools.read_raw(..., tag='model_name')
E = tools.add_total_energy(var_key_out='E', tag='model_name')

# You will also have the energy stored in the dataset if var_key_out is not None:
ds = tools['model_name']
E == ds.E

class gcm_toolkit.GCMT(p_unit='bar', time_unit='day', write='on')

The main gcm_toolkit class with which the user can interact.

Attributes:

modelsdict

Shortcut to return all GCMs in memory.

Methods

read_raw(self, gcm, data_path, tag=None):	Read in the raw data output from a GCM.
read_reduced(self, data_path, tag=None):	Read in the previously reduced GCM

add_horizontal_average(var_key, var_key_out=None, part='global', area_key='area_c', tag=None)

Calculate horizontal averaged quantities. Horizontal averages are calculated as area-weighted quantities q, so that:

bar{q} = int{q dA}/int{dA}

Parameters:

var_key: str, xarray.DataArray

The key or array of the variable quantity that should be averaged. If str, it will try to look up the key in the dataset. If DataArray, it will use this one instead.

var_key_out: str, optional

variable name used to store the outcome. If not provided, this script will just return the averages and not change the dataset inplace.

part: dict or str, optional

‘global’: global average ‘night’: only nightside (defined around +-180,0) ‘day’: only dayside (defined around 0,0) ‘morning’: morning terminator (average around lon=[-100,-80]) ‘evening’: evening terminator (average around lon=[80,100]) Alternatively you may specify a dict in the following way (example for morn. term.): part = {‘lon’: [-100,-80], ‘lat’:[-90,90], ‘inv’:False} The ‘lon’, ‘lat’ specify regions of lon and lat that should be used, whereas ‘inv’ (optional, default False) gives the option to invert the lon and lat regions (e.g., exclude instead of include for average)

area_key: str, optional

Variable key in the dataset for the area of grid cells

tagstr, optional

The tag of the dataset that should be used. If no tag is provided, and multiple datasets are available, an error is raised.

Returns:

avg: xarray.DataArray

Averaged data

add_rcb(tol=0.01, var_key_out=None, area_key='area_c', temp_key='T', part='global', tag=None)

Calculate the radiative convective boundary (rcb) by searching (from the bottom upwards) for the first occurance of a deviation from an adiabatic temperature profile.

Operates on and calculates the horizontal average.

Parameters:

tol: float

tolerance for the relative deviation from adiabat

var_key_out: str, optional

variable name used to store the outcome. If not provided, this script will just return the averages and not change the dataset inplace.

part: dict or str, optional

‘global’: global average ‘night’: only nightside (defined around +-180,0) ‘day’: only dayside (defined around 0,0) ‘morning’: morning terminator (average around lon=[-100,-80]) ‘evening’: evening terminator (average around lon=[80,100]) Alternatively you may specify a dict in the following way (example for morn. term.): part = {‘lon’: [-100,-80], ‘lat’:[-90,90], ‘inv’:False} The ‘lon’, ‘lat’ specify regions of lon and lat that should be used, whereas ‘inv’ (optional, default False) gives the option to invert the lon and lat regions (e.g., exclude instead of include for average)

area_key: str, optional

Variable key in the dataset for the area of grid cells

temp_key: str, optional

The key to look up the temperature (needed for density calculation)

tagstr, optional

The tag of the dataset that should be used. If no tag is provided, and multiple datasets are available, an error is raised.

Returns:

rcbxarray.DataArray

A dataArray with reduced dimensionality, containing the pressure of the rcb location.

add_theta(var_key_out=None, temp_key='T', tag=None)

Convert temperature to potential temperature with respect to model boundary.

Parameters:

var_key_out: str, optional

variable name used to store the outcome. If not provided, this script will just return theta and not change the dataset inplace.

temp_key: str, optional

The key to look up the temperature

tagstr, optional

The tag of the dataset that should be used. If no tag is provided, and multiple datasets are available, an error is raised.

Returns:

thetaxarray.DataArray

A dataArray with reduced dimensionality, containing the potential temperature

add_total_energy(var_key_out=None, area_key='area_c', temp_key='T', return_all=False, tag=None)

Calculate the total Energy of the GCM. See e.g., https://ui.adsabs.harvard.edu/abs/2014Icar..229..355P, Eq. 16

Parameters:

var_key_out: str, optional

variable name used to store the outcome. If not provided, this script will just return the averages and not change the dataset inplace.

area_key: str, optional

Variable key in the dataset for the area of grid cells

temp_key: str, optional

The key to look up the temperature

return_all: bool, optional

Also return the partial energies

tagstr, optional

The tag of the dataset that should be used. If no tag is provided, and multiple datasets are available, an error is raised.

Returns:

energyxarray.DataArray

A dataArray with reduced dimensionality, containing the total energy.

therm_energyxarray.DataArray, optional

A dataArray with reduced dimensionality, containing the thermal energy.

pot_energyxarray.DataArray, optional

A dataArray with reduced dimensionality, containing the potential energy.

kin_energyxarray.DataArray, optional

A dataArray with reduced dimensionality, containing the kinetic energy.

add_total_momentum(var_key_out=None, area_key='area_c', temp_key='T', tag=None)

Calculate the total angular momentum of the GCM. See e.g., https://ui.adsabs.harvard.edu/abs/2014Icar..229..355P, Eq. 17

Parameters:

var_key_out: str, optional

variable name used to store the outcome. If not provided, this script will just return the averages and not change the dataset inplace.

area_key: str, optional

Variable key in the dataset for the area of grid cells

temp_key: str, optional

The key to look up the temperature (needed for density calculation)

tagstr, optional

The tag of the dataset that should be used. If no tag is provided, and multiple datasets are available, an error is raised.

Returns:

momentumxarray.DataArray

A dataArray with reduced dimensionality, containing the total momentum.

add_meridional_overturning(v_data='V', var_key_out=None, tag=None)

Calculate meridional overturning streamfunction. This quantity psi is computed by integrating the zonal-mean meridional velocity ar V along pressure, and weighting with 2*pi*R_p / g times the cosine of latitude, where R_p is the planetary radius and g is the surface gravity: bar{psi} = 2 pi R_p / g cos{lat} int{bar{V} dp’} (see e.g. Carone et al. (2018), Eq. 7)

Parameters:

v_data: str

The key that holds the data to meridional velocity

var_key_out: str, optional

variable name used to store the outcome. If not provided, this script will just return the overturning circulation and not change the dataset inplace.

tagstr, optional

The tag of the dataset that should be used. If no tag is provided, and multiple datasets are available, an error is raised.

Returns:

psi: xarray.DataArray

Horizontal overturning

Plotting

gcm_toolkit includes some very simple (but yet powerful and highly customizable) plotting functions that can be used in two different ways. You can either use the gcm_toolkit class to invoke these functions or call these functions directly.

Example:

from gcm_toolkit import GCMT
import gcm_toolkit.gcm_plotting as gcmp

# Load data
tools = GCMT()
tools.read_raw(..., tag='model_name')
ds = tools['model_name']

# Plot data with GCMT
tools.isobaric_slice(tag = 'model_name', p=1e-2, lookup_method='nearest', var_key='T', wind_kwargs={'windstream':False, 'sample_one_in':2})

# Plot data with GCMTtools.gcm_plotting
tools.isobaric_slice(ds, p=1e-2, lookup_method='nearest', var_key='T', wind_kwargs={'windstream':False, 'sample_one_in':2})

Both ways are equally powerful and flexible.

gcm_toolkit.gcm_plotting.isobaric_slice(dsi, var_key, pres, time=-1, lookup_method='exact', ax=None, plot_windvectors=True, wind_kwargs=None, cbar_kwargs=None, add_colorbar=True, fs_labels=None, fs_ticks=None, title=None, xlabel='Longitude (deg)', ylabel='Latitude (deg)', contourf=False, **kwargs)

Plot an isobaric slice of the dataset and quantity at the given pressure level.

Parameters:

dsiDataSet

A gcm_toolkit-compatible dataset of a 3D climate simulation.

var_keystr

The key of the variable quantity that should be plotted.

presfloat

Pressure level for the isobaric slice to be plotted, expressed in the units specified in the dataset attributes (e.g., init of GCMT object).

timeint, optional

Timestamp that should be plotted. By default, the last time is selected.

lookup_methodstr, optional

The look-up method that is used to slice along pressure: ‘exact’ for exactly matching key look-up (default); ‘nearest’ to pick out the nearest neighbour Z coordinate; ‘interpolate’ for a linear interpolation along the Z axis.

axmatplotlib.axes.Axes, optional

The axis on which you want your plot to appear.

plot_windvectorsboolean, optional

If True (default), plot the horizontal wind vectors over the colour mesh.

wind_kwargsdict, optional

Additional keywords for the method plot_horizontal_wind.

cbar_kwargsdict, optional

Additional keywords for the colorbar.

add_colorbar: bool, optional,

Optionally decide if you want a colorbar or don’t

fs_labelsint, optional

Optionally set font size of the axis labels.

fs_ticksint, optional

Optionally set font size of the tick labels.

titlestr, optional

Title for the isobaric slice plot. By default, the selected pressure and time stamp of the slice are displayed.

xlabelstr, optional

X-axis label, longitude by default.

ylabelstr, optional

Y-axis label, latitude by default.

contourf: bool, optional

Decide if you want to do a contourplot or a pcolormesh plot

gcm_toolkit.gcm_plotting.zonal_mean(dsi, var_key, time=-1, ax=None, cbar_kwargs=None, fs_labels=None, xlabel='Latitude (deg)', ylabel='Z', add_ylabel_unit=True, title=None, add_colorbar=True, contourf=False, **kwargs)

Plot a zonal mean average of a quantity for the given dataset.

Parameters:

dsiDataSet

A gcm_toolkit-compatible dataset of a 3D climate simulation.

var_keystr

The key of the variable quantity that should be plotted.

timeint, optional

Timestamp that should be plotted. By default, the last time is selected.

axmatplotlib.axes.Axes, optional

The axis on which you want your plot to appear.

cbar_kwargsdict, optional

Additional keywords for the colorbar.

fs_labelsint, optional

Optionally set font size of the axis labels.

xlabel: str, optional

Label for x

ylabel: str, optional

Label for y

add_ylabel_unit: bool, optional

Optionally decide, if you want to add a unit to ylabel.

titlestr, optional

Title for the isobaric slice plot. By default, the selected pressure and time stamp of the slice are displayed.

add_colorbar: bool, optional

Optionally decide if you want a colorbar or don’t

contourf: bool, optional

Decide if you want to do a contourplot or a pcolormesh plot

gcm_toolkit.gcm_plotting.time_evol(dsi, var_key, ax=None, fs_labels=None, cbar_kwargs=None, add_colorbar=True, title=None, xlabel=None, ylabel='Z', add_ylabel_unit=True, **kwargs)

Function that plots the time evolution of a quantity in a 1D line collection plot, where the colorscale can be related to the time evolution. Note: var_key needs to contain data that is 2D in time and pressure.

Parameters:

dsiDataSet

A gcm_toolkit-compatible dataset of a 3D climate simulation.

var_keystr

The key of the variable quantity that should be plotted.

axmatplotlib.axes.Axes, optional

The axis on which you want your plot to appear.

fs_labelsint, optional

Optionally set font size of the axis labels.

cbar_kwargsdict, optional

Additional keywords for the colorbar.

add_colorbar: bool, optional

Decide if you want a colorbar

titlestr, optional

Title for the isobaric slice plot. By default, the selected pressure and time stamp of the slice are displayed.

xlabel: str, optional

Label for x

ylabel: str, optional

Label for y

add_ylabel_unit: bool, optional

Optionally decide, if you want to add a unit to ylabel.

Returns:

l: LineCollection

the collection of plotted lines

Interface

gcm_toolkit comes with some interfaces to other codes. We currently have support for petitRADTRANS to calculate spectra and phasecurves with prt_phasecure.

class gcm_toolkit.utils.interface.PrtInterface(tools, prt)

Interface with petitRADTRANS

__init__(tools, prt)

Constructs the Interface and links to pRT.

Parameters:

tools: GCMTools Object

A GCMTools object that is linked to the interface

pRT: petitRADTRANS.Radtrans

A pRT Radtrans object

set_data(time, tag=None, regrid_lowres=False, terminator_avg=False, lat_points=1, lon_resolution=10)

Set the data to be used for the interface

Parameters:

time: int

timestep to be used

tag: str

tag of the model to be used

regrid_lowres: bool, optional

Can be useful, if your GCMT uses a very detailed grid

terminator_avg: bool, optional

terminator avaraging. This requires lat_points and lon_resolution.

lat_points: int, optional

number of equally spaced latitude grid points for each terminator (morning and evening)

lon_resolution: float, optional

longitudinal opening angle for terminator avareging

calc_phase_spectrum(mmw, Rstar, Tstar, semimajoraxis, gravity=None, filename=None, normalize=True, **prt_args)

Calculate the spectrum for a phasecurve.

Parameters:

mmw: float or 1D-array

Mean molecular weight (in atomic units). Will be globally uniform if float or horizonatally uniform if 1D.

Rstar: float

stellar radius in !cm!.

Tstar: float

Temperature of the hoststar.

semimajoraxis: float

The distance between hoststar and planet in !cm!

gravity: float

surface gravity in !cgs!. Will default to the value provided by GCMT.

filename: str

path at which the output should be stored.

normalize: bool

Decide if you want the spectrum to be normalized to the star. The code will then: 1. correct intensity for the ratio between solar radius and

planetary radius

2. devide by stellar spectrum

prt_args:

All the args that should be parsed to calc_spectra. See the docs of prt_phasecurve for more info on the arguments.

Returns:

spectra: xr.DataArray

Dataarray containing the spectrum. The spectrum is normed to the stellar spectrum.

phase_curve(phases, spectra=None, filename=None)

Do the diskintegration of the spectrum to yield the phasecurve

Parameters:

phases (array(P)):

List of phases at which the phasecurve should be evaluated

spectra: str

The spectrum generated by calc_phase_spectrum

filename: str

Alternatively give a filename where the spectrum is saved

Returns:

phasecurve: xr.DataArray

DataArray containing the phasecurveinformation

chem_from_poorman(temp_key='T', co_ratio=0.55, feh_ratio=0.0)

Calculate equilibrium abundancies with poorman from pRT

Parameters:

temp_key: str, optional

The key to the temperature field used for the abundancies. Defaults to T.

co_ratio: float, optional

The C/O ratio. Currently only one global value allowed. Defaults to 0.55.

feh_ratio: float, optional

The metalicity ratio. Currently only one global value allowed.

Defaults to 0.0.