source: grass/trunk/lib/python/pygrass/vector/__init__.py

Last change on this file was 74506, checked in by annakrat, 5 years ago

pygrass: fixing doctests (only for Python 3)
  • Property svn:eol-style set to native
  • Property svn:mime-type set to text/x-python
File size: 32.0 KB
Line 
1# -*- coding: utf-8 -*-
2from __future__ import print_function
3
4from os.path import join, exists
5import grass.lib.gis as libgis
6libgis.G_gisinit('')
7import grass.lib.vector as libvect
8import ctypes
9
10#
11# import pygrass modules
12#
13from grass.pygrass.vector.vector_type import VTYPE
14from grass.pygrass.errors import GrassError, must_be_open
15from grass.pygrass.gis import Location
16
17from grass.pygrass.vector.geometry import GEOOBJ as _GEOOBJ
18from grass.pygrass.vector.geometry import read_line, read_next_line
19from grass.pygrass.vector.geometry import Area as _Area
20from grass.pygrass.vector.abstract import Info
21from grass.pygrass.vector.basic import Bbox, Cats, Ilist
22
23
24_NUMOF = {"areas": libvect.Vect_get_num_areas,
25 "dblinks": libvect.Vect_get_num_dblinks,
26 "faces": libvect.Vect_get_num_faces,
27 "holes": libvect.Vect_get_num_holes,
28 "islands": libvect.Vect_get_num_islands,
29 "kernels": libvect.Vect_get_num_kernels,
30 "points": (libvect.Vect_get_num_primitives, libvect.GV_POINT),
31 "lines": (libvect.Vect_get_num_primitives, libvect.GV_LINE),
32 "centroids": (libvect.Vect_get_num_primitives, libvect.GV_CENTROID),
33 "boundaries": (libvect.Vect_get_num_primitives, libvect.GV_BOUNDARY),
34 "nodes": libvect.Vect_get_num_nodes,
35 "updated_lines": libvect.Vect_get_num_updated_lines,
36 "updated_nodes": libvect.Vect_get_num_updated_nodes,
37 "volumes": libvect.Vect_get_num_volumes}
38
39# For test purposes
40test_vector_name = "vector_doctest_map"
41
42#=============================================
43# VECTOR
44#=============================================
45
46class Vector(Info):
47 """Vector class is the grass vector format without topology
48
49 >>> from grass.pygrass.vector import Vector
50 >>> test_vect = Vector(test_vector_name)
51 >>> test_vect.is_open()
52 False
53 >>> test_vect.mapset
54 ''
55 >>> test_vect.exist()
56 True
57 >>> test_vect.overwrite
58 False
59
60 """
61 def __init__(self, name, mapset='', *args, **kwargs):
62 # Set map name and mapset
63 super(Vector, self).__init__(name, mapset, *args, **kwargs)
64 self._topo_level = 1
65 self._class_name = 'Vector'
66 self.overwrite = False
67 self._cats = []
68
69 def __repr__(self):
70 if self.exist():
71 return "%s(%r, %r)" % (self._class_name, self.name, self.mapset)
72 else:
73 return "%s(%r)" % (self._class_name, self.name)
74
75 def __iter__(self):
76 """::
77
78 >>> test_vect = Vector(test_vector_name)
79 >>> test_vect.open(mode='r')
80 >>> features = [feature for feature in test_vect]
81 >>> features[:3]
82 [Point(10.000000, 6.000000), Point(12.000000, 6.000000), Point(14.000000, 6.000000)]
83 >>> test_vect.close()
84
85 ..
86 """
87 #return (self.read(f_id) for f_id in xrange(self.num_of_features()))
88 return self
89
90 @must_be_open
91 def __next__(self):
92 """::
93
94 >>> test_vect = Vector(test_vector_name)
95 >>> test_vect.open(mode='r')
96 >>> test_vect.next()
97 Point(10.000000, 6.000000)
98 >>> test_vect.next()
99 Point(12.000000, 6.000000)
100 >>> test_vect.close()
101
102 ..
103 """
104 return read_next_line(self.c_mapinfo, self.table, self.writeable,
105 is2D=not self.is_3D())
106
107 @must_be_open
108 def next(self):
109 return self.__next__()
110
111 @must_be_open
112 def rewind(self):
113 """Rewind vector map to cause reads to start at beginning."""
114 if libvect.Vect_rewind(self.c_mapinfo) == -1:
115 raise GrassError("Vect_rewind raise an error.")
116
117 @must_be_open
118 def write(self, geo_obj, cat=None, attrs=None):
119 """Write geometry features and attributes.
120
121 :param geo_obj: a geometry grass object define in
122 grass.pygrass.vector.geometry
123 :type geo_obj: geometry GRASS object
124 :param attrs: a list with the values that will be insert in the
125 attribute table.
126 :type attrs: list
127 :param cat: The category of the geometry feature, otherwise the
128 c_cats attribute of the geometry object will be used.
129 :type cat: integer
130
131 Open a new vector map ::
132
133 >>> new = VectorTopo('newvect')
134 >>> new.exist()
135 False
136
137 define the new columns of the attribute table ::
138
139 >>> cols = [(u'cat', 'INTEGER PRIMARY KEY'),
140 ... (u'name', 'TEXT')]
141
142 open the vector map in write mode
143
144 >>> new.open('w', tab_name='newvect', tab_cols=cols)
145
146 import a geometry feature ::
147
148 >>> from grass.pygrass.vector.geometry import Point
149
150 create two points ::
151
152 >>> point0 = Point(0, 0)
153 >>> point1 = Point(1, 1)
154
155 then write the two points on the map, with ::
156
157 >>> new.write(point0, cat=1, attrs=('pub',))
158 >>> new.write(point1, cat=2, attrs=('resturant',))
159
160 commit the db changes ::
161
162 >>> new.table.conn.commit()
163 >>> new.table.execute().fetchall()
164 [(1, 'pub'), (2, 'resturant')]
165
166 close the vector map ::
167
168 >>> new.close()
169 >>> new.exist()
170 True
171
172 then play with the map ::
173
174 >>> new.open(mode='r')
175 >>> new.read(1)
176 Point(0.000000, 0.000000)
177 >>> new.read(2)
178 Point(1.000000, 1.000000)
179 >>> new.read(1).attrs['name']
180 'pub'
181 >>> new.read(2).attrs['name']
182 'resturant'
183 >>> new.close()
184 >>> new.remove()
185
186 """
187 self.n_lines += 1
188 if not isinstance(cat, int) and not isinstance(cat, str):
189 # likely the case of using 7.0 API
190 import warnings
191 warnings.warn("Vector.write(geo_obj, attrs=(...)) is"
192 " depreciated, specify cat explicitly",
193 DeprecationWarning)
194 # try to accommodate
195 attrs = cat
196 cat = None
197 if attrs and cat is None:
198 # TODO: this does not work as expected when there are
199 # already features in the map when we opened it
200 cat = (self._cats[-1] if self._cats else 0) + 1
201
202 if cat is not None and cat not in self._cats:
203 self._cats.append(cat)
204 if self.table is not None and attrs is not None:
205 attr = [cat, ]
206 attr.extend(attrs)
207 cur = self.table.conn.cursor()
208 cur.execute(self.table.columns.insert_str, attr)
209 cur.close()
210
211 if cat is not None:
212 cats = Cats(geo_obj.c_cats)
213 cats.reset()
214 cats.set(cat, self.layer)
215
216 if geo_obj.gtype == _Area.gtype:
217 result = self._write_area(geo_obj)
218 result = libvect.Vect_write_line(self.c_mapinfo, geo_obj.gtype,
219 geo_obj.c_points, geo_obj.c_cats)
220 if result == -1:
221 raise GrassError("Not able to write the vector feature.")
222 if self._topo_level == 2:
223 # return new feature id (on level 2)
224 geo_obj.id = result
225 else:
226 # return offset into file where the feature starts (on level 1)
227 geo_obj.offset = result
228
229 @must_be_open
230 def has_color_table(self):
231 """Return if vector has color table associated in file system;
232 Color table stored in the vector's attribute table well be not checked
233
234 >>> test_vect = Vector(test_vector_name)
235 >>> test_vect.open(mode='r')
236 >>> test_vect.has_color_table()
237 False
238
239 >>> test_vect.close()
240 >>> from grass.pygrass.utils import copy, remove
241 >>> copy(test_vector_name,'mytest_vect','vect')
242 >>> from grass.pygrass.modules.shortcuts import vector as v
243 >>> v.colors(map='mytest_vect', color='population', column='value')
244 Module('v.colors')
245 >>> mytest_vect = Vector('mytest_vect')
246 >>> mytest_vect.open(mode='r')
247 >>> mytest_vect.has_color_table()
248 True
249 >>> mytest_vect.close()
250 >>> remove('mytest_vect', 'vect')
251 """
252 loc = Location()
253 path = join(loc.path(), self.mapset, 'vector', self.name, 'colr')
254 return True if exists(path) else False
255
256
257#=============================================
258# VECTOR WITH TOPOLOGY
259#=============================================
260
261class VectorTopo(Vector):
262 """Vector class with the support of the GRASS topology.
263
264 Open a vector map using the *with statement*: ::
265
266 >>> with VectorTopo(test_vector_name, mode='r') as test_vect:
267 ... for feature in test_vect[:7]:
268 ... print(feature.attrs['name'])
269 ...
270 point
271 point
272 point
273 line
274 line
275 line
276 >>> test_vect.is_open()
277 False
278
279 ..
280 """
281 def __init__(self, name, mapset='', *args, **kwargs):
282 super(VectorTopo, self).__init__(name, mapset, *args, **kwargs)
283 self._topo_level = 2
284 self._class_name = 'VectorTopo'
285
286 def __len__(self):
287 return libvect.Vect_get_num_lines(self.c_mapinfo)
288
289 def __getitem__(self, key):
290 """::
291
292 >>> test_vect = VectorTopo(test_vector_name)
293 >>> test_vect.open(mode='r')
294 >>> test_vect[:4]
295 [Point(10.000000, 6.000000), Point(12.000000, 6.000000), Point(14.000000, 6.000000)]
296 >>> test_vect.close()
297
298 ..
299 """
300 if isinstance(key, slice):
301 return [self.read(indx)
302 for indx in range(key.start if key.start else 1,
303 key.stop if key.stop else len(self),
304 key.step if key.step else 1)]
305 elif isinstance(key, int):
306 return self.read(key)
307 else:
308 raise ValueError("Invalid argument type: %r." % key)
309
310 @must_be_open
311 def num_primitive_of(self, primitive):
312 """Return the number of primitive
313
314 :param primitive: the name of primitive to query; the supported values are:
315
316 * *boundary*,
317 * *centroid*,
318 * *face*,
319 * *kernel*,
320 * *line*,
321 * *point*
322 * *area*
323 * *volume*
324
325 :type primitive: str
326
327 ::
328
329 >>> test_vect = VectorTopo(test_vector_name)
330 >>> test_vect.open(mode='r')
331 >>> test_vect.num_primitive_of('point')
332 3
333 >>> test_vect.num_primitive_of('line')
334 3
335 >>> test_vect.num_primitive_of('centroid')
336 4
337 >>> test_vect.num_primitive_of('boundary')
338 11
339 >>> test_vect.close()
340
341 ..
342 """
343 return libvect.Vect_get_num_primitives(self.c_mapinfo,
344 VTYPE[primitive])
345
346 @must_be_open
347 def number_of(self, vtype):
348 """Return the number of the chosen element type
349
350 :param vtype: the name of type to query; the supported values are:
351 *areas*, *dblinks*, *faces*, *holes*, *islands*,
352 *kernels*, *points*, *lines*, *centroids*, *boundaries*,
353 *nodes*, *line_points*, *update_lines*, *update_nodes*,
354 *volumes*
355 :type vtype: str
356
357 >>> test_vect = VectorTopo(test_vector_name)
358 >>> test_vect.open(mode='r')
359 >>> test_vect.number_of("areas")
360 4
361 >>> test_vect.number_of("islands")
362 2
363 >>> test_vect.number_of("holes")
364 0
365 >>> test_vect.number_of("lines")
366 3
367 >>> test_vect.number_of("nodes")
368 15
369 >>> test_vect.number_of("pizza")
370 ... # doctest: +ELLIPSIS +NORMALIZE_WHITESPACE
371 Traceback (most recent call last):
372 ...
373 ValueError: vtype not supported, use one of: 'areas', ...
374 >>> test_vect.close()
375
376
377 ..
378 """
379 if vtype in _NUMOF.keys():
380 if isinstance(_NUMOF[vtype], tuple):
381 fn, ptype = _NUMOF[vtype]
382 return fn(self.c_mapinfo, ptype)
383 else:
384 return _NUMOF[vtype](self.c_mapinfo)
385 else:
386 keys = "', '".join(sorted(_NUMOF.keys()))
387 raise ValueError("vtype not supported, use one of: '%s'" % keys)
388
389 @must_be_open
390 def num_primitives(self):
391 """Return dictionary with the number of all primitives
392 """
393 output = {}
394 for prim in VTYPE.keys():
395 output[prim] = self.num_primitive_of(prim)
396 return output
397
398 @must_be_open
399 def viter(self, vtype, idonly=False):
400 """Return an iterator of vector features
401
402 :param vtype: the name of type to query; the supported values are:
403 *areas*, *dblinks*, *faces*, *holes*, *islands*,
404 *kernels*, *line_points*, *lines*, *nodes*, *points*,
405 *update_lines*, *update_nodes*, *volumes*
406 :type vtype: str
407 :param idonly: variable to return only the id of features instead of
408 full features
409 :type idonly: bool
410
411 >>> test_vect = VectorTopo(test_vector_name, mode='r')
412 >>> test_vect.open(mode='r')
413 >>> areas = [area for area in test_vect.viter('areas')]
414 >>> areas[:3]
415 [Area(1), Area(2), Area(3)]
416
417
418 to sort the result in a efficient way, use: ::
419
420 >>> from operator import methodcaller as method
421 >>> areas.sort(key=method('area'), reverse=True) # sort the list
422 >>> for area in areas[:3]:
423 ... print(area, area.area())
424 Area(1) 12.0
425 Area(2) 8.0
426 Area(4) 8.0
427
428 >>> areas = [area for area in test_vect.viter('areas')]
429 >>> for area in areas:
430 ... print(area.centroid().cat)
431 3
432 3
433 3
434 3
435
436 >>> test_vect.close()
437 """
438 if vtype in _GEOOBJ.keys():
439 if _GEOOBJ[vtype] is not None:
440 ids = (indx for indx in range(1, self.number_of(vtype) + 1))
441 if idonly:
442 return ids
443 return (_GEOOBJ[vtype](v_id=indx, c_mapinfo=self.c_mapinfo,
444 table=self.table,
445 writeable=self.writeable)
446 for indx in ids)
447 else:
448 keys = "', '".join(sorted(_GEOOBJ.keys()))
449 raise ValueError("vtype not supported, use one of: '%s'" % keys)
450
451 @must_be_open
452 def rewind(self):
453 """Rewind vector map to cause reads to start at beginning. ::
454
455 >>> test_vect = VectorTopo(test_vector_name)
456 >>> test_vect.open(mode='r')
457 >>> test_vect.next()
458 Point(10.000000, 6.000000)
459 >>> test_vect.next()
460 Point(12.000000, 6.000000)
461 >>> test_vect.next()
462 Point(14.000000, 6.000000)
463 >>> test_vect.rewind()
464 >>> test_vect.next()
465 Point(10.000000, 6.000000)
466 >>> test_vect.close()
467
468 ..
469 """
470 libvect.Vect_rewind(self.c_mapinfo)
471
472 @must_be_open
473 def cat(self, cat_id, vtype, layer=None, generator=False, geo=None):
474 """Return the geometry features with category == cat_id.
475
476 :param cat_id: the category number
477 :type cat_id: int
478 :param vtype: the type of geometry feature that we are looking for
479 :type vtype: str
480 :param layer: the layer number that will be used
481 :type layer: int
482 :param generator: if True return a generator otherwise it return a
483 list of features
484 :type generator: bool
485 """
486 if geo is None and vtype not in _GEOOBJ:
487 keys = "', '".join(sorted(_GEOOBJ.keys()))
488 raise ValueError("vtype not supported, use one of: '%s'" % keys)
489 Obj = _GEOOBJ[vtype] if geo is None else geo
490 ilist = Ilist()
491 libvect.Vect_cidx_find_all(self.c_mapinfo,
492 layer if layer else self.layer,
493 Obj.gtype, cat_id, ilist.c_ilist)
494 is2D = not self.is_3D()
495 if generator:
496 return (read_line(feature_id=v_id, c_mapinfo=self.c_mapinfo,
497 table=self.table, writeable=self.writeable,
498 is2D=is2D)
499 for v_id in ilist)
500 else:
501 return [read_line(feature_id=v_id, c_mapinfo=self.c_mapinfo,
502 table=self.table, writeable=self.writeable,
503 is2D=is2D)
504 for v_id in ilist]
505
506 @must_be_open
507 def read(self, feature_id):
508 """Return a geometry object given the feature id.
509
510 :param int feature_id: the id of feature to obtain
511
512 >>> test_vect = VectorTopo(test_vector_name)
513 >>> test_vect.open(mode='r')
514 >>> feature1 = test_vect.read(0) #doctest: +ELLIPSIS
515 Traceback (most recent call last):
516 ...
517 ValueError: The index must be >0, 0 given.
518 >>> feature1 = test_vect.read(5)
519 >>> feature1
520 Line([Point(12.000000, 4.000000), Point(12.000000, 2.000000), Point(12.000000, 0.000000)])
521 >>> feature1.length()
522 4.0
523 >>> test_vect.read(-1)
524 Centroid(7.500000, 3.500000)
525 >>> len(test_vect)
526 21
527 >>> test_vect.read(21)
528 Centroid(7.500000, 3.500000)
529 >>> test_vect.read(22) #doctest: +ELLIPSIS
530 Traceback (most recent call last):
531 ...
532 IndexError: Index out of range
533 >>> test_vect.close()
534
535 """
536 return read_line(feature_id, self.c_mapinfo, self.table, self.writeable,
537 is2D=not self.is_3D())
538
539 @must_be_open
540 def is_empty(self):
541 """Return if a vector map is empty or not
542 """
543 primitives = self.num_primitives()
544 output = True
545 for v in primitives.values():
546 if v != 0:
547 output = False
548 break
549 return output
550
551 @must_be_open
552 def rewrite(self, geo_obj, cat, attrs=None, **kargs):
553 """Rewrite a geometry features
554
555 >>> cols = [(u'cat', 'INTEGER PRIMARY KEY'),
556 ... (u'name', 'TEXT')]
557
558 Generate a new vector map
559
560 >>> test_vect = VectorTopo('newvect_2')
561 >>> test_vect.open('w', tab_name='newvect_2', tab_cols=cols,
562 ... overwrite=True)
563
564 import a geometry feature ::
565
566 >>> from grass.pygrass.vector.geometry import Point
567
568 create two points ::
569
570 >>> point0 = Point(0, 0)
571 >>> point1 = Point(1, 1)
572 >>> point2 = Point(2, 2)
573
574 then write the two points on the map, with ::
575
576 >>> test_vect.write(point0, cat=1, attrs=('pub',))
577 >>> test_vect.write(point1, cat=2, attrs=('resturant',))
578 >>> test_vect.table.conn.commit() # save changes in the DB
579 >>> test_vect.table_to_dict()
580 {1: [1, 'pub'], 2: [2, 'resturant']}
581 >>> test_vect.close()
582
583 Now rewrite one point of the vector map: ::
584
585 >>> test_vect.open('rw')
586 >>> test_vect.rewrite(point2, cat=1, attrs=('Irish Pub',))
587 >>> test_vect.table.conn.commit() # save changes in the DB
588 >>> test_vect.close()
589
590 Check the output:
591
592 >>> test_vect.open('r')
593 >>> test_vect[1] == point2
594 True
595 >>> test_vect[1].attrs['name'] == 'Irish Pub'
596 True
597 >>> test_vect.close()
598 >>> test_vect.remove()
599 """
600 if self.table is not None and attrs:
601 self.table.update(key=cat, values=attrs)
602 elif self.table is None and attrs:
603 print("Table for vector {name} does not exist, attributes not"
604 " loaded".format(name=self.name))
605 libvect.Vect_cat_set(geo_obj.c_cats, self.layer, cat)
606 result = libvect.Vect_rewrite_line(self.c_mapinfo,
607 cat, geo_obj.gtype,
608 geo_obj.c_points,
609 geo_obj.c_cats)
610 if result == -1:
611 raise GrassError("Not able to write the vector feature.")
612
613 # return offset into file where the feature starts
614 geo_obj.offset = result
615
616 @must_be_open
617 def delete(self, feature_id):
618 """Remove a feature by its id
619
620 :param feature_id: the id of the feature
621 :type feature_id: int
622 """
623 if libvect.Vect_rewrite_line(self.c_mapinfo, feature_id) == -1:
624 raise GrassError("C function: Vect_rewrite_line.")
625
626 @must_be_open
627 def restore(self, geo_obj):
628 if hasattr(geo_obj, 'offset'):
629 if libvect.Vect_restore_line(self.c_mapinfo, geo_obj.offset,
630 geo_obj.id) == -1:
631 raise GrassError("C function: Vect_restore_line.")
632 else:
633 raise ValueError("The value have not an offset attribute.")
634
635 @must_be_open
636 def bbox(self):
637 """Return the BBox of the vecor map
638 """
639 bbox = Bbox()
640 if libvect.Vect_get_map_box(self.c_mapinfo, bbox.c_bbox) == 0:
641 raise GrassError("I can not find the Bbox.")
642 return bbox
643
644 def close(self, build=True, release=True):
645 """Close the VectorTopo map, if release is True, the memory
646 occupied by spatial index is released"""
647 if release:
648 libvect.Vect_set_release_support(self.c_mapinfo)
649 super(VectorTopo, self).close(build=build)
650
651 @must_be_open
652 def table_to_dict(self, where=None):
653 """Return the attribute table as a dictionary with the category as keys
654
655 The columns have the order of the self.table.columns.names() list.
656
657 Examples
658
659 >>> from grass.pygrass.vector import VectorTopo
660 >>> from grass.pygrass.vector.basic import Bbox
661 >>> test_vect = VectorTopo(test_vector_name)
662 >>> test_vect.open('r')
663
664 >>> test_vect.table_to_dict()
665 {1: [1, 'point', 1.0], 2: [2, 'line', 2.0], 3: [3, 'centroid', 3.0]}
666
667 >>> test_vect.table_to_dict(where="value > 2")
668 {3: [3, 'centroid', 3.0]}
669
670 >>> test_vect.table_to_dict(where="value > 0")
671 {1: [1, 'point', 1.0], 2: [2, 'line', 2.0], 3: [3, 'centroid', 3.0]}
672
673 >>> test_vect.table.filters.get_sql()
674 'SELECT cat,name,value FROM vector_doctest_map WHERE value > 0 ORDER BY cat;'
675
676 """
677
678 if self.table is not None:
679 table_dict = {}
680 # Get the category index
681 cat_index = self.table.columns.names().index("cat")
682 # Prepare a filter
683 if where is not None:
684 self.table.filters.where(where)
685
686 self.table.filters.order_by("cat")
687
688 self.table.filters.select(",".join(self.table.columns.names()))
689 # Execute the query and fetch the result
690 cur = self.table.execute()
691 l = cur.fetchall()
692 # Generate the dictionary
693 for entry in l:
694 table_dict[entry[cat_index]] = list(entry)
695
696 return(table_dict)
697
698 return None
699
700 @must_be_open
701 def features_to_wkb_list(self, bbox=None, feature_type="point", field=1):
702 """Return all features of type point, line, boundary or centroid
703 as a list of Well Known Binary representations (WKB)
704 (id, cat, wkb) triplets located in a specific
705 bounding box.
706
707 :param bbox: The boundingbox to search for features,
708 if bbox=None the boundingbox of the whole
709 vector map layer is used
710
711 :type bbox: grass.pygrass.vector.basic.Bbox
712
713 :param feature_type: The type of feature that should be converted to
714 the Well Known Binary (WKB) format. Supported are:
715 'point' -> libvect.GV_POINT 1
716 'line' -> libvect.GV_LINE 2
717 'boundary' -> libvect.GV_BOUNDARY 3
718 'centroid' -> libvect.GV_CENTROID 4
719 :type type: string
720
721 :param field: The category field
722 :type field: integer
723
724 :return: A list of triplets, or None if nothing was found
725
726 The well known binary are stored in byte arrays.
727
728 Examples:
729
730 >>> from grass.pygrass.vector import VectorTopo
731 >>> from grass.pygrass.vector.basic import Bbox
732 >>> test_vect = VectorTopo(test_vector_name)
733 >>> test_vect.open('r')
734
735 >>> bbox = Bbox(north=20, south=-1, east=20, west=-1)
736 >>> result = test_vect.features_to_wkb_list(bbox=bbox,
737 ... feature_type="point")
738 >>> len(result)
739 3
740 >>> for entry in result:
741 ... f_id, cat, wkb = entry
742 ... print((f_id, cat, len(wkb)))
743 (1, 1, 21)
744 (2, 1, 21)
745 (3, 1, 21)
746
747 >>> result = test_vect.features_to_wkb_list(bbox=None,
748 ... feature_type="line")
749 >>> len(result)
750 3
751 >>> for entry in result:
752 ... f_id, cat, wkb = entry
753 ... print((f_id, cat, len(wkb)))
754 (4, 2, 57)
755 (5, 2, 57)
756 (6, 2, 57)
757
758 >>> result = test_vect.features_to_wkb_list(bbox=bbox,
759 ... feature_type="boundary")
760 >>> len(result)
761 11
762
763 >>> result = test_vect.features_to_wkb_list(bbox=None,
764 ... feature_type="centroid")
765 >>> len(result)
766 4
767
768 >>> for entry in result:
769 ... f_id, cat, wkb = entry
770 ... print((f_id, cat, len(wkb)))
771 (19, 3, 21)
772 (18, 3, 21)
773 (20, 3, 21)
774 (21, 3, 21)
775
776 >>> result = test_vect.features_to_wkb_list(bbox=bbox,
777 ... feature_type="blub")
778 Traceback (most recent call last):
779 ...
780 GrassError: Unsupported feature type <blub>, supported are <point,line,boundary,centroid>
781
782 >>> test_vect.close()
783
784 """
785
786 supported = ['point', 'line', 'boundary', 'centroid']
787
788 if feature_type.lower() not in supported:
789 raise GrassError("Unsupported feature type <%s>, "\
790 "supported are <%s>"%(feature_type,
791 ",".join(supported)))
792
793 if bbox is None:
794 bbox = self.bbox()
795
796 bboxlist = self.find_by_bbox.geos(bbox, type=feature_type.lower(),
797 bboxlist_only = True)
798
799 if bboxlist is not None and len(bboxlist) > 0:
800
801 l = []
802 line_p = libvect.line_pnts()
803 line_c = libvect.line_cats()
804 size = ctypes.c_size_t()
805 cat = ctypes.c_int()
806 error = ctypes.c_int()
807
808 for f_id in bboxlist.ids:
809 barray = libvect.Vect_read_line_to_wkb(self.c_mapinfo,
810 ctypes.byref(line_p),
811 ctypes.byref(line_c),
812 f_id,
813 ctypes.byref(size),
814 ctypes.byref(error))
815 if not barray:
816 if error == -1:
817 raise GrassError(_("Unable to read line of feature %i"%(f_id)))
818 if error == -2:
819 print("Empty feature %i"%(f_id))
820 continue
821
822 ok = libvect.Vect_cat_get(ctypes.byref(line_c), field,
823 ctypes.byref(cat))
824 if ok < 1:
825 pcat = None
826 else:
827 pcat = cat.value
828
829 l.append((f_id, pcat, ctypes.string_at(barray, size.value)))
830 libgis.G_free(barray)
831
832 return l
833 return None
834
835 @must_be_open
836 def areas_to_wkb_list(self, bbox=None, field=1):
837 """Return all features of type point, line, boundary or centroid
838 as a list of Well Known Binary representations (WKB)
839 (id, cat, wkb) triplets located in a specific
840 bounding box.
841
842 :param bbox: The boundingbox to search for features,
843 if bbox=None the boundingbox of the whole
844 vector map layer is used
845
846 :type bbox: grass.pygrass.vector.basic.Bbox
847
848 :param field: The centroid category field
849 :type field: integer
850
851 :return: A list of triplets, or None if nothing was found
852
853 The well known binary are stored in byte arrays.
854
855 Examples:
856
857 >>> from grass.pygrass.vector import VectorTopo
858 >>> from grass.pygrass.vector.basic import Bbox
859 >>> test_vect = VectorTopo(test_vector_name)
860 >>> test_vect.open('r')
861
862 >>> bbox = Bbox(north=20, south=-1, east=20, west=-1)
863 >>> result = test_vect.areas_to_wkb_list(bbox=bbox)
864 >>> len(result)
865 4
866 >>> for entry in result:
867 ... a_id, cat, wkb = entry
868 ... print((a_id, cat, len(wkb)))
869 (1, 3, 225)
870 (2, 3, 141)
871 (3, 3, 93)
872 (4, 3, 141)
873
874 >>> result = test_vect.areas_to_wkb_list()
875 >>> len(result)
876 4
877 >>> for entry in result:
878 ... a_id, cat, wkb = entry
879 ... print((a_id, cat, len(wkb)))
880 (1, 3, 225)
881 (2, 3, 141)
882 (3, 3, 93)
883 (4, 3, 141)
884
885 >>> test_vect.close()
886
887
888 """
889 if bbox is None:
890 bbox = self.bbox()
891
892 bboxlist = self.find_by_bbox.areas(bbox, bboxlist_only = True)
893
894 if bboxlist is not None and len(bboxlist) > 0:
895
896 l = []
897 line_c = libvect.line_cats()
898 size = ctypes.c_size_t()
899 cat = ctypes.c_int()
900
901 for a_id in bboxlist.ids:
902 barray = libvect.Vect_read_area_to_wkb(self.c_mapinfo,
903 a_id,
904 ctypes.byref(size))
905 if not barray:
906 raise GrassError(_("Unable to read area with id %i"%(a_id)))
907
908 pcat = None
909 c_ok = libvect.Vect_get_area_cats(self.c_mapinfo, a_id,
910 ctypes.byref(line_c))
911 if c_ok == 0: # Centroid found
912
913 ok = libvect.Vect_cat_get(ctypes.byref(line_c), field,
914 ctypes.byref(cat))
915 if ok > 0:
916 pcat = cat.value
917
918 l.append((a_id, pcat, ctypes.string_at(barray, size.value)))
919 libgis.G_free(barray)
920
921 return l
922 return None
923
924if __name__ == "__main__":
925 import doctest
926 from grass.pygrass import utils
927 utils.create_test_vector_map(test_vector_name)
928 doctest.testmod()
929
930 """Remove the generated vector map, if exist"""
931 from grass.pygrass.utils import get_mapset_vector
932 from grass.script.core import run_command
933 mset = get_mapset_vector(test_vector_name, mapset='')
934 if mset:
935 run_command("g.remove", flags='f', type='vector', name=test_vector_name)
Note: See TracBrowser for help on using the repository browser.