SCM

SCM Repository

[rgeos] Annotation of /pkg/src/rgeos_poly2nb.c
ViewVC logotype

Annotation of /pkg/src/rgeos_poly2nb.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 466 - (view) (download) (as text)

1 : rsbivand 59 #include "rgeos.h"
2 :    
3 :     struct ud {
4 :     int count;
5 :     int *ids;
6 :     };
7 :    
8 :     void cb(void *item, void *userdata) {
9 :     struct ud *my_UD;
10 :     my_UD = userdata;
11 : rsbivand 258 my_UD->ids[my_UD->count] = *((int *) item);
12 :     // 110904 EJP
13 : rsbivand 59 my_UD->count++;
14 :     }
15 :    
16 :     static struct ud UD;
17 :    
18 : rsbivand 61 SEXP rgeos_poly_findInBox(SEXP env, SEXP pls, SEXP as_points) {
19 : rsbivand 59
20 :     GEOSGeom *bbs;
21 :     int npls, i, j, jj, pc=0;
22 :     GEOSGeom GC, bb;
23 :     SEXP pl, bblist;
24 :     GEOSSTRtree *str;
25 :     int *icard, *ids, *oids;
26 : rsbivand 61 int asPTS = LOGICAL_POINTER(as_points)[0];
27 : rsbivand 59
28 :     GEOSContextHandle_t GEOShandle = getContextHandle(env);
29 :    
30 :     str = (GEOSSTRtree *) GEOSSTRtree_create_r(GEOShandle, (size_t) 10);
31 :    
32 :     npls = length(pls);
33 : rsbivand 258 bbs = (GEOSGeom *) R_alloc((size_t) npls, sizeof(GEOSGeom));
34 : rsbivand 59 ids = (int *) R_alloc((size_t) npls, sizeof(int));
35 :     UD.ids = (int *) R_alloc((size_t) npls, sizeof(int));
36 :     oids = (int *) R_alloc((size_t) npls, sizeof(int));
37 :    
38 :     for (i=0; i<npls; i++) {
39 :     ids[i] = i;
40 :     pl = VECTOR_ELT(pls, i);
41 : rsbivand 61 if (asPTS) {
42 :     if ((GC = rgeos_Polygons2MP(env, pl)) == NULL) {
43 :     error("rgeos_poly2nb: MP GC[%d] not created", i);
44 :     }
45 :     } else {
46 : crundel 84 if ((GC = rgeos_Polygons2geospolygon(env, pl)) == NULL) {
47 : rsbivand 61 error("rgeos_poly2nb: GC[%d] not created", i);
48 :     }
49 : rsbivand 59 }
50 :     if ((bb = GEOSEnvelope_r(GEOShandle, GC)) == NULL) {
51 :     error("rgeos_poly2nb: envelope [%d] not created", i);
52 :     }
53 :     bbs[i] = bb;
54 : rsbivand 258 GEOSSTRtree_insert_r(GEOShandle, str, bb, &(ids[i]));
55 :     // 110904 EJP
56 :     GEOSGeom_destroy_r(GEOShandle, GC);
57 : rsbivand 59 }
58 :    
59 :     icard = (int *) R_alloc((size_t) npls, sizeof(int));
60 :     PROTECT(bblist = NEW_LIST(npls-1)); pc++;
61 :    
62 :     for (i=0; i<(npls-1); i++) {
63 :     UD.count = 0;
64 :     GEOSSTRtree_query_r(GEOShandle, str, bbs[i],
65 : crundel 256 (GEOSQueryCallback) cb, &UD);
66 : rsbivand 59 for (j=0, jj=0; j<UD.count; j++) if (UD.ids[j] > i) jj++;
67 :     icard[i] = jj;
68 :     if (icard[i] > 0) SET_VECTOR_ELT(bblist, i, NEW_INTEGER(icard[i]));
69 :    
70 :     for (j=0, jj=0; j<UD.count; j++) {
71 :     if (icard[i] > 0 && UD.ids[j] > i) {
72 :     oids[jj] = UD.ids[j] + R_OFFSET;
73 :     jj++;
74 :     }
75 :     }
76 :     R_isort(oids, jj);
77 :     for (j=0; j<jj; j++) {
78 :     INTEGER_POINTER(VECTOR_ELT(bblist, i))[j] = oids[j];
79 :     }
80 :     }
81 :    
82 : rsbivand 234 for (i=0; i<npls; i++) {
83 : rsbivand 258 GEOSSTRtree_remove_r(GEOShandle, str, bbs[i], &(ids[i]));
84 :     // 110904 EJP
85 : rsbivand 234 GEOSGeom_destroy_r(GEOShandle, bbs[i]);
86 :     }
87 : rsbivand 258 GEOSSTRtree_destroy_r(GEOShandle, str);
88 : rsbivand 234
89 : rsbivand 59 UNPROTECT(pc);
90 :     return(bblist);
91 :     }
92 :    
93 : rsbivand 446 /*GEOSSTRtree *rgeos_geom2tree(GEOSContextHandle_t GEOShandle, ) {
94 :     }*/
95 : rsbivand 444
96 : rsbivand 215 SEXP rgeos_binary_STRtree_query(SEXP env, SEXP obj1, SEXP obj2) {
97 :    
98 : rsbivand 321 GEOSGeom *bbs2;
99 : rsbivand 441 int nobj1, nobj2, i, j, pc=0, isPts=FALSE;
100 : rsbivand 466 GEOSGeom GC, GCpts=NULL, bb;
101 : rsbivand 214 SEXP pl, bblist;
102 :     GEOSSTRtree *str;
103 :     int *icard, *ids, *oids;
104 : rsbivand 215 char classbuf1[BUFSIZ], classbuf2[BUFSIZ];
105 : rsbivand 259 GEOSGeom (*rgeos_xx2MP)(SEXP, SEXP);
106 : rsbivand 59
107 : rsbivand 215 strcpy(classbuf1, CHAR(STRING_ELT(GET_CLASS(VECTOR_ELT(obj1, 0)), 0)));
108 : rsbivand 260 if (!strncmp(classbuf1, "Polygons", 8))
109 : rsbivand 259 rgeos_xx2MP = rgeos_Polygons2MP;
110 : rsbivand 260 else if (!strncmp(classbuf1, "Lines", 5))
111 : rsbivand 259 rgeos_xx2MP = rgeos_Lines2MP;
112 :     else
113 :     error("rgeos_binary_STRtree_query: object class %s unknown", classbuf1);
114 : rsbivand 215
115 : rsbivand 214 GEOSContextHandle_t GEOShandle = getContextHandle(env);
116 :    
117 :     str = (GEOSSTRtree *) GEOSSTRtree_create_r(GEOShandle, (size_t) 10);
118 :    
119 : rsbivand 215 nobj1 = length(obj1);
120 : rsbivand 441
121 :     SEXP cl2 = GET_CLASS(obj2);
122 :     if (cl2 == R_NilValue) strcpy(classbuf2, "\0");
123 :     else strcpy(classbuf2, CHAR(STRING_ELT(cl2, 0)));
124 :     if ( !strcmp( classbuf2, "SpatialPoints") ||
125 :     !strcmp(classbuf2, "SpatialPointsDataFrame")) {
126 :     isPts = TRUE;
127 :     SEXP crds = GET_SLOT(obj2, install("coords"));
128 :     SEXP dim = getAttrib(crds, install("dim"));
129 :     nobj2 = INTEGER_POINTER(dim)[0];
130 :     } else {
131 :     nobj2 = length(obj2);
132 :     }
133 : rsbivand 321 bbs2 = (GEOSGeom *) R_alloc((size_t) nobj2, sizeof(GEOSGeom));
134 : rsbivand 215 ids = (int *) R_alloc((size_t) nobj1, sizeof(int));
135 : rsbivand 214
136 : rsbivand 215 UD.ids = (int *) R_alloc((size_t) nobj1, sizeof(int));
137 :     oids = (int *) R_alloc((size_t) nobj1, sizeof(int));
138 : rsbivand 214
139 : rsbivand 215 for (i=0; i<nobj1; i++) {
140 : rsbivand 214 ids[i] = i;
141 : rsbivand 215 pl = VECTOR_ELT(obj1, i);
142 : rsbivand 259 GC = rgeos_xx2MP(env, pl);
143 : rsbivand 215 if (GC == NULL) {
144 : rsbivand 214 error("rgeos_binary_STRtree_query: MP GC[%d] not created", i);
145 :     }
146 :     if ((bb = GEOSEnvelope_r(GEOShandle, GC)) == NULL) {
147 :     error("rgeos_binary_STRtree_query: envelope [%d] not created", i);
148 :     }
149 : rsbivand 258 GEOSGeom_destroy_r(GEOShandle, GC);
150 :     GEOSSTRtree_insert_r(GEOShandle, str, bb, &(ids[i]));
151 : rsbivand 214 }
152 :    
153 : rsbivand 441 if (isPts) {
154 :     GCpts = rgeos_SpatialPoints2geospoint(env, obj2);
155 :     } else {
156 :     strcpy(classbuf2, CHAR(STRING_ELT(GET_CLASS(VECTOR_ELT(obj2, 0)), 0)));
157 :     if (!strncmp(classbuf2, "Polygons", 8))
158 :     rgeos_xx2MP = rgeos_Polygons2MP;
159 :     else if (!strncmp(classbuf2, "Lines", 5))
160 :     rgeos_xx2MP = rgeos_Lines2MP;
161 :     else
162 :     error("rgeos_binary_STRtree_query: object class %s unknown",
163 :     classbuf2);
164 :     }
165 : rsbivand 259
166 : rsbivand 215 for (i=0; i<nobj2; i++) {
167 : rsbivand 441 if (isPts) {
168 :     GC = (GEOSGeom) GEOSGetGeometryN_r(GEOShandle, GCpts, i);
169 :     } else {
170 :     pl = VECTOR_ELT(obj2, i);
171 :     GC = rgeos_xx2MP(env, pl);
172 :     }
173 : rsbivand 215 if (GC == NULL) {
174 : rsbivand 441 error("rgeos_binary_STRtree_query: GC[%d] not created", i);
175 : rsbivand 214 }
176 :     if ((bb = GEOSEnvelope_r(GEOShandle, GC)) == NULL) {
177 :     error("rgeos_binary_STRtree_query: envelope [%d] not created", i);
178 :     }
179 : rsbivand 258 GEOSGeom_destroy_r(GEOShandle, GC);
180 : rsbivand 441 // Rprintf("i: %d, bb %s\n", i, GEOSGeomType_r(GEOShandle, bb));
181 : rsbivand 321 bbs2[i] = bb;
182 :     }
183 :     // 110904 EJP
184 :     icard = (int *) R_alloc((size_t) nobj2, sizeof(int));
185 :     PROTECT(bblist = NEW_LIST(nobj2)); pc++;
186 :    
187 :     for (i=0; i<nobj2; i++) {
188 : rsbivand 214 UD.count = 0;
189 : rsbivand 321 GEOSSTRtree_query_r(GEOShandle, str, bbs2[i],
190 : crundel 256 (GEOSQueryCallback) cb, &UD);
191 : rsbivand 214
192 : rsbivand 215 icard[i] = UD.count;
193 :    
194 :     if (icard[i] > 0) {
195 :     SET_VECTOR_ELT(bblist, i, NEW_INTEGER(icard[i]));
196 :    
197 :     for (j=0; j<UD.count; j++) {
198 :     oids[j] = UD.ids[j] + R_OFFSET;
199 : rsbivand 214 }
200 : rsbivand 215 R_isort(oids, UD.count);
201 :     for (j=0; j<UD.count; j++) {
202 :     INTEGER_POINTER(VECTOR_ELT(bblist, i))[j] = oids[j];
203 :     }
204 : rsbivand 214 }
205 :     }
206 :    
207 : rsbivand 234 GEOSSTRtree_destroy_r(GEOShandle, str);
208 : rsbivand 321 for (i=0; i<nobj2; i++) {
209 : rsbivand 441 GEOSGeom_destroy_r(GEOShandle, bbs2[i]);
210 : rsbivand 321 }
211 : rsbivand 234
212 : rsbivand 214 UNPROTECT(pc);
213 :     return(bblist);
214 :     }
215 :    
216 : rsbivand 215 SEXP rgeos_unary_STRtree_query(SEXP env, SEXP obj) {
217 : rsbivand 214
218 :     GEOSGeom *bbs;
219 : rsbivand 215 int nobj, i, j, jj, pc=0;
220 : rsbivand 214 GEOSGeom GC, bb;
221 :     SEXP pl, bblist;
222 :     GEOSSTRtree *str;
223 :     int *icard, *ids, *oids;
224 : rsbivand 215 char classbuf[BUFSIZ];
225 : rsbivand 259 GEOSGeom (*rgeos_xx2MP)(SEXP, SEXP);
226 : rsbivand 214
227 : rsbivand 215 strcpy(classbuf, CHAR(STRING_ELT(GET_CLASS(VECTOR_ELT(obj, 0)), 0)));
228 : rsbivand 260 if (!strncmp(classbuf, "Polygons", 8))
229 : rsbivand 259 rgeos_xx2MP = rgeos_Polygons2MP;
230 : rsbivand 260 else if (!strncmp(classbuf, "Lines", 5))
231 : rsbivand 259 rgeos_xx2MP = rgeos_Lines2MP;
232 : rsbivand 260 else if (!strncmp(classbuf, "Polygon", 7))
233 : rsbivand 259 rgeos_xx2MP = rgeos_Polygon2MP;
234 :     else
235 :     error("rgeos_binary_STRtree_query: object class %s unknown", classbuf);
236 :    
237 : rsbivand 214 GEOSContextHandle_t GEOShandle = getContextHandle(env);
238 :    
239 :     str = (GEOSSTRtree *) GEOSSTRtree_create_r(GEOShandle, (size_t) 10);
240 :    
241 : rsbivand 215 nobj = length(obj);
242 : rsbivand 258 bbs = (GEOSGeom *) R_alloc((size_t) nobj, sizeof(GEOSGeom));
243 : rsbivand 215 ids = (int *) R_alloc((size_t) nobj, sizeof(int));
244 :     UD.ids = (int *) R_alloc((size_t) nobj, sizeof(int));
245 :     oids = (int *) R_alloc((size_t) nobj, sizeof(int));
246 : rsbivand 214
247 : rsbivand 215 for (i=0; i<nobj; i++) {
248 : rsbivand 214 ids[i] = i;
249 : rsbivand 215 pl = VECTOR_ELT(obj, i);
250 : rsbivand 259 GC = rgeos_xx2MP(env, pl);
251 : rsbivand 215 if (GC == NULL) {
252 :     error("rgeos_unary_STRtree_query: MP GC[%d] not created", i);
253 : rsbivand 214 }
254 :     if ((bb = GEOSEnvelope_r(GEOShandle, GC)) == NULL) {
255 : rsbivand 215 error("rgeos_unary_STRtree_query: envelope [%d] not created", i);
256 : rsbivand 214 }
257 :     bbs[i] = bb;
258 : rsbivand 258 GEOSSTRtree_insert_r(GEOShandle, str, bb, &(ids[i]));
259 :     GEOSGeom_destroy_r(GEOShandle, GC);
260 :     // 110904 EJP
261 : rsbivand 214 }
262 :    
263 : rsbivand 215 icard = (int *) R_alloc((size_t) nobj, sizeof(int));
264 :     PROTECT(bblist = NEW_LIST(nobj-1)); pc++;
265 : rsbivand 214
266 : rsbivand 215 for (i=0; i<(nobj-1); i++) {
267 : rsbivand 214 UD.count = 0;
268 :     GEOSSTRtree_query_r(GEOShandle, str, bbs[i],
269 : crundel 256 (GEOSQueryCallback) cb, &UD);
270 : rsbivand 214 for (j=0, jj=0; j<UD.count; j++) if (UD.ids[j] > i) jj++;
271 :     icard[i] = jj;
272 : rsbivand 215 if (icard[i] > 0) {
273 :     SET_VECTOR_ELT(bblist, i, NEW_INTEGER(icard[i]));
274 : rsbivand 214
275 : rsbivand 215 for (j=0, jj=0; j<UD.count; j++) {
276 :     if (UD.ids[j] > i) {
277 :     oids[jj] = UD.ids[j] + R_OFFSET;
278 :     jj++;
279 :     }
280 : rsbivand 214 }
281 : rsbivand 215 R_isort(oids, jj);
282 :     for (j=0; j<jj; j++) {
283 :     INTEGER_POINTER(VECTOR_ELT(bblist, i))[j] = oids[j];
284 :     }
285 : rsbivand 214 }
286 :     }
287 :    
288 : rsbivand 234 for (i=0; i<nobj; i++) {
289 : rsbivand 258 GEOSSTRtree_remove_r(GEOShandle, str, bbs[i], &(ids[i]));
290 :     // 110904 EJP
291 : rsbivand 234 GEOSGeom_destroy_r(GEOShandle, bbs[i]);
292 :     }
293 :     GEOSSTRtree_destroy_r(GEOShandle, str);
294 :    
295 : rsbivand 214 UNPROTECT(pc);
296 :     return(bblist);
297 :     }
298 :    
299 :    
300 : rsbivand 215

root@r-forge.r-project.org
ViewVC Help
Powered by ViewVC 1.0.0  
Thanks to:
Vienna University of Economics and Business Powered By FusionForge