How to use the rasterio.Env function in rasterio

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()

def test_reproject_nodata_nan(options, expected):

    with rasterio.Env(**options):
        params = uninvertable_reproject_params()
        source = np.ones((params.width, params.height), dtype=np.float32)
        out = np.zeros((params.dst_width, params.dst_height),
                       dtype=source.dtype)
        out.fill(120)  # Fill with arbitrary value

        reproject(
            source,
            out,
            src_transform=params.src_transform,
            src_crs=params.src_crs,
            src_nodata=np.nan,
            dst_transform=params.dst_transform,
            dst_crs=params.dst_crs,
            dst_nodata=np.nan
        )

def test_reproject_nodata(options, expected):
    # Older combinations of GDAL and PROJ might have got this transformation wrong.
    # Results look better with GDAL 3.
    nodata = 215

    with rasterio.Env(**options):
        params = uninvertable_reproject_params()
        source = np.ones((params.width, params.height), dtype=np.uint8)
        out = np.zeros((params.dst_width, params.dst_height), dtype=source.dtype)
        out.fill(120)  # Fill with arbitrary value

        reproject(
            source,
            out,
            src_transform=params.src_transform,
            src_crs=params.src_crs,
            src_nodata=nodata,
            dst_transform=params.dst_transform,
            dst_crs=params.dst_crs,
            dst_nodata=nodata,
        )

def test_reproject_nodata(options, expected):
    nodata = 215

    with rasterio.Env(**options):
        params = uninvertable_reproject_params()
        source = np.ones((params.width, params.height), dtype=np.uint8)
        out = np.zeros((params.dst_width, params.dst_height),
                       dtype=source.dtype)
        out.fill(120)  # Fill with arbitrary value

        reproject(
            source,
            out,
            src_transform=params.src_transform,
            src_crs=params.src_crs,
            src_nodata=nodata,
            dst_transform=params.dst_transform,
            dst_crs=params.dst_crs,
            dst_nodata=nodata
        )

template_mask = template_variable[:].mask

    template_y_name, template_x_name = template_variable.dimensions[-2:]
    template_coords = SpatialCoordinateVariables.from_dataset(
        template_ds,
        x_name=template_x_name,
        y_name=template_y_name,
        projection=template_prj
    )
    # template_geo_bbox = template_coords.bbox.project(ds_prj, edge_points=21)  # TODO: add when needing to subset

    ds_y_name, ds_x_name = ds.variables[variables[0]].dimensions[-2:]
    proj = Proj(init=ds_projection) if 'EPSG:' in ds_projection.upper() else Proj(str(ds_projection))
    ds_coords = SpatialCoordinateVariables.from_dataset(ds, x_name=ds_x_name, y_name=ds_y_name, projection=proj)

    with rasterio.Env():
        # Copy dimensions for variable across to output
        for dim_name in template_variable.dimensions:
            if not dim_name in out_ds.dimensions:
                if dim_name in template_ds.variables and not dim_name in out_ds.variables:
                    copy_variable(template_ds, out_ds, dim_name)
                else:
                    copy_dimension(template_ds, out_ds, dim_name)

        for variable_name in variables:
            click.echo('Processing: {0}'.format(variable_name))

            variable = ds.variables[variable_name]
            fill_value = getattr(variable, '_FillValue', variable[0, 0].fill_value)

            for dim_name in variable.dimensions[:-2]:
                if not dim_name in out_ds.dimensions:

# insert newline for nicer progress bar style
        click.echo('')

    sub_pbar_args = dict(
        disable=quiet,
        leave=False,
        bar_format='{l_bar}{bar}| {n_fmt}/{total_fmt}'
    )

    with contextlib.ExitStack() as outer_env:
        pbar = outer_env.enter_context(tqdm.tqdm(
            total=total_pixels, smoothing=0, disable=quiet,
            bar_format='{l_bar}{bar}| [{elapsed}<{remaining}{postfix}]',
            desc='Optimizing rasters'
        ))
        outer_env.enter_context(rasterio.Env(**GDAL_CONFIG))

        for input_file in raster_files_flat:
            if len(input_file.name) > 30:
                short_name = input_file.name[:13] + '...' + input_file.name[-13:]
            else:
                short_name = input_file.name

            pbar.set_postfix(file=short_name)

            output_file = output_folder / input_file.with_suffix('.tif').name

            if not overwrite and output_file.is_file():
                raise click.BadParameter(
                    f'Output file {output_file!s} exists (use --overwrite to ignore)'
                )

src,
        src_nodata=src_nodata,
        crs=bounds.crs,
        width=vrt_width,
        height=vrt_height,
        transform=vrt_transform,
        resampling=resampling,
        add_alpha=add_alpha,
    ) as vrt:
        dst_window = vrt.window(*bounds.bounds)

        data = vrt.read(out_shape=(vrt.count,) + target_shape, window=dst_window)

        mask = np.ma.nomask
        if source.mask:
            with rasterio.Env(OGR_ENABLE_PARTIAL_REPROJECTION=True):
                geom_mask = transform_geom(WGS84_CRS, bounds.crs, source.mask)

            mask_transform = from_bounds(
                *bounds.bounds, height=target_shape[0], width=target_shape[1]
            )
            mask = geometry_mask(
                [geom_mask], target_shape, transform=mask_transform, invert=True
            )

        if any([ColorInterp.alpha in vrt.colorinterp]):
            alpha_idx = vrt.colorinterp.index(ColorInterp.alpha)
            # TODO this coerces the alpha channel to a boolean mask; it should
            # eventually be used *as alpha*, setting transparency
            mask = [~data[alpha_idx].astype(np.bool) | mask] * (vrt.count - 1)
            bands = [data[i] for i in range(0, vrt.count) if i != alpha_idx]
            data = np.ma.masked_array(bands, mask=mask)

target_tile_height = tile_height = (256 + buffers[1] + buffers[3]) * scale

    window_width = (window[1][1] - window[1][0])
    window_height = (window[0][1] - window[0][0])

    # compare window to tile width to see how much we're up/downsampling
    scale_factor = math.ceil(256 / (window_width - buffers[0] - buffers[2]))

    if scale_factor > 1:
        # respect the window size so it can be interpolated cleanly
        tile_width = int(min(window_width, target_tile_width))
        tile_height = int(min(window_height, target_tile_height))

    buffers = map(lambda x: int(x * scale), buffers)

    with rasterio.Env(CPL_VSIL_CURL_ALLOWED_EXTENSIONS='.vrt,.tif,.ovr,.msk'):
        src = get_source(src_url)

        # TODO read the data and the mask in parallel
        if mask_url:
            data = src.read(
                out_shape=(src.count, tile_height, tile_width),
                window=window,
            )

            # handle masking ourselves (src.read(masked=True) doesn't use
            # overviews)
            mask = None
            if src.nodata is not None:
                # some datasets use the min value but report an alternate nodata value
                nodata_alt = None
                if np.issubdtype(data.dtype, int):

def main(infile, outfile, num_workers=4):
    """Process infile block-by-block and write to a new file

    The output is the same as the input, but with band order
    reversed.
    """

    with rasterio.Env():

        with rasterio.open(infile) as src:

            # Create a destination dataset based on source params. The
            # destination will be tiled, and we'll process the tiles
            # concurrently.
            profile = src.profile
            profile.update(blockxsize=128, blockysize=128, tiled=True)

            with rasterio.open(outfile, "w", **profile) as dst:

                # Materialize a list of destination block windows
                # that we will use in several statements below.
                windows = [window for ij, window in dst.block_windows()]

                # This generator comprehension gives us raster data

prj must be a pyproj.Proj object
    """

    if format != 'GTiff':
        raise NotImplementedError('Formats besides GTiff not yet supported')

    meta = {
        'width': data.shape[1],
        'height': data.shape[0],
        'dtype': data.dtype.name, # rasterio uses strings, not dtypes
        'transform': affine,
        'count': 1,
        'crs': projection.srs
    }

    with rasterio.Env():
        with rasterio.open(outfilename, 'w', driver=format, **meta) as out:
            out.write(data, 1)