ArcticDEM Mosaic Example¶

Purpose¶

Demonstrate how to sample the ArcticDEM at generated ATL06-SR points

Import Packages¶

In [1]:
import warnings
warnings.filterwarnings("ignore") # suppress warnings
In [2]:
import matplotlib.pyplot as plt
import matplotlib
import sliderule
from sliderule import icesat2

Initialize SlideRule Python Client¶

In [3]:
icesat2.init("slideruleearth.io", organization="bathy", verbose=True)

Make Processing Request to SlideRule¶

ATL06-SR request includes the samples parameter to specify that ArcticDEM Mosiac dataset should be sampled at each generated ATL06 elevation.

In [4]:
resource = "ATL03_20190314093716_11600203_005_01.h5"
region = sliderule.toregion("../data/dicksonfjord.geojson")
parms = { "poly": region['poly'],
          "cnf": "atl03_high",
          "ats": 5.0,
          "cnt": 5,
          "len": 20.0,
          "res": 10.0,
          "samples": {"mosaic": {"asset": "arcticdem-mosaic", "radius": 10.0, "zonal_stats": True}} }
gdf = icesat2.atl06p(parms, resources=[resource])

Display GeoDataFrame¶

Notice the columns that start with "mosaic"

In [5]:
gdf
Out[5]:
w_surface_window_final segment_id n_fit_photons spot gt dh_fit_dx rgt cycle region h_sigma ... mosaic.min mosaic.flags mosaic.median mosaic.value mosaic.stdev mosaic.mad mosaic.max mosaic.mean mosaic.count mosaic.time
time
2019-03-14 09:40:46.291445504 3.000000 405229 15 1 10 0.082324 1160 2 3 0.187902 ... 576.406250 0 584.289062 584.289062 4.606139 3.984568 592.898438 584.500000 81 1.358109e+12
2019-03-14 09:40:46.292858624 3.000000 405230 27 1 10 0.074614 1160 2 3 0.059746 ... 583.507812 0 592.257812 592.562500 3.477136 2.955249 596.804688 591.543306 81 1.358109e+12
2019-03-14 09:40:46.294273792 3.000000 405230 12 1 10 0.045681 1160 2 3 0.328370 ... 593.132812 0 596.843750 596.875000 2.326261 1.928062 602.914062 597.138600 81 1.358109e+12
2019-03-14 09:40:46.297095936 7.430891 405231 21 1 10 0.246724 1160 2 3 0.389433 ... 599.804688 0 605.265625 605.125000 2.346747 2.013160 608.570312 604.967496 81 1.358109e+12
2019-03-14 09:40:46.298507520 14.818856 405232 36 1 10 0.580637 1160 2 3 0.321795 ... 603.515625 0 608.179688 608.117188 2.340230 1.958574 612.867188 608.292921 81 1.358109e+12
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
2019-03-14 09:40:48.353273600 3.000000 405838 27 2 20 -0.021866 1160 2 3 0.053913 ... 1493.515625 0 1494.429688 1494.398438 0.497044 0.418148 1495.468750 1494.448110 81 1.358109e+12
2019-03-14 09:40:48.354683648 3.000000 405839 15 2 20 -0.028824 1160 2 3 0.114113 ... 1492.945312 0 1494.156250 1494.218750 0.559426 0.459164 1495.273438 1494.151427 81 1.358109e+12
2019-03-14 09:40:48.356088576 3.000000 405839 8 2 20 -0.112578 1160 2 3 0.046961 ... 1491.921875 0 1493.507812 1493.546875 0.750541 0.633445 1494.906250 1493.483507 81 1.358109e+12
2019-03-14 09:40:48.357494784 3.000000 405840 16 2 20 -0.106287 1160 2 3 0.032820 ... 1490.906250 0 1492.648438 1492.648438 0.844155 0.711853 1494.187500 1492.628954 81 1.358109e+12
2019-03-14 09:40:48.358897920 3.000000 405840 19 2 20 -0.078550 1160 2 3 0.028111 ... 1489.742188 0 1491.515625 1491.367188 0.878481 0.747692 1493.062500 1491.452932 81 1.358109e+12

1651 rows × 27 columns

Print Out File Directory¶

When a GeoDataFrame includes samples from rasters, each sample value has a file id that is used to look up the file name of the source raster for that value.

In [6]:
gdf.attrs['file_directory']
Out[6]:
{0: '/vsis3/pgc-opendata-dems/arcticdem/mosaics/v4.1/2m_dem_tiles.vrt'}

Demonstrate How To Access Source Raster Filename for Entry in GeoDataFrame¶

In [7]:
filedir = gdf.attrs['file_directory']
filedir[gdf['mosaic.file_id'][0]]
Out[7]:
'/vsis3/pgc-opendata-dems/arcticdem/mosaics/v4.1/2m_dem_tiles.vrt'

Difference the Sampled Value from ArcticDEM with SlideRule ATL06-SR¶

In [8]:
gdf["value_delta"] = gdf["h_mean"] - gdf["mosaic.value"]
gdf["value_delta"].describe()
Out[8]:
count    1651.000000
mean        1.110096
std         2.145814
min       -27.810532
25%         0.201006
50%         1.270788
75%         2.169877
max        14.924781
Name: value_delta, dtype: float64

Difference the Zonal Statistic Mean from ArcticDEM with SlideRule ATL06-SR¶

In [9]:
gdf["mean_delta"] = gdf["h_mean"] - gdf["mosaic.mean"]
gdf["mean_delta"].describe()
Out[9]:
count    1651.000000
mean        1.117867
std         2.146371
min       -28.240605
25%         0.171224
50%         1.285451
75%         2.180714
max        15.173141
Name: mean_delta, dtype: float64

Difference the Zonal Statistic Mdeian from ArcticDEM with SlideRule ATL06-SR¶

In [10]:
gdf["median_delta"] = gdf["h_mean"] - gdf["mosaic.median"]
gdf["median_delta"].describe()
Out[10]:
count    1651.000000
mean        1.115779
std         2.142610
min       -27.716782
25%         0.209473
50%         1.277488
75%         2.166412
max        14.924781
Name: median_delta, dtype: float64

Plot the Different ArcticDEM Values against the SlideRule ATL06-SR Values¶

In [11]:
# Setup Plot
fig,ax = plt.subplots(num=None, figsize=(10, 8))
fig.set_facecolor('white')
fig.canvas.header_visible = False
ax.set_title("SlideRule vs. ArcticDEM Elevations")
ax.set_xlabel('UTC')
ax.set_ylabel('height (m)')
legend_elements = []

# Plot SlideRule ATL06 Elevations
df = gdf[(gdf['rgt'] == 1160) & (gdf['gt'] == 10) & (gdf['cycle'] == 2)]
sc1 = ax.scatter(df.index.values, df["h_mean"].values, c='red', s=2.5)
legend_elements.append(matplotlib.lines.Line2D([0], [0], color='red', lw=6, label='ATL06-SR'))

# Plot ArcticDEM Elevations
sc2 = ax.scatter(df.index.values, df["mosaic.value"].values, c='blue', s=2.5)
legend_elements.append(matplotlib.lines.Line2D([0], [0], color='blue', lw=6, label='ArcticDEM'))

# Display Legend
lgd = ax.legend(handles=legend_elements, loc=3, frameon=True)
lgd.get_frame().set_alpha(1.0)
lgd.get_frame().set_edgecolor('white')

# Show Plot
plt.show()
No description has been provided for this image

Plot the Sampled Value and Zonal Statistic Mean Deltas to SlideRule ATL06-SR Values¶

In [12]:
# Setup Plot
fig,ax = plt.subplots(num=None, figsize=(10, 8))
fig.set_facecolor('white')
fig.canvas.header_visible = False
ax.set_title("Delta Elevations between SlideRule and ArcticDEM")
ax.set_xlabel('UTC')
ax.set_ylabel('height (m)')
ax.yaxis.grid(True)

# Plot Deltas
df1 = gdf[(gdf['rgt'] == 1160) & (gdf['gt'] == 10) & (gdf['cycle'] == 2)]
sc1 = ax.scatter(df1.index.values, df1["value_delta"].values, c='blue', s=2.5)

# Plot Deltas
df2 = gdf[(gdf['rgt'] == 1160) & (gdf['gt'] == 10) & (gdf['cycle'] == 2)]
sc2 = ax.scatter(df2.index.values, df2["mean_delta"].values, c='green', s=2.5)

# Show Plot
plt.show()
No description has been provided for this image
In [ ]: