Merge pull request #355 from azlinszkysinergise/planet_scripts

new kndvi script for planetscope

Author: blaz-pavlica

planet/planetscope/agriculture_growth_stage/README.md

@@ -0,0 +1,61 @@
+---
+title: Agricultural growth stage
+parent: PlanetScope
+grand_parent: Planet
+layout: script
+permalink: /planet/planetscope/agriculture_growth_stage
+nav_exclude: true
+scripts:
+  - [PlanetScope data, script.js]
+  - [Analysis-Ready PlanetScope data, analysis_ready_planetscope.js]
+
+examples:
+- zoom: '14'
+  lat: '47.71371'
+  lng: '23.09911'
+  datasetId: ccb1f8f0-e5bf-4c31-afe5-d8803bcbde2a
+  fromTime: '2023-05-01T00:00:00.000Z'
+  toTime: '2023-08-30T23:59:59.999Z'
+  platform:
+  - EOB
+  evalscripturl: https://custom-scripts.sentinel-hub.com/sentinel-2/agriculture_growth_stage/script.js
+  additionalQueryParams:
+  - - mosaickingOrder
+    - mostRecent
+---
+
+## Author of the script
+[@HarelDan](https://github.com/hareldunn/GIS_Repo/blob/master/Multi-Temporal%20NDVI%20for%20Sentinel%20Hub%20Custom%20Scripts){:target="_blank"} 
+
+Adapted to [PlanetScope](https://docs.planet.com/data/imagery/planetscope/) and [Analysis-Ready PlanetScope](https://docs.planet.com/data/imagery/arps/) by András Zlinszky and Github Copilot.
+
+## General description of the script
+Agricultural growth stage is a script visualizing the multi-temporal NDVI trends in Sentinel-2 imagery. It takes the current image as baseline and calculates the average NDVI for the previous 2 months.
+The script requires multi-temporal processing, so the parameter TEMPORAL=true should be added to the request.
+A simple stretching is applied to NDVI between 0.1 and 0.7 by default, then the mean NDVI from the first, second and third month is assigned to the Red, Green and Blue channels of the image respectively, creating a composite image. What you see on the composite is 
+- how dense and/or vigorous the vegetation is, represented by the brightness of the color from black (no vegetation at all) to white (dense green vegetation all year), with various shades of color in between
+- when the vegetation peak happens and how distinct it is. Vegetation with a single very distict peak will be one of the primary colors (Red, Green, Blue) while vegetation with a longer, more even growth season will be yellow (between Red and Green) or cyan (between Green and Blue). Purple color may indicate two vegetation peaks, one in the first month and another in the last, with a dry period or grassland mowing in between.
+
+## How to use
+
+- In Planet Insights Platform Browser, open the calendar panel dropdown (with the dropdown button on the right)
+- Select the time interval view (the calendar icon with arrows on the top right). You will now see two dates, labeled "from" and "until".
+- Select these dates to cover an interval of three months for the regular script.
+- Select your evalscript of choice from this script website and copy it from the code window above
+- In Planet Insights Browser, click the `Custom Script` in the list of layers and the `Custom Script` tab to open the evalscript code window
+- Select the full text inside the window (eg. with the Ctrl+A hotkey) and paste the evalscript code from the clipboard
+- Wait until the data loads - this may take some time for large areas.
+
+## Description of representative images
+
+The Agricultural growth stage script applied to the agricultural fields of the Nile Delta, just north of Cairo, Egypt. 
+
+![The Agricultural growth stage script applied to agricultural fields of the Nile Delta just north of Cairo.](fig/nile_delta_arps.jpg)
+
+## References
+Based on: 
+[source 1](https://twitter.com/sentinel_hub/status/922813457145221121){:target="_blank"}, 
+[source 2](https://twitter.com/sentinel_hub/status/1020755996359225344){:target="_blank"}
+
+
+

planet/planetscope/agriculture_growth_stage/analysis_ready_planetscope.js

@@ -0,0 +1,80 @@
+//VERSION=3 (auto-converted from 1)
+/*
+Source: @HarelDan - https://github.com/hareldunn/GIS_Repo/blob/master/Multi-Temporal%20NDVI%20for%20Sentinel%20Hub%20Custom%20Scripts
+Visualizing NDVI multi-temporal trends in Sentinel-2 imagery.
+will take the current image as baseline and calculate average NDVI for the previous 2 months
+Based on:
+https://twitter.com/sentinel_hub/status/922813457145221121
+https://twitter.com/sentinel_hub/status/1020755996359225344
+Script requires multi-temporal processing so parameter TEMPORAL=true should be added to the request Set the time interval to at least three months.
+Adapted to PlanetScope and Analysis-ready PlanetScope by @azlinszky.bsky.social
+This version is for Analysis-ready PlanetScope data, where cloudy pixels are already filtered out.
+*/
+
+function setup() {
+    return {
+        input: [
+            {
+                bands: ["red", "nir"],
+            },
+        ],
+        output: { bands: 3 },
+        mosaicking: "ORBIT",
+    };
+}
+
+function calcNDVI(sample) {
+    var denom = sample.red + sample.nir;
+    return denom != 0 ? (sample.nir - sample.red) / denom : 0.0;
+}
+function stretch(val, min, max) {
+    return (val - min) / (max - min);
+}
+
+function evaluatePixel(samples, scenes) {
+     if (!scenes || scenes.length === 0) {
+          return [0, 0.3, 0.3]; // just a color to indicate something is wrong
+    }
+    var avg1 = 0;
+    var count1 = 0;
+    var avg2 = 0;
+    var count2 = 0;
+    var avg3 = 0;
+    var count3 = 0;
+    var endMonth = scenes[0].date.getMonth();
+
+    for (var i = 0; i < samples.length; i++) {
+        var ndvi = calcNDVI(samples[i]);
+        if (scenes[i].date.getMonth() == endMonth) {
+            avg3 = avg3 + ndvi;
+            count3++;
+        } else if (scenes[i].date.getMonth() == endMonth - 1) {
+            avg2 = avg2 + ndvi;
+            count2++;
+        } else {
+            avg1 = avg1 + ndvi;
+            count1++;
+        }
+    }
+    avg1 = avg1 / count1;
+    avg2 = avg2 / count2;
+    avg3 = avg3 / count3;
+    avg1 = stretch(avg1, 0.1, 0.7);
+    avg2 = stretch(avg2, 0.1, 0.7);
+    avg3 = stretch(avg3, 0.1, 0.7);
+
+    return [avg1, avg2, avg3];
+}
+
+function preProcessScenes(collections) {
+    collections.scenes.orbits = collections.scenes.orbits.filter(function (
+        orbit
+    ) {
+        var orbitDateFrom = new Date(orbit.dateFrom);
+        return (
+            orbitDateFrom.getTime() >=
+            collections.to.getTime() - 3 * 31 * 24 * 3600 * 1000
+        );
+    });
+    return collections;
+}

planet/planetscope/agriculture_growth_stage/fig/nile_delta_arps.jpg

@@ -0,0 +1,82 @@
+//VERSION=3 (auto-converted from 1)
+/*
+Source: @HarelDan - https://github.com/hareldunn/GIS_Repo/blob/master/Multi-Temporal%20NDVI%20for%20Sentinel%20Hub%20Custom%20Scripts
+Visualizing NDVI multi-temporal trends in Sentinel-2 imagery.
+will take the current image as baseline and calculate average NDVI for the previous 2 months
+Based on:
+https://twitter.com/sentinel_hub/status/922813457145221121
+https://twitter.com/sentinel_hub/status/1020755996359225344
+Script requires multi-temporal processing so parameter TEMPORAL=true should be added to the request
+Adapted to PlanetScope and Analysis-ready PlanetScope by @azlinszky.bsky.social
+This version is for PlanetScope data, checking cloudiness with the "clear" band.
+*/
+
+function setup() {
+    return {
+        input: [
+            {
+                bands: ["red", "nir", "clear"],
+            },
+        ],
+        output: { bands: 3 },
+        mosaicking: "ORBIT",
+    };
+}
+
+function calcNDVI(sample) {
+    var denom = sample.red + sample.nir;
+    return denom != 0 ? (sample.nir - sample.red) / denom : 0.0;
+}
+function stretch(val, min, max) {
+    return (val - min) / (max - min);
+}
+
+function evaluatePixel(samples, scenes) {
+    if (!scenes || scenes.length === 0) {
+        return [0, 0.3, 0.3]; // just a color to indicate something is wrong
+    }
+    var avg1 = 0;
+    var count1 = 0;
+    var avg2 = 0;
+    var count2 = 0;
+    var avg3 = 0;
+    var count3 = 0;
+    var endMonth = scenes[0].date.getMonth();
+
+    for (var i = 0; i < samples.length; i++) {
+        // Only use clear pixels
+        if (!samples[i].clear) continue;
+        var ndvi = calcNDVI(samples[i]);
+        if (scenes[i].date.getMonth() == endMonth) {
+            avg3 = avg3 + ndvi;
+            count3++;
+        } else if (scenes[i].date.getMonth() == endMonth - 1) {
+            avg2 = avg2 + ndvi;
+            count2++;
+        } else {
+            avg1 = avg1 + ndvi;
+            count1++;
+        }
+    }
+    avg1 = count1 > 0 ? avg1 / count1 : 0;
+    avg2 = count2 > 0 ? avg2 / count2 : 0;
+    avg3 = count3 > 0 ? avg3 / count3 : 0;
+    avg1 = stretch(avg1, 0.1, 0.7);
+    avg2 = stretch(avg2, 0.1, 0.7);
+    avg3 = stretch(avg3, 0.1, 0.7);
+
+    return [avg1, avg2, avg3];
+}
+
+function preProcessScenes(collections) {
+    collections.scenes.orbits = collections.scenes.orbits.filter(function (
+        orbit
+    ) {
+        var orbitDateFrom = new Date(orbit.dateFrom);
+        return (
+            orbitDateFrom.getTime() >=
+            collections.to.getTime() - 3 * 31 * 24 * 3600 * 1000
+        );
+    });
+    return collections;
+}

planet/planetscope/agriculture_growth_stage/script.js

@@ -1,11 +1,21 @@
 //VERSION=3
 //Cloudless Mosaic with PlanetScope
+//adapted by András Zlinszky and Google Gemini to handle first quartile or median mosaicking
 
 function setup() {
     return {
         input: ["red", "green", "blue", "dataMask", "clear"],
         output: { bands: 4 },
-        mosaicking: "ORBIT",
+        mosaicking: "ORBIT", // 'ORBIT' or 'TILE' or 'NONE' - determines initial mosaicking behavior
+        // Define a custom parameter for the mosaicking method
+        processing: {
+            // 'mosaickingMethod' is the name of your custom parameter.
+            // Users can set this to 'median' or 'q1'.
+            mosaickingMethod: {
+                defaultValue: "q1", // Default to first quartile
+                validValues: ["median", "q1"],
+            },
+        },
     };
 }
 
@@ -26,52 +36,106 @@ function preProcessScenes(collections) {
 function getLastObservation(arr) {
     for (let i = arr.length - 1; i >= 0; i--) {
         if (arr[i] !== 0) {
-            // optional check if you are sure all invalid observations are filtered out
             return arr[i];
         }
     }
     return 0;
 }
 
 function getMedian(sortedValues) {
-    var index = Math.floor(sortedValues.length / 2);
-    return sortedValues[index];
+    const n = sortedValues.length;
+    if (n === 0) {
+        return undefined;
+    }
+
+    const mid = Math.floor(n / 2);
+
+    if (n % 2 === 1) {
+        return sortedValues[mid];
+    } else {
+        return (sortedValues[mid - 1] + sortedValues[mid]) / 2;
+    }
+}
+
+/**
+ * Calculates the first quartile (Q1) of a sorted array of numbers.
+ *
+ * @param {number[]} sortedValues An array of numbers sorted in ascending order.
+ * @returns {number} The first quartile (Q1) of the distribution, or undefined if array is empty.
+ */
+function getFirstQuartile(sortedValues) {
+    const n = sortedValues.length;
+    if (n === 0) {
+        return undefined;
+    }
+
+    const lowerHalfEndIndex = Math.floor(n / 2);
+    const lowerHalf = sortedValues.slice(0, lowerHalfEndIndex);
+
+    return getMedian(lowerHalf);
 }
 
-function evaluatePixel(samples, scenes) {
+// *** CRITICAL CHANGE HERE: ADD 'properties' AS THE THIRD ARGUMENT ***
+function evaluatePixel(samples, scenes, properties) {
     var reds = [];
     var greens = [];
-    var blues = []; //empty arrays for reds greens and blues
-    var a = 0; //incrementer
+    var blues = [];
 
+    // Collect clear samples
     for (var i = 0; i < samples.length; i++) {
-        //for each sample
-        var sample = samples[i]; //get current sample
-        var clear = sample.dataMask && sample.clear; //0 for clouds OR datamask, 1 for neither
+        var sample = samples[i];
+        var clear = sample.dataMask && sample.clear;
 
         if (clear === 1) {
-            //if not clouds nor datamask
-            reds[a] = sample.red; //assign values for that sample to the channel arrays
-            blues[a] = sample.blue;
-            greens[a] = sample.green;
-            a = a + 1; //increment a to represent that at this specific pixel, a value was detected
+            reds.push(sample.red);
+            blues.push(sample.blue);
+            greens.push(sample.green);
         }
     }
 
     var rValue;
     var gValue;
     var bValue;
+    var transparency;
+
+    if (reds.length > 0) {
+        // IMPORTANT: Sort the arrays by value before calculating statistics.
+        reds.sort((a, b) => a - b);
+        greens.sort((a, b) => a - b);
+        blues.sort((a, b) => a - b);
+
+        // *** Access the method from the 'properties' object ***
+        const method = properties.mosaickingMethod; // No 'processing' here, directly under properties
+
+        if (method === "median") {
+            rValue = getMedian(reds);
+            gValue = getMedian(greens);
+            bValue = getMedian(blues);
+        } else if (method === "q1") {
+            rValue = getFirstQuartile(reds);
+            gValue = getFirstQuartile(greens);
+            bValue = getFirstQuartile(blues);
+        } else {
+            // Fallback in case of an unexpected method value (shouldn't happen with validValues)
+            rValue = reds[0];
+            gValue = greens[0];
+            bValue = blues[0];
+        }
 
-    if (a > 0) {
-        rValue = getMedian(reds); // or call getLastObservation - which is less guaranteed to remove hazy images
-        gValue = getMedian(greens);
-        bValue = getMedian(blues);
         transparency = 1;
     } else {
-        rValue = 1;
-        gValue = 1;
-        bValue = 1;
+        // If no clear samples, default to black and fully transparent.
+        rValue = 0;
+        gValue = 0;
+        bValue = 0;
         transparency = 0;
     }
-    return [rValue / 3000, gValue / 3000, bValue / 3000, transparency];
-}
+
+    // Scale values for display (e.g., to 0-1 range).
+    return [
+        Math.min(1, Math.max(0, rValue / 3000)),
+        Math.min(1, Math.max(0, gValue / 3000)),
+        Math.min(1, Math.max(0, bValue / 3000)),
+        transparency
+    ];
+}

planet/planetscope/cloudless_mosaic/script.js

@@ -47,11 +47,14 @@ The spectral bands of PlanetScope data are the following if you order a 8-band p
  - [True Color]({% link planet/planetscope/true_color/index.md %})
  - [False Color]({% link planet/planetscope/false_color/index.md %})
  - [Cloudless Mosaic]({% link planet/planetscope/cloudless_mosaic/index.md %})
+ - [Agriculture Growth Stage]({% link planet/planetscope/agriculture_growth_stage/README.md %})
+ - [Green City Index]({% link planet/planetscope/green_city/index.md %})
+ - [Kernel NDVI]({% link planet/planetscope/kndvi/README.md %})
  - [NDVI]({% link planet/planetscope/ndvi/index.md %})
  - [Maximum NDVI]({% link planet/planetscope/max_ndvi/index.md %})
  - [NDVI Difference]({% link planet/planetscope/ndvi_difference/index.md %})
  - [NDWI]({% link planet/planetscope/ndwi/index.md %})
  - [NDCI - Normalized Difference Chlorophyll Index]({% link planet/planetscope/ndci/index.md %})
  - [NDRE - Normalized Difference Red Edge Index]({% link planet/planetscope/ndre/index.md %})
- - [Green City]({% link planet/planetscope/green_city/index.md %})
+ - [Penguin Locator]({% link planet/planetscope/penguin_locator/index.md %})
  - [UDM2 Cloud/Snow Classification]({% link planet/planetscope/cloud_classification/index.md %})