Note

Download full notebook here

Model Storage: Preservation and Reuse

Storing data is often tricky. You may be wondering what output to save, what output to throw away. You may also have questions surround what formats to use to store your data.

PICASO’s recommended model output is xarray. There are several xarray tutorials throughout the picaso documentation. The most thorough is in the data uniformity tutorial.

Here you will learn how to use a nifty picaso function that will store your the minimum data needed to reproduce your 1D spectrum. The one caveat is that we will not store opacity data. However, you can always get your reference opacity data using the auto citation tools.

For this tutorial we assume you already know the basics of running 1D models. Note that we have still not enabled this functionality for 3D and 4D runs. Please contact the developers if you are interested in this functionality.

[1]:

import picaso.justdoit as jdi
import picaso.justplotit as jpi
u = jdi.u #astropy units

WARNING: Failed to load Vega spectrum from /data/reference_data/picaso/ref4/stellar_grids/calspec/alpha_lyr_stis_011.fits; Functionality involving Vega will be severely limited: FileNotFoundError(2, 'No such file or directory') [stsynphot.spectrum]

Storing your run with xarray

Here we show a simple example of how to use the output_xarray function for preservation

[2]:

opa = jdi.opannection(wave_range=[0.3,1])
pl = jdi.inputs()#calculation='brown')
pl.gravity(radius=1, radius_unit=u.Unit('R_jup'),
           mass=1, mass_unit=u.Unit('M_jup'))
pl.atmosphere(filename=jdi.jupiter_pt(), sep='\s+')
pl.phase_angle(0)
pl.clouds(filename=jdi.jupiter_cld(), sep='\s+')
pl.star(opa, 5000,0,4, radius=1, radius_unit=u.Unit("R_sun"), semi_major=1, semi_major_unit=1*u.AU)
#MUST USE full output=True for this functionality
df= pl.spectrum(opa, calculation='reflected', full_output=True)

<>:5: SyntaxWarning: invalid escape sequence '\s'
<>:7: SyntaxWarning: invalid escape sequence '\s'
<>:5: SyntaxWarning: invalid escape sequence '\s'
<>:7: SyntaxWarning: invalid escape sequence '\s'
/tmp/ipykernel_104973/425976240.py:5: SyntaxWarning: invalid escape sequence '\s'
  pl.atmosphere(filename=jdi.jupiter_pt(), sep='\s+')
/tmp/ipykernel_104973/425976240.py:7: SyntaxWarning: invalid escape sequence '\s'
  pl.clouds(filename=jdi.jupiter_cld(), sep='\s+')

[3]:

preserve = jdi.output_xarray(df,pl,savefile='/data/picaso_dbs/model.nc')

[4]:

preserve

[4]:

<xarray.Dataset> Size: 724kB
Dimensions:         (pressure: 61, wavelength: 18060, pressure_layer: 60,
                     wavenumber_layer: 196)
Coordinates:
  * pressure        (pressure) float64 488B 1e-06 1.41e-06 2e-06 ... 708.0 1e+03
  * wavelength      (wavelength) float64 144kB 1.0 0.9999 0.9999 ... 0.3 0.3 0.3
Dimensions without coordinates: pressure_layer, wavenumber_layer
Data variables: (12/21)
    temperature     (pressure) float64 488B 150.9 149.7 ... 1.525e+03 1.675e+03
    albedo          (wavelength) float64 144kB 0.2477 0.3256 ... 0.7188 0.7187
    fpfs_reflected  (wavelength) float64 144kB 5.657e-08 7.436e-08 ... 1.642e-07
    e-              (pressure) float64 488B 4.5e-38 4.5e-38 ... 2.15e-11
    H2              (pressure) float64 488B 0.837 0.837 0.837 ... 0.836 0.836
    H               (pressure) float64 488B 4.5e-38 4.5e-38 ... 3.37e-06
    ...              ...
    NH3             (pressure) float64 488B 0.000137 0.000137 ... 0.000155
    N2              (pressure) float64 488B 5.42e-17 1.58e-17 ... 2.73e-05
    PH3             (pressure) float64 488B 4.5e-38 4.5e-38 ... 5.03e-07
    opd             (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
    ssa             (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
    asy             (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
Attributes:
    planet_params:   {"gravity": {"value": 2478.651947614914, "unit": "cm / s...
    stellar_params:  {"database": "ck04models", "steff": 5000, "feh": 0, "log...
    orbit_params:    {"sma": {"value": 14959787070000.0, "unit": "cm"}}

Reusing your run with xarray

We often revisit models time and time again. Maybe you want to change your wavelength range, or observing geometry (e.g. transit vs emission).

Here we show a simple example of how to use the input_xarray function for reuse. In this simple example, we will take our previous input but instead of computing reflected light, add thermal emission and transmission

[5]:

opa = jdi.opannection(wave_range=[0.3,14])
ds = jdi.xr.load_dataset('/data/picaso_dbs/model.nc')
reuse = jdi.input_xarray(ds, opa)
new_model = reuse.spectrum(opa, calculation='reflected+thermal', full_output=True)

[6]:

new_output= jdi.output_xarray(new_model, reuse)

[7]:

new_output

[7]:

<xarray.Dataset> Size: 2MB
Dimensions:        (pressure: 61, wavelength: 57646, pressure_layer: 60,
                    wavenumber_layer: 196)
Coordinates:
  * pressure       (pressure) float64 488B 1e-06 1.41e-06 2e-06 ... 708.0 1e+03
  * wavelength     (wavelength) float64 461kB 14.0 14.0 14.0 ... 0.3 0.3 0.3
Dimensions without coordinates: pressure_layer, wavenumber_layer
Data variables: (12/21)
    temperature    (pressure) float64 488B 150.9 149.7 ... 1.525e+03 1.675e+03
    flux_emission  (wavelength) float64 461kB 6.258e+06 6.259e+06 ... 1.238e-08
    albedo         (wavelength) float64 461kB 0.0195 0.0195 ... 0.7187 0.7187
    e-             (pressure) float64 488B 4.5e-38 4.5e-38 ... 4.46e-12 2.15e-11
    H2             (pressure) float64 488B 0.837 0.837 0.837 ... 0.836 0.836
    H              (pressure) float64 488B 4.5e-38 4.5e-38 ... 8.18e-07 3.37e-06
    ...             ...
    NH3            (pressure) float64 488B 0.000137 0.000137 ... 0.000155
    N2             (pressure) float64 488B 5.42e-17 1.58e-17 ... 2.73e-05
    PH3            (pressure) float64 488B 4.5e-38 4.5e-38 ... 5.11e-07 5.03e-07
    opd            (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
    ssa            (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
    asy            (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
Attributes:
    planet_params:   {"effective_temp": 88.91594769499913, "gravity": {"value...
    stellar_params:  {"database": "ck04models", "steff": 5000, "feh": 0, "log...
    orbit_params:    {"sma": {"value": 14959787070000.0, "unit": "cm"}}

Adding meta data

We have a tutorial that shows what some recommended meta data might be. output_xarray has a field add_output that allows you to add extra arguments to the xarray.Dataset.attrs. However, it has to be in a very specific format.

You can see the basic format by running this function:

[8]:

jdi.standard_metadata()

[8]:

{'author': 'optional',
 'contact': 'optional',
 'code': 'optional',
 'doi': 'optional',
 'planet_params': {'rp': 'usually taken from picaso',
  'mp': 'usually taken from picaso',
  'mh': 'optional',
  'cto': 'optional',
  'heat_redis': 'optional',
  'p_reference': 'usually taken from picaso',
  'tint': 'optional'},
 'stellar_params': {'logg': 'usually taken from picaso',
  'feh': 'usually taken from picaso',
  'steff': 'usually taken from picaso',
  'rs': 'usually taken from picaso',
  'ms': 'optional',
  'database': 'usually taken from picaso'},
 'orbit_params': {'sma': 'usually taken from picaso'}}

[9]:

opa = jdi.opannection(wave_range=[0.3,1])
pl = jdi.inputs()#calculation='brown')
pl.gravity(radius=1, radius_unit=u.Unit('R_jup'),
           mass=1, mass_unit=u.Unit('M_jup'))
mh = 0
cto = 0.55 #solar value

pl.atmosphere(filename=jdi.jupiter_pt(), sep='\s+')
pl.chemeq_visscher_2121(cto, mh)
pl.phase_angle(0)
pl.clouds(filename=jdi.jupiter_cld(), sep='\s+')
pl.star(opa, 5000,0,4, radius=1, radius_unit=u.Unit("R_sun"), semi_major=1, semi_major_unit=1*u.AU)
#MUST USE full output=True for this functionality
df= pl.spectrum(opa, calculation='reflected', full_output=True)

preserve = jdi.output_xarray(df,pl,add_output={
    'author':"Awesome Scientist",
    'contact' : "awesomescientist@universe.com",
    'code' : "picaso",
    'planet_params':{'mh':mh, 'co':cto},
    'cloud_params':{'fsed':3}
})

<>:8: SyntaxWarning: invalid escape sequence '\s'
<>:11: SyntaxWarning: invalid escape sequence '\s'
<>:8: SyntaxWarning: invalid escape sequence '\s'
<>:11: SyntaxWarning: invalid escape sequence '\s'
/tmp/ipykernel_104973/2498042645.py:8: SyntaxWarning: invalid escape sequence '\s'
  pl.atmosphere(filename=jdi.jupiter_pt(), sep='\s+')
/tmp/ipykernel_104973/2498042645.py:11: SyntaxWarning: invalid escape sequence '\s'
  pl.clouds(filename=jdi.jupiter_cld(), sep='\s+')

[10]:

preserve

[10]:

<xarray.Dataset> Size: 741kB
Dimensions:         (pressure: 61, wavelength: 18060, pressure_layer: 60,
                     wavenumber_layer: 196)
Coordinates:
  * pressure        (pressure) float64 488B 1e-06 1.41e-06 2e-06 ... 708.0 1e+03
  * wavelength      (wavelength) float64 144kB 1.0 0.9999 0.9999 ... 0.3 0.3 0.3
Dimensions without coordinates: pressure_layer, wavenumber_layer
Data variables: (12/56)
    temperature     (pressure) float64 488B 150.9 149.7 ... 1.525e+03 1.675e+03
    albedo          (wavelength) float64 144kB 0.2015 0.2819 ... 0.7187 0.7187
    fpfs_reflected  (wavelength) float64 144kB 4.601e-08 6.438e-08 ... 1.641e-07
    e-              (pressure) float64 488B 1e-50 1e-50 ... 4.203e-12 2.025e-11
    H2              (pressure) float64 488B 0.832 0.832 0.832 ... 0.832 0.832
    H               (pressure) float64 488B 1e-50 1e-50 ... 8.155e-07 3.361e-06
    ...              ...
    Ti              (pressure) float64 488B 1e-50 1e-50 ... 5.163e-18 2.011e-15
    Ti+             (pressure) float64 488B 1e-50 1e-50 ... 5.899e-30 4.525e-26
    C+              (pressure) float64 488B 1e-50 1e-50 ... 8.085e-50 1.678e-45
    opd             (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
    ssa             (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
    asy             (pressure_layer, wavenumber_layer) float64 94kB 0.0 ... 0.0
Attributes:
    author:          Awesome Scientist
    code:            picaso
    contact:         awesomescientist@universe.com
    planet_params:   {"gravity": {"value": 2478.651947614914, "unit": "cm / s...
    stellar_params:  {"database": "ck04models", "steff": 5000, "feh": 0, "log...
    orbit_params:    {"sma": {"value": 14959787070000.0, "unit": "cm"}}
    cloud_params:    {'fsed': 3}
[ ]: