|
5 | 5 | "metadata": {},
|
6 | 6 | "source": [
|
7 | 7 | "# Segment clouds and cloud shadows in Landsat images (L2A)\n",
|
8 |
| - "This notebook shows an example on how to use [ukis-csmask]() to segment clouds and cloud shadows in Level-2A images from Landsat-9, Landsat-8, Landsat-7 and Landsat-5 satellites. Images are acquired from [Planetary Computer]() and are prepared according to [ukis-csmask]() requirements.\n", |
| 8 | + "This notebook shows an example on how to use [ukis-csmask](https://github.com/dlr-eoc/ukis-csmask) to segment clouds and cloud shadows in Level-2A images from Landsat-9, Landsat-8, Landsat-7 and Landsat-5 satellites. Images are acquired from [Planetary Computer](https://planetarycomputer.microsoft.com) and are preprocessed.\n", |
9 | 9 | "\n",
|
10 | 10 | "> NOTE: to run this notebook, we first need to install some additional dependencies for work with Planetary Computer\n",
|
11 | 11 | "```shell\n",
|
12 |
| - "pip install planetary_computer pystac-client odc-stac tqdm\n", |
13 |
| - "pip install -e git+https://github.com/dlr-eoc/ukis-csmask/tree/l2a#egg=ukis-csmask[gpu]\n", |
| 12 | + "$ pip install planetary_computer rioxarray pystac-client odc-stac tqdm\n", |
14 | 13 | "```"
|
15 | 14 | ]
|
16 | 15 | },
|
17 | 16 | {
|
18 | 17 | "cell_type": "code",
|
19 |
| - "execution_count": 1, |
| 18 | + "execution_count": null, |
20 | 19 | "id": "703f3744-902d-470b-a80f-9a8d3ea08dfa",
|
21 | 20 | "metadata": {},
|
22 | 21 | "outputs": [],
|
23 | 22 | "source": [
|
24 | 23 | "import numpy as np\n",
|
25 | 24 | "import planetary_computer as pc\n",
|
| 25 | + "import rioxarray\n", |
26 | 26 | "\n",
|
27 | 27 | "from odc.stac import load\n",
|
| 28 | + "from pathlib import Path\n", |
28 | 29 | "from pystac_client import Client\n",
|
29 | 30 | "from tqdm import tqdm\n",
|
30 |
| - "from ukis_csmask.mask import CSmask" |
| 31 | + "from ukis_csmask.mask import CSmask\n", |
| 32 | + "from xarray import DataArray" |
31 | 33 | ]
|
32 | 34 | },
|
33 | 35 | {
|
34 | 36 | "cell_type": "code",
|
35 |
| - "execution_count": 13, |
| 37 | + "execution_count": null, |
36 | 38 | "id": "8ca03c78-1e24-479c-9786-a1b43206a08b",
|
37 | 39 | "metadata": {},
|
38 | 40 | "outputs": [],
|
39 | 41 | "source": [
|
| 42 | + "# user settings\n", |
40 | 43 | "stac_api_endpoint = \"https://planetarycomputer.microsoft.com/api/stac/v1\"\n",
|
41 | 44 | "collections = [\"landsat-c2-l2\"]\n",
|
42 | 45 | "ids = [\n",
|
43 |
| - " \"LC09_L2SP_194027_20240604_02_T1\",\n", |
44 |
| - " \"LC08_L2SP_191027_20240607_02_T1\",\n", |
45 |
| - " \"LE07_L2SP_192027_20050805_02_T1\",\n", |
46 |
| - " \"LT05_L2SP_192026_20100726_02_T1\",\n", |
| 46 | + " \"LC09_L2SP_192027_20240403_02_T2\",\n", |
| 47 | + " \"LC08_L2SP_192027_20240801_02_T1\",\n", |
| 48 | + " \"LE07_L2SP_192027_20000722_02_T1\",\n", |
| 49 | + " \"LT05_L2SP_192027_20100726_02_T1\",\n", |
47 | 50 | "]\n",
|
48 |
| - "bbox = [11.016, 47.666, 12.238, 48.472]\n", |
| 51 | + "bbox = [11.540, 47.463, 12.117, 47.872]\n", |
49 | 52 | "product_level = \"l2a\"\n",
|
50 | 53 | "band_order = [\"blue\", \"green\", \"red\", \"nir\", \"swir16\", \"swir22\"]\n",
|
51 |
| - "providers = [\"CUDAExecutionProvider\"]" |
| 54 | + "providers = [\"CUDAExecutionProvider\"]\n", |
| 55 | + "out_dir = \"ukis-csmask/examples\"" |
52 | 56 | ]
|
53 | 57 | },
|
54 | 58 | {
|
55 | 59 | "cell_type": "code",
|
56 |
| - "execution_count": 15, |
| 60 | + "execution_count": null, |
57 | 61 | "id": "7b568942-84e6-4baf-b490-a213b3787f80",
|
58 | 62 | "metadata": {},
|
59 |
| - "outputs": [ |
60 |
| - { |
61 |
| - "name": "stdout", |
62 |
| - "output_type": "stream", |
63 |
| - "text": [ |
64 |
| - "Search returned 4 item(s)\n" |
65 |
| - ] |
66 |
| - } |
67 |
| - ], |
| 63 | + "outputs": [], |
68 | 64 | "source": [
|
69 | 65 | "# search catalog by scene id\n",
|
70 | 66 | "catalog = Client.open(stac_api_endpoint)\n",
|
|
76 | 72 | },
|
77 | 73 | {
|
78 | 74 | "cell_type": "code",
|
79 |
| - "execution_count": 18, |
| 75 | + "execution_count": null, |
80 | 76 | "id": "68eb9c30-06f7-409e-914d-21e00f45de99",
|
81 | 77 | "metadata": {},
|
82 |
| - "outputs": [ |
83 |
| - { |
84 |
| - "name": "stderr", |
85 |
| - "output_type": "stream", |
86 |
| - "text": [ |
87 |
| - "Predict images: 0%| | 0/4 [00:00<?, ?it/s]" |
88 |
| - ] |
89 |
| - }, |
90 |
| - { |
91 |
| - "name": "stderr", |
92 |
| - "output_type": "stream", |
93 |
| - "text": [ |
94 |
| - "Predict images: 0%| | 0/4 [02:25<?, ?it/s]" |
95 |
| - ] |
96 |
| - }, |
97 |
| - { |
98 |
| - "name": "stdout", |
99 |
| - "output_type": "stream", |
100 |
| - "text": [ |
101 |
| - "(6, 6038, 7076) float32\n" |
102 |
| - ] |
103 |
| - }, |
104 |
| - { |
105 |
| - "name": "stderr", |
106 |
| - "output_type": "stream", |
107 |
| - "text": [ |
108 |
| - "\n" |
109 |
| - ] |
110 |
| - }, |
111 |
| - { |
112 |
| - "ename": "FileNotFoundError", |
113 |
| - "evalue": "[Errno 2] No such file or directory: '/mnt/data1/Workspace/python/ukis-csmask/src/ukis-csmask/ukis_csmask/model_6b_l2a.json'", |
114 |
| - "output_type": "error", |
115 |
| - "traceback": [ |
116 |
| - "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", |
117 |
| - "\u001b[0;31mFileNotFoundError\u001b[0m Traceback (most recent call last)", |
118 |
| - "Cell \u001b[0;32mIn[18], line 33\u001b[0m\n\u001b[1;32m 30\u001b[0m \u001b[38;5;28mprint\u001b[39m(arr\u001b[38;5;241m.\u001b[39mshape, arr\u001b[38;5;241m.\u001b[39mdtype)\n\u001b[1;32m 32\u001b[0m \u001b[38;5;66;03m# compute cloud and cloud shadow mask\u001b[39;00m\n\u001b[0;32m---> 33\u001b[0m csmask \u001b[38;5;241m=\u001b[39m \u001b[43mCSmask\u001b[49m\u001b[43m(\u001b[49m\n\u001b[1;32m 34\u001b[0m \u001b[43m \u001b[49m\u001b[43mimg\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mnp\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mmoveaxis\u001b[49m\u001b[43m(\u001b[49m\u001b[43marr\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mto_numpy\u001b[49m\u001b[43m(\u001b[49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m0\u001b[39;49m\u001b[43m,\u001b[49m\u001b[43m \u001b[49m\u001b[38;5;241;43m-\u001b[39;49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 35\u001b[0m \u001b[43m \u001b[49m\u001b[43mband_order\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mband_order\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 36\u001b[0m \u001b[43m \u001b[49m\u001b[43mproduct_level\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mproduct_level\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 37\u001b[0m \u001b[43m \u001b[49m\u001b[43mnodata_value\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;241;43m0\u001b[39;49m\u001b[43m,\u001b[49m\n\u001b[1;32m 38\u001b[0m \u001b[43m \u001b[49m\u001b[43minvalid_buffer\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;241;43m4\u001b[39;49m\u001b[43m,\u001b[49m\n\u001b[1;32m 39\u001b[0m \u001b[43m \u001b[49m\u001b[43mintra_op_num_threads\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;241;43m0\u001b[39;49m\u001b[43m,\u001b[49m\n\u001b[1;32m 40\u001b[0m \u001b[43m \u001b[49m\u001b[43minter_op_num_threads\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;241;43m0\u001b[39;49m\u001b[43m,\u001b[49m\n\u001b[1;32m 41\u001b[0m \u001b[43m \u001b[49m\u001b[43mproviders\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[43mproviders\u001b[49m\u001b[43m,\u001b[49m\n\u001b[1;32m 42\u001b[0m \u001b[43m \u001b[49m\u001b[43mbatch_size\u001b[49m\u001b[38;5;241;43m=\u001b[39;49m\u001b[38;5;241;43m1\u001b[39;49m\u001b[43m,\u001b[49m\n\u001b[1;32m 43\u001b[0m \u001b[43m\u001b[49m\u001b[43m)\u001b[49m\n\u001b[1;32m 45\u001b[0m \u001b[38;5;66;03m# access cloud and cloud shadow mask\u001b[39;00m\n\u001b[1;32m 46\u001b[0m csmask_csm \u001b[38;5;241m=\u001b[39m csmask\u001b[38;5;241m.\u001b[39mcsm\n", |
119 |
| - "File \u001b[0;32m/mnt/data1/Workspace/python/ukis-csmask/src/ukis-csmask/ukis_csmask/mask.py:90\u001b[0m, in \u001b[0;36mCSmask.__init__\u001b[0;34m(self, img, band_order, product_level, nodata_value, invalid_buffer, intra_op_num_threads, inter_op_num_threads, providers, batch_size)\u001b[0m\n\u001b[1;32m 85\u001b[0m \u001b[38;5;28;01melse\u001b[39;00m:\n\u001b[1;32m 86\u001b[0m \u001b[38;5;28;01mraise\u001b[39;00m \u001b[38;5;167;01mTypeError\u001b[39;00m(\n\u001b[1;32m 87\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mband_order must contain at least \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mblue\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m, \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mgreen\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m, \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mred\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m, \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mnir\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m \u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 88\u001b[0m \u001b[38;5;124mf\u001b[39m\u001b[38;5;124m\"\u001b[39m\u001b[38;5;124mand for better performance also \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mswir16\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m and \u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124mswir22\u001b[39m\u001b[38;5;124m'\u001b[39m\u001b[38;5;124m\"\u001b[39m\n\u001b[1;32m 89\u001b[0m )\n\u001b[0;32m---> 90\u001b[0m \u001b[38;5;28;01mwith\u001b[39;00m \u001b[38;5;28;43mopen\u001b[39;49m\u001b[43m(\u001b[49m\u001b[43mPath\u001b[49m\u001b[43m(\u001b[49m\u001b[43mmodel_file\u001b[49m\u001b[43m)\u001b[49m\u001b[38;5;241;43m.\u001b[39;49m\u001b[43mwith_suffix\u001b[49m\u001b[43m(\u001b[49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[38;5;124;43m.json\u001b[39;49m\u001b[38;5;124;43m\"\u001b[39;49m\u001b[43m)\u001b[49m\u001b[43m)\u001b[49m \u001b[38;5;28;01mas\u001b[39;00m f:\n\u001b[1;32m 91\u001b[0m model_dict \u001b[38;5;241m=\u001b[39m json\u001b[38;5;241m.\u001b[39mload(f)\n\u001b[1;32m 93\u001b[0m \u001b[38;5;66;03m# start onnx inference session and load model\u001b[39;00m\n", |
120 |
| - "\u001b[0;31mFileNotFoundError\u001b[0m: [Errno 2] No such file or directory: '/mnt/data1/Workspace/python/ukis-csmask/src/ukis-csmask/ukis_csmask/model_6b_l2a.json'" |
121 |
| - ] |
122 |
| - } |
123 |
| - ], |
| 78 | + "outputs": [], |
124 | 79 | "source": [
|
125 |
| - "import rioxarray\n", |
126 |
| - "\n", |
127 | 80 | "for item in tqdm(items, total=items_cnt, desc=\"Predict images\"):\n",
|
128 | 81 | " # near infrared band has different alias in landsat collections\n",
|
129 | 82 | " bands = [b.replace(\"nir\", \"nir08\") for b in band_order]\n",
|
|
143 | 96 | " .drop_vars(\"time\")\n",
|
144 | 97 | " )\n",
|
145 | 98 | " arr = arr.rename({\"variable\": \"band\"})\n",
|
146 |
| - " arr.rio.to_raster(raster_path=\"l8.tif\", driver=\"COG\")\n", |
147 | 99 | "\n",
|
148 | 100 | " # use band-specific rescale factors to convert DN to reflectance\n",
|
149 | 101 | " for idx, band_name in enumerate(bands):\n",
|
150 | 102 | " band_info = item.assets[band_name].extra_fields[\"raster:bands\"][0]\n",
|
151 | 103 | " arr[idx, :, :] = arr.sel(band=str(band_name)).astype(np.float32) * band_info[\"scale\"]\n",
|
152 | 104 | " arr[idx, :, :] += band_info[\"offset\"]\n",
|
153 | 105 | " arr[idx, :, :] = arr[idx, :, :].clip(min=0.0, max=1.0)\n",
|
154 |
| - " print(arr.shape, arr.dtype)\n", |
155 | 106 | "\n",
|
156 | 107 | " # compute cloud and cloud shadow mask\n",
|
157 | 108 | " csmask = CSmask(\n",
|
|
166 | 117 | " batch_size=1,\n",
|
167 | 118 | " )\n",
|
168 | 119 | "\n",
|
169 |
| - " # access cloud and cloud shadow mask\n", |
170 |
| - " csmask_csm = csmask.csm\n", |
171 |
| - "\n", |
172 |
| - " # access valid mask\n", |
173 |
| - " csmask_valid = csmask.valid" |
| 120 | + " # write image, csm and valid mask to file\n", |
| 121 | + " arr.rio.to_raster(raster_path=Path(out_dir) / Path(f\"{item.id}.tif\"), driver=\"COG\")\n", |
| 122 | + " DataArray(np.squeeze(np.moveaxis(csmask.csm, -1, 0)), coords=arr.sel(band=arr[\"band\"][0]).coords).rio.to_raster(\n", |
| 123 | + " raster_path=Path(out_dir) / Path(f\"{item.id}_csm.tif\"), driver=\"COG\"\n", |
| 124 | + " )\n", |
| 125 | + " DataArray(np.squeeze(np.moveaxis(csmask.valid, -1, 0)), coords=arr.sel(band=arr[\"band\"][0]).coords).rio.to_raster(\n", |
| 126 | + " raster_path=Path(out_dir) / Path(f\"{item.id}_valid.tif\"), driver=\"COG\"\n", |
| 127 | + " )" |
174 | 128 | ]
|
175 | 129 | }
|
176 | 130 | ],
|
|
0 commit comments