

翻訳は機械翻訳により提供されています。提供された翻訳内容と英語版の間で齟齬、不一致または矛盾がある場合、英語版が優先します。

# Sentinel-2 ラスターデータコレクションにアクセスし、土地セグメンテーションを実行する地球観測ジョブを作成します。
<a name="geospatial-demo"></a>

この Python ベースのチュートリアルでは、SDK for Python (Boto3) と Amazon SageMaker Studio Classic ノートブックを使用します。このデモを正常に完了するには、SageMaker 地理空間と Studio Classic を使用するために必要な AWS Identity and Access Management (IAM) アクセス許可があることを確認してください。SageMaker 地理空間では、Studio Classic にアクセスできるユーザー、グループ、またはロールが必要です。また、信頼ポリシーに、SageMaker 地理空間サービスプリンシパルである `sagemaker-geospatial.amazonaws.com` を指定する SageMaker AI 実行ロールも必要です。

これらの要件の詳細については、「[SageMaker geospatial IAM roles](sagemaker-geospatial-roles.md)」を参照してください。

このチュートリアルでは、SageMaker 地理空間 API を使用して次のタスクを完了する方法について説明します。
+ `list_raster_data_collections` で利用可能なラスターデータコレクションを検索する。
+ `search_raster_data_collection` を使用して、指定されたラスターデータコレクションを検索する。
+ `start_earth_observation_job` を使用して地球観測ジョブ (EOJ) を作成する。

## `list_raster_data_collections` を使用して利用可能なデータコレクションを検索する
<a name="demo-use-list-rdc"></a>

SageMaker 地理空間は、複数のラスターデータコレクションをサポートしています。利用可能なデータコレクションの詳細については、「[データコレクション](geospatial-data-collections.md)」を参照してください。

このデモでは、[Sentinel-2Cloud-Optimized GeoTIFF](https://registry.opendata.aws/sentinel-2-l2a-cogs/) 衛星から収集された衛星データを使用します。これらの衛星は、5 日ごとに世界中をカバーする地球の地表を提供します。地球の表面画像を収集することに加えて、Sentinel-2 衛星はさまざまなスペクトルバンドのデータも収集します。

対象エリア (AOI) を検索するには、Sentinel-2 衛星データに関連付けられている ARN が必要です。で使用可能なデータ収集と関連する ARNs を検索するには AWS リージョン、 `list_raster_data_collections` API オペレーションを使用します。

レスポンスはページ分割される場合があるため、 `get_paginator` オペレーションを使用して、関連するすべてのデータを返す必要があります。

```
import boto3
import sagemaker
import sagemaker_geospatial_map
import json 

## SageMaker Geospatial  is currently only avaialable in US-WEST-2  
session = boto3.Session(region_name='us-west-2')
execution_role = sagemaker.get_execution_role()

## Creates a SageMaker Geospatial client instance 
geospatial_client = session.client(service_name="sagemaker-geospatial")

# Creates a resusable Paginator for the list_raster_data_collections API operation 
paginator = geospatial_client.get_paginator("list_raster_data_collections")

# Create a PageIterator from the paginator class
page_iterator = paginator.paginate()

# Use the iterator to iterate throught the results of list_raster_data_collections
results = []
for page in page_iterator:
    results.append(page['RasterDataCollectionSummaries'])

print(results)
```

これは API `list_raster_data_collections` オペレーションからの JSON レスポンスサンプルです。コード例で使用されているデータコレクション (Sentinel-2) のみを含めるように切り捨てられています。特定のラスターデータコレクションの詳細を得るには、`get_raster_data_collection` を使用します。

```
{
    "Arn": "arn:aws:sagemaker-geospatial:us-west-2:378778860802:raster-data-collection/public/nmqj48dcu3g7ayw8",
    "Description": "Sentinel-2a and Sentinel-2b imagery, processed to Level 2A (Surface Reflectance) and converted to Cloud-Optimized GeoTIFFs",
    "DescriptionPageUrl": "https://registry.opendata.aws/sentinel-2-l2a-cogs",
    "Name": "Sentinel 2 L2A COGs",
    "SupportedFilters": [
        {
            "Maximum": 100,
            "Minimum": 0,
            "Name": "EoCloudCover",
            "Type": "number"
        },
        {
            "Maximum": 90,
            "Minimum": 0,
            "Name": "ViewOffNadir",
            "Type": "number"
        },
        {
            "Name": "Platform",
            "Type": "string"
        }
    ],
    "Tags": {},
    "Type": "PUBLIC"
}
```

前のコードサンプルを実行すると、Sentinel-2 ラスターデータコレクションの ARN、`arn:aws:sagemaker-geospatial:us-west-2:378778860802:raster-data-collection/public/nmqj48dcu3g7ayw8` が得られます。[次のセクション](#demo-search-raster-data) では、`search_raster_data_collection` API を使用して Sentinel-2 データコレクションのクエリを実行できます。

## `search_raster_data_collection` を使用して Sentinel-2 ラスターデータコレクションを検索する
<a name="demo-search-raster-data"></a>

前のセクションでは、 `list_raster_data_collections` を使用して Sentinel-2 データコレクションの ARN を取得しました。そのためその ARN を使用して、特定の対象エリア (AOI)、時間範囲、プロパティ、使用可能な UV バンドのデータコレクションを検索できるようになりました。

`search_raster_data_collection` API を呼び出すには、Python ディクショナリ内で `RasterDataCollectionQuery` パラメータに渡す必要があります。この例では、`AreaOfInterest`、`TimeRangeFilter`、`PropertyFilters`、`BandFilter` を使用します。簡単にするために、**search\$1rdc\$1query** 変数を使用して Python ディクショナリを指定して、検索クエリパラメータを保持できます。

```
search_rdc_query = {
    "AreaOfInterest": {
        "AreaOfInterestGeometry": {
            "PolygonGeometry": {
                "Coordinates": [
                    [
                        # coordinates are input as longitute followed by latitude 
                        [-114.529, 36.142],
                        [-114.373, 36.142],
                        [-114.373, 36.411],
                        [-114.529, 36.411],
                        [-114.529, 36.142],
                    ]
                ]
            }
        }
    },
    "TimeRangeFilter": {
        "StartTime": "2022-01-01T00:00:00Z",
        "EndTime": "2022-07-10T23:59:59Z"
    },
    "PropertyFilters": {
        "Properties": [
            {
                "Property": {
                    "EoCloudCover": {
                        "LowerBound": 0,
                        "UpperBound": 1
                    }
                }
            }
        ],
        "LogicalOperator": "AND"
    },
    "BandFilter": [
        "visual"
    ]
}
```

この例では、ユタ州の [Lake Mead](https://en.wikipedia.org/wiki/Lake_Mead) を含む `AreaOfInterest` のクエリを実行します。さらに、Sentinel-2 では複数の種類の画像バンドをサポートしています。水面の変化を測定するには、`visual` バンドのみが必要です。

クエリパラメータを作成したら、 `search_raster_data_collection` API を使用してリクエストを実行できます。

次のコードサンプルは、`search_raster_data_collection` API リクエストを実装します。この API は、`get_paginator` API を使用したページ分割をサポートしていません。完全な API レスポンスが収集されていることを確認するために、コードサンプルでは `while` ループを使用して、`NextToken` が存在することを確認します。次に、`.extend()` を使用して、衛星画像の URL およびその他のレスポンスメタデータを `items_list` に追加します。

`search_raster_data_collection` の詳細については、「*Amazon SageMaker AI API リファレンス*」の「[SearchRasterDataCollection](https://docs.aws.amazon.com/sagemaker/latest/APIReference/API_geospatial_SearchRasterDataCollection.html)」を参照してください。

```
search_rdc_response = sm_geo_client.search_raster_data_collection(
    Arn='arn:aws:sagemaker-geospatial:us-west-2:378778860802:raster-data-collection/public/nmqj48dcu3g7ayw8',
    RasterDataCollectionQuery=search_rdc_query
)


## items_list is the response from the API request. 
items_list = []

## Use the python .get() method to check that the 'NextToken' exists, if null returns None breaking the while loop 
while search_rdc_response.get('NextToken'):
    items_list.extend(search_rdc_response['Items'])
    search_rdc_response = sm_geo_client.search_raster_data_collection(
        Arn='arn:aws:sagemaker-geospatial:us-west-2:378778860802:raster-data-collection/public/nmqj48dcu3g7ayw8',
        RasterDataCollectionQuery=search_rdc_query, NextToken=search_rdc_response['NextToken']
    )

## Print the number of observation return based on the query
print (len(items_list))
```

以下に、 クエリからの JSON レスポンスを示します。わかりやすくするために切り捨てられています。リクエストで指定された **"BandFilter": ["visual"]** のみが `Assets` のキーと値のペアに返されます。

```
{
    'Assets': {
        'visual': {
            'Href': 'https://sentinel-cogs.s3.us-west-2.amazonaws.com/sentinel-s2-l2a-cogs/15/T/UH/2022/6/S2A_15TUH_20220623_0_L2A/TCI.tif'
        }
    },
    'DateTime': datetime.datetime(2022, 6, 23, 17, 22, 5, 926000, tzinfo = tzlocal()),
    'Geometry': {
        'Coordinates': [
            [
                [-114.529, 36.142],
                [-114.373, 36.142],
                [-114.373, 36.411],
                [-114.529, 36.411],
                [-114.529, 36.142],
            ]
        ],
        'Type': 'Polygon'
    },
    'Id': 'S2A_15TUH_20220623_0_L2A',
    'Properties': {
        'EoCloudCover': 0.046519,
        'Platform': 'sentinel-2a'
    }
}
```

クエリ結果が得られたら、次のセクションでは `matplotlib` を使用して結果を視覚化できます。これは、結果が正しい地域からもたらされたものであることを検証するためです。

## `matplotlib` を使用して `search_raster_data_collection` を視覚化する
<a name="demo-geospatial-visualize"></a>

地球観測ジョブ (EOJ) を開始する前に、`matplotlib` を使用してクエリの結果を視覚化できます。次のコードサンプルは、前のコードサンプルで作成された `items_list` 変数からの最初の項目 `items_list[0]["Assets"]["visual"]["Href"]` を取得し、`matplotlib` を使用して画像を印刷します。

```
# Visualize an example image.
import os
from urllib import request
import tifffile
import matplotlib.pyplot as plt

image_dir = "./images/lake_mead"
os.makedirs(image_dir, exist_ok=True)

image_dir = "./images/lake_mead"
os.makedirs(image_dir, exist_ok=True)

image_url = items_list[0]["Assets"]["visual"]["Href"]
img_id = image_url.split("/")[-2]
path_to_image = image_dir + "/" + img_id + "_TCI.tif"
response = request.urlretrieve(image_url, path_to_image)
print("Downloaded image: " + img_id)

tci = tifffile.imread(path_to_image)
plt.figure(figsize=(6, 6))
plt.imshow(tci)
plt.show()
```

結果が適切な地域にあることを確認したら、次のステップで地球観測ジョブ (EOJ) を開始できます。EOJ を使用し、土地セグメンテーションと呼ばれるプロセスを使用して衛星画像から水域を特定します。

## 一連の衛星画像で土地セグメンテーションを実行する地球観測ジョブ (EOJ) を開始する
<a name="demo-start-eoj"></a>

SageMaker 地理空間は、ラスターデータコレクションからの地理空間データを処理するために使用できる事前トレーニング済みの複数のモデルを提供します。利用可能な事前トレーニング済みモデルとカスタムオペレーションの詳細については、「[操作タイプ](geospatial-eoj-models.md)」を参照してください。

水面面積の変化を計算するには、画像内のどのピクセルが水に対応しているかを特定する必要があります。土地被覆セグメンテーションは、`start_earth_observation_job` API でサポートされているセマンティックセグメンテーションモデルです。セマンティックセグメンテーションモデルは、各画像のすべてのピクセルにラベルを関連付けます。その結果、各ピクセルには、モデルのクラスマップに基づくラベルが割り当てられます。土地セグメンテーションモデルのクラスマップを以下に示します。

```
{
    0: "No_data",
    1: "Saturated_or_defective",
    2: "Dark_area_pixels",
    3: "Cloud_shadows",
    4: "Vegetation",
    5: "Not_vegetated",
    6: "Water",
    7: "Unclassified",
    8: "Cloud_medium_probability",
    9: "Cloud_high_probability",
    10: "Thin_cirrus",
    11: "Snow_ice"
}
```

地球観測ジョブを開始するには、`start_earth_observation_job` API を使用します。リクエストを送信するときは、以下を指定する必要があります。
+ `InputConfig` (*dict*) – 検索するエリアの座標と、検索に関連付けられている他のメタデータを指定するために使用します。
+ `JobConfig` (*dict*) – データに対して実行した EOJ オペレーションのタイプを指定するために使用します。この例では **LandCoverSegmentationConfig** を使用します。
+ `ExecutionRoleArn` (*string*) – ジョブを実行するために必要なアクセス許可を持つ SageMaker AI 実行ロールの ARN。
+ `Name` (*string*) – 地球観測ジョブの名前。

`InputConfig` は Python ディクショナリです。次の変数 **eoj\$1input\$1config** を使用して、検索クエリパラメータを保持します。`start_earth_observation_job` API リクエストを行うときは、この変数を使用します。

```
# Perform land cover segmentation on images returned from the Sentinel-2 dataset.
eoj_input_config = {
    "RasterDataCollectionQuery": {
        "RasterDataCollectionArn": "arn:aws:sagemaker-geospatial:us-west-2:378778860802:raster-data-collection/public/nmqj48dcu3g7ayw8",
        "AreaOfInterest": {
            "AreaOfInterestGeometry": {
                "PolygonGeometry": {
                    "Coordinates":[
                        [
                            [-114.529, 36.142],
                            [-114.373, 36.142],
                            [-114.373, 36.411],
                            [-114.529, 36.411],
                            [-114.529, 36.142],
                        ]
                    ]
                }
            }
        },
        "TimeRangeFilter": {
            "StartTime": "2021-01-01T00:00:00Z",
            "EndTime": "2022-07-10T23:59:59Z",
        },
        "PropertyFilters": {
            "Properties": [{"Property": {"EoCloudCover": {"LowerBound": 0, "UpperBound": 1}}}],
            "LogicalOperator": "AND",
        },
    }
}
```

`JobConfig` は、データに対して実行する EOJ オペレーションを指定するために使用する Python ディクショナリです。

```
eoj_config = {"LandCoverSegmentationConfig": {}}
```

ディクショナリ要素を指定したら、次のコードサンプルを使用して `start_earth_observation_job` API リクエストを送信できます。

```
# Gets the execution role arn associated with current notebook instance 
execution_role_arn = sagemaker.get_execution_role()

# Starts an earth observation job
response = sm_geo_client.start_earth_observation_job(
    Name="lake-mead-landcover",
    InputConfig=eoj_input_config,
    JobConfig=eoj_config,
    ExecutionRoleArn=execution_role_arn,
)
            
print(response)
```

地球観測ジョブを開始すると、ARN と他のメタデータが返されます。

進行中のすべての地球観測ジョブと現在の地球観測ジョブのリストを取得するには、`list_earth_observation_jobs` API を使用します。単一の地球観測ジョブのステータスを監視するには、`get_earth_observation_job` API を使用します。このリクエストを行うには、EOJ リクエストの送信後に作成された ARN を使用します。詳細については、「*Amazon SageMaker AI API リファレンス*」の「[GetEarthObservationJob](https://docs.aws.amazon.com/sagemaker/latest/APIReference/API_geospatial_GetEarthObservationJob.html)」を参照してください。

EOJ に関連付けられている ARN を検索するには、`list_earth_observation_jobs` API オペレーションを使用します。詳細については、「*Amazon SageMaker AI API リファレンス*」の「[ListEarthObservationJobs](https://docs.aws.amazon.com//sagemaker/latest/APIReference/API_geospatial_ListEarthObservationJobs.html)」を参照してください。

```
# List all jobs in the account
sg_client.list_earth_observation_jobs()["EarthObservationJobSummaries"]
```

以下に、JSON レスポンスの例を示します。

```
{
    'Arn': 'arn:aws:sagemaker-geospatial:us-west-2:111122223333:earth-observation-job/futg3vuq935t',
    'CreationTime': datetime.datetime(2023, 10, 19, 4, 33, 54, 21481, tzinfo = tzlocal()),
    'DurationInSeconds': 3493,
    'Name': 'lake-mead-landcover',
    'OperationType': 'LAND_COVER_SEGMENTATION',
    'Status': 'COMPLETED',
    'Tags': {}
}, {
    'Arn': 'arn:aws:sagemaker-geospatial:us-west-2:111122223333:earth-observation-job/wu8j9x42zw3d',
    'CreationTime': datetime.datetime(2023, 10, 20, 0, 3, 27, 270920, tzinfo = tzlocal()),
    'DurationInSeconds': 1,
    'Name': 'mt-shasta-landcover',
    'OperationType': 'LAND_COVER_SEGMENTATION',
    'Status': 'INITIALIZING',
    'Tags': {}
}
```

EOJ ジョブのステータスが `COMPLETED` に変わったら、次のセクションに進み、Lake Mead's  の表面積の変化を計算します。

## Lake Mead の表面積の変化を計算する
<a name="demo-geospatial-calc"></a>

Lake Mead の表面積の変化を計算するには、まず `export_earth_observation_job` を使用して EOJ の結果を Amazon S3 にエクスポートします。

```
sagemaker_session = sagemaker.Session()
s3_bucket_name = sagemaker_session.default_bucket()  # Replace with your own bucket if needed
s3_bucket = session.resource("s3").Bucket(s3_bucket_name)
prefix = "export-lake-mead-eoj"  # Replace with the S3 prefix desired
export_bucket_and_key = f"s3://{s3_bucket_name}/{prefix}/"

eoj_output_config = {"S3Data": {"S3Uri": export_bucket_and_key}}
export_response = sm_geo_client.export_earth_observation_job(
    Arn="arn:aws:sagemaker-geospatial:us-west-2:111122223333:earth-observation-job/7xgwzijebynp",
    ExecutionRoleArn=execution_role_arn,
    OutputConfig=eoj_output_config,
    ExportSourceImages=False,
)
```

エクスポートジョブのステータスを確認するには、`get_earth_observation_job` を使用します。

```
export_job_details = sm_geo_client.get_earth_observation_job(Arn=export_response["Arn"])
```

Lake Mead の水位の変化を計算するには、土地被覆マスクをローカルの SageMaker ノートブックインスタンスにダウンロードし、前のクエリのソース画像をダウンロードします。土地セグメンテーションモデルのクラスマップでは、水のクラスインデックスは 6 です。

Sentinel-2 画像から水のマスクを抽出するには、次の手順に従います。まず、画像内で水 (クラスインデックス 6) としてマークされたピクセルの数をカウントします。次に、その数に各ピクセルがカバーする面積を掛けます。バンドの空間解像度は異なる場合があります。土地被覆セグメンテーションモデルでは、すべてのバンドが 60 メートルに等しい空間解像度までダウンサンプリングされます。

```
import os
from glob import glob
import cv2
import numpy as np
import tifffile
import matplotlib.pyplot as plt
from urllib.parse import urlparse
from botocore import UNSIGNED
from botocore.config import Config

# Download land cover masks
mask_dir = "./masks/lake_mead"
os.makedirs(mask_dir, exist_ok=True)
image_paths = []
for s3_object in s3_bucket.objects.filter(Prefix=prefix).all():
    path, filename = os.path.split(s3_object.key)
    if "output" in path:
        mask_name = mask_dir + "/" + filename
        s3_bucket.download_file(s3_object.key, mask_name)
        print("Downloaded mask: " + mask_name)

# Download source images for visualization
for tci_url in tci_urls:
    url_parts = urlparse(tci_url)
    img_id = url_parts.path.split("/")[-2]
    tci_download_path = image_dir + "/" + img_id + "_TCI.tif"
    cogs_bucket = session.resource(
        "s3", config=Config(signature_version=UNSIGNED, region_name="us-west-2")
    ).Bucket(url_parts.hostname.split(".")[0])
    cogs_bucket.download_file(url_parts.path[1:], tci_download_path)
    print("Downloaded image: " + img_id)

print("Downloads complete.")

image_files = glob("images/lake_mead/*.tif")
mask_files = glob("masks/lake_mead/*.tif")
image_files.sort(key=lambda x: x.split("SQA_")[1])
mask_files.sort(key=lambda x: x.split("SQA_")[1])
overlay_dir = "./masks/lake_mead_overlay"
os.makedirs(overlay_dir, exist_ok=True)
lake_areas = []
mask_dates = []

for image_file, mask_file in zip(image_files, mask_files):
    image_id = image_file.split("/")[-1].split("_TCI")[0]
    mask_id = mask_file.split("/")[-1].split(".tif")[0]
    mask_date = mask_id.split("_")[2]
    mask_dates.append(mask_date)
    assert image_id == mask_id
    image = tifffile.imread(image_file)
    image_ds = cv2.resize(image, (1830, 1830), interpolation=cv2.INTER_LINEAR)
    mask = tifffile.imread(mask_file)
    water_mask = np.isin(mask, [6]).astype(np.uint8)  # water has a class index 6
    lake_mask = water_mask[1000:, :1100]
    lake_area = lake_mask.sum() * 60 * 60 / (1000 * 1000)  # calculate the surface area
    lake_areas.append(lake_area)
    contour, _ = cv2.findContours(water_mask, cv2.RETR_TREE, cv2.CHAIN_APPROX_SIMPLE)
    combined = cv2.drawContours(image_ds, contour, -1, (255, 0, 0), 4)
    lake_crop = combined[1000:, :1100]
    cv2.putText(lake_crop, f"{mask_date}", (10,50), cv2.FONT_HERSHEY_SIMPLEX, 1.5, (0, 0, 0), 3, cv2.LINE_AA)
    cv2.putText(lake_crop, f"{lake_area} [sq km]", (10,100), cv2.FONT_HERSHEY_SIMPLEX, 1.5, (0, 0, 0), 3, cv2.LINE_AA)
    overlay_file = overlay_dir + '/' + mask_date + '.png'
    cv2.imwrite(overlay_file, cv2.cvtColor(lake_crop, cv2.COLOR_RGB2BGR))

# Plot water surface area vs. time.
plt.figure(figsize=(20,10))
plt.title('Lake Mead surface area for the 2021.02 - 2022.07 period.', fontsize=20)
plt.xticks(rotation=45)
plt.ylabel('Water surface area [sq km]', fontsize=14)
plt.plot(mask_dates, lake_areas, marker='o')
plt.grid('on')
plt.ylim(240, 320)
for i, v in enumerate(lake_areas):
    plt.text(i, v+2, "%d" %v, ha='center')
plt.show()
```

`matplotlib` を使用すると、結果をグラフで視覚化できます。このグラフは、Lake Mead の面積が 2021 年 1 月から 2022 年 7 月にかけて減少したことを示しています。

![\[Lake Mead の面積が 2021 年 1 月から 2022 年 7 月にかけて減少したことを示す線グラフ\]](http://docs.aws.amazon.com/ja_jp/sagemaker/latest/dg/images/lake-mead-decrease.png)
