(57)【特許請求の範囲】
【請求項１】
異なる第１時点及び第２時点で撮影された航空写真と、前記第１時点に対応する地図情報とを用いて土地利用変化を検出する変化領域検出装置であって、
変化領域検出を行う際に、入力された前記第１時点の航空写真の画素を特徴づける局所特徴ベクトル、及び前記第２時点の航空写真の画素を特徴づける局所特徴ベクトルを計算する局所特徴ベクトル計算手段と、
前記第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、前記画素ごとに変化量を計算する変化量計算手段と、
前記画素ごとの変化量から、変化領域候補を検出する変化領域候補決定手段と、
前記第１時点に対応する地図情報と前記第１時点の航空写真とから、前記第１時点の航空写真の小領域と、変化検出対象を含むか否かを示す情報とからなる検出対象識別器学習用データセットを作成する学習用データセット構築手段と、
前記検出対象識別器学習用データセットから、変化検出対象を含むか否かを識別するための検出対象識別器を学習する検出対象識別器構築手段と、
前記変化領域候補決定手段によって検出された変化領域候補の各々について、前記検出対象識別器を用いて、前記変化領域候補が前記変化検出対象を含むか否かを識別し、前記変化検出対象を含まない前記変化領域候補を除外する非検出対象領域除外手段と、
を含む変化領域検出装置。
【請求項２】
前記学習用データセット構築手段は、前記地図情報から、前記変化検出対象の物体の存在する画像座標を取得し、当該画像座標に基づいて前記検出対象識別器学習用データセットを作成する請求項１に記載の変化領域検出装置。
【請求項３】
前記局所特徴ベクトル計算手段は、注目画素及び前記注目画素の周辺の小領域の色又は輝度の情報も利用して、注目画素周辺の形又はテクスチャに関する情報を記述する前記局所特徴ベクトルを計算する請求項１又は２に記載の変化領域検出装置。
【請求項４】
前記検出対象識別器構築手段は、前記検出対象識別器学習用データセットに含まれる小領域の画素を特徴づける局所特徴ベクトルを計算し、前記計算された画素の前記局所特徴ベクトルをまとめ上げる処理によって、前記検出対象識別器学習用データセットに含まれる小領域の大域特徴ベクトルを計算し、前記大域特徴ベクトルに基づいて変化検出対象を含むか否かを識別するための前記検出対象識別器を学習し、
前記非検出対象領域除外手段は、前記変化領域候補決定手段によって検出された変化領域候補の各々について、前記第２時点の航空領域の前記変化領域候補の画素を特徴づける局所特徴ベクトルを計算し、前記計算された画素の前記局所特徴ベクトルをまとめ上げる処理によって、前記変化領域候補の大域特徴ベクトルを計算し、前記計算した大域特徴ベクトルに基づいて、前記検出対象識別器を用いて、前記変化領域候補が前記変化検出対象を含むか否かを識別する請求項１〜請求項３の何れか１項記載の変化領域検出装置。
【請求項５】
前記変化量計算手段は、前記第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、前記画素ごとに、利用する局所特徴ベクトルに合わせた距離関数から計算される距離を、前記変化量として計算する請求項１〜請求項４の何れか１項に記載の変化領域検出装置。
【請求項６】
局所特徴ベクトル計算手段、変化量計算手段、変化領域候補決定手段、学習用データセット構築手段、検出対象識別器構築手段、及び非検出対象領域除外手段を含み、異なる第１時点及び第２時点で撮影された航空写真と、前記第１時点に対応する地図情報とを用いて土地利用変化を検出する変化領域検出装置における変化領域検出方法であって、
前記局所特徴ベクトル計算手段が、変化領域検出を行う際に、入力された前記第１時点の航空写真の画素を特徴づける局所特徴ベクトル、及び前記第２時点の航空写真の画素を特徴づける局所特徴ベクトルを計算するステップと、
前記変化量計算手段が、前記第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、前記画素ごとに変化量を計算するステップと、
前記変化領域候補決定手段が、前記画素ごとの変化量から、変化領域候補を検出するステップと、
前記学習用データセット構築手段が、前記第１時点に対応する地図情報と前記第１時点の航空写真とから、前記第１時点の航空写真の小領域と、変化検出対象を含むか否かを示す情報とからなる検出対象識別器学習用データセットを作成するステップと、
前記検出対象識別器構築手段が、前記検出対象識別器学習用データセットから、変化検出対象を含むか否かを識別するための検出対象識別器を学習するステップと、
前記非検出対象領域除外手段が、前記変化領域候補決定手段によって検出された変化領域候補の各々について、前記検出対象識別器を用いて、前記変化領域候補が前記変化検出対象を含むか否かを識別し、前記変化検出対象を含まない前記変化領域候補を除外するステップと、
を含む変化領域検出方法。
【請求項７】
コンピュータを、請求項１〜請求項５のいずれか１項に記載の変化領域検出装置を構成する各手段として機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、変化領域検出装置、方法、及びプログラムに関する。
【背景技術】
【０００２】
航空写真や衛星画像からの変化検出はその有用性から、長い期間にわたりコンピュータビジョンやリモートセンシングの分野において研究されてきた。特に近年の急速な都市化における土地利用形態の変化のモニタリングや都市計画への利用のためにその重要性はますます増している。
【０００３】
異なる時点で撮影された航空写真や衛星画像から地上の変化を検出する技術は大きく２つに分けられる。１つは各時点での地上の三次元構造を復元し、その高さ方向の変化量から変化領域を検出しようとするものである。もう１つは画素間の地理的な対応関係が決定された複数の画像を利用するものである。これは各画素値を利用して、地上の物体の変化などを画素ごとに数値化し、このような数値の変化量によって表現される画像の変化から変化領域を検出するものである。
【０００４】
前者の技術として近年広く研究されているものとして、航空機によって撮影されたパノラマ画像を用いるものがある。これは、ステレオ画像を利用して撮影した領域の3次元構造（Digital Surface Modelと呼ばれる）を復元し、異なる時点における高さの変化が大きな領域を変化領域とするものである。三次元構造復元技術の進歩によって、地上に存在する物体の高さを精度よく得ることが出来るようになったため、これを利用して高い精度で変化領域を検出することが出来る（非特許文献１参照）。
【０００５】
最先端の研究成果によるとかなり高い精度で変化領域を推定することが出来るものの、航空機を利用してステレオ画像を取得することが難しいことに加え、通常の航空写真を撮影するよりもコストがかかるため技術の適用範囲が限定されてしまう。また、地上の変化は必ずしも高さ方向の変化を伴うものではない。例えば、２階建ての家屋が建て変わった場合などを想像してみればよい。したがってこのような手法によって検出できる変化の種類は限られている。
【０００６】
もう一方の三次元情報を利用しない技術について、ステレオ画像ではない通常の航空写真に画像処理の手法を適用することによって変化領域を推定する手法も数多く研究されている。こちらの技術の特徴としては、地上の物体の高さに依存したものではないため、より多くの種類の変化を検出できるということが挙げられる。また、ステレオ画像と比較すると通常の航空写真は取得しやすく、技術の適用範囲は広くなる。
【０００７】
そのようなもののひとつに建築物のスペクトル変化に注目したものがある。これは建築物が変化すると、これにともなって色だけではなく、反射に代表される光学的な性質が変化するという事実を利用し、このような変化を数値化する指標を用いて変化領域を検出する手法である。変化領域を決定する際にはファジー論理やデンプスター＝シェーファーの決定理論などが用いられる（非特許文献２参照）。
【０００８】
また、航空写真に対して画像処理の手法に加えてノンパラメトリックな密度推定手法とＥＭアルゴリズムと呼ばれる統計的機械学習の手法を適用することによって変化領域を検出する研究例もある（非特許文献３参照）。
【０００９】
２つの異なる時点において撮影された航空写真から確度の高い変化領域を発見し、当該領域を正例として教師あり学習を行い、識別器を構成するという研究例もある（非特許文献４参照）。なお、ここでいう識別器とは、数値ベクトルを入力としてうけとり、入力パターンが注目物体であるか否かを出力する関数一般を指す。
【００１０】
航空写真からの変化検出を行う際に従来から大きな障害となっているのが影の変化である。したがって、影の領域を推定し、その影の変化による影響を事前に除こうとする研究例が数多くある（非特許文献５参照）。
【００１１】
航空写真は撮影したタイミングによって影が生じる方向が変わったり、森林の植生が変わったりしてしまう。このような変化によっても航空機からの地上の見え方が変わってしまうため、本来検出したい土地利用変化以外の非本質的な変化を検出してしまうという難点がある。
【先行技術文献】
【非特許文献】
【００１２】
【非特許文献１】J.Tian, S.Cui, and P.Reinartz, “Building Change Detection Based on Satelite Stereo Imagery and Digital Surface Models,” IEEE Transactions on Geoscience and Remote Sensing, vol.52, No.1, 2014.
【非特許文献２】P.Du, S.Liu, P.Gamba, K.Tan, and J.Xia, “Fusion of Defference Images for Change Detection Over Urban Areas,” IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, Vol.5, No.4, 2012.
【非特許文献３】L.Bruzzone, and D.F.Prieto, “An Adaptive Semiparametric and Context-Based Approach to Unsupervised Change Detection in Multitemporal Remote-Sensing Images,” IEEE Transactions on Image Processing, Vol.11, No.4, 2002.
【非特許文献４】L.Castellana, A.D’addabbo, and G.Pasquariello, “A Composed supervised/unsupervised approach to improve change detection from remote sensing,” Pattern Recognition Letter, Vol.28, No.4, 2007.
【非特許文献５】T.T.Ngo, C.Collet, and V.Mazet, “MRF and Dempster-Shafer Theory for Simultaneous Shadow/Vegetation Detection on High Resolution Aerial Color Images,” In Proceedings of ICIP, 2014.
【発明の概要】
【発明が解決しようとする課題】
【００１３】
ステレオ画像を利用しない変化領域検出技術の多くは画素ごとの変化情報と統計的な推定を利用するものであり、パターンの意味づけや学習などの知的な情報処理を含まない。そのため、例えば赤色の建物が黒色の建物に変わったときに、それは黒い建物なのか影の影響によるものなのかを正確に判別することが出来ない。さらに、純粋に物理的な性質に基づいた指標のみを利用した場合、植生の変化など時間の経過による地上のテクスチャの変化を誤って土地利用形態の変化として検出してしまう可能性が高くなる。
【００１４】
一方、知的な情報処理を含む技術もこれまで提案されてきたが、学習のために自動的に生成された情報が必ずしも高品質ということができず、不正確なノイズの混入を防ぐことが出来ない。そのため、このようにして得られた識別器の効果は限定的であると考えられる。人力で学習のためのデータを生成することも考えられるが必要となるデータの量と労力を考えると現実的ではない。加えて、学習のための情報を必要としない技術も提案されてはきたものの、処理に時間がかかるため、大量のデータに適用するには不向きである。
【００１５】
本発明は、上記事情に鑑みてなされたものであり、変化領域を精度よく検出するための変化領域検出装置、方法、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【００１６】
上記目的を達成するために、本発明に係る変化領域検出装置は、異なる第１時点及び第２時点で撮影された航空写真と、前記第１時点に対応する地図情報とを用いて土地利用変化を検出する変化領域検出装置であって、変化領域検出を行う際に、入力された前記第１時点の航空写真の画素を特徴づける局所特徴ベクトル、及び前記第２時点の航空写真の画素を特徴づける局所特徴ベクトルを計算する局所特徴ベクトル計算手段と、前記第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、前記画素ごとに変化量を計算する変化量計算手段と、前記画素ごとの変化量から、変化領域候補を検出する変化領域候補決定手段と、前記第１時点に対応する地図情報と前記第１時点の航空写真とから、前記第１時点の航空写真の小領域と、変化検出対象を含むか否かを示す情報とからなる検出対象識別器学習用データセットを作成する学習用データセット構築手段と、前記検出対象識別器学習用データセットから、変化検出対象を含むか否かを識別するための検出対象識別器を学習する検出対象識別器構築手段と、前記変化領域候補決定手段によって検出された変化領域候補の各々について、前記検出対象識別器を用いて、前記変化領域候補が前記変化検出対象を含むか否かを識別し、前記変化検出対象を含まない前記変化領域候補を除外する非検出対象領域除外手段と、を含んで構成されている。
【００１７】
また、本発明に係る変化領域検出装置において、前記学習用データセット構築手段は、前記地図情報から、前記変化検出対象の物体の存在する画像座標を取得し、当該画像座標に基づいて前記学習用データセットを作成するようにすることができる。
【００１８】
また、本発明に係る変化領域検出装置において、前記局所特徴ベクトル計算手段は、注目画素及び前記注目画素の周辺の小領域の色又は輝度の情報も利用して、注目画素周辺の形又はテクスチャに関する情報を記述する前記局所特徴ベクトルを計算するようにすることができる。
【００１９】
また、本発明に係る変化領域検出装置において、前記検出対象識別器学習手段は、前記検出対象識別器学習用データセットに含まれる小領域の画素を特徴づける局所特徴ベクトルを計算し、前記計算された画素の前記局所特徴ベクトルをまとめ上げる処理によって、前記検出対象識別器学習用データセットに含まれる小領域の大域特徴ベクトルを計算し、前記大域特徴ベクトルに基づいて変化検出対象を含むか否かを識別するための前記検出対象識別器を学習し、前記非検出対象領域除外手段は、前記変化領域候補決定手段によって検出された変化領域候補の各々について、前記第２時点の航空領域の前記変化領域候補の画素を特徴づける局所特徴ベクトルを計算し、前記計算された画素の前記局所特徴ベクトルをまとめ上げる処理によって、前記変化領域候補の大域特徴ベクトルを計算し、前記計算した大域特徴ベクトルに基づいて、前記検出対象識別器を用いて、前記変化領域候補が前記変化検出対象を含むか否かを識別するようにすることができる。
【００２０】
また、本発明に係る変化領域検出装置において、前記変化量計算手段は、前記第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、前記画素ごとに、利用する局所特徴ベクトルに合わせた距離関数から計算される距離を、前記変化量として計算するようにすることができる。
【００２１】
本発明に係る変化領域検出方法は、局所特徴ベクトル計算手段、変化量計算手段、変化領域候補決定手段、学習用データセット構築手段、検出対象識別器構築手段、及び非検出対象領域除外手段を含み、異なる第１時点及び第２時点で撮影された航空写真と、前記第１時点に対応する地図情報とを用いて土地利用変化を検出する変化領域検出装置における変化領域検出方法であって、前記局所特徴ベクトル計算手段が、変化領域検出を行う際に、入力された前記第１時点の航空写真の画素を特徴づける局所特徴ベクトル、及び前記第２時点の航空写真の画素を特徴づける局所特徴ベクトルを計算するステップと、前記変化量計算手段が、前記第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、前記画素ごとに変化量を計算するステップと、前記変化領域候補決定手段が、前記画素ごとの変化量から、変化領域候補を検出するステップと、前記学習用データセット構築手段が、前記第１時点に対応する地図情報と前記第１時点の航空写真とから、前記第１時点の航空写真の小領域と、変化検出対象を含むか否かを示す情報とからなる検出対象識別器学習用データセットを作成するステップと、前記検出対象識別器構築手段が、前記検出対象識別器学習用データセットから、変化検出対象を含むか否かを識別するための検出対象識別器を学習するステップと、前記非検出対象領域除外手段が、前記変化領域候補決定手段によって検出された変化領域候補の各々について、前記検出対象識別器を用いて、前記変化領域候補が前記変化検出対象を含むか否かを識別し、前記変化検出対象を含まない前記変化領域候補を除外するステップと、を含む。
【００２２】
本発明に係るプログラムは、コンピュータを、上記の変化領域検出装置を構成する各手段として機能させるためのプログラムである。
【発明の効果】
【００２３】
本発明の変化領域検出装置、方法、及びプログラムによれば、第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから画素ごとに変化量を計算し、画素ごとの変化量から変化領域候補を検出し、第１時点に対応する地図情報と第１時点の航空写真とから検出対象識別器学習用データセットを作成し、検出対象識別器学習用データセットから、変化検出対象を含むか否かを識別するための検出対象識別器を学習し、検出された変化領域候補の各々について、検出対象識別器を用いて、変化領域候補が変化検出対象を含むか否かを識別し、変化検出対象を含まない変化領域候補を除外することにより、変化領域を精度よく検出することができる、という効果が得られる。
【図面の簡単な説明】
【００２４】
【図１】本発明の実施の形態に係る変化領域検出装置の機能的構成を示すブロック図である。
【図２】小領域を等間隔にサンプリングする場合における、格子点を中心とする小領域の取り出し方を説明するための図である。
【図３】局所特徴ベクトルの距離計算における格子点のペアの取り方を説明するための図である。
【図４】格子点以外の画素の距離決定の方法を説明するための図である。
【図５】小領域をフィルタリングする処理を説明するための図である。
【図６】本発明の実施の形態に係る変化領域検出装置における変化領域候補生成処理ルーチンのフローチャート図である。
【図７】本発明の実施の形態に係る変化領域検出装置における局所特徴ベクトル計算処理ルーチンのフローチャート図である。
【図８】本発明の実施の形態に係る変化領域検出装置における変化量計算処理ルーチンのフローチャート図である。
【図９】本発明の実施の形態に係る変化領域検出装置における変化領域候補決定処理ルーチンのフローチャート図である。
【図１０】本発明の実施の形態に係る変化領域検出装置における学習処理ルーチンのフローチャート図である。
【図１１】本発明の実施の形態に係る変化領域検出装置における学習用データセット構築処理ルーチンのフローチャート図である。
【図１２】本発明の実施の形態に係る変化領域検出装置における検出対象識別器構築処理ルーチンのフローチャート図である。
【図１３】本発明の実施の形態に係る変化領域検出装置における大域特徴ベクトル変換処理ルーチンのフローチャート図である。
【図１４】本発明の実施の形態に係る変化領域検出装置における変化領域検出処理ルーチンのフローチャート図である。
【図１５】本発明の実施の形態に係る変化領域検出装置における変化領域候補決定処理ルーチンのバリエーションを示すフローチャート図である。
【発明を実施するための形態】
【００２５】
以下、図面を参照して本発明の実施の形態を詳細に説明する。なお、実施の形態を説明するための全図において、同一機能を有するものは同一符号をつけ、その繰り返しの説明は省略する。
【００２６】
本発明の実施の形態のポイントは２点ある。
【００２７】
１点目は、過去に作成された地図を参照することによって得られる、変化検出対象が存在する画像座標から教師情報を生成することにより、人手で新たに教師情報を作成するという煩雑な作業を行うことなく、高品質な教師情報を得ることができるという点である。これを利用して機械学習手法により、注目対象を分類する識別器を構築し、この識別器を用いて本来検出すべきでない対象を最終的な結果から除外することにより、最終的に精度良く変化領域を検出することができる。
【００２８】
２点目は、一般物体認識において用いられている特徴量を変化検出対象の分類に用いることにより、画素ごとの情報ではなく小領域の見た目に基づく情報を利用しているという点である。これにより、従来利用されてきた画素ごとに注目していた特徴量では利用できなかった形やテクスチャに基づく情報を新たに利用できるようになり、地図から得た教師情報と組み合わせることにより、正確な分類を実現することができる。
【００２９】
本発明の実施の形態に係る変化領域検出装置は、撮影年月日が異なるが位置合わせがなされた同じ解像度を持つ２枚の航空写真と、電子的に記録された地図情報を組み合わせて変化検出対象を正確に識別できる識別器を学習し、学習された識別器を用いて変化検出対象を含まない変化領域候補を取り除き、最終的な変化領域を座標の集合として出力する装置である。また、本発明の実施の形態に係る変化領域検出装置は、航空写真からの変化領域検出に関し、特に統計的機械学習手法を用いて変化領域を検出する技術に関する。
【００３０】
＜本発明の実施の形態に係る変化領域検出装置の構成＞
【００３１】
図１に示すように、本発明の実施の形態に係る変化領域検出装置１は、ＣＰＵと、ＲＡＭと、後述する各処理ルーチンを実行するためのプログラムや各種データを記憶したＲＯＭと、を含むコンピュータで構成することが出来る。この変化領域検出装置１は、機能的には図１に示すように入力部１０と、演算部２０と、出力部４０を備えている。
【００３２】
本発明の実施の形態において、旧航空写真とは、入力された２枚の航空写真のうち、より過去において撮影されたものを指し、本発明の実施の形態では、旧航空写真を第１時点で撮影された航空写真と称する。また、新航空写真とは、より現在に近い時点において撮影されたものを指し、本発明の実施の形態では、新航空写真を第２時点で撮影された航空写真と称する。
【００３３】
本発明の実施の形態に係る変化領域検出装置１は、異なる第１時点及び第２時点で撮影された航空写真と、第１時点に対応する地図情報とを用いて土地利用変化を検出する。また、地図情報は、旧航空写真が撮影された当時のものを想定しているが、必ずしもその限りではない。
【００３４】
入力部１０は、第１時点で撮影された航空写真と第２時点で撮影された航空写真とを受け付ける。また、入力部１０は、後述する演算部２０によって用いられる、地図情報、検出対象情報、正例閾値、及び負例閾値を受け付ける。
【００３５】
演算部２０は、局所特徴ベクトル計算部２２と、変化量計算部２４と、変化領域候補決定部２６と、学習用データセット構築部２８と、検出対象識別器構築部３０と、非検出対象領域除外部３２とを備えている。
【００３６】
局所特徴ベクトル計算部２２は、変化領域検出を行う際に、入力部１０より入力された、第１時点の航空写真の画素を特徴づける局所特徴ベクトル、及び第２時点の航空写真の画素を特徴づける局所特徴ベクトルを計算する。
【００３７】
具体的には、局所特徴ベクトル計算部２２は、入力された航空写真上で等間隔に配置された格子点を中心とする小領域を、局所特徴ベクトルと呼ばれる数値ベクトルに変換する。
【００３８】
図２に、小領域を一定間隔で取り出す操作を説明するための図を示す。図２において黒点は画像上に等間隔で配置された格子点を表す。また、黒点を中心とした円は小領域を表す。図２においては見易さのために小領域を円で表したが、小領域は任意の形状でよい。局所特徴ベクトル計算部２２は、各小領域について、当該小領域に含まれる画素の色や輝度値に基づいて当該小領域を数値ベクトルに変換する。そして、局所特徴ベクトル計算部２２は、各格子点座標に対応する数値ベクトルの集合を出力する。
【００３９】
例えば、局所特徴ベクトル計算部２２は、第１時点の航空写真の注目画素及び注目画素の周辺の小領域の色又は輝度の情報も利用して、注目画素周辺の形又はテクスチャに関する情報を記述する局所特徴ベクトルを計算する。
【００４０】
変化量計算部２４は、局所特徴ベクトル計算部２２によって出力された、第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、画素ごとに変化量を計算する。
【００４１】
具体的には、変化量計算部２４は、第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、画素ごとに、利用する局所特徴ベクトルに合わせた距離関数から計算される距離を、変化量として計算し、航空写真の解像度と同じサイズで各要素に航空写真間の変化量が格納された二次元配列を出力する。
【００４２】
二次元配列を距離マップとすると、距離マップ［ｘ］［ｙ］には２つの航空写真上で座標が（ｘ，ｙ）であるような画素間の変化量が格納されている（画像左上をｘ＝０，ｙ＝０とする）。変化量を計算するに当たって、まず、格子点に相当する画素に対応した局所特徴ベクトル同士の距離を計算し、当該距離を各格子点同士の変化を表す量として、距離マップの当該箇所に格納する。
【００４３】
図３に、第１時点の航空写真と第２時点の航空写真との格子点同士での対応付けを表す図を示す。変化量計算部２４は、図３に示すように、格子点に相当する画素に対応した局所特徴ベクトル同士の距離を計算し、当該距離を各格子点同士の変化を表す量として、距離マップの当該箇所に格納する。
【００４４】
なお、ここでいう距離とは必ずしもユークリッド距離を指すものではなく、カーネル関数を利用するなど様々な形態が考えられる（例えば、非特許文献６参照）。ここで、ｖ１とｖ２をｄ次元ベクトル、
【００４５】
【数１】
0006397379.tif 000002
【００４６】
とすると、２つのベクトル間のユークリッド距離Ｌ（ｖ１，ｖ２）は、
【００４７】
【数２】
0006397379.tif 000003
【００４８】
と表される量である。
【００４９】
［非特許文献６］赤穂昭太郎，「カーネル多変量解析 --- 非線形データ解析の新しい展開」，岩波書店，2008．
【００５０】
次に、変化量計算部２４は、格子点以外の点の変化量を決定する。そのため、まず変化量計算部２４は、距離マップを、格子点を中心とする重複のない正方領域に分割し、各正方領域に含まれる座標の変化量を、正方領域の中心となっている座標の変化量と同じとする。
【００５１】
図４に、格子点を３画素ごとに配置して変化量を算出する場合の例を示す。図４に示すように、黒い四角は格子点に相当する配列要素を表し、太線で囲まれた正方領域に含まれる配列要素には中心となっている配列要素と同じ変化量が格納される。
【００５２】
変化量計算部２４は、以上のようにして作成された、入力画像と同じサイズをもち、かつ各要素が画素間の変化量であるような二次元配列を出力する。
【００５３】
変化領域候補決定部２６は、変化量計算部２４によって計算された画素ごとの変化量から、変化領域候補を検出する。
【００５４】
具体的には、変化領域候補決定部２６は、変化量計算部２４で計算された二次元配列である距離マップを入力として受け取り、当該距離マップに対して閾値処理を施す。閾値の決定方法としては人手で与えるほかに、変化量の分布から統計的に求める方法などが考えられる。本発明の実施の形態では、閾値が人手によって予め与えられている場合を例に説明する。
【００５５】
変化量計算部２４によって計算された画素ごとの変化量が、予め決定された閾値より大きな画素で、かつ、連続的につながっている領域を１つにまとめて得られる画像座標の集合をＳ１，Ｓ２，…，ＳＮとし、当該画像座標の集合を変化領域候補とする。ここでＮは自然数で得られた候補の数を表す。また、各Ｓｉ（ｉ＝１，…，Ｎ）は、具体的には、
【００５６】
【数３】
0006397379.tif 000004
【００５７】
のように表される（ＮｉはＳｉに含まれる画像座標の総数を表す）。そして、変化領域候補決定部２６は、変化領域候補の集合Ｓ＝｛Ｓ１，Ｓ２，…，ＳＮ｝を出力する。
【００５８】
学習用データセット構築部２８は、入力部１０によって受け付けられた第１時点に対応する地図情報と第１時点の航空写真とから、第１時点の航空写真の小領域と、変化検出対象を含むか否かを示す情報とからなる検出対象識別器学習用データセットを作成する。具体的には、学習用データセット構築部２８は、地図情報から、変化検出対象の物体の存在する画像座標を取得し、当該画像座標に基づいて学習用データセットを作成する。
なお、地図情報としては、例えば、第１時点に対応する地図を表す画像（例えば、複数の多角形領域を含んで表される。）を用いる。また、検出対象情報としては、例えば、地図情報の多角形領域の各々が何を表すのか（例えば、家屋又は道路など）を示す情報を用いる。
【００５９】
具体的に、まず、学習用データセット構築部２８は、第１時点の航空写真と、第１時点の航空写真を撮影した当時の地図情報及び検出対象情報を入力として、第１時点の航空写真からランダムかつ大量に小領域を取り出す。
【００６０】
次に、学習用データセット構築部２８は、小領域を抽出する小領域集合フィルタリング処理において、変化検出対象が、入力部１０により入力された正例閾値以上の割合で含まれている小領域を正例、入力部１０により入力された負例閾値以下の割合でしか含まれていない小領域を負例として検出対象識別器学習用データセットに追加していく。ランダムに取得された小領域に含まれる画素が変化検出対象であるかどうかは、入力された地図情報及び検出対象情報を参照することによって決定される。
【００６１】
図５に、小領域に含まれる変化検出対象の領域の比率を求める操作を表す。図５において点線で表される正方領域はランダムに選択された小領域を表し、灰色の長方形の領域は変化検出対象が当該領域に存在することを示す。
【００６２】
図５において、Ｐ１、Ｐ４、Ｐ６は変化検出対象の領域との重複が大きいため、正例として検出対象識別器学習用データセットに加えられる。Ｐ２、Ｐ３は変化検出対象の領域との重複が十分小さいため、負例として検出対象識別器学習用データセットに加えられる。Ｐ５は変化検出対象の領域との重複が、正例となるほど大きくもなく、負例となるほど小さくも無いため、どちらにも加えられず、検出対象識別器学習用データセットには含まれない。
【００６３】
そして、学習用データセット構築部２８は、検出対象識別器学習用データセットが所望のサイズに到達するまで繰り返し行い、最終的に得られた画像の集合と、各画像に変化検出対象が含まれるか否かという情報とを含む検出対象識別器学習用データセットを出力とする。
【００６４】
検出対象識別器構築部３０は、変化検出対象の種類毎に、学習用データセット構築部２８によって作成された検出対象識別器学習用データセットから、変化検出対象を含むか否かを識別するための検出対象識別器を学習する。
【００６５】
具体的には、検出対象識別器構築部３０は、検出対象識別器学習用データセットに含まれる小領域の画素を特徴づける局所特徴ベクトルを計算し、計算された画素の局所特徴ベクトルをまとめ上げる処理によって、検出対象識別器学習用データセットに含まれる小領域の大域特徴ベクトルを計算し、大域特徴ベクトルに基づいて変化検出対象を含むか否かを識別するための、各変化検出対象に対応した検出対象識別器を学習する。
【００６６】
ここで、変化検出対象をＯ１，Ｏ２，…，ＯＭとし、各変化検出対象に対応した検出対象識別器をｆ１，ｆ２，…，ｆＭとすると、ｆｊ（ｊ＝１、…、Ｍ）はＤ次元ベクトルｘを入力として受け取り数値を出力する関数で、ｆｊ（ｘ）＞０のとき入力ｘに対応する領域を変化検出対象物体Ｏｊであるとみなす。複数の識別器の出力が正であるときは最も大きな数値を出力した識別器に対応した変化検出対象物体であるとみなす。また、すべての識別器の出力が負であったときには、入力ｘに対応する領域には変化検出対象物体が含まれていないとみなす。
【００６７】
学習に当たって、検出対象識別器学習用データセットに含まれる画像や変化検出対象が含まれるか否かという教師情報は、学習しようとする検出対象識別器の種類に応じて適切に数値ベクトルなどに変換される。また、検出対象識別器学習用データセットに含まれる画像は大域特徴ベクトルと呼ばれるベクトルに変換する。
【００６８】
検出対象識別器構築部３０においては、入力領域として小領域全体を入力とする。また、局所特徴ベクトル統合処理としては、局所特徴ベクトルの平均を計算する、あるいは次元ごとに最大値を求めるなどの処理が考えられる。大域特徴ベクトルとしては一般物体認識などで用いられるBag-of-Visual-WordsやFisher Vectorといった様々な特徴ベクトルを利用することが考えられる。
【００６９】
検出対象識別器としては様々なものが考えられるが、具体的なものを挙げれば、パーセプトロンやニューラルネットワーク、サポートベクターマシンなどがある。
【００７０】
検出対象識別器構築部３０は、以上のようにして学習された検出対象識別器を出力する。
【００７１】
非検出対象領域除外部３２は、変化領域候補決定部２６によって検出された変化領域候補の各々について、検出対象識別器構築部３０によって出力された検出対象識別器を用いて、変化領域候補が変化検出対象を含むか否かを識別し、変化検出対象を含まない変化領域候補を除外する。
【００７２】
具体的には、非検出対象領域除外部３２は、変化領域候補決定部２６によって検出された変化領域候補の各々について、第２時点の航空領域の変化領域候補の画素を特徴づける局所特徴ベクトルを計算し、計算された画素の前記局所特徴ベクトルをまとめ上げる処理によって、変化領域候補の大域特徴ベクトルを計算し、計算した大域特徴ベクトルに基づいて、検出対象識別器を用いて、変化領域候補が変化検出対象を含むか否かを識別し、最終的に出力する変化領域を決定する。
【００７３】
まず、非検出対象領域除外部３２は、変化領域候補の集合から、各変化領域候補に相当する航空写真上の領域を検出対象識別器にかけるために、各変化領域候補に相当する航空写真上の領域を大域特徴ベクトルに変換する。変化領域候補に対応する画像の大域特徴ベクトル変換がこのプロセスにあたる。
【００７４】
検出対象識別器構築部３０においては入力領域として小領域全体を用いたが、本プロセスにおいては変化領域候補の集合Ｓに含まれるＳｉが入力となる。すなわちＳｉに含まれる画像領域から局所特徴ベクトルを計算し、局所特徴ベクトル統合処理によって抽出した局所特徴ベクトルを大域特徴ベクトルに変換する。ここで用いられる大域特徴ベクトルは検出対象識別器構築部３０において、検出対象識別器の学習に用いられたものと同じものである。
【００７５】
次に、非検出対象領域除外部３２は、大域特徴ベクトルを検出対象識別器にかけることにより、各変化領域候補に変化検出対象が含まれているかどうかを判別し、変化検出対象が含まれている変化領域候補のみを最終的に変化領域として出力する。変化領域候補Ｓｉから計算された大域特徴ベクトルをｘ（Ｓｉ）とすると、最終的に出力される領域集合Ｔは、
【００７６】
【数４】
0006397379.tif 000005
【００７７】
となる。
このとき、本質的ではない見た目の変化に対して頑健な大域特徴ベクトルを利用することにより、撮影時点の違いから生じる建物の影の変化や季節の移り変わりによる田畑や森林の変化、自動車などの動物体の移動によるテクスチャの変化といった、本来変化領域として検出すべきではない領域を正確に出力結果から除くことができ、最終的に精度よく変化領域を把握することが可能となる。
【００７８】
出力部４０は、非検出対象領域除外部３２によって出力された変化領域を結果として出力する。
【００７９】
＜本発明の実施の形態に係る変化領域検出装置の作用＞
【００８０】
次に、本発明の実施の形態に係る変化領域検出装置１の作用について説明する。第１時点で撮影された航空写真及び第２時点で撮影された航空写真が変化領域検出装置１に入力されると、変化領域検出装置１は、図６に示す変化領域候補生成処理ルーチンを実行する。
【００８１】
＜変化領域候補生成処理ルーチン＞
ステップＳ１００において、入力部１０は、第１時点で撮影された航空写真と第２時点で撮影された航空写真とを受け付ける。
【００８２】
次に、ステップＳ１０２において、局所特徴ベクトル計算部２２は、変化領域検出を行う際に、上記ステップＳ１００で受け付けた、第１時点の航空写真の画素を特徴づける局所特徴ベクトル、及び第２時点の航空写真の画素を特徴づける局所特徴ベクトルを計算する。ステップＳ１０２は、図７に示す局所特徴ベクトル計算処理ルーチンによって実現される。
【００８３】
＜局所特徴ベクトル計算処理ルーチン＞
ステップＳ２００において、上記ステップＳ１００で受け付けた、第１時点で撮影された航空写真と第２時点で撮影された航空写真とを取得する。
【００８４】
ステップＳ２０２において、上記ステップＳ２００で取得された第１時点で撮影された航空写真及び第２時点で撮影された航空写真について、当該航空写真上で等間隔に配置された格子点を中心とする小領域をサンプリングする。
【００８５】
ステップＳ２０４において、上記ステップＳ２０２でサンプリングされた小領域の各々について、当該小領域を局所特徴ベクトルに変換する。
【００８６】
ステップＳ２０６において、上記ステップＳ２０４で得られた局所特徴ベクトルの集合を結果として出力して、局所特徴ベクトル計算処理ルーチンを終了する。
【００８７】
次に、変化領域候補生成処理ルーチンに戻り、ステップＳ１０４において、変化量計算部２４は、上記ステップＳ１０２で出力された局所特徴ベクトルの集合について、第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから、画素ごとに変化量を計算する。ステップＳ１０４は、図８に示す変化量計算処理ルーチンによって実現される。
【００８８】
＜変化量計算処理ルーチン＞
ステップＳ３００において、上記ステップＳ１０２で出力された局所特徴ベクトルの集合を取得する。
【００８９】
ステップＳ３０２において、上記ステップＳ３００で取得した局所特徴ベクトルの集合に基づいて、各格子点のペアについて、当該格子点に相当する画素に対応した局所特徴ベクトル同士の距離を計算し、当該距離を各格子点同士の変化を表す量として、距離マップの当該箇所に格納する。
【００９０】
ステップＳ３０４において、上記ステップＳ３０２で格子点間の変化量が格納された距離マップを、格子点を中心とする重複のない正方領域に分割し、各正方領域に含まれる座標の変化量を正方領域の中心となっている格子点の変化量と同じとする。
【００９１】
ステップＳ３０６において、上記ステップＳ３０４で得られた距離マップを結果として出力して、変化量計算処理ルーチンを終了する。
【００９２】
次に、変化領域候補生成処理ルーチンに戻り、ステップＳ１０６において、変化領域候補決定部２６は、上記ステップＳ１０４で計算された画素ごとの変化量から、変化領域候補を検出する。ステップＳ１０６は、図９に示す変化領域候補決定処理ルーチンによって実現される。
【００９３】
＜変化領域候補決定処理ルーチン＞
ステップＳ４００において、上記ステップＳ１０４で出力された距離マップを取得する。
【００９４】
ステップＳ４０２において、上記ステップＳ４００で取得された距離マップに基づいて、距離マップの要素に対応する画素の各々について、当該画素の変化量が、予め決定された閾値より大きい画素を抽出する。
【００９５】
ステップＳ４０４において、上記ステップＳ４０２の閾値処理によって抽出された画素の各々について、連続的につながっている領域を１つに結合し変化領域候補とする。
【００９６】
ステップＳ４０６において、上記ステップＳ４０４で得られた変化領域候補の各々を変化領域候補の集合として出力し、変化領域候補決定処理ルーチンを終了する。
【００９７】
次に、変化領域候補生成処理ルーチンに戻り、ステップＳ１０８において、上記ステップＳ１０６で出力された変化領域候補の集合を結果として出力し、変化領域候補生成処理ルーチンを終了する。
【００９８】
次に、学習処理ルーチンについて説明する。
【００９９】
第１時点で撮影された航空写真と、地図情報、検出対象情報、正例閾値、及び負例閾値とが変化領域検出装置１に入力されると、変化領域検出装置１は、図１０に示す学習処理ルーチンを実行する。
【０１００】
＜学習処理ルーチン＞
ステップＳ５００において、入力部１０は、第１時点で撮影された航空写真と、地図情報、検出対象情報、正例閾値、及び負例閾値とを受け付ける。
【０１０１】
ステップＳ５０２において、学習用データセット構築部２８は、上記ステップＳ５００で受け付けられた第１時点に対応する地図情報と第１時点の航空写真とから、第１時点の航空写真の小領域と、変化検出対象を含むか否かを示す情報とからなる検出対象識別器学習用データセットを作成する。ステップＳ５０２は、図１１に示す学習用データセット構築処理ルーチンによって実現される。
【０１０２】
＜学習用データセット構築処理ルーチン＞
ステップＳ６００において、上記ステップＳ５００で受け付けられた、第１時点で撮影された航空写真と、地図情報、検出対象情報、正例閾値、及び負例閾値とを取得する。
【０１０３】
ステップＳ６０２において、上記ステップＳ６００で取得した第１時点の航空写真から、小領域をランダムに複数個取得する。
【０１０４】
ステップＳ６０４において、上記ステップＳ６００で取得された、地図情報、検出対象情報、正例閾値、及び負例閾値に基づいて、上記ステップＳ６０２で取得された小領域の各々について、変化検出対象が正例閾値以上の割合で含まれている当該小領域を正例、負例閾値以下の割合でしか含まれていない当該小領域を負例として検出対象識別器学習用データセットに追加する。
【０１０５】
ステップＳ６０６において、上記ステップＳ６０４のフィルタリング処理によって得られた検出対象識別器学習用データセットを結果として出力して、学習用データセット構築処理ルーチンを終了する。
【０１０６】
次に学習処理ルーチンに戻り、ステップＳ５０４において、検出対象識別器構築部３０は、変化検出対象の種類毎に、上記ステップＳ５０２で構築された検出対象識別器学習用データセットから、変化検出対象を含むか否かを識別するための検出対象識別器を学習する。ステップＳ５０４は、図１２に示す検出対象識別器構築処理ルーチンによって実現される。
【０１０７】
＜検出対象識別器構築処理ルーチン＞
ステップＳ７００において、上記ステップＳ５０２で出力された検出対象識別器学習用データセットを取得する。
【０１０８】
ステップＳ７０２において、検出対象識別器構築部３０は、上記ステップＳ７００で取得された検出対象識別器学習用データセットに含まれる小領域の画素を特徴づける局所特徴ベクトルを計算し、計算された画素の局所特徴ベクトルをまとめ上げる処理によって、検出対象識別器学習用データセットに含まれる小領域の大域特徴ベクトルを計算する。ステップＳ７０２は、図１３に示す大域特徴ベクトル変換処理ルーチンによって実現される。
【０１０９】
＜大域特徴ベクトル変換処理ルーチン＞
ステップＳ８００において、上記ステップＳ７００で取得された検出対象識別器学習用データセットに含まれる複数の小領域を取得する。
【０１１０】
ステップＳ８０２において、上記ステップＳ８００で取得された複数の小領域の各々について、当該小領域の画素を局所特徴ベクトルに変換する。
【０１１１】
ステップＳ８０４において、上記ステップＳ８００で取得された複数の小領域の各々について、上記ステップＳ８０２で得られた局所特徴ベクトルを統合して、大域特徴ベクトルを生成する。
【０１１２】
ステップＳ８０６において、上記ステップＳ８０４で得られた大域特徴ベクトルの各々を結果として出力して、大域特徴ベクトル変換処理ルーチンを終了する。
【０１１３】
次に検出対象識別器構築処理ルーチンに戻り、ステップＳ７０４において、検出対象識別器構築部３０は、上記ステップＳ７０２で出力された大域特徴ベクトルの各々に基づいて、各変化検出対象Ｏ１，Ｏ２，…，ＯＭに対応した検出対象識別器ｆ１，ｆ２，…，ｆＭを学習する。
【０１１４】
ステップＳ７０６において、検出対象識別器構築部３０は、上記ステップＳ７０４で学習された検出対象識別器ｆ１，ｆ２，…，ｆＭを結果として出力して、検出対象識別器構築処理ルーチンを終了する。
【０１１５】
次に、変化領域検出処理ルーチンについて説明する。
【０１１６】
変化領域候補生成処理ルーチンによって変化領域候補の集合が出力され、学習処理ルーチンによって検出対象識別器が出力されると、変化領域検出装置１は、図１４に示す変化領域検出処理ルーチンを実行する。
【０１１７】
＜変化領域検出処理ルーチン＞
ステップＳ８００において、非検出対象領域除外部３２は、変化領域候補生成処理ルーチンのステップＳ１００で受け付けた第２時点で撮影された航空写真と、変化領域候補生成処理ルーチンによって出力された変化領域候補の集合と、学習処理ルーチンによって出力された検出対象識別器とを取得する。
【０１１８】
ステップＳ８０２において、非検出対象領域除外部３２は、上記ステップＳ８００で取得された第２時点で撮影された航空写真と変化領域候補の集合とに基づいて、変化領域候補の集合に含まれる画像領域から局所特徴ベクトルを計算し、計算された局所特徴ベクトルを統合して大域特徴ベクトルに変換する。ここで用いられる大域特徴ベクトルは検出対象識別器構築部３０において、検出対象識別器の学習に用いられたものと同様のものである。
【０１１９】
ステップＳ８０４において、非検出対象領域除外部３２は、上記ステップＳ８００で取得された検出対象識別器と、上記ステップＳ８０２で得られた大域特徴ベクトルとに基づいて、変化領域候補の集合に含まれる変化領域候補の各々について、検出対象識別器を用いて、変化領域候補が変化検出対象を含むか否かを識別する。
【０１２０】
ステップＳ８０６において、非検出対象領域除外部３２は、上記ステップＳ８０４で得られた識別結果に基づいて、変化領域候補の集合から、変化検出対象を含まない変化領域候補を除外する。
【０１２１】
ステップＳ８０８において、出力部４０は、上記ステップＳ８０６によって変化検出対象を含まない変化領域候補が除外された変化領域候補の集合を結果として出力して、変化領域検出処理ルーチンを終了する。
【０１２２】
以上説明したように、本発明の実施の形態に係る変化領域検出装置によれば、第１時点及び第２時点の各々の航空写真で同じ座標を持つ画素において計算された局所特徴ベクトルから画素ごとに変化量を計算し、画素ごとの変化量から変化領域候補を検出し、第１時点に対応する地図情報と第１時点の航空写真とから検出対象識別器学習用データセットを作成し、検出対象識別器学習用データセットから、変化検出対象を含むか否かを識別するための検出対象識別器を学習し、検出された変化領域候補の各々について、検出対象識別器を用いて、変化領域候補が変化検出対象を含むか否かを識別し、変化検出対象を含まない変化領域候補を除外することにより、変化領域を精度よく検出することができる。
【０１２３】
また、過去にある地点において撮影された航空写真をもとに、同じ地点で新たに撮影された航空写真から、影や植生の変化がある場合でも、そのような変化に対して頑健に、かつ、新規に学習用のデータを作成するという煩雑な作業を介することなく変化領域を検出することが可能となる。
【０１２４】
また、電子的に利用可能となっている地図情報を識別器の学習に利用することにより変化領域の検出精度の改善を図ることができる。
【０１２５】
また、検出対象を識別するための識別器を学習するための高品質なデータセットを新たな労力をかけることなく構築し、構築されたデータセットを用いて学習した識別器を利用して、検出対象を含まない変化領域候補を最終的な出力結果から除くことができる。
【０１２６】
なお、本発明は、上述した実施の形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
【０１２７】
例えば、本発明の実施の形態では、変化領域候補の集合を生成する処理と、検出対象識別器を学習する処理と、変化領域を検出する処理とを１つの装置で行う場合を例に説明したが、変化領域候補の集合を生成する処理と、検出対象識別器を学習する処理と、変化領域を検出する処理とを、別々の装置で実行してもよい。その場合には、変化領域候補の集合を生成する装置は、局所特徴ベクトル計算部２２と、変化量計算部２４と、変化領域候補決定部２６とを含んで構成され、検出対象識別器を学習する装置は、学習用データセット構築部２８と、検出対象識別器構築部３０とを含んで構成され、変化領域を検出する処理装置は、非検出対象領域除外部３２を含んで構成される。
【０１２８】
また、本発明の実施の形態において、変化領域候補決定部２６は、予め定められた閾値に基づいて閾値処理を行う場合を例に説明したが、これに限定されるものではなく、例えば変化量の分布から閾値を統計的に求めてもよい。この場合、変化領域候補決定処理ルーチンは、図１５に示すように、閾値を決定するステップＳ４０１を更に含む。
【０１２９】
また、本発明の実施の形態において、学習用データセット構築部２８は、地図情報及び検出対象情報を用いて、学習用データセットを作成する場合を例に説明したが、これに限定されるものではなく、例えば、地図情報に検出対象情報が含まれており、学習用データセット構築部２８は、地図情報を用いて学習用データセットを作成してもよい。
【０１３０】
また、上述の変化領域検出装置は、内部にコンピュータシステムを有しているが、コンピュータシステムは、ＷＷＷシステムを利用している場合であれば、ホームページ提供環境（あるいは表示環境）も含むものとする。
【０１３１】
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読取り可能な記録媒体に格納して提供することも可能である。
【符号の説明】
【０１３２】
１０   入力部
２０   演算部
２２   局所特徴ベクトル計算部
２４   変化量計算部
２６   変化領域候補決定部
２８   学習用データセット構築部
３０   検出対象識別器構築部
３２   非検出対象領域除外部
４０   出力部
【図１】
0006397379.tif 000006
【図２】
0006397379.tif 000007
【図３】
0006397379.tif 000008
【図４】
0006397379.tif 000009
【図５】
0006397379.tif 000010
【図６】
0006397379.tif 000011
【図７】
0006397379.tif 000012
【図８】
0006397379.tif 000013
【図９】
0006397379.tif 000014
【図１０】
0006397379.tif 000015
【図１１】
0006397379.tif 000016
【図１２】
0006397379.tif 000017
【図１３】
0006397379.tif 000018
【図１４】
0006397379.tif 000019
【図１５】
0006397379.tif 000020