Tutorial: creating river concentrations file for regional domain in ROMS¶
The purpose of this package is to create gridded files of river biogeochemical properties. Such properties are measured locally or inferred at the mouth of the river and provided in spreadsheet-style databases. To force ocean model, we need gridded files of BGC concentrations, that can be multiplied by river and coastal runoff to create the riverine flux into the ocean. DICRIVERS will find the closest gridpoint to the true river mouth geographical location, then create a “plume” (i.e. a binary mask) on which to apply the concentration of each given river.
[1]:
import xarray as xr
import pandas as pd
from dicrivers import make_bgc_river_input
import cartopy.crs as ccrs
import matplotlib.pyplot as plt
[2]:
%matplotlib inline
Open the included example dataframe of the 20 largest rivers:
[3]:
river_df = pd.read_csv('../dicrivers/test/data/20major_rivers.csv')
[4]:
river_df
[4]:
| basinid | basinname | mouth_lon | mouth_lat | rspread | testvar | |
|---|---|---|---|---|---|---|
| 0 | 1 | Amazon | -51.75 | -1.25 | 5 | 99 |
| 1 | 2 | Nile | 31.25 | 31.25 | 5 | 98 |
| 2 | 3 | Zaire | 12.75 | -5.75 | 5 | 97 |
| 3 | 4 | Mississippi | -90.25 | 29.75 | 5 | 96 |
| 4 | 5 | Ob | 69.25 | 66.75 | 5 | 95 |
| 5 | 6 | Parana | -58.75 | -34.25 | 5 | 94 |
| 6 | 7 | Yenisei | 82.25 | 71.25 | 5 | 93 |
| 7 | 8 | Lena | 127.25 | 73.25 | 5 | 92 |
| 8 | 9 | Niger | 6.75 | 4.75 | 5 | 91 |
| 9 | 10 | Tamanrasett | -16.25 | 20.25 | 5 | 90 |
| 10 | 11 | Chang Jiang | 121.75 | 31.25 | 5 | 89 |
| 11 | 12 | Amur | 140.75 | 53.25 | 5 | 88 |
| 12 | 13 | Mackenzie | -134.75 | 69.25 | 5 | 87 |
| 13 | 14 | Ganges | 90.75 | 22.75 | 5 | 86 |
| 14 | 15 | CHARI BOUSSO/ Lake Chad | 14.25 | 12.25 | 5 | 85 |
| 15 | 16 | Volga | 48.25 | 46.25 | 5 | 84 |
| 16 | 17 | Zambezi | 36.25 | -18.75 | 5 | 83 |
| 17 | 18 | GreatArtesian/CopperCR | 136.75 | -28.75 | 5 | 82 |
| 18 | 19 | Tarim | 86.25 | 41.25 | 5 | 81 |
| 19 | 20 | Indus | 67.75 | 24.25 | 5 | 80 |
The example dataframe contains the longitude and latitude of the river mouth, the desired radius of spreading for each river plume and values for the various biogeochemical properties. Here the radius is of 5 grid points in each direction from the river mouth. This can be adjusted on an individual basis. As a rule of thumb, the higher the horizontal resolution the larger the spreading for a plume of similar extension.
The pandas framework allows you to add variables, select and sort entries very easily. Next we need the grid onto which to generate the plumes:
[5]:
roms_regional_grid = xr.open_dataset('../dicrivers/test/data/regional_ROMS_NWA_grid.nc')
[6]:
roms_regional_grid
[6]:
<xarray.Dataset>
Dimensions: (eta_rho: 362, xi_rho: 722)
Dimensions without coordinates: eta_rho, xi_rho
Data variables:
lat_rho (eta_rho, xi_rho) float64 ...
lon_rho (eta_rho, xi_rho) float64 ...
mask_rho (eta_rho, xi_rho) float64 ...
Attributes:
Description: ROMS grid
Author: pyroms.grid.write_grd
Created: 2013-04-01T09:31:11.293024
type: ROMS grid file
history: Wed Jun 5 14:35:38 2019: ncks -F -v lon_rho,lat_rho,mask_r...
NCO: netCDF Operators version 4.7.7 (Homepage = http://nco.sf.ne...
The xarray dataset contains the longitude, latitude of tracer points (look up Arakawa staggered grid for more details) and a land sea mask. The land sea mask is used to constrain the various plumes not to contaminate adjacent estuaries.
[7]:
fig = plt.figure(figsize=[16,8])
ax = plt.axes(projection=ccrs.PlateCarree())
landseamask = roms_regional_grid['mask_rho'].assign_coords(lon=roms_regional_grid['lon_rho'],
lat=roms_regional_grid['lat_rho'])
landseamask.plot(ax=ax, x='lon', y='lat',
transform=ccrs.PlateCarree())
[7]:
<matplotlib.collections.QuadMesh at 0xd1a79f198>
Let’s make the Mississippi plume bigger so we can see it. You probably don’t have to make it so big, here just for demo.
[8]:
river_df.loc[river_df.basinname == 'Mississippi', 'rspread'] = 50
[9]:
river_df
[9]:
| basinid | basinname | mouth_lon | mouth_lat | rspread | testvar | |
|---|---|---|---|---|---|---|
| 0 | 1 | Amazon | -51.75 | -1.25 | 5 | 99 |
| 1 | 2 | Nile | 31.25 | 31.25 | 5 | 98 |
| 2 | 3 | Zaire | 12.75 | -5.75 | 5 | 97 |
| 3 | 4 | Mississippi | -90.25 | 29.75 | 50 | 96 |
| 4 | 5 | Ob | 69.25 | 66.75 | 5 | 95 |
| 5 | 6 | Parana | -58.75 | -34.25 | 5 | 94 |
| 6 | 7 | Yenisei | 82.25 | 71.25 | 5 | 93 |
| 7 | 8 | Lena | 127.25 | 73.25 | 5 | 92 |
| 8 | 9 | Niger | 6.75 | 4.75 | 5 | 91 |
| 9 | 10 | Tamanrasett | -16.25 | 20.25 | 5 | 90 |
| 10 | 11 | Chang Jiang | 121.75 | 31.25 | 5 | 89 |
| 11 | 12 | Amur | 140.75 | 53.25 | 5 | 88 |
| 12 | 13 | Mackenzie | -134.75 | 69.25 | 5 | 87 |
| 13 | 14 | Ganges | 90.75 | 22.75 | 5 | 86 |
| 14 | 15 | CHARI BOUSSO/ Lake Chad | 14.25 | 12.25 | 5 | 85 |
| 15 | 16 | Volga | 48.25 | 46.25 | 5 | 84 |
| 16 | 17 | Zambezi | 36.25 | -18.75 | 5 | 83 |
| 17 | 18 | GreatArtesian/CopperCR | 136.75 | -28.75 | 5 | 82 |
| 18 | 19 | Tarim | 86.25 | 41.25 | 5 | 81 |
| 19 | 20 | Indus | 67.75 | 24.25 | 5 | 80 |
To generate the concentration arrays, we use the make_bgc_river_input function. The arguments are the river dataframe, the list of variables we want, lon/lat/mask of the ocean grid. Optional parameters control the method used to merge the plumes (default=‘average’), the proximity range prox to apply (default=200 km) and the max number of iterations allowed nitermax for the plume spreading algo (default=1000).
[10]:
river_conc = make_bgc_river_input(river_df, ['testvar'],
roms_regional_grid['lon_rho'].values,
roms_regional_grid['lat_rho'].values,
roms_regional_grid['mask_rho'].values,
prox=200.)
river mouth at (lon,lat) = (-51.750000,-1.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (31.250000,31.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (12.750000,-5.750000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (69.250000,66.750000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (-58.750000,-34.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (82.250000,71.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (127.250000,73.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (6.750000,4.750000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (-16.250000,20.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (121.750000,31.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (140.750000,53.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (-134.750000,69.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (90.750000,22.750000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (14.250000,12.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (48.250000,46.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (36.250000,-18.750000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (136.750000,-28.750000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (86.250000,41.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
river mouth at (lon,lat) = (67.750000,24.250000)
is too far from the ocean and cannot not used.
It can be out of the domain (regional case) or flowing into
an unresolved lake or sea. If you know that river should be
there, you may need to increase the proximity range (prox)
The obtained xarray dataset contain the concentrations for all selected variables:
[11]:
river_conc
[11]:
<xarray.Dataset>
Dimensions: (x: 722, y: 362)
Coordinates:
lat (y, x) float64 8.785 8.826 8.866 8.907 ... 53.74 53.77 53.8 53.83
lon (y, x) float64 275.5 275.6 275.6 275.7 ... 301.1 301.1 301.2 301.3
Dimensions without coordinates: x, y
Data variables:
testvar (y, x) float64 0.0 0.0 0.0 0.0 0.0 0.0 ... 0.0 0.0 0.0 0.0 0.0 0.0
Let’s make a visual check:
[12]:
fig = plt.figure(figsize=[16,8])
ax = plt.axes(projection=ccrs.PlateCarree())
river_conc['testvar'].where(river_conc['testvar'] != 0).plot(ax=ax, x='lon', y='lat', transform=ccrs.PlateCarree())
ax.coastlines()
[12]:
<cartopy.mpl.feature_artist.FeatureArtist at 0xd1bb13f28>
We can finally save the end result to a netcdf file:
[13]:
river_conc.to_netcdf('dicrivers_regional.nc')
In realistic applications, you can end up putting dozens if not hundreds of rivers, streams,… The rspread and prox parameters will help you tune the results. The final word is that your critical look at the results is something no computer yet can do, so open your eyes, be thorough and report any bugs you can find and propose improvements on the github page of DICRIVERS. Keep on rockin’ !!