1

我有 2 个形状文件:数据(空间点数据框)和多边形(多边形数据框)。我想做一个重叠,但似乎它不起作用。

这是数据和多边形:

> data
class       : SpatialPointsDataFrame 
features    : 12527 
extent      : 10.20075, 20.6108, 54.08669, 57.75905  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs 
variables   : 3
names       :    timestamp_pretty,     timestamp,     imo 
min values  : 01/04/2006 00:00:55, 1143849655232, 9048392 
max values  : 30/04/2006 23:59:36, 1146441576823, 9191541 
> polys
class       : SpatialPolygonsDataFrame 
features    : 436375 
extent      : 4210000, 5441000, 3395000, 4813000  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs 
variables   : 2
names       :    Id, Count 
min values  :     0,     0 
max values  : 99999,     9

为了管理重叠,我使用

proj4string(data) <- proj4string(polys)                            # to confirm the same reference system
inside <- !is.na(over(data, as(polys, "SpatialPolygons")))    # overlapping shape file and data

然后mean(inside)检查多边形中点的平均值。

但是什么也没发生,平均值始终为 0。我之前多次使用它并且它总是有效,我猜它不起作用不是因为 2 个 sph 文件的范围不同。有没有办法编辑这个?

谢谢!

4

1 回答 1

0

根据我的经验,问题可能出在您分配投影的方式上。

require(geosphere)
require(rgeos)
require(rgdal)
proj4string(data) <- spTransform(polys, "+proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs")   #this is the correct way to really project the spatial object on the same projection of polys
over(data, polys)   #this is just to check whether the two sp objects interact

尝试运行此代码并告诉我是否发生了什么事。

于 2018-01-16T10:20:02.420 回答