Derived Root-Zone Soil Water Content

//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: { bands: 4 },
        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];
}


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
const visualizer = new ColorRampVisualizer(cmap, vmin, vmax);


function evaluatePixel(samples, scenes) {
    // When there are no dates, return no data
    if (samples.length == 0) return [NaN, NaN, NaN, 0];

    // When there is no data for the last day, don't run calculation, return no data
    if (!samples[0].dataMask) return [NaN, NaN, NaN, 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 [...imgVals, samples[0].dataMask];
}
//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];
}

Evaluate and Visualize

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 observe 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 near-surface 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 Derived root-zone soil water content Sentinel-2 image

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 to 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]:

R(t)t=1T[S(t)R(t)].

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:

Rn=Rn1+Kn[SnRn1]

The gain K is defined as:

Kn=Kn1Kn1+eΔt/T

Δt is the time between the two consecutive observations. For the first time step, we set R1=S1 and Kn=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$$

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

  • 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

  • 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

  • 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