Add Derived Root-Zone Soil Water Content script

planetary-variables/soil-water-content/derived-root-zone-soil-water-content/eob.js

@@ -0,0 +1,140 @@
+//VERSION=3
+
+const nDays = 10;           // The number of days to load data for
+const scaleFactor = 1000;   // The scale factor for the SWC values
+const vmin = 0.1;           // The minimum value of the colormap
+const vmax = 0.4;           // The maximum value of the colormap
+
+function setup() {
+    return {
+        input: ["SWC", "dataMask"],
+        output: [
+            { id: "default", bands: 4 },
+            { id: "index", bands: 1, sampleType: "FLOAT32" },
+            { id: "eobrowserStats", bands: 1, sampleType: "FLOAT32" },
+            { id: "dataMask", bands: 1 },            
+        ],
+        mosaicking: "ORBIT"
+    };
+}
+
+
+function preProcessScenes(collections){
+    collections.scenes.orbits = collections.scenes.orbits.filter(function (orbit) {
+        var orbitDateFrom = new Date(orbit.dateFrom)
+        // Select all images within the last nDays
+        return orbitDateFrom.getTime() >= (collections.to.getTime() - (nDays * 24 * 3600 * 1000));
+    })
+    return collections
+}
+
+
+function expFilter(swc, dataMask, dates, timeConst) {
+    let deltaDays = 0;       // Initialize time difference with previous valid time [Unit days]
+    let gain = 1.0;          // Initialize gain
+    let referenceIndex = -1; // Initialize index of the first day with data
+    let referenceDate = NaN; // Initialize date of the first day with data
+    let referenceSWC = NaN;  // Initialize SWC value of the first day with data
+
+    // Find the first day with data in the time series
+    for (let i = 0; i < dates.length; i++) {
+        if (dataMask[i]) {
+            // Set index, date and SWC initial values to the first date with data
+            referenceIndex = i;
+            referenceDate = dates[i];
+            referenceSWC = swc[i];
+            break;
+        }
+    }
+
+    // Set first day of derived root-zone SWC equal to surface SWC
+    let drzswc = referenceSWC;
+
+    if (referenceIndex > -1) {
+        // Only apply the filter if there is at least one valid value in the time series
+        for (let i = referenceIndex; i < dates.length; i++) {
+            if (dataMask[i]) {
+                deltaDays = (dates[i] - referenceDate) / (24 * 3600 * 1000);
+                gain = gain / (gain + Math.exp(-deltaDays / timeConst));
+                drzswc = drzswc + gain * (swc[i] - drzswc);
+                referenceDate = dates[i]
+            }
+        }
+    }
+    return [drzswc];
+}
+
+
+function updateColormap(vmin, vmax) {
+    const numIntervals = cmap.length;
+    const intervalLength = (vmax - vmin) / (numIntervals - 1);
+    for (let i = 0; i < numIntervals; i++) {
+        cmap[i][0] = vmin + intervalLength * i;
+    }
+}
+
+
+const cmap = [
+    [0.0, 0xfff7ea],
+    [0.05, 0xfaedda],
+    [0.1, 0xede4cb],
+    [0.15, 0xdedcbd],
+    [0.2, 0xced3af],
+    [0.25, 0xbdcba3],
+    [0.3, 0xaac398],
+    [0.35, 0x96bc90],
+    [0.4, 0x80b48a],
+    [0.45, 0x68ac86],
+    [0.5, 0x4da484],
+    [0.55, 0x269c83],
+    [0.6, 0x009383],
+    [0.65, 0x008a85],
+    [0.7, 0x008186],
+    [0.75, 0x007788],
+    [0.8, 0x006d8a],
+    [0.85, 0x00618c],
+    [0.9, 0x00558d],
+    [0.95, 0x00478f],
+    [1.0, 0x003492],
+];
+
+// Prepare colormap based on provided min and max values
+updateColormap(vmin, vmax);
+const visualizer = new ColorRampVisualizer(cmap);
+
+
+function evaluatePixel(samples, scenes) {
+    // When there are no dates, return no data
+    if (samples.length == 0) return {
+        default: [NaN, NaN, NaN, 0],
+        index: [NaN],
+        eobrowserStats: [NaN],
+        dataMask: [0],
+    }  ;
+
+    // When there is no data for the last day, don't run calculation, return no data
+    if (!samples[0].dataMask) return {
+        default: [NaN, NaN, NaN, 0],
+        index: [NaN],
+        eobrowserStats: [NaN],
+        dataMask: [0],
+    }  ;
+
+    // Extract dates, data mask and scaled SWC values and sort them ascending (oldes to newest)
+    var datesASC = scenes.map(scene => new Date(scene.date)).reverse();
+    var swcASC = samples.map(sample => sample.SWC / scaleFactor).reverse();
+    var dataMaskASC = samples.map(sample => sample.dataMask).reverse();
+
+    // Calculate derived root-zone SWC by applying exponential filter
+    const drzswc = expFilter(swcASC, dataMaskASC, datesASC, nDays);
+
+    // Apply colormap
+    let imgVals = visualizer.process(drzswc);
+
+    return {
+        default: [...imgVals, samples[0].dataMask],
+        index: [drzswc],
+        eobrowserStats: [drzswc],
+        dataMask: [samples[0].dataMask],
+    }        
+}

planetary-variables/soil-water-content/derived-root-zone-soil-water-content/index.md

@@ -0,0 +1,79 @@
+---
+title: Derived Root-Zone Soil Water Content
+grand_parent: Planetary Variables
+parent: Soil Water Content
+layout: script
+nav_exclude: false
+scripts:
+    - [Visualization, script.js]
+    - [EO Browser, eob.js]
+    - [Raw Values, raw.js]
+examples:
+- zoom: '11'
+  lat: '41.1921'
+  lng: '-93.845'
+  datasetId: '65f7e4fb-a27a-4fae-8d79-06a59d7e6ede'
+  fromTime: '2022-05-01T00:00:00.000Z'
+  toTime: '2022-05-26T23:59:59.999Z'
+  platform:
+  - EOB
+  evalscripturl: https://custom-scripts.sentinel-hub.com/custom-scripts/planetary-variables/soil-water-content/derived-root-zone-soil-water-content/eob.js
+  additionalQueryParams: 
+  - - themeId
+    - PLANET_SANDBOX
+---
+## General description
+Here, we show how to compute and visualize derived root-zone soil water content (DRZSWC) using the Sentinel Hub EO Browser. DRZSWC is an estimate of the amount of water in the soil in the root zone: the depth range over which plant roots take up most of their water. The root zone depends on the type of vegetation, but typically covers the upper 100 cm of the soil [Stocker et al., 2023]. With satellites, we can observed soil water content (SWC) in the upper layer of the soil, typically covering the first 5 to 10 cm. To estimate SWC in the root zone from the upper-layed SWC observed by satellites, we can use an exponential filter.
+
+## Description of representative images
+
+| Near-surface soil water content (May 26, 2022) | Derived root-zone soil water content (May 26, 2022) | Sentinel-2 image (June 20, 2022) |
+|:----:|:----:|:----:| 
+| ![Near-surface soil water content](fig/sh_swc_22_05_26.png)  | ![Derived root-zone soil water content](fig/sh_drzswc_22_05_26.png) | ![Sentinel-2 image](fig/sh_opt_22_06_20.jpeg) |
+
+In the figure above, we show near-surface and root-zone soil water content in Iowa (USA), on a rainy day after a dry period. The surface has become wet (blue), while the deeper root zone is still dry (yellow) due to the long dry spell. 
+
+## Algorithm
+The main drivers of SWC changes (precipitation and evaporation) happen at the surface. As a result, near-surface soil water content tends to react more and faster on precipitation and evaporation than soil water content in the deeper root zone. Omitting lateral transport and drainage to deeper layers, we can describe the changes over time $$t$$ of DRZSWC $$R$$ as a simple diffusion process that is a function of near-surface soil water content $$S$$ [Wagner et al., 1998]:
+
+$$
+\frac{\partial R(t)}{\partial t} = \frac{1}{T} \left[S(t) - R(t)\right].
+$$
+
+Here, $$T$$ is a time constant that sets the typical time scale that determines how quickly DRZSWC reacts to changes in near-surface SWC. 
+
+Because we do not have continuous measurements of near-surface soil water content, we need to discretize this differential equation. Albergel et al. [2008] and Paulik et al. [2014] have proposed a recursive equation to solve this system, that computes DRZSWC for each time step $$n$$ for which we have an observation of $$S$$:
+
+$$
+R_n = R_{n-1} + K_n \cdot \left[S_n - R_{n-1} \right] 
+$$
+
+The gain $$K$$ is defined as:
+
+$$
+K_n = \frac{K_{n-1}}{K_{n-1} + e^{-{\Delta t / T}}}
+$$
+
+$$\Delta t$$ is the time between the two consecutive observations. For the first time step, we set $$R_1 = S_1$$ and $$K_n = 1$$. 
+
+The time constant $$T$$ sets the response time of DRZSWC to near-surface SWC changes, and depends to the soil type and the depth of the root zone. Here, $$T$$ has been set to 10 days, which follows Albergel et al. [2008] for a root-zone soil layer of about 1 meter.  
+
+The effect of setting $$T$$ can be seen in the figure below, which shows a typical near-surface SWC time series from the Netherlands, measured by SMAP:
+
+![Root-zone SWC for various values of $$T$$](fig/ts_varying_T.png "Root-zone SWC for various values of T")
+
+The longer the time constant $$T$$, the more the high-frequency variations in SWC are damped. Also note here that the procedure needs some spin-up days. 
+
+## References
+- Albergel, C., Rüdiger, C., Pellarin, T., Calvet, J.-C., Fritz, N., Froissard, F., Suquia, D., Petitpa, A., Piguet, B., & Martin, E. (2008). From near-surface to root-zone soil moisture using an exponential filter: An assessment of the method based on in-situ observations and model simulations. Hydrology and Earth System Sciences, 12(6), 1323–1337. [https://doi.org/10.5194/hess-12-1323-2008](https://doi.org/10.5194/hess-12-1323-2008)
+
+- Paulik, C., Dorigo, W., Wagner, W., & Kidd, R. (2014). Validation of the ASCAT Soil Water Index using in situ data from the International Soil Moisture Network. International Journal of Applied Earth Observation and Geoinformation, 30, 1–8. [https://doi.org/10.1016/j.jag.2014.01.007](https://doi.org/10.1016/j.jag.2014.01.007)
+
+- Stocker, B. D., Tumber-Dávila, S. J., Konings, A. G., Anderson, M. C., Hain, C., & Jackson, R. B. (2023). Global patterns of water storage in the rooting zones of vegetation. Nature Geoscience. [https://doi.org/10.1038/s41561-023-01125-2](https://doi.org/10.1038/s41561-023-01125-2)
+
+- Wagner, W., Lemoine, G., & Rott, H. (1999). A Method for Estimating Soil Moisture from ERS Scatterometer and Soil Data. Remote Sensing of Environment, 70(2), 191–207. [https://doi.org/10.1016/S0034-4257(99)00036-X](https://doi.org/10.1016/S0034-4257(99)00036-X)
+
+## Useful links
+-   [SWC Technical specifications](https://developers.planet.com/docs/planetary-variables/soil-water-content-technical-specification/)
+-   [SWC Data sheet](https://planet.widen.net/s/cv7bfjhhd5)
+-   [Sentinel Hub documentation about Soil Water Content](https://docs.sentinel-hub.com/api/latest/data/planetary-variables/soil-water-content/)

planetary-variables/soil-water-content/derived-root-zone-soil-water-content/raw.js

@@ -0,0 +1,76 @@
+//VERSION=3
+
+const nDays = 10;           // The number of days to load data for
+const scaleFactor = 1000;   // The scale factor for the SWC values
+
+function setup() {
+    return {
+        input: ["SWC", "dataMask"],
+        output: { bands: 1, sampleType: "FLOAT32" },
+        mosaicking: "ORBIT"
+    };
+}
+
+
+function preProcessScenes(collections){
+    collections.scenes.orbits = collections.scenes.orbits.filter(function (orbit) {
+        var orbitDateFrom = new Date(orbit.dateFrom)
+        // Select all images within the last nDays
+        return orbitDateFrom.getTime() >= (collections.to.getTime() - (nDays * 24 * 3600 * 1000));
+    })
+    return collections
+}
+
+
+function expFilter(swc, dataMask, dates, timeConst) {
+    let deltaDays = 0;       // Initialize time difference with previous valid time [Unit days]
+    let gain = 1.0;          // Initialize gain
+    let referenceIndex = -1; // Initialize index of the first day with data
+    let referenceDate = NaN; // Initialize date of the first day with data
+    let referenceSWC = NaN;  // Initialize SWC value of the first day with data
+
+    // Find the first day with data in the time series
+    for (let i = 0; i < dates.length; i++) {
+        if (dataMask[i]) {
+            // Set index, date and SWC initial values to the first date with data
+            referenceIndex = i;
+            referenceDate = dates[i];
+            referenceSWC = swc[i];
+            break;
+        }
+    }
+
+    // Set first day of derived root-zone SWC equal to surface SWC
+    let drzswc = referenceSWC;
+
+    if (referenceIndex > -1) {
+        // Only apply the filter if there is at least one valid value in the time series
+        for (let i = referenceIndex; i < dates.length; i++) {
+            if (dataMask[i]) {
+                deltaDays = (dates[i] - referenceDate) / (24 * 3600 * 1000);
+                gain = gain / (gain + Math.exp(-deltaDays / timeConst));
+                drzswc = drzswc + gain * (swc[i] - drzswc);
+                referenceDate = dates[i]
+            }
+        }
+    }
+    return [drzswc];
+}
+
+function evaluatePixel(samples, scenes) {
+    // When there are no dates, return no data
+    if (samples.length == 0) return [NaN];
+
+    // When there is no data for the last day, don't run calculation, return no data
+    if (!samples[0].dataMask) return [NaN];
+
+    // Extract dates, data mask and scaled SWC values and sort them ascending (oldes to newest)
+    var datesASC = scenes.map(scene => new Date(scene.date)).reverse();
+    var swcASC = samples.map(sample => sample.SWC / scaleFactor).reverse();
+    var dataMaskASC = samples.map(sample => sample.dataMask).reverse();
+
+    // Calculate derived root-zone SWC by applying exponential filter
+    const drzswc = expFilter(swcASC, dataMaskASC, datesASC, nDays);
+
+    return [drzswc];
+}