JP4058293B2 - レーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法及び高精度都市モデルの生成システム並びに高精度都市モデルの生成のプログラム - Google Patents
レーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法及び高精度都市モデルの生成システム並びに高精度都市モデルの生成のプログラム Download PDFInfo
- Publication number
- JP4058293B2 JP4058293B2 JP2002127672A JP2002127672A JP4058293B2 JP 4058293 B2 JP4058293 B2 JP 4058293B2 JP 2002127672 A JP2002127672 A JP 2002127672A JP 2002127672 A JP2002127672 A JP 2002127672A JP 4058293 B2 JP4058293 B2 JP 4058293B2
- Authority
- JP
- Japan
- Prior art keywords
- area
- dsm
- image
- stereo
- mask
- 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
- Instructional Devices (AREA)
- Processing Or Creating Images (AREA)
- Image Processing (AREA)
- Image Analysis (AREA)
Description
【０００１】
【発明の属する技術分野】
レーザスキャナデータと空中写真画像を用いて、空中写真画像の建物面を容易に精度よく３次元化できる高精度都市モデルの生成装置に関する。
【０００２】
【従来の技術】
高精度の都市モデルを生成する手法は、一般にレーザスキャナデータを用いる手法とステレオマッチングによる手法とがある。
【０００３】
例えば、レーザスキャナデータは、飛行機又はヘリコプター等の航空機から地表にレーザ光をスキャニングし、これらの受光データを解析することで、図２１に示すようなレーザデータを得る。この図２１は、図２１（ｂ）に対象エリア（航空写真で示している）を、図２１（ａ）にそのレーザデータを示している。図２０（ａ）に示すように、レーザデータは細かな点列になっており、それぞれが３次元座標（ｘ、ｙ、ｚ）を有している。
【０００４】
従って、レーザスキャナデータは、３次元座標を有するので航空写真のステレオマッチングなどに比べて、図２２に示す標高精度の高いＤＳＭ（Ｄｉｇｉｔａｌ Ｓｕｒｆａｃｅ Ｍｏｄｅｌの略で、建物などの地物を含む地表面の形状を表現するモデル）をとることができる。
【０００５】
一方、航空写真により自動標高計測を行うステレオマッチング手法（図２０を参照）もある。このステレオマッチングは、画像間の各画素の対応関係を求め、三角測量の原理で各対応点位置での３次元座標を求めるもので、高密度に標高データを計測することができる。
【０００６】
【発明が解決しようとする課題】
しかしながら、レーザスキャナーデータによる手法は標高に関する精度の高いデータが得られるが、３次元点群から都市モデルのような直線構造を多く持つデータを生成する場合には、境界部分を正確に取るために、データをさらに高密度化する必要がある。
【０００７】
つまり、レーザデータは細かな点列（３次元点群）であり、これ以上の解像度を得ることは容易ではない。
【０００８】
また、レーザスキャナデータを高密度化すると、データ取得にかかるコスト及びデータ容量が莫大になり、コスト面とデータハンドリングの両面で問題がある。
【０００９】
一方、ステレオマッチング手法に用いる画像は高分解能であり、人工構造物の境界や領域の識別に有効である。しかしながら、ステレオマッチングは隠れなどの影響により標高変化の激しい箇所での標高計測精度が低く、境界を直接３次元的に計測するのに問題が生じる。
【００１０】
本発明は、レーザスキャナデータの標高精度の利点と、航空写真の高分解能であることの利点、および既往の疎密度ＤＥＭ（Ｄｉｇｉｔａｌ Ｅｌｅｖａｔｉｏｎ Ｍｏｄｅｌの略で、地盤高を表現するモデル。概略地盤高がわかるので地盤高として利用できる）を組み合わせ、精度の高い都市３次元モデルを構築することを目的とする。
【００１１】
【課題を解決するための手段】
本発明のレーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法は、
所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶手段に記憶する工程と、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理工程と、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング工程と、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理工程と、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す工程と、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する工程と、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する工程と、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する工程と
を行うことを要旨とする。
【００１２】
本発明の高精度都市モデルの生成システムは、所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶した記憶手段とを備え、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理手段と、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング手段と、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理手段と、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す手段と、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する手段と、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する手段と、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する手段とを備えたことを要旨とする。
【００１３】
本発明の高精度都市モデルの生成のプログラムは、所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶した記憶手段とを備え、
コンピュータに、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理手段、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング手段、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理手段、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す手段、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する手段、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する手段、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する手段、
としての機能を実行させる。
【００１４】
【発明の実施の形態】
本実施の形態は、レーザスキャナデータにステレオ撮影された空中写真画像データから得られる領域情報・線構造情報を援用することによって、人工構造物の境界を精度よく抽出し、自動建物多角形化および３次元モデル化を行う。
【００１５】
まず、ステレオ空中写真画像の各々を領域分割する事により、色彩が均一な画像領域を得る。このような画像領域は、建物の屋上及び路面・グランド・その他広い平坦地に対応していると考えられる。
【００１６】
各色領域に重なるレーザスキャナの３次元点群の標高値を初期値として、ステレオマッチングを行う。その結果得られる標高値の計算結果が初期値とほぼ同じ高さになった場合はその領域に計算結果の標高値を与え、ＤＳＭ（Ｄｉｇｉｔａｌ Ｓｕｒｆａｃｅ Ｍｏｄｅｌの略で、建物等の地物を含む地表面の形状を表現するモデル）を生成する。これにより、レーザスキャナデータでは曖昧であった地物境界をより高精度に表現することができるようにしている。
【００１７】
続いて上記で得られたＤＳＭの高い方から順に水平な領域を領域分割によって認識し、各領域に対応するステレオ画像上の周辺領域で直線抽出を行い、その結果を多角形化する。既往のＤＥＭや地域の平均的地盤高に対するＤＳＭの差が任意に設定する閾値以下にある領域は、建物以外の領域であり、地盤高と推定する。これによって、建物領域を３次元多角形化することができる。更に多角形化された建物領域を角柱状に構成する事によって、３次元モデルを得る事が可能となるものである。
【００１８】
本発明を実現する３次元計測装置の構成を図１に示す。３次元計測システム１は、パーソナルコンピュータ本体２、ステレオ画像等を格納する外部記憶装置３、ディスプレイ４という一般的なパーソナルコンピュータからなる。また、次のことが前提として与えられている。
【００１９】
（１）入力データ
次の各種データが、外部記憶装置３に格納されているとする。
【００２０】
（ａ）ステレオ空中写真画像（ＩL およびＩR （図２０））とその標定要素と、
（ｂ）レーザスキャナデータが保存されている。このレーザスキャナデータは３次元座標の集まりであり、以下ではＬＰi （i=1...nLP）で表す（図２１（ａ））。
【００２１】
（ｃ）デジタル地形モデル（ＤＥＭ））。地盤高を表し、国土数値情報などの所与のデータを用いてよい（図２２）。
【００２２】
（２）計算モジュール
次の計算モジュールが与えられているものとする。この計算モジュールは、座標計算モジュール１０と、地盤高計算モジュール１１等からなる。
【００２３】
座標計算モジュール１０は、ステレオ空中写真画像の標定要素を読み取り、３次元座標と、ステレオ画像Ｉ1 およびＩ2 の座標との間の変換を行う。これは特願２０００−２５１４５６で座標変換関数・逆座標変換関数として与えられているものと同じとする。
【００２４】
つまり、輝度画像ＩＬ、ＩＲについて、地上の３次元座標が投影される各画像の画像座標算出のための座標変換関数ＦＸ、逆座標変換関数関数ＧＸ（総称して標定情報ともいう）を用いる。
【００２５】
例えば、３次元空間内の点（ｘ，ｙ，ｚ）が投影される左右画像のピクセル座標（ｘ１ｐｉｘ，ｙ１ｐｉｘ）および（ｘ２ｐｉｘ，ｙ２ｐｉｘ）を計算する。
【００２６】
すなわち、次の４つの座標変換関数ＦＸ１（ｘ，ｙ，ｚ）、ＦＹ１（ｘ，ｙ，ｚ）、ＦＸ２（ｘ，ｙ，ｚ）、ＦＹ２（ｘ，ｙ，ｚ）を解析的、もしくは手続き的に決定する。
【００２７】
ｘ１ｐｉｘ＝ＦＸ１（ｘ，ｙ，ｚ） 式（１）
ｙ１ｐｉｘ＝ＦＹ１（ｘ，ｙ，ｚ） 式（２）
ｘ２ｐｉｘ＝ＦＸ２（ｘ，ｙ，ｚ） 式（３）
ｙ２ｐｉｘ＝ＦＹ２（ｘ，ｙ，ｚ） 式（４）
また、ｚ＝一定として式（１），式（２）を解析的又は手続的に解くと、各画像上の点がある標高ｚにあると仮定した時の地上でのｘ座標とｙ座標を求めることができる。
【００２８】
ｘ＝ＧＸ１（ｘ１ｐｉｘ，ｙ１ｐｉｘ，ｚ） 式（５）
ｙ＝ＧＹ１（ｘ１ｐｉｘ，ｙ１ｐｉｘ，ｚ） 式（６）
同様に式（３），式（４）について解くと、
ｘ＝ＧＸ２（ｘ２ｐｉｘ，ｙ２ｐｉｘ，ｚ） 式（７）
ｙ＝ＧＹ２（ｘ２ｐｉｘ，ｙ２ｐｉｘ，ｚ） 式（８）
である。
【００２９】
このような関数は、例えば図２に示す不定形窓を用いているときに使用する。この図２は、航空機によって得た画像Ｉ１を画面に表示し、この画像上でマウス操作等で地物Ｂｕの輪郭に相当する閉領域Ｂ１を与える。そして、この閉領域に投影が一致する位置に不定形窓Ｐｌｈを発生させ、標高Ｈ（ｚ）を変化させる。
【００３０】
例えば３０ｍ、８０ｍ、１５０ｍの順に標高を変化させると、３０ｍ、８０ｍ、１５０ｍと変化したときの輝度画像Ｉ１、Ｉ２における不定形窓Ｐｌｈの輝度画像Ｉ１、Ｉ２の投影位置をその都度求めながら、この位置における不定形窓Ｐｌｈを投影したときの輝度画像Ｉ１、Ｉ２における輝度値を比較し、評価関数によって最も一致（近似）したと判定されるときの不定形窓Ｐ１ｈの高さ（標高値）としている。
【００３１】
地盤高計算モジュール１１は、ＤＥＭを読み取り、適当な内挿計算を用いることによって、与えられた水平座標位置における地盤高を計算する。
【００３２】
Ｈ＝ＧＥlv（ｘ，ｙ）
例えばＤＥＭが正方メッシュデータならば、バイリニア法等の内挿方法、ＤＥＭがランダムな配置の３次元点群ならば、三角メッシュを用いた内挿法を用いる。
【００３３】
すなわち、本実施の形態では、座標計算モジュール１０、地盤高計算モジュール１１等によって、レーザスキャナデータにステレオ撮影された空中写真画像データから得られる領域情報・線構造情報を援用して、例えば、人工構造物の境界を精度よく抽出し、自動建物多角形化および３次元モデル化を行う。
【００３４】
この３次元化モデルの生成においては、まず、ステレオ空中写真画像の各々を領域分割する事により、色彩が均一な画像領域を得る。
【００３５】
このような画像領域は、建物の屋上及び路面・グランド・その他広い平坦地に対応していると考えられる。
【００３６】
各色領域に重なるレーザスキャナの３次元点群の標高値を初期値として、ステレオマッチングを行い、その結果得られる標高値の計算結果が初期値とほぼ同じ高さになった場合はその領域に計算結果の標高値を与え、ＤＳＭを生成する。これにより、レーザスキャナデータでは曖昧であった地物境界をより高精度に表現することができる。
【００３７】
続いて上記で得られたＤＳＭの水平な領域を領域分割によって認識し、各領域に対応するステレオ画像上の周辺領域で直線抽出を行い、その結果を多角形化する。
【００３８】
既往のＤＥＭや地域の平均的地盤高に対するＤＳＭの差が任意に設定する閾値以下にある領域は、建物以外の領域であり、地盤高と推定する。これによって、建物領域を３次元多角形化することができる。更に多角形化された建物領域を角柱状に構成する事によって、３次元モデルを得る事ができる。
【００３９】
本発明の実施形態の概略を図１を用いて簡単に説明する。図１に示すように、次の５つの部分から成り立っている。
【００４０】
初めに、データベース（図示せず）に保存されているステレオ空中写真画像（ＩＲ、ＩＬ）を読み込んで、ステレオ空中写真画像の各々を領域分割（同じ色を一塊りとする領域：同じ色彩の画素が連続する領域）する事により、色彩が画像領域（建物の屋上及び路面・グランド・その他広い平坦地に対応していると考えられる）を得るステレオ画像の色彩による領域分割処理部１５と、各分割領域（分割色領域ともいう）に重なるレーザスキャナデータの３次元点群の標高値を初期値として、ステレオマッチングを行う水平面ステレオマッチング処理部１６と、この水平面ステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになった場合はその領域に計算結果の標高値を与え、ＤＳＭ（Ｄｉｇｉｔａｌ Ｓｕｒｆａｃｅ Ｍｏｄｅｌの略で、建物等の地物を含む地表面の形状を表現するモデル）を生成するＤＳＭの領域分割処理部１７と、多角形生成処理部１８と、３次元モデル化処理部１９等とからなる。
【００４１】
以下にこれらの各処理（ステレオ画像の色彩による領域分割処理、水平面ステレオマッチング処理、ＤＳＭの分割処理、多角形生成処理、３次元モデル化処理）の詳細をフローチャートを用いて説明する。
【００４２】
（ステレオ画像の色彩による領域分割処理部１５）
ステレオ画像の色彩による領域分割について図３のフローチャートを用いて詳細に説明する。
【００４３】
初めに、領域分割に必要なパラメータを設定する（Ｓ１０）。例えば、反復型領域拡張法によって領域分割を行う場合は、次の３つのパラメータがユーザにより指定されているものとする。
【００４４】
・初期しきい値
・しきい値の増分
・反復打ち切りのためのしきい値
次に、ステレオ空中写真画像の内、左画像ＩLを外部記憶装置３から内部に読み込む（Ｓ１１）。
【００４５】
次に、画像の輝度情報や色情報を用いて領域分割処理を行う（図４参照）。この領域分割処理では、まずメディアンフィルタなどで小さな孤立点を除去した後、除去画像の輝度や色が均一な連続領域ごとに識別番号を付ける。
【００４６】
なお、領域分割の方法としては、反復型領域拡張法など一般に広く行われている領域分割方法を用いてもよい。
【００４７】
そして、領域分割の結果（単色で塗られた領域）、各画素が何番目の識別番号の領域に属するかがわかるラベル番号が与えられるものとする。この左画像のラベル番号を格納する配列をＬＬ（ｉ，ｊ）とする。
【００４８】
すなわち、ＬＬ（ｉ，ｊ）＝ｎならば、画像ｉカラムｊラインの画素が、ラベル番号ｎの領域に属することを意味する。
【００４９】
次に、領域分割によって求められたラベル格納配列ＬＬ（ｉ，ｊ）を左画像のラベルとして外部記憶装置３に格納する領域分割結果保存処理を行う（Ｓ１３）。
【００５０】
そして、左画像読込（Ｓ１１）〜領域分割結果保存処理（Ｓ１３）と同様の処理を右画像ＩR について同様に行う右画像読込処理（Ｓ１４）と領域分割処理（Ｓ１５）とを行い、右画像のラベル格納配列ＬＲ（ｉ，ｊ）を外部記憶装置３に格納する（Ｓ１６）。
【００５１】
つまり、Ｓ１４〜Ｓ１６は右画像読み込み・領域分割・領域分割結果保存を行っている。
【００５２】
次に、水平面ステレオマッチング処理について説明する。
【００５３】
（水平面ステレオマッチング）
領域分割された各画像領域が水平であると仮定し、その拘束の元で領域分割された部分を「不定形窓を用いた画像間拡張イメージマッチング方法（特願２０００−２５１４５６）」を行い、デジタル地表面モデル（ＤＳＭ）を生成する。ここでのＤＳＭとは、地表面の標高をメッシュ点上で表現するものである。例えば、メッシュの原点を（ｘ０，ｙ０）、メッシュ間隔をｄとすると、ＤＳＭ配列ｈ（ｉｘ，ｉｙ）は、座標（ｘ０＋ｄ・ｉｘ、ｙ０＋ｄ・ｉｙ）での標高を与えるものとする。また、ＤＳＭ配列の原点や大きさは、ユーザが与える対象範囲にあわせて設定するものとする。前述の水平面ステレオマッチングは、各分割領域を不定形窓として、これを変化させて拡張ステレオモデルの両方の画像上の投影位置を求め、この投影位置における不定形窓で区切られた領域の輝度値をそれぞれ求めて輝度値を比較し、比較結果が最も一致するときの不定形窓の高さを拡張ステレオモデルにおける地物の高さとする拡張ステレオマッチングによる標高付与処理を行う。
【００５４】
図５は平面ステレオマッチングを説明するフローチャートである。初めに、左右画像をメモリ上に読み込む（Ｓ２０）。
【００５５】
次に、ＤＳＭ配列を格納するＤＳＭ配列を初期化する（Ｓ２１）。初期化するＤＳＭ配列は次の３つとする。
【００５６】
・ｈＬ（ｉｘ，ｉｙ）：左画像領域分割結果から生成されるＤＳＭ
・ｈＲ（ｉｘ，ｉｙ）：右画像領域分割結果から生成されるＤＳＭ
・ｈＬＲ（ｉｘ，ｉｙ）：左右画像領域分割結果から生成されるＤＳＭ
なお、ＤＳＭ配列の原点や範囲は、ユーザがあらかじめ与えるものとする。初期化にあたっては、すべての配列の値を一定値ｈｉｎｉｔとしておく。ｈｉｎｉｔは、実際にはありえない小さな値（−９９９９など）であり、ＤＳＭ配列がｈｉｎｉｔであるということは、平面ステレオマッチングによってＤＳＭが計算されていないことを示す。
【００５７】
次に、左右画像の画像領域の分割結果を読込む（Ｓ２２）。つまり、画像領域分割結果をメモリ上の配列Ｌ（ｉ，ｊ）に読み込む。
【００５８】
そして、各領域毎Ｌ（ｉ，ｊ）に記録されているラベル番号が付けられた各々の領域について、番号順に次の処理を行う。
【００５９】
この処理は、レーザスキャナデータの重ね合わせと概略標高の計算を行う（Ｓ２３）。
【００６０】
例えば、レーザスキャナデータＬpiの画像への投影位置を座標計算モジュールによって計算し、その中で処理中の領域に重なるものを検索する（図６を参照）。 次に、更に重なった点のヒストグラムを作り、ヒストグラムが最大になる標高Ｈmodeを、領域の概略標高として与える。
【００６１】
次に、平面ステレオマッチングによる標高値の計算を行う（Ｓ２４）。この計算は、領域形状を不定形窓として、特願２０００−２５１４５６の不定形窓法により画像領域の最適標高を計算する（図７参照）。
【００６２】
そして、計算結果の標高をＨopt とする。なお、Ｌ（ｉ，ｊ）にＬＬ（ｉ，ｊ）を読み込んだ場合は、画像ＩLをＩ1、ＩRをＩ2とする。また、ＬＲ（ｉ，ｊ）を読み込んだ場合は、画像ＩRをＩ1、ＩLをＩ2とする。
【００６３】
そして、標高値が次の条件を満たすかどうかを検定する（Ｓ２５）。
【００６４】
｜Ｈmode −Ｈopt ｜＜ΔＨ
ここでΔＨは標高検定のしきい値である。
【００６５】
次に、標高値の検定で条件を満たす場合は、座標計算モジュール１０を用いて画像領域に対応する地上座標を求め、処理画像が左画像ならＤＳＭ配列ｈＬ（ｉｘ，ｉｙ）に標高値Ｈopt を格納する。右画像ならＤＳＭは配列ｈＲ（ｉｘ，ｉｙ））に格納する。すべての領域について処理済みでなければ、Ｓ２３に戻る（Ｓ２７）。すべての領域について処理が終わった場合は、Ｓ２２に戻って次の画像を読み込み、同様の処理を行う（Ｓ２８）。
【００６６】
次に、ＤＳＭ配列ｈＬ（ｘ，ｙ）とｈＲ（ｘ，ｙ）を融合し、ｈＬＲ（ｘ，ｙ）を生成するＤＳＭの融合処理を行う（Ｓ２９：図８参照）。
【００６７】
【数１】
このＤＳＭ配列を外部記憶装置に格納する（Ｓ３０）。
【００６８】
（ＤＳＭの領域分割）
ＤＳＭを水平な領域に分割する。ＤＳＭの領域分割のフローを図９に示す。
【００６９】
ＤＳＭの領域分割は、図９に示すように、ＤＳＭｈＬＲ（ｘ，ｙ）をＤＳＭ配列ｈ（ｘ，ｙ）に読み込む（Ｓ３１）。
【００７０】
そして、この読み込んだＤＳＭ配列ｈ（ｘ，ｙ）、水平領域ごとに領域分割する（Ｓ３２）。
【００７１】
具体的には、高さが一定の範囲Ｄｈの間にある連続領域をひとつの領域とする。そして、各領域に識別番号をラベルとして付け、ラベル格納配列ＬＨ（ｘ，ｙ）に格納する。すなわち、ＬＨ（ｘ，ｙ）＝ｎならば、ＤＳＭｈ（ｘ，ｙ）が、識別番号ｎの領域に属することを意味する。識別番号がｎに対応する領域をＲＬＨ（ｎ）とする。
【００７２】
次に、このラベル格納配列ＬＨ（ｘ，ｙ）を外部記憶装置に格納する（Ｓ３３）。
【００７３】
（水平領域の多角形化）
次に、各水平領域の境界を直線によって多角形化する。多角形化にあたっては、ステレオ画像を援用する。多角形化の概略フローを図１０に示す。
【００７４】
初めに、ＤＳＭとその領域分割結果を読み込む（Ｓ４０）。具体的には、ラベルＬＨ（ｘ，ｙ）とそれに対応するＤＳＭｈ（ｘ，ｙ）を外部記憶装置から読み込む。
【００７５】
次に、水平領域ごとの処理として、ＬＨ（ｘ，ｙ）に格納されたラベル番号ｋ＝1〜ｎ（ｎはラベルの個数）の各々の領域について、以下の処理を行う。
【００７６】
まず、地盤高計算モジュール１１を用いて、各領域ＲＬＨ（ｋ）の位置での地盤高ＨＧを計算する。このとき、領域の標高Ｈ（ｋ）との差が閾値（２ｍ程度）より大きいかどうかを以下の式で判定する地表面からの高さの判定処理を行う（Ｓ４１）。
【００７７】
Ｈ（ｋ）−ＨＧ＞ΔＨＢ
そして、領域の標高Ｈ（ｋ）との差が閾値より、大きい場合は、この領域が建物であるとみなし、以下の処理を行う。
【００７８】
座標計算モジュール１０を用いてラベル領域ＲＬＨ（ｋ）に対応する画像領域を求め、その周辺を含めて切り出す「対応する画像領域の切り出し」処理を行う（Ｓ４２）。
【００７９】
このとき、図１１に示すように、切り出し原点を（ｘｃ０，ｙｃ０）、切り出した画像をＩｃ、切り出しサイズをｓｘｃ，ｓｙｃとする。また、Ｉｃと同じサイズをもつＤＳＭのマスク画像Ｉｃｍとする。
【００８０】
このＩｃｍのマスク画像は、Ｉｃ上でＲＬＨ（ｋ）に対応する領域でコードＴ（≠０）、それ以外で０となる２値画像である。
【００８１】
次に、画像Ｉｃ上で、Ｉｃｍで指定されるマスクの周辺部で直線抽出を行うために、線抽出のためのバッファＩｃｂを生成する線抽出域のバッファリング処理を行う（Ｓ４３：図１２参照）。
【００８２】
具体的には、Ｉｃｂは、図１２に示すように、ＩｃｍをＮｆａｔ回膨張させた膨張マスクＩｃｍｆａｔと、ＩｃｍｆａｔをＮｔｈｉｎ回収縮された収縮マスクＩｃｍｔｈｉｎとを生成し、この２つのマスクから次のような関係式を成立させる。
【００８３】
収縮マスクＩｃｍｔｈｉｎは、もとのマスクＩｃｍと同じかそれよりも少なくする。よって、Ｎｆａｔ＜Ｎｔｈｉｎとする。
【００８４】
次に、図１３に示すように、切り出し画像Ｉｃ上でエッジ抽出を行う（Ｓ４４）。エッジ抽出方法は例えば、次の仕様を満たすものとする。
【００８５】
・エッジ上の点とそれ以外で２値化される。
【００８６】
・エッジが細線化されている（エッジの幅が１ピクセルである）。
【００８７】
このようなエッジ抽出方法として、ＣａｎｎｙオペレータやＳｕｓａｎｓオペレータなどがある。また、Ｓｏｂｅｌ等のエッジ強調オペレータを施し、その後２値化・細線化を行っても良い。
【００８８】
ここで得たエッジをＩｃｅとする。これに線抽出用バッファのマスク画像Ｉｃｂの共通領域をとることによって、マスクの周辺部のみのエッジ画像をＩｃｅから得る。このエッジ画像をＩｃｅｂとする（図１４）。
【００８９】
次に、マスク周辺のエッジ抽出画像Ｉｃｅｂのエッジ上の点を使ってＨｏｕｇｈ変換の原理で直線の抽出を行う直線化処理を行う（Ｓ４５：図１４参照）。
【００９０】
Ｈｏｕｇｈ変換では、直線の方程式を法線方向θと原点からの距離ρを使ってｘcosθ＋ｙsinθ＝ρとおき、（θ−ρ）空間において直線を抽出する方法である。Ｈｏｕｇｈ変換は画像から直線や円等の図形を抽出するための一手法である。
【００９１】
また、抽出する直線の数ｎhoughとして、最もたくさんのエッジ上の点を通る直線のパラメータの対（θ1 ，ρ1 ）からｎhough番目にエッジの点を通る直線の対（θｎ houghρｎ hough）までを抽出する（図１５参照）。
【００９２】
次に、Ｈｏｕｇｈ変換で抽出されたすべての直線について交点を求め、線分化する直線の線分化処理を行う（Ｓ４６、図１４の線分化）。
【００９３】
そして、収縮マスクＩｃｍｔｈｉｎと重なる線分を削除するマスクとの重複部分削除処理を行う（図１４参照）。
【００９４】
次に、線分を領域毎に多角形化（点列化）し、更に、その多角形の中で、マスク画像Ｉｃｍと最も大きな共通領域を持つ多角形を、ラベル領域ＲＬＨ（ｋ）に対応する多角形とする（Ｓ４８）。この多角形は、画像座標列ＰｌＩ（ｋ）とする（図１６参照）。
【００９５】
次に、画像座標で与えられたＰｌＩ（ｋ）を、３次元空間の多角形ＰｌＴ（ｋ）に変換する３次元点列化処理を行う（Ｓ４９）。このとき、３次元空間での多角形のＺ座標は、領域の標高Ｈ（ｋ）とする（図１７参照）。
【００９６】
そして、得られた多角形座標列ＰｌＩ（ｋ）を外部記憶装置に保存した後に（Ｓ５０）、次の領域に移る（Ｓ５１）。
【００９７】
（３次元モデル化処理）
前述の３次元モデル化処理について説明を補充する。各多角形を角柱状の三次元モデルに変換する概略フローを図１８に示す。
【００９８】
３次元モデル化処理においては、図１６に示す多角形座標列ＰｌＩ（ｋ）を外部記憶装置から読み出し（Ｓ６０）、多角形の地表面位置の計算を行う（Ｓ６１）。
【００９９】
この計算は、地盤高計算モジュールにより、多角形座標列ＰｌＴ（ｋ）を構成する各点の水平座標位置において、地表面位置の高さを求め、ＰｌＴ（ｋ）に対応する底面多角形ＰｌＥ（ｋ）を計算する（図１９参照）。
【０１００】
次に、多角形ＰｌＴ（ｋ）と底面多角形ＰｌＥ（ｋ）の対応する辺と、これらの辺の両端を鉛直線で結んだ四角形を側面として、角柱３次元モデルＢＬＤを建物モデルとする多角形の建物モデル化処理を行う（Ｓ６２）。
【０１０１】
この角柱３次元モデルＢＬＤは、角柱の頂上多角形ＰｌＴ（ｋ）と底面多角形ＰｌＥ（ｋ）、および側面多角形群ＰｌＷi （ｉ＝１,...ｎPlT(k)）から構成される（図１９参照）。ここでｎPlT(k)は、頂上多角形ＰｌＴ（ｋ）の頂点数である。
【０１０２】
そして、得られた角柱３次元ＢＬＤを外部記憶装置に格納する（Ｓ６３）。格納にあたっては、ＶＲＭＬなどの３次元モデル表記用フォーマットを用いるのが好ましい。
【０１０３】
従って、航空写真内の建物形状の精密計算や、３次元モデル生成に利用可能である。
【０１０４】
なお、上記実施の形態では、水平面ステレオマッチングは不定形窓を用いた水平面ステレオマッチングとして説明したが特願２００２−１２７５１２に記載の平面ステレオマッチング方式を用いてもよい。
【０１０５】
この平面ステレオマッチングは、ステレオ写真に写る平面が中心投影の原理に基づき射影変換により、対応付けられることを利用し、その平面上にある３点（どの３点も一直線上にない）の対応関係についてエピポーラ条件を満たす拘束条件を設定することによって、その面の上での射影変換係数を自動的に決定する。
【発明の効果】
以上のように本発明によれば、データベースに保存されている所定地域の空中写真画像から連続して同じ色彩の画素の色領域を抽出し、該抽出後に、データベースに保存されている前記所定地域のレーザスキャナデータに基づいて得られたデジタル地表面モデルに対して前記色領域を用いたステレオマッチングを行いＤＳＭを生成する。
【０１０６】
そして、ＤＳＭ内の水平領域と画像を用いて、水平領域周辺部に対応する画像上の線構造情報から直線化を行い、この直線に囲まれた領域の画素列を３次元化することによって前記水平領域を３次元点列化する。
【０１０７】
これにより、レーザスキャナデータでは曖昧であった地物境界をより高精度に表現することができるという効果が得られている。
【図面の簡単な説明】
【図１】レーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成装置の概略構成図である。
【図２】不定形窓を用いた標高値の付与を説明する説明図である。
【図３】ステレオ画像の色彩による領域分割を説明するフローチャートである。
【図４】画像の領域分割を説明する説明図である。
【図５】水平面ステレオマッチングを説明するフローチャートである。
【図６】レーザスキャナーデータの重ね合わせを説明する説明図である。
【図７】不定形窓による水平ステレオマッチングを説明する説明図である。
【図８】ＤＳＭの合成を説明する説明図である。
【図９】ＤＳＭの領域分割を説明するフローチャートである。
【図１０】水平領域の多角形化を説明するフローチャートである。
【図１１】画像領域の切り出しを説明する説明図である。
【図１２】線抽出域のバッファリングを説明する説明図である。
【図１３】エッジ抽出を説明する説明図である。
【図１４】直線化および直線の線分化の説明図である。
【図１５】Ｈｏｕｇｈ変換の説明図である。
【図１６】多角形化の説明図である。
【図１７】３次元点列化を説明する説明図である。
【図１８】３次元モデル化を説明するフローチャートである。
【図１９】本実施の形態による建物の３次元モデル化を説明する説明図である。
【図２０】ステレオ空中写真を説明する説明図である。
【図２１】レーザスキャナーデータを説明する説明図である。
【図２２】デジタル地形モデルを説明する説明図である。
【符号の説明】
１０ 標計算モジュール
１１ 地盤高計算モジュール
１５ 領域分割処理部
１６ 水平面ステレオマッチング処理部
１７ 領域分割処理部
１８ 多角形生成処理部
１９ ３次元モデル化処理部
【発明の属する技術分野】
レーザスキャナデータと空中写真画像を用いて、空中写真画像の建物面を容易に精度よく３次元化できる高精度都市モデルの生成装置に関する。
【０００２】
【従来の技術】
高精度の都市モデルを生成する手法は、一般にレーザスキャナデータを用いる手法とステレオマッチングによる手法とがある。
【０００３】
例えば、レーザスキャナデータは、飛行機又はヘリコプター等の航空機から地表にレーザ光をスキャニングし、これらの受光データを解析することで、図２１に示すようなレーザデータを得る。この図２１は、図２１（ｂ）に対象エリア（航空写真で示している）を、図２１（ａ）にそのレーザデータを示している。図２０（ａ）に示すように、レーザデータは細かな点列になっており、それぞれが３次元座標（ｘ、ｙ、ｚ）を有している。
【０００４】
従って、レーザスキャナデータは、３次元座標を有するので航空写真のステレオマッチングなどに比べて、図２２に示す標高精度の高いＤＳＭ（Ｄｉｇｉｔａｌ Ｓｕｒｆａｃｅ Ｍｏｄｅｌの略で、建物などの地物を含む地表面の形状を表現するモデル）をとることができる。
【０００５】
一方、航空写真により自動標高計測を行うステレオマッチング手法（図２０を参照）もある。このステレオマッチングは、画像間の各画素の対応関係を求め、三角測量の原理で各対応点位置での３次元座標を求めるもので、高密度に標高データを計測することができる。
【０００６】
【発明が解決しようとする課題】
しかしながら、レーザスキャナーデータによる手法は標高に関する精度の高いデータが得られるが、３次元点群から都市モデルのような直線構造を多く持つデータを生成する場合には、境界部分を正確に取るために、データをさらに高密度化する必要がある。
【０００７】
つまり、レーザデータは細かな点列（３次元点群）であり、これ以上の解像度を得ることは容易ではない。
【０００８】
また、レーザスキャナデータを高密度化すると、データ取得にかかるコスト及びデータ容量が莫大になり、コスト面とデータハンドリングの両面で問題がある。
【０００９】
一方、ステレオマッチング手法に用いる画像は高分解能であり、人工構造物の境界や領域の識別に有効である。しかしながら、ステレオマッチングは隠れなどの影響により標高変化の激しい箇所での標高計測精度が低く、境界を直接３次元的に計測するのに問題が生じる。
【００１０】
本発明は、レーザスキャナデータの標高精度の利点と、航空写真の高分解能であることの利点、および既往の疎密度ＤＥＭ（Ｄｉｇｉｔａｌ Ｅｌｅｖａｔｉｏｎ Ｍｏｄｅｌの略で、地盤高を表現するモデル。概略地盤高がわかるので地盤高として利用できる）を組み合わせ、精度の高い都市３次元モデルを構築することを目的とする。
【００１１】
【課題を解決するための手段】
本発明のレーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法は、
所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶手段に記憶する工程と、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理工程と、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング工程と、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理工程と、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す工程と、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する工程と、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する工程と、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する工程と
を行うことを要旨とする。
【００１２】
本発明の高精度都市モデルの生成システムは、所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶した記憶手段とを備え、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理手段と、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング手段と、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理手段と、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す手段と、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する手段と、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する手段と、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する手段とを備えたことを要旨とする。
【００１３】
本発明の高精度都市モデルの生成のプログラムは、所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶した記憶手段とを備え、
コンピュータに、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理手段、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング手段、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理手段、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す手段、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する手段、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する手段、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する手段、
としての機能を実行させる。
【００１４】
【発明の実施の形態】
本実施の形態は、レーザスキャナデータにステレオ撮影された空中写真画像データから得られる領域情報・線構造情報を援用することによって、人工構造物の境界を精度よく抽出し、自動建物多角形化および３次元モデル化を行う。
【００１５】
まず、ステレオ空中写真画像の各々を領域分割する事により、色彩が均一な画像領域を得る。このような画像領域は、建物の屋上及び路面・グランド・その他広い平坦地に対応していると考えられる。
【００１６】
各色領域に重なるレーザスキャナの３次元点群の標高値を初期値として、ステレオマッチングを行う。その結果得られる標高値の計算結果が初期値とほぼ同じ高さになった場合はその領域に計算結果の標高値を与え、ＤＳＭ（Ｄｉｇｉｔａｌ Ｓｕｒｆａｃｅ Ｍｏｄｅｌの略で、建物等の地物を含む地表面の形状を表現するモデル）を生成する。これにより、レーザスキャナデータでは曖昧であった地物境界をより高精度に表現することができるようにしている。
【００１７】
続いて上記で得られたＤＳＭの高い方から順に水平な領域を領域分割によって認識し、各領域に対応するステレオ画像上の周辺領域で直線抽出を行い、その結果を多角形化する。既往のＤＥＭや地域の平均的地盤高に対するＤＳＭの差が任意に設定する閾値以下にある領域は、建物以外の領域であり、地盤高と推定する。これによって、建物領域を３次元多角形化することができる。更に多角形化された建物領域を角柱状に構成する事によって、３次元モデルを得る事が可能となるものである。
【００１８】
本発明を実現する３次元計測装置の構成を図１に示す。３次元計測システム１は、パーソナルコンピュータ本体２、ステレオ画像等を格納する外部記憶装置３、ディスプレイ４という一般的なパーソナルコンピュータからなる。また、次のことが前提として与えられている。
【００１９】
（１）入力データ
次の各種データが、外部記憶装置３に格納されているとする。
【００２０】
（ａ）ステレオ空中写真画像（ＩL およびＩR （図２０））とその標定要素と、
（ｂ）レーザスキャナデータが保存されている。このレーザスキャナデータは３次元座標の集まりであり、以下ではＬＰi （i=1...nLP）で表す（図２１（ａ））。
【００２１】
（ｃ）デジタル地形モデル（ＤＥＭ））。地盤高を表し、国土数値情報などの所与のデータを用いてよい（図２２）。
【００２２】
（２）計算モジュール
次の計算モジュールが与えられているものとする。この計算モジュールは、座標計算モジュール１０と、地盤高計算モジュール１１等からなる。
【００２３】
座標計算モジュール１０は、ステレオ空中写真画像の標定要素を読み取り、３次元座標と、ステレオ画像Ｉ1 およびＩ2 の座標との間の変換を行う。これは特願２０００−２５１４５６で座標変換関数・逆座標変換関数として与えられているものと同じとする。
【００２４】
つまり、輝度画像ＩＬ、ＩＲについて、地上の３次元座標が投影される各画像の画像座標算出のための座標変換関数ＦＸ、逆座標変換関数関数ＧＸ（総称して標定情報ともいう）を用いる。
【００２５】
例えば、３次元空間内の点（ｘ，ｙ，ｚ）が投影される左右画像のピクセル座標（ｘ１ｐｉｘ，ｙ１ｐｉｘ）および（ｘ２ｐｉｘ，ｙ２ｐｉｘ）を計算する。
【００２６】
すなわち、次の４つの座標変換関数ＦＸ１（ｘ，ｙ，ｚ）、ＦＹ１（ｘ，ｙ，ｚ）、ＦＸ２（ｘ，ｙ，ｚ）、ＦＹ２（ｘ，ｙ，ｚ）を解析的、もしくは手続き的に決定する。
【００２７】
ｘ１ｐｉｘ＝ＦＸ１（ｘ，ｙ，ｚ） 式（１）
ｙ１ｐｉｘ＝ＦＹ１（ｘ，ｙ，ｚ） 式（２）
ｘ２ｐｉｘ＝ＦＸ２（ｘ，ｙ，ｚ） 式（３）
ｙ２ｐｉｘ＝ＦＹ２（ｘ，ｙ，ｚ） 式（４）
また、ｚ＝一定として式（１），式（２）を解析的又は手続的に解くと、各画像上の点がある標高ｚにあると仮定した時の地上でのｘ座標とｙ座標を求めることができる。
【００２８】
ｘ＝ＧＸ１（ｘ１ｐｉｘ，ｙ１ｐｉｘ，ｚ） 式（５）
ｙ＝ＧＹ１（ｘ１ｐｉｘ，ｙ１ｐｉｘ，ｚ） 式（６）
同様に式（３），式（４）について解くと、
ｘ＝ＧＸ２（ｘ２ｐｉｘ，ｙ２ｐｉｘ，ｚ） 式（７）
ｙ＝ＧＹ２（ｘ２ｐｉｘ，ｙ２ｐｉｘ，ｚ） 式（８）
である。
【００２９】
このような関数は、例えば図２に示す不定形窓を用いているときに使用する。この図２は、航空機によって得た画像Ｉ１を画面に表示し、この画像上でマウス操作等で地物Ｂｕの輪郭に相当する閉領域Ｂ１を与える。そして、この閉領域に投影が一致する位置に不定形窓Ｐｌｈを発生させ、標高Ｈ（ｚ）を変化させる。
【００３０】
例えば３０ｍ、８０ｍ、１５０ｍの順に標高を変化させると、３０ｍ、８０ｍ、１５０ｍと変化したときの輝度画像Ｉ１、Ｉ２における不定形窓Ｐｌｈの輝度画像Ｉ１、Ｉ２の投影位置をその都度求めながら、この位置における不定形窓Ｐｌｈを投影したときの輝度画像Ｉ１、Ｉ２における輝度値を比較し、評価関数によって最も一致（近似）したと判定されるときの不定形窓Ｐ１ｈの高さ（標高値）としている。
【００３１】
地盤高計算モジュール１１は、ＤＥＭを読み取り、適当な内挿計算を用いることによって、与えられた水平座標位置における地盤高を計算する。
【００３２】
Ｈ＝ＧＥlv（ｘ，ｙ）
例えばＤＥＭが正方メッシュデータならば、バイリニア法等の内挿方法、ＤＥＭがランダムな配置の３次元点群ならば、三角メッシュを用いた内挿法を用いる。
【００３３】
すなわち、本実施の形態では、座標計算モジュール１０、地盤高計算モジュール１１等によって、レーザスキャナデータにステレオ撮影された空中写真画像データから得られる領域情報・線構造情報を援用して、例えば、人工構造物の境界を精度よく抽出し、自動建物多角形化および３次元モデル化を行う。
【００３４】
この３次元化モデルの生成においては、まず、ステレオ空中写真画像の各々を領域分割する事により、色彩が均一な画像領域を得る。
【００３５】
このような画像領域は、建物の屋上及び路面・グランド・その他広い平坦地に対応していると考えられる。
【００３６】
各色領域に重なるレーザスキャナの３次元点群の標高値を初期値として、ステレオマッチングを行い、その結果得られる標高値の計算結果が初期値とほぼ同じ高さになった場合はその領域に計算結果の標高値を与え、ＤＳＭを生成する。これにより、レーザスキャナデータでは曖昧であった地物境界をより高精度に表現することができる。
【００３７】
続いて上記で得られたＤＳＭの水平な領域を領域分割によって認識し、各領域に対応するステレオ画像上の周辺領域で直線抽出を行い、その結果を多角形化する。
【００３８】
既往のＤＥＭや地域の平均的地盤高に対するＤＳＭの差が任意に設定する閾値以下にある領域は、建物以外の領域であり、地盤高と推定する。これによって、建物領域を３次元多角形化することができる。更に多角形化された建物領域を角柱状に構成する事によって、３次元モデルを得る事ができる。
【００３９】
本発明の実施形態の概略を図１を用いて簡単に説明する。図１に示すように、次の５つの部分から成り立っている。
【００４０】
初めに、データベース（図示せず）に保存されているステレオ空中写真画像（ＩＲ、ＩＬ）を読み込んで、ステレオ空中写真画像の各々を領域分割（同じ色を一塊りとする領域：同じ色彩の画素が連続する領域）する事により、色彩が画像領域（建物の屋上及び路面・グランド・その他広い平坦地に対応していると考えられる）を得るステレオ画像の色彩による領域分割処理部１５と、各分割領域（分割色領域ともいう）に重なるレーザスキャナデータの３次元点群の標高値を初期値として、ステレオマッチングを行う水平面ステレオマッチング処理部１６と、この水平面ステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになった場合はその領域に計算結果の標高値を与え、ＤＳＭ（Ｄｉｇｉｔａｌ Ｓｕｒｆａｃｅ Ｍｏｄｅｌの略で、建物等の地物を含む地表面の形状を表現するモデル）を生成するＤＳＭの領域分割処理部１７と、多角形生成処理部１８と、３次元モデル化処理部１９等とからなる。
【００４１】
以下にこれらの各処理（ステレオ画像の色彩による領域分割処理、水平面ステレオマッチング処理、ＤＳＭの分割処理、多角形生成処理、３次元モデル化処理）の詳細をフローチャートを用いて説明する。
【００４２】
（ステレオ画像の色彩による領域分割処理部１５）
ステレオ画像の色彩による領域分割について図３のフローチャートを用いて詳細に説明する。
【００４３】
初めに、領域分割に必要なパラメータを設定する（Ｓ１０）。例えば、反復型領域拡張法によって領域分割を行う場合は、次の３つのパラメータがユーザにより指定されているものとする。
【００４４】
・初期しきい値
・しきい値の増分
・反復打ち切りのためのしきい値
次に、ステレオ空中写真画像の内、左画像ＩLを外部記憶装置３から内部に読み込む（Ｓ１１）。
【００４５】
次に、画像の輝度情報や色情報を用いて領域分割処理を行う（図４参照）。この領域分割処理では、まずメディアンフィルタなどで小さな孤立点を除去した後、除去画像の輝度や色が均一な連続領域ごとに識別番号を付ける。
【００４６】
なお、領域分割の方法としては、反復型領域拡張法など一般に広く行われている領域分割方法を用いてもよい。
【００４７】
そして、領域分割の結果（単色で塗られた領域）、各画素が何番目の識別番号の領域に属するかがわかるラベル番号が与えられるものとする。この左画像のラベル番号を格納する配列をＬＬ（ｉ，ｊ）とする。
【００４８】
すなわち、ＬＬ（ｉ，ｊ）＝ｎならば、画像ｉカラムｊラインの画素が、ラベル番号ｎの領域に属することを意味する。
【００４９】
次に、領域分割によって求められたラベル格納配列ＬＬ（ｉ，ｊ）を左画像のラベルとして外部記憶装置３に格納する領域分割結果保存処理を行う（Ｓ１３）。
【００５０】
そして、左画像読込（Ｓ１１）〜領域分割結果保存処理（Ｓ１３）と同様の処理を右画像ＩR について同様に行う右画像読込処理（Ｓ１４）と領域分割処理（Ｓ１５）とを行い、右画像のラベル格納配列ＬＲ（ｉ，ｊ）を外部記憶装置３に格納する（Ｓ１６）。
【００５１】
つまり、Ｓ１４〜Ｓ１６は右画像読み込み・領域分割・領域分割結果保存を行っている。
【００５２】
次に、水平面ステレオマッチング処理について説明する。
【００５３】
（水平面ステレオマッチング）
領域分割された各画像領域が水平であると仮定し、その拘束の元で領域分割された部分を「不定形窓を用いた画像間拡張イメージマッチング方法（特願２０００−２５１４５６）」を行い、デジタル地表面モデル（ＤＳＭ）を生成する。ここでのＤＳＭとは、地表面の標高をメッシュ点上で表現するものである。例えば、メッシュの原点を（ｘ０，ｙ０）、メッシュ間隔をｄとすると、ＤＳＭ配列ｈ（ｉｘ，ｉｙ）は、座標（ｘ０＋ｄ・ｉｘ、ｙ０＋ｄ・ｉｙ）での標高を与えるものとする。また、ＤＳＭ配列の原点や大きさは、ユーザが与える対象範囲にあわせて設定するものとする。前述の水平面ステレオマッチングは、各分割領域を不定形窓として、これを変化させて拡張ステレオモデルの両方の画像上の投影位置を求め、この投影位置における不定形窓で区切られた領域の輝度値をそれぞれ求めて輝度値を比較し、比較結果が最も一致するときの不定形窓の高さを拡張ステレオモデルにおける地物の高さとする拡張ステレオマッチングによる標高付与処理を行う。
【００５４】
図５は平面ステレオマッチングを説明するフローチャートである。初めに、左右画像をメモリ上に読み込む（Ｓ２０）。
【００５５】
次に、ＤＳＭ配列を格納するＤＳＭ配列を初期化する（Ｓ２１）。初期化するＤＳＭ配列は次の３つとする。
【００５６】
・ｈＬ（ｉｘ，ｉｙ）：左画像領域分割結果から生成されるＤＳＭ
・ｈＲ（ｉｘ，ｉｙ）：右画像領域分割結果から生成されるＤＳＭ
・ｈＬＲ（ｉｘ，ｉｙ）：左右画像領域分割結果から生成されるＤＳＭ
なお、ＤＳＭ配列の原点や範囲は、ユーザがあらかじめ与えるものとする。初期化にあたっては、すべての配列の値を一定値ｈｉｎｉｔとしておく。ｈｉｎｉｔは、実際にはありえない小さな値（−９９９９など）であり、ＤＳＭ配列がｈｉｎｉｔであるということは、平面ステレオマッチングによってＤＳＭが計算されていないことを示す。
【００５７】
次に、左右画像の画像領域の分割結果を読込む（Ｓ２２）。つまり、画像領域分割結果をメモリ上の配列Ｌ（ｉ，ｊ）に読み込む。
【００５８】
そして、各領域毎Ｌ（ｉ，ｊ）に記録されているラベル番号が付けられた各々の領域について、番号順に次の処理を行う。
【００５９】
この処理は、レーザスキャナデータの重ね合わせと概略標高の計算を行う（Ｓ２３）。
【００６０】
例えば、レーザスキャナデータＬpiの画像への投影位置を座標計算モジュールによって計算し、その中で処理中の領域に重なるものを検索する（図６を参照）。 次に、更に重なった点のヒストグラムを作り、ヒストグラムが最大になる標高Ｈmodeを、領域の概略標高として与える。
【００６１】
次に、平面ステレオマッチングによる標高値の計算を行う（Ｓ２４）。この計算は、領域形状を不定形窓として、特願２０００−２５１４５６の不定形窓法により画像領域の最適標高を計算する（図７参照）。
【００６２】
そして、計算結果の標高をＨopt とする。なお、Ｌ（ｉ，ｊ）にＬＬ（ｉ，ｊ）を読み込んだ場合は、画像ＩLをＩ1、ＩRをＩ2とする。また、ＬＲ（ｉ，ｊ）を読み込んだ場合は、画像ＩRをＩ1、ＩLをＩ2とする。
【００６３】
そして、標高値が次の条件を満たすかどうかを検定する（Ｓ２５）。
【００６４】
｜Ｈmode −Ｈopt ｜＜ΔＨ
ここでΔＨは標高検定のしきい値である。
【００６５】
次に、標高値の検定で条件を満たす場合は、座標計算モジュール１０を用いて画像領域に対応する地上座標を求め、処理画像が左画像ならＤＳＭ配列ｈＬ（ｉｘ，ｉｙ）に標高値Ｈopt を格納する。右画像ならＤＳＭは配列ｈＲ（ｉｘ，ｉｙ））に格納する。すべての領域について処理済みでなければ、Ｓ２３に戻る（Ｓ２７）。すべての領域について処理が終わった場合は、Ｓ２２に戻って次の画像を読み込み、同様の処理を行う（Ｓ２８）。
【００６６】
次に、ＤＳＭ配列ｈＬ（ｘ，ｙ）とｈＲ（ｘ，ｙ）を融合し、ｈＬＲ（ｘ，ｙ）を生成するＤＳＭの融合処理を行う（Ｓ２９：図８参照）。
【００６７】
【数１】
【００６８】
（ＤＳＭの領域分割）
ＤＳＭを水平な領域に分割する。ＤＳＭの領域分割のフローを図９に示す。
【００６９】
ＤＳＭの領域分割は、図９に示すように、ＤＳＭｈＬＲ（ｘ，ｙ）をＤＳＭ配列ｈ（ｘ，ｙ）に読み込む（Ｓ３１）。
【００７０】
そして、この読み込んだＤＳＭ配列ｈ（ｘ，ｙ）、水平領域ごとに領域分割する（Ｓ３２）。
【００７１】
具体的には、高さが一定の範囲Ｄｈの間にある連続領域をひとつの領域とする。そして、各領域に識別番号をラベルとして付け、ラベル格納配列ＬＨ（ｘ，ｙ）に格納する。すなわち、ＬＨ（ｘ，ｙ）＝ｎならば、ＤＳＭｈ（ｘ，ｙ）が、識別番号ｎの領域に属することを意味する。識別番号がｎに対応する領域をＲＬＨ（ｎ）とする。
【００７２】
次に、このラベル格納配列ＬＨ（ｘ，ｙ）を外部記憶装置に格納する（Ｓ３３）。
【００７３】
（水平領域の多角形化）
次に、各水平領域の境界を直線によって多角形化する。多角形化にあたっては、ステレオ画像を援用する。多角形化の概略フローを図１０に示す。
【００７４】
初めに、ＤＳＭとその領域分割結果を読み込む（Ｓ４０）。具体的には、ラベルＬＨ（ｘ，ｙ）とそれに対応するＤＳＭｈ（ｘ，ｙ）を外部記憶装置から読み込む。
【００７５】
次に、水平領域ごとの処理として、ＬＨ（ｘ，ｙ）に格納されたラベル番号ｋ＝1〜ｎ（ｎはラベルの個数）の各々の領域について、以下の処理を行う。
【００７６】
まず、地盤高計算モジュール１１を用いて、各領域ＲＬＨ（ｋ）の位置での地盤高ＨＧを計算する。このとき、領域の標高Ｈ（ｋ）との差が閾値（２ｍ程度）より大きいかどうかを以下の式で判定する地表面からの高さの判定処理を行う（Ｓ４１）。
【００７７】
Ｈ（ｋ）−ＨＧ＞ΔＨＢ
そして、領域の標高Ｈ（ｋ）との差が閾値より、大きい場合は、この領域が建物であるとみなし、以下の処理を行う。
【００７８】
座標計算モジュール１０を用いてラベル領域ＲＬＨ（ｋ）に対応する画像領域を求め、その周辺を含めて切り出す「対応する画像領域の切り出し」処理を行う（Ｓ４２）。
【００７９】
このとき、図１１に示すように、切り出し原点を（ｘｃ０，ｙｃ０）、切り出した画像をＩｃ、切り出しサイズをｓｘｃ，ｓｙｃとする。また、Ｉｃと同じサイズをもつＤＳＭのマスク画像Ｉｃｍとする。
【００８０】
このＩｃｍのマスク画像は、Ｉｃ上でＲＬＨ（ｋ）に対応する領域でコードＴ（≠０）、それ以外で０となる２値画像である。
【００８１】
次に、画像Ｉｃ上で、Ｉｃｍで指定されるマスクの周辺部で直線抽出を行うために、線抽出のためのバッファＩｃｂを生成する線抽出域のバッファリング処理を行う（Ｓ４３：図１２参照）。
【００８２】
具体的には、Ｉｃｂは、図１２に示すように、ＩｃｍをＮｆａｔ回膨張させた膨張マスクＩｃｍｆａｔと、ＩｃｍｆａｔをＮｔｈｉｎ回収縮された収縮マスクＩｃｍｔｈｉｎとを生成し、この２つのマスクから次のような関係式を成立させる。
【００８３】
【００８４】
次に、図１３に示すように、切り出し画像Ｉｃ上でエッジ抽出を行う（Ｓ４４）。エッジ抽出方法は例えば、次の仕様を満たすものとする。
【００８５】
・エッジ上の点とそれ以外で２値化される。
【００８６】
・エッジが細線化されている（エッジの幅が１ピクセルである）。
【００８７】
このようなエッジ抽出方法として、ＣａｎｎｙオペレータやＳｕｓａｎｓオペレータなどがある。また、Ｓｏｂｅｌ等のエッジ強調オペレータを施し、その後２値化・細線化を行っても良い。
【００８８】
ここで得たエッジをＩｃｅとする。これに線抽出用バッファのマスク画像Ｉｃｂの共通領域をとることによって、マスクの周辺部のみのエッジ画像をＩｃｅから得る。このエッジ画像をＩｃｅｂとする（図１４）。
【００８９】
次に、マスク周辺のエッジ抽出画像Ｉｃｅｂのエッジ上の点を使ってＨｏｕｇｈ変換の原理で直線の抽出を行う直線化処理を行う（Ｓ４５：図１４参照）。
【００９０】
Ｈｏｕｇｈ変換では、直線の方程式を法線方向θと原点からの距離ρを使ってｘcosθ＋ｙsinθ＝ρとおき、（θ−ρ）空間において直線を抽出する方法である。Ｈｏｕｇｈ変換は画像から直線や円等の図形を抽出するための一手法である。
【００９１】
また、抽出する直線の数ｎhoughとして、最もたくさんのエッジ上の点を通る直線のパラメータの対（θ1 ，ρ1 ）からｎhough番目にエッジの点を通る直線の対（θｎ houghρｎ hough）までを抽出する（図１５参照）。
【００９２】
次に、Ｈｏｕｇｈ変換で抽出されたすべての直線について交点を求め、線分化する直線の線分化処理を行う（Ｓ４６、図１４の線分化）。
【００９３】
そして、収縮マスクＩｃｍｔｈｉｎと重なる線分を削除するマスクとの重複部分削除処理を行う（図１４参照）。
【００９４】
次に、線分を領域毎に多角形化（点列化）し、更に、その多角形の中で、マスク画像Ｉｃｍと最も大きな共通領域を持つ多角形を、ラベル領域ＲＬＨ（ｋ）に対応する多角形とする（Ｓ４８）。この多角形は、画像座標列ＰｌＩ（ｋ）とする（図１６参照）。
【００９５】
次に、画像座標で与えられたＰｌＩ（ｋ）を、３次元空間の多角形ＰｌＴ（ｋ）に変換する３次元点列化処理を行う（Ｓ４９）。このとき、３次元空間での多角形のＺ座標は、領域の標高Ｈ（ｋ）とする（図１７参照）。
【００９６】
そして、得られた多角形座標列ＰｌＩ（ｋ）を外部記憶装置に保存した後に（Ｓ５０）、次の領域に移る（Ｓ５１）。
【００９７】
（３次元モデル化処理）
前述の３次元モデル化処理について説明を補充する。各多角形を角柱状の三次元モデルに変換する概略フローを図１８に示す。
【００９８】
３次元モデル化処理においては、図１６に示す多角形座標列ＰｌＩ（ｋ）を外部記憶装置から読み出し（Ｓ６０）、多角形の地表面位置の計算を行う（Ｓ６１）。
【００９９】
この計算は、地盤高計算モジュールにより、多角形座標列ＰｌＴ（ｋ）を構成する各点の水平座標位置において、地表面位置の高さを求め、ＰｌＴ（ｋ）に対応する底面多角形ＰｌＥ（ｋ）を計算する（図１９参照）。
【０１００】
次に、多角形ＰｌＴ（ｋ）と底面多角形ＰｌＥ（ｋ）の対応する辺と、これらの辺の両端を鉛直線で結んだ四角形を側面として、角柱３次元モデルＢＬＤを建物モデルとする多角形の建物モデル化処理を行う（Ｓ６２）。
【０１０１】
この角柱３次元モデルＢＬＤは、角柱の頂上多角形ＰｌＴ（ｋ）と底面多角形ＰｌＥ（ｋ）、および側面多角形群ＰｌＷi （ｉ＝１,...ｎPlT(k)）から構成される（図１９参照）。ここでｎPlT(k)は、頂上多角形ＰｌＴ（ｋ）の頂点数である。
【０１０２】
そして、得られた角柱３次元ＢＬＤを外部記憶装置に格納する（Ｓ６３）。格納にあたっては、ＶＲＭＬなどの３次元モデル表記用フォーマットを用いるのが好ましい。
【０１０３】
従って、航空写真内の建物形状の精密計算や、３次元モデル生成に利用可能である。
【０１０４】
なお、上記実施の形態では、水平面ステレオマッチングは不定形窓を用いた水平面ステレオマッチングとして説明したが特願２００２−１２７５１２に記載の平面ステレオマッチング方式を用いてもよい。
【０１０５】
この平面ステレオマッチングは、ステレオ写真に写る平面が中心投影の原理に基づき射影変換により、対応付けられることを利用し、その平面上にある３点（どの３点も一直線上にない）の対応関係についてエピポーラ条件を満たす拘束条件を設定することによって、その面の上での射影変換係数を自動的に決定する。
【発明の効果】
以上のように本発明によれば、データベースに保存されている所定地域の空中写真画像から連続して同じ色彩の画素の色領域を抽出し、該抽出後に、データベースに保存されている前記所定地域のレーザスキャナデータに基づいて得られたデジタル地表面モデルに対して前記色領域を用いたステレオマッチングを行いＤＳＭを生成する。
【０１０６】
そして、ＤＳＭ内の水平領域と画像を用いて、水平領域周辺部に対応する画像上の線構造情報から直線化を行い、この直線に囲まれた領域の画素列を３次元化することによって前記水平領域を３次元点列化する。
【０１０７】
これにより、レーザスキャナデータでは曖昧であった地物境界をより高精度に表現することができるという効果が得られている。
【図面の簡単な説明】
【図１】レーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成装置の概略構成図である。
【図２】不定形窓を用いた標高値の付与を説明する説明図である。
【図３】ステレオ画像の色彩による領域分割を説明するフローチャートである。
【図４】画像の領域分割を説明する説明図である。
【図５】水平面ステレオマッチングを説明するフローチャートである。
【図６】レーザスキャナーデータの重ね合わせを説明する説明図である。
【図７】不定形窓による水平ステレオマッチングを説明する説明図である。
【図８】ＤＳＭの合成を説明する説明図である。
【図９】ＤＳＭの領域分割を説明するフローチャートである。
【図１０】水平領域の多角形化を説明するフローチャートである。
【図１１】画像領域の切り出しを説明する説明図である。
【図１２】線抽出域のバッファリングを説明する説明図である。
【図１３】エッジ抽出を説明する説明図である。
【図１４】直線化および直線の線分化の説明図である。
【図１５】Ｈｏｕｇｈ変換の説明図である。
【図１６】多角形化の説明図である。
【図１７】３次元点列化を説明する説明図である。
【図１８】３次元モデル化を説明するフローチャートである。
【図１９】本実施の形態による建物の３次元モデル化を説明する説明図である。
【図２０】ステレオ空中写真を説明する説明図である。
【図２１】レーザスキャナーデータを説明する説明図である。
【図２２】デジタル地形モデルを説明する説明図である。
【符号の説明】
１０ 標計算モジュール
１１ 地盤高計算モジュール
１５ 領域分割処理部
１６ 水平面ステレオマッチング処理部
１７ 領域分割処理部
１８ 多角形生成処理部
１９ ３次元モデル化処理部
Claims (11)
- 所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶手段に記憶する工程と、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理工程と、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング工程と、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理工程と、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す工程と、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する工程と、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する工程と、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する工程と
を行うことを特徴とするレーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法。 - 前記直線化を施した後に、これらの直線に囲まれた各領域を抽出して、この領域が前記マスク画像に対して最も大きな共通領域を有しているときに前記ＤＳＭ領域の多角形とすることを特徴とする請求項１記載のレーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法。
- 前記多角形を構成する各点の水平座標によって前記デジタル地形モデル（ＤＳＭ）における前記多角形の底面高さを求める工程と、
前記多角形の底面の辺と、これらの辺の両端を鉛直線で結んだ四角形を側面とした角柱モデルを構築する工程と
を行うことを特徴とする請求項１又は２記載のレーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法。 - 前記ステレオマッチングは、前記ステレオ空中写真画像の色彩情報を利用して均一な色領域を抽出し、各色領域が水平であると仮定して行うことを特徴とする請求項１、又は３記載のレーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法。
- 所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶した記憶手段とを備え、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理手段と、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング手段と、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高 の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理手段と、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す手段と、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する手段と、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する手段と、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する手段と
を有することを特徴とする高精度都市モデルの生成システム。 - 前記直線化を施した後に、これらの直線に囲まれた各領域を抽出して、この領域が前記マスク画像に対して最も大きな共通領域を有しているときに前記ＤＳＭ領域の多角形とすることを特徴とする請求項５記載の高精度都市モデルの生成システム。
- 前記多角形を構成する各点の水平座標によって前記デジタル地形モデル（ＤＳＭ）における前記多角形の底面高さを求める手段と、
前記多角形の底面の辺と、これらの辺の両端を鉛直線で結んだ四角形を側面とした角柱モデルを構築する手段と
を有することを特徴とする請求項６又は７記載の高精度都市モデルの生成方法。 - 前記ステレオマッチングは、前記ステレオ空中写真画像の色彩情報を利用して均一な色領域を抽出し、各色領域が水平であると仮定して行うことを特徴とする請求項５、６又は７記載の高精度都市モデルの生成システム。
- 所定地域のステレオ空中写真画像（ＩＲ、ＩＬ）と、前記所定地域のレーザースキャナーデータと、デジタル地形モデル（ＤＳＭ）とを記憶した記憶手段とを備え、
コンピュータに、
前記記憶手段からステレオ空中写真画像を読み込んで、該ステレオ空中写真画像の各々の同じ色彩が連続する領域を分割領域として順次、抽出する領域分割処理手段、
前記ステレオ空中写真画像の前記分割領域が抽出される毎に、該分割領域に重なる前記レーザスキャナデータの３次元点群の標高値を初期値とし、かつ前記分割領域を不定形窓とし、前記初期値を用いてこれを水平方向に変化させてステレオマッチングを行うステレオマッチング手段、
このステレオマッチング処理後の結果得られる標高値の計算結果が初期値とほぼ同じ高さになる前記デジタル地形モデル（ＤＳＭ）の領域に計算結果の標高値を与え、この標高の高さにある連続領域を一つのＤＳＭ領域とするＤＳＭの領域分割処理手段、
前記ＤＳＭ領域と前記分割領域とを読み込み、該ＤＳＭ領域が地盤高さ以上の場合は、前記分割領域の周辺を含めた領域を、切出し領域として切出す手段、
前記切出し領域と同じ大きさを持つ前記ＤＳＭ領域のマスク画像を生成し、このマスク画像を、膨張させた膨張マスク及び収縮させた収縮マスクを合成して重なり合う部分以外を線抽出領域とした線抽出用バッファを生成する手段、
前記切出し領域内の画像のエッジの抽出を行い、この抽出されたエッジ画像に対して前記線抽出用バッファをかけて前記線抽出領域のエッジ画像を抽出する手段、
前記抽出したエッジ画像上の点を使ってＨｏｕｇｈ変換によって直線化を施して、これらの前記収縮マスクと重なる部分の直線を削除して、前記マスク画像に対して最も大きな共通領域を前記ＤＳＭ領域の多角形として出力する手段、
としての機能を実行させるための高精度都市モデルの生成のプログラム。 - 前記直線化を施した後に、これらの直線に囲まれた各領域を抽出して、この領域が前記マスク画像に対して最も大きな共通領域を有しているときに前記ＤＳＭ領域の多角形とすることを特徴とする請求項９記載の高精度都市モデルの生成のプログラム。
- 前記コンピュータに、
前記多角形を構成する各点の水平座標によって前記デジタル地形モデル（ＤＳＭ）における前記多角形の底面高さを求める手段、
前記多角形の底面の辺と、これらの辺の両端を鉛直線で結んだ四角形を側面とした角柱モデルを構築する手段、
としての機能を実行させるための請求項９又は１０記載のの高精度都市モデルの生成のプログラム。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2002127672A JP4058293B2 (ja) | 2002-04-26 | 2002-04-26 | レーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法及び高精度都市モデルの生成システム並びに高精度都市モデルの生成のプログラム |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2002127672A JP4058293B2 (ja) | 2002-04-26 | 2002-04-26 | レーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法及び高精度都市モデルの生成システム並びに高精度都市モデルの生成のプログラム |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2003323640A JP2003323640A (ja) | 2003-11-14 |
| JP4058293B2 true JP4058293B2 (ja) | 2008-03-05 |
Family
ID=29541675
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2002127672A Expired - Lifetime JP4058293B2 (ja) | 2002-04-26 | 2002-04-26 | レーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法及び高精度都市モデルの生成システム並びに高精度都市モデルの生成のプログラム |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP4058293B2 (ja) |
Families Citing this family (20)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| JP3854270B2 (ja) * | 2003-12-11 | 2006-12-06 | アジア航測株式会社 | 浸水深補正方法 |
| KR100715460B1 (ko) * | 2004-04-06 | 2007-05-07 | 이엔지정보기술 주식회사 | 레이저 스캐닝 데이터를 이용한 정밀 수치표고모델 제작장치 및 그 제작 방법 |
| JP4752409B2 (ja) * | 2005-09-09 | 2011-08-17 | ソニー株式会社 | 画像処理装置および方法、並びにプログラム |
| JP2008226189A (ja) * | 2007-03-15 | 2008-09-25 | Aisin Seiki Co Ltd | 特徴点検出装置及びプログラム |
| WO2009003529A1 (en) * | 2007-07-04 | 2009-01-08 | Saab Ab | Arrangement and method for providing a three dimensional map representation of an area |
| JP4852006B2 (ja) * | 2007-07-27 | 2012-01-11 | 株式会社パスコ | 空間情報データベース生成装置及び空間情報データベース生成プログラム |
| KR100912715B1 (ko) * | 2007-12-17 | 2009-08-19 | 한국전자통신연구원 | 이종 센서 통합 모델링에 의한 수치 사진 측량 방법 및장치 |
| EP2278553A4 (en) * | 2008-04-23 | 2013-10-23 | Pasco Corp | DEVICE, METHOD AND PROGRAM FOR DETECTING BUILDING ACOUSTIC CONTOURS |
| JP2010151577A (ja) * | 2008-12-25 | 2010-07-08 | Shimizu Corp | 配筋検査装置および配筋検査方法 |
| JP5356269B2 (ja) * | 2010-01-29 | 2013-12-04 | 株式会社パスコ | レーザデータのフィルタリング方法及び装置 |
| JP5305485B2 (ja) * | 2011-05-23 | 2013-10-02 | Ｎｅｃシステムテクノロジー株式会社 | 地盤高データ生成装置、地盤高データ生成方法、及びプログラム |
| JP2013072705A (ja) * | 2011-09-27 | 2013-04-22 | Kumagai Gumi Co Ltd | 岩盤斜面における不連続面の検出方法とその検出装置 |
| JP5986362B2 (ja) * | 2011-09-27 | 2016-09-06 | 株式会社熊谷組 | トンネル切羽の不連続面の抽出方法及びその装置 |
| JP5817422B2 (ja) * | 2011-10-18 | 2015-11-18 | 朝日航洋株式会社 | 建物抽出装置、方法及びプログラム |
| JP6296444B2 (ja) * | 2014-01-30 | 2018-03-20 | 国際航業株式会社 | 点群画像による図化方法、及び点群画像による図化装置 |
| US11657437B2 (en) | 2017-11-29 | 2023-05-23 | Angelswing Inc | Method and apparatus for providing drone data by matching user with provider |
| KR102058626B1 (ko) | 2017-11-29 | 2019-12-23 | (주)엔젤스윙 | 지형변화량 가시화 방법 및 장치 |
| CN108986207B (zh) * | 2018-06-29 | 2023-05-23 | 广东星舆科技有限公司 | 一种基于真实路面数据的道路及沿线建筑仿真建模方法 |
| CN113642065B (zh) * | 2021-06-28 | 2024-03-12 | 中国地质大学（武汉） | 一种基于dxf矢量平面图的室内半自动制图与建模方法 |
| CN113658336B (zh) * | 2021-07-28 | 2024-04-30 | 贵州图智信息技术有限公司 | 一种三维自动化建模方法及系统 |
-
2002
- 2002-04-26 JP JP2002127672A patent/JP4058293B2/ja not_active Expired - Lifetime
Also Published As
| Publication number | Publication date |
|---|---|
| JP2003323640A (ja) | 2003-11-14 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| JP4058293B2 (ja) | レーザスキャナデータと空中写真画像を用いた高精度都市モデルの生成方法及び高精度都市モデルの生成システム並びに高精度都市モデルの生成のプログラム | |
| CN109872397B (zh) | 一种基于多目立体视觉的飞机零件的三维重建方法 | |
| US7509241B2 (en) | Method and apparatus for automatically generating a site model | |
| Zhou et al. | Seamless fusion of LiDAR and aerial imagery for building extraction | |
| Chen et al. | Fusion of LIDAR data and optical imagery for building modeling | |
| JP4719753B2 (ja) | 異種センサ統合モデリングによるデジタル写真測量方法及びその装置 | |
| CN110866531A (zh) | 一种基于三维建模的建筑物特征提取方法、系统及存储介质 | |
| CN107767456A (zh) | 一种基于rgb‑d相机的物体三维重建方法 | |
| JP6534296B2 (ja) | ３次元モデル生成装置、３次元モデル生成方法、及びプログラム | |
| Chen et al. | Building reconstruction from LIDAR data and aerial imagery | |
| CN103839286B (zh) | 一种对象语义约束的真正射影像优化采样方法 | |
| JP4568845B2 (ja) | 変化領域認識装置 | |
| CN112767459A (zh) | 基于2d-3d转换的无人机激光点云与序列影像配准方法 | |
| Wang et al. | A method for detecting windows from mobile LiDAR data | |
| CN115690138A (zh) | 一种融合车载影像与点云的道路边界提取与矢量化方法 | |
| Verhoeven et al. | From 2D (to 3D) to 2.5 D: not all gridded digital surfaces are created equally | |
| Ebrahimikia et al. | True orthophoto generation based on unmanned aerial vehicle images using reconstructed edge points | |
| Zhao et al. | Alignment of continuous video onto 3D point clouds | |
| CN116051980B (zh) | 基于倾斜摄影的建筑识别方法、系统、电子设备及介质 | |
| Novacheva | Building roof reconstruction from LiDAR data and aerial images through plane extraction and colour edge detection | |
| JP3966419B2 (ja) | 変化領域認識装置および変化認識システム | |
| Habib et al. | Integration of lidar and airborne imagery for realistic visualization of 3d urban environments | |
| KR101103491B1 (ko) | 항공 레이저 측량 데이터를 이용한 도로 레이어 취득 시스템 및 그 방법 | |
| KR101079359B1 (ko) | 항공 사진과 항공 레이저 측량 데이터를 이용한 수치지도 구축 시스템 | |
| KR101083902B1 (ko) | 항공 레이저 측량 데이터의 필터링을 이용한 3차원 공간정보 구축 시스템 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20041117 |
|
| A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20070815 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20070821 |
|
| A521 | Request for written amendment filed |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20071022 |
|
| 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: 20071211 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20071217 |
|
| R150 | Certificate of patent or registration of utility model |
Ref document number: 4058293 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20101221 Year of fee payment: 3 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20111221 Year of fee payment: 4 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20121221 Year of fee payment: 5 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20131221 Year of fee payment: 6 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| EXPY | Cancellation because of completion of term |