Comparing NDVI mean with sentinel-hub request builder and the value from statistical API

[aswin.manohar]

I used the following modified evalscript to download statistics for NDVI and NDMI for an AOI. here is the evalscript:

ndvi_new_evalscript = """

//VERSION=3

function setup() {
  return {
    input: [{
      bands: ["B04","B05","B08", "B8A","B11", "B12", "CLM", "CLP", "dataMask"],
      units: "DN"
    }],
    output: [
      {
        id: "bands",
        bands: ["B04"],
        sampleType: "UINT16"
      },
      {
        id: "masks",
        bands: ["CLM"],
        sampleType: "UINT16"
      },
      {
        id: "indices",
        bands: ["NDVI", "NDVI_RE1", "NDMI","CLP"],
        sampleType: "UINT16"
      },
      {
        id: "dataMask",
        bands: 1
      }]
  }
}

function evaluatePixel(samples) {
    // Normalised Difference Vegetation Index and variation
    let NDVI = (samples.B08 - samples.B04) / (samples.B08 + samples.B04);
    let NDVI_RE1 = index(samples.B08, samples.B05);
    let NDMI = (samples.B08 - samples.B11) / (samples.B08 + samples.B11);

    // Bare Soil Index
   // let NBSI = index((samples.B11 + samples.B04), (samples.B08 + samples.B02));

    // cloud probability normalized to interval [0, 1]
    let CLP = samples.CLP / 255.0;

    // masking cloudy pixels
    let combinedMask = samples.dataMask
    if (samples.CLM > 0) {
        combinedMask = 0;
    }

    // masking water pixels 
    let waterMask = samples.dataMask
    if (samples.SCL == 6) {
        waterMask = 0;
    }



    const f = 5000;
    return {
        bands: [samples.B04],
        masks: [samples.CLM],
        indices: [toUINT(NDVI, f), toUINT(NDVI_RE1, f), toUINT(CLP, f), toUINT(NDMI, f) ],
        dataMask: [combinedMask*waterMask]
    };
}

function toUINT(product, constant){
  // Clamp the output to [-1, 10] and convert it to a UNIT16
  // value that can be converted back to float later.
  if (product < -1) {
    product = -1;
  } else if (product > 10) {
    product = 10;
  }
  return Math.round(product * constant) + constant;
}
"""

The python call was like this:

aggregation = SentinelHubStatistical.aggregation(
    evalscript=ndvi_new_evalscript, time_interval=yearly_time_interval, aggregation_interval="P1D", resolution=(10, 10)
)

#calculations = {"default": {"statistics": {"default": {"percentiles": {"k": [5, 50, 95]}}}}}
histogram_calculations = {"ndvi": {"histograms": {"default": {"nBins": 20, "lowEdge": -1.0, "highEdge": 1.0}}}}

features_requests = []
for geo_shape in AOI_polygon.geometry.values:
    request = SentinelHubStatistical(
        aggregation=aggregation,
        input_data=[SentinelHubStatistical.input_data(DataCollection.SENTINEL2_L2A)],
        geometry=Geometry(geo_shape, crs=CRS(AOI_polygon.crs)),
        calculations=histogram_calculations,
        config=config
    )

    features_requests.append(request)

I got an negative NDVI value and I wanted to check if it is the same Sentinel-2 image that is returned with just getting the least cloud cover image from that particular date or interval (2018-6-29 to 2018-6-30). Looks like the image from 2018-6-29 to 2018-6-30 returns a median 0.3 NDMI value as opposed to -1 that was returned by the statistical API.

It can also be that I am not using the exact evalscript and the masks that were used in the statistical API. Do you know if I can do a one-to-one comparison as further checks.

[william.ray]

Hi Aswin,

To use the image acquisition with the least cloud cover in your time interval you will need to explicitly define this as an argument in your Statistical API request; this will be part of the input_data section of your request payload. By default the API will obtain the most recent image, so this could be the reason why you are returning different NDMI values.

input_data=[
        SentinelHubStatistical.input_data(
            DataCollection.SENTINEL2_L2A,            
            other_args={"dataFilter": {"mosaickingOrder": "leastCC"}},            
      ),

If I have misinterpreted your question, then let us know

[aswin.manohar]

Thanks for the answers. Is there a way that I can also return RGB images for which the mean NDVI was calculated from the statistical API

[william.ray]

Hi Aswin, you can do this using a separate Process API request. Statistical API doesn’t return any pixels, only the underlying statistics.

[aswin.manohar]

Yes, I am wondering how I can adapt the same script to produce images since it has data masks and cloud masks.

[william.ray]

As I said, you need to use Process API to produce images; your current script uses Statistical API so you cannot produce images.

If you are looking for a True Color evalscript that uses a cloud mask, then I can point you to this example.

//VERSION=3
function setup() {
  return {
    input: ["B02", "B03", "B04", "CLM"],
    output: { bands: 3 }
  }
}

function evaluatePixel(sample) {
  if (sample.CLM == 1) {
    return [0.75 + sample.B04, sample.B03, sample.B02]
  }
  return [3.5*sample.B04, 3.5*sample.B03, 3.5*sample.B02];
}