R语言|在图像上点击一个点自动显示所在多边形(面)的轮廓

一、效果如下所示：在图像上点击几个点之后，自动显示所在多边形(面)的轮廓

二、思路：

1.读取要用到的shp(shapefile)文件

2.将其画在屏幕上

3.在屏幕上点几个点

4.使用判断点是否在面内的函数进行判断并画出相应面的轮廓

三、代码：

1.载入要用到的maptools包

library(maptools)#载入maptools包
#若为第一次使用，需要先下载该包，在联网的情况下，使用
#install.packages("maptools")
#一次下载完成后，不需要再次下载，只需载入即可

2.设置工作文件夹路径

setwd("D:/R/chapter_three")
getwd()#输出工作路径，用来检验是否和设置的相同
       #如果相同，就可以继续了

3.读入shp文件，并显示在屏幕上

LB<-readShapePoly("LondonBorough.shp",verbose = T)
#因为前面已经设置了工作路径，这里不再需要加上前面的完整路径了
plot(LB)

4.在图像上点几个点并显示出来

coords<-locator(type="p")
points(coords,pch="+",col="green")
#在这里有两个小问题：
#1.把点显示出来的这句代码其实不一定需要，只是给点涂上了颜色更明显而已
#2.使用locator函数点完点之后右键“停止”即可停止采点(本操作为在R软件中的方法)

5.画出所在面的轮廓

plids<-point.In.polygon(coords,LB)
p1<-slot(LB,"polygons")
for (i in 1:length(plids)) {
  lines(p1[[plids[i]]]@Polygons[[1]]@coords,col="red")
}

 6.其中，判断点是否在面内的函数point.In.polygon()如下：

1）首先初始化装着每个点所在面的面序号的"numeric"型变量为(0，0，0)，其实在这里假设了会在屏幕上点三个点；

2）所以在后面最外面一层循环也就确定了循环为1：3，在这一层里面，点击的点是确定的，通过内层循环确定它所在的面的面序号；

3）接下来的一层循环，用来依次遍历每一个多边形，判断该点是不是在该多边形内，循环的次数也就是多边形个数。

ps：

①注意使用@来层层拆解的用法

②在这层循环中，需要用到判断点是否在多边形内的算法，参考来源为：http://t.csdn.cn/BZsrn

③判断每一个多边形经过最内层的循环之后的标记逻辑值：

若为FALSE，则继续循环

若为TRUE，则表明已经找到该点所在多边形了，把该面序号赋值给（0，0，0）中对应的变量，并用break结束该次循环，进入到下一个点的循环，寻找下一个点所在面的面序号

④最后返回要求的plids(这些输出的yes no 和over都不是必要的，只是当时为了让我自己对整个运行过程更清晰写的)

point.In.polygon<-function(point,poly.point) #point是多个点 最后的结果是一个vector
{
  plids<-c(0,0,0)
  for (i in 1:3) {#在这一层循环里面 待测点是确定的
    for (j in 1:length(poly.point@polygons)) {#在这一层里面 多边形是确定的
      n<-length(poly.point@polygons[[j]]@Polygons[[1]]@coords)/2#边界点个数
      poly.dot<-slot(poly.point@polygons[[j]]@Polygons[[1]],"coords")#边界点
      poly.dot<-data.frame(x=poly.dot[,1],y=poly.dot[,2])
      #边界点x：poly.dot[,1][]  边界点y：poly.dot[,2][]
      oddnote<-FALSE
      x<-point$x[i]
      y<-point$y[i]
      for (m in 1:k) {
        if((poly.dot$y[m]<y&&poly.dot$y[n]>=y)|(poly.dot$y[m]>=y&&poly.dot$y[n]<y))
        {
          if(poly.dot$x[m]+(y-poly.dot$y[m])/(poly.dot$y[n]-poly.dot$y[m])*(poly.dot$x[n]-poly.dot$x[m])<x)
          {
            oddnote<-!oddnote
          }
        }
        n<-m
      }
      if(oddnote==TRUE)
      {
        plids[i]<-j
        print("yes")
        print(plids[i])
        break
      }
      else{
        print("no")
      }
    }
    print("over")
  }
  plids
}

ps：关于使用@读取数据的详细介绍http://t.csdn.cn/tw3Iq