Příklad vytvoření obalové zóny

Na tomto příkladu si ukážeme, jak otevřít vektorová data ve formátu Esri Shapefile, načíst datovou vrstvu, zobrazit atributy geoprvků a vytvořit obalovou zónu nad načtenými geoprvky.

Nejprve otevření souboru s daty:

from osgeo import ogr
ds = ogr.Open("data/chko.shp")
print(ds)

<osgeo.ogr.DataSource; proxy of <Swig Object of type 'OGRDataSourceShadow *' at 0x7f98d8152a50> >

print(ds.GetLayerCount())

1

Práce s vrstvou, její otevření:

l = ds.GetLayer(0)
print(l)

<osgeo.ogr.Layer; proxy of <Swig Object of type 'OGRLayerShadow *' at 0x7f98d80fa870> >

print(l.GetFeatureCount())

5626

Schéma vrstvy - definice typu geometrie a jednotlivých atributových polí:

print(l.GetGeomType())
print(l.GetGeomType() == ogr.wkbPolygon)

3
True

print(l.schema)
print(l.schema[4].name)

[<osgeo.ogr.FieldDefn; proxy of <Swig Object of type 'OGRFieldDefnShadow *' at 0x7f98d80fa9f0> >,
<osgeo.ogr.FieldDefn; proxy of <Swig Object of type 'OGRFieldDefnShadow *' at 0x7f98d80fa8...
'NAZEV'

Vypsání názvu geoprvku (atribut NAZEV):

features_nr = l.GetFeatureCount()
for i in range(features_nr):
    f = l.GetNextFeature()
    print(f.GetField('NAZEV'))

Český ráj
...

Vypsání vlastnosti geometrické složky popisu geoprvků (minimálního ohraničujícího obdélíku a centroidu polygonu, vytvoření obálky):

f = l.GetFeature(54)
f.GetField('NAZEV')
print (f.GetField('NAZEV'))

Český ráj

geom = f.GetGeometryRef()
print(geom.GetEnvelope())

(-683329.1875, -681265.625, -993228.75, -991528.0)

c = geom.Centroid()
print(c.GetPoint())

(-682407.4126500859, -992433.3498782327, 0.0)

buff = c.Buffer(100)
print(geom.Intersects(buff))

True