由于ARCGIS在国内的领先地位,在用R语言处理空间数据时和ESRI公司独有的以数据集的形式来存储的矢量和栅格数据是不可避免的问题。而R语言现有的raster包和rgdal 包好像都无法直接导入raster文件。好在ESRI公司发布了R-ArcGIS桥接工具并开发了arcgisbinding 包。
1 安装
安装部分直接参考虾神的这篇blog就可以了。需要注意的是这里用到的几个比较关键的函数是在arcgisbinding 的1.0.1.229版本中才添加进去的,需要到这里 下载比较新的版本。
2 导入
r.file <- system.file("pictures", "cea.tif", package="rgdal") r.file %>% arc.open
直接用包自带的案例吧,还要安装一个rgdal 包,需要看到这篇东西的人应该都会装吧。arc.open 是个通用函数,打开ESRI的shp文件也是用这个。看下导入进来的信息
dataset_type : RasterDatasetpath : F:/Administrator/Documents/R/win-library/3.4/rgdal/pictures/cea.tif format : TIFFpixel_type : U8 (8bit)compression_type: Nonenrow : 515ncol : 514extent : xmin=-28493.17, ymin=4224973, xmax=2358.212, ymax=4255885WKT : PROJCS["North_American_1927_Cylindrical_Equal_Area",GEOGCS["... bands : 1 ncol nrow Band_1 514 515
数据类型是RasterDataset ,这个格式raster包好像是不支持的。另外,ESRI公司的数据都是以WKT方式来标记参考坐标系的,而R语言里比较常用的是PROJ4STRING的方式,arcgisbinding 包提供了arc.fromWktToP4 和arc.fromP4ToWkt 两个函数进行互转。
3 转换
r.file %>% arc.raster %>% as.raster
转换的代码也是非常简单的,只不过需要有一个中间过程,把RasterDataset格式的数据转换成arc.raster格式,再转换成raster包支持的rasterlayer格式。
class : RasterLayer dimensions : 515, 514, 264710 (nrow, ncol, ncell) resolution : 60.02214, 60.02214 (x, y) extent : -28493.17, 2358.212, 4224973, 4255885 (xmin, xmax, ymin, ymax) coord. ref. : +proj=cea +lon_0=-117.333333333333 +lat_ts=33.75 +x_0=0 +y_0=0 +datum=NAD27 +units=m +no_defs +ellps=clrk66 +nadgrids=@conus,@alaska,@ntv2_0.gsb,@ntv1_can.dat data source : in memorynames : Band_1 values : 0, 255 (min, max)
到这里,就全部搞定了,我们看到参考坐标系也变成了proj4string 的方式。