JP3816825B2 - 市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体 - Google Patents
市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体 Download PDFInfo
- Publication number
- JP3816825B2 JP3816825B2 JP2002091060A JP2002091060A JP3816825B2 JP 3816825 B2 JP3816825 B2 JP 3816825B2 JP 2002091060 A JP2002091060 A JP 2002091060A JP 2002091060 A JP2002091060 A JP 2002091060A JP 3816825 B2 JP3816825 B2 JP 3816825B2
- Authority
- JP
- Japan
- Prior art keywords
- urban area
- image
- edge
- area extraction
- extraction processing
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Lifetime
Links
Images
Landscapes
- Image Processing (AREA)
- Image Analysis (AREA)
Description
【０００１】
【発明の属する技術分野】
本発明は，衛星写真画像や航空写真画像などの上空から撮影された画像の中から，そのテクスチャ情報およびマルチスペクトル情報にもとづき市街地領域を自動的に抽出する市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体に関するものである。
【０００２】
【従来の技術】
従来の市街地領域抽出技術では，特定の色や輝度条件を満たす閉領域に着目し，その閉領域の形状があらかじめ規定しておいたものに近い場合に，該閉領域を建造物領域であるとし，このような手法で抽出された建造物を統合して市街地であると認識していた。
【０００３】
【発明が解決しようとする課題】
しかし，この種の従来の方法は，画像がモノクロまたはグレースケールで色情報がない場合や，建造物の一部が木や別の建造物で隠れていて閉領域形状を取らない場合や，画像の解像度が低いために各建造物が十分な大きさを持たず閉領域の形状が崩れているような場合などに，安定して市街地領域抽出を行うことができないという問題があった。
【０００４】
本発明は，市街地領域でのテクスチャ特徴，すなわち市街地領域ではエッジのパワーが強く，かつエッジ接線方向が直交している傾向が強いという特徴に着目することにより，上記従来の技術では問題となるような場合でも，安定して市街地領域抽出を行うことができるようにすることを目的とする。
【０００５】
なお，ここでエッジパワーとは，画像各点の２次元エッジベクトルの長さをいい，またエッジ接線とは，画像各点における該点の２次元エッジベクトルを法線ベクトルとする直線をいう。
【０００６】
【課題を解決するための手段】
本発明は，前記課題を解決するため，市街地領域のエッジパワーとエッジ接線の直交性を定量的に表す特徴量にもとづき市街地領域判定を行うものであり，以下の方法を用いることを特徴とする。
１．可視光の赤色領域の波長をもつ電磁波をとらえた画像（Ｒ画像という）と近赤外領域の電磁波をとらえた画像（ＮＩＲ画像という）とを用い，水領域と森林領域を抽出する。
２．市街地領域抽出対象の画像に対し，上記１の結果を使用した適当な前処理を施す。
３．市街地領域抽出対象の画像の各画素について２次元のエッジベクトルを取得し，各画素を中心とするブロック内において，すべての２ベクトルの組み合わせについてベクトル外積の絶対値を計算し加算する。ベクトルＡ，Ｂの外積の絶対値は，この２ベクトルのなす角度をθとした場合，
｜Ａ×Ｂ｜＝｜（｜Ａ｜｜Ｂ｜ｓｉｎθ）｜
と計算されるため，各エッジベクトルのパワーが大きく，かつ２つのエッジベクトルの直交性が強いほど大きな値となる。この結果得られるベクトル外積和を，ブロック内の代表点に特徴量として付与し，この値にもとづき該点が市街地であるか否かを判定する。
４．市街地領域抽出対象の画像について各画素を中心とする任意サイズのブロックを設定し，各画素において横方向および縦方向のエッジ差分計算を行い，得られたエッジ差分値にもとづきエッジ接線の角度パラメータと距離パラメータとを求める。算出されたエッジ接線パラメータに対してハフ（Ｈｏｕｇｈ）変換を行い，角度パラメータと距離パラメータとを２つの軸とする２次元の投票空間に各エッジ接線パラメータを投票する。得られた投票空間から適当な特徴量を求め，これにもとづきブロックの代表点が市街地であるか否かを判定する。
５．上記３で抽出された市街地領域と上記４で抽出された市街地領域とを足し合わせ，そこに上記１で抽出される森林領域を除外して最終的な市街地領域を求める。
【０００７】
上記３の手法は，小さな家などが集合している市街地領域の抽出に有効であり，ブロックサイズは家程度の大きさにすると良い結果を示す。また，上記４の手法は，工場など大きな建造物の抽出に有効であり，ブロックサイズは上記３で用いるものより大きなものを使用すると良い結果を示す。
【０００８】
以上の方法により，市街地のテクスチャ特徴にもとづく安定した市街地領域抽出が可能となる。
【０００９】
【発明の実施の形態】
本発明の実施の形態を，図を用いて詳細に説明する。
【００１０】
図１に，本実施の形態による市街地領域抽出処理装置の構成例を示す。図中，１はＣＰＵおよびメモリなどからなる市街地領域抽出処理装置，１０は衛星写真画像や航空写真画像等，地表を可視光の赤色領域の波長をもつ電磁波でとらえたＲ画像，１１はＲ画像１０と同一地点を近赤外領域の電磁波でとらえたＮＩＲ画像，１２はＲ画像１０およびＮＩＲ画像１１と同一地点を可視光で撮影した画像（以下，Ｐ画像という），２０はＲ画像１０とＮＩＲ画像１１とから水領域を抽出する水領域抽出処理部，３０はＲ画像１０とＮＩＲ画像１１とから森林領域を抽出する森林領域抽出処理部，１３は抽出された森林領域，１４は抽出された水領域，４０はＰ画像１２に対して行うヒストグラム平坦化処理部，５０はベクトル外積和を特徴量とする市街地領域抽出処理部，１５は市街地領域抽出処理部５０において内積の符号チェックを行わない場合に抽出された市街地領域Ａ，１６は内積の符号チェックを行った場合に抽出された市街地領域Ｂ，６０はＨｏｕｇｈ特徴による市街地領域抽出処理部，１７は市街地領域抽出処理部６０の処理により抽出された市街地領域Ｃ，７０は市街地領域Ａ１５と市街地領域Ｃ１７とを統合するマージ処理部，１８はマージ処理部７０の結果得られるマージ結果，８０はマージ結果１８に対して森林領域の修正を行う森林領域修正処理部，１９は最終的に得られる市街地領域である。
【００１１】
Ｒ画像１０は，人工衛星や航空機を用いて，上空から地表を可視光の赤色領域の波長をもつ電磁波で撮影した画像である。
【００１２】
ＮＩＲ画像１１は，人工衛星や航空機を用いて，上空から地表を近赤外領域の波長をもつ電磁波で撮影した画像である。一般的に地表撮影用の人工衛星や航空機には，Ｒ画像やＮＩＲ画像を撮影するための機器が搭載されている。
【００１３】
Ｐ画像１２は，人工衛星や航空機を用いて，上空から地表を可視光で撮影した画像である。一般的にＰ画像は，Ｒ画像やＮＩＲ画像に比べて解像度が高いため，本実施の形態では，Ｐ画像から市街地領域を抽出する方法を説明する。Ｐ画像の代わりにＲ画像やＮＩＲ画像，または別の波長領域で撮影した画像を市街地領域抽出対象の画像とすることもできる。
【００１４】
水領域抽出処理部２０は，同一地点を撮影したＲ画像１０とＮＩＲ画像１１とから，水領域を抽出する処理手段である。はじめにアフィン変換などで２枚の画像の位置合わせが必要となる場合も考えられる。位置合わせの後，対応する各点において，
Ｘ＝（ＮＩＲ画像における輝度値）／（Ｒ画像における輝度値）
を計算する。水の反射率は可視光の赤色領域と近赤外領域とで大きく異なり，近赤外領域での反射率は可視光赤色領域の反射率よりたいへん小さい。したがって，水領域では上記Ｘは小さな値をとる。そこでＸの値が一定値以下である場合，該点を水領域であると判定する。この一定値は，照明条件やカメラ特性によって大きく値が変動するが，例えば，
ＮＩＲ／Ｒ＜０．７
の場合に，その画素を水領域と判定する。
【００１５】
また上記Ｘ以外にも，ＮＩＲ画像１１での輝度値とＲ画像１０での輝度値との差や，下記で説明するＮＤＶＩを水領域の判定尺度とする方法も考えられる。
【００１６】
森林領域抽出処理部３０は，水領域抽出処理部２０と同様に，同一地点を撮影したＲ画像１０とＮＩＲ画像１１とから森林領域を抽出する処理手段である。ここでもアフィン変換などで２枚の画像の位置合わせが必要となる場合も考えられる。その後，対応する各点において，
ＮＤＶＩ＝（ＮＩＲ画像における輝度値−Ｒ画像における輝度値）
／（ＮＩＲ画像における輝度値＋Ｒ画像における輝度値）
を計算する。森林領域では，クロロフィルの分光反射率特性より上記ＮＤＶＩの値が大きくなる（参照：“図解リモートセンシング”，日本測量協会）。ＮＤＶＩの値が一定値以上である場合，該点を森林領域と判定する。
【００１７】
森林領域１３は，森林領域抽出処理部３０により抽出された領域である。森林領域１３に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００１８】
水領域１４は，水領域抽出処理部２０により抽出された領域である。水領域１４に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００１９】
ヒストグラム平坦化処理部４０は，必要に応じてＰ画像に対して輝度変換を行う処理手段であり，本実施の形態では，ヒストグラム平坦化（参照：“コンピュータ画像処理入門”，総研出版，田村秀行監修）を行う。この処理により建物の屋根部分と影部分とのコントラストが大きくなり，市街地領域のテクスチャ特徴がより強調される。
【００２０】
また本実施の形態では，ヒストグラム平坦化の際に水領域１４を用いる方法について説明する。海や湖など，輝度値の小さい領域が画像中に大きく存在する場合，ヒストグラム平坦化処理では適切なコントラスト強調効果が現れない場合がある。そこで，ヒストグラム平坦化処理に必要な輝度ヒストグラムを，水領域を除外した部分で取得し，またヒストグラム平坦化処理自体も水領域を除外した部分にのみ適用する。このような方法を用いることにより，画像中に海や湖などが大きく存在している場合にも，市街地領域のコントラストが強調される。
【００２１】
ベクトル外積和による市街地領域抽出処理部５０は，ベクトル外積和を特徴量として市街地領域を抽出する処理手段である。図２に，ベクトル外積和による市街地領域抽出処理部５０による処理の流れを示す。
【００２２】
入力画像（Ｐ画像）５１は，ヒストグラム平坦化処理部４０で出力される画像である。エッジベクトル取得処理５２は，画像に差分計算を施してエッジベクトルを取得する処理である。差分計算を行う前に，前処理としてガウシアンフィルタなどでぼかし処理を施すことにより，ノイズの影響を小さくすることができる。
【００２３】
本実施の形態では，横方向と縦方向の差分計算を行う方法について説明する。図３は，そのエッジベクトルの算出方法を説明する図である。横方向の差分計算は，着目画素の左右の画素値を減算することより，縦方向の差分計算は上下の画素値を減算することにより行う。
【００２４】
図３（ａ）に示した例では，画素Ｃ（点Ｃ）における横方向の差分値をＣX ，縦方向の差分値をＣY として，エッジベクトルＣ（ＣX ，ＣY ）は，隣接画素Ａ，Ｂ，Ｄ，Ｅの画素値から，
ＣX ＝（Ｄ−Ｂ）／２
ＣY ＝（Ｅ−Ａ）／２
と計算される。
【００２５】
また，画像の端にある画素については，上下左右の隣接画素が存在しない場合もある。この場合には，隣接画素と着目画素自身との差分を計算して，エッジベクトルの成分とする。このような例を図３（ｂ）に示す。画素Ｇ（点Ｇ）におけるエッジベクトルＧ（ＧX ，ＧY ）は，隣接画素Ｆ，Ｈ，Ｉの画素値から，
ＧX ＝Ｈ−Ｇ
ＧY ＝（Ｉ−Ｆ）／２
と計算される。
【００２６】
エッジベクトル情報５３は，エッジベクトル取得処理５２の結果得られ，全画素に２次元のエッジベクトルが付与されている。図４に，エッジベクトル情報５３の例を示す。図４（ａ）は原画像の画素値，（ｂ）はエッジベクトル成分の値を示している。
【００２７】
ブロック取得処理５４は，着目画素を中心として一定サイズの領域を取得する処理である。図５は，領域のサイズを３画素×３画素とした場合のブロックの例を示している。
【００２８】
外積和計算処理５５は，ブロック取得処理５４で取得したブロック内において，すべての２つのエッジベクトルの組み合わせ間でのベクトル外積を計算し，その絶対値を加算していく処理である。ここでは，図６に示した例のように，９画素で構成されるブロック内のすべての２つのベクトルの組み合わせ，すなわち 9Ｃ2 ＝３６通りの組み合わせ間でベクトル外積の絶対値を計算して加算する。
【００２９】
また，図７（ａ）に示すように，白地に黒い領域と黒字に白い領域とを区別するため，外積を計算する２つのベクトルのうちのどちらか一方のエッジベクトルＶE と，該エッジベクトルＶE を付与されたエッジ点から他方のエッジ点に向かう方向ベクトルＶD との内積を計算し，図７（ｂ）に示すように内積値がプラス値の場合にのみベクトル外積の加算処理を行い，図７（ｃ）に示すように内積値がマイナス値の場合には，それら２つのエッジベクトルの組み合わせにおけるベクトル外積の加算処理をスキップする方法も考えられる（図７（ｄ））。
【００３０】
外積和５６は，外積和計算処理５５の結果得られ，画像各点に特徴量として与えられる。
【００３１】
市街地判定処理５７は，外積和計算処理５５で得られた外積和５６にもとづき，ブロックの代表点を市街地であるかどうかを判定する処理である。代表点はブロックの中心点，ブロックを構成する全点，ブロック中心付近の点など，任意に設定可能であるが，ここではブロックの中心点とする。外積和の値に対してしきい値処理を行って市街地と非市街地領域の二つに分別する。また，しきい値処理を行わずに，外積和の値をそのまま市街地らしさを表す値として出力することも考えられる。
【００３２】
また，海岸線などの誤抽出を抑えるため，水領域抽出処理部２０で抽出した水領域１４を用い水領域の近郊では市街地判定を取り消す方法も考えられる。
【００３３】
終了判定処理５８は，全画素についてブロック取得処理５４，外積和計算処理５５，市街地判定処理５７が行われたかどうかを判定する処理である。まだ未処理である画素が存在する場合には，ブロック取得処理５４に戻り，未処理画素について同様に市街地判定を行う。
【００３４】
市街地領域Ａ（内積チェック無）１５は，ベクトル外積和による市街地領域抽出処理部５０において内積の符号チェックを行わない場合に抽出された市街地領域である。市街地領域Ａ１５に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００３５】
また，市街地領域Ｂ（内積チェック有）１６は，ベクトル外積和による市街地領域抽出処理部５０において内積の符号チェックを行った場合に抽出された市街地領域である。市街地領域Ｂ１６に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００３６】
Ｈｏｕｇｈ特徴による市街地領域抽出処理部６０は，エッジ接線パラメータのＨｏｕｇｈ変換によって算出されるＨｏｕｇｈ特徴量にもとづいて市街地領域を抽出する処理手段である。図８に，Ｈｏｕｇｈ特徴による市街地領域抽出処理部６０の処理の流れを示す。
【００３７】
入力画像（Ｐ画像）６１は，ヒストグラム平坦化処理部４０で出力される画像である。ブロック取得処理６２は，着目画素を中心として一定サイズのブロックを取得する処理である。ブロックの形は正方形，円形など任意であるが，本実施の形態では，円形のブロックを取得するものとする。これはその後に行われるエッジ接線パラメータのＨｏｕｇｈ変換に都合が良いからである。円形ブロックの例を，図９に示す。
【００３８】
エッジ接線パラメータ計算処理６３は，ブロック内の各点におけるエッジ接線を求める処理である。ここでは，まず画像に横方向および縦方向の差分計算を施してエッジベクトルを取得し，取得したエッジベクトルにもとづきエッジ接線パラメータを求める手法を説明する。エッジベクトルを求める処理は，図２で説明したエッジベクトル取得処理５２と同様である。
【００３９】
次に，得られたエッジベクトルからエッジ接線パラメータを計算する。ここでブロックの中心を原点とおき，原点から接線までの距離ｒと，接線の法線ベクトルがＸ方向となす角度θを求めるべきエッジ接線パラメータとする。
【００４０】
図１０は，エッジ接線パラメータの例を説明する図である。図１０から明らかなように，点Ｐにおけるエッジベクトル成分ＰX とＰY とから，点Ｐにおけるエッジ接線のパラメータθP とｒP は簡単に計算される。
【００４１】
Ｈｏｕｇｈ変換処理６４は，算出されたエッジ接線パラメータに対してＨｏｕｇｈ変換を行い，２次元の投票空間に投票していく処理である。図１１は，その投票処理の例を説明する図である。図１０で示した点Ｐでのエッジ接線パラメータθP とｒP は，投票空間内の座標（θP ，ｒP ）に投票される。
【００４２】
投票される値は，＋１の場合もあるが，本実施の形態ではエッジパワー値そのもの，すなわち（ＰX 2 ＋ＰY 2 ）の平方根を投票値とする。この方法により，市街地領域に多く存在するパワーが強いエッジのパラメータが，より強調されることになる。
【００４３】
また，ノイズの悪影響の抑制等を目的として，エッジパワーが一定値以下の場合には投票処理を行わないとか，またはエッジパワーに何らかの重みをつけて投票値とするなどの方法も有効である。
【００４４】
投票空間６５は，Ｈｏｕｇｈ変換処理６４をブロック内の全画素について行った結果，得られる投票空間である。図１２にその例を示す。図１２（ａ）は原画像，図１２（ｂ）は投票空間の例を示している。図１２（ａ）の例では直径５０画素の円形ブロック６９内の全画素に対して投票を行っている。図１２（ｂ）の投票空間では，得票値の大小を濃度を変えて示している。投票空間の横軸は，角度パラメータθの軸（範囲０〜３６０°），縦軸は距離パラメータｒの軸（範囲０〜２５）である。画像内の直線が投票空間内では得票値の高い点として変換されていることが分かる。
【００４５】
ブロック内処理終了判定６６は，ブロック内の全画素についてエッジ接線パラメータ計算処理６３，Ｈｏｕｇｈ変換処理６４が終了したか否かを判定する処理である。未処理画素が存在する場合，再びエッジ接線パラメータ計算処理６３へ戻り，未処理画素の処理を行う。
【００４６】
市街地領域Ｂ（内積チェック有）１６は，図１のベクトル外積和による市街地領域抽出処理部５０で抽出された市街地領域Ｂである。水領域１４は，図１の水領域抽出処理部２０で抽出された水領域である。
【００４７】
市街地判定処理６７は，投票空間６５にもとづき，ブロックの代表点が市街地であるかどうかを判定する処理である。ブロックの代表点として，ブロックの中心画素，ブロック中心付近の画素，ブロックを構成する全画素，ブロック内の一定値以上のエッジパワーを持つ画素などが考えられる。本実施の形態では，２次元の投票空間をそのまま用いず，Ｘ軸方向（角度パラメータ軸方向）に射影をとり，１次元に縮退させたヒストグラム（以下，角度ヒストグラムと呼称する）を用いる方法を説明する。
【００４８】
図１３に，図１２（ｂ）に示す投票空間から作成した角度ヒストグラムを示す。角度ヒストグラムの最大値，最小値，平均値，分散値などの統計量から市街地を求める方法もあるが，ここでは（１）最大値／最小値，（２）エッジ接線直交性にもとづく判定処理について説明する。
【００４９】
（１）最大値／最小値
角度ヒストグラムの最大値と最小値を求め，最大値を最小値で除算して比率を計算する。この値が一定値以上である場合，ブロックの代表点を市街地と判定する。接線パラメータ計算時の原点をブロックの中心に設定しているため，見た目の傾きが同じ直線でもブロック中心に対して反対側に位置するとき，その角度パラメータθは１８０°の差を持つこととなる。そこで角度ヒストグラムの１８０〜３６０°部分を０〜１８０°部分に加算し，角度ヒストグラムの範囲を０〜１８０°とすることで，同一の傾きを持つ直線は同一の角度パラメータを持つように補正する方法も考えられる。
【００５０】
ランダムノイズのようにエッジ方向がばらばらである画像の場合，最大値／最小値は小さな値であり，この値が大きな値となるのはブロック内に強いエッジパワーを持つ直線が存在するなどの場合である。
【００５１】
（２）エッジ接線直交性
ブロック内に直交する直線が存在する場合，角度ヒストグラムには９０°の間隔を持つ２つのピークが発生する。９０°間隔の２つのピークがあるか否かは，次の方法により判定することができる。図１４に，その判定方法の説明図を示す。
【００５２】
図１４に示すヒストグラムは，上記（１）において説明したように１８０〜３６０°部分を０〜１８０°部分に重ねて角度ヒストグラムの角度範囲を０〜１８０゜としたものである。はじめにこのヒストグラムにおいて最大値を取る角度θ1 を求める。次に，θ1 から９０°の差を持つ角度θ2 を中心として任意の幅を持つ第２ピーク探索領域ＳA を設定し，その領域内での最大値（第２ピーク値）を求める。図１４の例では，θ2 を中心として±１０°の幅を第２ピーク探索領域ＳA の範囲と設定した。
【００５３】
第２ピーク値を求めた後，この値をヒストグラムの最小値で除算し，その比率が一定値以上の場合，該ブロックの代表点を市街地領域であると判定する。第２ピーク値と最小値との比率以外にも，第２ピーク値そのものの値，第２ピーク値と最小値との差分などを判定尺度としたり，最大値そのものの値も判定尺度に加えることも考えられる。
【００５４】
また，範囲を０〜１８０°とした角度ヒストグラムの，さらに９０〜１８０°部分を０〜９０°部分に重ね，角度範囲を０〜９０°としたヒストグラムを作成し，該ヒストグラムの最大値と最小値との比率を判定尺度とすることも考えられる。
【００５５】
ここでは，角度ヒストグラムにもとづき市街地判定する方法を説明したが，２次元の投票空間をそのまま使用して，適当な特徴量を算出する方法も考えられる。
【００５６】
また，水田や畑などが建造物として誤抽出されるのを抑制するため，市街地領域Ｂ（内積チェック有）１６を用い，上記判定と併せた再判定を行う方法もある。建造物は黒地に白い領域として存在する傾向が強いのに対し，水田や畑などは白地に黒い領域であるため，市街地領域Ｂ（内積チェック有）１６では抽出されない。そこで上記方法により市街地と判定されたブロック内に，市街地領域Ｂ（内積チェック有）１６が存在しない場合に，該ブロックは非市街地領域として前の判定結果をリジェクトする。
【００５７】
また，橋や埋め立て地などの誤抽出を抑えるため，水領域１４を用い水領域の近郊における市街地判定をリジェクトする方法も考えられる。
【００５８】
終了判定処理６８は，全画素について市街地判定処理が行われたか否かを判定する処理である。まだ未処理である画素が存在する場合，ブロック取得処理６２に再度戻り，未処理画素について同様に処理を行う。
【００５９】
市街地領域Ｃ１７は，Ｈｏｕｇｈ特徴による市街地領域抽出処理６０の結果抽出された市街地領域である。市街地領域Ｃ１７に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００６０】
マージ処理部７０は，ベクトル外積和を特徴量として抽出した市街地領域Ａ１５と，Ｈｏｕｇｈ特徴にもとづき抽出した市街地領域Ｃ１７とをマージする処理手段である。本実施の形態では，市街地領域Ａ１５と市街地領域Ｃ１７とを２値化して加算する単純な方法で行うこととする。本方法以外にも，市街地領域Ａ１５と市街地領域Ｃ１７の市街地らしさを多値で表し，その値に何らかの重みを付けて加算するなどの方法も考えられる。
【００６１】
マージ結果１８は，マージ処理部７０の結果，取得される結果である。マージ結果１８に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００６２】
森林領域修正処理部８０は，マージ結果１８に対して森林領域１３を除外し，最終的な市街地領域１９を抽出する処理手段である。除外する方法として，本実施の形態では，森林領域１３およびマージ結果１８を２値化して減算するといった単純な方法で行うこととする。
【００６３】
市街地領域１９は，森林領域修正処理部８０の結果取得される最終的な市街地領域の抽出結果である。この市街地領域１９に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００６４】
以上の市街地領域抽出処理装置１が行う処理は，コンピュータとソフトウェアプログラムとによって実現することができ，そのプログラムは，コンピュータが読み取り可能な可搬媒体メモリ，半導体メモリ，ハードディスク等の適当な記録媒体に格納して，そこから読み出すことによりコンピュータに実行させることができる。また，そのプログラムは通信回線を経由して他のコンピュータからダウンロードすることができ，それをインストールして実行させることもできる。
【００６５】
【発明の効果】
以上説明したように，本発明は，市街地領域でのテクスチャ特徴とマルチスペクトル特徴，特に市街地領域ではエッジのパワーが強く，かつエッジ接線方向が直交している傾向が強いというテクスチャ特徴に着目した市街地領域抽出を行う。家などの小さな建造物が集合している市街地領域の抽出に有効なベクトル外積和を特徴量とする市街地領域抽出と，工場などの大規模建造物の抽出に有効なＨｏｕｇｈ特徴にもとづく市街地領域抽出とを組み合わせることにより，建造物の大きさに依存しない安定した市街地領域抽出を行うことが可能である。また，マルチスペクトル情報を用いて森林領域と水領域を抽出することにより，より精度の高い市街地領域抽出が可能となった。
【図面の簡単な説明】
【図１】本実施の形態による市街地領域抽出処理装置の構成例を示す図である。
【図２】ベクトル外積和による市街地領域抽出処理部による処理の流れを示す図である。
【図３】エッジベクトルの算出方法を説明する図である。
【図４】エッジベクトル情報の例を示す図である。
【図５】ブロックの例を示す図である。
【図６】ベクトル外積和の計算の例を示す図である。
【図７】内積符号チェック方法の例を説明する図である。
【図８】Ｈｏｕｇｈ特徴による市街地領域抽出処理の流れを示す図である。
【図９】円形ブロックの例を示す図である。
【図１０】エッジ接線パラメータの例を説明する図である。
【図１１】投票処理の例を説明する図である。
【図１２】投票空間の例を示す図である。
【図１３】角度ヒストグラムの例を示す図である。
【図１４】第２ピーク探索処理を説明する図である。
【符号の説明】
１ 市街地領域抽出処理装置
１０ Ｒ画像
１１ ＮＩＲ画像
１２ Ｐ画像
１３ 森林領域
１４ 水領域
１５ 市街地領域Ａ（内積チェック無）
１６ 市街地領域Ｂ（内積チェック有）
１７ 市街地領域Ｃ
１８ マージ結果
１９ 市街地領域
２０ 水領域抽出処理部
３０ 森林領域抽出処理部
４０ ヒストグラム平坦化処理部
５０ ベクトル外積和による市街地領域抽出処理部
６０ Ｈｏｕｇｈ特徴による市街地領域抽出処理部
７０ マージ処理部
８０ 森林領域修正処理部
【発明の属する技術分野】
本発明は，衛星写真画像や航空写真画像などの上空から撮影された画像の中から，そのテクスチャ情報およびマルチスペクトル情報にもとづき市街地領域を自動的に抽出する市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体に関するものである。
【０００２】
【従来の技術】
従来の市街地領域抽出技術では，特定の色や輝度条件を満たす閉領域に着目し，その閉領域の形状があらかじめ規定しておいたものに近い場合に，該閉領域を建造物領域であるとし，このような手法で抽出された建造物を統合して市街地であると認識していた。
【０００３】
【発明が解決しようとする課題】
しかし，この種の従来の方法は，画像がモノクロまたはグレースケールで色情報がない場合や，建造物の一部が木や別の建造物で隠れていて閉領域形状を取らない場合や，画像の解像度が低いために各建造物が十分な大きさを持たず閉領域の形状が崩れているような場合などに，安定して市街地領域抽出を行うことができないという問題があった。
【０００４】
本発明は，市街地領域でのテクスチャ特徴，すなわち市街地領域ではエッジのパワーが強く，かつエッジ接線方向が直交している傾向が強いという特徴に着目することにより，上記従来の技術では問題となるような場合でも，安定して市街地領域抽出を行うことができるようにすることを目的とする。
【０００５】
なお，ここでエッジパワーとは，画像各点の２次元エッジベクトルの長さをいい，またエッジ接線とは，画像各点における該点の２次元エッジベクトルを法線ベクトルとする直線をいう。
【０００６】
【課題を解決するための手段】
本発明は，前記課題を解決するため，市街地領域のエッジパワーとエッジ接線の直交性を定量的に表す特徴量にもとづき市街地領域判定を行うものであり，以下の方法を用いることを特徴とする。
１．可視光の赤色領域の波長をもつ電磁波をとらえた画像（Ｒ画像という）と近赤外領域の電磁波をとらえた画像（ＮＩＲ画像という）とを用い，水領域と森林領域を抽出する。
２．市街地領域抽出対象の画像に対し，上記１の結果を使用した適当な前処理を施す。
３．市街地領域抽出対象の画像の各画素について２次元のエッジベクトルを取得し，各画素を中心とするブロック内において，すべての２ベクトルの組み合わせについてベクトル外積の絶対値を計算し加算する。ベクトルＡ，Ｂの外積の絶対値は，この２ベクトルのなす角度をθとした場合，
｜Ａ×Ｂ｜＝｜（｜Ａ｜｜Ｂ｜ｓｉｎθ）｜
と計算されるため，各エッジベクトルのパワーが大きく，かつ２つのエッジベクトルの直交性が強いほど大きな値となる。この結果得られるベクトル外積和を，ブロック内の代表点に特徴量として付与し，この値にもとづき該点が市街地であるか否かを判定する。
４．市街地領域抽出対象の画像について各画素を中心とする任意サイズのブロックを設定し，各画素において横方向および縦方向のエッジ差分計算を行い，得られたエッジ差分値にもとづきエッジ接線の角度パラメータと距離パラメータとを求める。算出されたエッジ接線パラメータに対してハフ（Ｈｏｕｇｈ）変換を行い，角度パラメータと距離パラメータとを２つの軸とする２次元の投票空間に各エッジ接線パラメータを投票する。得られた投票空間から適当な特徴量を求め，これにもとづきブロックの代表点が市街地であるか否かを判定する。
５．上記３で抽出された市街地領域と上記４で抽出された市街地領域とを足し合わせ，そこに上記１で抽出される森林領域を除外して最終的な市街地領域を求める。
【０００７】
上記３の手法は，小さな家などが集合している市街地領域の抽出に有効であり，ブロックサイズは家程度の大きさにすると良い結果を示す。また，上記４の手法は，工場など大きな建造物の抽出に有効であり，ブロックサイズは上記３で用いるものより大きなものを使用すると良い結果を示す。
【０００８】
以上の方法により，市街地のテクスチャ特徴にもとづく安定した市街地領域抽出が可能となる。
【０００９】
【発明の実施の形態】
本発明の実施の形態を，図を用いて詳細に説明する。
【００１０】
図１に，本実施の形態による市街地領域抽出処理装置の構成例を示す。図中，１はＣＰＵおよびメモリなどからなる市街地領域抽出処理装置，１０は衛星写真画像や航空写真画像等，地表を可視光の赤色領域の波長をもつ電磁波でとらえたＲ画像，１１はＲ画像１０と同一地点を近赤外領域の電磁波でとらえたＮＩＲ画像，１２はＲ画像１０およびＮＩＲ画像１１と同一地点を可視光で撮影した画像（以下，Ｐ画像という），２０はＲ画像１０とＮＩＲ画像１１とから水領域を抽出する水領域抽出処理部，３０はＲ画像１０とＮＩＲ画像１１とから森林領域を抽出する森林領域抽出処理部，１３は抽出された森林領域，１４は抽出された水領域，４０はＰ画像１２に対して行うヒストグラム平坦化処理部，５０はベクトル外積和を特徴量とする市街地領域抽出処理部，１５は市街地領域抽出処理部５０において内積の符号チェックを行わない場合に抽出された市街地領域Ａ，１６は内積の符号チェックを行った場合に抽出された市街地領域Ｂ，６０はＨｏｕｇｈ特徴による市街地領域抽出処理部，１７は市街地領域抽出処理部６０の処理により抽出された市街地領域Ｃ，７０は市街地領域Ａ１５と市街地領域Ｃ１７とを統合するマージ処理部，１８はマージ処理部７０の結果得られるマージ結果，８０はマージ結果１８に対して森林領域の修正を行う森林領域修正処理部，１９は最終的に得られる市街地領域である。
【００１１】
Ｒ画像１０は，人工衛星や航空機を用いて，上空から地表を可視光の赤色領域の波長をもつ電磁波で撮影した画像である。
【００１２】
ＮＩＲ画像１１は，人工衛星や航空機を用いて，上空から地表を近赤外領域の波長をもつ電磁波で撮影した画像である。一般的に地表撮影用の人工衛星や航空機には，Ｒ画像やＮＩＲ画像を撮影するための機器が搭載されている。
【００１３】
Ｐ画像１２は，人工衛星や航空機を用いて，上空から地表を可視光で撮影した画像である。一般的にＰ画像は，Ｒ画像やＮＩＲ画像に比べて解像度が高いため，本実施の形態では，Ｐ画像から市街地領域を抽出する方法を説明する。Ｐ画像の代わりにＲ画像やＮＩＲ画像，または別の波長領域で撮影した画像を市街地領域抽出対象の画像とすることもできる。
【００１４】
水領域抽出処理部２０は，同一地点を撮影したＲ画像１０とＮＩＲ画像１１とから，水領域を抽出する処理手段である。はじめにアフィン変換などで２枚の画像の位置合わせが必要となる場合も考えられる。位置合わせの後，対応する各点において，
Ｘ＝（ＮＩＲ画像における輝度値）／（Ｒ画像における輝度値）
を計算する。水の反射率は可視光の赤色領域と近赤外領域とで大きく異なり，近赤外領域での反射率は可視光赤色領域の反射率よりたいへん小さい。したがって，水領域では上記Ｘは小さな値をとる。そこでＸの値が一定値以下である場合，該点を水領域であると判定する。この一定値は，照明条件やカメラ特性によって大きく値が変動するが，例えば，
ＮＩＲ／Ｒ＜０．７
の場合に，その画素を水領域と判定する。
【００１５】
また上記Ｘ以外にも，ＮＩＲ画像１１での輝度値とＲ画像１０での輝度値との差や，下記で説明するＮＤＶＩを水領域の判定尺度とする方法も考えられる。
【００１６】
森林領域抽出処理部３０は，水領域抽出処理部２０と同様に，同一地点を撮影したＲ画像１０とＮＩＲ画像１１とから森林領域を抽出する処理手段である。ここでもアフィン変換などで２枚の画像の位置合わせが必要となる場合も考えられる。その後，対応する各点において，
ＮＤＶＩ＝（ＮＩＲ画像における輝度値−Ｒ画像における輝度値）
／（ＮＩＲ画像における輝度値＋Ｒ画像における輝度値）
を計算する。森林領域では，クロロフィルの分光反射率特性より上記ＮＤＶＩの値が大きくなる（参照：“図解リモートセンシング”，日本測量協会）。ＮＤＶＩの値が一定値以上である場合，該点を森林領域と判定する。
【００１７】
森林領域１３は，森林領域抽出処理部３０により抽出された領域である。森林領域１３に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００１８】
水領域１４は，水領域抽出処理部２０により抽出された領域である。水領域１４に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００１９】
ヒストグラム平坦化処理部４０は，必要に応じてＰ画像に対して輝度変換を行う処理手段であり，本実施の形態では，ヒストグラム平坦化（参照：“コンピュータ画像処理入門”，総研出版，田村秀行監修）を行う。この処理により建物の屋根部分と影部分とのコントラストが大きくなり，市街地領域のテクスチャ特徴がより強調される。
【００２０】
また本実施の形態では，ヒストグラム平坦化の際に水領域１４を用いる方法について説明する。海や湖など，輝度値の小さい領域が画像中に大きく存在する場合，ヒストグラム平坦化処理では適切なコントラスト強調効果が現れない場合がある。そこで，ヒストグラム平坦化処理に必要な輝度ヒストグラムを，水領域を除外した部分で取得し，またヒストグラム平坦化処理自体も水領域を除外した部分にのみ適用する。このような方法を用いることにより，画像中に海や湖などが大きく存在している場合にも，市街地領域のコントラストが強調される。
【００２１】
ベクトル外積和による市街地領域抽出処理部５０は，ベクトル外積和を特徴量として市街地領域を抽出する処理手段である。図２に，ベクトル外積和による市街地領域抽出処理部５０による処理の流れを示す。
【００２２】
入力画像（Ｐ画像）５１は，ヒストグラム平坦化処理部４０で出力される画像である。エッジベクトル取得処理５２は，画像に差分計算を施してエッジベクトルを取得する処理である。差分計算を行う前に，前処理としてガウシアンフィルタなどでぼかし処理を施すことにより，ノイズの影響を小さくすることができる。
【００２３】
本実施の形態では，横方向と縦方向の差分計算を行う方法について説明する。図３は，そのエッジベクトルの算出方法を説明する図である。横方向の差分計算は，着目画素の左右の画素値を減算することより，縦方向の差分計算は上下の画素値を減算することにより行う。
【００２４】
図３（ａ）に示した例では，画素Ｃ（点Ｃ）における横方向の差分値をＣX ，縦方向の差分値をＣY として，エッジベクトルＣ（ＣX ，ＣY ）は，隣接画素Ａ，Ｂ，Ｄ，Ｅの画素値から，
ＣX ＝（Ｄ−Ｂ）／２
ＣY ＝（Ｅ−Ａ）／２
と計算される。
【００２５】
また，画像の端にある画素については，上下左右の隣接画素が存在しない場合もある。この場合には，隣接画素と着目画素自身との差分を計算して，エッジベクトルの成分とする。このような例を図３（ｂ）に示す。画素Ｇ（点Ｇ）におけるエッジベクトルＧ（ＧX ，ＧY ）は，隣接画素Ｆ，Ｈ，Ｉの画素値から，
ＧX ＝Ｈ−Ｇ
ＧY ＝（Ｉ−Ｆ）／２
と計算される。
【００２６】
エッジベクトル情報５３は，エッジベクトル取得処理５２の結果得られ，全画素に２次元のエッジベクトルが付与されている。図４に，エッジベクトル情報５３の例を示す。図４（ａ）は原画像の画素値，（ｂ）はエッジベクトル成分の値を示している。
【００２７】
ブロック取得処理５４は，着目画素を中心として一定サイズの領域を取得する処理である。図５は，領域のサイズを３画素×３画素とした場合のブロックの例を示している。
【００２８】
外積和計算処理５５は，ブロック取得処理５４で取得したブロック内において，すべての２つのエッジベクトルの組み合わせ間でのベクトル外積を計算し，その絶対値を加算していく処理である。ここでは，図６に示した例のように，９画素で構成されるブロック内のすべての２つのベクトルの組み合わせ，すなわち 9Ｃ2 ＝３６通りの組み合わせ間でベクトル外積の絶対値を計算して加算する。
【００２９】
また，図７（ａ）に示すように，白地に黒い領域と黒字に白い領域とを区別するため，外積を計算する２つのベクトルのうちのどちらか一方のエッジベクトルＶE と，該エッジベクトルＶE を付与されたエッジ点から他方のエッジ点に向かう方向ベクトルＶD との内積を計算し，図７（ｂ）に示すように内積値がプラス値の場合にのみベクトル外積の加算処理を行い，図７（ｃ）に示すように内積値がマイナス値の場合には，それら２つのエッジベクトルの組み合わせにおけるベクトル外積の加算処理をスキップする方法も考えられる（図７（ｄ））。
【００３０】
外積和５６は，外積和計算処理５５の結果得られ，画像各点に特徴量として与えられる。
【００３１】
市街地判定処理５７は，外積和計算処理５５で得られた外積和５６にもとづき，ブロックの代表点を市街地であるかどうかを判定する処理である。代表点はブロックの中心点，ブロックを構成する全点，ブロック中心付近の点など，任意に設定可能であるが，ここではブロックの中心点とする。外積和の値に対してしきい値処理を行って市街地と非市街地領域の二つに分別する。また，しきい値処理を行わずに，外積和の値をそのまま市街地らしさを表す値として出力することも考えられる。
【００３２】
また，海岸線などの誤抽出を抑えるため，水領域抽出処理部２０で抽出した水領域１４を用い水領域の近郊では市街地判定を取り消す方法も考えられる。
【００３３】
終了判定処理５８は，全画素についてブロック取得処理５４，外積和計算処理５５，市街地判定処理５７が行われたかどうかを判定する処理である。まだ未処理である画素が存在する場合には，ブロック取得処理５４に戻り，未処理画素について同様に市街地判定を行う。
【００３４】
市街地領域Ａ（内積チェック無）１５は，ベクトル外積和による市街地領域抽出処理部５０において内積の符号チェックを行わない場合に抽出された市街地領域である。市街地領域Ａ１５に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００３５】
また，市街地領域Ｂ（内積チェック有）１６は，ベクトル外積和による市街地領域抽出処理部５０において内積の符号チェックを行った場合に抽出された市街地領域である。市街地領域Ｂ１６に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００３６】
Ｈｏｕｇｈ特徴による市街地領域抽出処理部６０は，エッジ接線パラメータのＨｏｕｇｈ変換によって算出されるＨｏｕｇｈ特徴量にもとづいて市街地領域を抽出する処理手段である。図８に，Ｈｏｕｇｈ特徴による市街地領域抽出処理部６０の処理の流れを示す。
【００３７】
入力画像（Ｐ画像）６１は，ヒストグラム平坦化処理部４０で出力される画像である。ブロック取得処理６２は，着目画素を中心として一定サイズのブロックを取得する処理である。ブロックの形は正方形，円形など任意であるが，本実施の形態では，円形のブロックを取得するものとする。これはその後に行われるエッジ接線パラメータのＨｏｕｇｈ変換に都合が良いからである。円形ブロックの例を，図９に示す。
【００３８】
エッジ接線パラメータ計算処理６３は，ブロック内の各点におけるエッジ接線を求める処理である。ここでは，まず画像に横方向および縦方向の差分計算を施してエッジベクトルを取得し，取得したエッジベクトルにもとづきエッジ接線パラメータを求める手法を説明する。エッジベクトルを求める処理は，図２で説明したエッジベクトル取得処理５２と同様である。
【００３９】
次に，得られたエッジベクトルからエッジ接線パラメータを計算する。ここでブロックの中心を原点とおき，原点から接線までの距離ｒと，接線の法線ベクトルがＸ方向となす角度θを求めるべきエッジ接線パラメータとする。
【００４０】
図１０は，エッジ接線パラメータの例を説明する図である。図１０から明らかなように，点Ｐにおけるエッジベクトル成分ＰX とＰY とから，点Ｐにおけるエッジ接線のパラメータθP とｒP は簡単に計算される。
【００４１】
Ｈｏｕｇｈ変換処理６４は，算出されたエッジ接線パラメータに対してＨｏｕｇｈ変換を行い，２次元の投票空間に投票していく処理である。図１１は，その投票処理の例を説明する図である。図１０で示した点Ｐでのエッジ接線パラメータθP とｒP は，投票空間内の座標（θP ，ｒP ）に投票される。
【００４２】
投票される値は，＋１の場合もあるが，本実施の形態ではエッジパワー値そのもの，すなわち（ＰX 2 ＋ＰY 2 ）の平方根を投票値とする。この方法により，市街地領域に多く存在するパワーが強いエッジのパラメータが，より強調されることになる。
【００４３】
また，ノイズの悪影響の抑制等を目的として，エッジパワーが一定値以下の場合には投票処理を行わないとか，またはエッジパワーに何らかの重みをつけて投票値とするなどの方法も有効である。
【００４４】
投票空間６５は，Ｈｏｕｇｈ変換処理６４をブロック内の全画素について行った結果，得られる投票空間である。図１２にその例を示す。図１２（ａ）は原画像，図１２（ｂ）は投票空間の例を示している。図１２（ａ）の例では直径５０画素の円形ブロック６９内の全画素に対して投票を行っている。図１２（ｂ）の投票空間では，得票値の大小を濃度を変えて示している。投票空間の横軸は，角度パラメータθの軸（範囲０〜３６０°），縦軸は距離パラメータｒの軸（範囲０〜２５）である。画像内の直線が投票空間内では得票値の高い点として変換されていることが分かる。
【００４５】
ブロック内処理終了判定６６は，ブロック内の全画素についてエッジ接線パラメータ計算処理６３，Ｈｏｕｇｈ変換処理６４が終了したか否かを判定する処理である。未処理画素が存在する場合，再びエッジ接線パラメータ計算処理６３へ戻り，未処理画素の処理を行う。
【００４６】
市街地領域Ｂ（内積チェック有）１６は，図１のベクトル外積和による市街地領域抽出処理部５０で抽出された市街地領域Ｂである。水領域１４は，図１の水領域抽出処理部２０で抽出された水領域である。
【００４７】
市街地判定処理６７は，投票空間６５にもとづき，ブロックの代表点が市街地であるかどうかを判定する処理である。ブロックの代表点として，ブロックの中心画素，ブロック中心付近の画素，ブロックを構成する全画素，ブロック内の一定値以上のエッジパワーを持つ画素などが考えられる。本実施の形態では，２次元の投票空間をそのまま用いず，Ｘ軸方向（角度パラメータ軸方向）に射影をとり，１次元に縮退させたヒストグラム（以下，角度ヒストグラムと呼称する）を用いる方法を説明する。
【００４８】
図１３に，図１２（ｂ）に示す投票空間から作成した角度ヒストグラムを示す。角度ヒストグラムの最大値，最小値，平均値，分散値などの統計量から市街地を求める方法もあるが，ここでは（１）最大値／最小値，（２）エッジ接線直交性にもとづく判定処理について説明する。
【００４９】
（１）最大値／最小値
角度ヒストグラムの最大値と最小値を求め，最大値を最小値で除算して比率を計算する。この値が一定値以上である場合，ブロックの代表点を市街地と判定する。接線パラメータ計算時の原点をブロックの中心に設定しているため，見た目の傾きが同じ直線でもブロック中心に対して反対側に位置するとき，その角度パラメータθは１８０°の差を持つこととなる。そこで角度ヒストグラムの１８０〜３６０°部分を０〜１８０°部分に加算し，角度ヒストグラムの範囲を０〜１８０°とすることで，同一の傾きを持つ直線は同一の角度パラメータを持つように補正する方法も考えられる。
【００５０】
ランダムノイズのようにエッジ方向がばらばらである画像の場合，最大値／最小値は小さな値であり，この値が大きな値となるのはブロック内に強いエッジパワーを持つ直線が存在するなどの場合である。
【００５１】
（２）エッジ接線直交性
ブロック内に直交する直線が存在する場合，角度ヒストグラムには９０°の間隔を持つ２つのピークが発生する。９０°間隔の２つのピークがあるか否かは，次の方法により判定することができる。図１４に，その判定方法の説明図を示す。
【００５２】
図１４に示すヒストグラムは，上記（１）において説明したように１８０〜３６０°部分を０〜１８０°部分に重ねて角度ヒストグラムの角度範囲を０〜１８０゜としたものである。はじめにこのヒストグラムにおいて最大値を取る角度θ1 を求める。次に，θ1 から９０°の差を持つ角度θ2 を中心として任意の幅を持つ第２ピーク探索領域ＳA を設定し，その領域内での最大値（第２ピーク値）を求める。図１４の例では，θ2 を中心として±１０°の幅を第２ピーク探索領域ＳA の範囲と設定した。
【００５３】
第２ピーク値を求めた後，この値をヒストグラムの最小値で除算し，その比率が一定値以上の場合，該ブロックの代表点を市街地領域であると判定する。第２ピーク値と最小値との比率以外にも，第２ピーク値そのものの値，第２ピーク値と最小値との差分などを判定尺度としたり，最大値そのものの値も判定尺度に加えることも考えられる。
【００５４】
また，範囲を０〜１８０°とした角度ヒストグラムの，さらに９０〜１８０°部分を０〜９０°部分に重ね，角度範囲を０〜９０°としたヒストグラムを作成し，該ヒストグラムの最大値と最小値との比率を判定尺度とすることも考えられる。
【００５５】
ここでは，角度ヒストグラムにもとづき市街地判定する方法を説明したが，２次元の投票空間をそのまま使用して，適当な特徴量を算出する方法も考えられる。
【００５６】
また，水田や畑などが建造物として誤抽出されるのを抑制するため，市街地領域Ｂ（内積チェック有）１６を用い，上記判定と併せた再判定を行う方法もある。建造物は黒地に白い領域として存在する傾向が強いのに対し，水田や畑などは白地に黒い領域であるため，市街地領域Ｂ（内積チェック有）１６では抽出されない。そこで上記方法により市街地と判定されたブロック内に，市街地領域Ｂ（内積チェック有）１６が存在しない場合に，該ブロックは非市街地領域として前の判定結果をリジェクトする。
【００５７】
また，橋や埋め立て地などの誤抽出を抑えるため，水領域１４を用い水領域の近郊における市街地判定をリジェクトする方法も考えられる。
【００５８】
終了判定処理６８は，全画素について市街地判定処理が行われたか否かを判定する処理である。まだ未処理である画素が存在する場合，ブロック取得処理６２に再度戻り，未処理画素について同様に処理を行う。
【００５９】
市街地領域Ｃ１７は，Ｈｏｕｇｈ特徴による市街地領域抽出処理６０の結果抽出された市街地領域である。市街地領域Ｃ１７に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００６０】
マージ処理部７０は，ベクトル外積和を特徴量として抽出した市街地領域Ａ１５と，Ｈｏｕｇｈ特徴にもとづき抽出した市街地領域Ｃ１７とをマージする処理手段である。本実施の形態では，市街地領域Ａ１５と市街地領域Ｃ１７とを２値化して加算する単純な方法で行うこととする。本方法以外にも，市街地領域Ａ１５と市街地領域Ｃ１７の市街地らしさを多値で表し，その値に何らかの重みを付けて加算するなどの方法も考えられる。
【００６１】
マージ結果１８は，マージ処理部７０の結果，取得される結果である。マージ結果１８に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００６２】
森林領域修正処理部８０は，マージ結果１８に対して森林領域１３を除外し，最終的な市街地領域１９を抽出する処理手段である。除外する方法として，本実施の形態では，森林領域１３およびマージ結果１８を２値化して減算するといった単純な方法で行うこととする。
【００６３】
市街地領域１９は，森林領域修正処理部８０の結果取得される最終的な市街地領域の抽出結果である。この市街地領域１９に対して，必要に応じてノイズ除去処理，ぼかし処理，最大値フィルタ処理，最小値フィルタ処理などの画像処理を施す方法もある。
【００６４】
以上の市街地領域抽出処理装置１が行う処理は，コンピュータとソフトウェアプログラムとによって実現することができ，そのプログラムは，コンピュータが読み取り可能な可搬媒体メモリ，半導体メモリ，ハードディスク等の適当な記録媒体に格納して，そこから読み出すことによりコンピュータに実行させることができる。また，そのプログラムは通信回線を経由して他のコンピュータからダウンロードすることができ，それをインストールして実行させることもできる。
【００６５】
【発明の効果】
以上説明したように，本発明は，市街地領域でのテクスチャ特徴とマルチスペクトル特徴，特に市街地領域ではエッジのパワーが強く，かつエッジ接線方向が直交している傾向が強いというテクスチャ特徴に着目した市街地領域抽出を行う。家などの小さな建造物が集合している市街地領域の抽出に有効なベクトル外積和を特徴量とする市街地領域抽出と，工場などの大規模建造物の抽出に有効なＨｏｕｇｈ特徴にもとづく市街地領域抽出とを組み合わせることにより，建造物の大きさに依存しない安定した市街地領域抽出を行うことが可能である。また，マルチスペクトル情報を用いて森林領域と水領域を抽出することにより，より精度の高い市街地領域抽出が可能となった。
【図面の簡単な説明】
【図１】本実施の形態による市街地領域抽出処理装置の構成例を示す図である。
【図２】ベクトル外積和による市街地領域抽出処理部による処理の流れを示す図である。
【図３】エッジベクトルの算出方法を説明する図である。
【図４】エッジベクトル情報の例を示す図である。
【図５】ブロックの例を示す図である。
【図６】ベクトル外積和の計算の例を示す図である。
【図７】内積符号チェック方法の例を説明する図である。
【図８】Ｈｏｕｇｈ特徴による市街地領域抽出処理の流れを示す図である。
【図９】円形ブロックの例を示す図である。
【図１０】エッジ接線パラメータの例を説明する図である。
【図１１】投票処理の例を説明する図である。
【図１２】投票空間の例を示す図である。
【図１３】角度ヒストグラムの例を示す図である。
【図１４】第２ピーク探索処理を説明する図である。
【符号の説明】
１ 市街地領域抽出処理装置
１０ Ｒ画像
１１ ＮＩＲ画像
１２ Ｐ画像
１３ 森林領域
１４ 水領域
１５ 市街地領域Ａ（内積チェック無）
１６ 市街地領域Ｂ（内積チェック有）
１７ 市街地領域Ｃ
１８ マージ結果
１９ 市街地領域
２０ 水領域抽出処理部
３０ 森林領域抽出処理部
４０ ヒストグラム平坦化処理部
５０ ベクトル外積和による市街地領域抽出処理部
６０ Ｈｏｕｇｈ特徴による市街地領域抽出処理部
７０ マージ処理部
８０ 森林領域修正処理部
Claims (7)
- 地表を撮影した画像の中から，コンピュータにより市街地領域を自動的に抽出する処理方法であって，
前記コンピュータが，
可視光の赤色領域の波長をもつ電磁波をとらえたＲ画像と，近赤外領域の電磁波をとらえたＮＩＲ画像とを用い，対応する各点においてＮＩＲ画像の輝度値とＲ画像の輝度値とから算出される値により，水領域または森林領域を抽出する第１の処理過程と，
前記第１の処理過程の結果を使用して，市街地領域抽出対象の画像のコントラストを大きくする輝度変換を行う第２の処理過程と，
輝度変換を行った画像各点について２次元のエッジベクトルを算出し，各画素を中心とする局所領域のブロック内において，全エッジベクトルのすべての２ベクトルの組み合わせについてベクトル外積の絶対値を計算して加算し，加算の結果得られるベクトル外積和にもとづきブロック内の代表点が市街地であるか否かを判定する第３の処理過程と，
市街地領域抽出対象の画像各点を中心とするブロックを設定し，該ブロック内の各点においてエッジ接線をそれぞれ求め，算出されたエッジ接線の角度パラメータと距離パラメータとに対してハフ変換を行い，投票空間用のメモリ内に該２つのパラメータを２軸とする２次元の投票空間を作成し，得られた投票空間にもとづき該ブロックの代表点が市街地であるか否かを判定する第４の処理過程と，
前記第１，第３および第４の処理過程による結果から最終的な市街地領域を抽出する処理過程とを有する
ことを特徴とする市街地領域抽出処理方法。 - 請求項１に記載の市街地領域抽出処理方法において，
前記第３の処理過程でエッジベクトルの外積を計算するとき，２つのベクトルのうちどちらか一方のエッジベクトルと，該エッジベクトルを付与されたエッジ点から他方のエッジ点に向かう方向ベクトルとの内積を計算し，
内積値がマイナス値の場合にはこのエッジベクトルの組み合わせにおけるベクトル外積の加算処理をスキップすること
を特徴とする市街地領域抽出処理方法。 - 請求項１または請求項２に記載の市街地領域抽出処理方法において，
前記第４の処理過程でエッジ接線パラメータに対してハフ変換を行う際には，エッジパワー値を投票値とする
ことを特徴とする市街地領域抽出処理方法。 - 請求項１，請求項２または請求項３に記載の市街地領域抽出処理方法において，
前記第４の処理過程において２次元の投票空間を角度パラメータ軸方向に射影をとり１次元のヒストグラムに縮退させ，
該ヒストグラムにもとづき該ブロックの代表点が市街地であるか否かを判定する
ことを特徴とする市街地領域抽出処理方法。 - 地表を撮影した画像の中から，コンピュータにより市街地領域を自動的に抽出する処理装置であって，
可視光の赤色領域の波長をもつ電磁波をとらえたＲ画像と，近赤外領域の電磁波をとらえたＮＩＲ画像とを用い，対応する各点においてＮＩＲ画像の輝度値とＲ画像の輝度値とから算出される値により，水領域または森林領域を抽出する第１の手段と，
前記第１の手段の結果を使用して，市街地領域抽出対象の画像のコントラストを大きくする輝度変換を行う第２の手段と，
輝度変換を行った画像各点について２次元のエッジベクトルを算出し，各画素を中心とする局所領域のブロック内において，全エッジベクトルのすべての２ベクトルの組み合わせについてベクトル外積の絶対値を計算して加算し，加算の結果得られるベクトル外積和にもとづきブロック内の代表点が市街地であるか否かを判定する第３の手段と，
市街地領域抽出対象の画像各点を中心とするブロックを設定し，該ブロック内の各点においてエッジ接線をそれぞれ求め，算出されたエッジ接線の角度パラメータと距離パラメータとに対してハフ変換を行い，該２つのパラメータを２軸とする２次元の投票空間を作成し，得られた投票空間にもとづき該ブロックの代表点が市街地であるか否かを判定する第４の手段と，
前記第１，第３および第４の手段による結果から最終的な市街地領域を抽出する手段とを備える
ことを特徴とする市街地領域抽出処理装置。 - 請求項１から請求項４までのいずれか一項に記載の市街地領域抽出処理方法を，コンピュータに実行させるための市街地領域抽出処理プログラム。
- 請求項１から請求項４までのいずれか一項に記載の市街地領域抽出処理方法を，コンピュータに実行させるためのプログラムを記録した
ことを特徴とする市街地領域抽出処理プログラムのコンピュータ読み取り可能な記録媒体。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2002091060A JP3816825B2 (ja) | 2002-03-28 | 2002-03-28 | 市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2002091060A JP3816825B2 (ja) | 2002-03-28 | 2002-03-28 | 市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2003288582A JP2003288582A (ja) | 2003-10-10 |
| JP3816825B2 true JP3816825B2 (ja) | 2006-08-30 |
Family
ID=29236240
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2002091060A Expired - Lifetime JP3816825B2 (ja) | 2002-03-28 | 2002-03-28 | 市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体 |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP3816825B2 (ja) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104881659A (zh) * | 2015-06-12 | 2015-09-02 | 中国地质大学（武汉） | 一种不透水层的提取方法及装置 |
Families Citing this family (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP4512778B2 (ja) * | 2004-04-26 | 2010-07-28 | 株式会社エヌ・ティ・ティ・データ | 不法投棄箇所検知装置、方法、およびプログラム |
| JP6727894B2 (ja) * | 2016-04-13 | 2020-07-22 | キヤノン株式会社 | 処理装置、処理方法、及びプログラム |
-
2002
- 2002-03-28 JP JP2002091060A patent/JP3816825B2/ja not_active Expired - Lifetime
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104881659A (zh) * | 2015-06-12 | 2015-09-02 | 中国地质大学（武汉） | 一种不透水层的提取方法及装置 |
| CN104881659B (zh) * | 2015-06-12 | 2018-05-01 | 中国地质大学（武汉） | 一种不透水层的提取方法及装置 |
Also Published As
| Publication number | Publication date |
|---|---|
| JP2003288582A (ja) | 2003-10-10 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Zhang et al. | Object-oriented shadow detection and removal from urban high-resolution remote sensing images | |
| CN106934803B (zh) | 电子器件表面缺陷的检测方法及装置 | |
| CN109961399B (zh) | 基于图像距离变换的最佳缝合线搜寻方法 | |
| CN109784344A (zh) | 一种用于地平面标识识别的图像非目标滤除方法 | |
| WO2009131108A1 (ja) | 建物屋根輪郭認識装置、建物屋根輪郭認識方法、及び建物屋根輪郭認識プログラム | |
| CN112184604B (zh) | 一种基于图像融合的彩色图像增强方法 | |
| US20130215234A1 (en) | Method and apparatus for stereo matching | |
| CN112614136A (zh) | 一种红外小目标实时实例分割方法及装置 | |
| CN116681636B (zh) | 基于卷积神经网络的轻量化红外与可见光图像融合方法 | |
| Zhao et al. | Infrared and visible image fusion algorithm based on saliency detection and adaptive double-channel spiking cortical model | |
| CN103310439A (zh) | 一种基于尺度空间的图像最稳极值区域检测方法 | |
| Palenichka et al. | Multiscale isotropic matched filtering for individual tree detection in LiDAR images | |
| Steinnocher | Adaptive fusion of multisource raster data applying filter techniques | |
| JP3816825B2 (ja) | 市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体 | |
| JP5648546B2 (ja) | 通路検出装置、方法、及びプログラム | |
| Aytekin et al. | Automatic and unsupervised building extraction in complex urban environments from multi spectral satellite imagery | |
| Xiao et al. | A Novel Image Fusion Method for Water Body Extraction Based on Optimal Band Combination. | |
| JP2005537562A (ja) | 皮膚紋イメージ処理方法 | |
| JP3816826B2 (ja) | 市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体 | |
| CN114219903A (zh) | 一种基于视频图像的河道网拓扑结构生成方法及其装置 | |
| Raja | Tree crown delineation from high resolution satellite images | |
| Saberi et al. | Multi-temporal satellite image correction based on CNN features and adaptive contrast enhancement | |
| JP3816827B2 (ja) | 市街地領域抽出処理方法，市街地領域抽出処理装置，市街地領域抽出処理プログラムおよびそのプログラムの記録媒体 | |
| Goud et al. | Evaluation of image fusion of multi focus images in spatial and frequency domain | |
| Wang et al. | Shadow Detection and Reconstruction of High-Resolution Remote Sensing Images in Mountainous and Hilly Environments |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20040121 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20060131 |
|
| A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20060317 |
|
| RD04 | Notification of resignation of power of attorney |
Free format text: JAPANESE INTERMEDIATE CODE: A7424 Effective date: 20060317 |
|
| TRDD | Decision of grant or rejection written | ||
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20060606 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20060608 |
|
| R150 | Certificate of patent or registration of utility model |
Free format text: JAPANESE INTERMEDIATE CODE: R150 Ref document number: 3816825 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
| RD02 | Notification of acceptance of power of attorney |
Free format text: JAPANESE INTERMEDIATE CODE: R3D02 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20090616 Year of fee payment: 3 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20100616 Year of fee payment: 4 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20100616 Year of fee payment: 4 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110616 Year of fee payment: 5 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20120616 Year of fee payment: 6 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20130616 Year of fee payment: 7 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20140616 Year of fee payment: 8 |
|
| S531 | Written request for registration of change of domicile |
Free format text: JAPANESE INTERMEDIATE CODE: R313531 |
|
| R350 | Written notification of registration of transfer |
Free format text: JAPANESE INTERMEDIATE CODE: R350 |
|
| EXPY | Cancellation because of completion of term |