土地被覆分類
ATBD とプログラム・テストデータ説明書（ ver. 0930 ）福江潔也
東海大学情報技術センター

2011 年 9 月 30 日

目

Ⅰ. ATBD次

1…………………………………………………………………………

1. 分類の原理（基本的な考え方）
1.1 分類対象データ…………………………………………………1………………………………………………………………1

1.2 分類に使用する特徴量
1.3 識別法（識別器）
1.4 分類結果

2………………………………………………………
………………………………………………………………4………………………………………………………………………5

2. ファジィ化

6………………………………………………………………………
2.1 空間的ファジィ化
6………………………………………………………………
2.2 時間的ファジィ化 ………………………………………………………………7
3. 分類対象画素フィルタリング
3.1 陸域・水域フィルタリング
3.2 雲・雪・氷フィルタリング……………………………………………………88………………………………………………………
………………………………………………………8

4. システム構成と入出力データ…………………………………………………10

Ⅱ. プログラムとテストデータ……………………………………………………11

1. フォルダ構造と記憶データ……………………………………………………

11

1.1 概要…………………………………………………………………………11
1.2 フォルダ“ 1.プログラム実行フォルダ”の内容
12………………………………
1.3 フォルダ“ 2.分類対象データ”の内容
13………………………………………
1.4 フォルダ“ 3.リストファイル”の内容 ……………………………………………16
1.5 フォルダ“ 4.トレーニングデータ”の内容
1.7 フォルダ“ 6.分類結果画像”の内容……………………………………1617………………………………………
…………………………………………17

1.8 フォルダ“ 7.ソースプログラム”の内容………………………………………17

1.6 フォルダ“ 5.カラーテーブル”の内容

2. 分類処理の入出力データ
2.1 入力データファイル
2.2 出力データファイル

18……………………………………………………
18…………………………………………………………
…………………………………………………………19

3. 分類プログラムの実行方法
4. 分類結果の表示方法

25…………………………………………………………

5. 分類プログラムのビルド方法
5.1 ソースプログラムとライブラリ
5.2 ビルド方法……………………………………………………2026…………………………………………………
…………………………………………………26……………………………………………………………………26－ ⅰ －Ⅰ. ATBD
1. 分類の原理（基本的な考え方）
分類プログラムの基本的な考え方を述べる．分類アルゴリズム（システム）は，どのような特徴に基づきどのように分類するか，すなわち使用する特徴量と識別器で特徴付けられる．そこで，特徴
量と識別器の２つの観点から，本研究で検討している分類アルゴリズムの骨子を以下に説明す
る．また，その理解がし易いように，想定する分類対象データと分類結果について簡単に述べる．
1.1 分類対象データ
極域を除く全域を対象範囲とする分光地表反射率の時系列データを分類対象データと想定す
る．極域は太陽が昇らない時期に対する地表反射率が得られないためである．対象観測時期の範囲は，過去１年間である．観測波長帯は，可視～中間赤外の複数バンドを想定する．
現在のところ次に示す SR と NBAR とよぶ MODIS の反射率プロダクトを使用してアルゴリズム
開発をしている．
・ MYD09A1.005 : Surface Reflectance 8-day L3（ SR ）・ MCD43A4.005 : Nadir BRDF Adjusted Reflectance 16-day L3（ NBAR ）

これらは，ともに 8 日単位に 1 年間あたり 46 回生産されており，表 1.1 に示すように可視～中間赤外まで 7 バンドから構成されている．座標系は SIN 投影座標系であり，公称 500m の地上分解能とされている．
表 1.1 MODIS 陸域バンドの仕様
バンド

波長域 [μm]

分解能 [m]

観測幅 [km]10.620 ～ 0.670

250

2330

2

0.841 ～ 0.876

30.459 ～ 0.479

4

0.545 ～ 0.565

51.230 ～ 1.250

6

1.628 ～ 1.652

7

2.105 ～ 2.155500

提案する分類アルゴリズムが上記の 2 種類の MODIS プロダクトの入力を要求しているわけではない．研究上，プロダクトの違いが分類結果に与える影響を把握するために用意しているだけであり，分類を行う際は，どちらか一方のプロダクトがあればよい．

-1-

1.2 分類に使用する特徴量(1)特徴量に関する従来の問題点
多時期データを対象としたときの大きな問題点の１つは，時間軸方向のデータ欠損と画素値の真値からの変動である．データ欠損は，衛星の軌道維持などの衛星・センサの維持管理のために
非観測時期が必ず発生することや，地上の受信・処理システムの何らかの問題により発生する．
また，雲・雪によって地表が観測されていない場合も，データ欠損と同じ扱いとなる．画素値の変動は，地表物体の BRDF の不均一性や起伏，また大気補正誤差などに起因し，それらは分類シ
ステム側から見ればノイズによる錯乱と見なせるのでノイズとよぶことにする．
従来，データ欠損については補間処理，またノイズに対しては平滑化処理などに多大な努力が払われている．補間処理がなされる理由の１つは，データ欠損があるような特徴量を識別器が受け付けないからである．一方，平滑化処理が行われる理由の１つは，分類に使用する特徴量の安定化を図り，識別器の性能の低下を抑制するためである．いずれも時間軸方向の連続性など
何らかの仮定を置いたアドホックな推定や近似を行うことになり，かつ対象データのデータ欠損やノイズの状況に応じた調整を行う必要がある．たとえば， MODIS データに対しては良好な補間や平滑化がなされたとしても SGLI データに対しては保証されないし，やはり時間と労力をかけて再調整することになろう．できれば，補間や平滑化を要請しない分類システムであることが望ましい．
(2)提案する特徴量
上記の入力データに関する問題を根本的に解決するには，データ欠損があっても識別器が動作し，ノイズがあっても識別性能に大きな影響を与えないような特徴量を使用すればよい．本研究では，時領域同時生起行列（ time domain co-occurrence matrix ）とよぶ特徴量を提案している．
従来の空間領域で定義される同時生起行列は，画像空間上で距離δ d 離れた 2 点の画素がそ
れぞれ画素値 i と j を持つ確率 P(i,j)で定義される．画素値のダイナミックレンジが L のとき L 行
L 列の行列となる．時領域同時生起行列は，これを時間δ t 離れた空間的には同一箇所の 2 つの画素に対して定義したものである．たとえば， Fig. 1(a)に示すような季節変化をある２つの画素（土地被覆）が示していた場合，それぞれ同図(b)に示す行列の要素の出現確率が 0 以外の値を持つ．すなわち，季節変化パターンの違いが行列に表れる．

pixel
value
i

piexel
value
j

pixel
value
i

d

t

piexelvaluej(a)空間領域で定義(b)時間領域で定義
図 1.1

同時生起行列の定義-2-

pixel value

pixel valueij(c)同時生起行列

occurrence
probability

pixel value
JAN
FEB
MAR
APR
MAY
JUN
JUL
AUG
SEP
OCT
NOV
DEC

(a)土地被覆の時系列変化
図 1.2(b)時領域同時生起行列

時領域同時生起行列の概念例

時領域同時生起行列では，ある時期の観測値にデータ欠損やノイズがあっても，それが出現確率に与える影響は小さく，ましてや識別器に入力できないなどということはない．さらに，たとえば農地は緯度方向で季節変化の違いによる作付け時期が異なるが，時期（位相）だけの違いで観測値の変化パターンに違いがあまりなければ，ほぼ同一の行列となる．すなわち，時領域同時生起行列は，データ欠損やノイズ，さらに位相変動に対して頑健な特徴量といえる．

-3-

1.3 識別法（識別器）(1)識別器に関する従来の問題点
従来たとえばツリー状分類法や最大尤度法など，最小値・最大値といった何らかのしきい値や平均・分散といった母集団パラメータを使用した分類法が多用されている．これらのパラメトリック
な分類法は，パラメータ（母数）を制御することできめ細かな制御が可能で処理も高速であるという利点を持つが，そのパラメータにしたがう分布を母集団がしているという前提を伴う．つまり処理対
象データに強く依存する．
このため，分類クラスの母集団がパラメータで示されるような分布をしているという前提を満足するように，分類クラスの設定やパラメータの設定に気を配らざるを得ず，分類対象データが変わるとそれに応じた調整を再度行うことになり，かなりの職人芸的な技術・経験そして時間が要求される．特に，本 GCOM-C1 のように，研究中には実際の分類対象データが存在せず，多時期のデータでスペクトル次元数が多く，処理対象データの種類も多い場合，あるいはまた異なるカテゴリ
ーセットごとに複数の土地被覆分類を行おうとすると，パラメータの種類や数が多くなるので，この問題はより大きなものとなる．
(2)提案する識別器
上記の識別器に関する問題を低減するため，本研究ではノンパラメトリック分類を識別器として使用することを提案している．ノンパラメトリック分類法として，現在のところノンパラメトリック最短距離法を使用している．
これは時領域同時生起行列の各要素を 1 次元のベクトルとして並べ直して得られる時領域同時生起行列の要素数だけの次元からなる特徴空間において，未知標本を最短距離にある分類クラスに分類する方法である．訓練エリアから計測された分類クラスの時領域同時生起行列（生起確率を表す）そのものを母集団分布（確率密度関数）であるとしてそのまま使用する．母集団の何らかの母数（パラメータ）を推定する必要はない．
ノンパラメトリック最短距離法を識別器として使用する場合，どのような距離を使用すべきかが問題となる．現在のところ，余弦距離（正規化相互相関ともよばれる）の使用を提案する．標本 xのチャネル b の時領域同時生起行列を Mx,b ，分類クラス c のチャネル b の時領域同時生起行列
を Mc,b としたとき， x と c との余弦距離 dn(x,c)は，次のように定義される．
d n  x, c   M i, j M i, j 
c ,b

x ,b

b

j

 M
b

j

i

i

2

  M 2 c,b i, j 

x ,b i, j

b

j

i

-4-

(1.1)

1.4 分類結果
多時期の分光地表反射率データから，分類結果として現在のところシングルプレーンの分類結果画像が１つ出力される．各画素は，決定された分類クラス ID 番号１つを画素値として持つ．カ
テゴリー ID 番号ではない．

-5-

2. ファジィ化
対象データである地表反射率時系列データには，様々な曖昧さが存在する．本研究では空間的な曖昧いさと時間的な曖昧さに注目して，この問題に対処している．
2.1 空間的ファジィ化
分類対象画像中のある特定の地表に対応する１画素の地表投影範囲は，時系列方向の複数
画像において一定ではない．そこで， 1 で述べた時領域同時生起行列では，時期差Δ t の画素対を考えるが，相手側の画素として空間的に拡げた複数画素を考える．このとき，出現頻度についてそれらしさを表す各画素に対する重みをメンバーシップ関数で表す．考える複数画素に対する重みの総和は 1.0 になるように設定する．

t

図 2.1

空間的ファイジィ化-6-

2.2 時間的ファジィ化
SR と NBAR はともに 8 日単位に生産されているが，それぞれ 8 日間と 16 日間の観測データから作成されている．したがって，ある１つの画素は，その 8 日間なり 16 日間のいずれかの日に観測されたものである．このため，時系列方向の地表の同一地点の画素値は， 8 日間隔で観測されたものではなく，ばらついている．
そこで，時領域同時生起行列を計測する時期差Δ t の画素対の相手側の画素として，Δ t の範囲を拡げて同一地点ではあるが生産日時が異なる複数の画素を考える．このとき，出現頻度についてそれらしさを表す各画素に対する重みをメンバーシップ関数で表す．考える複数画素に対する重みの総和は 1.0 になるように設定する．

t

図 2.2

時間的ファジィ化-7-3. 分類対象画素フィルタリング
分類対象画素とするか否か，あるいはそのバンドの画素値を使用するか否かについて，フィル
タリングを行う．
3.1 陸域・水域フィルタリング
土地被覆分類の対象は，陸域である．一方，地球表面上の 7 割は海とされているが，この海のように水域として時間的にあきらかに不変な領域の画素を分類対象画素から除く．これは，分類
処理時間を短縮するためである．
現在のところ， SR と NBAR それぞれのプロダクトごとに，異なったある方法で陸フラグを分類
対象画素値として設定している．すなわち，分類時にはその陸フラグは既に分類対象データに埋め込まれており，分類時はその陸フラグを見るだけであるので，分類処理部分のみ記述している
本ドキュメントでは陸フラグを設定するアルゴリズムについては説明しないことにする．ただし，分類時にはその陸フラグを見るので，その存在は知っておくべきであり，ここで言及した．

3.2 雲・雪・氷フィルタリング
分類対象プロダクトに付随している品質情報に基づき，地表反射率データを使用するか否かのフィルタリングを分類時に行う．これは， SR と NBAR のプロダクトそれぞれに対して異なるアルゴリズムを採用している．ここでは， NBAR プロダクトについてのみ以下に述べる．（ SR プロダクトについては，分類テストの実施に必要なデータ準備が今回間に合っておらず，説明しても意味をなさないためである．
(1)SR プロダクトに対する品質フィルタリング
上記の理由で，ここでは説明しない．
(2)NBAR プロダクトに対する品質フィルタリング
NBAR プ ロ ダ ク ト の 品 質 情 報 は ， MCD43A2 と い う 別 の プ ロ ダ ク ト で 生 産 さ れ て い る ．
MCD43A2 プロダクトは， 4 つのレイヤーで構成されているが，その内，・ Snow_BRDF_Albedo
・ BRDF_Albedo_Band_Qualityの２つのレイヤーの情報を利用する．
アルゴリズム上は，(1)フィルタリングしない，(2)ソフトフィルタリング，(3)ハードフィルタリング，の3 種類が選択できるようにしてある．その流れを図 3.1 に示す．

-8-

STARTno

Filteringyes
yes

Snow_BRDF_Albedo=1 (snow)nono
Hard Filteringyes
BRDF_Albedo_Band_Quality"best quality"であるband数Nを算出

noN<しきい値Tb
yes

分類対象画素から除く

図 3.1

分類対象画素とする

NBAR プロダクトに対する品質フィルタリング-9-4. システム構成と入出力データ
分類プログラムとその入出力データの論理構造を図 4.1 に示す．

処理パラメータ

前回使用された処理パラメータ

各種メッセージ

分類対象データ
ファイルのリスト

分類プログラム

使用された処理パラメータ

分類結果画像（分類クラスID
番号）

分類対象データ（多時期の分光地表
反射率データ）

訓練データ（分類クラスごとの時領域
同時生起行列）図 4.1

分類プログラムとその入出力データの論理構造- 10 -

Ⅱ. プログラム・テストデータ
1. フォルダ構造と記憶データ
1.1 概要HD のフォルダ構造を図 1.1 に示す．各フォルダに記憶してある内容の概要は次のとおりである．

ボリューム“提出用”
1.プログラム実行フォルダ
2.分類対象データ(1)“ 1.プログラム実行フォルダ”

3.リストファイル

分類プログラムを実行するときにカレントディレクトリ

4.トレーニングデータにするフォルダであり，実行形式のプログラムが格納さ

5.カラーテーブルれている．

6.分類結果画像
7.ソースプログラム(2)“ 2.分類対象データ”
分類プログラムをテスト実行する際に必要となる分

図 1.1

フォルダ構造概要

類対象データセットが格納されている．
(3)“ 3.リストファイル”
分類プログラムをテスト実行する際に，入力する処理パラメータの１つである分類対象データフ
ァイルのリストが格納されている．
(4)“ 4.トレーニングデータ”
分類プログラムをテスト実行する際に，入力する訓練データが格納されている．
(5)“ 5.カラーテーブル”
分類結果をカラーコーディングするためのカラーテーブルが格納されている．
(6)“ 6.分類結果画像”
分類結果例とカラーコーディング済みの分類結果画像例が格納されている．
(7)“ 7.ソースプログラム”
分類プログラムのソースプログラムおよびビルド時に必要となるライブラリが格納されている．

- 11 -

1.2 フォルダ“ 1.プログラム実行フォルダ”の内容
分類プログラムを実行するときにカレントディレクトリにすることを想定しているフォルダである．
内容は，表 1.1 に示すとおりである．
表 1.1

フォルダ“ 1.プログラム実行フォルダ”に格納されているファイル一覧

ファイル名

説

明

comCLimg2W_x64.exe

64 bit Windows 環境下で動く実行形式の分類プログラム

comCLimg2W_x32.exe

32 bit Windows 環境下で動く実行形式の分類プログラム

comCLimg2W.prm

前回対話入力された処理パラメータが格納されているテキストファ
イル（「処理パラメータファイル」とよんでいる）．
フォルダ“ 6.分類結果画像”に格納されている実行結果例を作成したとき使用したときの処理パラメータが格納されている．

実行画面例.txt

実行画面例（フォルダ“ 6.分類結果画像”に格納されている実行結果例を作成したとき使用したときのもの）[注意事項](a)分類プログラムは，コマンドプロンプトから実行するようにする．
(b)分類プログラムをコマンドプロンプトから実行するとき，フォルダ“ 1.プログラム実行フォルダ”をカレントディレクトリとする．その理由は，前回使用された処理パラメータがファイル
comCLimg2W.prm に記録されており，それを読み出す場所がカレントディレクトリを前提としているためである．もし， comCLimg2W.prm がカレントディレクトに存在しない場合，処理パラメータとして分類プログラム自身が持っている初期値が使用され，初心者にはどのような処理パラメータを設定してよいのか困ることになるからである．
(c)処理パラメータファイル comCLimg2W.prm は，分類プログラムを実行するたびに，内容が書き換えられる．したがって，現在，記録されているファイルの内容を書き換えられないように，現在のフ ァ イ ル comCLimg2W.prm を 一 旦 デ ス ク ト ッ プ か ど こ か に コ ピ ー し ， フ ァ イ ル 名 を
comCLimg2W.prm-原.txt などに変更する．そのうえで，このフォルダ“ 1.プログラム実行フォルダ”に改名したファイルを移動しておくようにすることを勧める．

- 12 -

1.3 フォルダ“ 2.分類対象データ”の内容
分類プログラムをテスト実行する際に必要となる分類対象データセットが格納されている．分類
対象データセットを構成する MODIS プロダクトは， MCD43A4, MCD43A2, MYD09A1 の 3 種類
であり，それぞれ"MCD43A4.005", "MCD43A2.005", "MYD09A1.005"というサブフォルダに格納されている．カクフォルダの末尾に".005"が付いているのは， MODIS プロジェクト側で，それらのプロダクトを生産したときのアルゴリズムが ver.5 だったためであり，当方でかってに付けたものではない．
プロダクト MCD43A4 と MYD09A1 は，それぞれ分光地表反射率プロダクトとナディア方向の分光地表反射率プロダクトである．また， MCD43A2 は MCD43A4 の品質情報プロダクトである．
MCD43A4 と MCD43A2 は，衛星 Terra と Aqua に搭載されている MODIS ２つの観測データから生産されているが， MYD09A1 は衛星 Aqua に搭載されている MODIS だけの観測データから生産されている．
各プロダクトは， 2007 年の１年間に 8 日間隔に生産されたもが本フォルダに格納されている．
各プロダクトのサブフォルダに， 2007.01.01, 2007.01.09, 2007.01.17,・・・, 2007.12.27 とあるのは，そのためである．各日付のサブフォルダの中にデータが格納されている．それを表 1.2 ，表 1.3 ，
表 1.4 に示す．

表 1.2

MCD43A2 プロダクトの各日付のサブフォルダに格納されているファイル一覧

ファイル名

説

明

M2007????-BRDF_Albedo_Band_Quality.raw

MCD43A2 プロダクトを構成する４M2007????-BRDF_Albedo_Band_Quality-Smp10_10.raw

つのレイヤーに対応するバイナリ

M2007????-Snow_BRDF_Albedo.raw

形式ファイル．

M2007????-Snow_BRDF_Albedo-Smp10_10.raw（上記のファイル名と同じ）.hdr

上記の４つのレイヤーに対応するファイルそれぞれのヘッ
ダーファイル．テキスト形式ファイル．ファイルサイズなどが記録されている．

注 1)文字“????”は， MMDD 形式でサブフォルダ名と同じ日付が入る．

- 13 -

表 1.3

MCD43A4 プロダクトの各日付のサブフォルダに格納されているファイル一覧

M2007????-C20070101-Nadir_Reflectance_Band1432567-land2-t3.03.06.00.00.00.00.5.raw
説明： MCD43A4 を構成する 7 レイヤーすなわち７つの波長バンドの地表反射率を 8 bit に非線形圧縮して BIP 形式にした画素間隔 0.005 度の経緯度座標系のバイナリ型ファイル．
M2007????-Nadir_Reflectance_Band1432567-land2-Smp10_10.raw
説明： 上記の画素間隔 0.005 度の画像データを横・縦方向それぞれ 1/10 サンプリングした画
素間隔 0.05 度のバイナリ型ファイル．
M2007????-Nadir_Reflectance_Band1432567-land2-Smp10_10-clst2000WD_BF5000.raw
説明： 上記の画素間隔 0.005 度の画像をスペクトルクラスタリングして得られたクラスタ画像．
バイナリ型ファイルであるが，今回の実行テストでは使用するようにドキュメントに記述していな
いので，このファイルの存在を無視してもらいたい．
（上記のファイル名と同じ）.hdr
説明： 上記の３つのファイルそれぞれのヘッダーファイル．テキスト形式ファイルである．ファイ
ルサイズなどが記録されている．
注 1)文字“????”は， MMDD 形式でサブフォルダ名と同じ日付が入る．

表 1.4

MYD09A1 プロダクトの各日付のサブフォルダに格納されているファイル一覧

C2007????-sur_refl_b1432567-land2-t3.03.06.00.00.00.00.0.raw
説明： MYD09A1 に含まれている７つの波長バンドの地表反射率を 8 bit に非線形圧縮してBIP 形式にした画素間隔 0.005 度の経緯度座標系のバイナリ型ファイル．
C2007????-sur_refl_b1432567-land2-Smp10_10.raw
説明： 上記の画素間隔 0.005 度の画像データを横・縦方向それぞれ 1/10 サンプリングした画
素間隔 0.05 度のバイナリ型ファイル．
C2007????-sur_refl_b1432567-land2-Smp10_10-clst2000WD.raw
説明： 上記の画素間隔 0.005 度の画像をスペクトルクラスタリングして得られたクラスタ画像．
バイナリ型ファイルであるが，今回の実行テストでは使用するようにドキュメントに記述していな
いので，このファイルの存在を無視してもらいたい．
M2007????-sur_refl_state_500m.raw
説明： MYD09A1 に含まれている品質情報の１つで，画素間隔 0.005 度の経緯度座標系のバ
イナリ型ファイル．
M2007????-sur_refl_state_500m-Smp10_10.raw
説明：

上記の品質情報ファイルを 1/10 サンプリングして画素間隔 0.05 度の経緯度座標系のバイナリ型ファイル．
（上記のファイル名と同じ）.hdr
説明： 上記の３つのファイルそれぞれのヘッダーファイル．テキスト形式ファイルである．ファイ
ルサイズなどが記録されている．
注 1)文字“????”は， MMDD 形式でサブフォルダ名と同じ日付が入る．

- 14 -

[注意](a)MCD ・・・と MYD ・・・の 2 種類の分光地表反射率データが本 HD に格納されているが，分類
ブログラムがその 2 種類の入力を必要としているわけではない．この正確の異なる地表反射率データをそれぞれ分類したときの分類性能を調べたいという当方の要望から用意してあるだけで，今回の負荷テストをするという意味では，どちらか一方を使用するだけで事は足りると思われる．
今回の実行テストでは， MCD ・・・を使用するという前提で，その他のデータや実行結果なども収録してある． MYD09A1 のデータがフォルダ“ 2.分類対象データ”に格納されているのは，もしJAXA のプロジェクト側で， MYD09A1 についても実行テストをしたいと思ったときに，困らないよ
うに念のために格納してあるだけである．
(b)スペクトルクラスタ画像が格納してあるが，今回の実行テストでは，実行するようにドキュメントに記述していないので，その存在を無視して頂きたい．

- 15 -

1.4 フォルダ“ 3.リストファイル”の内容
分類プログラムをテスト実行する際に，入力する処理パラメータの１つである分類対象データフ
ァイルのリストが格納されている．内容は，表 1.5 に示すとおりである．

表 1.5

フォルダ“ 3.リストファイル”に格納されているファイル一覧

List-MC_Band1432567-Smp10_10_CL 用.txt
説明： フォルダ“ 2.分類対象データ”にある 1/10 サンプリングされた画素間隔 0.05 度の BIP
形式 MCD43A4 地表反射率データファイルのリスト．テキスト形式ファイルである．
List-MC_Snow_BRDF_Albedo-Smp10_10_CL 用.txt
説明： フォルダ“ 2.分類対象データ”にある 1/10 サンプリングされた画素間隔 0.05 度の
MCD43A2 品質情報ファイルのリスト．テキスト形式ファイルである．
[注意](a)表 1.5 からわかるように， 1/10 サンプリングされた MCD ・・・・ファイルに対してのみのリストファ
イルが用意されている．つまり，それ以外のサンプリングされていない MCD ・・・ファイルや MYD・・・ファイルを対象に分類を行うには，それらのファイルのリストファイルをエディタ等で作成する必要がある．

1.5 フォルダ“ 4.トレーニングデータ”の内容
分類プログラムをテスト実行する際に，入力する訓練データが格納されている．内容は，表 1.5に示すとおりである．

表 1.6

フォルダ“ 4.トレーニングデータ”に格納されているファイル一覧

LCT1-MCimg-yama20110609-p15-hardF-test.TRNG
説明： MCD43A4 を分類対象データとした分類を行う際のバイナリ型訓練データファイルである． MCD43A4 であれば，画素間隔にかかわらずいずれのデータファイルに対しても使用できる．
（上記のファイル名と同じ）.hdr
説明： 上記のファイルのヘッダーファイル．テキスト形式ファイルである．メタデータが記録されている．

- 16 -

1.6 フォルダ“ 5.カラーテーブル”の内容
分類結果をカラーコーディングするためのカラーテーブルが格納されている．内容は，表 1.7に示すとおりである．
表 1.7

フォルダ“ 5.カラーテーブル”に格納されているファイル一覧

ファイル名
分類結果画像用.act

説明Adobe 社の Photoshop 形式のカラーテーブルでる．分類結果をイン
デックスカラー形式として，このカラーテーブルを読み込むと，カテゴ
リー別にカラーコーディングされた画像を得ることができる．

1.7 フォルダ“ 6.分類結果画像”の内容
分類結果例とカラーコーディング済みの分類結果画像例が格納されている．
表 1.8

フォルダ“ 6.分類結果画像”に格納されているファイル一覧

ファイル名
分類結果画像例.raw

説明フォルダ“ 1.プログラム実行フォルダ”に格納されている実行画面
例で分類プログラムを実行したときに得られた分類結果ファイル．
シングルチャネルのバイナリ型ファイルである．画素値にはクラス番号が記録されている．

分類結果画像例.raw.hdr

上記ファイルのヘッダーファイル．メタデータが記録されている．テ
キスト形式ファイルである．

分類結果画像例.psd

上記の分類結果ファイルをフォルダ“ 5.カラーテーブル”に格納されているカラーテーブルでカテゴライズした画像ファイルである．
Photoshop 形式のファイルである．

1.8 フォルダ“ 7.ソースプログラム”の内容
分類プログラムのソースプログラムおよびビルド時に必要となるライブラリが格納されている．プ
ログラムの開発には，プログラム統合環境 Microsoft Visual Studio 2008 Standard ，コンパイラ
Intel Visual Fortran 11.1 を使用している．分類プログラムのソースプログラムは，その開発環境ごとサブフォルダ“プログラム開発環境ごと”の中にサブフォルダ“ comCLimg2W ”として格納されている．また，使用しているライブラリは，フォルダ“ 7.ソースプログラム”の下のサブフォルダ“ LIB ”に格納されている．サブフォルダ“ LIB ”には， 32 bit 版と 64 bit 版のプログラムビルド向けのライ
ブラリが用意されている．関係するファイル一覧は，表 5.1 を参照されたい．

- 17 -

2. 分類処理の入出力データ
本 HD には，分類を行うプログラム comCLimg2W （ファイル名とは異なる）１つが格納されており，それに必要なデータや実行例が格納されている．そのプログラムと入出力データについて述
べる．
入力データは，対話（コンソール：標準出力装置から）入力される文字列と各種のデータファイ
ルに大別される．出力データも，コンソール（標準出力装置）に出力される文字列と各種のデータ
ファイルに大別される（．コンソール入出力以外の入力データファイルと出力データファイルの概要をそれぞれ 2.1 と 2.2 に示す．

処理パラメータ

処理パラメータ
ファイル

各種メッセージ

分類対象データ
ファイルのリスト
ファイル

分類プログラム

分類対象データ
ファイルとそのヘッダー
ファイル

処理パラメータ
ファイル

分類結果
ファイルと
そのヘッダー
ファイル

訓練データ
ファイルと
そのヘッダーファイル

図 2.1

分類プログラムの入出力データ

2.1 入力データファイル(1)処理パラメータファイル
前回分類プログラムを実行したとき使用された対話入力されたパラメータが格納されている．次に分類プログラムを実行するとき，これらの値が default 値として使用されるので，次に分類プログ
ラムを実行する際に， default 値そのままでよい場合は， Enter キーを押すだけで済む．
(2)処理対象データファイルのリストファイル
分類対象データファイルのリストである．分類対象データは，多時期の地表反射率データファイ- 18 -

ルであり， 8 日間隔で生産されて 1 年間 46 シーン存在する．分類プログラムは，とりあえず１年間のデータを対象に分類を行うことが想定されている．これを１つ１つ対話入力するのは大変なので，あらかじめそのリストを作成しておき，そのリストファイルを指示すればよいようにしてある．
(3)処理対象データファイル
MCD ・・・と MYD ・・・の 2 種類のプロダクト，それぞれ画素間隔 0.005 度と 0.05 度間隔の経緯度座標系の分類対象（分光地表反射率）データファイルが用意されている．
(4)処理対象データのヘッダーファイル
分類対象データファイルは地表反射率データのみが格納されたバイナリ型ファイルである．その分類対象データファイルを読み出すにあたっって必要な画像サイズなどのメタデータが格納されている．
(5)訓練データファイル
分類カテゴリーごとに複数のサブカテゴリーが設けられており，それを分類クラスとよぶが，その分類クラスごとの訓練（トレーニング）データすなわち時領域同時生起行列が格納されている．分類プログラムは，この訓練データに最も類似したクラスに各分類対象画素を分類する．
(6)訓練データファイルのヘッダーファイル
訓練データはバイナリ型ファイルであり，それを読み出すにあたって必要なメタデータが記録されている．テキスト形式のファイルである．

2.2 出力データファイル(1)処理パラメータファイル
対話（コンソール：標準入力装置）入力された処理パラメータが格納される．このデータは，次回，処理プログラムを実行する際に，対話入力の default 値として使用される．
(2)分類結果ファイル
分類対象画像の各画素は，いずれかの分類クラスに分類されるが，その分類ラス ID 番号が各画素に格納されている．シングルレイヤーの 8bit 画像ファイルである．
(3)分類結果ファイルのヘッダーファイル
分類結果ファイルの画像サイズなどのメタデータが格納されている．テキスト形式ファイルである．

- 19 -

3. 分類プログラムの実行方法
フォルダ“ 1.プログラム実行フォルダ”には，既にビルドされた実行形式のプログラムが格納されている． Win32 （ 32bit Windows ）版と x64 （ 64bit Windows ）版の内，どちらかを使用する． x64版の方が実行速度は速いので， 64 bit Windows 環境下で実行することを勧める．
comCLimg2W_x64.exe

：

x64 版

comCLimg2W_x32.exe：

Win32 版

実行は，以下のように行う．
(1)コマンドプロンプトの起動
まず最初に，“コマンドプロンプト”を起動して，カレントディレクトをフォルダ“ 1.プログラム実行
フォルダ”に移動する．プロンプトに対して分類プログラムファイル指定して，分類プログラムを起動する．（下記の実行画面例を参照）[注意]
上記のようなことをせずに，分類プログラムをエクスプローラ上でダブルクリックしても分類プログ
ラムは起動されるが，そうせずに上記のようにして起動することを勧める．これは，何か実行時にエラーが発生してプログラムが動作を停止したとき，もしダブルクリックして起動した場合，コマンドプ
ロンプトの画面が自動的にクローズされてしまい，何が発生したのか確認できなくなることを避けるためである．
(2)処理パラメータの指定
実行画面例（これはフォルダ“ 1.プログラム実行フォルダ”にファイル“実行画面例.txt ”として格納されている．）を以下に示す．この画面と等しくなるように処理パラメータを指定する．なお実行
画面例では，プログラム名“ comCLimg2W.exe ”を起動するようにしてあるが，実際には前述のどちらかのプログラムを指定する．また，ドライブ名はシステムごとに変わるのでそれを指定すること．
なお，このコマンドプロンプト画面でのパラメータの入力規則は，次のとおりである．
(a)“－”で始まるメッセージは，単なるメッセージであるが，“？”で始まるメッセージはプログラムがオペレータに対して何か入力を促しているので，オペレータは何か入力しなければならない．
(b)入力を促すとき，次のように前回入力されたパラメータが default 値として，コロン":"に挟まった形で表示される． dafault 値でよい場合は，単に Enter キーを押すだけでよい． default 値とは違うものを指定したいときのみキーインする．
[形式]? 質問内容： default 値：(c)入力待ち状態になっているとき，“^”（ハット），“^S ”，“^E ”をキーインすると，１つ前の質問に戻ったり，最初の質問にもどったり，さらにプログラムを直ちに終了することができる（表 3.1 参照）．
(d)ファイルやフォルダを指定するときは，エクスプローラ上からコマンドプロンプト画面にそのファイルやフォルダをドラッグ&ドロップすればよい（ Windows 7 や Windows XP の場合）．

- 20 -

表 3.1

質問に対する特別な入力

入力する文字 プログラムの動作^

ハットを入力すると，１つ前の質問に戻る^S

最初の質問に戻る^E

プログラムを直ちに終了する

Microsoft Windows [Version 6.1.7601]
Copyright (c) 2009 Microsoft Corporation. All rights reserved.
C:¥Users¥admin>q:Q:¥>cd Q:¥1.プ ロ グ ラ ム 実 行 フ ォ ル ダ
Q:¥1.プ ロ グ ラ ム 実 行 フ ォ ル ダ >Q:¥1.プ ロ グ ラ ム 実 行 フ ォ ル ダ ¥comCLimg2W.exe[同 時 生 起 行 列 を 使 用 し た 多 時 期 分 類 ： multi-ch image 専 用 ]- 同時生起行列のバイナリ型訓練ファイル? フ ァ イ ル 名 : Q:¥4.ト レ ー ニ ン グ デ ー タ ¥LCT1-MCimg-yama20110609-p15-hardF-test.TRNG :#classes, #channels= 69
7 ....detected.
? #class, #ch : 69 7 :? ok this training data file <y/n> : y :class No. classID #frequency
1
0
3077
2
1
0 .... warning !
3
2
2030
4
3
19
5
4
1417
6
10
4931
7
11
2829
（途中省略）
65
170
1084
66
171
465
67
172
0 .... warning !
68
94
358
69
151
12214
? classID の 確 認 ・ 再 設 定
Y. 別 途 classID table file か ら 読 み 込 ん で か ら ， 確 認 ・ 再 設 定 す る
N. 直 ち に キ ー ボ ー ド よ り キ ー イ ン し て 再 設 定 す る : n :? classID(分 類 出 力 値 ) : 0 1 2 3 4 10 11 20 21 22 23 30 31 40
41 42 51 52 54 60 61 70 72 73 74 75 76 80 93 82 83 84 91 92
100 101 102 110 111 112 120 121 103 123 124 125 126 127 130 131 132 133 134 135
140 141 142 160 161 162 163 164 165 166 170 171 172 94 151 :? 分 類 対 象 ク ラ ス 数 : 46 :? 分 類 対 象 か ら 除 く ク ラ ス No. : 2 4 7 8 9 16 19 22 25 29 35 38 39 47 49 50 5
2 53 55 61 66 67 68 :? 訓 練 デ ー タ の う ち ， 分 類 に 使 用 す る ch 数 : 7 :? 使 用 す る ch# : 1 2 3 4 5 6 7 :- リストファイル：分類対象画像ファイル? フ ァ イ ル 名 : Q:¥3.リ ス ト フ ァ イ ル ¥List-MC_Band1432567-Smp10_10_CL 用 .txt :- リスト中の指定する文字列を他の文字列に変更できます.
? 変 更 し ま す か <y/n> : n :46 scenes are detected.
open check...................................................finished#px/ln, #ln/fl, #ch/px = 7200 2800
7
ch. format=CIP, data type=int8u（次頁へつづく）- 21 -

（前ページよりつづき）? 対 象 ch 数 : 7 :? 処 理 対 象 ch# : 1 2 3 4 5 6 7 :? 処理対象入力画像データの下限，上限 :2254 :- 対象画像範囲? 開 始 ,終 了 画 素 No., 開 始 ,終 了 ラ イ ン No. :1 72001 2800 :? 画素方向，ライン方向の読み出しステップサイズ :55:? SIC state file を 利 用 し て 雪 ・ 氷 ・ 雲 を 検 出 し ま す か <y/n> : n :? プロダクト付属の品質情報を利用して分類に不適当な画素をフィルタリングしますか
N. No （ フ ィ ル タ リ ン グ し な い ）H. Hard Filtering （ 訓 練 デ ー タ 抽 出 時 と 同 じ く 強 い 条 件 で フ ィ ル タ リ ン グ す る
S. Soft Filtering （ 条 件 を 弱 め て フ ィ ル タ リ ン グ す る ） : S :- リ ス ト フ ァ イ ル ： MC-Snow_BRDF_Albedo フ ァ イ ル? フ ァ イ ル 名 : Q:¥3.リ ス ト フ ァ イ ル ¥List-MC_Snow_BRDF_Albedo-Smp10_10_CL 用 .txt :- リスト中の指定する文字列を他の文字列に変更できます.
? 変 更 し ま す か <y/n> : n :46 scenes are detected.
open check...................................................finished#px/ln, #ln/fl, #ch/px = 7200 2800
1
ch. format=CIP, data type=int8u? 時領域－同時生起行列を求める際の相手方の数 : 3 :? 相 手 方 の 相 対 位 置 : 14 15 16 :? 相 手 方 に 対 す る 重 み （ 総 和 =1.0 ） : 0.250000 0.500000 0.250000 :? 空間領域－相手方の画素に隣接する画素数
0. 隣 接 す る 画 素 は 考 え な い
2. 左 右 の ２ 画 素 を 考 え る
4. 4 近 傍 画 素 を 考 え る
8. 8 近 傍 画 素 を 考 え る:8:? 重 み － 中 心 画 素 ， 上 下 /左 右 画 素 ， 斜 め 画 素 : 0.320 0.100 0.070 :? 分類時に使用する距離
1.市 街 地 距 離
2.ユ ー ク リ ッ ド 距 離
3.余 弦 距 離
4.PLUS 類 似 度: 3:
5.相 互 相 関
? 分 類 結 果 画 像 の 背 景 画 素 値 : 255 :- 出力分類結果ファイル? フ ァ イ ル 名 : Q:¥7.分 類 結 果 画 像 ¥分 類 結 果 画 像 .raw :? サ ン プ リ ン グ し て い ま す が ， 出 力 画 像 フ ァ イ ル は パ ッ ク し ま す か <y/n> : y :? 分類対象を訓練領域のみに制限することができます
Y. 制 限 す る
N. 制 限 し な い : n :? debug の た め に 同 時 生 起 行 列 を text フ ァ イ ル に 出 力 し ま す か <y/n> : n :- 処理パラメータの設定は終了しました．
? 処 理 を 開 始 し て い い で す か <y.yes n.rekey-in parameters> : y :処 理 中 .................
1......... 1001......... 2001.......
経 過 時 間 [時 :分 :秒 ]= 1:02:29.37
CPU 時 間 [分 :秒 ]= 54:44.36- 処理終了しました.
? E. プ ロ グ ラ ム を 終 了
R. 最 初 か ら も う １ 度 : E :- 22 -

以上の操作で，分類結果ファイルが，指定したフォルダに指定したファイル名で作成される．本HD 内のフォルダ“分類結果画像”に同じファイル名で分類結果を出力すると，既に作成されている分類結果ファイル例が重ね書きされるので注意されたい．同じフォルダに出力するのであれ
ば，ファイル名を既存のものとは違う名前にすること．
(3)分類対象データの選択に関する注意
フォルダ“ 2.分類対象データ”には， 2 種類のプロダクト MCD43A4 （品質情報の MCD43A2が付随する）と MYD09A1 ごとに， 0.005 度の画素間隔の分光地表反射率データが 2007 年の 1
年間 46 シーン格納されている．ところが，その画像サイズは 72,000 画素× 24,000 ラインと非常に大きい．実行テストの目的からはこの全画素を処理する必要はないと思われる．そこで，その画像を横・縦両方向それぞれ 1/10 サンプリングした 7,200 画素× 2,400 ラインのデータも同フォルダに用意してある．上記の実行画面例では，その 1/10 サンプリングしてあるデータ（言い換えると， 0.05 度の画素間隔データ）を対象に分類処理を行っている．
さらに，その分類対象データは，フォルダ“ 3.リストファイル”に格納されている分類対象ファイ
ルのファイルリストのファイルで指示するようになっているが，“ 3.リストファイル”に実際に格納されているのは，上記の実行画面例で使用した MCD43A4 （と MCD43A2 ）の 1/10 サインプリングされたファイルリストのファイルだけである． 0.005 度の画素間隔のデータを対象に分類を行いたい
場合や， MYD09A1 を対象にする場合は，同様のファイルリストを作成する必要がある．
(4)画像サイズと処理対象画像のサンプリング読み出しに関する注意
上記に示した実行画面例では，水平方向（画素方向）と垂直方向（ライン方向）それぞれのステップサイズとして， 5, 5 を指定してある．画面例で示してある分類対象画像の画像サイズは，7,200 画素× 2,800 ラインであるので（これは横縦方向にそれぞれ 1/10 サンプリングされて， 0.01
度画素間隔の画像である）， 1,440 画素（=7,200/5 ）× 560 ライン（=2,800/5 ）のサイズの画像が分
類結果画像として出力される．これは実行テストの処理時間を短縮するためである．実行画面例に示されているように，そのときのスループットは約 1 時間であった．
もし，同じサイズの分類対象画像ファイルに対して，ステップサイズを 10, 10 と指定すると，出力分類結果画像のサイズはさらに半分の 720 画素（=7,200/10 ）× 280 ライン（=2,800/10 ）となり，
スループットは約 15 分（=60 分/4 ）となるなるはずである．このように，実行テストをする場合，実行画面例にある 1/10 にサンプリングされた分類対象画像を対象に，大きなステップサイズで読み出すようにすれば，実行テストを早く終了させることができる．
逆にいえば， 1/10 サンプリングしていない分類対象画像ですら， 0.005 度画素間隔で作成されており， SGLI の場合の 0.0025 度画素間隔の場合と比べて画素数は 1/4 のものである．さらに，実際にはサンプリングして分類処理をすることはないので，ステップサイズのべき乗をさらに10 のべき乗倍し，さらに 4 倍した処理時間が 0.0025 度画素間隔の全画素を分類したときのスループットである．
たとえば，上述の実行画面例の場合から 0.0025 度画素間隔の全画素を処理した場合のスループットを推定すると， 1 時間× 5 × 5 × 10 × 10 × 2 × 2= 10,000 時間であると予測される．
これは約 417 日間に相当する．すなわち，全画素をシーケンシャルに分類すると１年以上かかる．しかし，分類はシーケンシャルにする必要はないので，実際の SGLI の処理では，分類対象- 23 -

域を 6 コアの CPU ２つ搭載の PC １台で 12 地域に分割して処理させると， 35 日で処理を終える
ことができる．この PC を 3 台用意すれば，さらに 1/3 の 12 日となり，約半月で処理を終えることができる．ただし，これは単純計算である．実際には，分類対象データのコピーを複数のディスクドラ
イブに配置しておき，データ読み出し速度を低下させないようにする必要がある．
なお，上記画面例に示したプログラム実行時にデータは SATA-2 （ 3Gbps 規格）の HDD をローカル接続したものを使用している．これを USB2 規格の I/F で接続にすると読み出し速度は約半分になる（なった記憶がある．正確ではない）． NAS のように 1Gbps の LAN 接続にするとさらにその数分の１以下に低下するので注意が必要である．

- 24 -

4. 分類結果の表示方法
分類結果ファイルは， 8bit シングルチャネル（プレーン）の raw 形式ファイルとして出力される．
各画素には，土地被覆のカテゴリー ID 番号ではなく，分類クラス ID 番号が格納されている．これをカテゴリーごとにカラーコーディングして表示するが，それには， Adobe Photoshop を使用する．
分類結果画像ファイルと同名で最後に“.hdr ”が追加されたファイルがフォルダ“ 6.分類結果画像”に格納されている．これは，ヘッダーファイルとよんでいるもので，分類結果画像ファイルの画
像サイズなどが書き込まれたテキスト形式のファイルである．分類プログラムによって自動的に作成されるものである．これをエディタ等で開き，画像サイズなどを確認しておく．その後，
Photoshop で分類結果画像ファイルを RAW 形式を指定して開く．
次に，画像モードを「インデックスカラー」にした後，フォルダ“カラーテーブル”に用意されているカラーテーブル（拡張子が act のファイル）を読み出せば，土地被覆カテゴリーごとにカラーコー
ディングされた画像が表示される．
必要に応じて，画像モードをカラー画像にして， TIF 形式や jpg 形式など適当な形式で保存すればよい．

- 25 -

5. 分類プログラムのビルド方法
5.1 ソースプログラムとライブラリ
当方では， Windows 7 64bit 版の OS の元で，プログラム統合環境 Microsoft Visual Studio2008 Standard ，コンパイラ Intel Visual Fortran 11.1 を使用している．これからわかるように，ソー
スコードは Fortran で書かれている．使用している機能は， Fortran 90 以上である． HDD にはこの開発 環境ごとソースプログラムが フォルダ“ 7.ソースプログラム¥プログラム開発環境ごと¥comCLimg2W ”に格納されている．また，このプログラム comCLimg2W のビルドに必要なライブ
ラリが，フォルダ“ 7.ソースプログラム¥LIB ”に格納されている．
表 5.1
フォルダ

分類プログラムのソースプログラムとライブラリ

ファイル

説

明

¥7.ソースプログラム¥プログラム開発環境ごと¥comCLimg2W¥comCLimg2W
comCLimg2W.vfproj

上記の環境下で，このファイルをクリックすれば，統合環境が起動する．

comCLimg2W.F90

main ソースプログラム

comCLimg2W-prc.F90

関数ソースプログラム

changeString.F90

関数ソースプログラム

if_ckNpxClass.F90

関数ソースプログラム

if_setPrcPrm.F90

関数ソースプログラム

readListTargetFL.F90

関数ソースプログラム¥7.ソースプログラム¥LIB¥Win32-Release 用PaFI_LIB2010.lib

Win32 （ 32 bit ）版の実行プログラムをビルドするときに使用する 32b it 版ライブラリ¥7.ソースプログラム¥LIB¥X64-Release 用PaFI_LIB2010.lib

x64 （ 64 bit ）版の実行プログラムをビルドするときに使用する 64 bit 版ライブラリ

5.2 ビルド方法
適切なライブラリを指定してビルドすること以外，特に説明不要と思われる．

- 26 -