The example below is a full time-series data download, import and plot workflow in 4 steps : get the URL with the package getremotedata, download with the package httr, import and plot with the package raster. Packages purrr and sf are also used.
Data of interest are :
- Elevation (SRTM DEM) ;
- Daily rainfall (TAMSAT)
- Monthly nighttime lights (VIIRS_DNB_MONTH)
- Daily moonlight (MIRIADE)
- Hourly wind speed and direction (ERA5)
Region of interest :
Preliminary step : call the packages and set-up the region and the time range of interest
roi <- st_as_sf(data.frame(geom = "POLYGON ((-5.82 9.54, -5.42 9.55, -5.41 8.84, -5.81 8.84, -5.82 9.54))"), wkt = "geom", crs = 4326)
time_range <- as.Date(c("2017-01-01","2017-01-30"))1. Get the URLs of the data for a given ROI and time frame
roi <- st_as_sf(data.frame(geom = "POLYGON ((-5.82 9.54, -5.42 9.55, -5.41 8.84, -5.81 8.84, -5.82 9.54))"), wkt = "geom", crs = 4326)
time_range <- as.Date(c("2017-01-01","2017-01-30"))
# SRTM
strm_urls <- grd_get_url(collection = "SRTMGL1.003", roi = roi)
print(str(strm_urls))
#> 'data.frame': 2 obs. of 4 variables:
#> $ time_start: logi NA NA
#> $ name : chr "N08W006.SRTMGL1.hgt.zip" "N09W006.SRTMGL1.hgt.zip"
#> $ url : chr "http://e4ftl01.cr.usgs.gov/MEASURES/SRTMGL1.003/2000.02.11/N08W006.SRTMGL1.hgt.zip" "http://e4ftl01.cr.usgs.gov/MEASURES/SRTMGL1.003/2000.02.11/N09W006.SRTMGL1.hgt.zip"
#> $ destfile : chr "SRTMGL1.003/N08W006.SRTMGL1.hgt.zip" "SRTMGL1.003/N09W006.SRTMGL1.hgt.zip"
#> NULL
# TAMSAT
tamsat_urls <- grd_get_url(collection = "TAMSAT",variables = c("daily_rainfall_estimate"), time_range = time_range)
print(str(tamsat_urls))
#> 'data.frame': 30 obs. of 4 variables:
#> $ time_start: Date, format: "2017-01-30" "2017-01-29" ...
#> $ name : chr "rfe2017_01_30.v3.nc" "rfe2017_01_29.v3.nc" "rfe2017_01_28.v3.nc" "rfe2017_01_27.v3.nc" ...
#> $ url : chr "https://www.tamsat.org.uk/public_data/TAMSAT3/2017/01/rfe2017_01_30.v3.nc" "https://www.tamsat.org.uk/public_data/TAMSAT3/2017/01/rfe2017_01_29.v3.nc" "https://www.tamsat.org.uk/public_data/TAMSAT3/2017/01/rfe2017_01_28.v3.nc" "https://www.tamsat.org.uk/public_data/TAMSAT3/2017/01/rfe2017_01_27.v3.nc" ...
#> $ destfile : chr "TAMSAT/rfe2017_01_30.v3.nc" "TAMSAT/rfe2017_01_29.v3.nc" "TAMSAT/rfe2017_01_28.v3.nc" "TAMSAT/rfe2017_01_27.v3.nc" ...
#> NULL
# VIIRS_DNB_MONTH
viirsdnb_urls <- grd_get_url(collection = "VIIRS_DNB_MONTH",variables = c("Monthly_AvgRadiance"), roi = roi, time_range = time_range)
print(str(viirsdnb_urls))
#> 'data.frame': 1 obs. of 4 variables:
#> $ time_start: chr "2017-01-01"
#> $ name : chr "Monthly_AvgRadiance_201701.tif"
#> $ url : chr "https://gis.ngdc.noaa.gov/arcgis/rest/services/NPP_VIIRS_DNB/Monthly_AvgRadiance/ImageServer/exportImage?bbox=-"| __truncated__
#> $ destfile : chr "VIIRS_DNB_MONTH/Monthly_AvgRadiance_201701.tif"
#> NULL
# MIRIADE
imcce_urls <- grd_get_url(collection = "MIRIADE", roi = roi, time_range = time_range)
print(str(imcce_urls))
#> 'data.frame': 30 obs. of 4 variables:
#> $ time_start: Date, format: "2017-01-30" "2017-01-29" ...
#> $ name : chr "20170130.csv" "20170129.csv" "20170128.csv" "20170127.csv" ...
#> $ url : chr "http://vo.imcce.fr/webservices/miriade/ephemcc_query.php?-name=s:Moon&-type=Satellite&-ep=2017-01-30T23:30:00&-"| __truncated__ "http://vo.imcce.fr/webservices/miriade/ephemcc_query.php?-name=s:Moon&-type=Satellite&-ep=2017-01-29T23:30:00&-"| __truncated__ "http://vo.imcce.fr/webservices/miriade/ephemcc_query.php?-name=s:Moon&-type=Satellite&-ep=2017-01-28T23:30:00&-"| __truncated__ "http://vo.imcce.fr/webservices/miriade/ephemcc_query.php?-name=s:Moon&-type=Satellite&-ep=2017-01-27T23:30:00&-"| __truncated__ ...
#> $ destfile : chr "MIRIADE/20170130.csv" "MIRIADE/20170129.csv" "MIRIADE/20170128.csv" "MIRIADE/20170127.csv" ...
#> NULL
# ERA5
# ERA5 is an hourly database, so we keep only 1 day (i.e. not the whole month)
era5_urls <- grd_get_url(collection = "ERA5", variables = c("10m_u_component_of_wind","10m_v_component_of_wind"), roi = roi, time_range = as.Date(c("2017-01-01","2017-01-02")))2. Download the data
# Create directories if they do not exist
unique(dirname(data_to_dl$destfile)) %>% lapply(dir.create,recursive = TRUE)
# SRTM
# Login to Earthdata servers is needed to download STRM data. To create an account go to : https://urs.earthdata.nasa.gov/.
# Here we have stored our credentials in local environment variables
username <- Sys.getenv("earthdata_un")
password <- Sys.getenv("earthdata_pw")
srtm_dl <- map2(strm_urls$url,strm_urls$destfile,~GET(url = .x, write_disk(.y), progress() , authenticate(username,password)))
# TAMSAT
tamsat_dl <- map2(tamsat_urls$url,tamsat_urls$destfile,~GET(url = .x, write_disk(.y), progress()))
# VIIRS_DNB_MONTH
viirsdnb_dl <- map2(viirsdnb_urls$url,viirsdnb_urls$destfile,~GET(url = .x, write_disk(.y), progress()))
# MIRIADE
imcce_dl <- map2(imcce_urls$url,imcce_urls$destfile,~GET(url = .x, write_disk(.y), progress()))
# ERA5
# For ERA5 we must use a specific function to download the data
era5_dl <- grd_download_data_era5(era5_urls)3. Import the data in R
# SRTM
rast_srtm <- strm_urls$destfile %>%
map(~unzip(., exdir = dirname(.))) %>%
map(~raster(.)) %>%
do.call(merge,.) %>%
crop(roi)
# TAMSAT
rast_tamsat <- tamsat_urls$destfile %>%
map(~raster(.)) %>%
map(~crop(.,roi)) %>%
brick()
names(rast_tamsat) <- tamsat_urls$time_start
# VIIRS_DNB_MONTH
rasts_viirs <- viirsdnb_urls$destfile %>%
map(~raster(.)) %>%
brick(.)
names(rasts_viirs) <- viirsdnb_urls$time_start
# MIRIADE
miriade_dfs <- imcce_urls$destfile %>%
map(~read.csv(.,skip=10))
# ERA5
rasts_era5 <- era5_urls$destfile %>%
map(~raster(., varname = "u10")) %>%
brick(.)
names(rasts_era5) <- era5_urls$time_start
rast_srtm
#> class : RasterLayer
#> dimensions : 2556, 1475, 3770100 (nrow, ncol, ncell)
#> resolution : 0.0002777778, 0.0002777778 (x, y)
#> extent : -5.819861, -5.410139, 8.839861, 9.549861 (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
#> data source : in memory
#> names : layer
#> values : 268, 575 (min, max)
rast_tamsat
#> class : RasterBrick
#> dimensions : 19, 11, 209, 30 (nrow, ncol, ncell, nlayers)
#> resolution : 0.0375, 0.0375 (x, y)
#> extent : -5.83125, -5.41875, 8.83125, 9.54375 (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
#> data source : in memory
#> names : X2017.01.30, X2017.01.29, X2017.01.28, X2017.01.27, X2017.01.26, X2017.01.25, X2017.01.24, X2017.01.23, X2017.01.22, X2017.01.21, X2017.01.20, X2017.01.19, X2017.01.18, X2017.01.17, X2017.01.16, ...
#> min values : 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
#> max values : 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 5.8, 0.0, 0.0, 0.0, ...
rasts_viirs
#> class : RasterBrick
#> dimensions : 400, 400, 160000, 1 (nrow, ncol, ncell, nlayers)
#> resolution : 0.001775, 0.001775 (x, y)
#> extent : -5.97, -5.26, 8.84, 9.55 (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0
#> data source : in memory
#> names : X2017.01.01
#> min values : -0.08
#> max values : 63.86
print(str(miriade_dfs[1:3]))
#> List of 3
#> $ :'data.frame': 1 obs. of 10 variables:
#> ..$ X..Date.UTC..JD. : num 2457784
#> ..$ RA..h.m.s. : Factor w/ 1 level "23 12 23.31871": 1
#> ..$ DE..deg.arcmin.arcs : Factor w/ 1 level "-06 18 29.6643": 1
#> ..$ Distance..au. : num 0.00257
#> ..$ V.Mag : num -6.84
#> ..$ Phase..o. : num 144
#> ..$ Sun.elong..o. : num 35.6
#> ..$ muRAcosDE..arcsec.min.: num 43.1
#> ..$ muDE..arcsec.min. : num 11.2
#> ..$ Dist_dot..km.s. : num 0.248
#> $ :'data.frame': 1 obs. of 10 variables:
#> ..$ X..Date.UTC..JD. : num 2457783
#> ..$ RA..h.m.s. : Factor w/ 1 level "22 21 8.46964": 1
#> ..$ DE..deg.arcmin.arcs : Factor w/ 1 level "-10 10 5.8931": 1
#> ..$ Distance..au. : num 0.0026
#> ..$ V.Mag : num -5.26
#> ..$ Phase..o. : num 157
#> ..$ Sun.elong..o. : num 23.4
#> ..$ muRAcosDE..arcsec.min.: num 44.5
#> ..$ muDE..arcsec.min. : num 10.2
#> ..$ Dist_dot..km.s. : num 0.167
#> $ :'data.frame': 1 obs. of 10 variables:
#> ..$ X..Date.UTC..JD. : num 2457782
#> ..$ RA..h.m.s. : Factor w/ 1 level "21 29 36.14850": 1
#> ..$ DE..deg.arcmin.arcs : Factor w/ 1 level "-13 27 59.9202": 1
#> ..$ Distance..au. : num 0.00262
#> ..$ V.Mag : num -2.36
#> ..$ Phase..o. : num 168
#> ..$ Sun.elong..o. : num 11.5
#> ..$ muRAcosDE..arcsec.min.: num 45.3
#> ..$ muDE..arcsec.min. : num 8.32
#> ..$ Dist_dot..km.s. : num 0.079
#> NULL
rasts_era5
#> class : RasterBrick
#> dimensions : 3, 2, 6, 25 (nrow, ncol, ncell, nlayers)
#> resolution : 0.25, 0.25 (x, y)
#> extent : -5.945, -5.445, 8.715, 9.465 (xmin, xmax, ymin, ymax)
#> coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
#> data source : in memory
#> names : X2017.01.02.00.00.00, X2017.01.01.23.00.00, X2017.01.01.22.00.00, X2017.01.01.21.00.00, X2017.01.01.20.00.00, X2017.01.01.19.00.00, X2017.01.01.18.00.00, X2017.01.01.17.00.00, X2017.01.01.16.00.00, X2017.01.01.15.00.00, X2017.01.01.14.00.00, X2017.01.01.13.00.00, X2017.01.01.12.00.00, X2017.01.01.11.00.00, X2017.01.01.10.00.00, ...
#> min values : -0.8141584, -1.1385984, -1.2836151, -0.9246990, -0.8993572, -0.7861232, -0.7462901, -0.8951994, -1.0976095, -1.3974714, -1.8156195, -2.2225685, -2.5027256, -2.4738798, -2.3199224, ...
#> max values : -0.3400831, -0.6051512, -0.8189163, -0.6034481, -0.6186794, -0.4978389, -0.5015925, -0.7188536, -0.8169317, -1.1053572, -1.5913153, -2.0937271, -2.3445568, -2.1918898, -1.9231634, ...4. Display / Plot
# MIRIADE
print(str(miriade_dfs[1:3]))
#> List of 3
#> $ :'data.frame': 1 obs. of 10 variables:
#> ..$ X..Date.UTC..JD. : num 2457784
#> ..$ RA..h.m.s. : Factor w/ 1 level "23 12 23.31871": 1
#> ..$ DE..deg.arcmin.arcs : Factor w/ 1 level "-06 18 29.6643": 1
#> ..$ Distance..au. : num 0.00257
#> ..$ V.Mag : num -6.84
#> ..$ Phase..o. : num 144
#> ..$ Sun.elong..o. : num 35.6
#> ..$ muRAcosDE..arcsec.min.: num 43.1
#> ..$ muDE..arcsec.min. : num 11.2
#> ..$ Dist_dot..km.s. : num 0.248
#> $ :'data.frame': 1 obs. of 10 variables:
#> ..$ X..Date.UTC..JD. : num 2457783
#> ..$ RA..h.m.s. : Factor w/ 1 level "22 21 8.46964": 1
#> ..$ DE..deg.arcmin.arcs : Factor w/ 1 level "-10 10 5.8931": 1
#> ..$ Distance..au. : num 0.0026
#> ..$ V.Mag : num -5.26
#> ..$ Phase..o. : num 157
#> ..$ Sun.elong..o. : num 23.4
#> ..$ muRAcosDE..arcsec.min.: num 44.5
#> ..$ muDE..arcsec.min. : num 10.2
#> ..$ Dist_dot..km.s. : num 0.167
#> $ :'data.frame': 1 obs. of 10 variables:
#> ..$ X..Date.UTC..JD. : num 2457782
#> ..$ RA..h.m.s. : Factor w/ 1 level "21 29 36.14850": 1
#> ..$ DE..deg.arcmin.arcs : Factor w/ 1 level "-13 27 59.9202": 1
#> ..$ Distance..au. : num 0.00262
#> ..$ V.Mag : num -2.36
#> ..$ Phase..o. : num 168
#> ..$ Sun.elong..o. : num 11.5
#> ..$ muRAcosDE..arcsec.min.: num 45.3
#> ..$ muDE..arcsec.min. : num 8.32
#> ..$ Dist_dot..km.s. : num 0.079
#> NULL
# ERA5
plot(rasts_era5)