You can download the tar.gz archive of the source code here.

This archive corresponds to the version used by E. Marcq et al. in their 2017 paper available here.

Contents of README.txt

This code was used to produce the results discussed in Marcq et al. (2017) paper. It simulates a 1D atmospheric column of H2O and CO2 above a solid surface, and computes a low-resolution spectrum of the thermal emission at the top of the atmosphere. For more information about the physics of the model, please refer to the above cited paper:

It makes extensive use of the provided and publicly available DISORT radiative transfer code:


A working Fortran 90/95 compiler is required. Our code has only been tested using the following Fortran compilers and OS:

Although not strictly necessary, it is strongly recommended to have a working f2py installation, available in the NumPy library for Python 2.7+. This enables to use the model through any Python v2.7+ environment with Numpy installed.


Source code tree

 |   |
 |   ---steam.f90 
 |   |
 |   ---k_correle.f95
 |   |
 |   ---make_tau_grey.f95
     ---(all provided DISORT subroutines)

Auxiliary files

g-grid for k-correlated data
T-grid for k-correlated data
p-grid for k-correlated data
Q-grid for k-correlated data
H2O-H2O continuum LUT from MTCKD
CO2-CO2 continuum LUT from B. Bézard
k-correlated coefficients computed by KSPECTRUM
directory required for diagnostic Fortran outputs (created empty by unpacking the *.tar.gz file, so it should be OK).


  1. untar the archive file
  2. compile the source files
    1. using the script (if you have a working f2py installation) to produce a Python library named that you can import from Python — you may want to change the path to your python installation in this script.
    2. using the script for a Fortran-only installation. This will produce a radconv1d executable file. The main program is just a wrapper which provides the main subroutine the required inputs (see below for a description of these inputs).

Using the model with Python (recommended)

Example 1 (main subroutine)

from radconv1d import radconv
Ttop_out, flux_up, P0_v, P0_c, P0_o, spectrum = radconv(Tsurf, Ttop_in, Psurf_v, Psurf_c, Psurf_o, Ttop_flag,
                                                        cloud_flag, cont_flag, grey_flag, GP_flag, spher_flag,
                                                        alb_clear, alb_cloud, const_sol, g0, Rt, Nz)

Input parameters

Surface temperature in K
Starting mesospheric temperature in K. 200 is a recommanded value.
Total surface pressure of H2O in Pa (including liquid H2O if any)
Total surface pressure of CO2 in Pa
Total surface pressure of N2 in Pa
if == 1, Ttop will not be altered. Else, the model will try to adjust Ttop so that the TOA divergence of OLR is null. This results in unrealistically cold mesospheres, and is kept only in anticipation of stellar heating modeling!
if == 0, clouds won't be radiatively active. If == 1, Earth-like clouds are radiatively assumed. If == 2, Venus-like clouds are radiatively assumed (not realistic for H2O-rich atmospheres)
if == 0, no continuum opacity is included. If == 1, only CO2-CO2 opacity is included. If == 2, only H2O-H2O opacity is included. If == 3, H2O-H2O and CO2-CO2 opacities are included (recommended).
if == 0, line opacities are computed according to k-correlated data (slow, but more realistic). Else, line opacities will be averaged over the whole spectrum (much faster, but less accurate for the total integrated OLR, and plainly wrong for the OLR spectrum).
if ==1, H2O steam is treated as an ideal gas. Else, its EOS is given by the included steam tables.
future flag for sphericity corrections. Should be set equal to 0.
Prescribed albedo without clouds. Preliminary SW calculations yield about 0.2 around a G-star.
Prescribed albedo with clouds. Preliminary SW calculations yield about 0.7 around a G-star.
Averaged stellar constant in W/m². == 1/4 of its value at the substellar point.
Surface gravity in m/s²
Planetary radius in m
Number of computational layers. Typically around 200.


Final mesospheric temperature in K (differs from Ttop_in only if Ttop_flag != 1).
Total OLR in W/m²
Gaseous surface pressure of H2O in Pa (excluding surface liquid H2O if any).
Gaseous surface pressure of CO2 in Pa (== Psurf_c in this version)
Gaseous surface pressure of N2 in Pa (== Psurf_o in this version)
Spectrum of OLR in W/m²/cm-1 in each of the 36 bands. Boundary wavenumbers are given below:
sigmas = np.array([0., 40., 160., 280., 380., 500., 582., 600., 720.,
                    752., 800., 900., 1000., 1200., 1350., 1450., 1550.,
                    1650., 1750., 1850., 1950., 2050., 2200., 2500., 2800.,
                    3200., 3600., 4000., 4400., 4800., 5900.,6000., 6500.,
                    8000., 8300., 9300., 10100.]) # in inv-cm

Example 2 (make_profile only)

from radconv1d import make_profiles

Z, T, Pv, Pc, Po, rho_v, rho_c, rho_o, ibottom, itop = make_profiles(GP_flag, g0, Rt, Tsurf, Psurf_v, Psurf_c,
                                                                     Psurf_o, Ttop, Nz)


Input parameters

See above.


(Nz+1)-long 1D array of boundary altitudes for computational layers in m.
(Nz+1)-long 1D array of boundary temperatures for computational layers in K.
(Nz+1)-long 1D array of boundary H2O pressures for computational layers in Pa.
(Nz+1)-long 1D array of boundary CO2 pressures for computational layers in Pa.
(Nz+1)-long 1D array of boundary N2 pressures for computational layers in Pa.
(Nz)-long array of H2O densities within computational layers in kg/m3.
(Nz)-long array of CO2 densities within computational layers in kg/m3.
(Nz)-long array of N2 densities within computational layers in kg/m3.
index of the lowermost computational layer in the moist adiabat layer
index of the uppermost computational layer in the moist adiabat layer

Using the model with Fortran only

Starting from the provided test_radconv1d.f90 Fortran wrapper, the user is free to provide required inputs to the main subroutine named RADCONV, and process its ouputs afterwards. If no spectral output is required but merely vertical profiles in composition, density and temperature, you can call only MAKE_PROFILES instead of RADCONV.

For the physical meaning of the inputs and outputs of the aforementioned subroutines, please refer to the Python section hereabove.


If you find an error in the code or wish to submit a feature request, please contact the author by e-mail.


Dr Emmanuel Marcq, Associate Professor, LATMOS - Univ. Paris-Saclay Tel: +33 1 80 28 52 83


This wode has been released under Creative Common License BY-SA-NC

BY: Licensees may copy, distribute, display and perform the work and make derivative works and remixes based on it only if they give the author or licensor the credits (attribution) in the following manner:
"Author: E. Marcq/LATMOS-IPSL-UPSAY ("

SA: Licensees may distribute derivative works only under a license identical ("not more restrictive") to the license that governs the original work.

NC: Licensees may copy, distribute, display, and perform the work and make derivative works and remixes based on it only for non-commercial purposes.