13

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?

3
  • the thing is I have a raster with a bounding box. I don't want to create a bounding box using coordinates 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:)! Feb 29 '20 at 14:17
15

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:

enter image description here

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
1
  • gdaltindex is nice!
    – dorbodwolf
    Oct 20 '20 at 7:29
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': {}})
1
  • 1
    GeoPandas use Fiona
    – gene
    Feb 29 '20 at 17:20
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.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.