(57)【特許請求の範囲】
【請求項１】
処理部を備え、画像を解析して作物の特徴量を抽出する画像解析システムであって、
前記処理部は、
解析対象とする１つもしくは複数の作物の植付け間隔と、作物のスペクトルもしくは作物の半径を含む成長パターンをパラメータとして、解析対象とする前記画像の各領域から、解析対象とする前記作物の特徴量を抽出するものであり、
前記パラメータを用いて複数の擬似画像を生成し、解析対象の前記画像の各領域に対して、複数の前記擬似画像から最も適合度の高い擬似画像を選択し、選択した最も適合度の高い前記擬似画像の情報を用いて、解析対象とする前記作物の特徴量を抽出することを特徴とする画像解析システム。
【請求項２】
請求項１に記載の画像解析システムであって、
前記処理部は、
解析対象とする前記画像の撮影範囲もしくは、及び撮影時期によって、前記作物の特徴量を抽出するために利用する複数の前記パラメータを限定することを特徴とする画像解析システム。
【請求項３】
請求項１に記載の画像解析システムであって、
前記処理部は、
解析対象とする前記画像の解像度、バンド数に応じて、生成する前記擬似画像の解像度、バンド数を変更することを特徴とする画像解析システム。
【請求項４】
請求項１に記載の画像解析システムであって、
前記処理部は、
前記パラメータとして、前記作物の植付け状態が縦横方向に揃った格子状になっているかどうかを表す格子フラグを更に用いることを特徴とする画像解析システム。
【請求項５】
請求項１に記載の画像解析システムであって、
前記処理部は、
解析対象とする前記画像の解像度が所定の閾値以上の場合、前記作物の被覆率によるスペクトルの混合比率を前記擬似画像の情報とすることを特徴とする画像解析システム。
【請求項６】
請求項４に記載の画像解析システムであって、
前記処理部は、
解析対象とする前記画像の解像度が、所定の閾値以上の場合、前記作物の被覆率に基づく混合スペクトルを前記擬似画像の情報とすることを特徴とする画像解析システム。
【請求項７】
請求項４に記載の画像解析システムであって、
解析対象とする前記画像の解像度が、所定の閾値未満の場合、土壌スペクトルをバンドセットデータで平均化したスペクトルで、前記擬似画像の背景画像を生成することを特徴とする画像解析システム。
【請求項８】
請求項４に記載の画像解析システムであって、
解析対象とする前記画像の解像度が、所定の閾値未満の場合、前記格子フラグの値に基づき、前記擬似画像を格子画像あるいは直列画像として生成することを特徴とする画像解析システム。
【請求項９】
画像を解析して作物の特徴量を抽出する処理部による画像解析方法であって、
前記処理部は、
解析対象とする１つもしくは複数の作物の植付け間隔と、作物の成長パターンをパラメータとして、解析対象とする前記画像の各領域から、解析対象とする前記作物の特徴量を抽出する際に、
前記パラメータを用いて複数の擬似画像を生成し、解析対象の前記画像の各領域に対して、複数の前記擬似画像から最も適合度の高い擬似画像を選択し、選択した最も適合度の高い前記擬似画像の情報を用いて、解析対象とする前記作物の特徴量を抽出することを特徴とする画像解析方法。
【請求項１０】
請求項９に記載の画像解析方法であって、
前記処理部は、
解析対象とする前記画像の撮影範囲もしくは、及び撮影時期によって、前記作物の特徴量を抽出するために利用する複数の前記パラメータを限定することを特徴とする画像解析方法。
【請求項１１】
請求項９に記載の画像解析方法であって、
前記処理部は、
前記パラメータとして、前記作物の植付け状態が縦横方向に揃った格子状になっているかどうかを表す格子フラグを更に用いることを特徴とする画像解析方法。
【請求項１２】
請求項９に記載の画像解析方法であって、
前記処理部は、
解析対象とする前記画像の解像度、バンド数に応じて、生成する前記擬似画像の解像度、
バンド数を変更することを特徴とする画像解析方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、農業を行う圃場を撮影した画像を高精度に解析する画像解析技術に関するものである。
【背景技術】
【０００２】
農業統計や、営農支援を目的に、人工衛星や航空機によって撮影された画像が活用されている。特に衛星画像は広範囲の農業を撮影でき、面積当りの単価が航空画像に比べて安価なため、農業への応用が進んでいる。
【０００３】
具体的には、収穫前の小麦圃場を撮影した画像の赤バンドと近赤外線バンドのピクセル値の差と和の比によって定義される正規化植生指標ＮＤＶＩ（Ｎｏｒｍａｌｉｚｅｄ  Ｄｉｆｆｅｒｅｎｃｅ  Ｖｅｇｅｔａｔｉｏｎ  Ｉｎｄｅｘ）を用いて圃場ごとの小麦の水分量を推定し、収穫順を決定したり、同様に水稲圃場のたんぱく質含有量を推定し、酒造用の酒米の買い取り価格を決定するなどの活用例がある。
【０００４】
また、全国的な作付状況を把握するために、衛星画像を活用する取り組みも行われている。具体的には衛星に搭載した合成開口レーダーの画像を、田植え前の画像と、成長後の画像を比較することによって、水稲圃場を特定する技術や、光学画像を用いてピクセルごとのスペクトルの差を基準として分類することで圃場ごとの作付分類を行う技術が開発されている。
【０００５】
一方、民間で利用できる衛星画像の解像度は光学マルチスペクトルの場合で、最高２ｍとなっている。この解像度は、繁茂期の圃場中心部であれば１ピクセルが対象となる作物のスペクトルを表現していると考えられるレベルであるが、圃場周辺部では１ピクセル中に作物のスペクトルと土壌のスペクトルが混合されて撮影されることとなる。また、生育の初期の段階では、作物の被覆率が低いため、作物と土壌のスペクトルが混合されたピクセルの割合が多くなる。このようなピクセル内のスペクトル混合は、ＮＤＶＩの計算やスペクトルによる分類に悪影響を与える。
【０００６】
これを解決する手段として、ミクセル分解が提案されている（非特許文献１）。ミクセル分解は、各ピクセルのスペクトルが、複数の対象物のスペクトルの線形結合で決定されるというモデルに基づいている。複数の対象物のスペクトルが既知の場合は、混在比率のみを推定する問題となり、複数の対象物のスペクトルが未知の場合は、各スペクトルと混在比率の両者を推定する問題となる。ミクセル分解は、混在したスペクトルを複数の対象物のスペクトルに分解することで、ＮＤＶＩの計算やスペクトルによる分類を行う際の精度と安定性の向上につながる。
【０００７】
一方、作物の植付けパターンを画像解析に利用している例として特許文献１が挙げられる。この文献では、作物の形と作物の配置情報を含んだテンプレートを用いた画像を解析することによって、作物の位置を精度良く検出できるようにしている。
【先行技術文献】
【特許文献】
【０００８】
【特許文献１】特開２００３−１０２２７５号公報
【非特許文献】
【０００９】
【非特許文献１】秋山  侃他著「農業リモートセンシングハンドブック」
【発明の概要】
【発明が解決しようとする課題】
【００１０】
しかしながら、上述の非特許文献１においては、ミクセル分解ではピクセル内でのスペクトル・混合比率の推定に関して、農業に特化した情報を利用していないという問題がある。具体的には、作物ごとに異なる、作物の植付けパターン・植付け時期・成長パターンなどの情報を活用できていない。
【００１１】
一方、特許文献１では、複数の異なるテンプレートを同時に利用することが考慮されておらず、対象となる圃場に植えられた作物に対応したテンプレートを用いることで、作物の植付け位置を正確に検出することを目的としている。また、利用されるテンプレートは、植付け位置を検出することが目的のため、作物の生育による形状やスペクトルの変化が考慮されていない。
【００１２】
更に、衛星画像や航空画像の解析を行う場合、それぞれの圃場に作付された作物に関する知識が無い、圃場ごとに異なる作物が作付されている、作付時期が異なる、圃場の形状情報が無いなどの理由によって、特定のテンプレートを選択的に利用することが困難であり、各地点に対して最適なテンプレートを選択的に利用する必要があるという課題がある。
【００１３】
本発明の目的は、上記の課題を解決し、地域性、作物種類、作付時期等の差を考慮して、衛星画像や航空画像から作物の特徴を表す特徴量を精度良く抽出する画像解析システム、及び方法を提供することにある。
【課題を解決するための手段】
【００１４】
上記の目的を達成するため、本発明においては、処理部を備え、画像を解析して作物の特徴量を抽出する画像解析システムであって、処理部は、解析対象とする１つもしくは複数の作物の植付け間隔と、作物のスペクトルもしくは作物の半径のいずれかを含む成長パターンをパラメータとして、解析対象とする画像の各領域から、解析対象とする作物の特徴量を抽出する画像解析システムを提供する。
【００１５】
また、上記の目的を達成するため、本発明においては、画像を解析して作物の特徴量を抽出する処理部による画像解析方法であって、処理部は、解析対象とする１つもしくは複数の作物の植付け間隔と、成長パターンをパラメータとして、解析対象とする画像の各領域から、解析対象とする作物の特徴量を抽出する画像解析方法を提供する。
【００１６】
好適な態様の本発明の画像解析システムにあっては、対象とする作物と、その植付け間隔と、その成長パターンをパラメータセットとして、対象とする作物の画像から得られる特徴量を推定する。この特徴量の推定の手法の一つとして、対象とする作物の植付け間隔と成長パターンをパラメータセットとして、複数の擬似画像を生成し、解析対象の画像の各領域に対して、複数の擬似画像から最も適合度の高い擬似画像を選択し、最も適合度の高い擬似画像の情報を用いて、各領域から対象とする作物の特徴量を抽出する。
【発明の効果】
【００１７】
本発明の画像解析システムは、解析対象の画像から対象とする作物の特徴量を精度良く抽出することを可能にする。
【図面の簡単な説明】
【００１８】
【図１】第一の実施例に係る、画像解析システムの構成の一例を示す図である。
【図２】第一の実施例に係る、植生ＤＢの構成の一例を示す図である。
【図３】第一の実施例に係る、エリアテーブルの一例を示す図である。
【図４】第一の実施例に係る、地域作付テーブルの一例を示す図である。
【図５】第一の実施例に係る、作付テーブルの一例を示す図である。
【図６】第一の実施例に係る、作物テーブルの一例を示す図である。
【図７】第一の実施例に係る、作物のスペクトルパターンの一例を示す図である。
【図８】第一の実施例に係る、スペクトルテーブルの一例を示す図である。
【図９】第一の実施例に係る、バンドセットテーブルの一例を示す図である。
【図１０】第一の実施例に係る、成長パターンを表すグラフの一例を示す図である。
【図１１】第一の実施例に係る、成長テーブルの一例を示す図である。
【図１２】第一の実施例に係る、圃場ポリゴンテーブルの一例を示す図である。
【図１３】第一の実施例に係る、圃場作付テーブルの一例を示す図である。
【図１４】第一の実施例に係る、画像解析システムの動作全体フローの一例を示す図である。
【図１５】第一の実施例に係る、パターンＩＤリストの作成フローの一例を示す図である。
【図１６】第一の実施例に係る、擬似画像の作成フローの一例を示す図である。
【図１７】第一の実施例に係る、擬似画像の作成フローの一例（続き）を示す図である。
【図１８】第一の実施例に係る、最適テンプレート探索フローの一例を示す図である。
【図１９】第一の実施例に係る、最適テンプレート探索フローの一例（続き）を示す図である。
【図２０】第一の実施例に係る、格子画像の中央画像、端画像、角画像の一例を示す図である。
【図２１】第一の実施例に係る、土壌スペクトルテーブルの一例を示す図である。
【図２２】第一の実施例に係る、直列画像の一例を示す図である。
【図２３】第一の実施例に係る、特徴量抽出フローの一例を示す図である。
【図２４】第一の実施例に係る、格子状に配置されない植付けの一例を示す図である。
【発明を実施するための形態】
【００１９】
以下、図面を用いて、本発明の種々の実施例を説明する。なお、本明細書において、作物の成長パターンとは、農業の圃場等に植付けされた作物の成長の時間的な変化を示すパラメータであり、例えば作物の半径の経時変化、或いは作物のスペクトルの経時変化等を意味し、本発明においては、作物の成長パターンには、少なくとも作物のスペクトルもしくは作物の半径のいずれかを含むものとする。また、本明細書において、格子フラグとは、作物の植付け状態が縦横方向に揃った格子状になっているかどうかを表す情報を意味し、擬似画像の情報とは、解析対象画像の疑似画像の作物の空間的な配置情報や、作物の被覆率によるスペクトルの混合比率の情報等を意味する。
【実施例１】
【００２０】
図１は、第一の実施例の画像解析システムの一構成例を表した図である。同図のシステムにおいて、１００はコンピュータ構成を備えた画像解析サーバである。画像解析サーバ１００単体で画像解析システムを構成することができるが、同図に示すように、ネットワーク１０２を経由して端末１０１を画像解析サーバ１００と接続した構成とすることもできる。
【００２１】
画像解析サーバ１００は、データ通信バス１０３によって、制御部１０４、入力部１０５、表示部１０６、通信部１０７、植生データベース（Ｄａｔａ  Ｂａｓｅ：ＤＢ）編集部１０８、植生ＤＢ１０９、パターンリスト作成部１１０、擬似画像生成部１１１、最適テンプレート選択部１１２、特徴量抽出部１１３、画像ＤＢ１１４が接続された構成を備えることができる。制御部１０４は、画像解析サーバ１００上の処理全体を制御するものであり、中央処理部（Ｃｅｎｔｒａｌ  Ｐｒｏｃｅｓｓｉｎｇ  Ｕｎｉｔ：ＣＰＵ）等で実現される。この処理部であるＣＰＵは、省略されたメモリ等の記憶部に記憶されるプログラムを実行することにより、上述したパターンリスト作成部１１０、擬似画像生成部１１１、最適テンプレート選択部１１２、特徴量抽出部１１３、更には植生ＤＢ編集部１０８を実現する。図示を省略したメモリは記憶部を構成し、植生ＤＢ１０９や画像ＤＢ１１４が形成される。
【００２２】
入力部１０５は、利用者からのシステムに対しての情報提示の要求や、データ入力などを受付けるものであり、キーボード、マウス、タッチパネルなどによって実現される。表示部１０６は、利用者に対して各種の情報提示を行うものであり、液晶ディスプレイ等で実現される。通信部１０７は、ネットワーク１０２との接続インタフェース（Ｉｎｔｅｒｆａｃｅ：Ｉ／Ｆ）であり、ＬＡＮ等で実現される。植生ＤＢ編集部１０８は、後述する植生ＤＢ１０９の内容を編集するためのユーザインタフェースであり、ソフトウェアとして実現される。
【００２３】
植生ＤＢ１０９〜特徴量抽出部１１３が、本実施例の画像解析システムの主要構成要素である。植生ＤＢ１０９は、作付時期、作付作物、作付間隔、作物の成長パターン等の情報を保存するためのデータベースであり、リレーショナルデータベース等を用いて実現できる。ここで作物の成長パターンとは、上述の通り、作物の半径やスペクトル混合率等のスペクトルの経時変化を意味している。
【００２４】
図２に示すように、本実施例の植生ＤＢ１０９には、エリアテーブル２０１、地域作付テーブル２０２、作付テーブル２０３、作物テーブル２０４、スペクトルテーブル２０５、土壌スペクトルテーブル２０６、成長テーブル２０７、圃場ポリゴンテーブル２０８、圃場作付テーブル２０９が含まれている。
【００２５】
図３に一例を示すエリアテーブル２０１は、地域による栽培作物・作付時期・農法等の差を表現するための基本データとして、各地域を緯度経度情報を含むポリゴンで表現している。さらにエリアテーブル２０１は、各地域での土壌のタイプを表す土壌識別子（Ｉｄｅｎｔｉｆｉｅｒ：ＩＤ）を保持している。各地域での土壌タイプは複数存在しても良い。
【００２６】
図４に一例を示す地域作付テーブル２０２は、各地域でどのような農法で作付をしているかを表現している。パターンＩＤは、作付品種と農法を表現した作付テーブルへのインデックスである。作付時期・時期幅は、その地域でパターンＩＤで指定される農法で栽培される作物が、作付される時期とその時期幅を示している。収穫時期・時期幅も同様である。解析対象の画像の撮影時期と、作付時期・収穫時期を比較することで、撮影時期に、その地域で栽培されている作物を限定することができる。
【００２７】
図５に一例を示す作付テーブル２０３は、農法を表現している。どの作物（品種ＩＤ）を、どのような畝間・株間、即ち植付け間隔で植付けるかという植付け間隔の情報、植付け状態が縦横方向に揃った格子状になっているかどうかを表す情報である格子フラグを保有している。
【００２８】
図６に一例を示す作物テーブル２０４は、品種ＩＤと、品目名・品種名を対応づける情報を記録している。
【００２９】
図７の作物の反射スペクトル７００に示すように、各作物はそれぞれ異なったスペクトルのパターンを持ち、生育時期によってもスペクトルパターンは変化、即ち経時変化をする。そして、分光放射計等を用いて、詳細にこの作物のスペクトルパターンを計測することができる。
【００３０】
図８にその一例を示すスペクトルテーブル２０５は、スペクトルＩＤで示す各作物の生育ステージ毎に異なるスペクトルの情報を記録している。同図において、第一行目に波長(nm)を示している。ここでは、波長300nmから、波長1100nmまで5nm分解能で計測したスペクトルパターンを各列に記録しているが、波長域、分解能について限定するものではない。スペクトルパターンの単位としては、反射率を用いることができる。測定対象と反射率が既知の白色板をそれぞれ分光放射計で計測し、計測された分光強度と、白色板の反射率から、測定対象の反射率を計算することができる。スペクトルテーブル２０５は、周波数分解能の高い分光放射計等のセンサを用いた地上計測で作成することが可能である。
【００３１】
図２１にその一例を示す土壌スペクトルテーブル２０６は、土壌ＩＤで示す土壌タイプごとの土壌のスペクトルパターンを記録している。各列は各波長に対応するスペクトル値を示している。
【００３２】
図９にその一例を示すバンドセットテーブル２０７は、衛星画像や航空画像を撮影したセンサの周波数解像度を表現するテーブルであり、各センサＩＤに対応するセンサ名、バンド数、バンドセットデータを示している。上述の通り、スペクトルテーブル２０５は、周波数分解能の高い分光放射計等のセンサを用いた地上計測で作成することが可能であるが、衛星画像や航空画像では、波長域が４バンド、８バンド等と分解能が低い場合が多い。
【００３３】
図９に示すバンドセットテーブル２０７は、このような場合に、スペクトルテーブル２０５に記録されているスペクトルデータを波長域ごとに平均化して、衛星画像や航空画像に合わせるために用いる。ここでは平均化の処理を行っているが、撮影に用いたセンサの各波長に対する感度の情報が入手できる場合は、平均ではなく、感度の情報を用いて加重平均を行っても良い。
【００３４】
図１０は、作物の生育を経時変化で表現した生育パターンを示すグラフの一例を示す図である。図１０の(a)は、上から見た作物の大きさ（半径）の変化を表し、図１０の(b)は、生育に伴う作物のスペクトル混合率の変化を表している。両グラフとも、横軸には、例えば積算気温を用いることができるが時間の経過を表すものであれば、積算気温に限定されない。
【００３５】
図１１にその一例を示す成長テーブル２０８は、図１０で表現される生育パターンを記録したテーブルである。品種ＩＤで検索すると、複数のスペクトル/半径ＩＤが得られるように構成されている。このうちスペクトル/半径ＩＤ０は半径を表し、単位はcmで記録されている。０以外のスペクトル/半径ＩＤはスペクトルテーブル２０５のスペクトルＩＤに対応している。品種ＩＤ２で検索した場合、スペクトル/半径ＩＤは０，３，４が得られ、半径は積算温度が増えるにつれ、３ｃｍ→６ｃｍ→１３ｃｍ→１４ｃｍと大きくなることを示している。
【００３６】
一方、図１１のスペクトルパターンに関しては、植付け時（積算温度０）スペクトルＩＤが３のスペクトルパターンが１００％を占めているが、積算温度が２００を超えた時点でスペクトルパターン３と４の混合が始まり、積算気温３００になった時点で、１００％スペクトルパターン４になることを表現している。
【００３７】
図１２にその一例を示す圃場ポリゴンテーブル２０９は、各圃場の位置・形状の情報を緯度・経度形式のポリゴンとして記録し、さらに圃場ごとの土壌タイプをオプションで保持する。圃場ポリゴンテーブル２０９で土壌タイプが保持されている圃場に関しては、エリアテーブル２０１で指定された土壌タイプの代わりに、圃場ポリゴンテーブル２０９で指定された土壌タイプを用いて背景を生成する。
【００３８】
図１３にその一例を示す圃場作付テーブル２１０は、各圃場で作付されている作物と作付方法を現すパターンＩＤと、作付時期を記録している。全ての圃場ＩＤに対して、パターンＩＤ、作付時期が記録されている必要は無い。また、パターンＩＤだけが記録され、作付時期が記録されていない場合もある。
【００３９】
図１に示した画像解析サーバ１００のパターンリスト作成部１１０から特徴抽出部１１３までは、上述の通り、植生ＤＢ編集部１０８と同様に処理部で実行されるソフトウェアとして実現される。パターンリスト作成部１１０は、解析対象である画像の撮影範囲と撮影時期によって、マッチング対象とするパターンのリストを作成する。擬似画像生成部１１１は、上記パターンのリストに対応する農作物の各成長時期に応じた擬似画像を生成する。最適テンプレート選択部１１２は、解析対象画像の各領域に対して上記生成された擬似画像を用いてテンプレートマッチングを行い、最も適合度の高い擬似画像を選択する。特徴量抽出部１１３は、解析対象画像の各領域から、最も適合度の高い擬似画像の情報を用いて、作物、植生の特徴量を抽出する。
【００４０】
画像解析サーバ１０４中の画像ＤＢ１１４は、解析対象となる衛星画像や航空画像を記録するためのＤＢであり、リレーショナルデータベースやファイルシステムで構成することができる。図示していないＣＤ−ＲＯＭ・ＤＶＤ−ＲＯＭドライブなどを経由して画像を取り込んだり、通信部１０７を通じてネットワーク１０２経由で、図示していない外部の画像提供サービスから画像を取り込むことが出来る。
【００４１】
以下に本実施例の画像解析システムの動作について、図１４〜図１９に示すフローチャートを参照して説明する。
図１４は、本実施例の画像解析システム、画像解析サーバ１００の全体フローの一例を示す。まず解析対象となる画像の撮影範囲と撮影時期の情報を用いて、その地域でその時期に栽培されている可能性のある作物のパターンＩＤをリストアップする（Ｓ１４０１）。衛星画像は通常撮影範囲と撮影時期の情報が含まれた形で提供されているため、このステップはパターンＩＤを限定する上で重要である。航空画像に関しても、航空機にＧＰＳ（Ｇｌｏｂａｌ  Ｐｏｓｉｔｉｏｎｉｎｇ  Ｓｙｓｔｅｍ）やジャイロが搭載されている場合、航空機の位置と姿勢から、撮影された画像を地図上に投影できる。
【００４２】
次に、Ｓ１４０１でリストアップされた全パターンＩＤに対して擬似画像を生成する（Ｓ１４０２）。さらに、生成した擬似画像をテンプレートとして、解析対象に対するテンプレートマッチングを行い、各領域に対する最適テンプレートを探索し決定する（Ｓ１４０３）。このテンプレートマッチングにより決定された最適テンプレートには、作物の種類、植付け時期に関する情報が含まれており、さらに最適テンプレートが解析対象画像と最もマッチする位置と角度を得ることが出来るため、解析対象の各領域で、どの作物が、どのような位置関係で、どのような生育状態で育っているかを推定することができる。これらの情報は、圃場単位に決まる情報であるため、圃場単位特徴量と呼ぶ。
【００４３】
最後に、決定された最適テンプレートを用いて、解析対象画像からピクセル単位に作物の特徴量を抽出する（Ｓ１４０４）。抽出する特徴量は、スペクトル情報やＮＤＶＩ等の指標、ピクセルを中心としたテクスチャ指標等であり、これをピクセル単位特徴量と呼ぶ。本実施例の画像解析システムは、どの作物が、どのような位置関係で、どのような生育状態で育っているか推定できているために、従来の特徴量抽出手法よりも、より正確に作物のピクセル単位特徴量を抽出することが可能となる。
【００４４】
次に図１５を用いて、Ｓ１４０１の詳細なフローについて説明する。まず、画像の撮影範囲情報を用いて、エリアテーブル２０１に記録されている各ポリゴンと重なりのあるエリアを決定する（Ｓ１５０１）。撮影地域が複数地域にまたがる場合は、複数のエリアＩＤが決定される場合がある。次に、Ｓ１５０１で決定されたエリアＩＤと撮影時期によって、地域作付テーブルから該当するデータを抽出し、パターンＩＤリストを作成する（Ｓ１５０２）。撮影時期による該当判定は、撮影時期が、「時期作付−作付時期幅」と「収穫時期＋収穫時期幅」の間に入っているかどうかで判定する。解析対象画像に撮影時期の情報が含まれていない場合は、エリアＩＤの該当判定だけを行う。
【００４５】
次に図１６を用いて、Ｓ１４０２の詳細なフローについて説明する。Ｓ１４０１で作成されたパターンＩＤリストに含まれるパターンＩＤそれぞれに対して、次のステップを実行する。
まず、作付テーブル２０３から、パターンＩＤに該当する作付パラメータ（作物の品種ＩＤ、畝間や株間などの植付け間隔、格子フラグ）を取得する（Ｓ１６０１）。
【００４６】
次に、バンドセットテーブル２０７から、解析対象の画像を取得したセンサに対応したバンドセットを選択し、平均レンジ情報を取得する（Ｓ１６０２）。解析対象画像に対応したセンサＩＤは、画像ＤＢ１１４に画像を格納する際に付与しておく。次に成長テーブル２０８から、品種ＩＤに対応するデータを選択する（Ｓ１６０３）。
【００４７】
最後に、Ｓ１６０３で選択されたデータの各時点について図１７のフローを実行し、擬似画像を作成する（Ｓ１６０４）。ここで各時点とは、成長テーブル２０８の、０、１００、２００等の列ごとに記載されたパラメータを利用しても良いし、線形補間によって、ある列と次の列の間の時点のパラメータを生成し利用しても良い。
【００４８】
図１７を用いて引き続き、Ｓ１４０２の疑似画像の生成についての詳細なフローについて説明する。
擬似画像を生成する際には、解析対象と解像度を一致させて生成する必要がある。この際、解像度が大きい（粗い）と、一定間隔で植えられた作物１つ１つを判別することができなくなり、１ピクセルが作物と土壌が混合した状態を表現することになり、作物の空間的な配置情報は失われ、作物の被覆率によるスペクトルの混合比率が主な情報となる。一方で、航空画像、特に無人飛行機で撮影した航空画像については解像度が小さく（細かく）、作物の空間的な配置情報を含んだ形で擬似画像を生成することができる。そこで、まず解析対象の画像の解像度に基づいて処理を分ける（Ｓ１７０１）。ここでは閾値として、例えば５０ｃｍを設定する。これによって現在民間利用可能な衛星画像を解析対象とする場合は、Ｓ１７０２以降のステップが使われ、航空機や無人航空機で撮影された画像の場合は、ほとんどの場合Ｓ１７０４以降のステップが使われることになる。なお、上述の作物の空間的な配置情報や、作物の被覆率によるスペクトルの混合比率情報が、擬似画像の情報となる。
【００４９】
解像度が閾値以上だった場合、農作物の被覆率の算出を行う（Ｓ１７０２）。圃場の端では作物が植えられていない場合があるため、中央部の被覆率と周辺部の被覆率は異なる可能性がある。そこで、ここでは一例として、一様に作物が植えられている状態を中央画像、３分の１もしくは３分の２が土壌のみとなっている状態をそれぞれ1/3画像、2/3画像、全てが土壌となっている状態を土壌画像として、それぞれの被覆率を求める。中央画像の被覆率は、作物の半径が株間・畝間の半分以下である場合は、半径×半径×π/（畝間×株間）で求められる。半径がそれ以上になった場合も、重なり部分を考慮して計算すれば良い。土壌画像の被覆率は０であり、1/3画像の被覆率は、中央画像の被覆率×2/3であり、2/3画像の被覆率は、中央画像の被覆率×1/3である。
【００５０】
次に、求められた被覆率を用いて、混合スペクトルを算出する（Ｓ１７０３）。混合スペクトルの算出は、まず算出する時点の各スペクトルの混合比を成長テーブル２０８から得たパラメータより算出する。成長テーブル２０８に記載されている時点について算出する場合は、混合比のパラメータそのものを利用すれば良いし、そうでない場合はすでに述べたように線形補間をすれば良い。これによって得られた混合比で、バンドセットデータを用いて平均化した作物のスペクトルを加重平均し、その時点での作物スペクトルを算出し、さらにこの作物スペクトルと、バンドセットデータを用いて平均化した土壌スペクトルを被覆率で加重平均することで、擬似画像の情報としての混合スペクトルを得ることが出来る。
【００５１】
一方、解像度が閾値未満だった場合は、作成するパターンの格子フラグにより処理を振り分ける（Ｓ１７０４）。
格子フラグがｔｒｕｅ（真）の場合は、土壌スペクトルをバンドセットデータで平均化したスペクトルを用いて背景画像を生成する（Ｓ１７０５）。背景画像のサイズは例えば2m×2mとし、例えば1cm解像度で生成する。この場合、200pixel×200pixelの画像となる。生成される画像はバンドセットの数と同じレイヤー（チャンネル、バンドとも呼ぶ）を持ち、第Ｎレイヤーの全ての画素値は、上記土壌スペクトルの第Ｎ番目の値となる。
【００５２】
次に、格子画像を生成する（Ｓ１７０６）。格子画像は、ここでは中央画像・端画像・角画像の３種類を生成することとする。図２０に示すように、中央画像は、一様に植生が配置されており、端画像は左の一列が土壌になっており、角画像は左と上の一列が土壌になっている。
【００５３】
格子画像を生成するため、Ｓ１７０５で作成した背景の上に、Ｓ１７０３で説明したのと同様の方法で算出した各時点での作物スペクトルを用いて、畝間・株間のデータで指定された植付け間隔で、その生育時点での半径の円を描けば良い。背景より大きな範囲について、円を描画した上で背景の大きさでトリミングすることで、図２０のような周辺部で円が欠けた擬似画像を生成できる。端画像や角画像は、左や上を畝間・株間分空けてから円の描画を始めれば良い。
以上によって、格子フラグがｔｒｕｅだった場合について、格子状の疑似画像が作成される。
【００５４】
格子フラグがｆａｌｓｅ（偽）の場合は、作物の半径の２倍の幅を持ち、例えば長さ２ｍのサイズの背景を、土壌スペクトルをバンドセットデータで平均化したスペクトルを用いて作成する（Ｓ１７０７）。次に、直列画像を生成する（Ｓ１７０８）。この場合、図２４で表すように畝間・株間の植付け間隔は一定だが、隣接する畝の間で配置にズレがある可能性があることを意味している。
【００５５】
従って、図２０のような格子状配置の擬似画像をテンプレートとするのは適切でないため、図２２のような直列配置の擬似画像を生成する。この直列画像は、植生画像と土壌画像を生成し、植生画像は株間情報に従って１列に円を描画する。一方、土壌画像は、背景画像をそのまま使う。以上によって、格子フラグがｆａｌｓｅだった場合について、直列状の疑似画像が作成される。
【００５６】
格子フラグがｔｒｕｅだった場合、Ｓ１７０６で格子画像を生成した後に、図１７に示すように、さらにＳ１７０７とＳ１７０８を実行し、格子画像に加えて、直列画像を生成しても良い。これは、最適テンプレート探索を行う際に、格子画像ではなく、直列画像を用いて行う場合があるためである。
【００５７】
最後に、解析対象の画像相当の解像度に擬似画像を合わせて縮小することにより解像度の変更を行う（Ｓ１７０９）。例えば、10cm解像度の場合の格子画像は、200pixel×200pixelの画像から、20pixel×20pixelに縮小する。縮小のアルゴリズムは問わない。
【００５８】
次に、図１８を用いて、Ｓ１４０３の最適テンプレート探索フローの詳細を説明する。
まず、図１８のフローが開始される時点で、パターンＩＤリストには擬似画像を生成した全てのパターンＩＤが記録されていることを述べておく。
最適テンプレート探索は、ポリゴン情報の有無により大きく処理が異なるので有無を確認する（Ｓ１８０１）。ポリゴン情報が無い場合、対象画像全領域に対して処理を行う必要が出てくる。以下に処理内容を説明する。
【００５９】
まず解像度がＳ１７０１で用いたのと同じ閾値未満かどうかで処理を振り分ける（Ｓ１８０２）。
解像度が閾値より大きく、粗い場合、擬似画像は混合スペクトルを持つ１pixel画像となるため、解析対象画像と擬似画像との適合度はスペクトル距離のみで決定される。スペクトル距離を計算する上で、考慮しなければならないのは、擬似画像で用いられているスペクトルは、例えば反射率で記録されており、直接解析対象の画素値と比較することができないことである。
【００６０】
これを解決する方法として、解析対象画像の画素値を反射率に変換する方法がある。これには、アスファルトなどの反射率が既知の領域を画像から選択して、その領域の画素値とアスファルト反射率の変換係数を算出し、その変換係数を他の領域にも適用することによって実現できる。
【００６１】
別の手法として、擬似画像および解析対象画像の両者の画素値を正規化指標に変換する方法である。正規化指標はピクセルごとの（バンド１の値−バンド２の値）/(バンド１の値＋バンド２の値）で定義される。適当なバンドを基準バンドとして正規化指標の定義式のバンド２に適用し、残りの全てのバンドを用いて正規化指標を算出すると、当初のバンド数−１のバンド数の正規化指標画像を得ることができる。もしくは、バンド１、バンド２の全ての組み合わせを用いて正規化指標を計算することもできる。正規化指標はバンド間の比に基づくものであり、撮影条件に対してロバストな性質を持っているため、解析対象画像と擬似画像のスペクトル距離を安定して比較できるようになる。スペクトル距離としては、ユークリッド距離や、スペクトル同士の為す角度などが利用できる。以上のスペクトル距離に関する考え方は、本実施例の以降の記載に関しても同様である。
【００６２】
上記の考え方に基づき、解析対象の画像の全ピクセルに対して、パターンＩＤリストに含まれる全パターンの、全混合スペクトル（中央・1/3・2/3・土壌）とのスペクトル距離を算出し、最も近いスペクトルを選択する。
【００６３】
一方、解像度が閾値未満だった場合は、パターンＩＤリストに含まれる全パターンの中央画像のみを用いて、解析対象の画像全領域に対してテンプレートマッチングを行う。テンプレートマッチングの方法としては、正規化相互相関等の既存の手法を用いることができる。注目するピクセルを中心として、一定値範囲での平行移動と、回転を考慮して中央画像と解析対象の画像の正規化相互相関を計算し、最も相関の高い条件を求めることで、注目するピクセル周辺で最も適合度の高いパターンＩＤ，位置、角度の情報を得ることができる。注目するピクセルを一定距離ずらしてこの処理を行うことで、解析対象画像全体に対して、最も適合度の高いパターンＩＤ，位置、角度の情報を得ることができる。このテンプレートマッチングの手法は、本実施例の後の記述に関しても同様である。
【００６４】
これによって、解析対象画像の全ピクセルに対して、作物種類・植付け時期・角度等の情報を得ることができる。ここで、植付け角度とは、例えば図２０に示すような、作物が植付けられた畝等が基準となす角度を意味する。これら作物種類・植付け時期・角度等の情報のピクセル間での類似性に基づいてピクセルをグループ化することで、ピクセルのグループに対して圃場単位特徴量を算出することもできる。これは、次に説明するポリゴン情報を用いた圃場単位特徴量の算出の代用となるものである。ピクセルごとに得られた作物種類・植付け時期・角度等の情報はノイズが大きい場合が多いが、グループ化によってノイズを低減し、安定した推定が可能になる。
【００６５】
一方、Ｓ１８０１でポリゴン情報が存在するとなった場合は、全圃場に対して図１９に記載したフローを実施することになる（Ｓ１８０５）。ポリゴン情報を活用することで、より効率的で、精度の高い最適テンプレート探索が可能となる。近年、ＧＩＳ（Ｇｅｏｇｒａｐｈｉｃ  Ｉｎｆｏｒｍａｔｉｏｎ  Ｓｙｓｔｅｍ）の農業分野への応用が普及しつつあり、ポリゴン情報の活用は現実的になっている。また、農業向けＧＩＳでは、作付品目や作付時期等についても管理される場合が多い。図２に示した植生ＤＢの圃場ポリゴンテーブル２０９、圃場作付テーブル２１０に関して、このような農業向けＧＩＳで管理されている情報を利用することが可能である。
【００６６】
図１９を用いて、Ｓ１８０５の詳細なフローの一例を説明する。ポリゴン情報がある場合、さらに圃場の作付情報があれば、テンプレートマッチングに用いる擬似画像を限定することがでるため、圃場作付テーブル２１０を検索し、パターンＩＤの有無で処理を振り分け（Ｓ１９０１）、パターンＩＤリストの内容を圃場作付テーブル２１０に記載されたパターンＩＤのみに限定する（Ｓ１９０２）。
【００６７】
さらに圃場作付テーブル２１０に作付時期が記載されていた場合、作付時期と撮影時期の比較によってテンプレートとして利用する擬似画像を限定することができる（Ｓ１９０３）。具体的には、本実施例のように各擬似画像が積算温度と関連付けられている場合、ネットワーク１０２経由で、気象庁のアメダスから各地域の積算温度を取得できるので、作付時期からの積算温度と擬似画像に関連付けられた積算温度が閾値以下になる擬似画像のみに限定してテンプレートとすることができる。
【００６８】
次にＳ１８０２と同様に、解析対象の画像の解像度によって処理を振り分ける（Ｓ１９０４）。
解像度が閾値以上で粗い場合は、ポリゴン内で、ポリゴン境界から閾値以上距離がある領域（以下、領域Ｃと呼ぶ）の全ピクセルに対して、パターンＩＤリストに含まれる全パターンＩＤの中央画像スペクトルを用いてスペクトル距離を算出し、各ピクセルに対して最も距離が近いパターンＩＤを決定する（Ｓ１９０５）。ポリゴン境界からの距離を考慮しているのは、作物が一様に植えつけられている領域のみを処理対象にするためである。これによって、スペクトル距離を算出する対象を中央画像スペクトルに限定することができ、処理速度の向上と、安定したパターンＩＤ決定精度を得ることができる。
【００６９】
続いて、領域Ｃ内の各ピクセルに対して決定された最も適合度の高いパターンＩＤのうち、出現頻度が閾値以上のパターンＩＤを、その圃場内で栽培されている可能性の高い作物として、パターンＩＤリストに設定する（Ｓ１９０６）。
【００７０】
この際、最も出現頻度が高いパターンＩＤの出現頻度が別の閾値以上であった場合は、他のパターンＩＤの出現頻度に関わらずパターンＩＤリストに設定するパターンＩＤをこのパターンＩＤに設定することもできる。また、いずれのパターンＩＤの出現頻度も閾値を越えなかった場合、出現頻度の順にあらかじめ決めた個数のパターンＩＤをパターンＩＤリストに登録することができる。いずれにせよ、Ｓ１９０６が終了した時点で、パターンＩＤリストには、その圃場内で栽培されている可能性が高い作物のパターンＩＤが格納されていることになる。
【００７１】
次に、ポリゴン内で、ポリゴン境界から閾値以下の距離の領域（領域Ｄとする）の全ピクセルに対して、パターンＩＤリストに含まれる全パターンＩＤの、全画像スペクトルを用いて距離を算出し、各ピクセルに最も近いパターンＩＤと擬似画像の種類（中央・1/3・2/3・土壌）の組を決定する（Ｓ１９０７）。
【００７２】
最後に、領域Ｃ，領域Ｄを合わせて最も出現頻度の高いパターンＩＤを決定し、領域Ｃの全ピクセルに対して決定されたパターンＩＤの中央画像を対応付け，領域Ｄで決定されたパターンＩＤ以外が選択されている全ピクセルに対して、決定されたパターンＩＤの全画像スペクトルとのスペクトル距離を算出し、最も距離の短い擬似画像を対応付ける（Ｓ１９０８）。
これによって、圃場ポリゴン内の全てのピクセルがただ一つのパターンＩＤに対応付けられ、各ピクセルは、中央・1/3・2/3・土壌のいずれかの画像スペクトルに対応づけられる。
【００７３】
一方、Ｓ１９０４において、解析対象となる画像の解像度が閾値未満の場合は、作物配置を考慮したテンプレートマッチングを行う。
まず、格子フラグによって用いるテンプレートが変わるため、処理を振り分ける（Ｓ１９０９）。
【００７４】
格子フラグがｔｒｕｅの場合は、格子画像をテンプレートとして用いる。ポリゴン内で、ポリゴン境界から閾値以上距離がある領域（領域Ａとする）に対して、パターンＩＤリストに含まれる全パターンＩＤの中央画像を用いてテンプレートマッチングを行い、最も適合度の高い中央画像のパターンＩＤ，位置、角度を得る（Ｓ１９１０）。ここでも、植付けが一様な領域に中央画像を適用することで、処理量の低減と、処理結果の安定化という効果が得られる。
【００７５】
次に最も出現頻度の高いパターンＩＤを選択し、ポリゴン内でポリゴン境界から閾値以下の距離の領域（領域Ｂとする）に対して、選択されたパターンＩＤの全画像（中央・端・角）を用いてテンプレートマッチングを行い、最も高い擬似画像、位置、角度を得る（Ｓ１９１１）。ここでは最も出現頻度の高いパターンＩＤだけを利用して領域Ｂに対するテンプレートマッチングを行ったが、Ｓ１９０６からＳ１９０８で説明したように、出現頻度上位の複数のパターンＩＤに対する中央・端・角の画像を用いてテンプレートマッチングを行い、両機Ａ，Ｂを合わせて出現頻度の多いパターンＩＤを選択することもできる。
【００７６】
格子フラグがｆａｌｓｅの場合は、利用するテンプレートが直列画像である以外は、Ｓ１９１２、Ｓ１９１３ともに同様の処理となる。ただし、直列画像はパターンＩＤや成長の時期によってサイズが異なるため、得られた正規化相互相関の値をピクセル単位に正規化した上で、小さいテンプレートがマッチングしやすくなる傾向を防ぐために、テンプレートサイズによるペナルティ項を足すなどの方法で、サイズの差に影響されないテンプレートマッチングが実現できる。
【００７７】
Ｓ１９１０もしくはＳ１９１２のステップで得られた角度からヒストグラムを得て、そのモードを対象圃場の植付け角度とすることができる。この場合、Ｓ１９１１、Ｓ１９１３の処理の際に、角度を固定して探索することができ、処理量の低減と、処理結果の安定化という効果が得られる。ここでは最も出現頻度の高いパターンＩＤだけを利用して領域Ｂに対するテンプレートマッチングを行ったが、Ｓ１９０６からＳ１９０８で説明したように、出現頻度上位の複数のパターンＩＤに対する中央・端・角の画像を用いてテンプレートマッチングを行い、両機Ａ，Ｂを合わせて出現頻度の多いパターンＩＤを選択することもできる。
【００７８】
一方、Ｓ１９０１において、圃場情報にパターンＩＤが無い場合は、パターンＩＤリストの格子フラグを全てｆａｌｓｅに設定した後で、Ｓ１９０４移行のステップを実行する方法を用いることができる。これによって、解像度が閾値未満だった場合に、当初の格子フラグの値によらず、全てのパターンＩＤに対して、直列画像を用いてテンプレートマッチングを行うことになり、テンプレート間の適合度の比較が容易になる。Ｓ１７０６の後に、Ｓ１７０７以降のステップを実行して直列画像を生成したのはこのためである。
【００７９】
Ｓ１４０３を実行した段階で、各ポリゴン領域に対してただ１つのパターンＩＤが対応付けられているため、その領域での栽培品種や作付方法、作付時期等の情報が得られる。圃場情報がある場合は、圃場ポリゴン内の各領域に対応するパターンＩＤは同一で、作付時期や画像の種別が異なっている可能性がある。例えば作付時期の最頻値を求めることで、圃場ごとの作付時期を推定することができる。また、疑似画像をテンプレートマッチングした際の疑似画像の角度に関してもそれぞれ異なった値が得られるが、例えば角度の最頻値を求めることで、圃場の作付方向を算出することができる。
以上によって、各ポリゴン領域ごとに圃場特徴量を算出することができる。
【００８０】
最後に図２３を用いて、Ｓ１４０４の詳細フローについて説明する。Ｓ１４０４が実行される時点で、特定の領域に対して、ただ１つのテンプレートが選択されている。
まず解析対象の解像度で処理を振り分ける（Ｓ２３０１）。
【００８１】
解析対象の解像度が閾値未満の場合はテンプレートは１ピクセル画像であるためミクセル分解を行う。この際、選択されているテンプレートに対応した植物スペクトルと土壌スペクトルを初期エンドメンバとし、被覆率をこれらのエンドメンバの初期混合率として、ミクセル分解を行う（Ｓ２３０３）。上記の被覆率は、成長による変化だけでなく、圃場周辺での被覆率の変化も考慮しているため、初期値として適している。ミクセル分解の手法としては、非負値行列分解による方法を用いることができる。この手法は、局所最適解に収束することが知られているため、初期値の与え方によって現実とは異なる解が得られる可能性があるが、本実施例では上記の初期値を用いることによって、より最適な解に収束させることが可能となる。
【００８２】
解像度が閾値未満だった場合は、擬似画像の植物を表す円の内側のピクセルから、植物のピクセルの抽出を行う（Ｓ２３０３）。抽出方法としては、円内のピクセルを教師無し分類して最大のクラスを植物として抽出する方法や、テンプレートに対応する植物スペクトルとの距離が閾値以下のピクセルのみを植物として抽出する方法を用いることができる。前者は通常教師無し分類では分離が難しい植物のピクセルを対象範囲を限定することで実現できる。後者は、比較対照としてテンプレートに対応する植物スペクトルを用いることで、教師無し分類よりさらに精度良く、植物のピクセルを抽出することができる。
【００８３】
最後に、抽出された特徴量からＮＤＶＩ等の特徴量を計算する（Ｓ２３０４）。
以上によって、本発明の実施例によって、解析対象の画像から、精度良く植物の特徴量が抽出できることが示せた。
【００８４】
なお、本発明は上記した実施例に限定されるものではなく、様々な変形例が含まれる。例えば、上記した実施例は本発明を分かりやすく説明するために詳細に説明したのであり、必ずしも説明の全ての構成を備えるものに限定されものではない。また、実施例の構成の一部について、他の構成の追加・削除・置換をすることが可能である。
【００８５】
また、上記の各構成、機能、処理部、処理手段等は、それらの一部又は全部を、例えば集積回路で設計する等によりハードウェアで実現してもよい。また、上記の各構成、機能等は、それぞれの機能を実現するプログラムを実行することによりソフトウェアで実現する場合を例示して説明したが、各機能を実現するプログラム、テーブル、ファイル等の情報はメモリのみならず、ハードディスク、ＳＳＤ（Solid State Drive）等の記録装置、または、ＩＣカード、ＳＤカード、ＤＶＤ等の記録媒体におくことができるし、必要に応じてネットワーク等を介してダウンロード、インストールすることも可能である。
【産業上の利用可能性】
【００８６】
本発明によれば、解析対象の画像から対象とする作物の特徴量を精度良く抽出することが可能であり、農業を行う圃場を撮影した画像を高精度に解析する画像解析技術として極めて有用である。
【符号の説明】
【００８７】
１００    画像解析サーバ
１０１    端末
１０２    ネットワーク
１０３    データ通信バス
１０４    制御部
１０５    入力部
１０６    表示部
１０７    通信部
１０８    植生ＤＢ編集部
１０９    植生ＤＢ
１１０    パターンリスト生成部
１１１    擬似画像生成部
１１２    最適テンプレート選択部
１１３    特徴量抽出部
１１４    画像ＤＢ
１１５    制御部
１１６    入力部
１１７    表示部
１１８    通信部
２０１    エリアテーブル
２０２    地域作付テーブル
２０３    作付テーブル
２０４    作物テーブル
２０５    スペクトルテーブル
２０６    土壌スペクトルテーブル
２０７    バンドセットテーブル
２０８    成長テーブル
２０９    圃場ポリゴンテーブル
２１０    圃場作付テーブル
７００    スペクトルパターン
【図１】
0005886054.tif 000002
【図２】
0005886054.tif 000003
【図３】
0005886054.tif 000004
【図４】
0005886054.tif 000005
【図５】
0005886054.tif 000006
【図６】
0005886054.tif 000007
【図７】
0005886054.tif 000008
【図８】
0005886054.tif 000009
【図９】
0005886054.tif 000010
【図１０】
0005886054.tif 000011
【図１１】
0005886054.tif 000012
【図１２】
0005886054.tif 000013
【図１３】
0005886054.tif 000014
【図１４】
0005886054.tif 000015
【図１５】
0005886054.tif 000016
【図１６】
0005886054.tif 000017
【図１７】
0005886054.tif 000018
【図１８】
0005886054.tif 000019
【図１９】
0005886054.tif 000020
【図２０】
0005886054.tif 000021
【図２１】
0005886054.tif 000022
【図２２】
0005886054.tif 000023
【図２３】
0005886054.tif 000024
【図２４】
0005886054.tif 000025