I have a raster that I opened with rasterio. When I do raster.bounds I get
BoundingBox(left=399960.0, bottom=-3309780.0, right=609780.0, top=-3099960.0).
My question is now, if there is any simple way to make a shapefile of it?
Asked
Active
Viewed 1.2k times
16
Vince
- 20,017
- 15
- 45
- 64
Robin Kohrs
- 675
- 1
- 8
- 21
-
the thing is I have a raster with a bounding box. I don't want to create a bounding box using coordinates – Robin Kohrs Feb 29 '20 at 13:34
-
Please don't comment on your own question. Instead, [Edit] the Question to contain this information. One in mind that raster boundaries are sometimes computed from center-center of the corner pixels, so your polygon may be a half-pixel too small in each corner. – Vince Feb 29 '20 at 13:58
-
thank you very much for the advices:)! – Robin Kohrs Feb 29 '20 at 14:17
4 Answers
26
With shapely box and GeoPandas
import rasterio as rio
ra = rio.open("raster.tif")
bounds = ra.bounds
Convert bounds to shapely geometry
from shapely.geometry import box
geom = box(*bounds)
print(geom.wkt)
POLYGON ((205625.4144070312 88411.04799999669, 205625.4144070312 90534.35044407193, 202086.5770002506 90534.35044407193, 202086.5770002506 88411.04799999669, 205625.4144070312 88411.04799999669))
Create a shapefile with GeoPandas
df = gpd.GeoDataFrame({"id":1,"geometry":[geom]})
df.to_file("boundary.shp")
In one line:
df = gpd.GeoDataFrame({"id":1,"geometry":[box(*bounds)]})
Result:
gene
- 54,868
- 3
- 110
- 187
6
For scripting, other answers does the job. If you need generating shp from command line, I would recommand to use rasterio command line utility rio with
rio bounds /tmp/RGBA.byte.tif > tmp.geojson
ogr2ogr RGBA.shp tmp.geojson
rm tmp.geojson
You could also use gdaltindex e.g
gdaltindex RGBA.shp RGBA.byte.tif
ThomasG77
- 30,725
- 1
- 53
- 93
5
Here is a solution using shapely and fiona that takes into account the spatial reference of the raster and makes sure the new shapefile has it:
# imports
from shapely.geometry import box, mapping
import fiona
import rasterio
# create a Polygon from the raster bounds
bbox = box(*raster.bounds)
# create a schema with no properties
schema = {'geometry': 'Polygon', 'properties': {}}
# create shapefile
with fiona.open('bbox.shp', 'w', driver='ESRI Shapefile',
crs=raster.crs.to_dict(), schema=schema) as c:
c.write({'geometry': mapping(bbox), 'properties': {}})
Marcelo Villa
- 5,928
- 2
- 19
- 38
3
There is also a GDAL program gdaltindex https://gdal.org/programs/gdaltindex.html for this purpose.
Description
This program builds a shapefile with a record for each input raster file, an attribute containing the filename, and a polygon geometry outlining the raster.
user30184
- 65,331
- 4
- 65
- 118
