下面是利用這個數據進行的基本處理,尤其是白化操作(即將繪製區域之外的部分抹去),供大家參考!
# rm(list=ls())
# setwd('d:/zlh/zlh cs/Rmap/')
# library('geoR')
# library('maps')
# library('GISTools')
# library("akima")
# library("tseries")
# library('sp')
# library('ggplot2')
# library('ggmap')
# # library(animation)
# # library(XML)
#生成R色卡
# pdf.color="bar.pdf"
# pdf(pdf.color,height=120)
# par(mar=c(0,10,3,0)+0.1,yaxs='i')
# barplot(rep(1,length(colors())),
# col=rev(colors()),
# names.arg=rev(colors()),horiz=T,las=1,xaxt='n',
# main=expression("Bars of different colors in function"
# ~italic(colors())))
# dev.off()
# shell.exec(pdf.color)
#1---關於中國底圖繪製:填色底圖:polygon(lon,lat,col='red') #polygon {graphics}
#1-1 從世界地圖中提取,但存在邊界不準確問題 map {maps}
# map("world","China")
# map中提取經緯度數據方式
# borders=map("world","China")
# borders=data.frame(lon=borders$x,lat=borders$y)
# naloc=which(is.na(borders[,1])==T)
# borders=borders[-naloc,]
#1-2 從shapefile中提取,bou1_4l.shp、bou2_4l表示帶省界和不帶省界底圖readShapeLines {maptools}
# zg=readShapeLines('Rmap/bou1_4l.shp');# readShapePoly 讀取多邊形文件
# plot(c(70,140),c(18,55),type='n')
# plot(zg,add=T,col='blue',lty=1)
# map.axes()
# zg=readShapeLines('Rmap/bou2_4l.shp')
# plot(c(70,140),c(18,55),type='n')
# plot(zg,add=T,col='blue',lty=1)
#1-3 從surfer文件提取的邊界經緯度散點值,這種包含了右下角南海區域
# map_china=read.table("Rmap/china_map.txt",header=F)
# plot(map_china,col='black',type='l',lty=1,lwd=1.5,xlim=c(75,134),ylim=c(18.7,55),xlab='',ylab='')
# 116.863415373, 20.6994643144
# 135.062889024, 20.1632372666---參考坐標系(用於後期平移南海縮放圖到合適區域
#2---閾值場空間分布圖
#1-4 從shape文件中提取想要的區域邊界經緯度文件
#1-4-1 從輪廓文件*l.shp中提取邊界線
# zg=readShapeLines('Rmap/bou1_4l.shp');
# names(zg);
# unique(zg$GBCODE)
# plot(china_bj)
# plot(zg[zg$GBCODE==61020,],col='red')
# zg=readShapeLines('Rmap/bou2_4l.shp');
# names(zg)
# unique(zg$GBCODE)
# plot(china_bj)
# plot(zg[zg$GBCODE==61020,],col='red')
#1-4-2 從輪廓文件*p.shp中提取省區邊界線
# zg=readShapePoly('Rmap/bou2_4p.shp');
# names(zg)
# plot(china_bj)
# plot(zg[zg$GBCODE==61020,],col='red')
# plot(zg[zg$NAME=="河北省",])
#2---白化功能:關於在給定邊界情形下選取內部點,外部點設為NA值
# #2-1 均勻化格點產生方式
# #2-1-1 expand.grid {base} 類似於排列組合 combn函數
# gr=expand.grid(lon=seq(73.37,136,by=0.5),lat=seq(17.5,53.54,by=0.5))
# #2-1-2 pred_grid {geoR}
# gr=pred_grid(rplot[,1:2],by=0.5);names(gr)=c("lon","lat")
# #2-2 判斷位於邊界之內的格點,並給出對應位置
# #2-2-0 封閉邊界線:
# borders=read.table('Rmap/china_bj_mainland_simple.txt',header=F)
# #plot(borders,type='l')
# names(borders)=c("lon","lat")
# #2-2-1 locations.inside {geoR}
# gr.in1 <- locations.inside(gr,borders)
# loc1=paste(as.character(gr$lon),as.character(gr$lat),sep='+')
# loc2=paste(as.character(gr.in1$lon),as.character(gr.in1$lat),sep='+')
# loc=pmatch(loc2,loc1,dup=T)
# #2-2-2 polygrid {geoR} vec=T 控制輸出位於邊界之內的點所在位置
# gr.in2<-polygrid(seq(73.569,134.569,by=0.5),seq(20.33,53.33,by=0.5),borders,vec=F)
# loc=which(gr.in2$vec.inout==T)
# #2-2-3 over {sp} 克服上述不能用於多個封閉區域的缺點,通過形狀文件來進行
# # 2-2-3-1 通過提取的txt的surfer格式數據來進行
# mainland=read.table('Rmap/china_bj_mainland.txt',header=F,sep='')
# hainan=read.table('Rmap/china_bj_hainan.txt',header=F,sep='')
# taiwan=read.table('Rmap/china_bj_taiwan.txt',header=F,sep='')
# china_bj = SpatialPolygons(list(
# Polygons(list(Polygon(mainland)), ID="mainland"),
# Polygons(list(Polygon(hainan)), ID="hainan"),
# Polygons(list(Polygon(taiwan)), ID="taiwan")))
# # plot(china_bj)
# z=over(pts,china_bj);locout=which(is.na(z)==T)
# plot(map_china,col='black',type='l',lty=1,lwd=1.5,xlim=c(75,134),ylim=c(18.7,53.5),xlab='',ylab='')
# points(gr[-locout,],pch=21,bg='red',cex=0.25)
# 2-2-3-2 通過提取shp格式文件數據來進行
# china_bj=readShapePoly("Rmap/bou1_4p.shp")
# pts = SpatialPoints(gr);
# z=over(pts,china_bj);locout=which(is.na(z[,1])==T)
# plot(map_china,col='black',type='l',lty=1,lwd=1.5,xlim=c(75,134),ylim=c(18.7,53.5),xlab='',ylab='')
# points(gr,pch=21,bg='black',cex=0.005)
# points(gr[-locout,],pch=21,bg='red',cex=0.25)
附件說明:
連結:https://pan.baidu.com/s/1CmQaEY9XkUFYHhTkRHNUrg 密碼:lp2k
#1--- ArcGis中國區域shapefile文件,bou1_4*,bou2_4*分別為包含省界和不包含省界的底圖
[1] "bou1_4l.dbf"
[2] "bou1_4l.shp"
[3] "bou1_4l.shx"
[4] "bou1_4p.dbf"
[5] "bou1_4p.shp"
[6] "bou1_4p.shx"
[7] "bou2_4l.dbf"
[8] "bou2_4l.shp"
[9] "bou2_4l.shx"
[10] "bou2_4p.dbf"
[11] "bou2_4p.shp"
[12] "bou2_4p.shx"
#2--- R 軟體專用的中國底圖
[13] "Rmap_1_china_map_with_southsea_revision.txt" :
精細版中國底圖,包含南海縮小右放的底圖和長江、黃河(修正了surfer版本中長江、黃河不夠準確的情況)
[14] "Rmap_2_china_map_china_no_southsea_revision.txt" :
同上,但不包含南海縮放小圖
[15] "Rmap_3_southsea_simple_smallsize_box_revision.txt":
依據上述文件提取的單獨南海縮放小圖底圖,依據需要通過左右增加位移量放置不同位置,參考坐標為:【17.5,29.751】;【126.0148,136】
[16] "Rmap_4_changjiang_revision.txt":
精細版長江底圖
[17] "Rmap_5_huanghe_revision.txt":
精細版黃河底圖
#3--- R 軟體使用白化封閉區域,依據Rmap_1_china_map_with_southsea_revision.txt提取的各個區域單獨文件
[18] "Rmap_6_1_china_bj_mainland_revision.txt" :大陸底圖
[19] "Rmap_6_2_china_bj_hainan_revision.txt": 海南底圖
[20] "Rmap_6_3_china_bj_taiwan_revision.txt": 臺灣底圖
#4--- R 提取的簡單化邊界和同比例南海底圖
[21] "Rmap_7_southsea_comp_revision.txt":精細化同比例南海底圖,缺點是正常的九段線變為八段線,少了一個,因為同比例,需要特殊處理才能縮放用
[22] "Rmap_8_china_bj_mainland_simple_revision.txt":簡單化大陸底圖,因為只有大陸,所以白化使用起來還是需要結合前面海南、臺灣等封閉線進行