1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
|
point.in.polygon = function(point.x, point.y, pol.x, pol.y,
mode.checked=FALSE) {
if (mode.checked) res <- .Call("R_point_in_polygon_sp", point.x,
point.y, pol.x, pol.y, PACKAGE = "sp")
else res <- .Call("R_point_in_polygon_sp", as.numeric(point.x),
as.numeric(point.y), as.numeric(pol.x), as.numeric(pol.y),
PACKAGE = "sp")
res
}
pointsInPolygon = function(pts, Polygon,
mode.checked=FALSE) {
pts = coordinates(pts)
cc = slot(Polygon, "coords")
point.in.polygon(pts[,1], pts[,2], cc[,1], cc[,2],
mode.checked=mode.checked)
}
pointsInPolygons = function(pts, Polygons, #which = FALSE,
mode.checked=FALSE) {
rings = slot(Polygons, "Polygons")
res = matrix(unlist(lapply(rings, function(x, pts)
pointsInPolygon(pts, x, mode.checked=mode.checked),
pts = pts)), ncol=length(rings))
res <- res > 0
# holes <- sapply(rings, function(y) slot(y, "hole"))
# areas <- sapply(rings, function(x) slot(x, "area"))
# if (any(holes) && any(res[,holes])) {
# holerows <- which(res[,holes,drop=FALSE], arr.ind=TRUE)[,1]
# odd <- rowSums(res[holerows,,drop=FALSE])%%2 != 0
# for (i in seq(along = holerows)) {
# in_p <- which.min(areas[res[holerows[i],,drop=FALSE]])
# res[holerows[i],] <- FALSE
# if (odd[i]) res[holerows[i], in_p] <- TRUE
# }
# res[,holes] <- FALSE
# }
# revised 100716
ret <- rowSums(res) %% 2 != 0
# ret <- apply(res, 1, any)
# if (which) {
# reta <- integer(length(ret))
# for (i in seq(along = ret)) {
# if (ret[i]) reta[i] <- which(res[i,])
# else reta[i] <- as.integer(NA)
# }
# ret <- reta
# }
ret
}
#pointsInSpatialPolygons = function(pts, SpPolygons) {
# sr = slot(SpPolygons, "polygons")
# res = lapply(sr, function(x, pts) pointsInPolygons(pts, x), pts = pts)
# ret = rep(as.numeric(NA), nrow(coordinates(pts)))
# for (i in seq(along = res))
# ret[res[[i]] > 0] = i
# ret
#}
pointsInSpatialPolygons = function(pts, SpPolygons) {
pls = slot(SpPolygons, "polygons")
lb <- lapply(pls, function(x) as.double(bbox(x)))
cpts <- coordinates(pts)
storage.mode(cpts) <- "double"
mode.checked <- storage.mode(cpts) == "double"
cand0 <- .Call("pointsInBox", lb, cpts[,1], cpts[,2], PACKAGE="sp")
m <- length(pls)
cand <- .Call("tList", cand0, as.integer(m), PACKAGE="sp")
rm(cand0)
gc(verbose=FALSE)
res <- pointsInPolys2(pls, cand, cpts, mode.checked=mode.checked)
res
}
pointsInPolys2 <- function(pls, cand, pts, mode.checked=FALSE) {
n <- nrow(pts)
res <- rep(as.integer(NA), n)
for (i in seq(along=cand)) {
candi <- cand[[i]]
if (length(candi) > 0) {
ptsi <- pts[candi,,drop=FALSE]
ret <- pointsInPolygons(ptsi, pls[[i]], mode.checked=mode.checked)
for (j in seq(along=candi)) {
jj <- candi[j]
if (is.na(res[jj])) res[jj] <- ifelse(ret[j], i,
as.integer(NA))
}
}
}
res
}
|
