Loading the ECCOv4 native model grid parameters

Objectives

Briefly show how to open the ECCOv4 grid file as an xarray dataset, plot parameters, and create subsets of the grid dataset.

Introduction

The ECCOv4 model grid parameters are provided as a single NetCDF file. It can be downloaded using the ecco_podaac_download function as described in the download tutorial. The ShortName for the dataset is ECCO_L4_GEOMETRY_LLC0090GRID_V4R4. The grid parameters file has no time dimension, but ecco_podaac_download requires a StartDate and EndDate to be specified; any date in the range 1992-2017 can be used.

Load the ECCOv4 model grid parameter NetCDF file

Because the ECCOv4 model grid parameter data is provided in a single file you can use the open_dataset routine from xarray to open it.

First set up your environment.

[1]:
import numpy as np
import xarray as xr
import sys
import matplotlib.pyplot as plt
%matplotlib inline
[2]:
## Import the ecco_v4_py library into Python
## =========================================

## -- If ecco_v4_py is not installed in your local Python library,
##    tell Python where to find it.  For example, if your ecco_v4_py
##    files are in /Users/ifenty/ECCOv4-py/ecco_v4_py, then use:

from os.path import expanduser,join
import sys
user_home_dir = expanduser('~')
sys.path.append(join(user_home_dir,'ECCOv4-py'))
import ecco_v4_py as ecco
[3]:
## Set top-level file directory for the ECCO NetCDF files

## change ECCO_dir as needed
ECCO_dir = join(user_home_dir,'Downloads','ECCO_V4r4_PODAAC')
[4]:
# grid parameter file name and path (after it has been downloaded)
grid_params_shortname = "ECCO_L4_GEOMETRY_LLC0090GRID_V4R4"
grid_params_file = "GRID_GEOMETRY_ECCO_V4r4_native_llc0090.nc"
grid_params_file_path = join(ECCO_dir,grid_params_shortname,grid_params_file)

# open grid parameters file
grid = xr.open_dataset(grid_params_file_path)

# show contents of grid_dataset
grid
[4]:
<xarray.Dataset>
Dimensions:  (i: 90, i_g: 90, j: 90, j_g: 90, k: 50, k_u: 50, k_l: 50, k_p1: 51, tile: 13, nb: 4, nv: 2)
Coordinates: (12/20)
  * i        (i) int32 0 1 2 3 4 5 6 7 8 9 10 ... 80 81 82 83 84 85 86 87 88 89
  * i_g      (i_g) int32 0 1 2 3 4 5 6 7 8 9 ... 80 81 82 83 84 85 86 87 88 89
  * j        (j) int32 0 1 2 3 4 5 6 7 8 9 10 ... 80 81 82 83 84 85 86 87 88 89
  * j_g      (j_g) int32 0 1 2 3 4 5 6 7 8 9 ... 80 81 82 83 84 85 86 87 88 89
  * k        (k) int32 0 1 2 3 4 5 6 7 8 9 10 ... 40 41 42 43 44 45 46 47 48 49
  * k_u      (k_u) int32 0 1 2 3 4 5 6 7 8 9 ... 40 41 42 43 44 45 46 47 48 49
    ...       ...
    Zp1      (k_p1) float32 0.0 -10.0 -20.0 ... -5.244e+03 -5.678e+03 -6.134e+03
    Zu       (k_u) float32 -10.0 -20.0 -30.0 ... -5.678e+03 -6.134e+03
    Zl       (k_l) float32 0.0 -10.0 -20.0 ... -4.834e+03 -5.244e+03 -5.678e+03
    XC_bnds  (tile, j, i, nb) float32 ...
    YC_bnds  (tile, j, i, nb) float32 ...
    Z_bnds   (k, nv) float32 0.0 -10.0 -10.0 ... -5.678e+03 -6.134e+03
Dimensions without coordinates: nb, nv
Data variables: (12/21)
    CS       (tile, j, i) float32 ...
    SN       (tile, j, i) float32 ...
    rA       (tile, j, i) float32 ...
    dxG      (tile, j_g, i) float32 ...
    dyG      (tile, j, i_g) float32 ...
    Depth    (tile, j, i) float32 ...
    ...       ...
    hFacC    (k, tile, j, i) float32 ...
    hFacW    (k, tile, j, i_g) float32 ...
    hFacS    (k, tile, j_g, i) float32 ...
    maskC    (k, tile, j, i) bool ...
    maskW    (k, tile, j, i_g) bool ...
    maskS    (k, tile, j_g, i) bool ...
Attributes: (12/58)
    acknowledgement:                 This research was carried out by the Jet...
    author:                          Ian Fenty and Ou Wang
    cdm_data_type:                   Grid
    comment:                         Fields provided on the curvilinear lat-l...
    Conventions:                     CF-1.8, ACDD-1.3
    coordinates_comment:             Note: the global 'coordinates' attribute...
    ...                              ...
    references:                      ECCO Consortium, Fukumori, I., Wang, O.,...
    source:                          The ECCO V4r4 state estimate was produce...
    standard_name_vocabulary:        NetCDF Climate and Forecast (CF) Metadat...
    summary:                         This dataset provides geometric paramete...
    title:                           ECCO Geometry Parameters for the Lat-Lon...
    uuid:                            87ff7d24-86e5-11eb-9c5f-f8f21e2ee3e0

Plot grid parameters using plot_tiles function

Let’s plot two of the model grid parameter fields hFacC (tracer cell thickness factor) and rA (grid cell surface area)

First we plot hFac:

[5]:
ecco.plot_tiles(grid.hFacC.isel(k=0), cmap='gray', show_colorbar=True,);
_images/ECCO_v4_Loading_the_ECCOv4_native_model_grid_parameters_7_0.png
[6]:
ecco.plot_tiles(grid.rA, cmap='jet', show_colorbar=True);
plt.suptitle('Model grid cell surface area [m^2]')

plt.show()
_images/ECCO_v4_Loading_the_ECCOv4_native_model_grid_parameters_8_0.png

Plot subset of global domain

Once the file has been “opened” using open_dataset, we can use Dataset.isel to subset the file if we don’t want to plot the full global domain.

[7]:
grid_subset = grid.isel(tile=[1,10,12],k=[0,1,2,3])
grid_subset
[7]:
<xarray.Dataset>
Dimensions:  (i: 90, i_g: 90, j: 90, j_g: 90, k: 4, k_u: 50, k_l: 50, k_p1: 51, tile: 3, nb: 4, nv: 2)
Coordinates: (12/20)
  * i        (i) int32 0 1 2 3 4 5 6 7 8 9 10 ... 80 81 82 83 84 85 86 87 88 89
  * i_g      (i_g) int32 0 1 2 3 4 5 6 7 8 9 ... 80 81 82 83 84 85 86 87 88 89
  * j        (j) int32 0 1 2 3 4 5 6 7 8 9 10 ... 80 81 82 83 84 85 86 87 88 89
  * j_g      (j_g) int32 0 1 2 3 4 5 6 7 8 9 ... 80 81 82 83 84 85 86 87 88 89
  * k        (k) int32 0 1 2 3
  * k_u      (k_u) int32 0 1 2 3 4 5 6 7 8 9 ... 40 41 42 43 44 45 46 47 48 49
    ...       ...
    Zp1      (k_p1) float32 0.0 -10.0 -20.0 ... -5.244e+03 -5.678e+03 -6.134e+03
    Zu       (k_u) float32 -10.0 -20.0 -30.0 ... -5.678e+03 -6.134e+03
    Zl       (k_l) float32 0.0 -10.0 -20.0 ... -4.834e+03 -5.244e+03 -5.678e+03
    XC_bnds  (tile, j, i, nb) float32 -38.0 -37.0 -37.0 ... -115.0 -115.0 -108.5
    YC_bnds  (tile, j, i, nb) float32 -57.01 -57.01 -56.47 ... -88.18 -88.16
    Z_bnds   (k, nv) float32 0.0 -10.0 -10.0 -20.0 -20.0 -30.0 -30.0 -40.0
Dimensions without coordinates: nb, nv
Data variables: (12/21)
    CS       (tile, j, i) float32 1.0 1.0 1.0 1.0 ... -0.9601 -0.9854 -0.9984
    SN       (tile, j, i) float32 -0.0 6.524e-15 -6.524e-15 ... -0.1705 -0.05718
    rA       (tile, j, i) float32 3.625e+09 3.625e+09 ... 3.685e+08 3.611e+08
    dxG      (tile, j_g, i) float32 6.054e+04 6.054e+04 ... 2.36e+04 2.314e+04
    dyG      (tile, j, i_g) float32 5.944e+04 5.944e+04 ... 1.56e+04 1.558e+04
    Depth    (tile, j, i) float32 3.284e+03 3.486e+03 3.486e+03 ... 0.0 0.0 0.0
    ...       ...
    hFacC    (k, tile, j, i) float32 1.0 1.0 1.0 1.0 1.0 ... 0.0 0.0 0.0 0.0 0.0
    hFacW    (k, tile, j, i_g) float32 1.0 1.0 1.0 1.0 1.0 ... 0.0 0.0 0.0 0.0
    hFacS    (k, tile, j_g, i) float32 1.0 1.0 1.0 1.0 1.0 ... 0.0 0.0 0.0 0.0
    maskC    (k, tile, j, i) bool True True True True ... False False False
    maskW    (k, tile, j, i_g) bool True True True True ... False False False
    maskS    (k, tile, j_g, i) bool True True True True ... False False False
Attributes: (12/58)
    acknowledgement:                 This research was carried out by the Jet...
    author:                          Ian Fenty and Ou Wang
    cdm_data_type:                   Grid
    comment:                         Fields provided on the curvilinear lat-l...
    Conventions:                     CF-1.8, ACDD-1.3
    coordinates_comment:             Note: the global 'coordinates' attribute...
    ...                              ...
    references:                      ECCO Consortium, Fukumori, I., Wang, O.,...
    source:                          The ECCO V4r4 state estimate was produce...
    standard_name_vocabulary:        NetCDF Climate and Forecast (CF) Metadat...
    summary:                         This dataset provides geometric paramete...
    title:                           ECCO Geometry Parameters for the Lat-Lon...
    uuid:                            87ff7d24-86e5-11eb-9c5f-f8f21e2ee3e0

Notice that grid_subset only has 3 tiles (1, 10, 12) and 4 depth levels (0, 1, 2, 3), as expected. The subset indices can also be passed to isel as a Python dictionary, with the same result.

[8]:
subset_ind = {'tile':[1,10,12],'k':[0,1,2,3]}

grid_subset = grid.isel(subset_ind)
grid_subset
[8]:
<xarray.Dataset>
Dimensions:  (i: 90, i_g: 90, j: 90, j_g: 90, k: 4, k_u: 50, k_l: 50, k_p1: 51, tile: 3, nb: 4, nv: 2)
Coordinates: (12/20)
  * i        (i) int32 0 1 2 3 4 5 6 7 8 9 10 ... 80 81 82 83 84 85 86 87 88 89
  * i_g      (i_g) int32 0 1 2 3 4 5 6 7 8 9 ... 80 81 82 83 84 85 86 87 88 89
  * j        (j) int32 0 1 2 3 4 5 6 7 8 9 10 ... 80 81 82 83 84 85 86 87 88 89
  * j_g      (j_g) int32 0 1 2 3 4 5 6 7 8 9 ... 80 81 82 83 84 85 86 87 88 89
  * k        (k) int32 0 1 2 3
  * k_u      (k_u) int32 0 1 2 3 4 5 6 7 8 9 ... 40 41 42 43 44 45 46 47 48 49
    ...       ...
    Zp1      (k_p1) float32 0.0 -10.0 -20.0 ... -5.244e+03 -5.678e+03 -6.134e+03
    Zu       (k_u) float32 -10.0 -20.0 -30.0 ... -5.678e+03 -6.134e+03
    Zl       (k_l) float32 0.0 -10.0 -20.0 ... -4.834e+03 -5.244e+03 -5.678e+03
    XC_bnds  (tile, j, i, nb) float32 -38.0 -37.0 -37.0 ... -115.0 -115.0 -108.5
    YC_bnds  (tile, j, i, nb) float32 -57.01 -57.01 -56.47 ... -88.18 -88.16
    Z_bnds   (k, nv) float32 0.0 -10.0 -10.0 -20.0 -20.0 -30.0 -30.0 -40.0
Dimensions without coordinates: nb, nv
Data variables: (12/21)
    CS       (tile, j, i) float32 1.0 1.0 1.0 1.0 ... -0.9601 -0.9854 -0.9984
    SN       (tile, j, i) float32 -0.0 6.524e-15 -6.524e-15 ... -0.1705 -0.05718
    rA       (tile, j, i) float32 3.625e+09 3.625e+09 ... 3.685e+08 3.611e+08
    dxG      (tile, j_g, i) float32 6.054e+04 6.054e+04 ... 2.36e+04 2.314e+04
    dyG      (tile, j, i_g) float32 5.944e+04 5.944e+04 ... 1.56e+04 1.558e+04
    Depth    (tile, j, i) float32 3.284e+03 3.486e+03 3.486e+03 ... 0.0 0.0 0.0
    ...       ...
    hFacC    (k, tile, j, i) float32 1.0 1.0 1.0 1.0 1.0 ... 0.0 0.0 0.0 0.0 0.0
    hFacW    (k, tile, j, i_g) float32 1.0 1.0 1.0 1.0 1.0 ... 0.0 0.0 0.0 0.0
    hFacS    (k, tile, j_g, i) float32 1.0 1.0 1.0 1.0 1.0 ... 0.0 0.0 0.0 0.0
    maskC    (k, tile, j, i) bool True True True True ... False False False
    maskW    (k, tile, j, i_g) bool True True True True ... False False False
    maskS    (k, tile, j_g, i) bool True True True True ... False False False
Attributes: (12/58)
    acknowledgement:                 This research was carried out by the Jet...
    author:                          Ian Fenty and Ou Wang
    cdm_data_type:                   Grid
    comment:                         Fields provided on the curvilinear lat-l...
    Conventions:                     CF-1.8, ACDD-1.3
    coordinates_comment:             Note: the global 'coordinates' attribute...
    ...                              ...
    references:                      ECCO Consortium, Fukumori, I., Wang, O.,...
    source:                          The ECCO V4r4 state estimate was produce...
    standard_name_vocabulary:        NetCDF Climate and Forecast (CF) Metadat...
    summary:                         This dataset provides geometric paramete...
    title:                           ECCO Geometry Parameters for the Lat-Lon...
    uuid:                            87ff7d24-86e5-11eb-9c5f-f8f21e2ee3e0

Let’s plot hFacC and rA again

[9]:
ecco.plot_tiles(grid_subset.hFacC.isel(k=0), cmap='gray', show_colorbar=True,);
'Model grid cell surface area [m^2] in tiles 1, 10, and 12 '
[9]:
'Model grid cell surface area [m^2] in tiles 1, 10, and 12 '
_images/ECCO_v4_Loading_the_ECCOv4_native_model_grid_parameters_14_1.png

Notice that 10 of the 13 tiles are blank because they were not loaded.

[10]:
ecco.plot_tiles(grid_subset.rA, cmap='jet', show_colorbar=True);
'Model grid cell surface area [m^2]'
[10]:
'Model grid cell surface area [m^2]'
_images/ECCO_v4_Loading_the_ECCOv4_native_model_grid_parameters_16_1.png