Skip to content

Commit ab39960

Browse files
timothypratleytimothypratley
authored
Merge pull request #87 from larzeitlin/water-ndwi
notebook: Remote sensing on water
2 parents aaeefad + 566692e commit ab39960

File tree

2 files changed

+202
-1
lines changed

2 files changed

+202
-1
lines changed

‎site/db.edn‎

Lines changed: 4 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -121,7 +121,10 @@
121121
:affiliation [:cim]}
122122
{:id :tanvi
123123
:name "Tanvi Nagar"
124-
:affiliation [:cim]}]
124+
:affiliation [:cim]}
125+
{:id :luke-zeitlin
126+
:name "Luke Zeitlin"
127+
:affiliation []}]
125128

126129
:affiliation
127130
[{:id :clojure.core

‎src/earth_observation/waterquality.clj‎

Lines changed: 198 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,198 @@
1+
^{:kindly/hide-code true
2+
:clay {:title "Remote Sensing - Water Quality"
3+
:quarto {:author :luke-zeitlin
4+
:type :post
5+
:date "2025-09-26"
6+
:category :clojure
7+
:tags [:clay :reagent :scittle
8+
:data-science :remote-sensing]}}}
9+
10+
11+
(ns earth-observation.waterquality
12+
(:require [scicloj.kindly.v4.kind :as kind]))
13+
14+
;; # Remote sensing water quality analysis
15+
16+
;; This notebook walks through the steps required to create an
17+
;; interactive map widget displaying simple band-ratio algorithms for
18+
;; water quality analysis using cloud-optimized GeoTIFFs.
19+
;; This example will show how to use
20+
;; [Normalized Difference Water Index](https://en.wikipedia.org/wiki/Normalized_difference_water_index)
21+
;; to differentiate water bodies from land.
22+
23+
24+
25+
;; ## Fetching OpenLayers packages
26+
;; The [OpenLayers](https://openlayers.org/) map widget will allow us to fetch
27+
;; multi-spectral satellite imagery in GeoTIFF format and manipulate the data
28+
;; with WebGL.
29+
30+
;; First fetch the packages from CDNs:
31+
(kind/hiccup
32+
[:div
33+
[:link {:rel "stylesheet"
34+
:href "https://cdn.jsdelivr.net/npm/ol@v10.6.0/ol.css"}]
35+
[:script {:src "https://cdn.jsdelivr.net/npm/geotiff"}]
36+
[:script {:src "https://cdn.jsdelivr.net/npm/ol@v10.6.0/dist/ol.js"}]])
37+
38+
;; ## Handling the GeoTIFF files
39+
;; We will create an OpenLayers WebGL layer that uses a GeoTIFF source.
40+
;; Here we are using a
41+
;; [Sentinel-2](https://dataspace.copernicus.eu/data-collections/copernicus-sentinel-data/sentinel-2)
42+
;; 10m resolution GeoTIFF image that
43+
;; contains blue, green, red and near infra-red bands.
44+
45+
(kind/scittle
46+
'(def geotiff-sources
47+
(js/ol.source.GeoTIFF.
48+
(clj->js
49+
{:sources [{:min 0
50+
:nodata 0
51+
:max 10000
52+
:bands [1 ;; B02 blue (490nm) -> band 1
53+
2 ;; B03 green (560nm) -> band 2
54+
3 ;; B04 red (665nm) -> band 3
55+
4 ;; B08 NIR (841nm) -> band 4
56+
]
57+
:url "https://s2downloads.eox.at/demo/Sentinel-2/3857/R10m.tif"}]}))))
58+
59+
;; ## OpenLayers WebGL DSL helpers
60+
;; Openlayers has it's own [lisp-oid DSL](https://openlayers.org/en/latest/apidoc/module-ol_expr_expression.html#~ExpressionValue)
61+
;; that gets converted to
62+
;; GLSL. Here we set up a few helper functions to allow us to
63+
;; write it more clearly.
64+
65+
(kind/scittle
66+
;; helper fns for writing in OpenLayers GL DSL
67+
'(do
68+
(defn gl-var [kw]
69+
["var" kw])
70+
(defn gl-band [band-key]
71+
["band" (band-key {:blue 1
72+
:green 2
73+
:red 3
74+
:NIR 4})])))
75+
76+
;; ## Normalized Difference Water Index
77+
^:kindly/hide-code
78+
(kind/tex "\\text{NDWI} = \\frac{\\text{Green} - \\text{NIR}}{\\text{Green} + \\text{NIR}}")
79+
;; NDWI is a band ratio Index that allows us to detect the presence of
80+
;; bodies of water.
81+
;; Water bodies absorb most NIR and reflect plenty of green light,
82+
;; so they result in a positive NDWI.
83+
;; Vegetation and soil strongly reflect NIR and and moderately reflect
84+
;; green light so they result in a negative NDWI.
85+
86+
87+
(kind/scittle
88+
'(def NDWI
89+
["/"
90+
["-" (gl-band :green) (gl-band :NIR)]
91+
["+" (gl-band :green) (gl-band :NIR)]]))
92+
93+
;; ## WebGL instructions
94+
;; Here we compose the WebGL instructions that will get converted
95+
;; into a GLSL fragment shader. We can access variables that will
96+
;; be passed in a separate variables map, and also the values from
97+
;; the bands selected above.
98+
99+
;; In this example we show the NDWI value in a gradient from blue
100+
;; to green. NDWI values below the `land_threshold` are displayed in RGB
101+
;; with a brightness variable so we can dim the land pixels dynamically.
102+
103+
(kind/scittle
104+
'(def color
105+
(let [opaque-blue [1 0 0 255]
106+
opaque-green [0 1 0 255]]
107+
["case"
108+
["<" NDWI (gl-var :land_threshold)]
109+
["color"
110+
["*" (gl-var :land_brightness) (gl-band :red)]
111+
["*" (gl-var :land_brightness) (gl-band :green)]
112+
["*" (gl-var :land_brightness) (gl-band :blue)]
113+
225]
114+
115+
["interpolate"
116+
["linear"]
117+
NDWI
118+
(gl-var :ndwi_upper) opaque-blue ; (lower bound)
119+
(gl-var :ndwi_lower) opaque-green ; (upper bound)
120+
]])))
121+
122+
123+
(kind/scittle
124+
;; default variable values for the map vis
125+
'(def param-defaults {:land_threshold 0.0
126+
:land_brightness 3.0
127+
:ndwi_upper 1.0
128+
:ndwi_lower 0.0}))
129+
130+
;; ## Openlayers map widget setup
131+
;; Below we create a map widget with a WebGL tile and a
132+
;; view based on the bounds of our GeoTIFF file.
133+
134+
(kind/scittle
135+
'(def webgl-tile-layer
136+
(js/ol.layer.WebGLTile.
137+
(clj->js {:style {:variables param-defaults
138+
:color color}
139+
:source geotiff-sources}))))
140+
141+
(kind/scittle
142+
'(let [osm-layer (->> {:source (js/ol.source.OSM.)}
143+
clj->js
144+
js/ol.layer.Tile.)
145+
146+
view (->> {:center [0 0]
147+
:zoom 2}
148+
clj->js
149+
js/ol.View.)]
150+
(->> {:target "map"
151+
:layers [webgl-tile-layer]
152+
:view (.getView geotiff-sources)}
153+
clj->js
154+
js/ol.Map.)))
155+
156+
;; ## Dynamic controls
157+
;; Using reagent we can create some sliders to control the
158+
;; WebGL variables, triggering a redraw on change.
159+
160+
(kind/reagent
161+
['(fn []
162+
(let [*map-params (reagent.core/atom param-defaults)
163+
slider
164+
(fn [label param minv maxv]
165+
(let [value (param @*map-params)]
166+
[:div
167+
[:label {:for param} label] " : "
168+
[:input
169+
{:id param
170+
:type "range"
171+
:min minv
172+
:max maxv
173+
:value value
174+
:step 0.01
175+
:on-change (fn [e]
176+
(let [newv (js/parseFloat (.. e -target -value))]
177+
(swap! *map-params #(assoc % param newv))
178+
(.updateStyleVariables
179+
webgl-tile-layer
180+
(clj->js @*map-params))))}]]))]
181+
182+
(fn []
183+
[:div {:style {:text-align "right"}}
184+
[slider "Land Threshold" :land_threshold 0.0 1.0]
185+
[slider "Land Brightness" :land_brightness 0.0 10.0]
186+
[slider "NDWI Upper" :ndwi_upper 0.0 1.0]
187+
[slider "NDWI Lower" :ndwi_lower 0.0 1.0]])))])
188+
189+
;; ## Results
190+
;; We can see the river clearly outlined by the NDWI index. The sliders
191+
;; allow us to dim the land to see the water outline more clearly, and
192+
;; adjust the thresholds to control how the water is displayed.
193+
194+
^:kindly/hide-code
195+
(kind/hiccup
196+
;; a div to render the map in
197+
[:div#map {:style {:width "500px"
198+
:height "500px"}}])

0 commit comments

Comments
 (0)