TransWikia.com

Changing coordinate system with rasterio affects dataset opening

Geographic Information Systems Asked by SaskiaG on January 1, 2021

I have troubles changing coordinate system of a raster with rasterio. Here is my code:

  1. Here I create my .tif (it’s a RGB image)
import numpy as np
import rasterio
from rasterio.warp import calculate_default_transform, reproject, Resampling


    with rasterio.open(result_path + "" + "test_infra_{}.tif".format(product_group), 'w', **meta) as dst:
        dst.write(bandstack_infra)
        dst.close()

  1. Then I re-open it and change the meta and create a new .tif with the new coordinate system, according to the piece of code given in the docs (https://rasterio.readthedocs.io/en/latest/topics/reproject.html).
    dst_crs = 'EPSG:4326'
    with rasterio.open(result_path + "" + "test_infra_{}.tif".format(product_group)) as src:
        transform, width, height = calculate_default_transform(src.crs, dst_crs, src.width, src.height, *src.bounds)
        kwargs = src.meta.copy()
        kwargs.update({'crs': dst_crs, 'transform': transform, 'width': width, 'height': height})

    with rasterio.open(result_path + "" + "test_infra_{}_wgs.tif".format(product_group), 'w', **kwargs) as dst:
        for i in range(1, src.count + 1):
            reproject(
                source=rasterio.band(src, i),
                destination=rasterio.band(dst, i),
                src_transform=src.transform,
                src_crs=src.crs,
                dst_transform=transform,
                dst_crs=dst_crs,
                resampling=Resampling.nearest)

Then I get the error:

rasterio.errors.RasterioIOError: Dataset is closed: --> path to the original data, open as src.

I understand that the new dataset cannot be created if the original one (open as src) is not open. However, I didn’t close it (with src.close())!

Why is it closed and how can I solve this problem?

One Answer

In Python the with statement ensures that the file will be closed when leaving the processing therefore your raster is closed ( With statement in Python)

Therefore try:

dst_crs = 'EPSG:4326'
with rasterio.open(result_path + "" + "test_infra_{}.tif".format(product_group)) as src:
    transform, width, height = calculate_default_transform(src.crs, dst_crs, src.width, src.height, *src.bounds)
    kwargs = src.meta.copy()
    kwargs.update({'crs': dst_crs, 'transform': transform, 'width': width, 'height': height})
    with rasterio.open(result_path + "" + "test_infra_{}_wgs.tif".format(product_group), 'w', **kwargs) as dst:
        for i in range(1, src.count + 1):
            reproject(
                source=rasterio.band(src, i),
                destination=rasterio.band(dst, i),
                src_transform=src.transform,
                src_crs=src.crs,
                dst_transform=transform,
                dst_crs=dst_crs,
                resampling=Resampling.nearest)

and the src raster will not be closed in the second with

Answered by gene on January 1, 2021

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP