{ "cells": [ { "cell_type": "markdown", "id": "dc6a17b0-c9a5-4e6c-885f-4950e41df86c", "metadata": {}, "source": [ "# Common raster operations" ] }, { "cell_type": "markdown", "id": "0103b09a-4c6f-4b39-9b45-1ba790e0085e", "metadata": {}, "source": [ "## Creating a raster mosaic with `rioxarray`\n", "\n", "Quite often you need to merge multiple raster files together and create a `raster mosaic`. This can be done easily with the `merge_datasets()` -function in `rioxarray`.\n", "Here, we will create a mosaic based on DEM files (altogether 4 files) covering Kilimanjaro region in Tanzania. First we will read elevation data from a S3 bucket for Kilimanjaro region in Africa." ] }, { "cell_type": "code", "execution_count": 1, "id": "8b42d329", "metadata": {}, "outputs": [], "source": [ "import xarray as xr\n", "import os\n", "import rioxarray as rxr\n", "from rioxarray.merge import merge_datasets\n", "\n", "# S3 bucket containing the data\n", "bucket = \"https://a3s.fi/swift/v1/AUTH_0914d8aff9684df589041a759b549fc2/PythonGIS\"\n", "\n", "# Generate urls for the elevation files\n", "urls = [\n", " os.path.join(bucket, \"elevation/kilimanjaro/ASTGTMV003_S03E036_dem.tif\"),\n", " os.path.join(bucket, \"elevation/kilimanjaro/ASTGTMV003_S03E037_dem.tif\"),\n", " os.path.join(bucket, \"elevation/kilimanjaro/ASTGTMV003_S04E036_dem.tif\"),\n", " os.path.join(bucket, \"elevation/kilimanjaro/ASTGTMV003_S04E037_dem.tif\"),\n", "]\n", "\n", "# Read the files\n", "datasets = [\n", " xr.open_dataset(url, engine=\"rasterio\").squeeze(\"band\", drop=True) for url in urls\n", "]" ] }, { "cell_type": "markdown", "id": "d0f28dcf", "metadata": {}, "source": [ "Investigate how our data looks like:" ] }, { "cell_type": "code", "execution_count": 2, "id": "a4a24ab9", "metadata": {}, "outputs": [ { "data": { "text/html": [ "
<xarray.Dataset>\n", "Dimensions: (x: 3601, y: 3601)\n", "Coordinates:\n", " * x (x) float64 36.0 36.0 36.0 36.0 36.0 ... 37.0 37.0 37.0 37.0\n", " * y (y) float64 -2.0 -2.0 -2.001 -2.001 ... -2.999 -2.999 -3.0 -3.0\n", " spatial_ref int64 0\n", "Data variables:\n", " band_data (y, x) float32 ...
array([36. , 36.000278, 36.000556, ..., 36.999444, 36.999722, 37. ])
array([-2. , -2.000278, -2.000556, ..., -2.999444, -2.999722, -3. ])
array(0)
[12967201 values with dtype=float32]