COS_Spatial_Analyses/plot_bounds.py at master · Timothy-W-Hilton/COS_Spatial_Analyses · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
"""plot color map image of Jul/Aug climatological boundary [COS]

horizontal axis: boundary ring cell number (`see models-3 I/O API
documentation
<https://www.cmascenter.org/ioapi/documentation/3.1/html/DATATYPES.html#bndary>`_)
vertical axis: meters above ground level

Uses pre-computed July/August climatological boundary [COS] calculated
by climatological_bounds_vertprofile.py.  These [COS] data are assumed
to be in
./ClimatologicalBounds/climatological_COS_bdy_22levs_124x124.nc
"""

import matplotlib
matplotlib.use('AGG')

import os
import os.path

import matplotlib.pyplot as plt
import netCDF4
import numpy as np
from timutils import colormap_nlevs
from stem_pytools.na_map import NAMapFigure
from stem_pytools import domain

if __name__ == "__main__":

    # plot the boundary ring cell number vs. vertical cell number
    nc = netCDF4.Dataset(
        os.path.join(os.getcwd(),
                     'ClimatologicalBounds',
                     'climatological_COS_bdy_22levs_124x124.nc'))
    ppb_2_ppt = 1e3  # unit conversion factor
    cos = nc.variables['CO2_TRACER1'][:].squeeze() * ppb_2_ppt

    cmap, norm = colormap_nlevs.setup_colormap(vmin=cos.min() - 1,
                                               vmax=cos.max() + 1,
                                               nlevs=20,
                                               cmap=plt.get_cmap('Blues'),
                                               extend='neither')

    d = domain.STEM_Domain()
    d.get_STEMZ_height(wrfheight_fname=os.path.join(os.getenv('SARIKA_INPUT'),
                                                    'wrfheight-124x124-22levs.nc'))
    agl_perim = np.array([domain.get_2d_perimeter(d.agl[z, ...]).mean()
                          for z in range(22)])

    fontsz = 14
    matplotlib.rcParams.update({'font.size': fontsz})
    fig, ax = plt.subplots(nrows=1, ncols=1, figsize=(7, 6))
    cm = ax.pcolormesh(np.arange(500), agl_perim, cos,
                       cmap=cmap, norm=norm,
                       linewidth=0, rasterized=True)
    ax.set_xlim([0, 500])
    ax.set_ylim([agl_perim.min(), agl_perim.max()])
    ax.set_ylabel('meters above ground', fontdict={'fontsize': fontsz})
    ax.set_xlabel('lateral boundary index', fontdict={'fontsize': fontsz})
    ax_cb = plt.colorbar(cm, cmap=cmap, norm=norm, ax=ax)
    ax_cb.set_label('[COS] (ppt)', fontdict={'fontsize': fontsz})
    ax_cb.solids.set_rasterized(True)

    fig.tight_layout()
    fig.savefig(os.path.join(os.getenv('HOME'),
                             'plots',
                             'climatological_bounds.pdf'))
    plt.close(fig)