Transcript 搬送波位相測定値による 精密測位の理論及び解析処理 Precise

GPS/GNSSシンポジウム2005
搬送波位相測定値による
精密測位の理論及び解析処理
Precise positioning theory and analysis with
carrier-phase measurements
技術コンサルタント 高須 知二
Tomoji TAKASU
内容 Contents
• A.4.1 はじめに
• A.4.2 精密測位のモデル
• A.4.3 パラメータ推定の手法
• A.4.4 相対測位
• A.4.5 精密単独測位(PPP)
• A.4.6 精密解析の実際
• A.4.付 座標系の定義と変換
http://gpspp.sakura.ne.jp
A.4.1 はじめに
精密測位
Precise Positioning
GPS/GNSS搬送波位相測定値に
よる高精度測位技術
http://gpspp.sakura.ne.jp
A.4.2.1 擬似距離
測位コードにより測定された信号伝搬時間
×光速
s
P  c(tr  t )   P
観測時刻
(受信機時計)
送信時刻
(衛星時計)
測位コード(PRNコード)
tr
http://gpspp.sakura.ne.jp
A.4.2.1 擬似距離
測位コードにより測定された信号伝搬時間
×光速
s
P  c(tr  t )   P 
  c(dt  dT)  I  T   P
幾何学距離
電離層
観測誤差
衛星 遅延 対流圏 (10cm
受信機
遅延 ～数m)
時計誤差 時計誤差
http://gpspp.sakura.ne.jp
A.4.2.1 搬送波位相
GPS/GNSS搬送波位相角の連続測定値
L  Φ
搬送波波長(m) 搬送波位相(cycle)
搬送波
tr
http://gpspp.sakura.ne.jp
A.4.2.1 搬送波位相
GPS/GNSS搬送波位相角の連続測定値
L  Φ 
  c(dt  dT)  I  T  N   L
幾何学距離
電離層
衛星 遅延 対流圏
受信機
遅延
時計誤差 時計誤差
観測誤差
(1mm
～数cm)
http://gpspp.sakura.ne.jp
A.4.2.1 搬送波位相
GPS/GNSS搬送波位相角の連続測定値
L  Φ 
  c(dt  dT)  I  T  N   L
N
搬送波位相バイアス
s
 0r 0  n （整数不定性）
http://gpspp.sakura.ne.jp
A.4.2.2 幾何学距離
y
Satellite
Receiver
z
x
r s (t  )
送信時
http://gpspp.sakura.ne.jp
A.4.2.2 幾何学距離
y
Satellite
Receiver
rr (t)

z
x
RZ (E )r s (t  )
E
受信時
http://gpspp.sakura.ne.jp
A.4.2.2 幾何学距離
光差方程式(Light-Time Equation)を解く
  rr (t)  RZ (E )r (t  )
s
観測点位置(受信時)
衛星位置(送信時)
   / c, t  tr  dt
電波伝搬時間
http://gpspp.sakura.ne.jp
A.4.2.2 幾何学距離
光差方程式(Light-Time Equation)を解く
  rr (t)  RZ (E )r (t  )
s
   / c, t  tr  dt
受信時刻
受信機
観測時刻
(受信機時計) 時計誤差
http://gpspp.sakura.ne.jp
A.4.2.2 衛星方位・仰角
• 衛星方位・仰角 Az, El
z (Up) Satellite
r '  Er ( Rz (E )r (t  )  rr (t ))
s
rs'
s
 ( x s , y s , z s )T
Az  ATAN2( x s , y s )
El  arcsin
z
s
Receiver
El
Az
rs'
x (East)
y (North)
http://gpspp.sakura.ne.jp
A.4.2.3 衛星軌道・衛星時計
s
r : 衛星位置
dT : 衛星時計誤差
• 放送暦 (Broadcast Ephemeris)
航法メッセージ中に含まれる
精度 : 軌道～数m, 時計～10nsec
• 精密暦 (Precise Ephemeris)
世界中の観測局データにより決定した
高精度衛星軌道・時計情報
http://gpspp.sakura.ne.jp
A.4.2.3 IGS精密暦
種類
決定/予報値
Final
Rapid Ultra-Rapid
(最終暦) (速報暦) (超速報暦)
決定値
決定値 決定 予報
時間遅れ
13 日
提供頻度
時間 軌道
間隔 時計
軌道
精度
時計
1 週間毎
15 分
5分
<5cm
<0.1ns
17 時間 3 時間 即時
1 日毎
15 分
5分
<5cm
0.1ns
6 時間毎
15 分
15 分
<5cm 10cm
0.2ns 5ns
http://gpspp.sakura.ne.jp
A.4.2.3 精密暦
• 衛星軌道
一定間隔の衛星位置(地球固定座標)
→多次多項式補間(10次以上)
• 衛星時計
一定間隔の衛星時計誤差(バイアス)
→補間誤差、時計の定義に注意
http://gpspp.sakura.ne.jp
A.4.2.4 電離層遅延
40.3
40.310
I  2 Ndl 
TEC
2
f
f

16
Total Electron Content
総電子数 (TECU=1016el/m2)
f
: 搬送波周波数 (Hz)
http://gpspp.sakura.ne.jp
A.4.2.4 電離層フリーLC
二周波搬送波位相による電離層遅延項消去
LC  C1L1  C2 L2 
  c(dt  dT)  T  N LC   LC
2
f1
2
f2
C1  2
,
C

2
2
2
2
f1  f 2
f1  f 2
http://gpspp.sakura.ne.jp
A.4.2.4 電離層モデル
Single Layerモデル
Ionospheric
Pierce Point
z’
Satellite
1 40.31016
I
2
cos z'
f
 TEC(t, pp,  pp )
H
z
Ionosphere
Receiver
Earth
RE
http://gpspp.sakura.ne.jp
GIM(Global Ionos Maps)
100
90
80
70
60
50
40
30
20
10
http://gpspp.sakura.ne.jp
0
A.4.2.5 対流圏遅延
静水圧遅延
湿潤遅延
T  Mdry(El)ZHD Mwet (El)ZWD
天頂静水圧遅延
(2.0～2.4m)
天頂湿潤遅延
(0～0.4m)
M dry, M wet : マッピング関数
http://gpspp.sakura.ne.jp
A.4.2.5 対流圏遅延
• 天頂静水圧遅延(ZHD)
地上気圧(hPa)
0.0022768P0
ZHD 
7
1  0.00266 cos 2  2.8 10 h
緯度
• 天頂湿潤遅延(ZWD)
未知パラメータとして推定
高度(m)
http://gpspp.sakura.ne.jp
A.4.2.5 マッピング関数
NMF(Niell Mapping Function)
1
M ( El) 
仰角 sin El 
a
Zenith
Satellite
b
1
1 c
a
b
sin El 
sin El  c
a,b,c : 係数(緯度、高度、通算日)
El
Receiver
http://gpspp.sakura.ne.jp
A.4.2.6 受信アンテナ位相中心
Satellite
Antenna
Phase Center
Geodetic
Reference
Point
Antenna
Phase Center
Offset
Antenna
Reference Point
(ARP)
http://gpspp.sakura.ne.jp
A.4.2.6 受信アンテナ位相中心
• アンテナ位相中心パラメータ
L1,L2 アンテナ位相中心オフセット
L1,L2 アンテナ位相中心変動(PCV)
• アンテナ機種毎テーブル
地上精密計測データベース
http://gpspp.sakura.ne.jp
A.4.2.7 衛星アンテナ位相中心
• 精密暦 : 衛星重心位置
• 衛星アンテナ位相中心オフセット補正
衛星
アンテナ位相中心(IGS)
x(m)
y(m)
z(m)
Block II/IIA
0.279
0.000
1.023
Block IIR
0.000
0.000
0.000
http://gpspp.sakura.ne.jp
A.4.2.7 衛星固定座標系
Sun
Satellite
y
x
z
Earth
http://gpspp.sakura.ne.jp
A.4.2.8 局位置変動
• 地球潮汐：
月・太陽の引力による弾性変形
固体地球潮汐(Solid Earth Tide)
海洋荷重(Ocean Loading)
極運動潮汐(Pole Tide) etc
• 局位置変動(Site Displacement)
最大上下変動：～20cm
http://gpspp.sakura.ne.jp
A.4.2.8 局位置変動
Up (m)
North (m)
East (m)
TSKB : 2004/1/1 0:00-1/2 0:00
0.2
0.1
0
-0.1
-0.2
0.2
0.1
0
-0.1
-0.2
0.2
0.1
0
-0.1
-0.2
1/1 0:00
1/1 6:00
1/1 12:00
1/1 18:00
1/2 0:00
http://gpspp.sakura.ne.jp
A.4.2.9 Phase Wind-up効果
• GPS搬送波
右旋円偏波
• 衛星ー受信機アンテナの相対的回転
→搬送波位相の進行、後退
http://gpspp.sakura.ne.jp
A.4.2.10 相対論効果
• 特殊相対論効果：衛星時計の遅れ
• 衛星軌道離心率 ≠ 0
→衛星時計誤差の周期変動
→最大で10m以上
• 放送暦・精密暦時計：相対論効果を除去
2r  v
dT  dT0 
2
c
衛星時計誤差
s
s
衛星位置・速度
(慣性座標系)
http://gpspp.sakura.ne.jp
A.4.2.11 マルチパス・観測誤差
• マルチパス
• 電離層遅延・対流圏遅延補正残差
• 受信機雑音
• カットオフ角
• 仰角によるWeighting
http://gpspp.sakura.ne.jp
A.4.2.12 線形結合
記号
線形結合
L1
L2
LC
LG
WL
NL
MW
MP1
MP2
(L1 搬送波位相)
(L2 搬送波位相)
電離層フリー
幾何学フリー
ワイドレーン
ナローレーン
Melbourne-Wübbena
L1 マルチパス
L2 マルチパス
波長 電離 雑音
(cm) 層 (cm)
19.0 1.0 0.3
24.4 1.6 0.3
0.0 0.9
0.6 0.4
86.2 1.3 1.7
10.7 1.3 1.7
86.2 0.0
21
0.0
30
0.0
30
http://gpspp.sakura.ne.jp
A.4.3.1 最小二乗法
観測方程式
z  h( x)  ε
観測量
観測誤差
観測モデル パラメータ
z  (z1, z2 , z3 ,..., zm )T
x  (x1, x2 , x3 ,..., xn )T
http://gpspp.sakura.ne.jp
A.4.3.1 最小二乗法
非線形最小二乗法
初期値 x0 の周りで線形化
z  h( x)  ε
 h( x0 )  H ( x  x0 )  ε
h( x)
H
x
x  x0
 h1 / x1
 h / x
  2 1
 h / x
 m 1
h1 / x2
h2 / x2

hm / x2
n
 h1 / xn 
 h2 / xn 




 hm / xn 
m
http://gpspp.sakura.ne.jp
A.4.3.1 最小二乗法
非線形最小二乗法
初期値 x0 の周りで線形化
z  h( x)  ε
 h( x0 )  H ( x  x0 )  ε
最小二乗推定値
1
xˆ  x0  (H WH) H W (z  h( x0 ))
T
T
http://gpspp.sakura.ne.jp
A.4.3.2 カルマンフィルタ
逐次推定法
推定値
観測データ
z1, z2 ,...
観測更新
xˆ1, xˆ 2 ,...
時間更新
http://gpspp.sakura.ne.jp
A.4.3.2 カルマンフィルタ
観測更新
K k  Pk () H k ( H k Pk () H k  Rk ) 1
xˆ k ()  xˆ k ()  K k ( zk  h( xˆ k ()))
T
zk
T
x̂k
Pk ()  ( I  K k H k ) Pk ()
時間更新
xˆ k 1 ()  xˆ k () 

tk 1
tk
f ( xˆ k (), )d
Pk 1 ()  Φ(tk 1, tk ) Pk ()Φ(tk 1, tk )T  Qk
http://gpspp.sakura.ne.jp
A.4.4 相対測位
• 最も一般的な精密測位
• 基本観測量：搬送波位相二重差
• 基準観測点からの相対位置
(基線ベクトル)
• 干渉測位、基線解析
http://gpspp.sakura.ne.jp
A.4.4.1 二重差
LbA
Satellite A
LaA
Receiver a
AB
Lab
LBa
Satellite B
LBb
Receiver b
A
B
A
B
 (La  La )  (Lb  Lb )
http://gpspp.sakura.ne.jp
A.4.4.1 二重差
AB
Lab
AB
AB
AB
AB
AB
AB
 ab  c (dtab  dTab )  I ab  Tab  Nab
AB
  ab
AB
dtab
 (dtaA  dtaB )  (dtbA  dtbB )  0
AB
dTab
 (dTaA  dTaB )  (dTbA  dTbB )  0
AB
AB
Nab
 (NaA  NaB )  (NbA  NbB )  nab
AB
AB
AB
AB
AB
AB
Lab  ab  I ab  Tab  nab   ab
AB
AB
AB
AB
Lab  ab  nab   ab （短基
線）
http://gpspp.sakura.ne.jp
A.4.4.2 相対測位の原理
Satellite 3
Satellite 2
Satellite 1
Receiver r
Satellite n
Receiver R
http://gpspp.sakura.ne.jp
A.4.4.2 相対測位の原理
未知パラメータ
12
13
1n T
x  (rr T , NrR
, NrR
,..., NrR
)
x0  (r0T ,0,0,...,0)T
観測値
12 13
1n T
zk  (LrR , LrR ,..., LrR )
1
2
  12





/

r



/

r

0

0
rR
r
0
r
0




1
3

  13



/

r



r
0
r / r0 0   0 
rR H 
hk ( x0 )  

 k 

   

  
 1 / r   n / r 0 0   
  1n 
r
0
 r 0

 rR 
http://gpspp.sakura.ne.jp
A.4.4.2 相対測位の原理
全観測値
z
T
T
T T
 (z1 , z2 ,..., zm )
h( x0 )  (h1( x0 ) , h2 ( x0 ) ,...,hm ( x0 ) )
T
T
T T
H  ( H1 , H n ,..., H m )
T
最小二乗推定値
T
1
T T
xˆ  x0  (H WH) H W (z  h( x0 ))
T
T
http://gpspp.sakura.ne.jp
A.4.4.3 整数不定性決定
最小二乗推定値
12 ˆ 13
1n T
ˆ
ˆx  (rˆr T , Nˆ rR
, NrR ,..., NrR )
• 搬送波位相バイアス：実数解
→FLOAT解
• 搬送波位相バイアス：
整数条件を利用した整数化
→FIX解
http://gpspp.sakura.ne.jp
A.4.4.4 相対測位の応用
• RTK(リアルタイムキネマティック)
• 仮想基準点
http://gpspp.sakura.ne.jp
A.4.5 単独測位
軌道/時計
放送暦
～3m/～10 ns
電離層
補正モデル
1～10 m
雑音+マルチパス
コード
0.1～1 m
対流圏
補正モデル
～0.3 m
測位精度
数m～数10m
http://gpspp.sakura.ne.jp
A.4.5 精密単独測位(PPP)
軌道/時計
精密暦
～3cm/～0.1ns
電離層
２周波線形結合
～0.5 cm?
雑音+マルチパス
搬送波位相
0.3～1 cm
対流圏
モデル＋推定
～0.5 cm
測位精度
1cm～数cm
http://gpspp.sakura.ne.jp
A.4.5.1 精密単独測位の原理
Satellite 3
Satellite 2
Satellite n
Satellite 1
L2r
L3r
Lnr
L1r
Receiver r
http://gpspp.sakura.ne.jp
A.4.5.1 精密単独測位の原理
未知パラメータ
1
2
n T
x  (rr , dt, ZWD, NLCr , NLCr ,..., NLCr )
T
初期推定値・共分散行列
xˆ1()  (r0T , dt0 , ZWD0 , NLC1r 0 , NLC2r 0 ,..., NLCnr0 )T
P1 ()  diag( x 2 ,  y 2 ,  z 2 ,  dt2 ,  ZWD2 ,  N 2 ,..., N 2 )
http://gpspp.sakura.ne.jp
A.4.5.1 精密単独測位の原理
観測更新(1)
1
2
n T
zk  (LCr , LCr ,..., LCr )
 1r  c(dtˆ  dT1 )  Tr1  Nˆ LC1 
r 

 r2  c(dtˆ  dT 2 )  Tr2  Nˆ LC 2r 
h( xˆ k ())  




  n  c(dtˆ  dT n )  T n  Nˆ n 
r
LC r 
 r
http://gpspp.sakura.ne.jp
A.4.5.1 精密単独測位の原理
観測更新(2)
 1r / rˆr

 r2 / rˆr
Hk  


  n / rˆ
r
 r
c M wet (Elr1 )
c M wet ( Elr2 )


c M wet (Elrn )
1
0

0
0
1

0




0 
0


1 
2
2
2
Rk  diag(1 , 2 ,..., n )
推定値
xˆ k ()  xˆ k ()  Kk (zk  h( xˆ k ()))
http://gpspp.sakura.ne.jp
A.4.5.1 精密単独測位の原理
時間更新
ˆxk 1()  (rˆr T , dt0 , ZWˆD, Nˆ LC1r , Nˆ LC2r ,..., Nˆ LCnr )T
Pk 1()  Φ(tk 1, tk )Pk ()Φ(tk 1, tk )T  QK
Φ(tk 1, tk )  diag(1,1,1,0,1,1,1,...,1)
Qk  diag(0,0,0,  dt2 , ZWD2t,0,0,...,0)
http://gpspp.sakura.ne.jp
Up (m)
North (m)
East (m)
A.4.5.1 精密単独測位の例
1
TSKB 2004/10/3 0:00-3:00
0
-1
1
0
-1
1
0
-1
10/3 0:00
10/3 1:00
10/3 2:00
10/3 3:00
http://gpspp.sakura.ne.jp
A.4.5.2 精密単独測位の応用
• キネマティックPPP
• GDGPS(Global Differential GPS)
• 精密時刻同期
• GPS可降水量
http://gpspp.sakura.ne.jp
A.4.6 精密解析の実際
• 解析データ
• 解析データ編集
• 単独測位による概略値計算
• 受信機時計飛びの検出・修正
• サイクルスリップの検出・修正
• パラメータ推定
• アウトライアの検出・除去
http://gpspp.sakura.ne.jp
まとめ
• 精密測位のモデル
幾何学距離、電離層遅延、対流圏遅延、
アンテナ位相中心 ...
• パラメータ推定の手法
• 相対測位
• 精密単独測位(PPP)
• 精密解析の実際
http://gpspp.sakura.ne.jp