Downloader DTM 5m Regione Calabria

downloader_dtm_rc.R

download.dtm.calabria <- function(comune, qu) {
  require("RCurl", quietly = TRUE)
  require(downloader, quietly = TRUE)
  require(rgdal, quietly = TRUE)
  require(gdalUtils, quietly = TRUE)
  
  if (missing(comune))
    stop("Error: the 'comune' argument is mandatory.")
  if (missing(qu))
    stop("Error: the 'qu' argument is mandatory.")
  
  dir.out <- paste0('download/')
  dir.create(dir.out, showWarnings = FALSE)
  url <- "ftp://ftpopendata:OPENDATA2013@geoportale.regione.calabria.it/DTM5X5/"
  wfs.url <- 'WFS:http://wms.pcn.minambiente.it/ogc?map=/ms_ogc/wfs/LimitiAmministrativi_2011.map&version=1.0.0&request=GetFeature&typeName=UA.UNITAAMMINISTRATIVE.COMUNI'
  wfs.filter <- paste0(
    '&filter=<Filter><PropertyIsEqualTo><PropertyName>comune</PropertyName><Literal>',
    comune,
    '</Literal></PropertyIsEqualTo></Filter>'
    )
  
  grid <- readOGR(qu)
  proj4string(grid) <- CRS('+init=epsg:32633')
  
  ogr2ogr(
    paste0(wfs.url, wfs.filter),
    file.path(getwd(), dir.out),
    "UA.UNITAAMMINISTRATIVE.COMUNI",
    f = "ESRI Shapefile"
  )
  shp <- readOGR(paste0(dir.out, "UA.UNITAAMMINISTRATIVE.COMUNI.shp"), stringsAsFactors = FALSE)
  shp <- spTransform(shp, CRS('+init=epsg:32633'))
  
  grid.intersection <- gIntersection(grid, shp, byid = TRUE)
  grid <- grid[grid.intersection, ]
  
  cnt <- 1
  for (i in seq(dim(grid@data)[1])) {
    res <- try(download(paste0(url,
                               "grid", grid@data$elemento[i], '.img'),
                        paste0(dir.out, "grid", grid@data$elemento[i], '.img'),
                        mode = "wb"))
    if (class(res) == 'try-error') {
      print('ERROR')
    }
    print(paste0(cnt, ' of ', dim(grid@data)[1]))
    cnt <- cnt + 1
  }
  
  raster.files <- list.files(
    path = dir.out,
    full.names = TRUE,
    pattern = "^grid",
    all.files = FALSE
  )
  
  mosaic_rasters(
    gdalfile = raster.files,
    dst_dataset = paste0(dir.out, "DTM5M.tif"),
    of = "GTiff",
    a_srs = '+init=epsg:32633 +proj=utm +zone=33 +datum=WGS84 +units=m +no_defs +ellps=WGS84 '
  )
  tmp.file <-
    list.files(path = dir.out,
               pattern = "UA.UNITAAMMINISTRATIVE.COMUNI|.img",
               full.names = T)
  file.remove(tmp.file)
}