空间数据导入导出

空间数据导入导出.....................................................金百锁/////////////////////////////////////////////////////2017年2月24日1.坐标参考系CRS.....................................................ESPG清单1ESPG是欧洲石油调查组织，现在是石油和天然气生产商测量及定位委员会。在1986年开始收集大地测量参数的数据集。2rgdal包提供了ESPG清单的一个副本，有PROJ.4样式描述，允许基准和投影转换。基准转换基于WGS84座标系。2grep函数从rgdal中读取ESPG清单。3/23空间数据导入导出ÍÏ使用ESPG清单library(rgdal)EPSG-make_EPSG()EPSG[grep(^#ED50$,EPSG$note),]codenote1594230#ED50prj4159+proj=longlat+ellps=intl+towgs84=-87,-98,-121,0,0,0,0+no_defs4/23空间数据导入导出ÍÏPROJ.4规范CRS(+init=epsg:4230)CRSarguments:+init=epsg:4230+proj=longlat+ellps=intl+towgs84=-87,-98,-121,0,0,0,0+no_defsED50-CRS(+init=epsg:4230+towgs84=-87,-96,-120,0,0,0,0)5/23空间数据导入导出ÍÏ投影和座标变换IJ.east-as(char2dms(4d31\’00\E),numeric)IJ.north-as(char2dms(52d28\’00\N),numeric)IJ.ED50-SpatialPoints(cbind(x=IJ.east,y=IJ.north),proj4string=ED50)res-spTransform(IJ.ED50,CRS(+proj=longlat+datum=WGS84))6/23空间数据导入导出ÍÏ座标距离（椭圆上大弧距离）coordinates(IJ.ED50)xy[1,]4.51666752.46667numeric)coordinates(res)xy[1,]4.51535952.46589spDistsN1(coordinates(IJ.ED50),coordinates(res),longlat=TRUE)*1000[1]124.13727/23空间数据导入导出ÍÏ坐标方位library(maptools)gzAzimuth(coordinates(IJ.ED50),coordinates(res))x-134.36748/23空间数据导入导出ÍÏ美国国家地图投影EPSG[grep(Atlas,EPSG$note),1:2]codenote6392163#USNationalAtlasEqualArea23413978#NAD83/CanadaAtlasLambert23423979#NAD83(CSRS)/CanadaAtlasLambertCRS(+init=epsg:2163)CRSarguments:+init=epsg:2163+proj=laea+lat_0=45+lon_0=-100+x_0=0+y_0=0+a=6370997+b=6370997+units=m+no_defs9/23空间数据导入导出ÍÏ中国投影EPSG[grep(China,EPSG$note),1:2]codenote2544490#ChinaGeodeticCoordinateSystem200018823415#WGS72BE/SouthChinaSeaLambert48364479#ChinaGeodeticCoordinateSystem200010/23空间数据导入导出ÍÏ兰伯特等积方位投影proj-projInfo(proj)proj[proj$name==laea,]namedescription53laeaLambertAzimuthalEqualAreaellps-projInfo(ellps)ellps[grep(a=6370997,ellps$major),]namemajorell43spherea=6370997.0b=6370997.0descriptionNormalSphere(r=6370997)11/23空间数据导入导出ÍÏDMS类（度、分、秒）IJ.dms.E-4d31\’00\EIJ.dms.N-52d28\’00\NIJ_east-char2dms(IJ.dms.E)IJ_north-char2dms(IJ.dms.N)IJ_east[1]4d31’EIJ_north[1]52d28’NgetSlots(DMS)WSdegminsecNSlogicalnumericnumericnumericlogical12/23空间数据导入导出ÍÏ2.矢量文件格式.....................................................shapefile格式1shapefile格式由ESRI规定，至少包含三个文件.shp,.shx,.dbf.2.shp为几何文件。3.shx为指向几何文件的索引文件。4.dbf为存储属性数据的DBF文件。14/23空间数据导入导出ÍÏreadOGR函数1readOGR包含两个参数。数据源名称（dsn）和图层（layer）2dsn通常包含导入的三个shapefile文件目录名（当前工作目录用’.’表示）3layer没有.shp扩展名的shapefile文件名。15/23空间数据导入导出ÍÏ苏格兰shapefile文件scot_LL-readOGR(dsn=.,layer=scot,integer64=allow.loss)OGRdatasourcewithdriver:ESRIShapefileSource:.,layer:scotwith56featuresIthas2fieldsInteger64fieldsreadassigned32-bitintegers:IDproj4string(scot_LL)[1]NAproj4string(scot_LL)-CRS(+proj=longlat+ellps=WGS84)16/23空间数据导入导出ÍÏ苏格兰唇癌数据library(gcmr)data(scotland)scot_dat-scotlandscot_dat$District=1:nrow(scot_dat)names(scot_dat)observedexpectedAFFlatitudelongitudeDistrict17/23空间数据导入导出ÍÏ唇癌数据导入多边形ID_D-match(scot_LL$ID,scot_dat$District)scot_dat1-scot_dat[ID_D,]row.names(scot_dat1)-row.names(scot_LL)library(maptools)scot_LLa-spCbind(scot_LL,scot_dat1)all.equal(scot_LLa$ID,scot_LLa$District)[1]TRUE18/23空间数据导入导出ÍÏ相对风险与经验贝叶斯稳定相对风险library(spdep)O-scot_LLa$observedE-scot_LLa$expectedscot_LLa$SMR-probmap(O,E)$relRisk/100library(DCluster)scot_LLa$smth-empbaysmooth(O,E)$smthrr19/23空间数据导入导出ÍÏ画图scot_BNG-spTransform(scot_LLa,CRS(+init=epsg:27700))library(RColorBrewer)spplot(scot_BNG,c(SMR,smth),+at=c(0,0.25,0.5,0.8,1,1.5,2.5,4.5,7),+col.regions=rev(brewer.pal(8,RdBu)))20/23空间数据导入导出ÍÏ比较图21/23空间数据导入导出ÍÏ导出Shapefile文件数据drv-ESRIShapefilewriteOGR(scot_BNG,dsn=.,layer=scot_BNG,driver=drv,overwrite_layer=TRUE)22/23空间数据导入导出ÍÏ导出KLM文件数据KLM(KeyholeMarkupLanguage)文件，可以覆盖GoogleEarthwriteOGR(scot_LLa[ID],dsn=scot_district.kml,layer=borders,driver=KML)llCRS-CRS(+proj=longlat+ellps=WGS84)scot_SP_LL-SpatialPointsDataFrame(coordinates(scot_LLa),proj4string=llCRS,data=as(scot_LLa,data.frame))writeOGR(scot_SP_LL,dsn=scot_rates.kml,layer=rates,driver=KML)23/23空间数据导入导出ÍÏ