ibm-granite/granite-geospatial-biomassのメモ①

[ibm-granite/granite-geospatial-biomass](https://huggingface.co/ibm-granite/granite-geospatial-biomass
)のZero-shot for all biomesまでをGoogle Colabで実施してみたので以下メモ。

```python
# -*- coding: utf-8 -*-
# Commented out IPython magic to ensure Python compatibility.
# %load_ext autoreload
# %autoreload 2

from google.colab import drive
drive.mount('/content/drive')

# Clone the ibm-granite/granite-geospatial-biomass GitHub
!git clone https://github.com/ibm-granite/granite-geospatial-biomass.git

# Commented out IPython magic to ensure Python compatibility.
# Change directory. Move inside the tsfm repo.
# %cd granite-geospatial-biomass/

# Install the package
!pip install -e .

# Commented out IPython magic to ensure Python compatibility.
# Change directory. Standardize so below paths work with local notebook and colab.
# %cd notebooks

# Download dataset from Zenodo
dataset_path = '../granite-geospatial-biomass-datasets.tar.gz'
!wget "https://zenodo.org/records/12356481/files/granite-geospatial-biomass-datasets.tar.gz?download=1" -O {dataset_path}

# Commented out IPython magic to ensure Python compatibility.
# %mv ../granite-geospatial-biomass-datasets.tar.gz /content/drive/MyDrive/

# Commented out IPython magic to ensure Python compatibility.
# %ls /content/drive/MyDrive/granite-geospatial-biomass-datasets.tar.gz

# Unpack compressed dataset
dataset_path = '/content/drive/MyDrive/granite-geospatial-biomass-datasets.tar.gz'
target_directory = '../'
!tar -xvf {dataset_path} --directory {target_directory}

# Commented out IPython magic to ensure Python compatibility.
# %pwd

# Standard
import os
import shutil
import glob
import matplotlib.pyplot as plt
import yaml
import pandas as pd
import numpy as np

# Third Party
import rioxarray as rio
from lightning.pytorch import Trainer
from terratorch.cli_tools import LightningInferenceModel
from huggingface_hub import hf_hub_download

# First Party
from utils.binwise_rmse import calc_binwise_rmse
from utils.dataset_scalers import calc_mean_std
from utils.plotting import plot_rgb_agb_gedi

ckpt_path = hf_hub_download(repo_id="ibm-granite/granite-geospatial-biomass", filename="biomass_model.ckpt")

### Provide all necessary fils, paths, and hyperparameter

# Path to configuration file which contains all hyperparameters
config_path = '../configs/config.yaml'

# Path to directory with geotiff test images
predict_input_dir = '../granite-geospatial-biomass-datasets/all_ecos_datasplit/test_images/'

# List to define the bands in the input images. As currently configured, the model looks for following
# HLS bands: BLUE, GREEN, RED, NIR_NARROW, SWIR_1, SWIR_2
# The line below names all the bands in the input, so the ones above can be extracted. we use -1 for placeholders, as we dont care about those
UNUSED_BAND = "-1"
predict_dataset_bands = [UNUSED_BAND,"BLUE","GREEN","RED","NIR_NARROW","SWIR_1","SWIR_2",UNUSED_BAND,UNUSED_BAND,UNUSED_BAND,UNUSED_BAND]

# Path to directory with GEDI test labels - these will be used for plotting and evaluation further below
test_label_dir = '../granite-geospatial-biomass-datasets/all_ecos_datasplit/test_labels/'

# Commented out IPython magic to ensure Python compatibility.
# %cat ../configs/config.yaml

# Commented out IPython magic to ensure Python compatibility.
# %rm -rf inference_images/  inference_labels/

# create subset of images for inference
images_for_inference = ['T10SFF_144_tile_img.tif', 'T33TUM_34_tile_img.tif', 'T47PRR_79_tile_img.tif']

labels_for_inference = ['T10SFF_144_tile_label.tif', 'T33TUM_34_tile_label.tif', 'T47PRR_79_tile_label.tif']

if not os.path.isdir("inference_images"):
    os.mkdir("inference_images")
    for inference_image in images_for_inference:
        shutil.copy(os.path.join(predict_input_dir, inference_image), os.path.join("inference_images", inference_image))

if not os.path.isdir("inference_labels"):
    os.mkdir("inference_labels")
    for inference_label in labels_for_inference:
        shutil.copy(os.path.join(test_label_dir, inference_label), os.path.join("inference_labels", inference_label))

# Commented out IPython magic to ensure Python compatibility.
# %ls inference_images inference_labels

model = LightningInferenceModel.from_config(config_path, ckpt_path, predict_dataset_bands)

inference_results, input_file_names = model.inference_on_dir("inference_images")
inference_results
inference_results.shape
input_file_names

# 推定結果の可視化
i = 0
tmp = inference_results.numpy()[i]
plt.imshow(tmp, cmap='Greens')
plt.title(label=input_file_names[i])
plt.colorbar()

i = 1
tmp = inference_results.numpy()[i]
plt.imshow(tmp, cmap='Greens')
plt.title(label=input_file_names[i])
plt.colorbar()

i = 2
tmp = inference_results.numpy()[i]
plt.imshow(tmp, cmap='Greens')
plt.title(label=input_file_names[i])
plt.colorbar()

# 結果の保存
import torch
torch.save(inference_results, "inference_results.pt")

# Commented out IPython magic to ensure Python compatibility.
# %ls inference_results.pt

with open('input_file_names.txt', 'w') as f:
    for file_name in input_file_names:
        f.write(file_name + '\n')

# Commented out IPython magic to ensure Python compatibility.
# %cat input_file_names.txt

# Commented out IPython magic to ensure Python compatibility.
# %mv inference_results.pt input_file_names.txt /content/drive/MyDrive/

# Commented out IPython magic to ensure Python compatibility.
# %cp -irp ./inference_images /content/drive/MyDrive/
```

入力データの仕様と推定結果の画像
T10SFF_144_tile_img.tif 

```bash
$ gdalinfo T10SFF_144_tile_img.tif 
Driver: GTiff/GeoTIFF
Files: T10SFF_144_tile_img.tif
Size is 224, 224
Coordinate System is:
PROJCRS["WGS 84 / UTM zone 10N",
    BASEGEOGCRS["WGS 84",
        DATUM["World Geodetic System 1984",
            ELLIPSOID["WGS 84",6378137,298.257223563,
                LENGTHUNIT["metre",1]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4326]],
    CONVERSION["UTM zone 10N",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",-123,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",0.9996,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",500000,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["(E)",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["(N)",north,
            ORDER[2],
            LENGTHUNIT["metre",1]],
    USAGE[
        SCOPE["Navigation and medium accuracy spatial referencing."],
        AREA["Between 126°W and 120°W, northern hemisphere between equator and 84°N, onshore and offshore. Canada - British Columbia (BC); Northwest Territories (NWT); Nunavut; Yukon. United States (USA) - Alaska (AK)."],
        BBOX[0,-126,84,-120]],
    ID["EPSG",32610]]
Data axis to CRS axis mapping: 1,2
Origin = (600000.000000000000000,4039560.000000000000000)
Pixel Size = (30.000000000000000,-30.000000000000000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  600000.000, 4039560.000) (121d53' 0.41"W, 36d29'46.12"N)
Lower Left  (  600000.000, 4032840.000) (121d53' 3.54"W, 36d26' 8.07"N)
Upper Right (  606720.000, 4039560.000) (121d48'30.34"W, 36d29'43.51"N)
Lower Right (  606720.000, 4032840.000) (121d48'33.68"W, 36d26' 5.46"N)
Center      (  603360.000, 4036200.000) (121d50'46.99"W, 36d27'55.81"N)
Band 1 Block=224x1 Type=Float32, ColorInterp=Gray
  NoData Value=-9999
Band 2 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 3 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 4 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 5 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 6 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 7 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 8 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 9 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 10 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3756172/a42620c1-2dfb-6ceb-cb8f-5e01551ce9dd.png)

T33TUM_34_tile_img.tif 

```bash
% gdalinfo T33TUM_34_tile_img.tif 
Driver: GTiff/GeoTIFF
Files: T33TUM_34_tile_img.tif
Size is 224, 224
Coordinate System is:
PROJCRS["WGS 84 / UTM zone 33N",
    BASEGEOGCRS["WGS 84",
        DATUM["World Geodetic System 1984",
            ELLIPSOID["WGS 84",6378137,298.257223563,
                LENGTHUNIT["metre",1]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4326]],
    CONVERSION["UTM zone 33N",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",15,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",0.9996,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",500000,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["(E)",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["(N)",north,
            ORDER[2],
            LENGTHUNIT["metre",1]],
    USAGE[
        SCOPE["Navigation and medium accuracy spatial referencing."],
        AREA["Between 12°E and 18°E, northern hemisphere between equator and 84°N, onshore and offshore. Austria. Bosnia and Herzegovina. Cameroon. Central African Republic. Chad. Congo. Croatia. Czechia. Democratic Republic of the Congo (Zaire). Gabon. Germany. Hungary. Italy. Libya. Malta. Niger. Nigeria. Norway. Poland. San Marino. Slovakia. Slovenia. Svalbard. Sweden. Vatican City State."],
        BBOX[0,12,84,18]],
    ID["EPSG",32633]]
Data axis to CRS axis mapping: 1,2
Origin = (313440.000000000000000,5186580.000000000000000)
Pixel Size = (30.000000000000000,-30.000000000000000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  313440.000, 5186580.000) ( 12d33'17.46"E, 46d48'23.90"N)
Lower Left  (  313440.000, 5179860.000) ( 12d33'27.31"E, 46d44'46.39"N)
Upper Right (  320160.000, 5186580.000) ( 12d38'34.25"E, 46d48'30.55"N)
Lower Right (  320160.000, 5179860.000) ( 12d38'43.75"E, 46d44'53.03"N)
Center      (  316800.000, 5183220.000) ( 12d36' 0.69"E, 46d46'38.50"N)
Band 1 Block=224x1 Type=Float32, ColorInterp=Gray
  NoData Value=-9999
Band 2 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 3 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 4 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 5 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 6 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 7 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 8 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 9 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 10 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3756172/7f334c97-c1a8-d0ab-3eef-3e3f6231d7de.png)

T47PRR_79_tile_img.tif

```bash
$ gdalinfo T47PRR_79_tile_img.tif 
Driver: GTiff/GeoTIFF
Files: T47PRR_79_tile_img.tif
Size is 224, 224
Coordinate System is:
PROJCRS["WGS 84 / UTM zone 47N",
    BASEGEOGCRS["WGS 84",
        DATUM["World Geodetic System 1984",
            ELLIPSOID["WGS 84",6378137,298.257223563,
                LENGTHUNIT["metre",1]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4326]],
    CONVERSION["UTM zone 47N",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",99,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",0.9996,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",500000,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["(E)",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["(N)",north,
            ORDER[2],
            LENGTHUNIT["metre",1]],
    USAGE[
        SCOPE["Navigation and medium accuracy spatial referencing."],
        AREA["Between 96°E and 102°E, northern hemisphere between equator and 84°N, onshore and offshore. China. Indonesia. Laos. Malaysia - West Malaysia. Mongolia. Myanmar (Burma). Russian Federation. Thailand."],
        BBOX[0,96,84,102]],
    ID["EPSG",32647]]
Data axis to CRS axis mapping: 1,2
Origin = (900780.000000000000000,1573140.000000000000000)
Pixel Size = (30.000000000000000,-30.000000000000000)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  COMPRESSION=LZW
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (  900780.000, 1573140.000) (102d42'43.89"E, 14d12' 3.46"N)
Lower Left  (  900780.000, 1566420.000) (102d42'40.34"E, 14d 8'25.18"N)
Upper Right (  907500.000, 1573140.000) (102d46'27.62"E, 14d11'59.95"N)
Lower Right (  907500.000, 1566420.000) (102d46'24.01"E, 14d 8'21.69"N)
Center      (  904140.000, 1569780.000) (102d44'33.96"E, 14d10'12.57"N)
Band 1 Block=224x1 Type=Float32, ColorInterp=Gray
  NoData Value=-9999
Band 2 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 3 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 4 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 5 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 6 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 7 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 8 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 9 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
Band 10 Block=224x1 Type=Float32, ColorInterp=Undefined
  NoData Value=-9999
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/3756172/0b6956eb-b83e-6d67-1156-c9d642ce856b.png)

入力データのオリジナルは[hlsl30v002](https://lpdaac.usgs.gov/products/hlsl30v002/)。