Topics

Combining rasters using rasterio


James David Smith
 

Hello,

Thanks for the help last week with the OVR files. I think I've solved that now.

My next question is about combining rasters of shared space, but different extents. Each raster has one band. I want to combine them into one raster, keeping the values of the highest cell value where they overlap. Hopefully the attached image explains better than I can in writing. If anyone can point me in the right direction or provide some examples I'd be very grateful. A previous colleague completed this work using ArcPY, and I'm doing my best to re-write the code in open-source. Thanks!

Capture.PNG


vincent.sarago@...
 

On Wed, Jul 3, 2019 at 06:52 AM, James David Smith wrote:
g. If anyone can point me in the right direction or provide some examples I'd be very grateful. A previous colleague completed this work using ArcPY, and I'm doing my best to re-write the code in open-source. Thanks!
Hi James,
Rasterio has a cli command called `merge` (`rio merge`) which is the equivalent of `gdal_merge.py` with rasterio. If you are looking for a way to do it within python you can check https://github.com/mapbox/rasterio/blob/master/rasterio/rio/merge.py

Vincent 


James David Smith
 

Hi Vincent,

Thanks for the reply. I looked at merge (
https://rasterio.readthedocs.io/en/stable/api/rasterio.merge.html ),
but was a bit confused about what to do for the bounds. The
documentation says:

"Geospatial bounds and resolution of a new output file in the units of
the input file coordinate reference system may be provided and are
otherwise taken from the first input file."

So in my example, the bounds would be set from raster_one . What I
actually need to do is set the bounds as the min_x of all 3 rasters,
the min_y of all 3 rasters, the max_x of all 3 rasters, and the max_y
of all 3 rasters.

I couldn't figure out how to do that.

Thanks, James

On Fri, 5 Jul 2019 at 17:06, <vincent.sarago@gmail.com> wrote:

On Wed, Jul 3, 2019 at 06:52 AM, James David Smith wrote:

g. If anyone can point me in the right direction or provide some examples I'd be very grateful. A previous colleague completed this work using ArcPY, and I'm doing my best to re-write the code in open-source. Thanks!

Hi James,
Rasterio has a cli command called `merge` (`rio merge`) which is the equivalent of `gdal_merge.py` with rasterio. If you are looking for a way to do it within python you can check https://github.com/mapbox/rasterio/blob/master/rasterio/rio/merge.py

Vincent


vincent.sarago@...
 

It seems the docs it wrong, looking at the code it take the max extent of all the dataset https://github.com/mapbox/rasterio/blob/b9f34ee559039239c7c0c97bd911b466701a39cd/rasterio/merge.py#L77-L92


James David Smith
 

Oh right. Well that sounds promising then. So for a set of single band
geotiffs, would I be doing something like this .... ?

file_one = rasterio.open('raster_one.tif', 'r')
file_two = rasterio.open('raster_two.tif', 'r')
file_three = rasterio.open('raster_three.tif', 'r')

new_raster = rasterio.merge.merge([file_one, file_two, file_three], indexes=1)

Thanks again, James

On Fri, 5 Jul 2019 at 17:33, <vincent.sarago@gmail.com> wrote:

It seems the docs it wrong, looking at the code it take the max extent of all the dataset https://github.com/mapbox/rasterio/blob/b9f34ee559039239c7c0c97bd911b466701a39cd/rasterio/merge.py#L77-L92