Read Sentinel-2 bands into R retrieves high values

Viewed 35

When I directly read into R the jp2 band files I get unusual high values compared to when I read the files in SNAP (version 9). To read the bands into R I use terra package (you can also use raster package) and the values range from 0 to 18000 more or less. I was wondering if SNAP is doing some conversion that I am not aware of to show values that range from 0 to 0.15 more or less.

> r_10
class       : SpatRaster 
dimensions  : 10980, 10980, 4  (nrow, ncol, nlyr)
resolution  : 10, 10  (x, y)
extent      : 499980, 609780, 6690240, 6800040  (xmin, xmax, ymin, ymax)
coord. ref. : WGS 84 / UTM zone 32N (EPSG:32632) 
sources     : T32VNN_20181018T105031_B02_10m.jp2  
              T32VNN_20181018T105031_B03_10m.jp2  
              T32VNN_20181018T105031_B04_10m.jp2  
              ... and 1 more source(s)
names       : B02_10m_m10_2018, B03_10m_m10_2018, B04_10m_m10_2018, B08_10m_m10_2018 
min values  :                0,                0,                0,                0 
max values  :            18815,            17880,            17023,            15608 
> 

I have tried to export the bands from SNAP into TIF to see if it is a problem of format but it takes forever. I was hoping that there is a convesion factor to show the actual values that I need for my analysis.

3 Answers

Legacy formats ([Geo]TIFF, JPEG) have incomplete metadata support. You may lose missing value codes, offset and scale factors, processing history, etc. NetCDF4-CF has good metadata support for applications that actually use what is available. The R terra package can import NetCDF4-CF format, but is selective about what metadata are imported. For the files I have tested, missing data, scale, and offset values are used by the rast() function, but other metadata are lost.

With terra you can also set a scale/offset if you want:

library(terra)
#terra 1.6.21
f <- system.file("ex/elev.tif", package="terra")
r <- rast(f)
r
#class       : SpatRaster 
#dimensions  : 90, 95, 1  (nrow, ncol, nlyr)
#resolution  : 0.008333333, 0.008333333  (x, y)
#extent      : 5.741667, 6.533333, 49.44167, 50.19167  (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 (EPSG:4326) 
#source      : elev.tif 
#name        : elevation 
#min value   :       141 
#max value   :       547 

scoff(r) <- cbind(10, 0)
r
#class       : SpatRaster 
#dimensions  : 90, 95, 1  (nrow, ncol, nlyr)
#resolution  : 0.008333333, 0.008333333  (x, y)
#extent      : 5.741667, 6.533333, 49.44167, 50.19167  (xmin, xmax, ymin, ymax)
#coord. ref. : lon/lat WGS 84 (EPSG:4326) 
#source      : elev.tif 
#name        : elevation 
#min value   :      1410 
#max value   :      5470 

Essentially this is way to delay the evaluation of value * scale + offset, but now it may need to be done multiple times (each time that r is used), so it is not something I would generally recommend doing.

Related