(57)【要約】【課題】異なる時刻に撮影された２つの画像間の位置ずれ量を効率よく自動的に求めることを目的とする。【解決手段】マッチング処理部２０２は、異なる時刻に撮影された２つの画像（例えば、地震が発生する直前及び直後の画像である発災前画像及び発災後画像）をメモリ２５２から読み取り、一方の画像（例えば、発災後画像）にて所定数の画素からなる第１の領域を１画素単位で移動し、第１の領域を移動する度に、第１の領域と、他方の画像（例えば、発災前画像）にて第１の領域と同じ位置にある第２の領域との間の相関係数をプロセッサ２５１により算出する。位置ずれ量算出処理部２０３は、マッチング処理部２０２が算出した相関係数に基づいて、上記２つの画像間の位置ずれ量を１画素未満の単位でプロセッサ２５１により推定する。【選択図】図２
2009295062.tif 000002
【特許請求の範囲】
【請求項１】
異なる時刻に撮影された２つの画像を予めメモリに記憶する画像処理装置であって、
前記２つの画像をメモリから読み取り、一方の画像にて所定数の画素からなる第１の領域を１画素以上の単位で移動し、前記第１の領域を移動する度に、前記第１の領域と、他方の画像にて前記第１の領域と同じ位置にある第２の領域との間の相関係数をプロセッサにより算出するマッチング処理部と、
前記マッチング処理部が算出した相関係数に基づいて、前記２つの画像間の位置ずれ量を１画素未満の単位でプロセッサにより推定する位置ずれ量算出処理部とを備えることを特徴とする画像処理装置。
【請求項２】
前記マッチング処理部は、前記第１の領域を移動する度に、前記第１の領域の位置と、前記第１の領域と前記第２の領域との間の相関係数とを対応付けてメモリに記憶し、
前記位置ずれ量算出処理部は、前記第１の領域の位置と前記相関係数とをメモリから読み取り、前記第１の領域の位置と前記相関係数との関係を示す近似式をプロセッサにより算出し、前記近似式を用いて、前記２つの画像間の位置ずれ量を１画素未満の単位で推定することを特徴とする請求項１に記載の画像処理装置。
【請求項３】
前記画像処理装置は、さらに、
前記一方の画像を所定数の画素からなる複数の領域に分割し、各領域における画素値の分散とエッジ量との少なくともいずれかに基づいて、前記第１の領域を選定するマッチング小領域選定処理部を備えることを特徴とする請求項１又は２に記載の画像処理装置。
【請求項４】
前記画像処理装置は、さらに、
前記一方の画像を前記位置ずれ量算出処理部が推定した位置ずれ量ずらした画像と前記他方の画像との特徴量をプロセッサにより比較し、当該比較結果に基づいて、前記２つの画像間の変化を抽出する変化抽出処理部を備えることを特徴とする請求項１から３までのいずれかに記載の画像処理装置。
【請求項５】
前記２つの画像は、上空から撮影された地表面の画像であり、
前記画像処理装置は、さらに、
前記一方の画像を前記位置ずれ量算出処理部が推定した位置ずれ量ずらした画像と前記他方の画像とで異なる位置に映っている物体をプロセッサにより推定し、前記異なる位置の間の距離が所定の閾値以上である場合、前記物体が雲であると判断する雲抽出処理部を備えることを特徴とする請求項１から４までのいずれかに記載の画像処理装置。
【請求項６】
前記２つの画像は、上空から撮影された地表面の画像であり、
前記画像処理装置は、さらに、
前記一方の画像を前記位置ずれ量算出処理部が推定した位置ずれ量ずらした画像と前記他方の画像との少なくともいずれかの各画素の赤波長値とＮＤＶＩ（Ｎｏｒｍａｌｉｚｅｄ・Ｄｉｆｆｅｒｅｎｃｅ・Ｖｅｇｅｔａｔｉｏｎ・Ｉｎｄｅｘ）値とＳ（Ｓａｔｕｒａｔｉｏｎ）値とＩ（Ｉｎｔｅｎｓｉｔｙ）値とをプロセッサにより算出し、各画素の赤波長値とＮＤＶＩ値とＳ値とＩ値とに基づいて、各画素が植物と水域と陸地と山肌と雲とのいずれであるか判断するテクスチャ分類処理部を備えることを特徴とする請求項１から５までのいずれかに記載の画像処理装置。
【請求項７】
異なる時刻に撮影された２つの画像を予めメモリに記憶する画像処理装置を用いた画像処理方法であって、
前記画像処理装置のマッチング処理部が、前記２つの画像をメモリから読み取り、一方の画像にて所定数の画素からなる第１の領域を１画素以上の単位で移動し、前記第１の領域を移動する度に、前記第１の領域と、他方の画像にて前記第１の領域と同じ位置にある第２の領域との間の相関係数をプロセッサにより算出する第１ステップと、
前記画像処理装置の位置ずれ量算出処理部が、前記第１ステップで算出された相関係数に基づいて、前記２つの画像間の位置ずれ量を１画素未満の単位でプロセッサにより推定する第２ステップとを備えることを特徴とする画像処理方法。
【請求項８】
異なる時刻に撮影された２つの画像を予めメモリに記憶するコンピュータで実行される画像処理プログラムであって、
前記２つの画像をメモリから読み取り、一方の画像にて所定数の画素からなる第１の領域を１画素以上の単位で移動し、前記第１の領域を移動する度に、前記第１の領域と、他方の画像にて前記第１の領域と同じ位置にある第２の領域との間の相関係数をプロセッサにより算出するマッチング処理と、
前記マッチング処理が算出した相関係数に基づいて、前記２つの画像間の位置ずれ量を１画素未満の単位でプロセッサにより推定する位置ずれ量算出処理とをコンピュータに実行させることを特徴とする画像処理プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、画像処理装置及び画像処理方法及び画像処理プログラムに関するものである。本発明は、特に、人工衛星や航空機から撮像された複数の画像データを解析し、画像の変化を抽出する技術（変化抽出システム）に関するものである。
【背景技術】
【０００２】
特許文献１又は特許文献２又は特許文献３には、異なる時刻に撮像された２つの画像データ（ペア画像）から、地表面の変化を抽出する方法及び装置が開示されている。
【特許文献１】特開２００５−７８５６６号公報
【特許文献２】特開２００５−２４１８８６号公報
【特許文献３】特開２００７−１１４９９４号公報
【発明の開示】
【発明が解決しようとする課題】
【０００３】
画像の変化を抽出する処理は、人手に頼らないで自動的に行えることが望ましい。しかし、特許文献１の装置では、ペア画像の間で季節や時刻といった時期的な撮像条件が異なる場合には、事前にどちらか一方の画像をもう一方の画像の撮像条件に合わせるような補正をしておく必要があった。さらに、特許文献１又は特許文献２の方法及び装置では、ペア画像の位置合わせ（レジストレーション）処理を実施する手段を具備していないため、事前にレジストレーションを手動で行っておく必要があった。さらに、特許文献１又は特許文献２又は特許文献３の装置では、画像から雲の抽出処理を実施する手段を具備していないため、抽出された変化の原因が雲であるか否かを人間が判読する必要があった。さらに、特許文献１又は特許文献２又は特許文献３の装置では、テクスチャの分類処理を実施する手段を具備していないため、どのような地表面で変化があったのかを人間が判読する必要があった。即ち、特許文献１又は特許文献２又は特許文献３の装置では、人手による作業を必要とするため、自動的に変化を抽出することが困難であった。
【０００４】
本発明は、例えば、異なる時刻に撮影された２つの画像間の位置ずれ量を効率よく自動的に求めることを目的とする。また、本発明は、例えば、人手による作業を必要とせずに、異なる時刻に撮影された２つの画像間の変化を迅速に抽出することを目的とする。
【課題を解決するための手段】
【０００５】
本発明の一の態様に係る画像処理装置は、
異なる時刻に撮影された２つの画像を予めメモリに記憶する画像処理装置であって、
前記２つの画像をメモリから読み取り、一方の画像にて所定数の画素からなる第１の領域を１画素以上の単位で移動し、前記第１の領域を移動する度に、前記第１の領域と、他方の画像にて前記第１の領域と同じ位置にある第２の領域との間の相関係数をプロセッサにより算出するマッチング処理部と、
前記マッチング処理部が算出した相関係数に基づいて、前記２つの画像間の位置ずれ量を１画素未満の単位でプロセッサにより推定する位置ずれ量算出処理部とを備えることを特徴とする。
【発明の効果】
【０００６】
本発明の一の態様によれば、画像処理装置において、マッチング処理部が、異なる時刻に撮影された２つの画像のうち、一方の画像にて第１の領域を１画素以上の単位で移動し、前記第１の領域を移動する度に、前記第１の領域と、他方の画像にて前記第１の領域と同じ位置にある第２の領域との間の相関係数を算出し、位置ずれ量算出処理部が、前記相関係数に基づいて、前記２つの画像間の位置ずれ量を１画素未満の単位で推定することにより、前記２つの画像間の位置ずれ量を効率よく自動的に求めることが可能となる。
【発明を実施するための最良の形態】
【０００７】
以下、本発明の実施の形態について、図を用いて説明する。
【０００８】
実施の形態１．
図１は、本実施の形態に係る変化抽出システム１００の構成の一例を示すブロック図である。
【０００９】
変化抽出システム１００は、静止光学衛星１０１（飛翔体の一例）と衛星画像受信センタ１０２と衛星画像データベース１０３と画像処理装置１０４と端末装置１０５とから構成される。
【００１０】
静止光学衛星１０１は、赤道上空の高度約３６，０００ｋｍの円軌道を、地球の自転と同じ周期で公転し、光学センサ（撮像装置の一例）によって地表面の画像（衛星画像）を撮影する人工衛星である。静止光学衛星１０１は、撮影した衛星画像を、無線通信により衛星画像受信センタ１０２へ送信する。なお、静止光学衛星１０１は、周回衛星、航空機など、他の飛翔体で置き換えてもよい。また、光学センサは、合成開口レーダなど、他の撮像装置で置き換えてもよい。
【００１１】
衛星画像受信センタ１０２は、通信用アンテナなどの装置を具備し、２４時間体制で静止光学衛星１０１と交信する設備である。衛星画像受信センタ１０２は、静止光学衛星１０１より受信した衛星画像を、衛星画像データベース１０３に格納する。
【００１２】
衛星画像データベース１０３は、３〜４ヶ月分の衛星画像の格納が可能な大容量のハードディスクドライブ（ＨＤＤ）などを搭載したファイルサーバである。衛星画像データベース１０３は、衛星画像受信センタ１０２及び画像処理装置１０４及び端末装置１０５と、ネットワーク１０６を介して接続されている。
【００１３】
画像処理装置１０４は、ＣＰＵ（Ｃｅｎｔｒａｌ・Ｐｒｏｃｅｓｓｉｎｇ・Ｕｎｉｔ）やメモリなどを搭載し、端末装置１０５からの処理要求に応答するアプリケーションサーバである。画像処理装置１０４は、その処理機能として、レジストレーション処理部１４１と、変化抽出処理部１４２と、雲抽出処理部１４３と、テクスチャ分類処理部１４４とを実装する。なお、図では画像処理装置１０４を１台のみ示しているが、複数台あってもよい。
【００１４】
端末装置１０５は、モニタやキーボードやマウスなどを具備したワークステーションであり、画像処理装置１０４に対して処理（サービス）の依頼を行い、その提供を受けるクライアントである。なお、図では端末装置１０５を１台のみ示しているが、複数台あってもよい。
【００１５】
このように、本実施の形態では、人工衛星や航空機から撮像された複数の画像データを解析し、画像の変化を抽出する変化抽出システム１００にて、画像処理装置１０４が画像のレジストレーション処理と変化抽出処理と雲抽出処理とテクスチャ分類処理とを自動的に実施することで、画像判読者の負担を軽減することができる。
【００１６】
図２は、画像処理装置１０４の構成の一例を示すブロック図である。
【００１７】
前述したように、画像処理装置１０４は、レジストレーション処理部１４１、変化抽出処理部１４２、雲抽出処理部１４３、テクスチャ分類処理部１４４を備える。また、画像処理装置１０４は、プロセッサ２５１、メモリ２５２、入力装置２５３、出力装置２５４などのハードウェアを備える。ハードウェアは画像処理装置１０４の各部によって利用される。画像処理装置１０４は、異なる時刻に撮影された少なくとも２つの画像を予めメモリ２５２に記憶している。画像処理装置１０４は、具体的には、上空から撮影された地表面の画像、即ち、前述した衛星画像をメモリ２５２に記憶している。
【００１８】
レジストレーション処理部１４１は、衛星画像のうち、異なる時刻に撮像された２つの画像データ（ペア画像）をメモリ２５２から読み取り、プロセッサ２５１を用いて自動的に画像の位置合わせをする。後述するように、レジストレーション処理部１４１は、マッチングのための小領域を選定する指標として、分散値とエッジ量とを使用する。レジストレーション処理部１４１は、マッチング小領域選定処理部２０１、マッチング処理部２０２、位置ずれ量算出処理部２０３、リサンプリング処理部２０４を備える。各部の動作については後述する。
【００１９】
変化抽出処理部１４２は、ペア画像における変化を、プロセッサ２５１を用いて自動的に抽出する。後述するように、変化抽出処理部１４２は、変化抽出の指標として、分散値とエッジ量と相関係数とを使用する。変化抽出処理部１４２は、特徴量算出処理部２０５を備える。特徴量算出処理部２０５の動作については後述する。
【００２０】
雲抽出処理部１４３は、ペア画像における雲の移動ベクトル情報から、ペア画像に映っている物体が雲であるか否かを、プロセッサ２５１を用いて自動的に判定する。後述するように、雲抽出処理部１４３は、雲の移動ベクトルから雲であるか否かを判定する際、ベクトルの信頼性の指標として、分散値とエッジ量と相関係数とを使用する。雲抽出処理部１４３は、移動ベクトル抽出処理部２０６、移動量処理部２０７、品質指標算出処理部２０８を備える。各部の動作については後述する。
【００２１】
テクスチャ分類処理部１４４は、画像の地表面を、植物及び水域及び陸地及び山肌及び雲のいずれかに、プロセッサ２５１を用いて自動的に分類する。後述するように、テクスチャ分類処理部１４４は、植物及び水域及び陸地及び山肌及び雲のいずれかにテクスチャ分類するための指標として、赤波長の値とＮＤＶＩ（Ｎｏｒｍａｌｉｚｅｄ・Ｄｉｆｆｅｒｅｎｃｅ・Ｖｅｇｅｔａｔｉｏｎ・Ｉｎｄｅｘ）値とＨＳＩ（Ｈｕｅ・Ｓａｔｕｒａｔｉｏｎ・Ｉｎｔｅｎｓｉｔｙ）変換後の彩度（Ｓａｔｕｒａｔｉｏｎ）と明度（Ｉｎｔｅｎｓｉｔｙ）とを使用する。テクスチャ分類処理部１４４は、ＮＤＶＩ算出処理部２０９、ＨＳＩ変換処理部２１０を備える。各部の動作については後述する。
【００２２】
図３は、画像処理装置１０４のハードウェア構成の一例を示す図である。
【００２３】
図３において、画像処理装置１０４は、コンピュータであり、ＣＲＴ（Ｃａｔｈｏｄｅ・Ｒａｙ・Ｔｕｂｅ）やＬＣＤ（液晶ディスプレイ）の表示画面を有する表示装置９０１、キーボード９０２（Ｋ／Ｂ）、マウス９０３、ＦＤＤ９０４（Ｆｌｅｘｉｂｌｅ・Ｄｉｓｋ・Ｄｒｉｖｅ）、ＣＤＤ９０５（Ｃｏｍｐａｃｔ・Ｄｉｓｃ・Ｄｒｉｖｅ）、プリンタ装置９０６などのハードウェアを備え、これらはケーブルや信号線で接続されている。
【００２４】
画像処理装置１０４は、プログラムを実行するＣＰＵ９１１を備えている。ＣＰＵ９１１は、プロセッサ２５１の一例である。ＣＰＵ９１１は、バス９１２を介してＲＯＭ９１３（Ｒｅａｄ・Ｏｎｌｙ・Ｍｅｍｏｒｙ）、ＲＡＭ９１４（Ｒａｎｄｏｍ・Ａｃｃｅｓｓ・Ｍｅｍｏｒｙ）、通信ボード９１５、表示装置９０１、キーボード９０２、マウス９０３、ＦＤＤ９０４、ＣＤＤ９０５、プリンタ装置９０６、磁気ディスク装置９２０と接続され、これらのハードウェアデバイスを制御する。磁気ディスク装置９２０の代わりに、光ディスク装置、メモリカードリーダライタなどの記憶媒体が用いられてもよい。
【００２５】
ＲＡＭ９１４は、揮発性メモリの一例である。ＲＯＭ９１３、ＦＤＤ９０４、ＣＤＤ９０５、磁気ディスク装置９２０の記憶媒体は、不揮発性メモリの一例である。これらは、メモリ２５２の一例である。通信ボード９１５、キーボード９０２、マウス９０３、ＦＤＤ９０４、ＣＤＤ９０５などは、入力装置２５３の一例である。また、通信ボード９１５、表示装置９０１、プリンタ装置９０６などは、出力装置２５４の一例である。
【００２６】
通信ボード９１５は、ＬＡＮ（ローカルエリアネットワーク）などに接続されている。通信ボード９１５は、ＬＡＮに限らず、インターネット、あるいは、ＩＰ−ＶＰＮ（Ｉｎｔｅｒｎｅｔ・Ｐｒｏｔｏｃｏｌ・Ｖｉｒｔｕａｌ・Ｐｒｉｖａｔｅ・Ｎｅｔｗｏｒｋ）、広域ＬＡＮ、ＡＴＭ（Ａｓｙｎｃｈｒｏｎｏｕｓ・Ｔｒａｎｓｆｅｒ・Ｍｏｄｅ）ネットワークなどのＷＡＮ（ワイドエリアネットワーク）などに接続されていても構わない。ＬＡＮ、インターネット、ＷＡＮは、ネットワーク１０６の一例である。
【００２７】
磁気ディスク装置９２０には、オペレーティングシステム９２１（ＯＳ）、ウィンドウシステム９２２、プログラム群９２３、ファイル群９２４が記憶されている。プログラム群９２３のプログラムは、ＣＰＵ９１１、オペレーティングシステム９２１、ウィンドウシステム９２２により実行される。プログラム群９２３には、本実施の形態の説明において「〜部」として説明する機能を実行するプログラムが記憶されている。プログラムは、ＣＰＵ９１１により読み出され実行される。また、ファイル群９２４には、本実施の形態の説明において、「〜データ」、「〜情報」、「〜ＩＤ（識別子）」、「〜フラグ」、「〜結果」として説明するデータや情報や信号値や変数値やパラメータが、「〜ファイル」や「〜データベース」や「〜テーブル」の各項目として記憶されている。「〜ファイル」や「〜データベース」や「〜テーブル」は、ディスクやメモリなどの記憶媒体に記憶される。ディスクやメモリなどの記憶媒体に記憶されたデータや情報や信号値や変数値やパラメータは、読み書き回路を介してＣＰＵ９１１によりメインメモリやキャッシュメモリに読み出され、抽出・検索・参照・比較・演算・計算・制御・出力・印刷・表示などのＣＰＵ９１１の処理（動作）に用いられる。抽出・検索・参照・比較・演算・計算・制御・出力・印刷・表示などのＣＰＵ９１１の処理中、データや情報や信号値や変数値やパラメータは、メインメモリやキャッシュメモリやバッファメモリに一時的に記憶される。
【００２８】
また、本実施の形態の説明において用いるブロック図やフローチャートの矢印の部分は主としてデータや信号の入出力を示し、データや信号は、ＲＡＭ９１４などのメモリ、ＦＤＤ９０４のフレキシブルディスク（ＦＤ）、ＣＤＤ９０５のコンパクトディスク（ＣＤ）、磁気ディスク装置９２０の磁気ディスク、その他光ディスク、ミニディスク（ＭＤ）、ＤＶＤ（Ｄｉｇｉｔａｌ・Ｖｅｒｓａｔｉｌｅ・Ｄｉｓｃ）などの記録媒体に記録される。また、データや信号は、バス９１２や信号線やケーブルその他の伝送媒体により伝送される。
【００２９】
また、本実施の形態の説明において「〜部」として説明するものは、「〜回路」、「〜装置」、「〜機器」であってもよく、また、「〜ステップ」、「〜工程」、「〜手順」、「〜処理」であってもよい。即ち、「〜部」として説明するものは、ＲＯＭ９１３に記憶されたファームウェアで実現されていても構わない。あるいは、ソフトウェアのみ、あるいは、素子・デバイス・基板・配線などのハードウェアのみ、あるいは、ソフトウェアとハードウェアとの組み合わせ、さらには、ファームウェアとの組み合わせで実現されていても構わない。ファームウェアとソフトウェアは、プログラムとして、磁気ディスク、フレキシブルディスク、光ディスク、コンパクトディスク、ミニディスク、ＤＶＤなどの記録媒体に記憶される。このプログラムはＣＰＵ９１１により読み出され、ＣＰＵ９１１により実行される。即ち、プログラムは、本実施の形態の説明で述べる「〜部」としてコンピュータを機能させるものである。あるいは、本実施の形態の説明で述べる「〜部」の手順や方法をコンピュータに実行させるものである。
【００３０】
以下、変化抽出システム１００の動作を、地震が発生した直後の例を挙げ説明する。
【００３１】
地震発生の一報を受けた画像判読者は、端末装置１０５のキーボード及びマウスなどを操作して発災時刻と推定される震源の位置とを入力し、画像処理装置１０４に対して処理の依頼を行う。
【００３２】
次に、画像処理装置１０４は入力された情報を基に、衛星画像データベース１０３から発災前後のペア画像を、ネットワーク１０６を介して取得し、レジストレーション処理と、変化抽出処理と、雲抽出処理と、テクスチャ分類処理とを実施し、処理結果を端末装置１０５へ提供する。それぞれの処理内容については、後述する。
【００３３】
次に、端末装置１０５は、画像処理装置１０４より提供された処理結果をモニタに表示する。ここで、処理結果とは、ペア画像のレジストレーション処理結果と変化抽出処理結果と雲抽出処理結果とテクスチャ分類処理結果である。なお、処理結果をモニタに表示するときは、画像をセル画のように重ねて表示（レイヤ表示）する。
【００３４】
最後に、画像判読者はモニタに表示された処理結果を基に、被災地及び被災状況を確認する。例えば、変化抽出結果で抽出された場所が何処であるかを確認するためには、レイヤをレジストレーション結果に重ねて表示する。また、変化抽出結果で抽出された原因が雲であるかどうかを確認するためには、レイヤを雲抽出結果に重ねて表示する。また、変化抽出結果で抽出された場所がどのような地表面であるかを確認するためには、レイヤをテクスチャ分類結果に重ねて表示する。そして、画像判読者は、確認した被災地及び被災状況を災害対策部門に報告する。必要であれば、救急隊の出動を要請する。
【００３５】
このようにして、災害の発生直後に変化抽出システム１００を利用し、被災地及び被災状況を早期に特定することにより、迅速な救助活動が可能となる。
【００３６】
なお、上述の例では地震発生の場合について説明したが、津波などの他の災害についても同様である。
【００３７】
以下では、画像処理装置１０４のレジストレーション処理部１４１の動作、即ち、レジストレーション処理について、内容を説明する。処理の概要を示すフローチャートを図４に示す。図４に示す通り、レジストレーション処理部１４１は、ペア画像を入力装置２５３により入力し、メモリ２５２に記憶する。その後、レジストレーション処理部１４１は、メモリ２５２に記憶したペア画像に対し、マッチング小領域選定処理（ステップＳ４１１）、マッチング処理（ステップＳ４１２）、位置ずれ量算出処理（ステップＳ４１３）、リサンプリング処理（ステップＳ４１４）をプロセッサ２５１により実施する。そして、レジストレーション処理部１４１は、位置合わせしたペア画像をメモリ２５２に記憶するとともに、出力装置２５４により出力する。詳細を以下に示す。
【００３８】
まず、レジストレーション処理部１４１のマッチング小領域選定処理部２０１は、マッチング小領域選定処理（ステップＳ４１１）を実施する。なお、この処理は入力されたペア画像のどちらか一方に対してのみ行う。ここでは、撮像された時刻が発災後の画像（発災後画像）に対して処理を実施するものとする。その具体的な手順は次の通りである。まず、マッチング小領域選定処理部２０１は、発災後画像をＮ×Ｎ（ＮはＮ＞０となる整数）の小領域に分割する。次に、マッチング小領域選定処理部２０１は、分割したそれぞれの小領域内で分散値（画素値の分散）及びエッジ量を下記の式（１）及び式（２）より算出する。最後に、マッチング小領域選定処理部２０１は、分散値及び／又はエッジ量が最も大きい小領域をマッチングのための小領域として選定する。これにより、後述する処理において、最も特徴的な小領域について相関係数分布が求められることになり、レジストレーション処理の精度が向上する。
【００３９】
なお、分散値は次式で求められる。
【００４０】
【数１】
2009295062.tif 000003
【００４１】
なお、エッジ量は次式で求められる。
【００４２】
【数２】
2009295062.tif 000004
【００４３】
即ち、ＳＯＢＥＬなどのエッジ抽出オペレータによりエッジ抽出した小領域画像の輝度値を合計したものをエッジ量とする。
【００４４】
上記のように、マッチング小領域選定処理部２０１は、異なる時刻に撮影された２つの画像（例えば、発災前画像及び発災後画像）をメモリ２５２から読み取り、一方の画像（例えば、発災後画像）を所定数の画素からなる複数（例えば、Ｎ×Ｎ）の領域に分割し、各領域における画素値の分散（分散値）とエッジ量との少なくともいずれかに基づいて、マッチング処理部２０２の処理対象とする第１の領域（例えば、分散値及び／又はエッジ量が最も大きい領域）を選定する。
【００４５】
次に、レジストレーション処理部１４１のマッチング処理部２０２は、マッチング処理（ステップＳ４１２）を実施する。この処理は、正規化相関法（面積相関法）を利用して行う。正規化相関法は、２枚の画像における同一地点（対応点）の検出のために用いられる手法で、両画像から取り出した小領域の一致の程度を定量的に表すものであり、求めた相関係数の値が大きいほど２つの小領域が一致していることを示唆する。その具体的な手順は次の通りである。まず、マッチング処理部２０２は、マッチング小領域選定処理（ステップＳ４１１）にて選定されたマッチングのための小領域を、発災前画像に重ね合わせる。次に、マッチング処理部２０２は、正規化相関法により相関係数を下記の式（３）より算出する。以降、マッチング処理部２０２は、重ね合わせの位置を上下左右に少なくとも１画素ずつ移動させ、それぞれの位置で相関係数を算出して相関係数分布（相関係数マップ）を作成する。
【００４６】
なお、正規化相関法の相関係数は次式により求められる。
【００４７】
【数３】
2009295062.tif 000005
【００４８】
次に、レジストレーション処理部１４１の位置ずれ量算出処理部２０３は、位置ずれ量算出処理（ステップＳ４１３）を実施する。この処理は、ペア画像間の位置ずれの程度を画素単位以下の精度で算出する処理である。その具体的な手順は次の通りである。まず、位置ずれ量算出処理部２０３は、マッチング処理（ステップＳ４１２）にて作成された相関係数分布から、相関係数の最大値を抽出する。次に、位置ずれ量算出処理部２０３は、最大値の両隣の相関係数を抽出する。最後に、位置ずれ量算出処理部２０３は、抽出した３点を２次多項式（２次関数）に近似し、２次多項式のピーク位置（２次関数の放物線の頂点）を求める。このピーク位置の座標を、発災前画像と発災後画像との位置ずれ量とする。なお、位置ずれ量はＸ方向（横方向）及びＹ方向（縦方向）それぞれにおいて算出する。
【００４９】
上記のように、マッチング処理部２０２は、異なる時刻に撮影された２つの画像（例えば、発災前画像及び発災後画像）をメモリ２５２から読み取り、一方の画像（例えば、発災後画像）にて所定数の画素からなる第１の領域（例えば、マッチング小領域選定処理部２０１が選定した領域）を１画素以上の単位（好ましくは、１画素単位）で移動し、第１の領域を移動する度に、第１の領域と、他方の画像（例えば、発災前画像）にて第１の領域と同じ位置にある第２の領域との間の相関係数をプロセッサ２５１により算出する。位置ずれ量算出処理部２０３は、マッチング処理部２０２が算出した相関係数に基づいて、上記２つの画像間の位置ずれ量を１画素未満の単位でプロセッサ２５１により推定する。
【００５０】
なお、マッチング処理部２０２は、第１の領域を移動する度に、第１の領域の位置（座標）と、第１の領域と第２の領域との間の相関係数とを対応付けてメモリ２５２に記憶する。このときメモリ２５２に記憶されるデータ（第１の領域の位置と相関係数との組み合わせ）が相関係数マップに相当する。位置ずれ量算出処理部２０３は、第１の領域の位置と相関係数とをメモリ２５２から読み取り、第１の領域の位置と相関係数との関係を示す近似式をプロセッサ２５１により算出し、この近似式を用いて、上記２つの画像間の位置ずれ量を１画素未満の単位で推定する。上記のように、例えば、位置ずれ量算出処理部２０３は、第１の領域にて相関係数の最大値に対応する点のＸ座標とその左右の点のＸ座標、及び、それぞれの点に対応する相関係数をメモリ２５２から読み取り、当該３点のＸ座標及び対応する相関係数から得られる近似２次曲線を算出する。そして、位置ずれ量算出処理部２０３は、この近似２次曲線の頂点（相関係数の最大値の近似値に対応する点）のＸ座標を、横方向における位置ずれ量と推定する。同様に、位置ずれ量算出処理部２０３は、第１の領域にて相関係数の最大値に対応する点のＹ座標とその上下の点のＹ座標、及び、それぞれの点に対応する相関係数をメモリ２５２から読み取り、当該３点のＹ座標及び対応する相関係数から得られる近似２次曲線を算出する。そして、位置ずれ量算出処理部２０３は、この近似２次曲線の頂点（相関係数の最大値の近似値に対応する点）のＹ座標を、縦方向における位置ずれ量と推定する。このような処理を行うことにより、マッチング処理部２０２が第１の領域を画素単位で移動しても、位置ずれ量算出処理部２０３が位置ずれ量を１画素未満の単位で推定できる（例えば、位置ずれ量が縦方向に−０．１画素、横方向に＋０．２画素であるといった結果が得られる）。
【００５１】
次に、レジストレーション処理部１４１のリサンプリング処理部２０４は、リサンプリング処理（ステップＳ４１４）を実施する。この処理は、ペア画像のどちらか一方をリサンプリングして、もう片方に位置合わせする処理である。ここでは、発災後画像に対して処理を実施するものとする。その具体的な手順は次の通りである。まず、リサンプリング処理部２０４は、発災後画像を位置ずれ量算出処理（ステップＳ４１３）にて算出された位置ずれ量だけずらした画像を生成する。ここで、サンプリングされていない画素の値は、周囲の値から計算により内挿する。内挿アルゴリズムとしてはＮｅａｒｅｓｔＮｅｉｇｈｂｏｒやＢｉｌｉｎｅａｒ、Ｂｉｃｕｂｉｃ、Ｌａｎｃｚｏｓなどを利用することができる。なお、位置ずれ量算出処理（ステップＳ４１３）にて位置ずれ量は１画素未満の単位で推定されるものの、推定された位置ずれ量が１画素の整数倍（例えば、＋２．０画素）となることもある。この場合には、位置ずれ量だけずらした画像の各画素は既にサンプリングされているため、画素値の内挿処理は必要ない。一方、推定された位置ずれ量が１画素の小数倍（例えば、＋０．２画素）となった場合、即ち、推定された位置ずれ量の単位が１画素未満であった場合には、位置ずれ量だけずらした画像の各画素は未だサンプリングされていないものとなるため、上記のような画素値の内挿処理が必要となる。
【００５２】
最後に、レジストレーション処理部１４１は、位置合わせした発災前画像と発災後画像とを、同じ大きさに切り出して出力する。
【００５３】
以下では、画像処理装置１０４の変化抽出処理部１４２の動作、即ち、変化抽出処理について、内容を説明する。処理の概要を示すフローチャートを図５に示す。図５に示す通り、変化抽出処理部１４２は、レジストレーション済みペア画像をメモリ２５２から読み取って入力した後、特徴量算出処理（ステップＳ４２１）をプロセッサ２５１により実施する。そして、変化抽出処理部１４２は、得られた特徴量を比較して変化の有無をプロセッサ２５１により判定し、変化抽出画像をメモリ２５２に記憶するとともに、出力装置２５４により出力する。詳細を以下に示す。
【００５４】
まず、変化抽出処理部１４２の特徴量算出処理部２０５は、特徴量算出処理（ステップＳ４２１）を実施する。この処理はペア画像の両方に対して行う。なお、様々な特徴量が考えられるが、ここでは分散値とエッジ量と相関係数を特徴量と呼ぶことにする。その具体的な手順は次の通りである。まず、特徴量算出処理部２０５は、ペア画像をそれぞれＮ×Ｎの小領域に分割する。次に、特徴量算出処理部２０５は、それぞれ分割した小領域内で分散値及びエッジ量を前述した式（１）式及び式（２）より算出する。次に、特徴量算出処理部２０５は、ペア画像間で同じ位置の小領域の相関係数を前述した式（３）より算出する。
【００５５】
次に、変化抽出処理部１４２は、特徴量を比較する。この処理は、ペア画像から得られた特徴量を比較して、変化の有無を判定する処理である。その具体的な手順は次の通りである。まず、変化抽出処理部１４２は、ペア画像の同じ小領域における分散値の差と、エッジ量の差と、相関係数とを、それぞれの閾値と比較する。それぞれの閾値は予めメモリ２５２に記憶されており、変化抽出処理部１４２によって読み出されるものとする。次に、変化抽出処理部１４２は、比較した結果全てが閾値以上の場合にその小領域を変化ありと判断し、そうでない場合は変化なしとする。以降、変化抽出処理部１４２は、全ての小領域について変化の有無を判定する。
【００５６】
最後に、変化抽出処理部１４２は、変化の有無を画像化し、変化抽出画像を作成する。なお、様々な画像化の方法が考えられるが、変化の有無を２値化して表示する方法や、変化の大きさに応じて色分けして表示する方法などがある。
【００５７】
なお、以上で説明した変化抽出処理は、前述したレジストレーション処理が高い精度で実施されていることが前提となっているが、そうでない場合も考慮して次のような工夫をしてもよい。図５のフローチャートにおいて、マッチング処理（ステップＳ４１２）と特徴量算出処理（ステップＳ４２１）とを、特徴量比較の直前に追加する。追加したマッチング処理（ステップＳ４１２）では、変化抽出処理部１４２は、前述した方法で相関係数分布を作成する。追加した特徴量算出処理（ステップＳ４２１）では、変化抽出処理部１４２は、得られた相関係数分布の中から最も相関が高い小領域において、発災前画像の特徴量を算出する。このようにすることで、レジストレーション処理が高い精度で実施されていない場合においても、位置ずれの影響を抑制することができる。
【００５８】
上記のように、変化抽出処理部１４２は、異なる時刻に撮影された２つの画像（例えば、発災前画像及び発災後画像）のうち、一方の画像を位置ずれ量算出処理部２０３が推定した位置ずれ量ずらした画像と他方の画像と（レジストレーション済みペア画像）の特徴量をプロセッサ２５１により比較し、当該比較結果に基づいて、上記２つの画像間の変化を抽出する。
【００５９】
以下では、画像処理装置１０４の雲抽出処理部１４３の動作、即ち、雲抽出処理について、内容を説明する。処理の概要を示すフローチャートを図６に示す。図６に示す通り、雲抽出処理部１４３は、レジストレーション済みペア画像をメモリ２５２から読み取って入力した後、移動ベクトル抽出処理（ステップＳ４３１）、移動量算出処理（ステップＳ４３２）、品質指標算出処理（ステップＳ４３３）をプロセッサ２５１により実施する。そして、雲抽出処理部１４３は、得られた品質指標と移動量とを判断材料として画像中の対象物が雲であるか否かをプロセッサ２５１により判定し、雲抽出画像をメモリ２５２に記憶するとともに、出力装置２５４により出力する。詳細を以下に示す。
【００６０】
まず、雲抽出処理部１４３の移動ベクトル抽出処理部２０６は、移動ベクトル抽出処理（ステップＳ４３１）を実施する。この処理は、ペア画像から雲の移動ベクトルを抽出する処理である。その具体的な手順は次の通りである。まず、移動ベクトル抽出処理部２０６は、発災後画像をＮ×Ｎの小領域に分割する。次に、移動ベクトル抽出処理部２０６は、前述したマッチング処理（ステップＳ４１２）と同様に、それぞれの小領域と発災前画像とのマッチングを前述した式（３）の正規化相関法により行い、相関係数分布を作成する。次に、移動ベクトル抽出処理部２０６は、それぞれの相関係数分布から相関係数が最大となる位置を抽出する。最後に、移動ベクトル抽出処理部２０６は、それぞれの相関係数分布から抽出した位置、即ち、それぞれの小領域についてマッチングを行った際の相関係数が最大値となる位置の発災前画像におけるＸＹ座標（始点）と発災後画像におけるＸＹ座標（終点）とを求め、画像中の物体の移動ベクトルとする。
【００６１】
次に、雲抽出処理部１４３の移動量処理部２０７は、移動量算出処理（ステップＳ４３２）を実施する。この処理は、前述した移動ベクトル抽出処理（ステップＳ４３１）により得られた移動ベクトルの大きさを算出する処理である。その具体的な手順として、移動量処理部２０７は、移動ベクトルの始点から終点までのユークリッド距離を求め、移動ベクトルの大きさとする。
【００６２】
次に、雲抽出処理部１４３の品質指標算出処理部２０８は、品質指標算出処理（ステップＳ４３３）を実施する。この処理は、前述した移動ベクトル抽出処理（ステップＳ４３１）により得られた移動ベクトルの、品質の良し悪しを判断するための指標を算出する処理である。なお、様々な品質指標が考えられるが、ここでは分散値とエッジ量と相関係数を品質指標とすることにする。その具体的な手順は次の通りである。まず、品質指標算出処理部２０８は、発災後画像をＮ×Ｎの小領域に分割する。次に、品質指標算出処理部２０８は、それぞれの小領域内で分散値及びエッジ量を前述した式（１）及び式（２）より算出する。最後に、品質指標算出処理部２０８は、前述したマッチング処理（ステップＳ４１２）と同様に、それぞれの小領域と発災前画像とのマッチングを前述した式（３）の正規化相関法により行い、相関係数分布より最大値を抽出する。なお、品質指標算出処理部２０８は、移動ベクトル抽出処理（ステップＳ４３１）で作成された相関係数分布を利用して、それぞれの相関係数分布より最大値を抽出してもよいし、移動ベクトル抽出処理（ステップＳ４３１）でそれぞれの相関係数分布から抽出された相関係数の最大値をそのまま利用してもよい。
【００６３】
次に、雲抽出処理部１４３は、品質指標と移動量とを判断材料として、雲であるか否かを判定する。その具体的な手順は次の通りである。まず、雲抽出処理部１４３は、１つの小領域における品質指標（分散値、エッジ量、相関係数）を、それぞれの閾値と比較する。それぞれの閾値は予めメモリ２５２に記憶されており、雲抽出処理部１４３によって読み出されるものとする。次に、雲抽出処理部１４３は、比較した結果全てが閾値以上の場合に、移動量と閾値とを比較する。次に、雲抽出処理部１４３は、比較した移動量が閾値以上の場合に、その小領域を雲であると判断する。以降、雲抽出処理部１４３は、全ての小領域について雲の判定をする。このように、品質指標となる属性が一定の閾値以上である小領域のみについて雲の判定をすることで、判定の誤りの発生率を抑えることができる。また、地表面の画像において雲以外に移動するものはないことが想定されるため、移動量が一定の閾値以上である小領域のみを雲と判定することで、さらに、判定の誤りの発生率を抑えることができる。
【００６４】
最後に、雲抽出処理部１４３は、それぞれの小領域が雲であるか否かを画像化し、雲抽出画像を作成する。なお、様々な画像化の方法が考えられるが、雲の有無を２値化して表示する方法や、雲であることの確からしさに応じて色分けして表示する方法などがある。
【００６５】
上記のように、雲抽出処理部１４３は、異なる時刻に撮影された２つの画像（例えば、発災前画像及び発災後画像）のうち、一方の画像を位置ずれ量算出処理部２０３が推定した位置ずれ量ずらした画像と他方の画像と（レジストレーション済みペア画像）で異なる位置に映っている物体（例えば、発災前画像及び発災後画像の間でマッチングをした際に相関係数が最大値となる小領域）をプロセッサ２５１により推定し、上記異なる位置の間の距離（移動量）が所定の閾値以上である場合、その物体が雲であると判断する。
【００６６】
以下では、画像処理装置１０４のテクスチャ分類処理部１４４の動作、即ち、テクスチャ分類処理について、内容を説明する。処理の概要を示すフローチャートを図７に示す。図７に示す通り、テクスチャ分類処理部１４４は、レジストレーション済み画像をメモリ２５２から読み取って入力した後、ＮＤＶＩ算出処理（ステップＳ４４１）、ＨＳＩ変換処理（ステップＳ４４２）をプロセッサ２５１により実施する。そして、テクスチャ分類処理部１４４は、得られたＮＤＶＩ値とＨＳＩ値とを判断材料としてテクスチャ（植物、水域、陸地、山肌、雲）をプロセッサ２５１により分類し、テクスチャ分類画像をメモリ２５２に記憶するとともに、出力装置２５４により出力する。詳細を以下に示す。
【００６７】
まず、テクスチャ分類処理部１４４のＮＤＶＩ算出処理部２０９は、ＮＤＶＩ算出処理（ステップＳ４４１）を実施する。この処理は、植物の緑葉は青領域と赤領域の波長を吸収し近赤外線領域の波長を強く反射するという特徴を活かし、植生指標を算出する処理である。その具体的な手順として、ＮＤＶＩ算出処理部２０９は、下記の式（４）よりＮＤＶＩ値を算出する。
【００６８】
ＮＤＶＩ値は次式で求められる。
【００６９】
【数４】
2009295062.tif 000006
【００７０】
次に、テクスチャ分類処理部１４４のＨＳＩ変換処理部２１０は、ＨＳＩ変換処理（ステップＳ４４２）を実施する。この処理は、ＲＧＢ（Ｒｅｄ・Ｇｒｅｅｎ・Ｂｌｕｅ）の輝度値を持つ画像を顕色系のＨＳＩ空間に変換し、色相（Ｈｕｅ）と彩度（Ｓａｔｕｒａｔｉｏｎ）と明度（Ｉｎｔｅｎｓｉｔｙ）とを算出する処理である。その具体的な手順として、ＨＳＩ変換処理部２１０は、下記の式（５）及び式（６）より明度及び彩度を算出する。なお、色相については後続の処理で使用しないため、ここでは算出しない。
【００７１】
ＨＳＩ変換における明度及び彩度はそれぞれ次式で求められる。
【００７２】
【数５】
2009295062.tif 000007
【００７３】
次に、テクスチャ分類処理部１４４は、画素をテクスチャに分類する。この処理は、得られたＮＤＶＩ値とＨＳＩ値とを判断材料として、植物及び水域及び陸地及び山肌及び雲のいずれかに分類する処理である。その具体的な手順は次の通りである。まず、テクスチャ分類処理部１４４は、ＮＤＶＩ値を閾値と比較し、結果が閾値以上の場合にその画素を植物に分類する。次に、テクスチャ分類処理部１４４は、その画素が植物に分類されなかった場合に、ＮＤＶＩ値を閾値と比較し、結果が閾値未満の場合にその画素を水域に分類する。次に、テクスチャ分類処理部１４４は、その画素が水域に分類されなかった場合に、赤波長の値を閾値と比較し、結果が閾値未満の場合にその画素を陸地に分類する。次に、テクスチャ分類処理部１４４は、その画素が陸地に分類されなかった場合に、Ｓ値とＩ値をそれぞれの閾値と比較し、Ｓ値が閾値以上かつＩ値が閾値未満の場合にその画素を山肌に分類する。最後に、テクスチャ分類処理部１４４は、山肌に分類されなかった画素を雲に分類する。以降、テクスチャ分類処理部１４４は、全ての画素について植物又は水域又は陸地又は山肌又は雲に分類する。それぞれの閾値は予めメモリ２５２に記憶されており、テクスチャ分類処理部１４４によって読み出されるものとする。
【００７４】
最後に、テクスチャ分類処理部１４４は、分類結果を画像化し、テクスチャ分類画像を作成する。なお、様々な画像化の方法が考えられるが、分類に応じて色分けして表示する方法や、分類に応じて異なる模様を表示する方法などがある。
【００７５】
上記のように、テクスチャ分類処理部１４４は、異なる時刻に撮影された２つの画像（例えば、発災前画像及び発災後画像）のうち、一方の画像を位置ずれ量算出処理部２０３が推定した位置ずれ量ずらした画像と他方の画像と（レジストレーション済みペア画像）の少なくともいずれかの各画素の赤波長値とＮＤＶＩ値とＳ値とＩ値とをプロセッサ２５１により算出し、各画素の赤波長値とＮＤＶＩ値とＳ値とＩ値とに基づいて、各画素が植物と水域と陸地と山肌と雲とのいずれであるか判断する。
【００７６】
以下、本実施の形態が奏する効果を説明する。
【００７７】
本実施の形態によれば、変化抽出システム１００の画像処理装置１０４が、ペア画像のレジストレーション処理を実施するレジストレーション処理部１４１を具備しているため、撮像された位置にずれがあるペア画像間でも変化抽出が可能である。したがって、事前にペア画像のレジストレーションを行っておくといった、人手による作業が不要となる効果が得られる。
【００７８】
さらに、本実施の形態によれば、画像処理装置１０４が具備する変化抽出処理部１４２が、変化抽出のための特徴量として、分散値とエッジ量と相関係数を用いているため、ペア画像の間で季節や時刻といった時期的な撮像条件が異なる場合においても、それに大きく影響されることなく変化抽出が可能である。したがって、事前にどちらか一方の画像をもう一方の画像の撮像条件に合わせるような補正をするといった、人手による作業が不要となる効果が得られる。
【００７９】
さらに、本実施の形態によれば、画像処理装置１０４が、画像から雲の抽出処理を実施する雲抽出処理部１４３を具備しているため、雲を含む画像でもそれに大きく影響されることなく変化抽出が可能である。したがって、変化抽出後にその変化の原因が雲の移動によるものであるか否かを人間が判読するといった、人手による作業が不要となる効果が得られる。
【００８０】
さらに、本実施の形態によれば、画像処理装置１０４が、テクスチャの分類処理を実施するテクスチャ分類処理部１４４を具備しているため、自動的にテクスチャの分類をすることが可能である。したがって、変化抽出後にその変化が発生した場所がどのようなところなのかを人間が判読するといった、人手による作業が不要となる効果が得られる。
【図面の簡単な説明】
【００８１】
【図１】実施の形態１に係る変化抽出システムの構成の一例を示すブロック図である。
【図２】実施の形態１に係る画像処理装置の構成の一例を示すブロック図である。
【図３】実施の形態１に係る画像処理装置のハードウェア構成の一例を示す図である。
【図４】実施の形態１に係るレジストレーション処理の概要を示すフローチャートである。
【図５】実施の形態１に係る変化抽出処理の概要を示すフローチャートである。
【図６】実施の形態１に係る雲抽出処理の概要を示すフローチャートである。
【図７】実施の形態１に係るテクスチャ分類処理の概要を示すフローチャートである。
【符号の説明】
【００８２】
１００  変化抽出システム、１０１  静止光学衛星、１０２  衛星画像受信センタ、１０３  衛星画像データベース、１０４  画像処理装置、１０５  端末装置、１０６  ネットワーク、１４１  レジストレーション処理部、１４２  変化抽出処理部、１４３  雲抽出処理部、１４４  テクスチャ分類処理部、２０１  マッチング小領域選定処理部、２０２  マッチング処理部、２０３  位置ずれ量算出処理部、２０４  リサンプリング処理部、２０５  特徴量算出処理部、２０６  移動ベクトル抽出処理部、２０７  移動量処理部、２０８  品質指標算出処理部、２０９  ＮＤＶＩ算出処理部、２１０  ＨＳＩ変換処理部、２５１  プロセッサ、２５２  メモリ、２５３  入力装置、２５４  出力装置、９０１  表示装置、９０２  キーボード、９０３  マウス、９０４  ＦＤＤ、９０５  ＣＤＤ、９０６  プリンタ装置、９１１  ＣＰＵ、９１２  バス、９１３  ＲＯＭ、９１４  ＲＡＭ、９１５  通信ボード、９２０  磁気ディスク装置、９２１  オペレーティングシステム、９２２  ウィンドウシステム、９２３  プログラム群、９２４  ファイル群。
【図１】
2009295062.tif 000008
【図２】
2009295062.tif 000009
【図３】
2009295062.tif 000010
【図４】
2009295062.tif 000011
【図５】
2009295062.tif 000012
【図６】
2009295062.tif 000013
【図７】
2009295062.tif 000014