Opened 6 years ago

Closed 6 years ago

#6858 closed defect (invalid)

Shapefile extent not updated when a geometry is updated

Reported by: Mike Taves Owned by: warmerdam
Priority: normal Milestone:
Component: default Version: svn-trunk
Severity: normal Keywords:
Cc:

Description

I can't seem to update a Shapefile's geometry and get a correct extent for a layer. For example:

from osgeo import gdal, ogr
print(gdal.VersionInfo())
ogr.UseExceptions()

# Create a new Shapefile with one simple geometry
shp_drv = ogr.GetDriverByName('Esri Shapefile')
shp_fname = 'tmp.shp'
ds = shp_drv.CreateDataSource(shp_fname)

ly = ds.CreateLayer('', None, ogr.wkbPolygon)
fd = ly.GetLayerDefn()
ft = ogr.Feature(fd)
g = ogr.CreateGeometryFromWkt('POLYGON((0 0, 0 1, 1 0, 0 0))')
assert ft.SetGeometry(g) == 0
assert ly.CreateFeature(ft) == 0
ft = g = None
ds = ly = fd = None

# Confirm the previous extent, and edit the geometry
ds = ogr.Open(shp_fname, 1)
ly = ds.GetLayer()
print(ly.GetExtent())  # (0.0, 1.0, 0.0, 1.0))
ft = ly.GetFeature(0)
g = ogr.CreateGeometryFromWkt('POLYGON((100 100, 100 101, 101 100, 100 100))')
assert ft.SetGeometry(g) == 0
assert ly.SetFeature(ft) == 0
ft = g = None

# This is wrong
print(ly.GetExtent(force=1))  # (0.0, 101.0, 0.0, 101.0)
ly = ds = None

# Open again, and confirm things
ds = ogr.Open(shp_fname, 1)
ly = ds.GetLayer()
assert ly.GetFeatureCount() == 1
ft = ly.GetFeature(0)
g = ft.GetGeometryRef()
print(g.GetEnvelope())  # (100.0, 101.0, 100.0, 101.0)
print(ly.GetExtent(force=1))  # (0.0, 101.0, 0.0, 101.0)
ly = ds = ft = g = None

As you can see, the geometry envelope was correctly updated, but the layer extent was not.

I've tried a similar exercise with the memory driver (but keeping the dataset open), but I see ly.GetExtent(force=1) return the correct extent (100.0, 101.0, 100.0, 101.0). Therefore I might assume this bug is limited to the Shapefile driver.

I've tested this with version '2010300' on Windows 10 x64 via Anaconda2, and CentOS6 with Python 2.6 compiled from source. Same result with '1100100' on Ubuntu with python-gdal, so this behavior not new.

Change History (3)

comment:1 by Mike Taves, 6 years ago

I should have noted the extent was updated to be larger (thus it's a misleading bug title), but the minimum X and Y extents were not contracted.

Reading the docs for OGRLayer::GetExtent:

If bForce is TRUE then some implementations will actually scan the entire layer once to compute the MBR of all the features in the layer.

It seems that the extent may have been a union of the previous extent and the new extent(s), rather than an extent based on a single layer scan.

comment:2 by Jukka Rahkonen, 6 years ago

The source code of shape driver is in https://trac.osgeo.org/gdal/browser/trunk/gdal/ogr/ogrsf_frmts/shape/ogrshapelayer.cpp.

There seems to be a comment:

// /* GetExtent() */ /* */ /* Fetch extent of the data currently stored in the dataset. */ /* The bForce flag has no effect on SHP files since that value */ /* is always in the header. */ /* */ /* Returns OGRERR_NONE/OGRRERR_FAILURE. */

I can also see this:

// /* RecomputeExtent() */ /* */ /* Force recomputation of the extent of the .SHP file */ //

While I can read the comments I can't really read the code but perhaps this is something you should do.

comment:3 by Even Rouault, 6 years ago

Resolution: invalid
Status: newclosed

Quoting http://gdal.org/drv_shapefile.html :

Spatial extent

Shapefiles store the layer spatial extent in the .SHP file. The layer spatial extent is automatically updated when inserting a new feature in a shapefile. However when updating an existing feature, if its previous shape was touching the bounding box of the layer extent but the updated shape does not touch the new extent, the computed extent will not be correct. It is then necessary to force a recomputation by invoking the SQL 'RECOMPUTE EXTENT ON <tablename>' via the datasource ExecuteSQL() method. The same applies for the deletion of a shape.
Note: See TracTickets for help on using tickets.