In this tutorial, we will read the Active Layer Thickness from Airborne L- and P- band SAR dataset using direct access method from NASA Earthdata Cloud, and calculate soil volumetric water content using the UAVSAR data.
This dataset provides active layer thickness (ALT) and the vertical soil moisture profile retrievals imultaneously use L- and P-band synthetic aperture radar (SAR) data acquired by the NASA/JPL Uninhabited Aerial Vehicle Synthetic Aperture Radar (UAVSAR) instruments across the ABoVE domain. The data are provided in NetCDF Version 4 format.
import earthaccess
import s3fs
import xarray as xr
import numpy as np
import warnings
import matplotlib.pyplot as plt
# suppress future warnings
warnings.filterwarnings('ignore', category=FutureWarning)Authentication¶
We recommend authenticating your Earthdata Login (EDL) information using the earthaccess python library as follows:
# works if the EDL login already been persisted to a netrc
try:
auth = earthaccess.login(strategy="netrc")
except FileNotFoundError:
# ask for EDL credentials and persist them in a .netrc file
auth = earthaccess.login(strategy="interactive", persist=True)Search granules¶
We will using the earthaccess module to search the granules within the UAVSAR dataset. W
doi = '10.3334/ORNLDAAC/2004' # uavsar data
granules = earthaccess.search_data(
granule_name = f"*.nc4", # retrieve only netcdfs
cloud_hosted=True, # make sure dataset is in cloud for direct access
doi=doi
)
print(f'Granules found: {len(granules)}')Granules found: 51
Directly access a granule¶
Let’s use earthaccess.open to open the first granule as S3FileSystem object, which mimics the standard Python file protocol, allowing you to read and write data to S3 object.
s3fs = earthaccess.open(granules[:1])
s3fs[<File-like object S3FileSystem, ornl-cumulus-prod-protected/above/ABoVE_ReSALT_InSAR_PolSAR_V3/data/PDO_ReSALT_scoaoi_2017_03.nc4>]Let’s open the above S3FileSystem object using xarray module.
ds = xr.open_dataset(s3fs[0], engine='h5netcdf', chunks={})
dsPlots of selected variables¶
Let’s plot three key variables we will use later to compute soil volumetric moisture content.
Soil water saturation fraction at the surface (z=0)¶
ds.Sw0.plot()
plt.show()
Water table depth¶
ds.wtd.plot()
plt.show()
Active layer thickness¶
ds.alt.plot()
plt.show()
Calculate soil volumetric water content (VWC)¶
Now, that we have opened a single granule above, we can compute the soil volumetric water content (VWC). The calculation will be based on this code provided with the dataset. For this tutorial, let’s compute the volumetric water content (VWC) at 0.3m soil depth.
# soil depth in meters
davg = 0.3Let’s instantiate porosity profile at various soil depths assuming a 15-cm surface organic layer.
# soil depths upto 1.5m at 0.01m interval
depth = np.arange(0, 1.51, 0.01)
# porosity profile
poros = np.zeros(depth.shape)
poros[0:60] = [0.898544856, 0.898499062, 0.898423619, 0.898299395,
0.898095018, 0.897759235, 0.897208803, 0.896309857,
0.894850646, 0.892505296, 0.888794675, 0.883063585,
0.874497473, 0.862106338, 0.844293206, 0.817876684,
0.779170531, 0.729407351, 0.676689093, 0.629853558,
0.593348919, 0.567343140, 0.549886512, 0.538609107,
0.531497790, 0.527080555, 0.524362119, 0.522698640,
0.521684243, 0.521066965, 0.520691823, 0.520464014,
0.520325740, 0.520241835, 0.520190930, 0.520160050,
0.520141318, 0.520129956, 0.520123065, 0.520118885,
0.520116349, 0.520114812, 0.520113879, 0.520113313,
0.520112970, 0.520112762, 0.520112636, 0.520112559,
0.520112513, 0.520112484, 0.520112467, 0.520112457,
0.520112451, 0.520112447, 0.520112445, 0.520112443,
0.520112442, 0.520112442, 0.520112441, 0.520112441]
poros[60:] = 0.520112441Now, compute soil volumetric water content averaged to depth davg (0.3m) and plot.
alt = ds.alt.values
Sw0 = ds.Sw0.values
wtd = ds.wtd.values
# masking the pixels where wtd is not retrieved and where alt < 0.3 soil depth
mv_avg = np.where((np.isnan(wtd) | (alt < davg)), np.nan, 0)
# calculate vwc averaged to dvag
for idx, z in enumerate(depth):
if z < davg:
Sw = np.where(wtd > z, ((Sw0 - 1) / wtd**2) * (z - wtd)**2 + 1, 1)
mv = Sw * poros[idx]
mv_avg = mv_avg + mv * 0.01
else:
break
mv_avg = mv_avg / davg
# create vwc as a new variable in dataset
ds = ds.assign(mv_avg = ds.Sw0)
ds.mv_avg[:] = np.nan
ds.mv_avg[:] = mv_avg
ds.mv_avg.attrs['long_name'] = f'Soil volumetric water content at {davg} m'
# plot
ds.mv_avg.plot()
plt.show()
- Chen, R. H., Michaelides, R. J., Chen, J., Chen, A. C., Clayton, L. K., Bakian-Dogaheh, K., Huang, L., Jafarov, E., Liu, L., Moghaddam, M., Parsekian, A. D., Sullivan, T. D., Tabatabaeenejad, A., Wig, E., Zebker, H. A., & Zhao, Y. (2022). ABoVE: Active Layer Thickness from Airborne L- and P- band SAR, Alaska, 2017, Ver. 3. ORNL Distributed Active Archive Center. 10.3334/ORNLDAAC/2004