Process on Cube

Mattia Rossi

Import Libraries

library(CubeR)
library(tibble)
library(magrittr)
library(raster)
library(tidyr)
library(lubridate)
library(purrr)
library(MonalisR)
library(dplyr)
library(ggplot2)

The norm_diff_pixel() function allows to calculate any Normalized Difference (e.g. NDVI or NDSI) between two bands. The following example extracts the NDVI from one Pixel over the complete Temporal Extent of the Data Cube.

coverage<-getCapability()[11]
bands <- coverage_get_bands(coverage = coverage)
coords<-c("620812.02", "5171499.36") %>% as.numeric
NIR <- bands[5]
Red <- bands[4]

nd1.tab <- norm_diff_pixel(coverage, coords, band1=NIR, band2=Red, date = NULL,plot = F) %>% as.tibble
## Warning: `as.tibble()` is deprecated, use `as_tibble()` (but mind the new semantics).
## This warning is displayed once per session.
## # A tibble: 171 x 2
##    times                        res
##    <fct>                      <dbl>
##  1 2015-06-27T00:00:00.000Z 0.498  
##  2 2015-07-04T00:00:00.000Z 0.459  
##  3 2015-07-24T00:00:00.000Z 0.799  
##  4 2015-08-03T00:00:00.000Z 0.924  
##  5 2015-08-06T00:00:00.000Z 0.908  
##  6 2015-08-13T00:00:00.000Z 0.905  
##  7 2015-08-16T00:00:00.000Z 0.0717 
##  8 2015-08-23T00:00:00.000Z 0.00936
##  9 2015-08-26T00:00:00.000Z 0.905  
## 10 2015-09-02T00:00:00.000Z 0.0261 
## # ... with 161 more rows

The norm_diff_raster() function allows to calculate any Normalized Difference (e.g. NDVI or NDSI) between two bands. It is similar to the previous function but works spatially at one date. It is possible to extract the whole image or only a subset

coverage<-getCapability()[11]
NIR <- bands[5]
Red <- bands[4]
bands <- coverage_get_bands(coverage = coverage)
coords<-c("620812.02", "5171499.36") %>% as.numeric
slice_E <- as.character(c(coords[1]-5000, coords[1]+5000))
slice_N <- as.character(c(coords[2]-5000, coords[2]+5000))

times <- coverage_get_timestamps(coverage = coverage)
date<-times[8]

nd1.tiff <- norm_diff_raster(coverage, slice_E, slice_N, date=date, band1=NIR, band2=Red)

nd1.tiff
## class      : RasterLayer 
## dimensions : 1001, 1001, 1002001  (nrow, ncol, ncell)
## resolution : 9.99001, 9.99001  (x, y)
## extent     : 615812, 625812, 5166499, 5176499  (xmin, xmax, ymin, ymax)
## crs        : +init=epsg:32632 +proj=utm +zone=32 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0 
## source     : memory
## names      : layer 
## values     : 0.4352941, 0.8980392  (min, max)

In the last examples we combined the functionality of NDVI computation provided by this package with an in-situ database implemented at EURAC. This so called MONALISA Database holds NDVI measurements in 15-minute temporal resolution and can be accessed with a second R-package developed at EURAC Research.

# Download Sentinel Data
coverage<-getCapability()[11] # S2 10m Reflectance
bands <- coverage_get_bands(coverage = coverage)
coords<-c("620812.02", "5171499.36") %>% as.numeric

nd1.tab <- norm_diff_pixel(coverage, coords, band1= bands[5], band2= bands[4], plot = F) %>% as.tibble

# Download InSitu Data
mnls_down<-downloadMonalisa(foi="vimes1500",
                            datestart = "2016-03-01 00:00",
                            dateend = "2016-12-31 00:00",
                            property = "Normalized Difference Vegetation Index - average",
                            procedure = "SpectralReflectanceSensor_vimes1500") %>% .[[1]] %>% as.tibble

# Format the Datasets
ndvi.MONA<-mnls_down %>%
  tidyr::separate(TimeStamp,c("Date","Time"),sep=" ") %>% 
  dplyr::filter(.,Time=="10:30:00") %>% 
  dplyr::select(-Time) %>% 
  dplyr::mutate(Date=as_date(Date)) %>% 
  tibble::add_column(Sensor="MONALISA NDVI")

ndvi.S2<-nd1.tab %>% 
  dplyr::mutate(Date=as_date(times)) %>% 
  tibble::add_column(Value=nd1.tab$res) %>% 
  tibble::add_column(Sensor="Raw S2 NDVI 10m") %>% 
  dplyr::select(-c(times,res)) %>% 
  dplyr::filter(year(Date)=="2016" & month(Date)>2)

# Bind
allndvi<-rbind(ndvi.MONA,ndvi.S2) %>% arrange(Date)
allndvi
## # A tibble: 364 x 3
##    Date         Value Sensor         
##    <date>       <dbl> <chr>          
##  1 2016-03-01 -0.038  MONALISA NDVI  
##  2 2016-03-02  0.307  MONALISA NDVI  
##  3 2016-03-03 -0.069  MONALISA NDVI  
##  4 2016-03-03  0.0186 Raw S2 NDVI 10m
##  5 2016-03-04 -0.038  MONALISA NDVI  
##  6 2016-03-05 -0.074  MONALISA NDVI  
##  7 2016-03-06 -0.099  MONALISA NDVI  
##  8 2016-03-07 -0.076  MONALISA NDVI  
##  9 2016-03-08 -0.036  MONALISA NDVI  
## 10 2016-03-09 -0.019  MONALISA NDVI  
## # ... with 354 more rows