首页
学习
活动
专区
圈层
工具
发布
社区首页 >问答首页 >用栅格软件包下载SRTM数据?

用栅格软件包下载SRTM数据?
EN

Stack Overflow用户
提问于 2017-12-28 03:03:31
回答 3查看 2.3K关注 0票数 6

我尝试在R中使用"raster“包获取SRTM数据,但是一旦我在getData命令中选择SRTM,我就会得到以下错误:

代码语言:javascript
复制
library(raster)

srtm <- getData('SRTM', lon=16, lat=48)
trying URL 'ftp://xftp.jrc.it/pub/srtmV4/tiff/srtm_40_03.zip'
trying URL 'http://hypersphere.telascience.org/elevation/cgiar_srtm_v4/tiff/zip/srtm_40_03.ZIP'
downloaded 572 bytes

Error in .SRTM(..., download = download, path = path) : file not found
In addition: Warning messages:
1: In utils::download.file(url = aurl, destfile = fn, method = "auto",  :
  URL 'ftp://xftp.jrc.it/pub/srtmV4/tiff/srtm_40_03.zip': status was 'Couldn't resolve host name'
2: In utils::unzip(zipfilename, exdir = dirname(zipfilename)) :
  error 1 in extracting from zip file

你知道这个错误是为了什么吗?

EN

回答 3

Stack Overflow用户

发布于 2018-06-13 23:25:33

我也有同样的问题,似乎是个bug。raster包中的getData函数检查三个不同url中的栅格文件的可用性。

代码语言:javascript
复制
1. ftp://xftp.jrc.it/pub/srtmV4/tiff/FILENAME
2. http://hypersphere.telascience.org/elevation/cgiar_srtm_v4/tiff/zip/FILENAME
3. http://srtm.csi.cgiar.org/SRT-ZIP/SRTM_V41/SRTM_Data_GeoTiff/FILENAME

它们中的前两个(截至今天)不工作或无法访问。然而,由于某些原因,有一小部分数据正在通过服务器传输,因此该包假定它是一个可用文件,但最终导致utils出错。然而,第三个url是三个url中最可信的一个。

我做了一些挖掘,在稍微修改了raster包本身,使其使用第三个url之后,我想出了以下函数。您可以在此处输入LongitudeLatitude值。请注意,只有当您想要下载基于经纬度的文件时,此选项才有用。

代码语言:javascript
复制
SRTM<-function(lon, lat) {
  stopifnot(lon >= -180 & lon <= 180)
  stopifnot(lat >= -60 & lat <= 60)
  rs <- raster(nrows=24, ncols=72, xmn=-180, xmx=180, ymn=-60, ymx=60 )
  rowTile <- rowFromY(rs, lat)
  colTile <- colFromX(rs, lon)
  if (rowTile < 10) { rowTile <- paste('0', rowTile, sep='') }
  if (colTile < 10) { colTile <- paste('0', colTile, sep='') }

  f <- paste('srtm_', colTile, '_', rowTile, sep="")
  theurl <- paste("http://srtm.csi.cgiar.org/wp-content/uploads/files/srtm_5x5/TIFF/", f, ".ZIP", sep="")
  utils::download.file(url=theurl, destfile='srtm_40_0.zip', method="auto", quiet = FALSE, mode = "wb", cacheOK = TRUE)
}

示例:

代码语言:javascript
复制
SRTM(lon=16, lat=48)

这将在您的文件夹中生成一个名为srtm_40_03.zip的文件,其中通常包含一个tif、一个hdr和一个同名的tfw文件。像往常一样使用它们进行进一步的处理。

EDIT DATE 22-JAN-19:更改的srtm链接(以及),上面的代码已经适应以反映这一点。

票数 4
EN

Stack Overflow用户

发布于 2018-06-18 22:33:34

当我尝试使用raster::getData()下载SRTM数据时,我遇到了完全相同的错误。正如samAct所提到的,这个问题与不能再访问的网址有关。

我的解决方案是编辑getData函数,使其首先尝试http://srtm.csi.cgiar.org/SRT-ZIP/SRTM_V41/SRTM_Data_GeoTiff/FILENAME URL。

要使用我编辑过的光栅包版本,只需使用`devtools::install_github(“pokyah”)即可。这将覆盖您当前的版本,并使getData与有效的URL一起工作。

代码语言:javascript
复制
devtools::install_github("pokyah/raster")
library(raster)

然后,您可以使用此函数下载您感兴趣区域的SRTM数据:

代码语言:javascript
复制
# @param country_code.chr a character specifying the ISO contrycode. Ex : BE for belgium
# @param NAME_1.chr a character specifying the NAME_1 value for lower than country level information
# @param aggregation_factor.num a numeric specifying the aggregation factor to get the desired spatial resolution
# @param EPSG.chr a character specifying the EPSG code of the desired Coordiante Reference System (CRS)
# @param path.chr a character specifying the path where to dowload the SRTM data

build_highRes_terrain_rasters.fun <- function(country_code.chr, NAME_1.chr=NULL, aggregation_factor.num=NULL, EPSG.chr=NULL, path.chr) {
  # Path to downloaded SRTM Tiles refs
  srtm.tiles.ref <- raster::shapefile("<PATH_TO_DOWNLOADED_TILES_REFS")

  # Get country geometry first
  extent.sp <- raster::getData('GADM', country=country_code.chr, level=1)

  if(!is.null(NAME_1.chr)){
    extent.sp <- subset(extent.sp, NAME_1 == NAME_1.chr)
  }

  # Intersect extent geometry and tile grid
  intersects <- rgeos::gIntersects(extent.sp, srtm.tiles.ref, byid=T)
  tiles      <- srtm.tiles.ref[intersects[,1],]

  # Download tiles using getData
  # inspired from https://www.gis-blog.com/download-srtm-for-an-entire-country/
  srtm_list  <- list()
  for(i in 1:length(tiles)) {
    lon <- extent(tiles[i,])[1]  + (extent(tiles[i,])[2] - extent(tiles[i,])[1]) / 2
    lat <- extent(tiles[i,])[3]  + (extent(tiles[i,])[4] - extent(tiles[i,])[3]) / 2

    tile <- getData('SRTM', #data are downloaded from http://www.cgiar-csi.org/. See getData do of pokyah/raster repo on github
                    lon=lon,
                    lat=lat,
                    download = FALSE,
                    path = path.chr)

    srtm_list[[i]] <- tile
  }

  # Mosaic tiles
  srtm_list$fun <- mean
  srtm_mosaic.ras <- do.call(raster::mosaic, srtm_list)

  # Crop tiles to extent borders
  extent.elevation.ras <- raster::crop(srtm_mosaic.ras, extent.sp)
  extent.elevation.ras <- raster::mask(extent.elevation.ras, extent.sp)

  # transform to desired CRS
  if(!is.null(EPSG.chr)){
    raster::projectRaster(extent.elevation.ras, crs = toString((dplyr::filter(rgdal::make_EPSG(), code==EPSG.chr))$prj4))
  }

  # aggregate to lower resolution
  # inspired from https://stackoverflow.com/questions/32278825/how-to-change-the-resolution-of-a-raster-layer-in-r
  if(!is.null(aggregation_factor.num)){
    extent.elevation.ras <- raster::aggregate(extent.elevation.ras, fact=aggregation_factor.num)
  }

  # compute the slope from the elevation
  # inspired from https://rpubs.com/etiennebr/visualraster
  extent.slope.ras <- raster::terrain(extent.elevation.ras, opt="slope", unit="degrees")
  extent.aspect.ras <- raster::terrain(extent.elevation.ras, opt="aspect", unit="degrees")
  extent.roughness.ras <- raster::terrain(extent.elevation.ras, opt="roughness")

  # compute the aspect from the elevation
  extent.terrain.ras = raster::stack(
    extent.elevation.ras,
    extent.slope.ras,
    extent.aspect.ras,
    extent.roughness.ras)
}

希望它能帮上忙!

票数 0
EN

Stack Overflow用户

发布于 2019-01-22 06:29:33

今天,2019年1月21日,链接仍然断开,即使是被pokyah更正的链接。新的工作版本位于:

代码语言:javascript
复制
devtools::install_github("fedefyco/raster")
library(raster)
票数 0
EN
页面原文内容由Stack Overflow提供。腾讯云小微IT领域专用引擎提供翻译支持
原文链接:

https://stackoverflow.com/questions/47997347

复制
相关文章

相似问题

领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档