(57)【要約】【課題】  精度の高い擬似距離を選択して測位計算できるようにする。【解決手段】  測位装置１０は、ＧＰＳ受信機２と赤外線カメラ７と画像処理機８と慣性装置４と観測データフィルタ５と航法計算装置６とを備えて測位計算を行う。観測データフィルタ５は、ＧＰＳ受信機２から各衛星に対する擬似距離と衛星位置とを入力し、さらに、画像処理機８から赤外線カメラ７で撮像した全周囲の建物を仰角・方位角で表した撮像データを取得する。次に、自己位置と衛星位置と建物との位置関係に基づいて各衛星からの測位信号が直接波か反射波か判定する。そして、直接波である測位信号から求めた擬似距離を測位計算に使用する擬似距離として出力する。また、反射波に対して測位計算における重みを設定して出力する。航法計算装置６は、観測データフィルタ５が出力した擬似距離および擬似距離の重みに基づいて測位計算を行う。【選択図】        図２
2007093483.tif 000002
【特許請求の範囲】
【請求項１】
測位衛星から測位信号を受信して測位衛星との擬似距離と測位衛星の位置とを算出する受信機と、
撮影して撮像データを生成する撮像機と、
前記撮像機が生成した撮像データから測位信号を遮蔽・反射する地物を仰角と方位角とで表す撮像データに変換する撮像データ変換部と、
前記撮像機が撮像したときの前記受信機の概算位置を算出する概算部と、
前記受信機が算出した測位衛星の位置と前記概算部が算出した前記受信機の概算位置とに基づいて前記受信機から測位衛星への仰角と所定の方向から測位衛星への方位角とを算出し、前記撮像データ変換部の変換した撮像データが示す測位衛星への方位角に位置する地物の仰角と測位衛星への仰角とを比較し、地物の仰角が測位衛星への仰角より大きい場合に前記受信機の算出した擬似距離は反射波に基づいて算出された擬似距離であると判定する擬似距離判定部と、
前記擬似距離判定部の判定結果に基づいて前記受信機が算出した擬似距離を使用して測位する測位部と
を備えたことを特徴とする測位装置。
【請求項２】
前記撮像データ変換部は、
前記撮像機が生成した撮像データを二値化し、二値化した撮像データから地物のエッジを抽出し、各方位角において抽出したエッジが示す地物の仰角を算出することで、前記撮像機が生成した撮像データから地物を仰角と方位角とで表す撮像データに変換する
ことを特徴とする請求項１記載の測位装置。
【請求項３】
前記概算部は、概算位置と共に、
前記撮像機が撮像したときの前記撮像機のロール角とピッチ角とを算出し、
前記撮像データ変換部は、
抽出したエッジが示す各方位角に位置する地物の仰角を前記概算部が算出したロール角とピッチ角とで補正することで、前記撮像機が生成した撮像データから地物を仰角と方位角とで表す撮像データに変換する
ことを特徴とする請求項２記載の測位装置。
【請求項４】
前記測位部は、
前記擬似距離判定部が反射波に基づいて算出されたと判定しなかった擬似距離を使用して測位する
ことを特徴とする請求項１〜請求項３いずれかに記載の測位装置。
【請求項５】
前記測位部は、重み付けした各擬似距離を使用して重みつき最小二乗法により測位し
重みつき最小二乗法における重みについて、前記擬似距離判定部が反射波に基づいて算出されたと判定した擬似距離の重み付けを前記擬似距離判定部が反射波に基づいて算出されたと判定しなかった擬似距離より軽くする
ことを特徴とする請求項１〜請求項３いずれかに記載の測位装置。
【請求項６】
前記擬似距離判定部は、
前記受信機の算出した擬似距離は反射波に基づいて算出された擬似距離であると判定しなかった場合、前記撮像データ変換部の変換した撮像データの示す測位衛星への方位角と反対の方位角に位置する地物の仰角が前記撮像データ変換部の変換した撮像データの示す測位衛星への方位角に位置する地物の仰角に対する比率を所定の比率と比較し、地物の仰角の比率が所定の比率より大きい場合に前記受信機の算出した擬似距離は直接波と反射波との合成波に基づいて算出された擬似距離であると判定し、地物の仰角の比率が所定の比率以下である場合に前記受信機の算出した擬似距離は直接波に基づいて算出された擬似距離であると判定する
ことを特徴とする請求項１〜請求項３いずれかに記載の測位装置。
【請求項７】
前記擬似距離判定部は、
前記受信機の算出した擬似距離は反射波に基づいて算出された擬似距離であると判定しなかった場合、測位衛星への仰角と所定の仰角とを比較し、測位衛星への仰角が所定の仰角より大きい場合に前記受信機の算出した擬似距離は直接波に基づいて算出された擬似距離であると判定し、測位衛星への仰角が所定の仰角以下である場合に前記撮像データ変換部の変換した撮像データの示す測位衛星への方位角と反対の方位角に位置する地物の仰角が前記撮像データ変換部の変換した撮像データの示す測位衛星への方位角に位置する地物の仰角に対する比率を所定の比率と比較し、地物の仰角の比率が所定の比率より大きい場合に前記受信機の算出した擬似距離は直接波と反射波との合成波に基づいて算出された擬似距離であると判定し、地物の仰角の比率が所定の比率以下である場合に前記受信機の算出した擬似距離は直接波に基づいて算出された擬似距離であると判定する
ことを特徴とする請求項１〜請求項３いずれかに記載の測位装置。
【請求項８】
前記測位部は、
前記擬似距離判定部が直接波に基づいて算出されたと判定した擬似距離を使用して測位する
ことを特徴とする請求項６または請求項７記載の測位装置。
【請求項９】
前記測位部は、重み付けした各擬似距離を使用して重みつき最小二乗法により測位し、
重みつき最小二乗法における重みについて、前記擬似距離判定部が反射波に基づいて算出されたと判定した擬似距離、前記擬似距離判定部が合成波に基づいて算出されたと判定した擬似距離、前記擬似距離判定部が直接波に基づいて算出されたと判定した擬似距離の順で重み付けを軽くする
ことを特徴とする請求項６または請求項７記載の測位装置。
【請求項１０】
受信機が測位衛星から測位信号を受信して測位衛星との擬似距離と測位衛星の位置とを算出する受信処理を行い、
撮像機が撮影して撮像データを生成する撮像処理を行い、
撮像データ変換部が前記撮像機の生成した撮像データから測位信号を遮蔽・反射する地物を仰角と方位角とで表す撮像データに変換する撮像データ変換処理を行い、
概算部が前記撮像機の撮像したときの前記受信機の概算位置を算出する概算処理を行い、
擬似距離判定部が前記受信機の算出した測位衛星の位置と前記概算部の算出した前記受信機の概算位置とに基づいて前記受信機から測位衛星への仰角と所定の方向から測位衛星への方位角とを算出し、前記撮像データ変換部の変換した撮像データが示す測位衛星への方位角に位置する地物の仰角と測位衛星への仰角とを比較し、地物の仰角が測位衛星への仰角より大きい場合に前記受信機の算出した擬似距離は反射波に基づいて算出された擬似距離であると判定する擬似距離判定処理を行い、
測位部が前記擬似距離判定部の判定結果に基づいて前記受信機の算出した擬似距離を使用して測位する測位処理を行う
ことを特徴とする測位方法。
【請求項１１】
請求項１０記載の測位方法をコンピュータに実行させることを特徴とする測位プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、米国のＧＰＳ（Ｇｌｏｂａｌ  Ｐｏｓｉｔｉｏｎｉｎｇ  Ｓｙｓｔｅｍ）をはじめとするＧＮＳＳ（Ｇｌｏｂａｌ  Ｎａｖｉｇａｔｉｏｎ  Ｓａｔｅｌｌｉｔｅ  Ｓｙｓｔｅｍ）測位システムを利用する測位装置、測位方法および測位プログラムに関するものである。
【背景技術】
【０００２】
ＧＰＳ測位は、理想的な環境下ではスペックどおりの精度で測位が可能であるが、マルチパスが多く発生する環境下や、建物による測位信号の遮蔽がおきる場合は測位精度が大きく劣化する。このためＧＰＳ受信機の中にはマルチパス低減技術を実装し、擬似距離の測位精度劣化を一定範囲に抑えるように工夫しているものもある。しかし、ハードウェアの追加により受信機のサイズに制約が生じることやコストがよりかかるという問題が生じる。
【０００３】
ここで、これまでのＧＰＳ受信機に一般的に用いられているマルチパス低減技術について説明する。
マルチパスエラーには疑似距離マルチパス（Ｐｓｅｕｄｏｒａｎｇｅ  ｍｕｌｔｉｐａｔｈ）と搬送波マルチパス（Ｃａｒｒｉｅｒ  ｐｈａｓｅ  ｍｕｌｔｉｐａｔｈ）とがあるが、ここでは疑似距離マルチパスの低減技術について説明する。
【０００４】
良く知られているように、マルチパスによる誤差はＧＰＳ衛星からの直接波の他に、地面や建物からの反射波を同時に受信することによって発生する。マルチパスは常に存在し、観測精度を悪くするだけでなく、それが原因となって測位信号の捕捉を失うことすらある。また、高精度のＤＧＰＳ（Ｄｉｆｆｅｒｅｎｔｉａｌ  ＧＰＳ）受信機ではマルチパスが最も大きな誤差要因となることが多い。ＤＧＰＳについては後述で説明する。
ＧＰＳ受信機が疑似距離を算出する際、ＧＰＳ受信機内部で生成したコード（Ｌ１測位信号のＣ／Ａコードなど）と受信信号の相互相関を取ることによって測位信号の到着時刻を知り、疑似距離に換算する。その際、マルチパスが存在せず直接波だけ受信していれば、その相互相関関数の形状が自己相関関数の形状に一致し、正確にピークの位置を知ることができるため測位信号の正確な時刻同定が可能となる。
ところが、マルチパスが存在する場合、受信する測位信号の相互相関関数の形状は、直接波（図１（ａ））によるものとマルチパス波（図１（ｂ））によるものとの重ねあわせ（図１（ｃ））となる。このため、受信する測位信号の相互相関関数の形状は、直接波（図１（ａ））と比べ形状が崩れてしまいピーク位置もずれてしまう。これにより、そのずれ分が受信時刻誤差となり距離に換算されて疑似距離誤差となる。
【０００５】
ＧＰＳ受信機の多くは、複数の相関器（ｃｏｒｒｅｌａｔｏｒ）を使うＤＬＬ（Ｄｅｌａｙ  Ｌｏｃｋ  Ｌｏｏｐ）を実装してコードタイミングを測定しており、１．０  ｃｈｉｐ  ｅａｒｌｙ−ｌａｔｅ  ｃｏｒｒｅｌａｔｏｒ  ｓｐａｃｉｎｇと呼ばれる技術が一般的に使われていた。
【０００６】
この最も一般的であった技術に対してマルチパス誤差を減らす手法として以下のものが挙げられる。
［１］Ｎａｒｒｏｗ−Ｃｏｒｒｅｌａｔｏｒ  ｔｅｃｈｎｏｌｏｇｙ
［２］Ｌｅａｄｉｎｇ−Ｅｄｇｅ  Ｔｅｃｈｎｉｑｕｅｓ
［３］Ｃｏｒｒｅｌａｔｉｏｎ  Ｆｕｎｃｔｉｏｎ  Ｓｈａｐｅ−Ｂａｓｅｄ  Ｍｅｔｈｏｄｓ（ＥＬＳ  ｍｅｔｈｏｄ，  ＭＥＤＬＬなど）
［４］Ｍｏｄｉｆｉｅｄ  Ｃｏｒｒｅｌａｔｏｒ  Ｒｅｆｅｒｅｎｃｅ  Ｗａｖｅｆｏｒｍｓ（ｓｔｒｏｂｅ  ｃｏｒｒｅｌａｔｏｒ）
［５］ＭＭＴ  Ｔｅｃｈｎｏｌｏｇｙ
【０００７】
例えば［１］Ｎａｒｒｏｗ  Ｃｏｒｒｅｌａｔｏｒは、通常ｅａｒｌｙ  ｃｏｒｒｅｌａｔｏｒとｌａｔｅ  ｃｏｒｒｅｌａｔｏｒのｓｐａｃｉｎｇを１ｃｈｉｐ分にしていたものを、もっと少ない時間差（０．１ｃｈｉｐ分など）にしたもので、非常に簡単な原理で疑似距離誤差を低減できる。ただし、帯域幅が狭いと相関関数の形状が崩れて鈍化することから、Ｎａｒｒｏｗ  Ｃｏｒｒｅｌａｔｏｒを使用するには帯域幅を確保し、それに見合った計算処理能力をＧＰＳ受信機が備える必要がある。
【０００８】
また、ＧＰＳ受信機にマルチパス低減技術が使われるか使われないかにかかわらず、ＧＰＳ受信機で得た擬似距離から精度の悪い擬似距離を特定する技術もある。その一つが、ＲＡＩＭ（Ｒｅｃｅｉｖｅｒ  Ａｕｔｏｎｏｍｏｕｓ  Ｉｎｔｅｇｒｉｔｙ  Ｍｏｎｉｔｏｒｉｎｇ）で使われているアルゴリズムをベースとした手法である。しかし、この手法には５衛星あるいはそれ以上の衛星の観測データがなければ使えないという制約がある。
【０００９】
今後、ＧＰＳ受信機の性能向上のため、マルチパス低減技術の向上のほかに、ＧＰＳ受信機の感度向上と信号補足時間（Ａｃｑｕｉｓｉｔｉｏｎ  ｔｉｍｅ）の短縮が見込まれる。これにより、ＧＰＳ衛星からの直接波が届かない場所でも、反射波だけでその信号を捕捉することが可能になる。このようなＧＰＳ受信機は高感度ＧＰＳ受信機と呼ばれ、すでに実用段階に入っているものも存在する。また、欧州のＧＡＬＩＬＥＯ計画やロシアが運用するＧＬＯＮＡＳＳ（ＧＬＯｂａｌ  ＮＡｖｉｇａｔｉｏｎ  Ｓａｔｅｌｌｉｔｅ  Ｓｙｓｔｅｍ）が計画どおりに執行されれば、２８衛星（ＧＰＳ）＋３０衛星（Ｇａｌｉｌｅｏ）＋２４衛星（ＧＬＯＮＡＳＳ）の合計８２機の測位衛星が利用可能になる。このような状況においては可観測衛星数が大幅に増えることから、精度が悪い観測データを使用しなくても測位計算が行える機会が増えることが予想される。
ＧＰＳ受信機に関する技術として特許文献１・２がある。
【特許文献１】特開２００３−１３９８４３号公報
【特許文献２】特開２００３−００４８３８号公報
【発明の開示】
【発明が解決しようとする課題】
【００１０】
従来の技術では、ＲＡＩＭベースのアルゴリズムを用いて精度の悪い擬似距離を除去する技術のほかに、測位信号のキャリアノイズ比ＣＮ０を基準に、測位に使う信号と使わない信号を区別する技術（特許文献１）、最小二乗法による測位計算にキャリアノイズ比を使って重み付けをする手法などがある。
ここで、ＲＡＩＭのアルゴリズムは先に述べたように５衛星以上の観測データが必要であるので、観測衛星数が４以下の場合は擬似距離精度の悪い衛星を特定できない欠点がある。また、マルチパスが存在しても必ずしも測位信号のキャリアノイズ比が悪くなるとは限らない。よって、キャリアノイズ比が悪ければ擬似距離精度が低下していることが予想されるが、キャリアノイズ比がよいからといって必ずしも擬似距離精度が良いとは限らない。一般にマルチパスの存在による擬似距離精度の低下を、１エポックデータのキャリアノイズ比だけで判断するのは難しい。よって、キャリアノイズ比だけを判断基準とする手法が有効に働かない瞬間もあると予想される。
【００１１】
本発明は、上記の課題などを解決するためになされたもので、キャリアノイズ比によらず、１衛星の観測データに基づいて、その衛星との擬似距離が精度の悪いものか判定して測位計算できるようにする、つまり、精度の高い擬似距離を選択して測位計算できるようにすることを目的とする。
【課題を解決するための手段】
【００１２】
本発明の測位装置は、測位衛星から測位信号を受信して測位衛星との擬似距離と測位衛星の位置とを算出する受信機と、撮影して撮像データを生成する撮像機と、前記撮像機が生成した撮像データから測位信号を遮蔽・反射する地物を仰角と方位角とで表す撮像データに変換する撮像データ変換部と、前記撮像機が撮像したときの前記受信機の概算位置を算出する概算部と、前記受信機が算出した測位衛星の位置と前記概算部が算出した前記受信機の概算位置とに基づいて前記受信機から測位衛星への仰角と所定の方向から測位衛星への方位角とを算出し、前記撮像データ変換部の変換した撮像データが示す測位衛星への方位角に位置する地物の仰角と測位衛星への仰角とを比較し、地物の仰角が測位衛星への仰角より大きい場合に前記受信機の算出した擬似距離は反射波に基づいて算出された擬似距離であると判定する擬似距離判定部と、前記擬似距離判定部の判定結果に基づいて前記受信機が算出した擬似距離を使用して測位する測位部とを備えたことを特徴とする。
【発明の効果】
【００１３】
本発明によれば、撮像データを利用して各擬似距離が反射波に基づいて算出されたかを判定できることにより、精度の高い擬似距離を選択して測位することができる。これにより、都市部などのマルチパス条件下において測位精度の劣化を防ぐことができる。
【発明を実施するための最良の形態】
【００１４】
実施の形態１．
ＧＰＳの測位方式は大きく分類するとコード測位とＧＰＳ測量（干渉測位とも呼ばれる）とに分けられる。コード測位ではＣ／Ａ（Ｃｏａｒｓｅ／Ａｃｑｕｉｓｉｔｉｏｎ）コードなどを利用し、ＧＰＳ測量では搬送波位相積算値を利用する（ＧＰＳではＣ／Ａコードの他にＰ（Ｐｒｅｃｉｓｉｏｎ又はＰｒｏｔｅｃｔｅｄ）コードもあるが秘匿操作によってＹコードになっており民生用には開放されていない）。また、コード測位で使用する擬似距離（Ｐｓｅｕｄｏｒａｎｇｅ）を搬送波位相積算値で平滑化する手法などもあるが、これはコード測位に分類されると考える。実施の形態ではコード測位用のＧＰＳ受信機の使用を前提として説明を行う。
【００１５】
ＧＰＳによるコード測位では複数の衛星から直接波を受信して、直接波の到達時間をもとに各衛星から自分の位置までの擬似距離を算出し、算出した擬似距離をもとに自己位置を計算するのが基本である。この際、地表やビルによって反射した反射波を受信して測位を行ってしまうと擬似距離の精度が劣化し、測位精度が悪くなってしまう。そこで実施の形態１では、大きな誤差を持つ擬似距離データを撮像データを用いて予測する手法を提案する。
【００１６】
擬似距離データの誤差要因には、マルチパス、樹木やポールなどの障害物による回折、散乱などがある。通常は直接波にマルチパス（反射波）が加わって誤差が増大するが、高感度ＧＰＳの場合はマルチパスだけで擬似距離を算出することもあるため、その場合は直接波との距離差分がそのまま測距誤差となる。
直接波とマルチパスを同時に受信するケースでは、マルチパス低減技術で誤差を軽減できる。しかし、前記の通り、必ずしも最適なマルチパス低減技術がＧＰＳ受信機に組み込まれているとは限らない。そのような場合は特に、マルチパスの存在を予測して、観測データの扱いをそれに応じて変えることが望ましい。
【００１７】
図２は、実施の形態１における測位装置１０の構成図である。
実施の形態１における測位装置１０の構成について、図２に基づいて以下に説明する。但し、測位装置１０の構成は複数種類考えられ、図２は測位装置１０の代表的な構成を示す。
【００１８】
測位装置１０は以下のものを備える。
ＧＰＳ衛星が発信したＧＰＳ測位を行うための測位信号を受信するＧＰＳアンテナ１とＧＰＳ受信機２。
ＧＰＳとの複合航法システムを構成するために、ジャイロ・加速度計を有し、ジャイロで測定した姿勢と加速度計で測定した加速度から算出する移動距離とに基づいて、位置（緯度，経度，高度）・姿勢（ロール角，ピッチ角，方位角）を算出する慣性装置４。
全周（３６０度）を撮像して画像データ（以下、天空画像とする）を出力する赤外線カメラ７。
天空画像を仰角（Ｅｌｅｖａｔｉｏｎ）と方位角（Ａｚｉｍｕｔｈ）とを表す天空画像（以下、天空画像投影データとする）に変換する画像処理機８。
ＧＰＳ受信機２が受信した測位信号に基づく擬似距離を測位計算などに使用するかどうかを判断する観測データフィルタ５。
観測データフィルタ５の判断結果に基づいて測位を行う航法計算装置６。
【００１９】
測位装置１０は測位対象である移動体に搭載される。そして、慣性装置４は移動体と姿勢軸が一致し、移動体の位置・姿勢を算出する。また、赤外線カメラ７も移動体と姿勢軸が一致し、撮像した画像は慣性装置４が算出した位置・姿勢に対応する。
【００２０】
ＧＰＳ受信機２は、ＧＰＳアンテナ１からのＲＦ（Ｒａｄｉｏ  Ｆｒｅｑｕｅｎｃｙ）を受信してデジタルデータに変換する回路を備え、そのデジタルデータを使用して疑似距離、疑似距離変化率（または搬送波位相、ドップラー観測値）、測位結果、衛星軌道データ（ＧＰＳ航法メッセージ中のエフェメリスデータあるいはアルマナックデータ）や、衛星軌道データに基づく各ＧＰＳ衛星の位置などの測位関連情報を出力するものである。但し、ＧＰＳ受信機２は、上記ＧＡＬＩＬＥＯやＧＬＯＮＡＳＳで使用する受信機であっても構わない。
【００２１】
ＧＰＳ衛星から発信される測位信号は、電離層や対流圏などの影響を受け伝搬速度が変化し、ＧＰＳ受信機２が算出する擬似距離には誤差が含まれる。そこで、ＧＰＳ受信機２は誤差情報を通信回線経由で取得することで擬似距離の精度を向上して測位精度を向上させることができる。このような誤差情報は、あらかじめＧＰＳアンテナ１の位置が分かっているＧＰＳ受信機が受信した観測データに基づいて作成される。このような誤差情報に基づいてＧＰＳ測位精度を向上するための手法を一般にＤＧＰＳと呼ぶ。
ＤＧＰＳには、測位結果に誤差情報を適用して補正するもの、衛星との擬似距離に誤差情報を適用して補正するものなどがある。本実施の形態では、擬似距離に誤差情報を適用して測位精度を上げるＤＧＰＳを使用する場合を説明する。但し、このＤＧＰＳはあくまで精度向上のための補助的手段であり、擬似距離に誤差情報を適用する以外の方法で誤差情報を利用しても構わないし、ＤＧＰＳを実施しなくてもよい。
【００２２】
赤外線カメラ７は、赤外線用の魚眼レンズや以下に説明する反射型レンズを上方に向け全周を撮像する。
図３は、実施の形態１における反射型レンズを用いた赤外線カメラ７の平面図である。
図４は、実施の形態１における反射型レンズを用いた赤外線カメラ７の正面図である。
図５は、実施の形態１における反射型レンズを用いた赤外線カメラ７の右側面図である。
図６は、実施の形態１における反射型レンズを用いた赤外線カメラ７のＡ−Ａ線断面図である。
図３において、符号「１０３」は“フレーム”、「１１３」、「１１５」は“ネジ”である。
図４において、符号「１０７」は“ミラー”、「１１２」、「１１３」、「１１５」は“ネジ”である。
図５において、符号「１０１」は“カメラヘッド”、「１０３」「１０５」「１０６」は“フレーム”である。
図６において、符号「１０２」は“レンズ”、「１０８」は“ミラー”、「１１２」、「１１７」は“ネジ”である。
反射型レンズを用いた赤外線カメラ７は、図６に示すように、フレーム１０５を透過させて入射した光をミラー１０８で前方に反射し、ミラー１０８で反射した光をミラー１０７で後方に反射する。ミラー１０７で反射した光は、ミラー１０８中央の開口穴を通過して赤外線用のレンズ１０２に入光し、レンズ１０２に入光した光は赤外線カメラ７内部の結像面（図示しない）に集光されて、赤外線カメラ７による撮像が行なわれる。図３〜図６に示すような反射型レンズを用いた赤外線カメラ７は、仰角“１５度〜７０度”程の全周画像を撮像することができる。
【００２３】
建物等によるマルチパスが発生するときには上方に向けたカメラで撮像した天空画像に必ず建物等が写ることに着目した点はポイントの１つである。そこで、処理した天空図と捕捉衛星の可視／不可視を判定することで、不可視衛星を棄却することができ、マルチパス対策が行える。
【００２４】
また、撮像するカメラに赤外線カメラを使用する点もポイントの１つである。
可視カメラを用いた場合、可視光の波長（４００ｎｍ〜８００ｎｍ）が太陽光の中心波長（約５５０ｎｍ）と重なるために、天頂にカメラを向けることによる太陽光外乱による画像飽和やダイナミックレンジの不足が起こり、建物エッジの抽出が困難であるという課題がある。
赤外線の大部分は大気中の水蒸気や二酸化炭素によって吸収されるが、赤外線の領域でも一部は大気をほぼ素通りするものもある。特に中赤帯３〜５μｍ、遠赤帯８〜１２μｍは“大気の窓”と言われ、大気の透過率が極めて高い。このため、この中赤帯または遠赤帯に感度を持つ赤外線カメラを使うことにより、屋外において、太陽光によるカメラの画像飽和や白ボケ／スミアを起すことなく、空の部分と建物や木々との境界エッジ（建物エッジ）を簡単なアルゴリズムの画像処理により抽出が可能になり、可視カメラを用いた場合の課題を解決できる。
【００２５】
航法計算装置６は、ＧＰＳ受信機２から擬似距離、搬送波位相、ドップラー周波数、キャリアノイズ比などの観測生データを取得し、慣性装置４からは加速度・レートなどの慣性データを取得し、航法フィルタを通して現在の自己位置、速度、方向などのデータを予測し出力する。通常、航法フィルタにはカルマンフィルタ（Ｋａｌｍａｎ  ｆｉｌｔｅｒ）が使われ、ＧＰＳの測位結果を利用する疎結合方式（ｌｏｏｓｅ  ｃｏｕｐｌｉｎｇ）と、擬似距離などの観測生データを使用する密結合方式（ｔｉｇｈｔ  ｃｏｕｐｌｉｎｇ）とがある。コード測位（ＧＰＳのＬ１測位信号のＣ／Ａ  ｃｏｄｅによる疑似距離を使って自己位置を算出する測位方式）において、疑似距離を利用して測位計算を行う際、誤差分散値が小さな疑似距離観測値をより信頼して解を算出する方が測位精度を高くするため望ましい。そこで、航法計算装置６は、重みつき最小二乗法により測位計算を行う。
【００２６】
図７は、実施の形態１における観測データフィルタ５の構成図である。
実施の形態１における観測データフィルタ５の構成について図７に基づいて以下に説明する。
観測データフィルタ５は以下のものを備える。
ＧＰＳ受信機２が出力する測位関連情報と画像処理機８が出力する天空画像を入力する入力部５１。
自己位置（赤外線カメラおよびＧＰＳアンテナ１の位置）を推定する自己位置推定部５７。
入力部５１が入力した測位関連情報と天空画像および自己位置推定部５７が推定した自己位置・姿勢に基づいて各測位衛星とＧＰＳアンテナ１との間に測位信号を遮蔽する遮蔽物が有るかの判定を行う遮蔽判定部５３。
入力部５１が入力した測位関連情報と天空画像とに基づいて各測位衛星が発信した測位信号をＧＰＳアンテナ１に反射する反射物が有るかの判定を行う反射判定部５４。
遮蔽判定部５３と反射判定部５４とを備え、各測位衛星との擬似距離がそれぞれ、直接波に基づいて算出された擬似距離か、反射波に基づいて算出された擬似距離か、直接波と反射波との合成波に基づいて算出された擬似距離かを判定する擬似距離判定部５２。
擬似距離判定部５２が反射波または合成波に基づいて算出されたと判定した擬似距離に対して航法計算装置６の測位計算における重要度を示す重みを誤差分散値を使用して設定する重み設定部５５。
擬似距離判定部５２が判定した直接波に基づいて算出されたと判定した擬似距離と反射波に基づいて算出されたと判定した擬似距離と合成波に基づいて算出されたと判定した擬似距離と重み設定部５５が設定した重みとを航法計算装置６に出力する出力部５６。
【００２７】
自己位置推定部５７は、直前の自己位置・姿勢を現在位置・姿勢と推定する方法、カルマンフィルタを使ってダイナミクスに応じた自己位置・姿勢を推定する方法、ジャイロなどを用いたデッドレコニングにより自己位置・姿勢を推定する方法、航法計算装置６から自己位置・姿勢推定値を取得する方法、ＧＰＳ受信機２が測位した位置を用いる方法など、どのような方法で自己位置・姿勢の推定を行っても構わない。
遮蔽判定部５３は、ＧＰＳアンテナ１との間に測位信号を遮蔽する遮蔽物が無い測位衛星との擬似距離を直接波と合成波のいずれかに基づいて算出された擬似距離と判定し、ＧＰＳアンテナ１との間に測位信号を遮蔽する遮蔽物が有る測位衛星との擬似距離を反射波に基づいて算出された擬似距離と判定する。
反射判定部５４は、測位信号をＧＰＳアンテナ１に反射する反射物が無い測位衛星との擬似距離を直接波に基づいて算出された擬似距離と判定し、測位信号をＧＰＳアンテナ１に反射する反射物が有る測位衛星との擬似距離を反射波と合成波のいずれかに基づいて算出された擬似距離と判定する。
【００２８】
図８は、実施の形態１における測位装置１０のハードウェア構成図である。
図８において、測位装置１０は、プログラムを実行するＣＰＵ（Ｃｅｎｔｒａｌ  Ｐｒｏｃｅｓｓｉｎｇ  Ｕｎｉｔ）９１１を備えている。ＣＰＵ９１１は、バス９１２を介してＲＯＭ９１３、ＲＡＭ９１４、通信ボード９１５、磁気ディスク装置９２０と接続されている。
ＲＡＭ９１４は、揮発性メモリの一例である。ＲＯＭ９１３、磁気ディスク装置９２０は、不揮発性メモリの一例である。これらは、記憶装置あるいは記憶部の一例である。
通信ボード９１５は、ＦＡＸ機、電話器、ＬＡＮ、インターネット等に接続されている。
例えば、通信ボード９１５は情報入力部の一例であり、出力部の一例である。
【００２９】
磁気ディスク装置９２０には、オペレーティングシステム（ＯＳ）９２１、プログラム群９２３、ファイル群９２４が記憶されている。プログラム群９２３は、ＣＰＵ９１１、ＯＳ９２１により実行される。
【００３０】
上記プログラム群９２３には、実施の形態の説明において「〜部」、「〜機」、「〜装置」、「〜フィルタ」として説明する機能を実行するプログラムが記憶されている。プログラムは、ＣＰＵ９１１により読み出され実行される。
ファイル群９２４には、以下に述べる実施の形態の説明において、「〜を判定し」、「〜を判定した結果」、「〜を計算し」、「〜を計算した結果」、「〜を処理し」、「〜を処理した結果」のような表現で説明する結果情報が、「〜ファイル」として記憶されている。
また、以下に述べる実施の形態の説明において説明するフローチャートの矢印の部分は主としてデータの入出力を示し、そのデータの入出力のためにデータは、磁気ディスク装置９２０、ＦＤ（Ｆｌｅｘｉｂｌｅ  Ｄｉｓｋ  ｃａｒｔｒｉｄｇｅ）、光ディスク、ＣＤ（コンパクトディスク）、ＭＤ（ミニディスク）、ＤＶＤ（Ｄｉｇｉｔａｌ  Ｖｅｒｓａｔｉｌｅ  Ｄｉｓｋ）、その他の記録媒体に記録される。あるいは、信号線やその他の伝送媒体により伝送される。
【００３１】
また、実施の形態の説明において「〜部」、「〜機」、「〜装置」、「〜フィルタ」として説明するものは、ＲＯＭ９１３に記憶されたファームウェアで実現されていても構わない。或いは、ソフトウェアのみ、或いは、ハードウェアのみ、或いは、ソフトウェアとハードウェアとの組み合わせ、さらには、ファームウェアとの組み合わせで実施されても構わない。
【００３２】
また、実施の形態を実施するプログラムは、磁気ディスク装置９２０、ＦＤ、光ディスク、ＣＤ、ＭＤ、ＤＶＤ、その他の記録媒体による記録装置を用いて記憶されても構わない。
【００３３】
以下に、観測データフィルタ５の処理について説明する。
観測データフィルタ５は、航法計算装置６がＧＰＳ受信機２の受信した観測生データ（擬似距離）を利用するのかどうかを判定する。また、観測生データ（擬似距離）を利用する場合、航法計算装置６の測位計算における重要度の設定である重み付けを行う。
【００３４】
図９は、実施の形態１における観測生データ利用判定パターンを示すフローチャートである。
まず、測位装置１０の利用形態に基づく観測生データの利用を判定する４つのパターンについて図９に基づいて以下に説明する。
【００３５】
判定パターンを判断する条件の一つは、「観測生データ（擬似距離）の数を極力４以下にしない（Ｓ１０１）」かどうかである。
この条件は、「とにかくそのエポックのデータで測位計算をできるようにしたいかどうか」を問うものである。そのエポックの擬似距離を４つ以上揃えて測位計算しなくても、デッドレコニング（推測航法）やカルマンフィルタを用いて位置を推測することができるため判断条件の一つとなる。
例えば、疎結合ＧＰＳ／ＩＮＳ複合航法を使う場合では、ＧＰＳデータの誤差分散値が比較的小さくかつ正確に見積もれる状況であれば「ＹＥＳ」を選択して４つ以上の擬似距離に基づいて測位すれば良いし、誤差分散値が大きいと同時に不正確であれば「ＮＯ」を選択してＩＮＳのデータに基づいて位置を推測した方が良いこともある。
【００３６】
判定パターンを判断する条件のもう一つは、「マルチパス判定をする（Ｓ１０２、Ｓ１０５）」かどうかである。
受信した測位信号が直接波であるか反射波であるかを判定する方法には、衛星と測位装置１０との間に建物などの遮蔽物が存在する場合に測位装置１０の受信したその衛星の測位信号は直接波ではなく反射波であると判定する方法（遮蔽判定）がある。
但し、遮蔽物が存在しない場合に測位装置１０の受信した測位信号には直接波の他に反射波（マルチパス）を含むこともある。この反射波（マルチパス）を含むかを判定するのがマルチパス判定である。
つまり、遮蔽判定のみを行うか、遮蔽判定とマルチパス判定とを共に行うかを判断する条件である。
【００３７】
上記２つの判断条件に基づく４つの判定パターンについて以下に説明する。
図１０は、実施の形態１における判定パターン１の処理を示すフローチャートである。
まず、「擬似距離の数を極力４以下にしない」且つ「マルチパス判定をしない」場合（Ｓ１０３）の擬似距離の利用判定処理を図１０に基づいて以下に説明する。
【００３８】
観測データフィルタ５は、以下の処理を行い、航法計算装置６がどの衛星の擬似距離を利用して測位処理を行うかを判定する。
【００３９】
まず、入力部５１は、ＧＰＳ受信機２が記憶部に記憶した観測生データ（擬似距離、衛星の位置など）、衛星軌道データ、受信時刻などを取得する（Ｓ２０１）。
次に、入力部５１は、画像処理機８が赤外線カメラ７の撮像した天空画像から遮蔽物・反射物を仰角と方位角とで表した天空画像投影データを取得する（Ｓ２０２）。
【００４０】
次に、自己位置推定部５７は自己位置の推定を行う（Ｓ２０３）。
【００４１】
次に、擬似距離判定部５２は、前記観測生データから観測衛星情報を抽出し衛星の選択を行う（Ｓ２０４）。
以下、擬似距離判定部５２は、Ｓ２０５〜Ｓ２０８の処理を全衛星に対して行う。
【００４２】
次に、遮蔽判定部５３は、Ｓ２０１で取得した衛星の位置、Ｓ２０３で推定した自己位置、Ｓ２０２で取得した天空画像投影データを用いて遮蔽判定を行う（Ｓ２０５）。
【００４３】
次に、擬似距離判定部５２は、Ｓ２０５で遮蔽判定した結果に基づいて処理を分岐する（Ｓ２０６）。
擬似距離判定部５２は、遮蔽物が有ると判定した場合は、擬似距離の利用有無を示すフラグＦｕに、擬似距離を利用しないことを示す「０」を設定する（Ｓ２０７）。
遮蔽物が有る場合は、擬似距離が直接波ではなく反射波に基づくデータであると判断できるためである。つまり、擬似距離の精度が悪いため航法計算装置６の測位処理に利用しない。
また、擬似距離判定部５２は、遮蔽物が無いと判定した場合は、フラグＦｕに擬似距離を利用することを示す「１」を設定する（Ｓ２０８）。
【００４４】
ここで、擬似距離判定部５２は、全衛星に対する処理が終わっていなければＳ２０４に処理を分岐してＳ２０５〜Ｓ２０８の処理を繰り返す。また、全衛星に対して処理済みである場合はＳ２１０に処理を分岐する（Ｓ２０９）。
【００４５】
全衛星に対するＳ２０５〜Ｓ２０８の処理を終了後、重み設定部５５は、擬似距離を利用することを示す「１」を設定したフラグＦｕの個数を判定する（Ｓ２１０）。
「１」を設定したフラグＦｕの個数が４以上である場合、出力部５６はフラグＦｕに「１」を設定した擬似距離を航法計算装置６に出力して処理を終了する。
【００４６】
「１」を設定したフラグＦｕの個数が４未満である場合、重み設定部５５は、まず「０」を設定したフラグＦｕに「１」を設定する（Ｓ２１１）。
そして、重み設定部５５は、Ｓ２１１の処理でフラグＦｕに「１」を設定した擬似距離の重みを算出するための誤差分散σＧＩＳを計算する（Ｓ２１２）。
Ｓ２１１の処理でフラグＦｕに「１」を設定した擬似距離は精度が悪いので、擬似距離の誤差分散にσＧＩＳを加えることで、Ｓ２０８の処理でフラグＦｕに「１」を設定した精度の良い擬似距離に対して、測位処理における重みを低くすることができる。
但し、Ｓ２１１において、「０」を設定したフラグＦｕの一部（例えば必要数分）に「１」をするようにしてもよい。
誤差分散σＧＩＳを計算後、出力部５６はフラグＦｕに「１」を設定した擬似距離と誤差分散σＧＩＳを航法計算装置６に出力して処理を終了する。
【００４７】
図１１は、実施の形態１における判定パターン２の処理を示すフローチャートである。
「擬似距離の数を極力４以下にしない」且つ「マルチパス判定をする」場合（Ｓ１０４）の擬似距離の利用判定処理を図１１に基づいて以下に説明する。
【００４８】
観測データフィルタ５は、以下の処理を行い、航法計算装置６がどの衛星の擬似距離を利用して測位処理を行うかを判定する。
【００４９】
図１１においてＳ３０１〜Ｓ３０５は判定パターン１（図１０）のＳ２０１〜Ｓ２０５と同じ処理である。
【００５０】
次に、Ｓ３０６〜Ｓ３１７の処理を説明する。
遮蔽判定（Ｓ３０５）後、まず、反射判定部５４は、衛星の位置、自己位置・姿勢、Ｓ２０２で取得した天空画像投影データを用いてマルチパス判定を行う（Ｓ３０６）。
次に、擬似距離判定部５２は、遮蔽判定（Ｓ３０５）とマルチパス判定（Ｓ３０６）の結果に基づいて処理を分岐する（Ｓ３０７、Ｓ３０８、Ｓ３１１）。
【００５１】
擬似距離判定部５２は、「遮蔽物が有り且つマルチパスが有る」と判定した場合、この衛星との擬似距離は反射波に基づくデータであると判断できる為、フラグＦｕに擬似距離を利用しないことを示す「０」を設定する。また、遮蔽判定とマルチパス判定の結果を示すフラグＦｋに「遮蔽物が有り且つマルチパスが有る」ことを示す「１」を設定する（Ｓ３０９）。
【００５２】
擬似距離判定部５２は、「遮蔽物が有り且つマルチパスが無い」と判定した場合、直接波、反射波共に受信できず擬似距離は算出されない為、フラグＦｕに「０」を設定する。また、フラグＦｋに「遮蔽物が有り且つマルチパスが無い」ことを示す「３」を設定する（Ｓ３１０）。
【００５３】
擬似距離判定部５２は、「遮蔽物が無く且つマルチパスが有る」と判定した場合、この衛星との擬似距離は直接波の他に反射波にも基づくデータであると判断できる為、フラグＦｕに「０」を設定する。また、フラグＦｋに「遮蔽物が無く且つマルチパスが有る」ことを示す「４」を設定する（Ｓ３１２）。
【００５４】
擬似距離判定部５２は、「遮蔽物が無く且つマルチパスが無い」と判定した場合、この衛星との擬似距離は直接波に基づくデータであると判断できる為、フラグＦｕに擬似距離を利用することを示す「１」を設定する。また、フラグＦｋに「遮蔽物が無く且つマルチパスが無い」ことを示す「２」を設定する（Ｓ３１３）。
【００５５】
ここで、擬似距離判定部５２は、全衛星に対する処理が終わっていなければＳ３０４に処理を分岐してＳ３０５〜Ｓ３１３の処理を繰り返す。また、全衛星に対して処理済みである場合はＳ３１５に処理を分岐する（Ｓ３１４）。
【００５６】
全衛星に対するＳ３０５〜Ｓ３１３の処理を終了後、重み設定部５５は、擬似距離を利用することを示す「１」を設定したフラグＦｕの個数を判定する（Ｓ３１５）。
「１」を設定したフラグＦｕの個数が４以上である場合、出力部５６はフラグＦｕに「１」を設定した擬似距離を航法計算装置６に出力して処理を終了する。
【００５７】
「１」を設定したフラグＦｕの個数が４未満であれば、重み設定部５５は、「遮蔽物が有り且つマルチパスが有る」ことを示す「１」または「遮蔽物が無く且つマルチパスが有る」ことを示す「４」を設定したフラグＦｋに対応するフラグＦｕに「１」を設定する（Ｓ３１６）。
そして、重み設定部５５は、Ｓ３１６の処理でフラグＦｕに「１」を設定した擬似距離の重みを算出するための誤差分散σＧＩＳを計算する（Ｓ３１７）。
Ｓ３１６の処理でフラグＦｕに「１」を設定した擬似距離は精度が悪いので、擬似距離の誤差分散にσＧＩＳを加えることで、Ｓ３１３の処理でフラグＦｕに「１」を設定した精度の良い擬似距離に対して、測位処理における重みを低くすることができる。
但し、Ｓ３１６において、「１」または「４」を設定したフラグＦｋに対応するフラグＦｕの一部（例えば、フラグＦｋに「４」を設定したフラグＦｕを優先した必要数分）に「１」を設定するようにしてもよい。
誤差分散σＧＩＳを計算後、出力部５６は、フラグＦｕに１を設定した擬似距離と誤差分散σＧＩＳを航法計算装置６に出力して処理を終了する。
【００５８】
次に、図９に示す判定パターン３（Ｓ１０７：「擬似距離の数が４以下でも構わない」且つ「マルチパス判定をしない」）の擬似距離の利用判定処理について説明する。
判定パターン３では、図１０に示す判定パターン１のフローチャートにおいてＳ２１０〜Ｓ２１２の処理（利用する擬似距離の個数が４未満の場合に利用しないと判定した擬似距離を利用するようにフラグを設定する）を実行しない。その他の処理は図１０に示す判定パターン１と同じである。
【００５９】
次に、図９に示す判定パターン４（Ｓ１０６：「観測生データの数が４以下でも構わない」且つ「マルチパス判定をする」）の擬似距離の利用判定処理について説明する。
判定パターン４では、図１１に示す判定パターン２のフローチャートにおいてＳ３１５〜Ｓ３１７の処理（利用する擬似距離の個数が４未満の場合に利用しないと判定した擬似距離を利用するようにフラグを設定する）を実行しない。その他の処理は図１１に示す判定パターン２と同じである。
【００６０】
判定パターン３や判定パターン４の使用は、具体的には密結合ＧＰＳ／ＩＮＳ複合航法フィルタを使う場合を想定している。この場合には、計算に使える擬似距離の数（観測衛星数）が４を下回っても、航法フィルタ（カルマンフィルタ）は観測更新を行えるため、ＩＮＳによる累積誤差を小さくすることができる。そのため、疎結合方式のＧＰＳ／ＩＮＳよりも精度が良くなるのが普通である。従って密結合方式の場合は、航法フィルタが測位結果を利用する疎結合方式とは違い、建物等に遮蔽されていると判断された衛星を無理に使用する必要がない。
【００６１】
以下に、天空画像取得処理（Ｓ２０２・Ｓ３０２）の際の画像処理機８の“天空画像変換処理”と、遮蔽判定計算処理（Ｓ２０５・Ｓ３０５）と、マルチパス判定計算処理（Ｓ３０６）と、誤差分散σＧＩＳ計算処理（Ｓ２１２・Ｓ３１７）の詳細について説明する。
【００６２】
まず、画像処理機８が赤外線カメラ７の撮像した天空画像から遮蔽物・反射物の仰角と方位角とを表す天空画像に変換する天空画像変換処理（Ｓ２０２・Ｓ３０２）の流れについて、図１２に基づいて以下に説明する。図１２は、実施の形態１における画像処理機８が実行する天空画像変換処理を示すフローチャートである。
【００６３】
まず、画像処理機８は赤外線カメラ７の出力した天空画像を入力する。天空画像が表す画像の例を図１３に示す。図１３に示すように、赤外線カメラで撮像した画像は空が黒く建物が明るく写る（Ｓ４０１）。
次に、画像処理機８は入力した天空画像を所定の明度を閾値として二値化する。図１３に示した天空画像を二値化した図を図１４に示す。図１４に示すように、画像処理機８は、天空画像の二値化により、赤外線カメラ７の撮像した画像を建物（および樹木など他の地物）と建物以外とに区別する（Ｓ４０２）。
次に、画像処理機８は二値化した天空画像から建物のエッジ（遮蔽物・反射物となる地物の輪郭）を抽出する。図１４に示した二値化した天空画像から建物のエッジを抽出した図を図１５に示す。このとき、赤外線カメラ７の撮像した画像の縁も抽出される（Ｓ４０３）。
次に、画像処理機８は抽出したエッジの外周（画像の縁）を削除する。図１５に示した建物シルエット以外の外周を削除した図を図１６に示す。図１６に示すように、画像処理機８は、エッジの外周を削除することにより、画像の縁を除いて建物のエッジを抽出する（Ｓ４０４）。
次に、画像処理機８はエッジの外周を削除して残った建物シルエットそれぞれのピクセルの仰角・方位角を算出し、仰角・方位角面に投影した天空画像投影データを生成する。図１６に示した建物のエッジを仰角・方位角面に投影した天空画像投影データが表すグラフを図１７に示す。図１７のグラフは、方位角を「０〜３６０度」とし、分解能を「１／２度」とした、測位点の各方位角分解能に位置する建物のエッジの赤外線カメラ７からの仰角を表す（Ｓ４０５）。
次に、画像処理機８は、赤外線カメラ７が撮像時の赤外線カメラ７の姿勢角（ロール角・ピッチ角・方位角）を慣性装置４から取得し、姿勢角に基づいて天空画像投影データを補正する（Ｓ４０６）。
そして、画像処理機８は生成した天空画像投影データを出力する（Ｓ４０７）。
【００６４】
ここで、画像処理機８がＳ４０５において実行する天空画像投影データの生成処理について、図１８に基づいて以下に説明する。Ｓ４０５において行う仰角・方位角平面への投影は、外周を削除したエッジ画像を定量的に取り扱うため、建物エッジのｘ、ｙ座標を仰角、方位角に変換する。図１８は、実施の形態１における天空画像の関係図である。
図１８に示すように、天空画像の左上を原点としたＸ−Ｙ軸において、建物のエッジを示す各点を（ｘ，ｙ）、天空画像の中心で表される天頂座標を（ｘｃ，ｙｃ）、天頂座標（ｘｃ，ｙｃ）から建物のエッジを示す各点（ｘ，ｙ）までの長さで表される像高を「ｈ」とする。また、天空画像を撮像した赤外線カメラの焦点距離を「ｆ」とする。この関係において以下の式１〜式４により、画像処理機８は各方位角θにおける天頂角φ（＝９０°−“仰角”）を算出して天空画像投影データを生成する。以下の式１〜式４は、撮像機が正射影方式である場合の一例であり、画像処理機８は撮像機の射影方式（立体射影・等距離射影・等立体角射影など）に応じた式を用いて天空画像投影データを生成する。
【００６５】
【数１】
2007093483.tif 000003
【００６６】
次に、画像処理機８が、Ｓ４０６において実行する、天空画像投影データを姿勢角で補正する処理について、図１９に基づいて以下に説明する。赤外線カメラ７の仰角・方位角平面の方位０度方向と慣性装置４の姿勢角検出軸は一致しているものとする。図１９は、実施の形態１における天空画像投影データに対する姿勢角補正処理を示すフローチャートである。
【００６７】
画像処理機８は、慣性装置４から取得したピッチ角φｐｉｔｃｈ分の仰角補正角φｐｉｔｃｈ’を各方位角θ（０〜３６０度）について式５により算出する。そして、仰角補正角φｐｉｔｃｈ’を加算して天空画像投影データの示す各方位角θにおける仰角を補正する。ピッチ角φｐｉｔｃｈは赤外線カメラ７の傾きを示す（Ｓ５０１）。
【００６８】
また、画像処理機８は、慣性装置４から取得したロール角φｒｏｌｌ分の仰角補正角φｒｏｌｌ’を各方位角θ（０〜３６０度）について式６により算出する。そして、仰角補正角φｒｏｌｌ’を加算して天空画像投影データの示す各方位角θにおける仰角を補正する。ロール角φｒｏｌｌは赤外線カメラ７の傾きを示す（Ｓ５０２）。
【００６９】
画像処理機８は、慣性装置４から取得した方位角θｈｅａｄｉｎｇを減算して天空画像投影データの示す各方位角θ（０〜３６０度）を補正する。方位角θｈｅａｄｉｎｇは赤外線カメラの方位を示す（Ｓ５０３）。
【００７０】
【数２】
2007093483.tif 000004
【００７１】
姿勢角補正処理で行う“ピッチ角分の補正”“ロール角分の補正”“方位角分の補正”は、上記の順序に限らず、任意の順序で構わない。
【００７２】
上記の姿勢角補正処理により、画像処理機８は赤外線カメラ７の姿勢角（オイラー角）に対応した天空画像投影データを生成することができる。例えば、真北を方位角０度、ロール角・ピッチ角が０度である水平面を仰角０度とする天空画像投影データを生成することができる。
このように、測位装置１０の移動体搭載を前提にして、車体の姿勢角（ロール角，ピッチ角，方位角）を慣性センサ（ジャイロ，加速度計）を赤外画像と同期して計測し、これらの情報により車体動揺によらず光軸が重力方向と常に平行に天頂を向くように画像を切り出せるようにした点はポイントの１つである。
【００７３】
投影した天空画像投影データを生成する。図１７に示した天空画像投影データが表すグラフに対する姿勢角補正処理後のグラフを図２０に示す。
【００７４】
上記のように、画像処理機８は、カメラで撮像した天空画像を、カメラの射影方式に応じて、カメラ座標（カメラの姿勢軸を基準とする座標）からＧＰＳ座標系と同じ真北を方位角零度、水平面を基準とする極座標系の仰角・方位角平面に変換できる。
そして、擬似距離判定部５２は、画像処理機８が天空画像を変換して生成した生成した仰角・方位角平面（天空画像投影データ）と衛星軌道データから計算される衛星の存在する方位角・仰角とを利用して、その衛星の可視・不可視を判断する。
【００７５】
次に、遮蔽判定部５３が行う遮蔽判定計算処理（Ｓ２０５・Ｓ３０５）について説明する。このとき、遮蔽判定部５３は、自己位置と衛星の位置とに基づいて、ＧＰＳアンテナ１を基準位置とする衛星の存在する方位角（ＧＰＳアンテナ１から衛星への方位角）と仰角（ＧＰＳアンテナ１から衛星への仰角）とを算出する。そして、遮蔽判定部５３は、天空画像投影データが示す衛星の存在する方位角に位置する遮蔽物の仰角と衛星の存在する仰角とを比較し、遮蔽物の仰角が衛星の存在する仰角より大きい角度であれば、当該衛星は不可視であり、測位信号の遮蔽有りと判定する。つまり、図２０に示した仰角・方位角平面グラフに、自己位置から衛星への方位角と仰角をプロットした場合に衛星仰角がグラフ線よりも上方であれば移動体から見てその衛星は可視であり、衛星仰角が低ければその衛星は不可視となる。
【００７６】
次に、反射判定部５４が衛星の位置、自己位置、天空画像投影データを用いて反射波の有無を判定するマルチパス判定計算処理（Ｓ３０６）の流れについて以下に説明する。
遮蔽判定部５３は仰角・方位角平面と衛星の仰角の関係だけを使って、衛星の建物遮断による可視・不可視判断を行ったが、反射判定部５４は仰角・方位角平面図上で建物面と衛星の位置関係から反射波発生の可能性判定を行う。
ここで、反射判定部５４は、遮蔽判定部５３により、測位信号の遮蔽無しと判定された衛星について以下のマルチパス判定計算処理を行う。また、反射判定部５４は、遮蔽判定部５３により、測位信号の遮蔽有りと判定された衛星について、測位信号のマルチパス有りと判定する。
【００７７】
図２１は、実施の形態１における直接波と反射波の混入を示す図である。
マルチパス判定計算処理（Ｓ３０６）では、図２１に示すように、反射を起こす建物面（反射面）は地面から垂直で曲面を持たない比較的単純な面を想定し、反射波の建物面への入射角と反射角はほぼ等しいものとする。また、衛星とＧＰＳアンテナ１との距離は十分に離れているため直接波と建物面に入射前の反射波とは平行に進行するものとする。
【００７８】
図２２は、実施の形態１における反射判定部５４が実行するマルチパス判定計算処理を示すフローチャートである。
【００７９】
まず、「反射波の建物面への入射角と反射角はほぼ等しい」という前提から、高い仰角に存在する衛星からの測位信号は入射角・反射角が浅いため反射波としてＧＰＳアンテナ１には届かないと考え、反射波の届かない衛星の仰角を判定する。そこで、反射判定部５４は、所定の仰角閾値φｔｈｒｅｄと、自己位置と衛星の位置とに基づいて算出したＧＰＳアンテナ１を基準位置とする衛星の存在する仰角とを比較し、衛星の存在する仰角が所定の仰角閾値φｔｈｒｅｄ以下である衛星を抽出する。また、反射判定部５４は、衛星の存在する仰角が所定の仰角閾値φｔｈｒｅｄより大きい衛星について、測位信号のマルチパス無しと判定する。例えば、仰角閾値φｔｈｒｅｄを４０度または２０度〜６０度の間の値にする（Ｓ６０１）。
【００８０】
反射判定部５４は、マルチパス無しと判定した場合、マルチパス判定計算処理を終了し、
マルチパス無しと判定しなかった場合、Ｓ６０３〜Ｓ６０４によりマルチパス判定を行う（Ｓ６０２）。
【００８１】
マルチパス無しと判定しなかった場合は、まず、抽出した衛星からの反射波が到来する方位角を算出する。ここで、信号強度の強い反射波は直接波と反対の方位角（衛星の存在する方位角からおよそ１８０度の方向）から到来すると想定する（これ以外の反射ケースは信号強度が下がり、一般的な受信機の信号強度棄却処理で対処可能と考える）が、反射を起こす建物面がいくらかの曲面であったり直接波の伝搬方向に対するいくらかの角度を持っていたりする場合を考えて、反射波が到来する方位角を算出する。そこで、反射判定部５４は、自己位置と衛星の位置とに基づいて算出したＧＰＳアンテナ１を基準位置とする衛星の存在する方位角θｃａｎｄｉｄａｔｅ１［ｄｅｇ］に対する反射波の到来する方位角θｃａｎｄｉｄａｔｅ２を以下の式により算出する（Ｓ６０３）。
以下の式においてθｏｆｆｓｅｔ（ｊ）、ｊ＝１〜ｎ、但しｎはθｏｆｆｓｅｔの個数。このθｏｆｆｓｅｔ（ｊ）は反射を起こす建物面が、いくらかの曲面である場合や、直接波の伝搬方向に対するいくらかの角度を持つ場合を考慮したパラメータである。例えば、ｎ＝４の場合、θｏｆｆｓｅｔ（１）＝５［ｄｅｇ］、θｏｆｆｓｅｔ（２）＝−５［ｄｅｇ］、θｏｆｆｓｅｔ（３）＝１０［ｄｅｇ］、θｏｆｆｓｅｔ（４）＝−１０［ｄｅｇ］とする。
【００８２】
θｃａｎｄｉｄａｔｅ２＝θｃａｎｄｉｄａｔｅ１＋１８０＋θｏｆｆｓｅｔ（ｊ）            式７
【００８３】
次に、図２３に示すような関係にある、衛星の存在する方位角θｃａｎｄｉｄａｔｅ１に位置する遮蔽物の仰角と、反射波の到来する方位角θｃａｎｄｉｄａｔｅ２に位置する反射物の仰角とを天空画像投影データから取得し、反射物の仰角判定を行う。そこで、反射判定部４５は、以下の条件式８を満たす場合、当該衛星について測位信号のマルチパス有りと判定し、以下の条件式８を満たさない場合、当該衛星について測位信号のマルチパス無しと判定する（Ｓ６０４）。
条件式８において、天空画像投影データが示す方位角θにおいて遮蔽物・反射物の仰角をφ（θ）で表す。また、Ｐｔｈｒｅａｄは所定の比率閾値である。Ｐｔｈｒｅａｄは、道路の中央にＧＰＳアンテナ１が位置する場合に、例えば“１”とする。また、複数車線を走行する自動車などを移動体として測位する場合、ＧＰＳアンテナ１は道路の中央に位置しないため、Ｐｔｈｒｅａｄは、“道路幅”“通行車線（右側通行または左側通行）”“車線数”などを考慮して設定する。例えば、ＧＰＳアンテナ１が道路中央から道路脇に偏よったところに位置する場合、Ｐｔｈｒｅａｄは道路中央のときの値よりも大きい値を設定し、判定するのが良い（例えば、“１”よりも大きい値）。また、車線数が少ない例えば“６ｍ”程度の幅の“２車線”道路であれば、道路中央からの偏りが少ないため、Ｐｔｈｒｅａｄを“１”として判定して構わない。当然の如く、測定誤差を考慮して適宜適切な値を設定しても良い。
また、Ｓ６０３において複数の方位角θｃａｎｄｉｄａｔｅ２を算出した場合、任意の条件によりマルチパスの有無を判定する。例えば、複数の方位角θｃａｎｄｉｄａｔｅ２のうち、１つでも条件式８を満たす場合にマルチパス有りと判定し、１つも条件式８を満たさない場合にマルチパス無しと判定する。
【００８４】
【数３】
2007093483.tif 000005
【００８５】
そして、反射判定部４５は、全てのθｏｆｆｓｅｔ（ｊ）についてＳ６０３とＳ６０４を実施済みか判定する。判定の結果、全てのθｏｆｆｓｅｔ（ｊ）についてＳ６０３とＳ６０４を実施済みの場合はマルチパス判定計算処理を終了し、実施済みでない場合はＳ６０３に処理を移行する（Ｓ６０５）。
【００８６】
図２１に示すように、建物に遮蔽はされていないが、衛星仰角が低く、直接波の入射方向と反対側の方向に高い建物がある場合において、直接波と建物に１回反射した反射波とを同時に計測する場合が都市部環境においては多く発生する。そこで、反射判定部５４が上記のマルチパス判定計算処理を行うことにより、都市部の高いビルの間を走行しながら測位する場合に発生しやすい事象に対応してマルチパス判定することができる。
【００８７】
次に、Ｓ２１２で行う誤差分散σＧＩＳの計算処理について説明する。
擬似距離の精度に対する重みを測位結果に反映するには、重みつき最小二乗法によって測位計算を行う。この重みを下げるには誤差分散値を大きくすればよい。航法計算装置６が測位処理に使用する重みつき最小二乗法（後述で説明する）における重みは、重み行列Ｃ（Ｗｅｉｇｈｔ  ｍａｔｒｉｘ）によって決定される。重み行列は疑似距離観測値誤差を衛星毎に見積もることによって得られる。疑似距離は衛星と観測者間の測位信号伝播時間を計測することによって得られるが、伝搬遅延などの計測誤差要因は複数あり、全計測誤差要因は式９で表される。
【００８８】
【数４】
2007093483.tif 000006
【００８９】
従って、全計測誤差の分散σＤ２は式１０で計算される。
【００９０】
【数５】
2007093483.tif 000007
【００９１】
ＧＰＳ受信機のキャリアノイズ比ＣＮ０（Ｃａｒｒｉｅｒ  ｎｏｉｓｅ  ｒａｔｉｏ）を利用して重み行列を計算する手法があるが、これはトラッキングエラーの分散をＧＰＳ受信機の設計内容およびＣＮ０に応じて計算して、全遅延量の分散σＤ２に反映させる方法である。
【００９２】
これに対し実施の形態１では、カメラによる撮像データを使用して判定した遮蔽・マルチパスに対する誤差分散値σＧＩＳ２を新たに追加する。遮蔽・マルチパスに対する誤差分散値σＧＩＳ２を追加後の全遅延量の分散σＤ２を式１１に示す。
【００９３】
【数６】
2007093483.tif 000008
【００９４】
σＧＩＳ２は各衛星の視線方向（ＬＯＳ：Ｌｉｎｅ  Ｏｆ  Ｓｉｇｈｔ）に建造物などの障害物があると判断された場合は所与の値を設定する。例えば式１２の計算によりσＧＩＳ２の値を設定する。
【００９５】
【数７】
2007093483.tif 000009
【００９６】
σＧ２はあらかじめ決めた一定値で、Ｆｕはその観測生データ（擬似距離）を測位に利用するかどうかを示すフラグである。σＧ２の値は、例えば以下のように設定することができる。
まず、直接波と反射波の両方を受信していると考えられる場合は、測距誤差ＵＥＲＥ（Ｕｓｅｒ  Ｅｑｕｉｖａｌｅｎｔ  Ｒａｎｇｅ  Ｅｒｒｏｒ）が５メートル程度であることを考慮してσＧ＝５（すなわち、σＧ２＝２５）とする。
また、反射波しか受信していないと考えられる場合は、擬似距離誤差が数十メートルから数百メートルになることを考慮してσＧ＝１００（すなわち、σＧ２＝１００００）とする。
【００９７】
図２４は、実施の形態１における判定パターン１および３を示す図である。
上記の通り、判定パターン１、判定パターン３はマルチパス判定を行わず、遮蔽判定で擬似距離を利用するか判定する。
図２４において、遮蔽がある場合は、受信波が反射波であること、擬似距離を利用しない（判定パターン３）こと、擬似距離を利用する場合は誤差分散σＧＩＳに１００（ｍ）を設定する（判定パターン１）ことを示している（図２４（１））。
また、遮蔽がない場合は、受信波が直接波および反射波であること、擬似距離を利用する（判定パターン１および３）ことを示している（図２４（２））。
【００９８】
図２５は、実施の形態１における判定パターン２および４を示す図である。
上記の通り、判定パターン２、判定パターン４は遮蔽判定とマルチパス判定とで擬似距離を利用するか判定する。
図２５において、遮蔽がありマルチパスがある場合は、受信波が反射波であること、擬似距離を利用しない（判定パターン４）こと、擬似距離を利用する場合は誤差分散σＧＩＳに１００（ｍ）を設定する（判定パターン２）ことを示している（図２５（１））。
また、遮蔽がなくマルチパスがない場合は、受信波が直接波であること、擬似距離を利用する（判定パターン２および４）ことを示している（図２５（２））。
また、遮蔽がありマルチパスがない場合は、受信波が無いこと、擬似距離を利用できない（判定パターン２および４）ことを示している（図２５（３））。
また、遮蔽がなくマルチパスがある場合は、受信波が直接波および反射波であること、擬似距離を利用しない（判定パターン４）こと、擬似距離を利用する場合は誤差分散σＧＩＳに５（ｍ）を設定する（判定パターン２）ことを示している（図２５（４））。
【００９９】
次に、観測データフィルタ５が利用すると判定した擬似距離と誤差分散σＧＩＳを使用して算出した重みとに基づいて航法計算装置６が行う重みつき最小二乗法による測位処理について説明する。
【０１００】
まず、疑似距離の観測値は式１３のように表される。
【０１０１】
【数８】
2007093483.tif 000010
【０１０２】
ここで、ｎ衛星の観測値が得られている時、式１３を点（Ｘ０，Ｙ０，Ｚ０）の周りでテイラー展開して式１４を得る。
【０１０３】
【数９】
2007093483.tif 000011
【０１０４】
ｈｏｔ’ｓは繰り返し計算によって小さくすることができるので無視するものとする。また誤差の項Ｘを疑似距離観測値に含めて考えると、疑似距離残差をδρとして式１５を得る。ここで式１６〜式１８を使用するものとする。
【０１０５】
【数１０】
2007093483.tif 000012
【０１０６】
式１５を解く際、δρｉの誤差分散値σｉ２を見積もることができていれば、δρの疑似距離観測値誤差の共分散行列Ｐが求められる。通常、それぞれの疑似距離観測値誤差の間には相関がないものとして、式１９のように共分散行列Ｐを決めるからである。
【０１０７】
【数１１】
2007093483.tif 000013
【０１０８】
共分散行列Ｐを用いて重み行列Ｃを式２０で表す。
【０１０９】
【数１２】
2007093483.tif 000014
【０１１０】
そして、式１５を変形して式２１により最小二乗解を得る。
【０１１１】
【数１３】
2007093483.tif 000015
【０１１２】
σ０２は分散ファクター（Ｖａｒｉａｎｃｅ  ｆａｃｔｏｒ，  ｔｈｅ  ｖａｒｉａｎｃｅ  ｏｆ  ｕｎｉｔ  ｗｅｉｇｈｔ）である。
なお、この計算は最初に（Ｘ０，Ｙ０，Ｚ０）を適当に仮定して計算し、
【０１１３】
【数１４】
2007093483.tif 000016
【０１１４】
のｃΔｔｒ以外の要素（（Ｘ−Ｘ０）など）の値が十分小さな値に収束するまで繰り返し計算する。
分散ファクターの値は適当に決めて計算するが、その値が適切な値であれば、
分散ファクターの見積もり値は１に近くなる。分散ファクターの見積もり値は式２２、式２３で見積もられる。
【０１１５】
【数１５】
2007093483.tif 000017
【０１１６】
【数１６】
2007093483.tif 000018
【０１１７】
ただし行列Ｈがｍ×ｎ行列の場合とする。（ｍ−ｎ）は方程式の自由度となる。
【０１１８】
【数１７】
2007093483.tif 000019
【０１１９】
を使えば、解
【０１２０】
【数１８】
2007093483.tif 000020
【０１２１】
の分散値は式２４のように見積もられる。
【０１２２】
【数１９】
2007093483.tif 000021
【０１２３】
以上から、重みつき最小二乗法では式２１、式２４によって解とその分散値を得る。
ここで、観測データフィルタ５が算出したσＧＩＳに基づく分散σＤ２（式１１）は航法計算装置６で算出するものとする。
【０１２４】
実施の形態１において、遮蔽判定または遮蔽判定とマルチパス判定により測位計算に使用する擬似距離の判定を行ったが、遮蔽判定を行わずマルチパス判定で測位計算に使用する擬似距離の判定を行っても構わない。
【０１２５】
実施の形態１では、以下のようなことを説明した。
例えば、高感度ＧＰＳ受信機においてマルチパスのみで測距を行っている場合は、擬似距離誤差は大きくなることから測位精度は悪くなる。この高感度ＧＰＳ受信機は、測位結果を全く得られないよりは精度が悪くても自己位置を知ることができたほうが良い場面で使用される。
ここで、観測によって得た複数の擬似距離の中には精度のよいものもあるため、マルチパスのみで測距を行っている精度の悪いデータなのかそうではないのかを判断する基準が必要となる。
最も代表的な判断基準はキャリアノイズ比であるが、実施の形態１では撮像データを利用した判断基準を以下のように提供した。
撮像データを利用して、直接波の経路、すなわち衛星への視線方向（ＬＯＳ）に建物がないかどうかを判断する。また、マルチパスを発生させる建物がないかどうかを判断する。もし、建物が存在すると判断される場合は、その衛星の観測データは反射波によるものである可能性が高いと判断する。
つまり、撮像データを利用して精度の低い擬似距離を判断できることが、実施の形態１のポイントの１つである。
【０１２６】
上記の実施の形態１によれば、天空画像の撮影に天空に向けた魚眼レンズまたは反射光学系を備えた超広角の赤外線カメラを用いることにより、赤外線透過率の高い空と建物等その他の部分とを容易に区別することができ、これにより太陽光による外乱に影響を受けず、天候に左右されない簡易なアルゴリズムで建物エッジの抽出が可能となる。また、カメラの射影方式によらず、建物エッジを自己位置を中心とする極座標系の仰角・方位角グラフに変換できる。また、本センサ類を搭載する車体のローリング、ピッチング、ヨーイングといった動揺によらず、建物エッジをＧＰＳ座標系と同じ真北を方位角零度、水平面を基準とする極座標系の仰角・方位角グラフに変換でき、その位置におけるＧＰＳ衛星の存在する方位角（Ａｚｉｍｕｔｈ）と仰角（Ｅｌｅｖａｔｉｏｎ）の方位を同一グラフ上でプロットすることで衛星の可視・不可視判定が容易に可能である。また、本マルチパス判定方式によるれば、衛星からの反射波到達点の方位とその方位に位置する建物高さから反射波混入の可能性を判定し、適宜その衛星の擬似距離の重みを下げることで、マルチパスの影響を軽減することができる。
【０１２７】
そして、マルチパス（反射波）に基づく擬似距離を測位計算に使用する場合に、誤差分散値を大きく設定して、重みつき最小二乗法による測位計算を行った。これにより、より精度の高い測位計算を行うことができる。これも実施の形態１のポイントの１つである。
【０１２８】
実施の形態１で示した処理方法は、マルチパスである反射波を判断し、測位計算に使用しない又は測位計算における重みを下げることができることから、マルチパス除去技術ということもできる。
【図面の簡単な説明】
【０１２９】
【図１】測位信号の相互相関関数を示す図。
【図２】実施の形態１における測位装置１０の構成図。
【図３】実施の形態１における反射型レンズを用いた赤外線カメラ７の平面図。
【図４】実施の形態１における反射型レンズを用いた赤外線カメラ７の正面図。
【図５】実施の形態１における反射型レンズを用いた赤外線カメラ７の右側面図。
【図６】実施の形態１における反射型レンズを用いた赤外線カメラ７のＡ−Ａ線断面図。
【図７】実施の形態１における観測データフィルタ５の構成図。
【図８】実施の形態１における測位装置１０のハードウェア構成図。
【図９】実施の形態１における観測生データ利用判定パターンを示すフローチャート。
【図１０】実施の形態１における判定パターン１の処理を示すフローチャート。
【図１１】実施の形態１における判定パターン２の処理を示すフローチャート。
【図１２】実施の形態１における画像処理機８が実行する天空画像変換処理を示すフローチャート。
【図１３】実施の形態１における天空画像が表す画像の例。
【図１４】実施の形態１における二値化画像の例。
【図１５】実施の形態１におけるエッジ抽出画像の例。
【図１６】実施の形態１における外周削除画像の例。
【図１７】実施の形態１における天空画像投影データが表す仰角・方位角平面グラフの例。
【図１８】実施の形態１における天空画像の関係図。
【図１９】実施の形態１における天空画像投影データに対する姿勢角補正処理を示すフローチャート。
【図２０】実施の形態１における天空画像投影データが表す仰角・方位角平面グラフの例。
【図２１】実施の形態１における直接波と反射波の混入を示す図。
【図２２】実施の形態１における反射判定部５４が実行するマルチパス判定計算処理を示すフローチャート。
【図２３】実施の形態１における直接波と反射波の方位角を示す仰角・方位角平面グラフの例。
【図２４】実施の形態１における判定パターン１および３を示す図。
【図２５】実施の形態１における判定パターン２および４を示す図。
【符号の説明】
【０１３０】
１  ＧＰＳアンテナ、２  ＧＰＳ受信機、４  慣性装置、５  観測データフィルタ、６  航法計算装置、７  赤外線カメラ、８  画像処理機、１０  測位装置、５１  入力部、５２  擬似距離判定部、５３  遮蔽判定部、５４  反射判定部、５５  重み設定部、５６  出力部、５７  自己位置推定部、１０１  カメラヘッド、１０２  レンズ、１０３，１０５，１０６  フレーム、１０７，１０８  ミラー、１１２，１１３，１１５，１１７  ネジ、９１１  ＣＰＵ、９１２  バス、９１３  ＲＯＭ、９１４  ＲＡＭ、９１５  通信ボード、９２０  磁気ディスク装置、９２１  ＯＳ、９２３  プログラム群、９２４  ファイル群。
【図１】
2007093483.tif 000022
【図２】
2007093483.tif 000023
【図３】
2007093483.tif 000024
【図４】
2007093483.tif 000025
【図５】
2007093483.tif 000026
【図６】
2007093483.tif 000027
【図７】
2007093483.tif 000028
【図８】
2007093483.tif 000029
【図９】
2007093483.tif 000030
【図１０】
2007093483.tif 000031
【図１１】
2007093483.tif 000032
【図１２】
2007093483.tif 000033
【図１４】
2007093483.tif 000034
【図１５】
2007093483.tif 000035
【図１６】
2007093483.tif 000036
【図１７】
2007093483.tif 000037
【図１８】
2007093483.tif 000038
【図１９】
2007093483.tif 000039
【図２０】
2007093483.tif 000040
【図２１】
2007093483.tif 000041
【図２２】
2007093483.tif 000042
【図２３】
2007093483.tif 000043
【図２４】
2007093483.tif 000044
【図２５】
2007093483.tif 000045
【図１３】
2007093483.tif 000046