SENTINEL

From GRASS-Wiki
Revision as of 11:21, 29 January 2018 by ⚠️Sbl (talk | contribs) (→‎Data download: link Martins new addon(s))
Jump to navigation Jump to search

Sentinel 2

Sentinel-2 is an Earth observation mission developed by ESA as part of the Copernicus Programme to perform terrestrial observations in support of services such as forest monitoring, land cover changes detection, and natural disaster management. It consists of two identical satellites, Sentinel-2A and Sentinel-2B. Details about the mission, user and technical guides as well as other documentation can be found at the ESA dedicated website [1] and a number of other sites.

Overview of Sentinel-2 mission capabilities

Sentinel-2 mission has the following features:

  • Multi-spectral data with 13 bands in the visible, near infrared (NIR), and short wave infrared (SWIR).
  • Systematic global coverage of land surfaces from 56° S to 84° N, coastal waters, and all of the Mediterranean Sea
  • Revisiting every 5 days under the same viewing angles. At high latitudes, Sentinel-2 swath overlaps and some regions are observed twice or more every 5 days, but with different viewing angles.
  • Spatial resolutions of 10 m (Bands 2, 3, 4 and 8 - Blue, Green, Red and NIR), 20 m (Bands 5, 6 and 7 - Red Edge bands, Band 8A - Narrow infrared and Bands 11 and 12 - SWIR) and 60 m (Bands 1 - Coastal aerosol, Band 9 - Water vapour and Band 10 - SWIR/Cirrus)
  • 290 km field of view
  • Free and open data policy

Sentinel-2 products

  • Level-1B product provides radiometrically corrected imagery in Top-Of-Atmosphere (TOA) radiance values and in sensor geometry. Level-1B product is composed of an ensemble of granules that are 25 km across track (AC) by 23 km along track (AL). Additionally, this product includes the refined geometry which is used to generate the Level-1C product. Level-1B pixel coordinates refer to the centre of each pixel.
  • Level-1C product is composed of 100 km² single tiles (ortho-images in UTM/WGS84 projection). The Level-1C product results from using a Digital Elevation Model (DEM) to project the image in cartographic coordinates. Per-pixel radiometric measurements are provided in TOA reflectances with all parameters to transform them into radiances. Level-1C products are resampled with a constant Ground Sampling Distance (GSD) of 10, 20 and 60 m depending on the native resolution of the different spectral bands. In Level-1C products, pixel coordinates refer to the upper left corner of the pixel. Before September 2016, Level-1C was delivered as multi-tile product.
  • Level-2A product provides Bottom Of Atmosphere (BOA) reflectance images derived from the associated Level-1C products. Therefore, each Level-2A product is also a 100 km² tile in cartographic geometry (UTM/WGS84 projection). BOA reflectances in Level-2A products are scaled by a factor of 10000.

Data download

Sentinel-2 data might be obtained from the Copernicus Open Access Hub, either by using the online interactive interface and browsing to your area of interest or by means of the API Hub. There are several other sites from which it is possible to download Sentinel-2 imagery as well. Here, some of them:

For those preferring or needing to download data via the API, there also exist several options:

  • sentinelsat: Python utilities to access the API of Copernicus Sentinels Scientific Data Hub (useful for all Sentinel satellite imagery)
  • Sentinel-download: Automated download of Sentinel-2 L1C data from ESA (through wget)

The probably most convenient way to download (and import) Sentinel data is to use the dedicated AddOn:

It is build on top of the sentinelsat library (and also requires python-pandas).

Data pre-processing

As from May 2017, ESA started to provide Level 2A products for Europe and they will reprocess the archive during this year (2017). However, it is (still) possible to download Level-1C products and customize the atmospheric correction by using different algorithms (i.e.: 6S, 6SV, arcsi) and/or a different DEM (Sen2Cor uses STRM 90m by default).

Sentinel-2 in GRASS GIS

To exemplify the handling of Sentinel-2 data with GRASS GIS, we will use S2 L1C scene from North Carolina (NC) from November 2016. Here are the details: 

S2A_OPER_PRD_MSIL1C_PDMC_20161111T123605_R097_V20161110T160502_20161110T160502
830f2a5a-9914-427c-b4e0-69c6c14c71d4
2016-11-10T15:05:02Z

In the zipped folder delivered by ESA, the data itself is hidden under: 

S2A_OPER_PRD_MSIL1C_PDMC_20161111T123605_R097_V20161110T160502_20161110T160502.SAFE/GRANULE/S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_N02.04/IMG_DATA/

S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B01.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B02.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B03.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B04.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B05.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B06.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B07.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B08.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B09.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B10.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B11.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B12.jp2
S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B8A.jp2


Import data

Sentinel-2 data comes in JPEG2000 format that can be easily imported into GRASS GIS either by means of r.in.gdal (in a UTM location) or using r.import that provides reprojection on the fly. In this example, we will use r.import to reproject and import into the NC sample location delivered with GRASS GIS. 

for map in `ls *.jp2` ; do
 outmap=`basename $map .jp2`
 r.import input=$map output=$outmap
done

If you prefer shorter map names, then it is possible to use g.rename or the add-on g.rename.many. The other option is editing the file names before importing with sed, for example. Let us keep the S2A, date and band number parts only: 

for map in `ls *.jp2` ; do
 outmap1=`basename $map .jp2`
 outmap2=`echo $outmap1 | sed s/OPER_MSI_L1C_TL_MTI__//g | sed s/T210803_A007244_T17SPV//g`
 r.import input=$map output=$outmap2
done

The output resolution estimated by r.import is slightly different than the nominal resolution of Sentinel-2 bands. 

r.info S2A_20161110_B08 -g | grep res
nsres=9.80178836938829
ewres=9.80224757766693

We need to set the proper resolution at import by means of the resolution and resolution_value options in r.import. 

for i in 02 03 04 08 ; do
 g.message message="importing band ${i} with R 10m"
 r.import S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B${i}.jp2 out=S2A_20161110_B${i} resolution=value resolution_value=10 extent=input
done

for i in 05 06 07 8A 11 12 ; do 
  g.message message="importing band ${i} with R 20m"
  r.import S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B${i}.jp2 out=S2A_20161110_B${i} resolution=value resolution_value=20 extent=input
done

for i in 01 09 10 ; do 
 g.message message="importing band ${i} with R 60m"
 r.import S2A_OPER_MSI_L1C_TL_MTI__20161110T210803_A007244_T17SPV_B${i}.jp2 out=S2A_20161110_B${i} resolution=value resolution_value=60 extent=input
done

Unfortunately, nodata are not encoded in the Sentinel data and need to be set manually by replacing 0 with NULL. This should be done in the original projection, before reprojecting the data.

In L2A products, there are artefacts along edges of swath paths, cell values are too low, caused by a bug in sen2cor. All cells that are bordering a NULL cell must be set to NULL. This step might need to be repeated twice.

Because neighboring Sentinel tiles are overlapping each other, file sizes can be reduced by patching together neighboring tiles with the same datatake sensing start time.

Before the actual reprojection, the computational region must be set properly. In the original projection, four 10 meter cells fall exactly into one 20 meter cell, and nine 20 meter cells fall exactly into one 60 meter cell. This can be achieved by first setting the region to the input extents (obtained with r.proj -g), then aligning the region to a resolution of 60 meter, after that setting the region's resolution to the actual resolution of the input band, this time without alignment. Alternatively, a reference grid can be created and the region for the current input can then be aligned to that reference grid.

Now the data can be reprojected, preferably with method=bicubic_f or method=lanczos_f to avoid sawtooth patterns of linear features which can appear with method=nearest (default).

Atmospheric correction with i.atcorr

Example here

Alternative: Atmospheric and Radiometric Correction of Satellite Imagery (ARCSI)


Estimate diverse indices

Given that Sentinel-2 MSI comes with 3 red edge bands, it allows to estimate several other vegetation indices, aside from the most classic ones. There is a quite exhaustive list at the Sentinel Hub website.

  • EVI = 2.5*(B08 - B04) / (B08 + 6*B04 - 7.5*B02 + 1)
  • ARVI = (B08 - (2*B04 - B02)) / (B08 + (2*B04 - B02))
  • NDII = (B08 - B11) / (B08 + B11)
  • MSI = B11 / B08
  • PSRI-NIR = (B04 - B02) / B08
  • NDVI-GREEN = B03 * (B08 - B04) / (B08 + B04)
  • MTCI = (B06 - B05) / (B05 - B04)
  • IRECI = (B07 - B04) * B06 / B05
  • GRVI1 = (B04 - B03) / (B04 + B03)
  • GNDVI = (B08 - B03) / (B08 + B03)
  • EVI2 = 2.5 * (B08 - B04) / (B08 + 2.4 * B04 + 1)
  • ChlRE = B05 / B08

Classification

Example here

Further reading

IMAGICO.DE blog

Retrieved from "https://grasswiki.osgeo.org/w/index.php?title=SENTINEL&oldid=25005"