Get_Data_ita.Rmd
E’ possibile scaricare e/o estrarre i Datacube utilizzando la funzione della libreria get_coverage()
.
Tale funzione utilizza la funzione GetCoverage dei servizi WCS (Vedere l’immagine GetCoverage.png per comprendere l’analogia con l’interfaccia grafica di Rasdaman); per cui è possibile scaricare un raster o parte di esso dichiarando diverse opzioni.
NB: La funzione implementata in questa libreria non permette lo scaricamento contemporaneo di più raster, ovvero non è possibile inserire un intervallo di date o più bande. Ogni richiesta deve essere fatta per una singola data/ora e per una singola banda; il che corriponde allo scaricamento di un singolo raster.
Per scaricare invece raster multibanda utilizzare la funzione image_from_coverage()
.
get_coverage()
Attraverso l’utilizzo della funzione get_coverage()
è possibile scaricare parte di una Coverage (effettuando uno “slicing” del DataCube) in diversi formati:
La funzione è così costruita:
get_coverage(coverage, DATA, FORMAT, SUBSET_E=NULL, SUBSET_N=NULL, BAND=NULL, CRS_Extension=NULL, filename=NULL, others_opt=NULL)
E’ dunque possibile specificare:
Per maggiori in formazioni sulle opzioni della funzione GetCoverage: OGC WCS - Rasdaman Doc
1. Scaricamento di un raster (intero) specificandone solo la data
coverage="rh_ana"
data="2020-05-18T12:00:00.000Z"
formato="image/tiff"
r=get_coverage(coverage=coverage, DATA=data, FORMAT=formato)
#[1] "Messaggio informativo: Il sistema di riferimento di default è: http://localhost:8080/def/crs/EPSG/0/32632"
r
#class : RasterLayer
#dimensions : 249, 251, 62499 (nrow, ncol, ncell)
#resolution : 1000, 1000 (x, y)
#extent : 436000, 687000, 4918000, 5167000 (xmin, xmax, ymin, ymax)
#crs : +proj=utm +zone=32 +datum=WGS84 +units=m +no_defs
#source : C:/Users/sgrasso/AppData/Local/Temp/RtmpkzJWGw/file33e4114ebfe.tiff
#names : file33e4114ebfe
plot(r)
2. Scaricamento di un raster, nel sistema di riferimento MonteMario (ESPG:3003) - diverso dal sistema di riferimento della coverage, salvandolo con nome
coverage="rh_ana"
data="2020-05-18T12:00:00.000Z"
formato="image/tiff"
CRSext="http://localhost:8080/def/crs/EPSG/0/3003"
r2=get_coverage(coverage=coverage, DATA=data, FORMAT=formato, CRS_Extension=CRSext)
r2
#class : RasterLayer
#dimensions : 249, 251, 62499 (nrow, ncol, ncell)
#resolution : 1000.019, 1000.019 (x, y)
#extent : 1436026, 1687031, 4918018, 5167023 (xmin, xmax, ymin, ymax)
#crs : +proj=tmerc +lat_0=0 +lon_0=9 +k=0.9996 +x_0=1500000 +y_0=0 +ellps=intl +units=m +no_defs
#source : C:/Users/sgrasso/AppData/Local/Temp/RtmpkzJWGw/file33e472e214d.tiff
#names : file33e472e214d
plot(r2)
3. Scaricamento di un raster con opzione aggiuntiva Clipping Extension Scaricamento in una porzione del raster rh_ana (un quadrato di circa 20Km intorno a Milano Via Rosellini - EPSG 3003 (Roma40 GB Fuso Ovest): 1515230 , 5037450)
coverage="rh_ana"
data="2020-05-18T12:00:00.000Z"
formato="image/tiff"
CRSext="http://localhost:8080/def/crs/EPSG/0/3003"
more="clip=POLYGON((1505230 5027450,1525230 5027450, 1525230 5047450, 1505230 5047450, 1505230 5027450))"
r3=get_coverage(coverage=coverage, DATA=data, FORMAT=formato, CRS_Extension=CRSext, others_opt=more)
plot(r2, col='black', legend=FALSE)
plot(r3, box=FALSE, axes=FALSE, add = TRUE)
4. Scaricamento e salvataggio come immagine (jpeg) di una parte di un raster determinata da un subset lungo l’asse delle ascisse (E/X)
coverage="rh_ana"
data="2020-05-18T12:00:00.000Z"
formato="image/jpeg"
nomefile="C:\\Users\\sgrasso\\Documents\\test.jpg"
subsetE="E(436000,550000)"
get_coverage(coverage=coverage, DATA=data, FORMAT=formato, filename=nomefile, SUBSET_E=subsetE)
#[1] "Il sistema di riferimento di default è: http://localhost:8080/def/crs/EPSG/0/32632"
#[1] "L'immagine è stata salvata: C:\\Users\\sgrasso\\Documents\\test.jpg"
## Per plottare una immagine in R si può ad esempio utilizzare la libreria 'imager':
install.packages("imager")
library(imager)
library(devtools)
img <- load.image(nomefile)
# img
# Image. Width: 114 pix Height: 249 pix Depth: 1 Colour channels: 1
plot(img)
image_from_coverage()
Questa funziona permette in particolare di scaricare un raster multibanda espicitando nel parametro bands lista delle bande che si vogliono esportare. Nel caso di raster multibanda viene creato un RasterStack. Il raster viene sempre scaricato in formato image/tiff.
coords<-c("515200", "5037430") %>% as.numeric
slice_E <- as.character(c(coords[1]-10000, coords[1]+10000))
slice_N <- as.character(c(coords[2]-10000, coords[2]+10000))
data="2020-10-01"
outfile="C:\\Users\\sgrasso\\Desktop\\prova_rh_ana.tiff"
#coords
#[1] 515200 5037430
#slice_E
#[1] "505200" "525200"
#slice_N
#[1] "5027430" "5047430"
a=image_from_coverage(coverage, slice_E, slice_N, DATA=data, bands=NULL,filename=outfile, query_url=NULL)
#[1] "Il raster ha una sola banda. Verrà restituito un raster con una sola banda"
#[1] "File temporaneoC:\\Users\\sgrasso\\AppData\\Local\\Temp\\RtmpkzJWGw/rh_ana_2020-10-01.tiff"
plot(a)
Nel caso ad oggi aggiornato di Rasdaman di ARPA non vi sono salvati raster multibanda. Immaginiamo per esempio che il raster multibanda scaricato sia "myCubeR_RGB_Ortho.tif" (Immagine collezionata dalla Piattaforma NEON Airborne Observation. Ogni immagine RGB image è un raster di 3bande.
Vediamo qualche manipolazione di raster stack come proposto qui https://datacarpentry.org/r-raster-vector-geospatial/05-raster-multi-band-in-r/.
file.path="C:\\WSL_folder\\myCubeR\\vignettes\\test\\HARV_RGB_Ortho.tif"
RGB_stack_HARV <- stack(file.path)
RGB_stack_HARV
#class : RasterStack
#dimensions : 2317, 3073, 7120141, 3 (nrow, ncol, ncell, nlayers)
#resolution : 0.25, 0.25 (x, y)
#extent : 731998.5, 732766.8, 4712956, 4713536 (xmin, xmax, ymin, ymax)
#crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
#names : HARV_RGB_Ortho.1, HARV_RGB_Ortho.2, HARV_RGB_Ortho.3
#min values : 0, 0, 0
#max values : 255, 255, 255
## Per estrarre o lavorare con le diverse bande:
banda1=RGB_stack_HARV[[1]] #oppure: banda1=raster(file.path, band = 2)
banda1=RGB_stack_HARV[[2]]
banda3=RGB_stack_HARV[[3]]
banda1
#class : RasterLayer
#band : 1 (of 3 bands)
#dimensions : 2317, 3073, 7120141 (nrow, ncol, ncell)
#resolution : 0.25, 0.25 (x, y)
#extent : 731998.5, 732766.8, 4712956, 4713536 (xmin, xmax, ymin, ymax)
#crs : +proj=utm +zone=18 +datum=WGS84 +units=m +no_defs
#source : C:/WSL_folder/myCubeR/vignettes/test/HARV_RGB_Ortho.tif
#names : HARV_RGB_Ortho.1
#values : 0, 255 (min, max)
plot(banda1)
## Per creare una immagine
plotRGB(RGB_stack_HARV, r = 1, g = 2, b = 3)