diff --git a/doc/source/conf.py b/doc/source/conf.py
index c8d5dd7a..19df486e 100644
--- a/doc/source/conf.py
+++ b/doc/source/conf.py
@@ -79,6 +79,7 @@
     "geopandas": ("https://geopandas.org/en/stable", None),
     "xarray": ("https://docs.xarray.dev/en/stable/", None),
     "rioxarray": ("https://corteva.github.io/rioxarray/stable/", None),
+    "sliderule": ("https://slideruleearth.io/web/rtd/", None),
 }
 
 sphinx_gallery_conf = {
diff --git a/doc/source/sliderule.md b/doc/source/sliderule.md
index 45459b85..e38eb365 100644
--- a/doc/source/sliderule.md
+++ b/doc/source/sliderule.md
@@ -14,7 +14,7 @@ kernelspec:
 ---
 (cheatsheet)=
 
-# Elevation reference with SlideRule
+# Pair with reference data
 
 Most analysis of **xDEM relies on independent, high-precision elevation data to use as reference**, whether for
 coregistration, bias-corrections or uncertainty analysis.
@@ -22,7 +22,7 @@ coregistration, bias-corrections or uncertainty analysis.
 [SlideRule](https://slideruleearth.io/) provides the ideal way to retrieve such data by accessing big data archives
 of high-precision elevations, such as ICESat-2 or GEDI, efficiently in the cloud.
 
-Below, a short example to coregister and perform uncertainty analysis of a DEM with ICESat-2 ATL06.
+Below, a short example to coregister our example DEM with ICESat-2 ATL06 reference data.
 
 ```{code-cell} ipython3
 :tags: [remove-cell]
@@ -34,6 +34,10 @@ pyplot.rcParams['savefig.dpi'] = 600
 pyplot.rcParams['font.size'] = 9  # Default 10 is a bit too big for coregistration plots
 ```
 
+## Retrieving reference data
+
+We load our example DEM data, and define the region of interest from its footprint.
+
 ```{code-cell} ipython3
 from sliderule import sliderule, icesat2
 import geoutils as gu
@@ -52,7 +56,11 @@ inlier_mask = ~glacier_outlines.create_mask(dem)
 # Define region of interest as DEM footprint
 bounds = list(dem.get_bounds_projected(4326))
 region = sliderule.toregion(bounds)["poly"]
+```
 
+We can then initialize the SlideRule client, and fetch the ICESat-2 ATL06 reference data using {func}`icesat2.atl06sp`.
+
+```{code-cell} ipython3
 # Initialize SlideRule client
 sliderule.init("slideruleearth.io")
 
@@ -71,15 +79,22 @@ gdf = gdf[np.isfinite(gdf["h_li"])] # Keep valid data
 gdf = gdf[gdf["atl06_quality_summary"]==0]  # Keep very high-confidence data
 ```
 
+## Coregistration
+
+Running the coregistration is as simple as passing the geodataframe to {func}`~xdem.coreg.Coreg.fit_and_apply`,
+specifying the name of the column to use.
+
 ```{code-cell} ipython3
 # Run a translation coregistration
 nk = xdem.coreg.NuthKaab()
-aligned_dem = nk.fit_and_apply(reference_elev=gdf, to_be_aligned_elev=dem, inlier_mask=inlier_mask, z_name="h_li")
+aligned_dem = nk.fit_and_apply(reference_elev=gdf, to_be_aligned_elev=dem, inlier_mask=inlier_mask, z_name="h_li", random_state=42)
 
 # Print the estimated translation parameters
 print([k+f': {nk.meta["outputs"]["affine"][k]:.2f} meters' for k in ["shift_x", "shift_y", "shift_z"]])
 ```
 
+We can plot the improvement due to the coregistration.
+
 ```{code-cell} ipython3
 :tags: [hide-input]
 :mystnb:
@@ -96,9 +111,13 @@ pc_mask = inlier_mask.astype("uint8").interp_points((gdf.geometry.x.values, gdf.
 vect.ds = vect.ds[pc_mask == 1]
 # Get point differences before and after
 z_pc = dem.interp_points((vect.ds.geometry.x.values, vect.ds.geometry.y.values))
-vect.ds["dh_bef"] = vect.ds["h_li"] - z_pc
 z_pc_aligned = aligned_dem.interp_points((vect.ds.geometry.x.values, vect.ds.geometry.y.values))
-vect.ds["dh_aft"] = vect.ds["h_li"] - z_pc_aligned
+import warnings
+# GeoPandas raising unexpected warning
+with warnings.catch_warnings():
+    warnings.simplefilter("ignore")
+    vect["dh_bef"] = vect["h_li"] - z_pc
+    vect["dh_aft"] = vect["h_li"] - z_pc_aligned
 
 # Plot before and after
 import matplotlib.pyplot as plt
@@ -113,6 +132,8 @@ _ = ax[1].set_yticklabels([])
 plt.tight_layout()
 ```
 
+And the full information around the coregistration using {func}`~xdem.coreg.Coreg.info`.
+
 ```{code-cell} ipython3
 # Show full coregistration summary
 nk.info()