Examples¶
aospy comes with some sample data files, which can be used to illustrate some of its basic features. These files contain monthly mean time series output of two variables (the large-scale and convective components of the total precipitation rate) from an idealized aquaplanet climate model. A simple computation one could seek to do from this model output would be to compute some statistics of the total precipitation rate (large-scale plus convective).
Here’s a quick summary of the included data:
In [1]: import xarray as xr
In [2]: xr.open_mfdataset('../aospy/test/data/netcdf/000[4-6]0101.precip_monthly.nc',
...: decode_times=False)
...:
Out[2]:
<xarray.Dataset>
Dimensions: (lat: 64, latb: 65, lon: 128, lonb: 129, nv: 2, time: 36)
Coordinates:
* lonb (lonb) float64 -1.406 1.406 4.219 7.031 9.844 12.66 ...
* lon (lon) float64 0.0 2.812 5.625 8.438 11.25 14.06 16.88 ...
* latb (latb) float64 -90.0 -86.58 -83.76 -80.96 -78.16 ...
* lat (lat) float64 -87.86 -85.1 -82.31 -79.53 -76.74 ...
* nv (nv) float64 1.0 2.0
* time (time) float64 1.111e+03 1.139e+03 1.17e+03 1.2e+03 ...
Data variables:
time_bounds (time, nv) float64 1.095e+03 1.126e+03 1.126e+03 ...
condensation_rain (time, lat, lon) float64 5.768e-06 5.784e-06 ...
convection_rain (time, lat, lon) float64 0.0 0.0 0.0 0.0 0.0 0.0 0.0 ...
average_DT (time) float64 31.0 28.0 31.0 30.0 31.0 30.0 31.0 ...
In this particular model, the large-scale component of the precipitation rate is called “condensation_rain” and the convective component is called “convection_rain.”
Defining the simulation metadata¶
The first step is to create an aospy.Run
object that
stores metadata about this simulation. This includes giving it a
name, a description, and specifying where files are located through a
DataLoader.
In [3]: from datetime import datetime
In [4]: from aospy import Run
In [5]: from aospy.data_loader import DictDataLoader
In [6]: file_map = {'monthly': '../aospy/test/data/netcdf/000[4-6]0101.precip_monthly.nc'}
In [7]: example_run = Run(
...: name='example_run',
...: description=(
...: 'Control simulation of the idealized moist model'
...: ),
...: data_loader=DictDataLoader(file_map)
...: )
...:
We then need to associate this Run
with an aospy.Model
object:
In [8]: from aospy import Model
In [9]: example_model = Model(
...: name='example_model',
...: grid_file_paths=(
...: '../aospy/test/data/netcdf/00040101.precip_monthly.nc',
...: '../aospy/test/data/netcdf/im.landmask.nc'
...: ),
...: runs=[example_run]
...: )
...:
Finally, we need to associate the Model
object with an
aospy.Proj
object. Here we can specify the location that
aospy will save its output files.
In [10]: from aospy import Proj
In [11]: example_proj = Proj(
....: 'example_proj',
....: direc_out='example-output',
....: tar_direc_out='example-tar-output',
....: models=(example_model,)
....: )
....:
Now the metadata associated with this simulation is fully defined. We can move on to computing the total precipitation.
Computing the annual mean total precipitation rate¶
We can start by defining a simple python function that computes the total precipitation from condensation and convection rain arguments:
In [12]: def total_precipitation(condensation_rain, convection_rain):
....: return condensation_rain + convection_rain
....:
To hook this function into the aospy framework, we need to connect it
to an aospy.Var
object, as well as define the Var
objects it depends on (variables that are natively stored in model
output files).
In [13]: from aospy import Var
In [14]: condensation_rain = Var(
....: name='condensation_rain',
....: alt_names=('prec_ls',),
....: def_time=True,
....: description=('condensation rain'),
....: )
....:
In [15]: convection_rain = Var(
....: name='convection_rain',
....: alt_names=('prec_conv',),
....: def_time=True,
....: description=('convection rain'),
....: )
....:
In [16]: precip = Var(
....: name='total_precipitation',
....: def_time=True,
....: description=('total precipitation rate'),
....: func=total_precipitation,
....: variables=(condensation_rain, convection_rain)
....: )
....:
Here the func attribute of the precip Var
object is the function
we defined, and the variables attribute is a tuple containing the
Var
objects the function depends on, in the order of the
function’s call signature.
If we’d like to compute the time-mean total precipitation rate from
year four to year six using aospy, we can create an
aospy.Calc
object. This is currently done through passing
an aospy.CalcInterface
object to a Calc
object; once
created, the computation can be submitted by simply calling the
compute function of Calc
.
In [17]: from aospy import CalcInterface, Calc
In [18]: calc_int = CalcInterface(
....: proj=example_proj,
....: model=example_model,
....: run=example_run,
....: var=precip,
....: date_range=(datetime(4, 1, 1), datetime(6, 12, 31)),
....: intvl_in='monthly',
....: dtype_in_time='ts',
....: intvl_out='ann',
....: dtype_out_time='av'
....: )
....:
In [19]: Calc(calc_int).compute()
The result is stored in a netcdf file, whose path and filename contains metadata about where it came from:
In [20]: calc_int.path_out['av']
Out[20]: 'example-output/example_proj/example_model/example_run/total_precipitation/total_precipitation.ann.av.from_monthly_ts.example_model.example_run.0004-0006.nc'
Using xarray we can open and plot the results of the calculation:
In [21]: xr.open_dataset(calc_int.path_out['av']).total_precipitation.plot()
Out[21]: <matplotlib.collections.QuadMesh at 0x7fc18a4d2e10>
Computing the global annual mean total precipitation rate¶
Not only does aospy enable reductions along the time dimension, it
also enables area weighted regional averages. As a simple
introduction, we’ll show how to compute the global mean total
precipitation rate from this Run
. To do so, we’ll make use of the
infrastructure defined above, and also define an
aospy.Region
object:
In [22]: from aospy import Region
In [23]: globe = Region(
....: name='globe',
....: description='Entire globe',
....: lat_bounds=(-90, 90),
....: lon_bounds=(0, 360),
....: do_land_mask=False
....: )
....:
To compute the global annual mean total precipitation rate, we can now create
another Calc
object:
In [24]: calc_int = CalcInterface(
....: proj=example_proj,
....: model=example_model,
....: run=example_run,
....: var=precip,
....: date_range=(datetime(4, 1, 1), datetime(6, 12, 31)),
....: intvl_in='monthly',
....: dtype_in_time='ts',
....: intvl_out='ann',
....: dtype_out_time='reg.av',
....: region={'globe': globe}
....: )
....:
In [25]: Calc(calc_int).compute()
This produces a new file, located in:
In [26]: calc_int.path_out['reg.av']
Out[26]: 'example-output/example_proj/example_model/example_run/total_precipitation/total_precipitation.ann.reg.av.from_monthly_ts.example_model.example_run.0004-0006.nc'
We find that the global annual mean total precipitation rate for this run (converting to units of mm per day) is:
In [27]: xr.open_dataset(calc_int.path_out['reg.av']).globe * 86400.
Out[27]:
<xarray.DataArray 'globe' ()>
array(3.025191320965487)
Coordinates:
avg_start_date datetime64[ns] 1678-01-01
avg_end_date datetime64[ns] 1681-01-01