Write chunk array as GeoTIFF¶

Or: How to safe a single layer chunk array as georeferenced raster with MultiRasterIO?

[4]:

import pandas as pd
from pathlib import Path
import rasterio

from eobox.raster import MultiRasterIO
from eobox import sampledata


year = 2008
dataset = sampledata.get_dataset("lsts")
layers_paths = [Path(p) for p in dataset["raster_files"]]
layers_df = pd.Series([p.stem for p in layers_paths]).str.split("_", expand=True) \
.rename({0: "sceneid", 1:"band"}, axis=1)

layers_df["date"] = pd.to_datetime(layers_df.sceneid.str[9:16], format="%Y%j")
layers_df["uname"] = layers_df.sceneid.str[:3] + "_" + layers_df.date.dt.strftime("%Y-%m-%d") + "_" + layers_df.band.str[::]
layers_df["path"] = layers_paths

layers_df = layers_df.sort_values(["date", "band"])
layers_df = layers_df.reset_index(drop=True)

df_layers = layers_df[(layers_df.date >= str(year)) & (layers_df.date < str(year+1))]
df_layers = df_layers.reset_index(drop=True)
df_layers.head()

[4]:

	sceneid	band	date	uname	path
0	LT50350322008110PAC01	b3	2008-04-19	LT5_2008-04-19_b3	/home/ben/anaconda/envs/eocube/lib/python3.6/s...
1	LT50350322008110PAC01	b4	2008-04-19	LT5_2008-04-19_b4	/home/ben/anaconda/envs/eocube/lib/python3.6/s...
2	LT50350322008110PAC01	b5	2008-04-19	LT5_2008-04-19_b5	/home/ben/anaconda/envs/eocube/lib/python3.6/s...
3	LT50350322008110PAC01	fmask	2008-04-19	LT5_2008-04-19_fmask	/home/ben/anaconda/envs/eocube/lib/python3.6/s...
4	LE70350322008118EDC00	b3	2008-04-27	LE7_2008-04-27_b3	/home/ben/anaconda/envs/eocube/lib/python3.6/s...

[11]:

mrio = MultiRasterIO(df_layers.path.values) \
    .windows_from_blocksize(2**5)
n_chunks = len(mrio.windows)
print("Number of chunks : ", n_chunks)

write_to_disk = False

with rasterio.open(mrio._get_template_for_given_resolution(mrio.dst_res, "path")) as src_layer:
    pass # later we need src_layer for src_layer.window_transform(win)

for ji in range(n_chunks):
    dst_path = f"./xxx_uncontrolled_99_chunk_ji{ji:02}.tif"
    chunk_arrays_ji = mrio.get_arrays(ji)
    chunk_layer_ji = chunk_arrays_ji[:,:,[0]]
    print("shape of chunk_layer_ji : ", chunk_layer_ji.shape)
    win = mrio.windows[ji]
    kwargs = mrio._get_template_for_given_resolution(
        res=mrio.dst_res, return_="meta").copy()
    kwargs.update({"height": win.height,
                   "width": win.width,
                   "transform": src_layer.window_transform(win)})
    kwargs["dtype"] = chunk_layer_ji.dtype
    with rasterio.open(dst_path, "w", **kwargs) as dst:
        dst.write(chunk_layer_ji[:,:,0], 1)
        print(f"Chunk written to path: {dst_path}")

Number of chunks :  4
shape of chunk_layer_ji :  (32, 32, 1)
Chunk written to path: ./xxx_uncontrolled_99_chunk_ji00.tif
shape of chunk_layer_ji :  (32, 29, 1)
Chunk written to path: ./xxx_uncontrolled_99_chunk_ji01.tif
shape of chunk_layer_ji :  (29, 32, 1)
Chunk written to path: ./xxx_uncontrolled_99_chunk_ji02.tif
shape of chunk_layer_ji :  (29, 29, 1)
Chunk written to path: ./xxx_uncontrolled_99_chunk_ji03.tif