JP2014186567A - ３次元点群解析方法 - Google Patents

Abstract

【課題】 複数の地物の外形を表す３次元点群を解析して、地物単位の点群に分離する。
【解決手段】 複数の地物の外形を表す３次元点群データを用意する。まず、３次元点群データを高さ方向に階層に区分する。そして、各階層において、点間の距離に基づいて、点群を同じ地物に属するクラスタに分離する。また、各クラスタについて地物の上面か側面かの属性を付与する。そして、上面のクラスタから順に下の階層に対し、相互に重なっているクラスタを関連づけることによってツリー構造を生成する。全地物のクラスタが関連づけられている下位のクラスタを地表面と判断し、その上の階層に属するクラスタに基づいて、３次元点群を地物単位に分離する。
こうすることによって、分岐などがある複雑な形状の地物が含まれている場合でも、３次元点群を適切に地物単位に分離することができる。
【選択図】 図１

Description

本発明は、複数の地物等の形状を表した３次元点群を解析する３次元点群解析方法に関する。

地物の形状を３次元的に表した３次元地図においては、地物の３次元モデルを精度良く、かつ低い処理負荷で生成することが求められる。このような技術の一つとして、３次元復元の利用が検討されている。
３次元復元処理とは、ＳＩＦＴ（Scale-Invariant Feature Transform）などのように、複数の地点から撮影した複数の画像間の特徴点の対応関係に基づいて、撮影された地物の３次元点群を得る処理である。３次元復元を利用すれば、処理の対象範囲内の画像を航空機等から複数撮影し、この画像を解析することで、対象範囲を表す３次元の点群を得ることができる。かかる点群が得られれば、それを解析することによって３次元モデルを生成することも可能となる。
３次元モデルを生成するためには、３次元点群を地物単位に分離する必要がある。
特許文献１、２は共に、３次元点群をレイヤに区分し、上のレイヤの凸包が、下のレイヤの点を内包しないという条件を付した上で、レイヤごとの点群の並びを順次併合しながら２次元形状を求める技術を開示している。

特開２００４−２４０５９９号公報 特開２００４−３０２２６号公報

３次元点群は、独立した点の集合に過ぎないため、それぞれの点がいかなる面を構成し、どの地物を構成するかを特定することは容易ではない。特許文献１、２に開示された方法も、１つの解析方法ではあるが、例えば、ツインタワーのように上方が分岐した建物や、中間の階層だけ分岐しているような建物のような特殊な形状の建物については、必ずしも適切な２次元形状を得ることができず、３次元点群を地物単位に分離できない場合があった。本発明は、かかる課題に基づき、３次元点群を解析して、地物単位に３次元点群を分離する技術を提供することを目的とする。

本発明は、
複数の地物の外形を表す複数の点からなる３次元点群をコンピュータによって解析する３次元点群解析方法であって、
前記３次元点群を構成する３次元点群データは、前記地物の外形を表す複数の点に対し、点ごとに位置座標が格納されており、
前記３次元点群解析方法は、前記コンピュータが実行するステップとして、
（ａ） 前記３次元点群データを読み込むステップと、
（ｂ） 前記３次元点群データを、高さ方向に複数の階層に区分するステップと、
（ｃ） 前記各階層において、前記３次元点群を点間の最短距離が所定値以上開いている部分で複数のクラスタに分離するステップと、
（ｄ） 前記上下に隣接する階層間で、各階層に属する前記クラスタのうち、相互に少なくとも一部が重なるクラスタ同士を同一地物に属するクラスタとして関連づけることにより、前記３次元点群を地物単位に分離するステップとを備える
３次元点群解析方法として構成することができる。
ここで、本発明で扱う３次元点群は、例えばＳＩＦＴ、レーザー計測など種々の方法で得られたものを対象とすることができる。

本発明によれば、階層に区分するため、特殊な形状の地物であっても、階層内でみれば、同じ地物を構成するクラスタに点群を分離することができる。本発明は、各地物の外形形状を決めることが目的ではなく、３次元点群を地物単位に分離することが目的であるから、階層内の処理についても、点群を塊として分離するだけでよいため、比較的容易かつ精度よく分離することができる。
このように階層に分けて点群を分離した後、本発明では、異なる階層に属するクラスタ同士の位置関係に基づき、同一地物に属するクラスタを集約する。このように、階層内での処理と、階層間にわたる集約とを行うことにより、本発明によれば、地物単位のクラスタ、ひいては地物単位の３次元点群を、比較的容易かつ精度良く特定することができる。

異なる階層間のクラスタの一部が重なるか否かは、種々の方法によって判断することができる。例えば、異なる階層のクラスタに属する点群を、同一の水平面に投影した上で、異なるクラスタに属する点群間の距離の最小値を求め、これが所定値以下の場合には、重なっているものと判断する方法をとってもよい。また、各クラスタの点群を包含する凸形状を特定し、この凸形状同士が重なるか否かで判断する方法や、凸形状の重心位置の距離が所定値以下か否かに基づいて判断する方法をとってもよい。これらの種々の方法は、いずれかを選択して用いても良いし、併用してもよい。

本発明においては、
前記ステップ（ｄ）は、階層間で１対多または多対１の関連付けも許容して前記関連付けを行うものとしてもよい。
こうすることにより、上の階層または中間の階層で分岐しているような特殊な形状の地物も、関連付けを行うことができ、地物ごとの３次元点群の分離を実現することができる。１対多とは、上階層の１つのクラスタが、下階層の２以上のクラスタに関連づけられていることを言い、例えば、単一の屋上部分を複数本の柱で支える形状の地物などでこうした関連付けが生じ得る。また、多対１とは、逆に、上位層で複数のクラスタが、下位層の１つのクラスタに関連づけられていることを言い、例えば、ツインタワーのように上方で分岐している地物などでこうした関連付けが生じ得る。
こうした関連付けは、いわゆるツリー構造として表すことが可能である。

上述の態様のように関連付けを行う場合において、
前記３次元点群は、地表面を表す点群も含んでいるときは、
前記ステップ（ｄ）は、前記関連付けの結果、所定割合以上のクラスタが集約して関連づけられている低階層のクラスタを地表面と判断し、該クラスタが存在する階層の上に位置する階層におけるクラスタによって前記地物単位への分離を行うものとしてもよい。

３次元点群が地表面も含んでいる場合には、地表面を含む階層も生成されることになる。そして、全ての地物は、最終的には地表面につながっているはずであるから、かかる点群の欠損等を考慮したとしても、クラスタ間の関連付けを表すツリー構造において、所定割合以上のクラスタが集約しているノードが地表面ということになる。地表面の上の階層に属するクラスタは、それぞれの地物の最下位の階層を表していると考えられるから、この階層で分離することにより、上位の層で複雑な形状をなしている地物であっても、精度良く地物単位に３次元点群を分離することが可能となる。

本発明においては、さらに、
前記ステップ（ｃ）で得られた各階層のクラスタに対し、該クラスタ内に前記地物の上面を構成する点が含まれる場合には地物の上面を表す旨の属性を付すステップを有し、
前記ステップ（ｄ）における関連付けは、上面の属性が付されたクラスタから、下側の階層に向かって順次行うものとしてもよい。
上面か否かの判断は種々の方法で行うことができる。例えば、３次元点群を構成する各点について、位置座標と併せて、その点が属すると推測される面の法線ベクトルが格納されている場合には、法線ベクトルが重力方向に沿っているか否か、即ち重力方向と法線ベクトルとのなす角が所定値以下であるか否かに基づいて上面か否かを判断するものとしてもよい。また、別の方法として、クラスタの点群を水平面に投影し、輪郭部分だけでなく、その内面にも多数の点群が存在する場合には上面であると判断するようにしてもよい。
このように上面を特定した上で関連付けを行うことにより、地物の構造上の端点が明確にした状態で関連づけを行うことができるため、精度向上を図ることができる。

本発明は、上述した特徴を必ずしも全てを備えている必要はなく、適宜、その一部を省略したり、組み合わせたりしてもよい。上述した３次元点群解析方法を実行する３次元点群解析装置として構成してもよいし、かかる解析をコンピュータに行わせるためのコンピュータプログラムとして構成してもよい。さらに、かかるコンピュータプログラムを記録したコンピュータが読み取り可能な記録媒体として構成してもよい。

３次元モデル生成システムの構成を示す説明図である。 対象領域の空撮方法を示す説明図である。 飛行経路決定処理のフローチャートである。 ３次元点群生成処理のフローチャートである。 ３次元モデル生成処理のフローチャートである。 重力方向解析処理のフローチャート（１）である。 重力方向解析処理のフローチャート（２）である。 重力方向解析処理のフローチャート（３）である。 重力方向解析処理のフローチャート（４）である。 点群分離処理のフローチャート（１）である。 点群分離処理のフローチャート（２）である。 点群分離処理のフローチャート（３）である。 点群分離処理のフローチャート（４）である。 建物枠生成処理のフローチャート（１）である。 建物枠生成処理のフローチャート（２）である。 位置合わせ処理のフローチャートである。 位置合わせ処理における評価関数を示す説明図である。 ポリゴン生成処理例を示す説明図である。 テクスチャ生成処理のフローチャート（１）である。 テクスチャ生成処理のフローチャート（２）である。

Ａ．システム構成：
図１は、３次元モデル生成システムの構成を示す説明図である。３次元モデル生成システムは、処理対象となる地域の写真を航空機から撮影し、この写真を解析して地物ごとの３次元形状を表す３次元モデルを生成するシステムである。航空機に搭載され、地域の写真の画像データを取得する装置を撮影管理装置１００と呼ぶ。画像データから３次元モデルを生成する装置を、３次元モデル生成装置２００と呼ぶ。以下、それぞれの構成について順に説明する。

（１）撮影管理装置１００の構成：
撮影管理装置１００は、航空機に搭載されている。航空機は、固定翼、回転翼いずれでもよく、有人／無人のいずれであってもよい。
撮影管理装置１００には、撮影用のカメラ１０１、位置情報を取得するためのセンサであるＧＰＳ（Global Positioning System）１０２、および図示する各機能ブロックが備えられている。機能ブロックは、それぞれの機能を実現するためのコンピュータプログラムを、コンピュータにインストールすることによって構成することもできるし、ハードウェア的に構成してもよい。
カメラ１０１は、３次元モデルを生成するための解析に使用される写真を撮影するためのものである。解析精度を向上させるためには、高解像度の写真が望ましい。本実施例では高解像度での撮影が可能な静止画用のディジタルカメラを用いた。十分な解像度が得られる場合には、動画撮影用のカメラを用いてもよい。

撮影シーケンス記憶部１１２は、撮影を行うための航空機の飛行経路および撮影ポイントを記憶している。撮影ポイントは、カメラ１０１の画角も考慮して設定されており、処理対象となる領域について、２０枚以上の写真が、それぞれ６０％程度以上は重複した状態で、順次撮影されるように設定されている。飛行経路の設定方法については後述する。飛行経路等は、予め別のコンピュータ等によって設定したものを、撮影シーケンス記憶部１１２に記憶させる方法をとってもよいし、撮影シーケンス記憶部１１２で設定するようにしてもよい。
有人の航空機を用いる場合、飛行経路に従って航空機を飛行させるのはパイロットの役目となる。かかる場合には、撮影シーケンス記憶部１１２は、パイロットに対し、適宜、設定された飛行経路に沿った飛行を支援するための情報を提供してもよい。
無人の航空機を用いる場合には、撮影シーケンス記憶部１１２は、飛行を制御する制御装置に、飛行経路に関する情報を送信するように構成してもよい。

撮影制御部１１０は、カメラ１０１を制御し、撮影シーケンス記憶部１１２で設定された撮影ポイントにおいて処理対象領域の写真を撮影する。航空機の現実の飛行経路と、撮影シーケンス記憶部１１２に記憶された経路とは誤差があるため、撮影制御部１１０は、かかる誤差も踏まえて、予め設定された撮影ポイント付近に到達したと判断される時点で撮影を行うようにしてもよい。
また、撮影シーケンスは、必ずしも撮影ポイントを特定する方法だけでなく、一定の時間周期で撮影を行うように設定してもよい。かかる場合には、撮影制御部１１０は、設定された時間間隔ごとに撮影を行うようカメラ１０１を制御することになる。
撮影制御部１１０は、いずれの方法で撮影を行う場合においても、カメラ１０１による撮影と同期して、ＧＰＳ１０２から撮影位置の位置情報を取得し、両者を関連づけて画像記録部１２０に記録する。
撮影制御部１１０は、自動的に撮影を行う態様として構成する他、航空機に搭乗した撮影者に対して、撮影タイミングを知らせるものとしてもよい。

（２）３次元モデル生成装置２００の構成：
３次元モデル生成装置２００には、図示する機能ブロックが備えられている。これらの機能ブロックは、それぞれの機能を実現するためのコンピュータプログラムを、コンピュータにインストールすることによって構成されるし、ハードウェア的に構成してもよい。
以下、各機能ブロックの機能について説明する。

画像データ記憶部２０１は、撮影管理装置１００で撮影、収集された画像データおよび撮影位置の位置情報を読込み、記憶する。撮影管理装置１００から画像データ記憶部２０１への画像データの受け渡しは、有線、無線または記録媒体を介して行うことができる。
３次元点群生成部２０２は、画像データ記憶部２０１の画像データおよび撮影位置を用いて、３次元点群を生成する。３次元点群とは、それぞれ３次元の座標値および法線ベクトルが付与された多数の点である。法線ベクトルとは、各点が存在すると推定される面の法線ベクトルを意味する。
本実施例では、画像データを撮影した時系列に並べると、前後の画像データ同士は、約６０％程度以上の面積が重複した状態で撮影されている。つまり、対象領域内の建物は、撮影ポイントが異なる複数の画像データに撮影されていることになる。従って、複数の画像データ内で相互に対応する点を特定できれば、いわゆる三角測量と同様の原理によってその点の３次元座標を特定することが可能となる。複数の画像データ内で対応する点を特定する技術としては、種々の技術が知られているが、本実施例では、ＳＩＦＴ（Scale-Invariant Feature Transform）と呼ばれる公知の技術を適用した。
３次元点群データ記憶部２０３は、３次元点群生成部２０２によって得られた３次元点群を記憶する。画像データに写っている部分を解析して得られた点群であるから、３次元点群は、主として建物等の外壁、屋根などの外面部分を表すものとなる。また、画像内には、地面、道路、街路樹、車両など建物以外も多数写っており、これらについても画像内の対応関係が特定される限り３次元点群は生成可能であるから、３次元点群には、こうした建物以外を表すものも多数含まれている。３次元点群は、それぞれ３次元の座標値と法線ベクトルを有しているが、これらの座標値、法線ベクトルは、地球座標系で定義されたものではなく、解析のために設定された架空の３次元空間における座標値である。
このように地球座標系での位置関係が特定されておらず、建物単位に区分もされていない３次元点群から、建物ごとに地図として利用可能な精度で位置座標を有する３次元ポリゴンを生成するのが、以下で示す各機能ブロックの機能である。

重力方向解析部２０４は、点群の重力方向を解析する。通常、直交座標系で地理的な位置が表されている場合には、重力方向はＺ軸方向と一致するが、３次元点群の座標系は地球座標系で定義されたものではないため、必ずしもＺ軸が重力方向になるとは限らない。重力方向解析部２０４は、解析によって重力方向を求めることによって、点群の座標系を地球座標系に対応づける第1段階としての機能を奏する。
点群分離処理部２０５は、点群を建物単位に分離する。３次元点群は、画像に写っている点の３次元座標を求めたものであり、建物単位で分離もされていなければ、建物の地面等との分離もされていない。点群分離処理部２０５は、かかる３次元点群を解析して、建物ごとの点群に分離する機能を奏する。
建物枠生成部２０６は、建物ごとに分離された点群に基づいて２次元的な外形形状を表す建物枠を生成する。
位置合わせ処理部２０７は、２次元地図データ２１０に記憶された２次元地図に複数の建物枠を適合させることによって、重力方向以外の座標系を地球座標系に対応づける機能を奏する。位置合わせ処理部２０７は、建物ごとに分離された点群を1つずつ２次元地図に適合させていくのではない。点群で表されている複数の建物全体を、その相対的な関係を概ね保った状態で、全体として２次元地図に適合させるのである。このように複数の建物全体を適合させることによって、精度良く位置合わせを行うことができる。
ポリゴン生成部２０８は、以上の処理によって２次元的な形状が定まった建物ごとに、点群に基づき側面および上面のポリゴンを生成する。
テクスチャ生成部２０９は、生成されたポリゴンの表面に貼り付けるテクスチャを、画像データ記憶部２０１に記憶された画像データによって生成する。
以上の処理で生成された３次元ポリゴンおよびテクスチャは、３次元モデルデータ２１１に格納される。

Ｂ．対象領域の空撮：
次に、解析に利用する画像を撮影する飛行経路および撮影ポイントの設定方法について説明する。
図２は、対象領域の空撮方法を示す説明図である。本実施例では、目標となる建物の周辺の所定の領域に含まれる複数の建物について３次元モデルを生成する。一般に街並みは、駅などのランドマークを中心に展開していることが多いため、こうしたランドマークに相当する建物を目標として、その周辺を撮影し、３次元モデルを生成することにより、効率的に３次元モデルを生成することができる。

図２（ａ）に示すように、本実施例では、航空機を用いて３次元モデルの生成対象となる対象領域の画像を撮影する。航空機は、固定翼、回転翼のいずれでもよい。飛行経路は、目標となる建物の周囲を周回する円軌道を描く。円を描くことによって、生成対象の領域を多方面から均等に撮影することができるため、他の建物の陰になって視認できない部分を減らし、３次元モデルの生成精度を向上させることができる。飛行経路は、目標建物の周囲を周回する軌道であれば、必ずしも円に限られるものではなく、楕円その他の曲線状であってもよい。

図２（ｂ）には、円軌道とする場合の飛行経路の決定方法を示した。３次元復元では、複数の画像内に共通に撮影されている特徴点同士の対応関係を特定することで，各特徴点の３次元位置を求める。かかる対応関係を精度良く求めるためには、各特徴点が多くの画像に重複して撮影されるように、飛行経路および撮影ポイントを設定する必要がある。本実施例では、撮影された画像同士が、約６０％以上の面積で重複するように、飛行経路および撮影ポイントを設定した。つまり、図２（ｂ）に示すように、ある撮影ポイントＰ１から撮影した場合の画像ＩＭＧ１と、次の撮影ポイントＰ２で撮影した場合の画像ＩＭＧ２とは、約６０％以上の範囲で重複する。図の例では、画像ＩＭＧ１、ＩＭＧ２は、目標を中心に回転移動した状態にある例を示したが、目標が中心からずれていても構わない。
撮影に用いるカメラの画角および飛行高度を決定すれば、目標を含む所定範囲を撮影するための飛行経路の概ねの半径を求めることができる。

目標周辺を１周する間の撮影枚数は、２０枚以上とすることが望ましい。従って、撮影ポイントＰ１、Ｐ２の間の角度θは、１８度以下とすることが望ましい。撮影枚数があまりに多数であっても、解析に利用できない無駄な画像が増えるだけであるため、動画像を撮影するカメラを用いて多数のフレームを撮影するよりは、静止画撮影用のカメラを用いて、高解像度の画像を必要枚数撮影する方が好ましい。

図３は、飛行経路決定処理のフローチャートである。図２で説明した飛行経路を決定するための処理である。この処理は、飛行前に地上のコンピュータで予め行っておくものとしてもよいし、飛行中に航空機に搭載した撮影管理装置１００で実行するものとしてもよい。
処理を開始すると、コンピュータはまず、撮影目標地点を入力する（ステップＳ１）。また、カメラの取付角度ｒおよび画角を入力する（ステップＳ２）。図中にカメラの取付角度ｒを示した。本実施例のように、目標物の周囲を周回しながら撮影する場合には、カメラの画角の中心軸を鉛直下方に向けて取り付けるよりも、飛行経路の中心側に向けて傾けて取り付けた方が効率的な場合がある。航空機の鉛直軸と、カメラの画角の中心軸とのなす角度が、取付角度ｒである。取付角度ｒは、カメラの撮影範囲に影響することになる。
コンピュータは、撮影高度Ｈを入力する（ステップＳ３）。そして、以上で入力した値を用いて、飛行経路の半径Ｒを算出する（ステップＳ４）。図中に半径Ｒの算出方法を示した。カメラの画角および取付角度ｒを考慮した上で、撮影目標が撮影範囲に入るよう、飛行高度Ｈに応じた半径Ｒを算出する。例えば、撮影目標が撮影範囲の中心に来るように撮影する場合には、飛行経路の半径Ｒは、Ｒ＝Ｈ×ｔａｎ（ｒ）で求められることになる。
撮影ポイントは、飛行経路に対して、図２（ｂ）に示したように、所定の角度間隔θを設定するようにしてもよいし、所定の時間間隔で撮影するように設定してもよい。

こうして設定された飛行経路および撮影ポイントの情報は、撮影管理装置１００の撮影シーケンス記憶部１１２に記憶され、撮影時に活用される。例えば、航空機が無人の場合には、航空機およびカメラの自動制御に活用することができる。航空機が有人の場合には、飛行経路を維持するのはパイロットの役目となるから、パイロットに適宜、飛行経路および撮影ポイントを示す情報を提供する方法で活用することができる。

以上で説明した本実施例の空撮方法によれば、目標の周囲を周回する飛行経路をとることにより、次のような利点がある。
本実施例では、画像間で特徴点の対応関係を解析するため、画像同士が重複して撮影されている必要があるところ、本実施例のような周回軌道であれば、目標周辺の領域の画像を、相互に十分重複させながら、多方面からの画像を、無駄なく効率的に撮影することができる。
また、画像間の特徴点の対応関係の認識精度が向上する利点もある。つまり、画像間で特徴点を認識する際には、画像の色成分を利用することになるため、撮影時のライティングの条件が異なると、対応関係の認識精度が低下する。これに対し、本実施例のような周回軌道であれば、目標周辺の領域を、比較的短時間に集中的に撮影することが可能となるため、ライティングの条件の変化が小さくなり、認識精度を向上させることができる。
さらに、目標および対象領域を、ほぼ全方位から漏れなく撮影することが可能となる。

Ｃ．３次元点群生成処理：
飛行経路および撮影ポイントが決まると、これに従って、目標の周辺を飛行しながら撮影が行われる。撮影された画像データは、図１で説明した通り、逐次、撮影ポイントの位置座標とともに、撮影管理装置１００の画像記録部１２０に記録される。撮影が完了すると、画像データおよび位置座標は、３次元モデル生成装置２００に入力され、３次元点群が生成される。画像データの３次元モデル生成装置２００への読込みは、種々の方法で行うことができるが、本実施例では、記録媒体を介して読み込ませるものとした。

図４は、３次元点群生成処理のフローチャートである。３次元モデル生成装置２００の３次元点群生成部２０２が実行する処理であり、ハードウェア的には、３次元モデル生成装置２００のＣＰＵが実行する処理である。
処理を開始すると、３次元モデル生成装置２００は、画像データおよび位置情報を読み込む（ステップＳ１０）。そして、ＳＩＦＴ（Scale-Invariant Feature Transform）マッチングにより、画像間の特徴点の対応付けを行う（ステップＳ１１）。図中に特徴点の対応付けの様子を示した。図示するように画像１、画像２のそれぞれから、特徴点（図中の４つの×印）を抽出し、それぞれの特徴点同士の対応関係を認識するのである。ＳＩＦＴマッチングは、画像内の色成分に基づいて特徴点を抽出し、画像間の対応関係を求める解析技術であり、解析原理およびアルゴリズムともに周知であるため、詳細な説明は省略する。ステップＳ１１では、画像間の特徴点同士を対応づけるものであれば、ＳＩＦＴマッチングに限らず、種々の手法を用いることができる。

こうして特徴点の対応関係が求まると、各特徴点の３次元座標および法線ベクトルを算出することができる（ステップＳ１２）。図中に算出原理を図示した。実施例で使用する各画像には、それぞれ撮影ポイントの位置情報が付されているため、複数の画像間で対応関係が得られた特徴点については、２つ以上の撮影ポイント１，２が対応づけられることになる。撮影ポイント１，２からの画角等に基づいて撮影した画像内の特徴点の位置を解析すれば、それぞれの撮影ポイント１，２から特徴点を視た方向を決定することができ、撮影ポイント１、２の座標値に対して特徴点の相対的な位置関係を特定することができる。ここで得られる特徴点の位置座標は、地球座標系での絶対的な位置座標ではなく、解析用の空間に設定された架空の３次元座標（ｘａ、ｙａ、ｚａ）上での座標値である。
法線ベクトルとは、各特徴点が存在すると推定される推定面の法線方向のベクトルである。実施例では、対応する画像において多数の特徴点が抽出され、各画像内の特徴点同士の相対的な位置関係が維持されることを前提として、画像間での特徴点の対応関係を認識する。このように特徴点同士の相対的な位置関係から、その特徴点が存在する推定面が特定され、その法線ベクトルを求めることができる。
これらの画像間の相対的位置関係、３次元座標および法線ベクトルの算出も、ＳｆＭ（Structure from Motion）技術、およびそれを利用した手法として周知であるため、さらなる詳細な説明は省略する。
３次元座標および法線ベクトルが算出された特徴点の集合を、以下、３次元点群と呼ぶ。３次元モデル生成装置２００は、得られた３次元点群のデータを格納し（ステップＳ１３）、３次元点群生成処理を終了する。
本実施例では、撮影した画像データから３次元点群を取得しているが、以下の処理は、レーザ計測など、他の方法によって得られた３次元点群を用いて行うこともできる。

Ｄ．３次元モデル生成処理：
図５は、３次元モデル生成処理のフローチャートである。先に説明した３次元点群データに基づいて、建物ごとの３次元形状を表す３次元モデルを生成するための処理である。この処理は、図１に示した重力方向解析部２０４、点群分離処理部２０５、建物枠生成部２０６、位置合わせ処理部２０７、ポリゴン生成部２０８、テクスチャ生成部２０９によって、それぞれ実現される処理である。ハードウェア的には、３次元モデル生成装置２００のＣＰＵが実行する処理となる。

図５には、全体の処理の流れを示した。それぞれの処理内容については、後で詳述する。
３次元モデル生成処理では、まず重力方向解析処理が行われる（ステップＳ２０）。これは、重力方向解析部２０４の機能に相当する処理であり、３次元点群を解析して重力方向を特定する処理である。
次に、点群分離処理が行われる（ステップＳ２１）。これは、点群分離処理部２０５の機能に相当する処理であり、３次元点群を建物単位に分離する処理である。
次に、建物枠生成処理が行われる（ステップＳ２２）。これは、建物枠生成部２０６の機能に相当する処理であり、建物の２次元的な形状を特定する処理である。
次に、位置合わせ処理が行われる（ステップＳ２３）。これは、位置合わせ処理部２０７の機能に相当する処理であり、３次元点群を２次元の地図に適合させる処理である。
次に、ポリゴン生成処理が行われる（ステップＳ２４）。これは、ポリゴン生成部２０８の機能に相当する処理であり、建物ごとの３次元点群に基づいて、側面および上面のポリゴンを生成する処理である。
最後に、テクスチャ生成処理が行われる（ステップＳ２５）。これは、テクスチャ生成部２０９の機能に相当する処理であり、撮影された画像データを用いてポリゴンの表面に貼るテクスチャを生成する処理である。
以上の処理によって、３次元地図として利用可能な位置精度を有する建物ごとの３次元モデルが生成される。
以下、それぞれの処理について順に内容を詳述する。

Ｄ１．重力方向解析処理：
図６〜９は、重力方向解析処理のフローチャートである。３次元モデル生成処理（図５）のステップＳ２０における処理に相当する。
処理を開始すると、３次元モデル生成装置２００は、３次元点群データを読み込む（ステップＳ３１）。
そして、３次元点群全体の主成分分析を行い、第３主成分を仮の重力方向Ｇｔとして設定する（ステップＳ３２）。図中に主成分分析の概要を示した。本実施例の３次元点群は、図２に示すように目標を中心とする周辺の領域を撮影して得られるものである。撮影される建物が高層ビルであっても、その高さは、撮影される領域の水平方向の広がりに比べれば小さい。従って、３次元点群は、全体的に見れば、図示するように、比較的平たい領域に分布することとなる。
主成分分析とは、３次元空間内における３次元点群の分布を特徴づける方向を求めるものであり、３次元点群の分布が大きい方向から順に第１主成分、第２主成分、第３主成分となる。平たい領域に分布した３次元点群については、水平方向の分布に基づいて第１主成分、第２主成分が求まり、それに直交する方向が第３主成分となる。従って、第３主成分は、３次元点群の重力方向Ｇｔを表すことになる。ただし、点群の分布に基づく誤差が含まれているため、第３主成分は、あくまでも仮の重力方向を示すに過ぎない。
本実施例では、このように平たい領域に分布する３次元点群を利用するため、垂直方向と推測される面のみを抽出するなどといった複雑な前処理を行うまでなく、主成分分析により仮の重力方向を求めることが可能となる。

３次元モデル生成装置２００は、次に、仮の重力方向Ｇｔとの角度θが閾値ＴＨａ以下の点群を除外する（ステップＳ３３）。実施例の３次元点群は、３次元の位置座標とともに、法線ベクトルを有している。この法線ベクトルと重力方向Ｇｔとの角度θを０〜１８０度の範囲で求め、この値が閾値ＴＨａ以下の点群、即ち、重力方向に比較的沿う方向の法線ベクトルを有する点群を除外すればよい。閾値ＴＨａは、点群の法線ベクトルの誤差および仮の重力方向の誤差を考慮して、法線ベクトルがほぼ重力方向を向いている点群を除外できる値に設定すればよい。
図中に処理の様子を示した。図中の×および○は、それぞれ３次元点群を示している。実施例の点群は、建物を撮影した画像から得られるものであるから、建物の外壁に沿って分布する。点Ｐ１は、建物の水平面Ｓ１に含まれる点であり、法線ベクトルｎ１は上向きとなっている。法線ベクトルｎ１と、仮の重力方向Ｇｔとの角度θ１は、閾値ＴＨａよりも小さくなり、ステップＳ３３の処理において除外される。同様に、面Ｓ１に存在する×で示した他の点群も除外される。
点Ｐ２は、建物の側面Ｓ２に含まれる点であり、法線ベクトルｎ２は、水平に近い方向を向いている。法線ベクトルｎ２と仮の重力方向Ｇｔとの角度θ２は、閾値ＴＨａよりも大きくなり、ステップＳ３３の処理において除外されない。同様に、面Ｓ２に存在する○で示した他の点群も除外されない。
点Ｐ３は、建物の水平面Ｓ３に含まれる点であり、法線ベクトルｎ３は下向きとなっている。このような面Ｓ３は、建物の床面等に相当するものであり、外部から撮影した画像に写ることはないため、実施例における処理では、このような点群は得られない。ただし、ここでは説明の便宜上、法線ベクトルｎ３が下向きとなっている場合についての例も示したものである。点Ｐ３についても、法線ベクトルｎ３と、仮の重力方向Ｇｔとの角度θ１は、閾値ＴＨａよりも小さくなるため、ステップＳ３３の処理において除外される。同様に、面Ｓ３に存在する×で示した他の点群も除外される。
以上の通り、ステップＳ３３の処理では、建物を構成する点群のうち、水平面を構成する点群を除外することができる。

図７に移り、３次元モデル生成装置２００は、抽出された点群について、点間の距離によるクラスタリングを行う（ステップＳ３４）。つまり、比較的近くにある点群同士をクラスタ化することによって、点群を建物単位で分離するのである。
図中に処理の様子を示した。左側の図は、クラスタリングの方法を模式的に示したものである。３次元モデル生成装置２００は、処理対象となる点群の中から、任意に1つの点Ｐ１を選択し、固有のラベル（図の例では「Ａ」）を付与する。そして、この点Ｐ１から所定の距離ｒａ内にある点を探索し、これらに同じラベルを付与する。図の例では、点Ｐ２，Ｐ３に対してラベル「Ａ」が付与されている。距離ｒａは、クラスタリングするか否かの判断基準となる距離である。クラスタリングの目的は、点群を建物ごとに分離することであるから、距離ｒａは、処理対象となる領域における建物の間隔よりも小さい範囲で任意に設定することができる。
点Ｐ２、Ｐ３にラベルを付すと、３次元モデル生成装置２００は、これらの点からさらに距離ｒａ内の点を探索し、これらの同じラベルを付与する。図の例では、点Ｐ４、Ｐ５にラベル「Ａ」が付与されている。かかる処理を繰り返し実行することによって、相互に距離ｒａ内にある点を探索することができ、それらに同じラベルを付与することができる。距離ｒａ内にある点が探索できなくなった時点で、ラベルの付与されていない未処理の点を選択し、新たなラベル（例えば、「Ｂ」）を付与しながら、同様の処理を実行することによって、点間の距離によるクラスタリングを行うことができる。
図の右側に、クラスタリングの結果を例示した。図示するように、ビル１に属する点群（図中の×）、ビル２に属する点群（図中の○）というように、建物単位に点群を分離することができる。

次に、３次元モデル生成装置２００は、建物ごとに点群の法線ベクトルの水平方向分布を作成する（ステップＳ３５）。ステップＳ３３の処理で建物の水平面を除く点群、即ち側面を構成する点群が抽出されており、ステップＳ３４の処理で建物ごとに分離されているから、その法線ベクトル分布を作成することにより、側面の法線ベクトルを求めることができる。
図中に処理の様子を示した。左側には、法線ベクトルの方向を決める方法を示している。３次元モデル生成装置２００は、処理対象の建物の任意の点Ｐ１の法線ｎ１に対して、その水平方向の成分ｎｈ１を求める。ここに言う水平方向の成分とは、仮の重力方向Ｇｔに直交する成分という意味である。そして、水平成分ｎｈ１の方向を、基準軸からの角度Ａｇ１で表す。この角度Ａｇ１は、図の右側に示すように０〜１８０度の範囲で求めるものとした。基準軸からの角度は、反時計回りを正、時計回りを負として表されるのが通常であるが、ここでは分布を求めるという処理の性質上、０〜１８０度の範囲で絶対値を求めるのである。点Ｐ２に対しても同様に、法線ベクトルｎ２の水平成分ｎｈ２の方向を表す角度Ａｇ２が求められる。

水平成分の分布が求まると、３次元モデル生成装置２００は、この分布に基づいて面を分離する（図８のステップＳ３６）。点群がそれぞれ建物の側面を構成していると考えれば、法線ベクトルの分布に基づき、そのピークを求めれば、それらが側面の法線ベクトルを表すと考えられる。従って、ピークごとに点群をクラスタ化することによって側面ごとに点群を分離することができる。
図中に処理の様子を示した。下側の図に示すように、法線ベクトルの方向ごとに、それに該当する点の数を計数し、図示することによって法線ベクトルの分布を得ることができる。上側の図に示すような平面的な側面からなる建物の場合には、法線ベクトルに比較的明瞭なピークが生じる。図の例では、最大のピーク１は、面積が大きい側の側面に対応し、２番目のピーク２は、面積が小さい側の側面に対応する。
従って、法線ベクトルの水平成分が、それぞれのピーク１、ピーク２に該当する点群を抽出することによって、建物の点群を、側面単位に分離することができる。ここで、図の例では、ピーク１、ピーク２に所定の幅を持たせた範囲（図中の枠で囲んだ部分）で、点群を抽出している。この幅は、ピーク１、ピーク２を区別可能な範囲で任意に設定可能である。ただし、実施例では、水平成分を０〜１８０度の範囲で求めているため、面積が最大となる側の表面／裏面の双方がクラスタとして分離されることになり、面積が小さい側の側面についても表面／裏面の双方がクラスタとして分離される。

次に、３次元モデル生成装置２００は、分離された点群ごとに点間の距離によるクラスタリングを行う（ステップＳ３７）。クラスタリングの方法は、ステップＳ３４の処理と同様である。ここでは、クラスタリングの判定基準となる距離は、建物の表面／裏面を分離可能な程度の値とすればよい。
図中に処理の様子を示した。左側に示すように、ステップＳ３６の処理では、建物の表面Ｓ１、裏面Ｓ２が１つのクラスタとして認識されている。これを点間の距離によるクラスタリングを行うと、右側に示すように、距離が離れている表面Ｓ１、裏面Ｓ２が分類される。図中では、表面Ｓ１の点群を×、裏面Ｓ２の点群を○で示した。
以上の処理によって、建物の側面ごとに点群が分離されたことになる。

図９に移り、３次元モデル生成装置２００は、各面ごとに主成分分析を行い、第３主成分Ｖｉ、および第３固有値ｅｉを算出する（ステップＳ３８）。第３主成分は、点群の分布のうち、３番目に分散が大きくなる方向である。第３固有値は、第３主成分の分散値である。
図中に処理の様子を示した。処理対象となる点群は、建物の側面を構成する点群であるから、主成分分析を行うと、面内の分散に応じて第１主成分、第２主成分が求まり、これらに直交する方向、即ち面の法線方向が第３主成分となる。そして、その分散値が第３固有値となる。精度よく点群が得られていれば、点群は精度良く側面上に分布するはずであるから法線方向の分散値は非常に小さくなるはずである。従って、第３固有値が大きいということは、法線方向に点群のばらつきが大きいことを示しており、法線の精度が低いことを表している。この処理を、各面について実行することにより、その面の法線ベクトルを得るとともに、その精度を表す評価値を得ることができる。

３次元モデル生成装置２００は、こうして得られた第３主成分を用いて重力方向Ｇを算出する(ステップＳ３９)。
図中に示すように、任意の２つの側面の法線ベクトルＶｉ、Ｖｊを用いて、これらの外積Ｖｉ×Ｖｊを求めれば、重力方向Ｇが得られるはずである。ただし、法線ベクトルＶｉ、Ｖｊには、それぞれ誤差が含まれているから、こうして得られる重力方向Ｇにも誤差が含まれている。そこで、実施例では、得られている法線ベクトルから、平行ではない任意の組み合わせを選択して、外積により個別の重力方向を求め、その加重平均によって最終的な重力方向Ｇを算出するものとした。
ここで加重平均は、法線ベクトルの精度を表す第３固有値の和（ｅｉ＋ｅｊ）の逆数とした。第３固有値が大きいほど法線ベクトルの精度が低いことを表しているから、第３固有値の逆数を重み値として乗じることにより、法線ベクトルの精度が高い項ほど重み値が大きくなることになる。以上の考え方に基づき、本実施例では、重力方向は、（Ｖｉ×Ｖｊ）／（ｅｉ＋ｅｊ）の平均によって求める。
重力方向Ｇの算出は、上述の方法に限られるものではない。例えば、重み値は、法線ベクトルの精度を反映できればよく、ｅｉ＋ｅｊの他、ｅｉとｅｊの積などとしてもよい。また、全ての法線ベクトルを用いるのではなく、第３固有値が所定の閾値以下の精度が高いと見なすことができる法線ベクトルのみを用いるようにしてもよい。

以上の処理により、３次元点群に基づいて重力方向を決定することができる。本実施例によれば、仮の重力方向の決定（図６のステップＳ３２）や、側面の法線ベクトルの決定（図９のステップＳ３８）などに主成分分析を適用する。このように統計的な処理を施すことにより、それぞれの点の位置座標および法線ベクトルに含まれる誤差の影響を抑制し、高い精度で重力方向を決定することが可能となる。

Ｄ２．点群分離処理：
図１０〜１３は、点群分離処理のフローチャートである。３次元モデル生成処理（図５）のステップＳ２１における処理に相当する。
先に実行した重力方向解析処理において、点間の距離によるクラスタリングを利用して、点群を建物単位、さらには建物の側面単位に分離する処理を実行した（図７のステップＳ３４および図８のステップＳ３７）。しかし、建物には、いわゆるツインタワーのように高層階のみが複数に分岐して建築されているものや、中間の階層のみ複数に分岐して建築されているものなどが存在する。重力方向解析処理におけるクラスタリングでは、こうした特殊な形状の建物について、点群を建物単位に分離することはできない。点群分離処理は、こうした建物の存在も考慮し、点群を建物単位に分離する処理を実行する。点群分離処理は、特殊な形状の建物のみを対象とすることも可能ではあるが、３次元点群の中から、かかる形状の建物を抽出すること自体が、容易ではないため、本実施例では、形状が特殊か否かに関わらず、全建物を対象として処理を行うものとしている。

３次元モデル生成装置２００は、処理を開始すると、処理対象となる点群データを読み込む（ステップＳ５０）。図中に処理対象となる点群を例示した。点群の密度は十分に高いため、表示スケールによっては、このように点群のみでも、既に建物の形状が完成しているように視認できるが、あくまでも点群同士は何ら関連づけられておらず、建物ごとに分離もされていない点の集合に過ぎない。

３次元モデル生成装置２００は、読み込んだ点群を階層ごとに分離する（ステップ５１）。階層の高さＨは、任意に決定可能である。高層階のみが複数に分岐して建築されているなど、特殊な形状の建物であっても、階と階との境で分岐させているのが通常であることを考えると、階層の高さＨは、建物一般の各階の高さ相当の値とすることが好ましい。
図中に階層ごとに分離する様子を示した。図中の×は、３次元点群を表している。階層ごとへの分離とは、３次元点群が存在する３次元空間内に、所定の高さＨごとに境界平面を設定することによって行う。図の例では、左側の建物の屋上面に階層の境界平面が位置するように示してあるが、境界面は、建物形状とは無関係に設定すればよい。

次に、３次元モデル生成装置２００は、階層内で点間の距離によるクラスタリングを行う（ステップＳ５２）。クラスタリングの方法は、重力方向解析処理で説明したのと同様である（図７のステップＳ３４参照）。
図中に処理の様子を示した。処理対象となる階層内に×で示す点群が存在する場合を考える。距離によってクラスタリングを行うと、点群は、左右の各クラスタＣＬ１、ＣＬ２に分離される。屋上に相当する面の点群が、クラスタＣＬ１にのみ含まれているのは、クラスタＣＬ２側の建物は、処理対象となった階層よりも高く、この階層には側面の点群しか含まれないからである。

図１１に移り、３次元モデル生成装置２００は、分離された各クラスタについて、上向き法線の分布を計算する（ステップＳ５３）。上向き法線は、例えば、図中の左側に示すように、点の法線ベクトルｎと重力方向Ｇとの内積の絶対値｜ｎ・Ｇ｜（以下、この処理の説明において、「評価値」と呼ぶこともある）によって評価することができる。この評価値が所定値以上であれば、重力方向Ｇに近く、上向きの法線であると判断することができるのである。
図の右側に、上向き法線の分布を示した。側面のみに点群が存在するクラスタＣＬ２の結果は、実線で示すように、全点群のうち、評価値｜ｎ・Ｇ｜が比較的小さい値となる点の割合が高くなる。一方、屋上面の点群も含まれるクラスタＣＬ１の結果は、破線で示すように、全点群のうち、評価値｜ｎ・Ｇ｜が比較的小さい値となる点の割合および値１付近の点の割合の２カ所が高くなる。

３次元モデル生成装置２００は、こうして得られた分布に従い、各クラスタに属性を付与する（ステップＳ５４）。即ち、上向き法線の割合が所定の閾値ＴＨＶを超える場合には上面と定義し、そうでない場合には側面と定義するのである。上面というのは、その点群が必ずしも上面のみから構成されるという意味ではなく、上面を含む点群であることを意味する。ステップＳ５３に示した例によれば、クラスタＣＬ１は上面と定義され、クラスタＣＬ２は側面と定義されることになる。定義に用いる閾値ＴＨＶは、上面が存在するか否かの判断基準となる値であり、点群の法線ベクトルに含まれる誤差や、通常の建物の上面の面積に基づいて任意に設定することができる。
ステップＳ５４の図中に、属性の付与結果を示した。ここでは、独立した直方体形状の建物と、高層階のみ２つに分岐している建物を例示している。箱を積み重ねたように描かれているのは、階層への分離結果を模式的に表すものである。図中では、上面をハッチングで、側面を白抜きで示した。図に示すように、独立した建物では、屋上面を含む階層Ｒ１が上面と判断される。分岐した建物では、分岐部分の各屋上面Ｒ２、Ｒ３、および分岐の根元に当たる中間階層Ｒ４が上面と判断される。
さらに、本実施例では、地表面も上面と判断される。３次元点群においては、建物と地表面との分離もなされておらず、地表面は、上向き法線の割合が高い面に相当するからである。

図１２に移り、３次元モデル生成装置２００は、隣接階層のクラスタ間の距離を算出する（ステップＳ５５）。先に付与したクラスタの属性（上面、側面）は、この処理では考慮する必要はない。
図中に処理の様子を示した。ここでは、上階層に属するクラスタＣＬＵと、その下階層に属するクラスタＣＬＬ１、ＣＬＬ２との距離の算出について説明する。人間が見ると、クラスタＣＬＵと、クラスタＣＬＬ２が、同じ建物の上下に位置するクラスタであることは容易に識別できるが、コンピュータが３次元点群の解析をする過程においては、かかる判断を行うためのアルゴリズムが必要である。隣接階層のクラスタ間の距離の算出は、同じ建物の上下に配置するクラスタを特定するために行う処理であるから、上階層に属するクラスタＣＬＵは、その下階層に属する全クラスタＣＬＬ１、ＣＬＬ２を対象に、クラスタ間の距離を算出することになる。
３次元モデル生成装置２００は、まず、処理の対象となる各クラスタＣＬＵ、ＣＬＬ１、ＣＬＬ２の点群を、水平面に投影する。そして、投影面内において、クラスタＣＬＵに属する点とクラスタＣＬＬ１に属する点の間の距離を算出し、その距離の最小値を両クラスタ間の距離とする。つまり、クラスタＣＬＵに属する点（図中の〇）およびクラスタＣＬＬ１に属する点（図中の×）は多数存在するため、それぞれから１点ずつを任意に選出し、点間の距離を算出する。全ての点の組み合わせについて距離を算出し、その中の最小値を求めることにより、クラスタＣＬＵとクラスタＣＬＬ１との距離を求めることができる。図の例では、クラスタＣＬＵ内の点Ｐ１と、クラスタＣＬＬ１の点Ｐ３間の距離Ｄ１１が、両クラスタ間の距離となる。
同様の方法で、３次元モデル生成装置２００は、クラスタＣＬＵと、クラスタＣＬＬ２との間の距離を算出する。図の例では、クラスタＣＬＵ内の点Ｐ１と、クラスタＣＬＬ１の点Ｐ２間の距離Ｄ１２が、両クラスタ間の距離となる。
クラスタＣＬＬ１、ＣＬＬ２は、同階層に属するクラスタであるため、クラスタ間の距離を算出する必要はない。

次に、３次元モデル生成装置２００は、算出されたクラスタ間距離＜Ｔｈｄとなるクラスタ間を、同じ建物で上下に重なっているものと判断し、これらを関連づけする（ステップＳ５６）。
図中に処理の様子を示した。クラスタＣＬＵから見た場合、クラスタＣＬＬ１との間のクラスタ間距離Ｄ１１は閾値Ｔｈｄよりも大きく、クラスタＣＬＬ２との間のクラスタ間距離Ｄ１２は閾値Ｔｈｄよりも小さいとする。この結果、クラスタＣＬＵに対しては、クラスタＣＬＬ２が同じ建物で上下関係にあるものと関連づけられ、クラスタＣＬＬ１は、関連づけられない。
閾値Ｔｈｄは、このように上下関係を判断するための基準となる値である。その値は、建築物の分岐部分の間隔よりも小さい範囲で任意に設定することができる。
実施例では、最も下位層には、地表面が存在しており、建物と地表面との間のクラスタ間距離は小さくなるから、全ての建物は、地表面に関連づけられることになる。

図１３に移り、３次元モデル生成装置２００は、全クラスタについて関連付けを行い、全体の関連付けを表すツリー構造を生成する（ステップＳ５７）。
図中にツリー構造の例を示した。左側に示す建物を例にとり、右側にそのツリー構造を示してある。独立して建築されている建物ＬＶ３については、特に分岐や合流することなく、一本の筋として各クラスタが関連づけられ、地表面ＬＶ５に到達する。
分岐している建物については、高層階の分岐部分ＬＶ１、ＬＶ２から低層階に向けて、それぞれ分岐した状態で関連づけが行われた後、中間の階層ＬＶ４で合流し、その後、低層階では一本の筋として各クラスタが関連づけられ、地表面ＬＶ５に到達する。他の建物についても、同様の方法によりツリー構造が生成される。
図中では、ツリー構造を図示しているが、ステップＳ５７におけるツリー構造生成とは、必ずしもこのような図を描くことを意味するものではなく、全クラスタ間の関連付けを特定することを意味している。

次に、３次元モデル生成装置２００は、ツリー構造を参照し、最下位層のクラスタを地表面と判定する（ステップＳ５８）。３次元点群において地表面が含まれていることが既知である場合には、前処理として、最下位層を地表面と判断する方法をとってもよい。本実施例のように上下関係の関連づけにおいて最下位層を判断する方法をとれば、地表面の存在が確認できていない場合でも、地表面を特定することが可能となる。

次に、３次元モデル生成装置２００は、地表面に隣接するクラスタ単位で点群を分離する（ステップＳ５９）。全建物は最終的には地表面に関連づけられていると考えると、地表面に隣接するクラスタ、即ち、地表面のすぐ上の階層のクラスタは、建物の基礎部分を表すものと考えられる。従って、基礎部分の階層ごとに点群を分離すれば、建物ごとに点群を分離することが可能となる。
図中に処理の様子を示した。最下位層の地表面（図中の破線）を除き、その直上の階層ＣＬＲ１、ＣＬＲ２が建物の基礎部分に相当すると判断する。従って、階層ＣＬＲ１に関連づけられている一覧のクラスタは建物１を表すものと判断される。また、階層ＣＬＲ２に関連づけられている一連のクラスタは、途中で分岐する部分も含めて、建物２を表すものと判断される。
このようにクラスタの関連付けおよびツリー構造を利用することにより、特殊な形状の建物も含め、３次元点群を建物ごとに分離することが可能となる。

Ｄ３．建物枠生成処理：
図１４、１５は、建物枠生成処理のフローチャートである。３次元モデル生成処理（図５）のステップＳ２２における処理に相当する。建物ごとに分離された３次元点群に基づいて、建物の２次元形状を特定する処理である。
３次元モデル生成装置２００は、処理対象となる建造物の３次元点群データを読み込む（ステップＳ７０）。
そして、点群の法線ベクトルに基づいて、上面の点群、側面の点群を分離する（ステップＳ７１）。
図中に処理の様子を示した。上面の点群は、点Ｐ１のように法線ベクトルが上を向いているものとなり、側面の点群は、点Ｐ２のように法線ベクトルがほぼ水平方向を向いているものとなる。上向きの法線ベクトルを有する点群は、点群分離処理と同様の方法で特定することができる（図１１のステップＳ５３参照）。
図中の例では、建造物について、上面１、上面２という２つの高さの異なる上面（図中にハッチングを付した部分）と、その他の側面とが分離される。建物枠生成処理は、通常、建物の２次元形状は上面を外包する形状となっているであろうとの考えに基づいて行われるものである。従って、ステップＳ７０で分離する上面は、図中に示した上面１、上面２のように高さの異なる面であってもよい。

３次元モデル生成装置２００は、分離された側面に沿わせてグリッドを配置する（ステップＳ７２）。
図中に処理の様子を示した。分離された側面には、多数の点群が存在する。これらの点群を水平面に投影すると建造物の平面形状、即ち建物枠に沿って分布する点群が得られる。これらの多数の点群を通過する線分を解析によって求めることも可能ではあるが、本実施例では、より簡易な方法として、グリッドを利用することで、建物枠の決定に際し、考慮すべき点群を減らすものとした。
図の右側にグリッドを配置した例を示した。グリッドは建造物に対して任意の状態で配置することが可能であるが、建物枠の形状を精度良く表すために、グリッドＧＬ１のように建物枠に対して斜めに配置するよりは、グリッドＧＬ２のように建物枠の最も長い辺に沿って配置することが好ましい。
かかる配置を行う方法としては、例えば、建物枠の点群を解析して最も長い辺の方向を求め、グリッドＧＬ２の配置方向を決定する方法をとることができる。
また、別の方法として、次の２段階の方法で求めても良い。第１に、グリッドを仮配置した上で、点群が属している格子を抽出し、この格子がグリッドの縦または横方向に並ぶ最長の長さを評価値として求める。第２に、グリッドの配置方向を回転させながら、同様にして評価値を求め、この評価値が最大となる方向を求める。こうすることによっても、グリッドを建物枠に沿わせて配置することができる。
グリッドの格子サイズは任意に設定可能である。格子サイズが小さければ、建物形状を精度良く表すことができる一方、処理負荷が高くなる。格子サイズが大きければ、建物形状を表す精度が粗くなるが、処理負荷が軽減される。これらを考慮しながら、格子サイズは決定すればよい。

図１５に移り、３次元モデル生成装置２００は、配置されたグリッドに対して上面の３次元点群を投影する（ステップＳ７３）。
図中に点群を投影した状態を示した。点は×で表してある。各グリッド内に１個または複数の点群が存在する格子、全く点が存在しない格子が表れる。
上面に代えて、側面の点群を投影するものとしてもよいが、途中から分岐している建物（図１１のステップＳ５４参照）や、側面よりも上面が張り出している建物では、建物枠を正確に特定できない場合が生じ得るため、上面を用いることが好ましい。

次に、３次元モデル生成装置２００は、点群が存在するグリッドを抽出することで、建物枠を決定する（ステップＳ７４）。
図中に処理の様子を示した。ステップＳ７３に示した点群の投影結果に基づき、１又は複数の点が存在する格子を抽出する。ステップＳ７４内の図では、このような格子にハッチングを付して示した。３次元モデル生成装置２００は、ハッチングを付した各格子の代表点を結ぶことによって建物枠を特定する。本実施例では、格子の重心（図中の●で示した点）を代表点として用いた。
この処理によって、処理対象となるべき点数は、建物枠に対応する位置にある格子点の数にまで減少される。

３次元モデル生成装置２００は、さらに、頂点を間引く処理を行う（ステップＳ７５）。
本実施例では、図中の左側に示すように、ほぼ直線的に配列されていると判断される部分について、中央に位置する頂点を削除する方法を採用した。例えば、頂点ａ、ｂ、ｃが存在する場合、線分ａｂ、ｂｃの外角αを求める。そして、α＜閾値の場合には、これらの３点がほぼ直線的に並んでいるものと判断し、頂点ｂを削除するのである。閾値は、頂点ｂを削除するか否かの判断基準となる値であり、頂点ｂを削除することによって建物枠に有意な形状変化を与えるか否かという観点から任意に設定することができる。
図の右側には、頂点を削除した結果を例示した。この処理によって、ステップＳ７４で示した形状に対し、ステップＳ７５の図中に示すように、頂点Ｐ１〜Ｐ４の４点に削除することができる。
以上の処理によって、３次元モデル生成装置２００は、建物枠を決定することができる。

Ｄ４．位置合わせ処理：
図１６は、位置合わせ処理のフローチャートである。３次元モデル生成処理（図５）のステップＳ２３における処理に相当する。
先に説明した処理によって、３次元点群は建物単位に分割され、しかも建物枠までは求められている。しかし、こうして得られた建物枠は、点群の解析によるものであり、点群の位置誤差および解析誤差などが含まれているため、必ずしも３次元の地図モデルとして使用できるほどの精度は有していない。位置合わせ処理では、２次元地図データを参照することによって、建物枠の形状や位置の精度の向上を図る。
ここで、最初に本実施例での処理の考え方を説明しておく。本実施例における位置合わせ処理は、建物ごとに２次元地図データと対比して、形状や位置の微修正を行うものではない。本実施例では、目標の建物周辺の領域について３次元点群が得られており、これらの３次元点群は、建物ごとに上面、側面を再現するのと同様の精度で、建物間の相対的な位置関係も表しているはずである。かかる考え方の下、本実施例では、これまでの処理で得られている複数の建物枠の全体を，２次元地図に最も精度良く適合させ得る形状および位置の変換方法を見いだすのである。
建物ごとに位置、形状を合わせる方法では、建物の相対的な位置関係について３次元点群が有する情報を無視することになる。これに対し、本実施例の方法によれば、建物ごとの形状、および建物の相対的な位置関係について、３次元点群が有する情報を十分に活用することができるため、精度良く位置合わせを行うことが可能となるのである。以下、その処理方法について説明する。

処理を開始すると、３次元モデル生成装置２００は、３次元点群データおよび２次元地図データを読み込む（ステップＳ８０）。
そして、３次元点群に基づき、側面の端点位置を算出する（ステップＳ８１）。
図中に処理の様子を示した。側面の点群に基づき、側面を表す法線ベクトルｎを求める。法線ベクトルｎは、重力方向解析処理と同様の方法によって求めることができる（図９のステップＳ３８参照）。
また、重力方向ｇと法線ベクトルｎとの外積ｇ×ｎによって、これらに直交する方向、即ち側面の幅方向を求める。こうして求められた幅方向において側面の点群が最も遠方に位置する点Ｗ１、Ｗ２を求め、これを側面の端点位置とするのである。ここで、３次元点群は散在しているため、幅方向の線分上に点Ｗ１、Ｗ２が存在するとは限らない。そこで、側面を構成する点群を、幅方向の線分上に投影した上で、点Ｗ１、Ｗ２を求める方法をとることができる。

次に、３次元モデル生成装置２００は、側面の端点位置を、２次元地図に合わせるよう位置合わせを行う（ステップＳ８２）。
図中に処理の様子を示した。左側に破線で示したのが３次元点群の分布によって得られる建物枠の形状を表し、右側の実線および点が２次元地図上の建物枠およびその頂点を表している。
この例では、３次元点群のうち、２つの建物から求められた側面の端点Ｗ１、Ｗ２、Ｗ３、Ｗ４を２次元地図に合わせる例を示した。３次元モデル生成装置２００は、まず、端点Ｗ１〜Ｗ４について、２次元地図上の頂点と任意に対応づけを行う。３次元点群と２次元地図とを人間が見れば、端点Ｗ１〜Ｗ４と２次元地図上の頂点との対応関係は容易に特定することができるが、コンピュータが処理を行う際には、両者を対応づける情報が何も存在しないため、双方から任意に選択したものを対応づけることになるのである。
例えば、端点Ｗ１、Ｗ２、Ｗ３、Ｗ４を頂点ｆ１、ｆ２、ｆ５、ｆ６にそれぞれ対応づけたものとする。次に、３次元モデル生成装置２００は、この対応づけに適合するような、３次元点群の位置を変換するための変換式を求める。本実施例では、ホモグラフィ変換の変換行列を求めるものとした。そして、かかる変換により３次元点群の形状および位置がどの程度変形するのかを表す評価値を求める。この評価値の算出方法については後述する。
３次元モデル生成装置２００は、端点と頂点との対応関係を変更しながら、上述の変換行列の算出、評価値の算出を繰り返す。そして、最も小さい変形で、３次元点群を２次元地図に適合させられる対応関係を特定するのである。
図の例では、端点Ｗ１、Ｗ２、Ｗ３、Ｗ４を頂点ｆ１、ｆ２、ｆ５、ｆ６に対応づけた場合には、これを実現する変換では、３次元点群の形状、位置が大きく変形する。これに対し、端点Ｗ１、Ｗ２、Ｗ３、Ｗ４を頂点ｆ１、ｆ２、ｆ３、ｆ４に対応づけると、これを実現する変換では、３次元点群の形状、位置の変形は小さい。この結果、端点Ｗ１、Ｗ２、Ｗ３、Ｗ４は頂点ｆ１、ｆ２、ｆ３、ｆ４に対応するものと特定することができる。
かかる考え方により、最適な対応関係を特定する方法を焼きなまし法、またはシミューレーテッドアニーリング法と呼ぶ。対応関係の特定は、他の方法を用いるものとしてもよい。

図１７は、位置合わせ処理における評価値を示す説明図である。
先に説明した通り、位置合わせ処理では、３次元点群を２次元地図に合わせた場合の変形の大きさを評価する評価値を用いることになる。
図１７（ａ）は、正方形ＡＢＣＤを用いた評価値の算出方法を示している。図１６のステップＳ８２において、３次元点群の端点から２次元地図の頂点への変換方法が求められると、この変換方法によって正方形ＡＢＣＤを変換することによって四角形Ａ１、Ｂ１、Ｃ１、Ｄ１が得られる。変換方法による変形が小さいほど、変形後の四角形Ａ１、Ｂ１、Ｃ１、Ｄ１は元の正方形ＡＢＣＤに近い形状をなすはずである。従って、正方形からの形状の変化を、この変換方法の評価値として用いることができる。かかる評価値としては例えば、次の評価値Ｅａ１が考えられる。
評価値Ｅａ１＝｜Ａ１／（π／２）−１｜＋｜Ｂ１／（π／２）−１｜＋｜Ｃ１／（π／２）−１｜＋｜Ｄ１／（π／２）−１｜；
ここで、Ａ１、Ｂ１、Ｃ１、Ｄ１は、各頂角の大きさ（ｒａｄ）を表している。
評価値Ｅａ１は、変形後の四角形Ａ１、Ｂ１、Ｃ１、Ｄ１の頂角の９０度（π／２（ｒａｄ））からのずれに着目した評価値である。この例では、頂角Ａ１、Ｂ１、Ｃ１、Ｄ１の全てを用いて評価値を求めているが、いずれかの頂角で代表させてもよい。

また、別の評価値として、次の評価値Ｅａ２を用いても良い。
評価値Ｅａ２＝｜１−（Ａ１Ｂ１）／（Ｃ１Ｄ１）｜；
ここで、Ａ１Ｂ１は辺Ａ１Ｂ１の長さ、Ｃ１Ｄ１は辺Ｃ１Ｄ１の長さを表している。
評価値Ｅａ２は、変形後の四角形Ａ１、Ｂ１、Ｃ１、Ｄ１の向かい合う辺の長さの変化に着目した評価値である。この例では、辺Ａ１Ｂ１と辺Ｃ１Ｄ１を用いているが、他の対辺、即ち辺Ａ１Ｄ１と辺Ｂ１Ｃ１とを併せて用いるようにしてもよい。

図１７（ｂ）は、３次元点群の端点と、それに対応する２次元地図の頂点との相対的な位置関係に基づく評価値の算出方法を示している。図の左側には、３次元点群における端点Ｗ１〜Ｗ５を示した。これらの端点Ｗ１〜Ｗ５の重心はＧＷと求められ、各端点の相対的な位置関係は、端点から重心に向かう線分のなす角によって表すことができる。例えば、端点Ｗ１、Ｗ２の位置関係は、線分Ｗ１、ＧＷと線分Ｗ２、ＧＷとのなす角θＷ１２で表される。
右側には、２次元地図において端点Ｗ１〜Ｗ５に対応づけられている頂点ｆ１〜ｆ５を示した。これらの頂点ｆ１〜ｆ５の重心はＧｆと求められ、各頂点の相対的な位置関係も端点と同様に定義できる。例えば、頂点ｆ１、ｆ２の位置関係は、線分ｆ１、Ｇｆと線分ｆ２、Ｇｆのなす角θｆ１２で表される。
このとき、変形の程度を表す評価値Ｅｂ１としては、次の値が考えられる。
評価値Ｅｂ１＝（θｆｉｊ−θＷｉｊ）の分散または標準偏差；
評価値Ｅｂ１は、端点と頂点の相対的な位置関係のずれの大きさに着目した評価値である。
また別の評価値として、次の評価値Ｅｂ２を用いても良い。
評価値Ｅｂ２＝｜（ｆｉＧｆ）／（ＷｉＧＷ）｜の分散または標準偏差；
ここで、ｆｉＧｆは線分ｆｉ、Ｇｆの長さを表し、ＷｉＧＷは線分Ｗｉ、ＧＷの長さを表している。
評価値Ｅｂ２は、端点および頂点とそれぞれの重心位置との相対的な長さの変化に着目した評価値である。
上述した評価値Ｅａ１、Ｅａ２、Ｅｂ１、Ｅｂ２は、いずれか一つを選択して用いても良いし、これらを併用してもよい。また、これらを用いた計算式によって新たな評価値を求めるものとしてもよい。

以上の処理により、３次元点群を、２次元地図に整合させることができ、建物を表す点群と地図上の地物との対応関係を特定することができる。また、本実施例の方法によれば、建物ごとの形状、および建物の相対的な位置関係について、３次元点群が有する情報を十分に活用することができるため、精度良く位置合わせを行うことが可能となる。

Ｄ５．ポリゴン生成処理：
図１８は、ポリゴン生成処理例を示す説明図である。３次元モデル生成処理（図５）のステップＳ２４における処理に相当する。
先に説明した処理により、３次元点群は建物ごとに分離されており、また建物枠も求められている。しかも、３次元点群と２次元地図との対応関係も求められているため、建物枠の形状も、対応関係で用いられたホモグラフィ変換により、さらに精度が向上している。
３次元モデル生成装置２００は、こうして得られた情報に基づき、３次元点群から、各建物のポリゴンを生成する。
具体的には、建物枠に基づいて上面を定義し、次にこれを下方に平行移動することで、各辺に対して建物の側面のポリゴンを設定する。特殊な形状の建物（図１３のステップＳ５７参照）についても同様に、分岐している各上面を下方に平行移動することによって分岐部分の側面を生成し、また低階層の上面を下方に平行移動することによって低階層の側面ポリゴンを生成する。このように上面の形状を利用する他、建物枠生成処理で分離された側面の点群（図１４のステップＳ７１参照）に基づいてポリゴンを定義してもよい。
図１８（ａ）には、生成前の３次元点群の様子を示し、図１８（ｂ）には生成されたポリゴンの様子を示した。ポリゴンによって３次元点群が精度良く再現されていることが分かる。

Ｄ６．テクスチャ生成処理：
図１９、２０は、テクスチャ生成処理のフローチャートである。３次元モデル生成処理（図５）のステップＳ２５における処理に相当する。
処理を開始すると、３次元モデル生成装置２００は、ポリゴン生成処理で生成されたポリゴンデータおよび画像データ、撮影条件を入力する（ステップＳ９０）。撮影条件とは、撮影位置、撮影時のカメラの方向、画角などである。
３次元モデル生成装置２００は、次に、ポリゴンが画像内に描かれる画像データを抽出する（ステップＳ９１）。これは、建物に対するテクスチャを切り出すのに適した画像を選択するための処理である。かかる画像を選択する方法としては、撮影した画像の画像解析によって目的となる建物が写っているか否かを判断する方法も可能ではあるが、多大な処理負荷がかかる。そこで、本実施例では、既に建物のポリゴンが得られていることから、このポリゴンを、ステップＳ９０で読み込んだ撮影位置、撮影条件に相当する投影条件で透視投影し、得られた投影画像に基づいて目的の建物が写っているか否かを判断するのである。
図中に処理の様子を例示した。左側の画像１では、撮影位置、撮影条件に相当する投影条件で建物を描いたところ、画像外にはみ出す部分が生じている。右側の画像２では、画像内に建物全体が描かれている。本実施例では、このようにポリゴンを描画した結果に基づき、建物全体が描かれている画像２を、テクスチャの切り出し用に抽出する。

次に、３次元モデル生成装置２００は、抽出された画像のうち、処理対象のポリゴンが背面となっている画像データを除外する（ステップＳ９２）。この処理も、画像解析ではなく、建物のポリゴンを利用して行うものとした。
図中に処理の様子を示した。左側の画像では、処理対象となるべきポリゴンが背面として描かれている。一方、右側の画像では、視認できる範囲に描かれている。従って、左側の画像データは、テクスチャを切り出すのに不適切と判断し、除外されることになる。例えば、処理対象となるべきポリゴンの法線ベクトルと、視線方向を表すベクトルとのなす角度が９０度以下である場合には、目標のポリゴンが背面となると判断することができる。

３次元モデル生成装置２００は、さらに、処理対象のポリゴンが他のポリゴンで遮蔽されている画像データを除外する（ステップＳ９３）。例えば、処理対象となるポリゴンの辺が、他のポリゴンの辺と交差する場合に、遮蔽されていると判断する方法をとることができる。

以上の処理によって、テクスチャを切り出す候補画像としては、画像内に建物が写っているもの、かつ、処理対象のポリゴンが背面とはなっていないもの、かつ、他のポリゴンによって遮蔽されていないものが選択されているはずである。即ち、処理対象となるポリゴンの全体が、隠れたり切れたりすることなく写っている画像が候補画像となる。
３次元モデル生成装置２００は、これらの候補画像データの中から、処理対象ポリゴンの面積が最大となるものを選択する（ステップＳ９４）。
図中に処理の様子を示した。図示するように、候補画像が３通りある場合、処理対象ポリゴン（図中のハッチングを付した部分）を正面から撮影した中央の画像が、面積最大となるため、これを選択することになる。この処理は、処理対象となるポリゴンの面積を算出することによって行うことができる。また、ポリゴンを正面から撮影した画像が面積最大になることが通常であるため、ポリゴンの法線ベクトルと視線方向とのなす角が最小となる画像を候補として選択するようにしてもよい。
候補の選択にあたっては、処理対象となるポリゴンの面積だけでなく、さらに画像データの画質などを考慮するようにしてもよい。

３次元モデル生成装置２００は、こうしてテクスチャを切り出す候補画像を決めると、候補画像から、処理対象ポリゴンに対応する領域を切り出し、テクスチャとして格納する（ステップＳ９５）。
本実施例のテクスチャ生成処理によれば、処理対象となるポリゴンを正面付近から撮影した画像を利用することにより、画質に優れたテクスチャを生成することができる。
また、かかる画像を選択するに当たり、シミュレーション、即ち３次元点群から得られた建物のポリゴンを用いて、各画像データに対応する撮影位置、撮影条件で透視投影を行うことで得られる画像データを利用することにより、画像解析を利用する場合に比較して、軽い処理負荷で精度良く処理を行うことができる。

以上、本発明の実施例について説明した。本発明は、必ずしも上述した実施例の全ての機能を備えている必要はなく、一部のみを実現するようにしてもよい。また、上述した内容に追加の機能を設けてもよい。
本発明は上述の実施例に限定されず、その趣旨を逸脱しない範囲で種々の構成を採ることができることはいうまでもない。例えば、実施例においてハードウェア的に構成されている部分は、ソフトウェア的に構成することもでき、その逆も可能である。

本発明は、３次元点群から３次元モデルを生成するために利用可能である。

１００…撮影管理装置
１０１…カメラ
１０２…ＧＰＳ
１１０…撮影制御部
１１２…撮影シーケンス記憶部
１２０…画像記録部
２００…３次元モデル生成装置
２０１…画像データ記憶部
２０２…３次元点群生成部
２０３…３次元点群データ記憶部
２０４…重力方向解析部
２０５…点群分離処理部
２０６…建物枠生成部
２０７…位置合わせ処理部
２０８…ポリゴン生成部
２０９…テクスチャ生成部
２１０…２次元地図データ
２１１…３次元モデルデータ

Claims (6)

1. 複数の地物の外形を表す複数の点からなる３次元点群をコンピュータによって解析する３次元点群解析方法であって、
   前記３次元点群を構成する３次元点群データは、前記地物の外形を表す複数の点に対し、点ごとに位置座標が格納されており、
   前記３次元点群解析方法は、前記コンピュータが実行するステップとして、
   （ａ） 前記３次元点群データを読み込むステップと、
   （ｂ） 前記３次元点群データを、高さ方向に複数の階層に区分するステップと、
   （ｃ） 前記複数の階層のそれぞれにおいて、前記３次元点群を点間の最短距離が所定値以上開いている部分で複数のクラスタに分離するステップと、
   （ｄ） 前記上下に隣接する階層間で、各階層に属する前記クラスタのうち、相互に少なくとも一部が重なるクラスタ同士を同一地物に属するクラスタとして関連づけることにより、前記３次元点群を地物単位に分離するステップとを備える
   ３次元点群解析方法。
2. 請求項１記載の３次元点群解析方法であって、
   前記ステップ（ｄ）は、階層間で１対多または多対１の関連付けも許容して前記関連付けを行う
   ３次元点群解析方法。
3. 請求項２記載の３次元点群解析方法であって、
   前記３次元点群は、地表面を表す点群も含んでおり、
   前記ステップ（ｄ）は、前記関連付けの結果、所定割合以上のクラスタが集約して関連づけられている低階層のクラスタを地表面と判断し、該クラスタが存在する階層の上に位置する階層におけるクラスタによって前記地物単位への分離を行う
   ３次元点群解析方法。
4. 請求項１〜３いずれか記載の３次元点群解析方法であって、さらに、
   前記ステップ（ｃ）で得られた各階層のクラスタに対し、該クラスタ内に前記地物の上面を構成する点が含まれる場合には地物の上面を表す旨の属性を付すステップを有し、
   前記ステップ（ｄ）における関連付けは、上面の属性が付されたクラスタから、下側の階層に向かって順次行う
   ３次元点群解析方法。
5. 複数の地物の外形を表す複数の点からなる３次元点群を解析する３次元点群解析装置であって、
   前記３次元点群を構成する３次元点群データとして、前記地物の外形を表す複数の点に対し、点ごとに位置座標を格納する３次元点群データ格納部と、
   前記３次元点群を前記地物単位の点群に分離する点群分離処理部とを備え、
   前記点群分離処理部は、
   前記３次元点群データを読み込み、
   前記３次元点群データを、高さ方向に複数の階層に区分し、
   前記複数の階層のそれぞれにおいて、前記３次元点群を点間の最短距離が所定値以上開いている部分で複数のクラスタに分離し、
   前記上下に隣接する階層間で、各階層に属する前記クラスタのうち、相互に少なくとも一部が重なるクラスタ同士を同一地物に属するクラスタとして関連づけることにより、前記３次元点群を地物単位に分離する
   ３次元点群解析装置。
6. 複数の地物の外形を表す複数の点からなる３次元点群を解析するためのコンピュータプログラムであって、
   前記３次元点群を構成する３次元点群データは、前記地物の外形を表す複数の点に対し、点ごとに位置座標が格納されており、
   前記コンピュータプログラムは、
   （ａ） 前記３次元点群データを読み込む機能と、
   （ｂ） 前記３次元点群データを、高さ方向に複数の階層に区分する機能と、
   （ｃ） 前記複数の階層のそれぞれにおいて、前記３次元点群を点間の最短距離が所定値以上開いている部分で複数のクラスタに分離する機能と、
   （ｄ） 前記上下に隣接する階層間で、各階層に属する前記クラスタのうち、相互に少なくとも一部が重なるクラスタ同士を同一地物に属するクラスタとして関連づけることにより、前記３次元点群を地物単位に分離する機能とをコンピュータに実現させるコンピュータプログラム。