How to use rasterio - 10 common examples

source2[rows // 2:, cols // 2:] = 100

    kwargs = {
        "count": 1,
        "width": params.width,
        "height": params.height,
        "dtype": np.uint8,
        "driver": "GTiff",
        "crs": params.dst_crs,
        "transform": params.dst_transform,
    }

    with rasterio.open(tiffname, "w", **kwargs) as dst:
        reproject(
            source1,
            rasterio.band(dst, 1),
            src_transform=params.src_transform,
            src_crs=params.src_crs,
            src_nodata=0.0,
            dst_transform=params.dst_transform,
            dst_crs=params.dst_crs,
            dst_nodata=0.0,
        )

        reproject(
            source2,
            rasterio.band(dst, 1),
            src_transform=params.src_transform,
            src_crs=params.src_crs,
            src_nodata=0.0,
            dst_transform=params.dst_transform,
            dst_crs=params.dst_crs,

window = evaluate(window, height=height, width=width)

    height_shape = block_shapes[0][0]
    width_shape = block_shapes[0][1]

    (row_start, row_stop), (col_start, col_stop) = window.toranges()

    row_min = int(row_start // height_shape) * height_shape
    row_max = int(row_stop // height_shape) * height_shape + \
        (height_shape if row_stop % height_shape != 0 else 0)

    col_min = int(col_start // width_shape) * width_shape
    col_max = int(col_stop // width_shape) * width_shape + \
        (width_shape if col_stop % width_shape != 0 else 0)

    return Window(col_min, row_min, col_max - col_min, row_max - row_min)

"""
    Create Cloud Optimized Geotiff.
    Parameters
    ----------
    src_path : str or PathLike object
        A dataset path or URL. Will be opened in "r" mode.
    This script is the rasterio equivalent of
    https://svn.osgeo.org/gdal/trunk/gdal/swig/python/samples/validate_cloud_optimized_geotiff.py
    """
    errors = []
    details = {}

    if not GDALVersion.runtime().at_least("2.2"):
        raise Exception("GDAL 2.2 or above required")

    with rasterio.open(src_path) as src:
        if not src.driver == "GTiff":
            raise Exception("The file is not a GeoTIFF")

        filelist = [os.path.basename(f) for f in src.files]
        src_bname = os.path.basename(src_path)
        if len(filelist) > 1 and src_bname + ".ovr" in filelist:
            errors.append(
                "Overviews found in external .ovr file. They should be internal"
            )

        if src.width >= 512 or src.height >= 512:
            if not src.is_tiled:
                errors.append(
                    "The file is greater than 512xH or 512xW, but is not tiled"
                )

def test_rd_internals_crs():
    from rasterio.crs import CRS as RioCRS

    assert _dc_crs(None) is None
    assert _dc_crs(RioCRS()) is None
    assert _dc_crs(RioCRS.from_epsg(3857)).epsg == 3857
    assert _dc_crs(RioCRS.from_wkt(SAMPLE_WKT_WITHOUT_AUTHORITY)).epsg is None

def test_reproject_invalid_src_nodata():
    """src_nodata must be in range for data type."""
    params = default_reproject_params()

    source = np.ones((params.width, params.height), dtype=np.uint8)
    out = source.copy()

    with pytest.raises(ValueError):
        reproject(
            source,
            out,
            src_transform=params.src_transform,
            src_crs=params.src_crs,
            src_nodata=999999999,
            dst_transform=params.dst_transform,
            dst_crs=params.dst_crs,
            dst_nodata=215,
        )

def test_reproject_epsg():
    with rasterio.open('tests/data/RGB.byte.tif') as src:
        source = src.read(1)

    dst_crs = {'init': 'EPSG:3857'}
    out = np.empty(src.shape, dtype=np.uint8)
    reproject(
        source,
        out,
        src_transform=src.transform,
        src_crs=src.crs,
        dst_transform=DST_TRANSFORM,
        dst_crs=dst_crs,
        resampling=Resampling.nearest)
    assert (out > 0).sum() == 438113

profile = src.profile.copy()

        dst_crs = {'init': 'EPSG:32619'}

        # Calculate the ideal dimensions and transformation in the new crs
        dst_affine, dst_width, dst_height = calculate_default_transform(
            src.crs, dst_crs, src.width, src.height, *src.bounds)

        profile['height'] = dst_height
        profile['width'] = dst_width

        out = np.empty(shape=(dst_height, dst_width), dtype=np.uint8)

        # see #614, some resampling methods succeed but produce blank images
        out = np.empty(src.shape, dtype=np.uint8)
        reproject(
            source,
            out,
            src_transform=src.transform,
            src_crs=src.crs,
            dst_transform=dst_affine,
            dst_crs=dst_crs,
            resampling=method)

        assert out.mean() > 0

def test_reproject_init_dest_nodata():
    """No pixels should transfer over"""
    crs = CRS.from_epsg(4326)
    transform = Affine.identity()
    source = np.zeros((1, 100, 100))
    destination = np.ones((1, 100, 100))
    reproject(
        source, destination, src_crs=crs, src_transform=transform,
        dst_crs=crs, dst_transform=transform,
        src_nodata=0, init_dest_nodata=False
    )
    assert destination.all()

'transform': params.dst_transform
    }

    with rasterio.open(tiffname, 'w', **kwargs) as dst:
        reproject(
            source1,
            rasterio.band(dst, 1),
            src_transform=params.src_transform,
            src_crs=params.src_crs,
            src_nodata=0.0,
            dst_transform=params.dst_transform,
            dst_crs=params.dst_crs,
            dst_nodata=0.0
        )

        reproject(
            source2,
            rasterio.band(dst, 1),
            src_transform=params.src_transform,
            src_crs=params.src_crs,
            src_nodata=0.0,
            dst_transform=params.dst_transform,
            dst_crs=params.dst_crs,
            dst_nodata=0.0,
            init_dest_nodata=False
        )

    # 200s should remain along with 100s
    with rasterio.open(tiffname) as src:
        data = src.read()

    assert data.max() == 200

def test_threads_main_env():
    """Get raster data using ThreadPoolExecutor with main thread Env"""
    with rasterio.Env(), ThreadPoolExecutor(4) as pool:
        for res in pool.map(get_data, ['tests/data/RGB.byte.tif'] * 10):
            assert res.any()