Sentinel-3 Synergy Vegetation Product.#
Products used: s3_syn_2_vg1
Keywords data used; Sentinel-3 Level-2 VG1 SYN
Background#
The Sentinel-3 SYN Level-2 VG1 product is a 1 km VEGETATION-like daily synthesis product generated from the Sentinel-3 Synergy processing chain. It combines observations from the Ocean and Land Colour Instrument (OLCI) and the Sea and Land Surface Temperature Radiometer (SLSTR) onboard Sentinel-3A and Sentinel-3B to provide harmonised land surface observations for vegetation monitoring. The product provides Top-of-Canopy (TOC) surface reflectance and NDVI on a regular latitude–longitude grid with an approximate 1 km spatial resolution at the equator.
The SY_2_VG1 product is designed as a SPOT/VEGETATION continuity product supporting long-term monitoring of vegetation dynamics across large areas. It provides daily synthesis observations using the best available pixel selection approach, based on the maximum NDVI compositing method, to improve consistency in vegetation monitoring applications. The product includes four main surface reflectance channels: blue (B0), red (B2), near-infrared (B3), and mid-infrared/SWIR (MIR). NDVI is also provided and is derived from the B3 near-infrared and B2 red channels.
The product delivers consistent measurements of vegetation condition and land surface reflectance at 1 km resolution, with a daily synthesis frequency, making it valuable for broad-scale environmental and agricultural applications, including:
Vegetation health and greenness monitoring
Drought early warning and food security assessments
Crop and rangeland condition monitoring
Land productivity and land degradation assessment
Phenology and seasonal vegetation dynamics monitoring
Long-term vegetation trend analysis
Climate, ecosystem and land surface process studies
Regional to continental-scale agricultural monitoring
In addition to reflectance and NDVI, SY_2_VG1 includes supporting atmospheric and geometric information such as aerosol optical thickness, total column ozone, total column water vapour and sun-viewing geometry. These variables are useful for interpreting atmospheric effects, filtering poor-quality observations and improving the reliability of vegetation time-series analysis.
More information on the product can be found in the Sentinel-3 Synergy Products documentation and the Sentinel-3 SYN Land User Handbook.
Description#
This notebook will cover following topics:
Inspecting the products and measurements available in the datacube
Loading Sentinel-3 SYN Level-2 VG1 Data.
Plotting the results
Getting started#
To run this analysis, run all the cells in the notebook, starting with the “Load packages” cell.
Load packages#
Import Python packages that are used for the analysis.
[1]:
%matplotlib inline
import datacube
import geopandas as gpd
import numpy as np
import xarray as xr
import pandas as pd
import matplotlib.pyplot as plt
from odc.geo.geom import Geometry
from deafrica_tools.plotting import display_map
from deafrica_tools.areaofinterest import define_area
Connect to the datacube#
Connect to the datacube so we can access DE Africa data.
[2]:
dc = datacube.Datacube(app="Sentinel-3 SYN SY_2_VG1")
List products#
We can use datacube’s list_products functionality to inspect DE Africa’s products that are available in the datacube. The table below shows the product names that we will use to load the data, a brief description of the data, and the satellite instrument that acquired the data.
[3]:
dc.list_products().loc[dc.list_products()['description'].str.contains('Sentinel-3')]
[3]:
| name | description | license | default_crs | default_resolution | |
|---|---|---|---|---|---|
| name | |||||
| s3_ol_2_wfr_nrt | s3_ol_2_wfr_nrt | Sentinel-3 Level 2 Water Full Resolution (WFR)... | CC-BY-4.0 | EPSG:4326 | Resolution(x=0.003, y=-0.003) |
| s3_olci_l2_lfr | s3_olci_l2_lfr | Sentinel-3 OLCI L2 LAND | CC-BY-4.0 | EPSG:4326 | Resolution(x=0.003, y=-0.003) |
| s3_olci_l2_wfr | s3_olci_l2_wfr | Sentinel-3 OLCI L2 WATER | CC-BY-4.0 | EPSG:4326 | Resolution(x=0.003, y=-0.003) |
| s3_slstr_l2_lst | s3_slstr_l2_lst | Sentinel-3 LST | CC-BY-4.0 | EPSG:4326 | Resolution(x=0.003, y=-0.003) |
| s3_syn_2_vg1 | s3_syn_2_vg1 | Sentinel-3 Level-2 VG1 SYN product | CC-BY-4.0 | EPSG:4326 | Resolution(x=0.01, y=-0.01) |
Please choose the product by specifying its name, with the product s3_syn_2_vg1
[4]:
product = "s3_syn_2_vg1"
List measurements#
We can further inspect the data available for the Sentinel-3 Land s3_olci_l2_lfr product using datacube’s list_measurements functionality. The table below lists each of the measurements available in the data.
[5]:
measurements = dc.list_measurements()
measurements.loc[product]
[5]:
| name | dtype | units | nodata | aliases | flags_definition | add_offset | scale_factor | |
|---|---|---|---|---|---|---|---|---|
| measurement | ||||||||
| AG_VG1 | AG_VG1 | uint8 | 1 | 255.0 | [aerosol_optical_thickness] | NaN | 0.0 | 1.0000 |
| B0_VG1 | B0_VG1 | int16 | 1 | 32767.0 | NaN | NaN | 0.0 | 0.0001 |
| B2_VG1 | B2_VG1 | int16 | 1 | 32767.0 | NaN | NaN | 0.0 | 0.0001 |
| B3_VG1 | B3_VG1 | int16 | 1 | 32767.0 | NaN | NaN | 0.0 | 0.0001 |
| MIR_VG1 | MIR_VG1 | int16 | 1 | 32767.0 | [middle_infrared] | NaN | 0.0 | 0.0001 |
| NDVI_VG1 | NDVI_VG1 | int16 | 1 | 32767.0 | NaN | NaN | 0.0 | 0.0001 |
| OG_VG1 | OG_VG1 | uint8 | atm.cm | 255.0 | [total_column_ozone] | NaN | 0.0 | 1.0000 |
| SAA_VG1 | SAA_VG1 | int16 | degrees | 32767.0 | [solar_azimuth_angle] | NaN | 0.0 | 1.0000 |
| SZA_VG1 | SZA_VG1 | uint8 | degrees | 255.0 | [solar_zenith_angle] | {'mask': {'bits': [0, 1, 2, 3, 4, 8, 16, 32, 6... | 0.0 | 1.0000 |
| TOA_NDVI_VG1 | TOA_NDVI_VG1 | int16 | 1 | 32767.0 | NaN | NaN | 0.0 | 0.0001 |
| VAA_VG1 | VAA_VG1 | int16 | degrees | 32767.0 | [view_azimuth_angle] | NaN | 0.0 | 1.0000 |
| VZA_VG1 | VZA_VG1 | uint8 | degrees | 255.0 | [view_zenith_angle] | NaN | 0.0 | 1.0000 |
| WVG_VG1 | WVG_VG1 | uint8 | g/cm^2 | 255.0 | [total_column_water_vapour] | NaN | 0.0 | 1.0000 |
The Sentinel-3 SYN SY_2_VG1 product has thirteen measurements:
Measurement |
Description |
Application |
|---|---|---|
B0_VG1 |
Surface reflectance in the blue band at about 450 nm. |
Atmospheric correction support, surface reflectance analysis, vegetation/background separation. |
B2_VG1 |
Surface reflectance in the red band at about 645 nm. |
Used with NIR to assess vegetation greenness and photosynthetic activity. |
B3_VG1 |
Surface reflectance in the near-infrared band at about 835 nm. |
Vegetation condition monitoring, biomass sensitivity, NDVI calculation and land productivity analysis. |
MIR_VG1 |
Surface reflectance in the mid-infrared / SWIR band at about 1665 nm. |
Vegetation moisture stress, drought monitoring, burned/dry vegetation detection, land surface condition analysis. |
NDVI_VG1 |
Normalized Difference Vegetation Index computed from surface reflectance. |
Main variable for vegetation productivity, drought early warning, crop condition, rangeland monitoring and land degradation assessment. |
TOA_NDVI_VG1 |
NDVI computed from Top-of-Atmosphere reflectance before full surface correction. |
Useful for comparison with surface NDVI, quality checking and understanding atmospheric influence. |
AG_VG1 |
Aerosol optical thickness at 550 nm. |
Helps assess atmospheric contamination and supports quality screening of vegetation observations. |
OG_VG1 |
Total column ozone. |
Atmospheric correction support and quality interpretation of reflectance/NDVI values. |
WVG_VG1 |
Total column water vapour. |
Atmospheric correction support, climate-related analysis and interpretation of vegetation stress conditions. |
SZA_VG1 |
Solar zenith angle. |
Helps assess illumination conditions and filter poor-quality observations. |
SAA_VG1 |
Solar azimuth angle. |
Used with viewing geometry to understand directional reflectance effects. |
VZA_VG1 |
Viewing or observation zenith angle. |
Helps identify observations affected by extreme viewing angles. |
VAA_VG1 |
Viewing or observation azimuth angle. |
Supports bidirectional reflectance interpretation and quality filtering. |
Analysis parameters#
The following cell sets the parameters, which define the area of interest to conduct the analysis over. #### Select location To define the area of interest, there are two methods available:
By specifying the latitude, longitude, and buffer, or separate latitude and longitude buffers, this method allows you to define an area of interest around a central point. You can input the central latitude, central longitude, and a buffer value in degrees to create a square area around the center point. For example,
lat = 10.338,lon = -1.055, andbuffer = 0.1will select an area with a radius of 0.1 square degrees around the point with coordinates(10.338, -1.055).Alternatively, you can provide separate buffer values for latitude and longitude for a rectangular area. For example,
lat = 10.338,lon = -1.055, andlat_buffer = 0.1andlon_buffer = 0.08will select a rectangular area extending 0.1 degrees north and south, and 0.08 degrees east and west from the point(10.338, -1.055).For reasonable loading times, set the buffer as
0.1or lower.By uploading a polygon as a
GeoJSON or Esri Shapefile. If you choose this option, you will need to upload the geojson or ESRI shapefile into the Sandbox using Upload Files button
in the top left corner of the Jupyter Notebook interface. ESRI shapefiles must be uploaded with all the related files (.cpg, .dbf, .shp, .shx). Once uploaded, you can use the shapefile or geojson to define the area of interest. Remember to update the code to call the file you have uploaded.
To use one of these methods, you can uncomment the relevant line of code and comment out the other one. To comment out a line, add the "#" symbol before the code you want to comment out. By default, the first option which defines the location using latitude, longitude, and buffer is being used.
If running the notebook for the first time, keep the default settings below. This will demonstrate how the analysis works and provide meaningful results. The example covers a central Congo Basin belt (evergreen lowland rainforest) in Democratic Republic of the Congo (DRC), northeastern part of the country—within Tshopo Province. This region is characterised by frequent, heavy cloud cover so is often difficult to monitor effectively with optical missions that have lower revisit frequency, such as the Sentinel and Landsat missions.
To run the notebook for a different area, make sure Sentinel-3 SYN SY_2_VG1 data is available for the chosen area using the DEAfrica Explorer.
[6]:
# Method 1: Specify the latitude, longitude, and buffer)
aoi = define_area(lat=1.066, lon=26.455, buffer=0.9)
# Method 2: Use a polygon as a GeoJSON or Esri Shapefile.
# aoi = define_area(vector_path='aoi.shp')
#Create a geopolygon and geodataframe of the area of interest
geopolygon = Geometry(aoi["features"][0]["geometry"], crs="epsg:4326")
geopolygon_gdf = gpd.GeoDataFrame(geometry=[geopolygon], crs=geopolygon.crs)
# Get the latitude and longitude range of the geopolygon
lat_range = (geopolygon_gdf.total_bounds[1], geopolygon_gdf.total_bounds[3])
lon_range = (geopolygon_gdf.total_bounds[0], geopolygon_gdf.total_bounds[2])
[7]:
display_map(x=lon_range, y=lat_range)
[7]:
Load Sentinel SYN SY_2_VG1 dataset using dc.load()#
Now that we know what products and measurements are available for the product, we can load data from the datacube using dc.load. We will load data from spectral satellite bands. By specifying output_crs='EPSG:4326' and resolution=(-0.01, 0.01), we request that datacube reproject our data to the global geographic coordinate reference system (CRS), with 1000 x 1000 m pixels. Finally, group_by='solar_day' ensures that overlapping images taken within seconds of each other as the
satellite passes over are combined into a single time step in the data.
[8]:
query = {
'x': (lon_range),
'y': (lat_range),
'time':('2025-04-01', '2025-04-30'),
'output_crs': 'EPSG:4326',
'resolution': (-0.01, 0.01)}
ds_S3_syn_2_vg1 = dc.load(product=product,
group_by="solar_day",
**query)
ds_S3_syn_2_vg1
[8]:
<xarray.Dataset> Size: 21MB
Dimensions: (time: 30, latitude: 181, longitude: 181)
Coordinates:
* time (time) datetime64[ns] 240B 2025-04-01 ... 2025-04-30
* latitude (latitude) float64 1kB 1.965 1.955 1.945 ... 0.185 0.175 0.165
* longitude (longitude) float64 1kB 25.55 25.57 25.57 ... 27.34 27.36
spatial_ref int32 4B 4326
Data variables: (12/13)
AG_VG1 (time, latitude, longitude) uint8 983kB 0 0 0 0 0 ... 0 0 0 0
B0_VG1 (time, latitude, longitude) int16 2MB 429 437 388 ... 0 0 0
B2_VG1 (time, latitude, longitude) int16 2MB 783 783 723 ... 0 0 0
B3_VG1 (time, latitude, longitude) int16 2MB 4358 4391 4419 ... 0 0 0
MIR_VG1 (time, latitude, longitude) int16 2MB 2463 2441 2353 ... 0 0 0
NDVI_VG1 (time, latitude, longitude) int16 2MB 6960 6960 7200 ... 0 0 0
... ...
SAA_VG1 (time, latitude, longitude) int16 2MB 84 84 84 84 ... 0 0 0 0
SZA_VG1 (time, latitude, longitude) uint8 983kB 35 35 35 35 ... 0 0 0
TOA_NDVI_VG1 (time, latitude, longitude) int16 2MB 5520 5520 5680 ... 0 0 0
VAA_VG1 (time, latitude, longitude) int16 2MB 102 102 102 ... 0 0 0
VZA_VG1 (time, latitude, longitude) uint8 983kB 32 32 32 32 ... 0 0 0
WVG_VG1 (time, latitude, longitude) uint8 983kB 4 4 4 4 4 ... 0 0 0 0
Attributes:
crs: EPSG:4326
grid_mapping: spatial_ref⚠️ Note: Some Sentinel-3 SYN (SY_2_VG1) measurements are stored as integers with a scale factor (e.g., 0.0001) or an offset (e.g., 1000) applied. This approach reduces storage requirements while preserving precision. To obtain the correct physical values, the scale factor and/or offset must be reversed during processing.
The example below demonstrates how to apply this correction. In this dataset, several Sentinel-3 measurements use a scale factor of 10,000, as shown in the output of the List Measurements section above. You can recover the original values by either dividing by 10,000, multiplying by 0.0001, or, more generally, by retrieving the scale factor from the measurement attributes using
.attrs.get("scale_factor"), as demonstrated in the example below.
The cell below visualizes the NDVI from surface reflectance at the time of observation.
Without applying the scale factor, pay close attention to the units.
[9]:
(ds_S3_syn_2_vg1['NDVI_VG1'].isel(time=[0, 15, 20, -2])).plot(robust=True, col="time", col_wrap=4);
After applying the scale factor, pay close attention to the units.
[10]:
# Retrieve scale factor using the .attrs.get('scale_factor') measurement attribute
scale_factor = ds_S3_syn_2_vg1['NDVI_VG1'].attrs.get('scale_factor')
[11]:
(ds_S3_syn_2_vg1['NDVI_VG1'].isel(time=[0, 15, 20, -2])*scale_factor).plot(robust=True, col="time", col_wrap=4);
Before using any measurement band, always check whether a scale factor or offset needs to be applied before using the data for further analysis. After confirming this, the same process described above can be used to visualize all the other measurement bands available for the product.
Understanding BOA and TOA NDVI#
The Sentinel-3 SYN SY_2_VG1 product provides both BOA/surface NDVI (NDVI_VG1) and Top-of-Atmosphere NDVI (TOA_NDVI_VG1). BOA NDVI represents surface-level vegetation conditions after atmospheric correction, while TOA NDVI is calculated before full atmospheric correction and still includes atmospheric effects such as aerosols, water vapour, and clouds.
Comparing BOA NDVI and TOA NDVI helps users understand how atmospheric correction affects the NDVI values. In this notebook, NDVI_VG1 is used as the main vegetation layer for analysis, while TOA_NDVI_VG1 is included as a supporting layer to help users understand the difference between surface-corrected and atmosphere-influenced NDVI.
Most Sentinel-3 products include a dataMask or quality layer that helps users remove pixels affected by invalid observations, cloud, missing data, or other quality issues. However, because a readily available dataMask is not available for Sentinel-3 SYN SY_2_VG1, we need to apply a simple quality-masking approach ourselves by removing unreliable values, such as 0 NDVI values and invalid reflectance pixels, before producing maps, plots, and time-series analysis. The steps below provide the necessary guidance for applying quality masking. Users can choose the method that best suits your analysis and data requirements.
Applying the scale factor#
[12]:
ds_S3_syn_2_vg1_scaled = ds_S3_syn_2_vg1*scale_factor
The cell below extracts the surface/BOA NDVI layer (NDVI_VG1) and the Top-of-Atmosphere NDVI layer (TOA_NDVI_VG1) from the scaled Sentinel-3 SYN SY_2_VG1 dataset. It then calculates and prints the minimum and maximum values for each NDVI layer while ignoring missing values. This step helps confirm that the NDVI values are within the expected range of -1 to 1 and supports the quality-checking process before plotting or further analysis.
[13]:
boa_ndvi = ds_S3_syn_2_vg1_scaled["NDVI_VG1"]
toa_ndvi = ds_S3_syn_2_vg1_scaled["TOA_NDVI_VG1"]
print("BOA NDVI min:", float(boa_ndvi.min(skipna=True)))
print("BOA NDVI max:", float(boa_ndvi.max(skipna=True)))
print("TOA NDVI min:", float(toa_ndvi.min(skipna=True)))
print("TOA NDVI max:", float(toa_ndvi.max(skipna=True)))
BOA NDVI min: -0.08
BOA NDVI max: 0.92
TOA NDVI min: -0.08
TOA NDVI max: 0.912
This code identifies pixels where BOA NDVI_VG1 is equal to 0 and calculates how common these 0 values are for each observation date. It counts the total number of valid pixels, the number of zero-value NDVI pixels, and the percentage of valid pixels affected by 0 values. This summary helps assess whether 0 NDVI values are isolated or widespread, and supports the decision on whether they should be masked before plotting and time-series analysis.
[14]:
zero_mask = boa_ndvi == 0
zero_count = zero_mask.sum(dim=["latitude", "longitude"])
valid_count = boa_ndvi.notnull().sum(dim=["latitude", "longitude"])
zero_percent = (zero_count / valid_count) * 100
zero_summary = pd.DataFrame({
"time": boa_ndvi.time.values,
"valid_pixels": valid_count.values,
"zero_ndvi_pixels": zero_count.values,
"zero_ndvi_percent": zero_percent.values,
})
zero_summary
[14]:
| time | valid_pixels | zero_ndvi_pixels | zero_ndvi_percent | |
|---|---|---|---|---|
| 0 | 2025-04-01 | 32761 | 11609 | 35.435426 |
| 1 | 2025-04-02 | 32761 | 31803 | 97.075791 |
| 2 | 2025-04-03 | 32761 | 32761 | 100.000000 |
| 3 | 2025-04-04 | 32761 | 32729 | 99.902323 |
| 4 | 2025-04-05 | 32761 | 32761 | 100.000000 |
| 5 | 2025-04-06 | 32761 | 32761 | 100.000000 |
| 6 | 2025-04-07 | 32761 | 32761 | 100.000000 |
| 7 | 2025-04-08 | 32761 | 24863 | 75.892067 |
| 8 | 2025-04-09 | 32761 | 32761 | 100.000000 |
| 9 | 2025-04-10 | 32761 | 32761 | 100.000000 |
| 10 | 2025-04-11 | 32761 | 32761 | 100.000000 |
| 11 | 2025-04-12 | 32761 | 32761 | 100.000000 |
| 12 | 2025-04-13 | 32761 | 4462 | 13.619853 |
| 13 | 2025-04-14 | 32761 | 32761 | 100.000000 |
| 14 | 2025-04-15 | 32761 | 32761 | 100.000000 |
| 15 | 2025-04-16 | 32761 | 5325 | 16.254083 |
| 16 | 2025-04-17 | 32761 | 30805 | 94.029486 |
| 17 | 2025-04-18 | 32761 | 32761 | 100.000000 |
| 18 | 2025-04-19 | 32761 | 32761 | 100.000000 |
| 19 | 2025-04-20 | 32761 | 2429 | 7.414304 |
| 20 | 2025-04-21 | 32761 | 2298 | 7.014438 |
| 21 | 2025-04-22 | 32761 | 32761 | 100.000000 |
| 22 | 2025-04-23 | 32761 | 32761 | 100.000000 |
| 23 | 2025-04-24 | 32761 | 23801 | 72.650407 |
| 24 | 2025-04-25 | 32761 | 32714 | 99.856537 |
| 25 | 2025-04-26 | 32761 | 32761 | 100.000000 |
| 26 | 2025-04-27 | 32761 | 27359 | 83.510882 |
| 27 | 2025-04-28 | 32761 | 32761 | 100.000000 |
| 28 | 2025-04-29 | 32761 | 29783 | 90.909923 |
| 29 | 2025-04-30 | 32761 | 32761 | 100.000000 |
The next cells identifies observation dates where more than 10% of the valid BOA NDVI_VG1 pixels have values equal to 0. These dates are important because a high percentage of 0 NDVI values may indicate poor-quality observations, missing data, failed retrievals, or cloud-affected pixels. By identifying these dates, users can decide whether to mask the affected pixels, remove poor-quality dates, or investigate them further before producing NDVI plots and time-series analysis.
[15]:
high_zero_dates = zero_summary[zero_summary["zero_ndvi_percent"] > 10]
high_zero_dates
[15]:
| time | valid_pixels | zero_ndvi_pixels | zero_ndvi_percent | |
|---|---|---|---|---|
| 0 | 2025-04-01 | 32761 | 11609 | 35.435426 |
| 1 | 2025-04-02 | 32761 | 31803 | 97.075791 |
| 2 | 2025-04-03 | 32761 | 32761 | 100.000000 |
| 3 | 2025-04-04 | 32761 | 32729 | 99.902323 |
| 4 | 2025-04-05 | 32761 | 32761 | 100.000000 |
| 5 | 2025-04-06 | 32761 | 32761 | 100.000000 |
| 6 | 2025-04-07 | 32761 | 32761 | 100.000000 |
| 7 | 2025-04-08 | 32761 | 24863 | 75.892067 |
| 8 | 2025-04-09 | 32761 | 32761 | 100.000000 |
| 9 | 2025-04-10 | 32761 | 32761 | 100.000000 |
| 10 | 2025-04-11 | 32761 | 32761 | 100.000000 |
| 11 | 2025-04-12 | 32761 | 32761 | 100.000000 |
| 12 | 2025-04-13 | 32761 | 4462 | 13.619853 |
| 13 | 2025-04-14 | 32761 | 32761 | 100.000000 |
| 14 | 2025-04-15 | 32761 | 32761 | 100.000000 |
| 15 | 2025-04-16 | 32761 | 5325 | 16.254083 |
| 16 | 2025-04-17 | 32761 | 30805 | 94.029486 |
| 17 | 2025-04-18 | 32761 | 32761 | 100.000000 |
| 18 | 2025-04-19 | 32761 | 32761 | 100.000000 |
| 21 | 2025-04-22 | 32761 | 32761 | 100.000000 |
| 22 | 2025-04-23 | 32761 | 32761 | 100.000000 |
| 23 | 2025-04-24 | 32761 | 23801 | 72.650407 |
| 24 | 2025-04-25 | 32761 | 32714 | 99.856537 |
| 25 | 2025-04-26 | 32761 | 32761 | 100.000000 |
| 26 | 2025-04-27 | 32761 | 27359 | 83.510882 |
| 27 | 2025-04-28 | 32761 | 32761 | 100.000000 |
| 28 | 2025-04-29 | 32761 | 29783 | 90.909923 |
| 29 | 2025-04-30 | 32761 | 32761 | 100.000000 |
Quality Masking for Reliable NDVI Analysis#
Most Sentinel-3 products include a dataMask or quality layer that helps users remove pixels affected by invalid observations, cloud, missing data, or other quality issues. However, because a readily available dataMask is not exposed for Sentinel-3 SYN SY_2_VG1 in our current workflow, we apply a simple quality-masking approach ourselves by removing unreliable values, such as 0 NDVI values and invalid reflectance pixels, before producing maps, plots, and time-series analysis.
This code removes 0 values from the BOA NDVI_VG1 layer and compares the original NDVI image with the masked NDVI image. The purpose is to assess whether 0 values are influencing the visual interpretation of vegetation condition. In the masked image, pixels with NDVI equal to 0 are treated as no-data, allowing users to see the vegetation pattern without potentially unreliable low-value pixels affecting the plot.
[16]:
#keeps only pixels where boa_ndvi value is not equal to 0.
boa_ndvi_masked_zero = boa_ndvi.where(boa_ndvi != 0)
[17]:
#selects the first time step
time_index = 0
#plot the original against the masked data
fig, axes = plt.subplots(1, 2, figsize=(14, 5))
#original
boa_ndvi.isel(time=time_index).plot(
ax=axes[0],
robust=True,
cmap="RdYlGn"
)
axes[0].set_title("Original BOA NDVI_VG1")
#masked
boa_ndvi_masked_zero.isel(time=time_index).plot(
ax=axes[1],
robust=True,
cmap="RdYlGn"
)
axes[1].set_title("BOA NDVI_VG1 with 0 values masked")
plt.tight_layout()
plt.show()
Compare BOA NDVI VG1 with TOA NDVI VG1#
This code compares TOA_NDVI_VG1 and BOA NDVI_VG1 for a selected observation date using a scatter plot. Each point represents one pixel, with TOA NDVI shown on the x-axis and BOA NDVI shown on the y-axis. We can see that the BOA values are generally higher than TOA values for the same acquisition as values are to the top-left of the 1:1 line. This is because atmospheric corrections remove the effect of aerosols, haze, or other interference and lead to a higher, ‘cleaner’, NDVI value.
[18]:
# Select the first time step for the comparison
time_index = 0
boa_flat = boa_ndvi.isel(time=time_index).values.flatten()
toa_flat = toa_ndvi.isel(time=time_index).values.flatten()
valid = np.isfinite(boa_flat) & np.isfinite(toa_flat)
plt.figure(figsize=(6, 6))
plt.scatter(toa_flat[valid], boa_flat[valid], s=2, alpha=0.3)
plt.axhline(0, linestyle="--")
plt.axvline(0, linestyle="--")
plt.axline((0, 0), slope=1, color="black", linestyle="--", label="1:1 Line")
plt.xlabel("TOA_NDVI_VG1")
plt.ylabel("BOA NDVI_VG1")
plt.title("TOA NDVI vs BOA NDVI")
plt.grid(True)
plt.show()
Additional information#
License The code in this notebook is licensed under the Apache License, Version 2.0.
Digital Earth Africa data is licensed under the Creative Commons by Attribution 4.0 license.
Contact If you need assistance, please post a question on the DE Africa Slack channel or on the GIS Stack Exchange using the open-data-cube tag (you can view previously asked questions here).
If you would like to report an issue with this notebook, you can file one on Github.
Compatible datacube version
[19]:
print(datacube.__version__)
1.9.13
Last Tested:
[20]:
from datetime import datetime
datetime.today().strftime('%Y-%m-%d')
[20]:
'2026-09-01'