第2章位置の計算

Size: px

Start display at page:

Download "第2章位置の計算"

あおいおなか
6 years ago
Views:

1 1 第章位置の計算測量において地球上の水平位置は経度と緯度平面座標又は極座標の距離と方向角 ( 球面距離と方位角 ) で表される以下では経緯度と平面座標の関係等について計算方法について例題 ( プログラム ) を用いて示すことにしよう以下の計算での数値は Bessel 楕円体から GRS80 楕円体に変更してあるので注意のこと.1 経緯度および方位角図.1.1 に示すように与点 Aの緯度と経度 1 (φ, λ) は既知とするいま点 Aと新点 Bとの球面距離 (S) と点 Aにおける点 Bの方位角 (α) を測定したとすれば点 Bの緯度と経度 (φ, λ ) と点 BにおけるAの方位角 ( 反方位角 ) は次に示す式で求められる φ = φ + W3 1 e u 3W4 e (1 e ) u sin φ cos φ W4 tan φ (1 e ) v W5 e (1 e ) 3 [1 sin φ e sin φ(5 6 sin φ)]u 3 W5 [ (1 e ) tan φ e tan φ(13 10 sin φ)]uv (.1.1) λ λ cosφ = Wv + W uv + W tan φ + e 1 e cos φ u v 1 3 (W3 tan φ)v 3 (.1.) Np dλ 90-φ 90-φ B α S A α φ =φ+δφ-δφ φ λ 赤道 λ 図.1.1 緯度経度と方位角 α = α + (Wtanφ)v + W (1 + tan φ + e 1 e cos φ)uv W3 tanφ tan φ + e 1 e cos φ (1 e ) cos4 φ u v 4e 1 6 W3 tanφ 1 + tan φ + 1 e cos φ v o e (.1.3) 1latitude( 緯度 ) longitude( 経度 ): 楕円体上の地点の位置 azimuthal angle( 方位角 ): 楕円体上ある位置の子午線の北を基準に他の位置まで右回りに測った角

2 u = S a cosα, v = S a sinα W = 1 e sin φ e = a b a [ 例題.1] 新点 Bの緯度経度 (φ B, λ B ) および反方位角 (α ) をプログラム (SOK1.BAS) を用いて求めてみようただし与点 Aの緯度経度 (φ A, λ A ) AB 間の球面距離 (S) および点 Aにおける点 Bの方位角 (αは次の値とする φa =34 o , λa =135 o S=14, m, α=134 o ( 解 ) 緯度 B = 34º 経度 L = 135º 方位角 α = 314º 点の経緯度より球面距離および方位角の計算点の緯度経度が与えられている場合前節の式.1.1~ 式.1.3 を用い繰り返し反復計算によって (u, v) を求めるとこれらの点間の球面距離 (S) および方位角 ( α) は次式から解ける S = a u + v α = tan 1 v u (..1) ( 第 1 近似 ) u f 1 v g 1 ( 第近似 ) u f 1 + f u 1 + f 3 v 1 v g 1 g u v 1 ( 第 3 近似 ) u f 1 + f u + f 3 v + f 4 u 3 + f 5 u v v g 1 g uv g 3 u v + g 4 v 3...(..) ここに f 1 = (φ φ) 1 e W 3 f = 3We (1 e ) sinφcosφ f 3 = Wtanφ f 4 = W e (1 e ) [1 sin φ e sin φ(5 6sin φ)]

3 3 f 5 = W 6(1 e ) [1 + 3tan φ e tan φ(13 10sin φ)] g 1 = (λ λ) cosφ W g = Wtanφ g 3 = W 3 (1 + 3tan φ + e 1 e cos φ) g 4 = W 3 tan φ (..3) [ 例題.] 点 A,Bの緯度経度が与えられているものとするこれらの値からプログラム (SOK.BAS) によって球面距離 (S) と点 Aにおける点 Bの方位角 (α) を求めてみようただし点 Aの緯度経度 (φ A, λ A ) および点 Bの緯度経度 (φ B, λ B ) は次のとおりとする φa=34º φb=34º λa=135º λb=135º ( 解 ) 方位 α =134º 球面距離 S= 14, m.3 点の平面座標による球面距離および方位角の計算点 A,Bの平面直角座標 1 はそれぞれA(x A, y A ), B(x B, y A ) とすれば AB 間の平面距離 (s) および点 Aにおける点 Bの方向角 (t) は次式から計算できる ( 平面距離 ) s = (x B x A ) + (y B y A ) (.3.1) ( 方向角 ) t = tan 1 y B y A x B x A ( 球面距離 ) s ( s S ) (.3.) (.3.3) ( 縮尺係数 ) 縮尺係数はその点の縮率であり日本の場合楕円体に横円筒をかぶせたときの ds/ds は次式で求められる y m = m o 1 + M 1 N 1 m + y 4 o 4M 1 N 1 m 4 o...(.3.4) 平面直角座標系の場合原点の縮尺係数 m 0 = UTM の場合 m 0 = とする 1plane rectangular coordinates: 地球楕円体に横円筒をかぶせて楕円体の中心を投影原点とし地点を円筒面に投影しこの円筒面を切り開いて平面にした座標 direction angle: この円筒面の真上である座標の北を基準として他の点まで右回りに測った角

4 4 ( 距離の補正 )s/s は平面距離 s を求めたり球面距離 S を求めたりする係数である S s = m o 6r m m (y 1 + y 1 y + y ) + o η t 6r 3 m o 3 (x x 1 )(y y 1 ) = m o 8r m m (y 1 + y ) 1 o 4r m m (y y 1 ) + o 6r 3 3 m (x x 1 )(y y 1 ) o η t ( 方位角 ) (.3.5) α = t + γ + T 1 t 1 γ:p における子午線収差 T 1 t 1 は方向補正 α = t + γ + T t (.3.6) ( 方向補正 ) T 1 t 1 = x x 1 (y 6r 1 + y ) η t 1 9r 3 (x x 1 ) (y y 1 ) + η t 6r 3 (y y 1 )(3y 1 + y 1 y + y ) T t = x x 1 (y 6r 1 + y ) + η t 1 9r 3 (x x 1 ) (y y 1 ) η t 6r 3 (y y 1 )(y 1 + y 1 y + 3y ) R = MN (.3.7) M = a(1 e ) W 3 N = a W W = 1 e sin φ e = a b a [ 例題.3] 平面直角座標系における点 A,Bの値が与えられているものとするこれらの値を用いてプログラム (SO K3.BAS) により平面距離 (s) 球面距離(S) 方向角(t) および点 A における点 B の方位角 ( α) を求めてみようただし点 A(-144, m,107, m), 緯度 φa=34º 点 B(-155,04.18m,118, m) および点 Aにおける子午線収差は γ=+0º である ( 解 ) 平面距離 s= 15, 方向角 t= 球面距離 S= 14, m 方位角 α=

5 5.4 経緯度から平面座標の計算地球の基準楕円体 1 に横円筒をかぶせて各地点をその円筒上に投影しておきその円筒を切り開けば楕円体から平面への写像が直接的に得ることができるこれは横メルカトール図法と呼ばれている測量で多く用いられる横メルカトール図法としての平面直角座標系では投影原点の縮率を mo= として距離の精度が 1/10,000 以内に収まるように配慮されている地点の緯度経度 (φ, λ) から平面座標 (x,y) および子午線収差 (γ) を求める式は次のとおりである x=m o [ φ + λ 4 Δλ Nsinφcosφ + Nsinφ 4 cos3 φ(5 t + 9η + 4η 4 ) + λ6 70 Nsinφcos5 φ(61 58t + t η 330t η )] (.4.1) y = m o ΔλNcosφ + Δλ3 N 6 cos3 φ(1 t + η + λ 5 10 Ncos5 φ(5 18t + t η 58t η ] (.4.) γ = Δλsinφ + Δλ3 3 sinφcos φ(1 + 3η + η 4 ) + Δλ 5 15 sinφcos4 φ( t ) (.4.3) Np x Δλ=λ-λo 原子午線子午線 y γ P φo φ Δφ δφ x y φ λ 赤道 λo 平行圏図.4.1 緯度経度と平面座標 Δφ = a(1 e )[k 1 φ φ o 1 k sinφ φ o + 1 k 4 3 sin4φ 4φ o 1 6 k 4 sin6φ 6φ o k 5 sin8φ 8φ o 1 10 k 6 sin10φ 10φ o + ] (.4.4) k1= k= k3= k4= reference ellipsoid: 準拠楕円体 transverse mercator projection: 横メルカトル投影

6 6 k5= k6= Δλ=λ-λ o λ o = 座標系原点の経度 t=tanφ η=e cosφ e = a b b [ 例題.4] 点 A の緯度経度は (φ=34º , λ=135º ) とするこの点の平面直角座標 (x, y) の値および子午線収差 (γ) をプログラム (SOK4.BAS) を用いて求めてみようただし平面座標系は第 5 系 (φo=36º,λo=134º0 ) とする ( 解 ) = Y= γ = 平面座標から経緯度への変換ある地点における平面直角座標 (x,y) から緯度経度 (φ, λ) および子午線収差 (γ) を求めるには次の式で行える ) MN φ = φ 1 (y m o λ = λ o + γ = y mo y mo 3 Ncosφ y mo N t y mo 3 3N 3 t + y mo 4 t(5 + 4MN 3 3t + η 9t η ) y mo 6 t( MN 5 90t + 45t 4 )...(.5.1) y mo 5 (1 + 6N 3 cosφ t + η ) + (5 + 10N 5 cosφ 8t + 4t 4 )...(.5.) φ 1 = φ1 (i 1) S (i 1) D M (i 1) (1 + t η ) + y mo 5 t( + 5t + 3t 4 )...(.5.3) 15N 5 D = a(1 e ) k 1 φ o k sinφ o + k 3 4 sin4φ o k 4 6 sin6φ o + k 5 8 sin8φ o k 6 10 sin10φ x o + + m o S = a(1 e )(k 1 φ 1 k sinφ 1 M = a(1 e ) (1 e sin φ 1 ) 3/ + k 3sin4φ 1 4 k 4sin6φ k 5sin8φ 1 8 k 6sin10φ ) [ 例題.5] 点 A の平面直角座標 xy の値は第 5 系 (-144, m,107, m) とするこの点の緯度経度 (φ,λ) 及び子午線収差 (γ) をプログラム (SOK5.BAS) を用いて求めてみよう ( 解 ) 緯度 B= 34º 経度 L= 135º

7 7 子午線収差 (γ)= UTM 図法 UTM 投影 (Universal Transverse Mercator Projection) による平面座標は世界的に利用されているわが国では縮尺 1/1 万地形図から 1/5000 基本図,1/5 万編集図における投影に採用されている UTM 座標系は赤道を 6 ずつ東回りに 60 個のゾーンに分割して横メルカトールで投影する図法であるつまり最初のゾーン 1 は東周りに西経 (-)180 ~174 範囲 ( 中央経線西経 177 ) ゾーンは西経 174º~168º としゾーン 30 は西経 6 ~0 ( 中央経線西経 3 ) ゾーン 31 は東経 0 ~6 ゾーン 60 は東経 174 ~180 ( 中央経線東経 177 ) である因みに大阪や沖ノ鳥島はゾーン 53 で東経 13 ~138 中央経線東経 135 ( 明石市 ) に該当する ( この計算は utm_zone.xlsm で作成している ) 日本はゾーン 51~56( 中央経線 13º( 尖閣 ) 19º( 沖縄 ) 135º( 大阪 ) 141º( 東京 ) 147º 153º( 北方四島 )) を採用する各経度帯内ではその中央経線と赤道との交点を座標原点として横メルカトル図法で投影する中央経線から東西方向に離れるに従って距離ひずみが増大するので平面距離 s と球面距離 S の点の投影誤差 ds/ds の比率は ±4/10000 で収まるように考えるつまり中央経線上の縮尺係数はに設定するので中央経線より東西に 180km 離れたところの縮率はそして同様に 70km 離れた地点の縮率はとなる UTM 投影では中央経線より 70km までは距離の相対精度が ±4/10000 に収まるので各 6 幅の経度帯内で地球を平面とみなして投影することができる ( ただし原点から 3 離れると 334kmなのでを超える )UTM 座標においてx 座標 (N 軸 ) は中央経線 ( 縦軸 ) に沿って北半球では 0 より増加し南半球では 10,000km( 擬北距 =F.N.) より減少させる y 座標 (E 軸 ) は赤道 ( 横軸 ) を 500km( 擬東距 =F.E.) として E 軸に沿い東側へは増加西側へは減少する Bessel 緯度経度 ( 以下の Ki は GRS80 での値である ) と xy(xn,ye) 座標との関係式は以下のように表すことができる.6.1 緯度経度から UTM 座標の変換 x N = F. N. +m o [Δφ + Δλ Δλ4 Nsinφcosφ + 4 Nsinφcos3 φ(5 t + 9η + 4η 4 ) + Δλ 6 70 Nsinφcos5 φ(61 58t + t η 330η t )] (.6.1) y E = F. E. +m o [Δλ Ncosφ + Δλ3 6 Ncos3 φ(1 t + η ) + Δλ 5 10 Ncos5 φ(5 18t + t η 58t η )] (.6.) γ = Δλsinφ + Δλ3 3 sinφcos φ 1 + 3η + η 4 + Δλ 5 15 sinφcos4 φ( t ) (.6.3) Δφ = a(1 e )[k 1 φ φ o 1 k sinφ φ o k 3 sin4φ 4φ o 1 k 6 4 sin6φ 6φ o + 1 k 8 5 sin8φ 8φ o 1 k 10 6 sin10φ 10φ o + ] (.4.4) K1=

8 8 K = D -3 K3 = D -5 K4 = D -8 K5 = D -11 K6 = D -14 Δλ=λ-λo λo= 座標系原点の経度, t=tanφ η=e cosφ e = a b b mo= a=6,378,137m 1/f= 南北半球では F.E.=500km, 北半球では北距 F.N.=0m 南半球 F.N.=10,000km [ 例題.6.1] ある三角点の GRS80 による世界測地系の緯度経度は次のとおりとするときこの点の UTM 座標および子午線収差 ( 角 ) をプログラム (SOK61.TT) により求めてみよう φ=35º λ=139º またこの点の UTM 座標の原点は次の値である φo=0º λo=141º ( 解 ) = 3,946, Y= 386, 子午線収差 γ =- 0 o UTM から緯度経度への変換 φ = φ 1 y E NM t + y E 4 4MN 3 t(5 + 3t + η 9t η ) y E 6 70MN 5 t( t + 45t 4 ) (.6.6)...(..1) λ = λ o + [ y E Ncosφ y E 3 6N 3 cosφ (1 + t + η ) + y E 5 10N 5 cosφ (5 + 8t + 4t 4 )] (.6.7) γ = y E N t y E 3 3N 3 t(1 + t η ) + y E 5 15N 5 t( + t + 3t 4 )...(.6.8) 南北半球では y E = (y E 500,000[m])/m o

9 9 m o = ただし γ は子午線収差角である φ 1 = φ1 (i 1) S (i 1) D M (i 1) D = a(1 e ) k 1 φ o k sinφ o + k 3 4 sin4φ o k 4 6 sin6φ o + k 5 8 sin8φ o k 6 10 sin10φ x o + + m o S = a(1 e )(k 1 φ 1 k sinφ 1 + k 3sin4φ 1 4 k 4sin6φ k 5sin8φ 1 8 k 6sin10φ ) M = a(1 e ) (1 e sin φ 1 ) 3/ (..4) 北半球では x = x N F. N. = x N 南半球では x = x S F. N. = x S 10,000,000[m] (..5) である [ 例題.6.] ここでは例題.6.1の逆計算をプログラム (SOK6.TT) を用いて行ってみようただし緯度経度を求めるために用いる UTM 座標は次のとおりである xn=3,946,757.90m yn=386, m ( 解 ) Input data(utm x,y) n= 3,946,757.90m Ye= 386, m LatitudeB= LongitudeL= Meridian Convergence γ= GPS 観測値の平面直角座標への計算 GPS 観測値であるWGS84における緯度経度楕円体高 ( ϕ, λ, h) から Bessel の (φ, λ, h) に変換する式は以下のとおりである (1) GPS 観測値参照点または求点のDGPS 観測から得られたWGS84 1 の経緯度および楕円体高は次のとおりとする φ GPS = λ...(3..1) h WGS84 1 World Geodetic System(1984):1984 年に人工衛星を利用して決められた回転楕円体

10 10 () GPS 観測値から WGS84 の YZ への変換 WGS84 楕円体における 3 次元直交座標 YZ は次式で算出する WGS84 = (N 84 + h 84 )cosφ 84 cosλ 84 Y WGS84 = (N 84 + h 84 )cosφ 84 sinλ 84 (3..) Z WGS84 = [N 84 (1 e 84 ) + h 84 ]sinφ 84 なお以上で用いた要素は次のとおりである卯酉線曲率半径 N 84 = a 84 W 84 W 84 = 1 e 84 sin φ 84 扁平率の逆数 1/f 84 = a 84 a 84 b 84 短半径 b 84 = a 84 (1 f 84 ) 第一離心率 e 84 = a 84 b 84 a...(3..3) 84 また WGS84 の地球楕円体の値は次のとおりである a 84 = 6,378,137m 1 f 84 = b 84 = 6,356,75.314m (3..4) (3) WGS84 から ITRF89 系座標への変換 Y = R T Y (3..5) Z ITRF Z (4)ITRF89 系から筑波座標系 9(Bessel) 座標への変換 Y = R T Y (3..6) Z Tsukuba Z ITRF (5) Bessel 経緯度への変換 Y = Y (3..7) Z Bessel Z Tsukuba Bessel 経度緯度および正規標高 (Orthometric Height) はこれらの YZ 値を用いて以下のように求めることができる tanλ B = Y/...(3..8) ( 近似式 ) tanφ B = tanθ = ( 厳密式 ) Z+(e ) bsin 3 θ +Y e acosθ Z a b +Y...(3..9)

11 11 tanφ B = Z+N Be B sinφ B +Y H B = +Y N cosφ B B h B = H B H G...(3..10) 上に示した Bessel 楕円体上の標高は次式で求められる hb=( 水準標高 H)+( ジオイド高 HG)...(3..11) またジオイド高 1 はジオイド図を参照して求められる Bessel 楕円体は次のとおりの値である a B = 6,377, m b B = 6,356, m 1 f B = (3..1) GRS80(geodetic reference system 1980) では長半径 a= m 及び扁平率 f= 1/ b=a(1-f)= b = m WGS84(GPS 衛星で使用 ) a= m, f=1/ ロシアの GLONASS の使用する楕円体 (PZ90) GLObal naya NAvigationnaya Sputonikovaya Sistema= 英語 Global Navigation Satellite System ロシア宇宙軍の運用する GPS 衛星が GLONASS であり GLONASS では座標系として PZ-90 という北極点の位置を 1900 年から 1905 年までの測定位置の平均としたものを採用しているこれに対して GPS では WGS84 という 1984 年時点での北極点の測定位置を元にした測地系を採用している PZ-90 (Parametry Zemli 1990, PE-90; Parameters of the Earth 1990) は測地系の基準楕円体に SGS-90 を用いた ECEF 座標系であり原点を地球重心にとり Z 軸が IERS の国際基準座標系の 3 軸を慣用極方向として軸はグリニジ子午線と赤道面との交点と地球重心を結ぶ軸 Y 軸は赤道面上で軸に直角な軸で -Y-Z は右手系を採用 PE-90 楕円体は長軸半径 a= m 扁平率の逆数 1/f= である 1 基準楕円体表面の標高をゼロとしたときの平均海面高つまり楕円体面からのジオイド面の高さ

Microsoft Word - 03基準点成果表

Microsoft Word - 03基準点成果表基準点成果表 ( 情報 ) < 試験合格へのポイント > 基準点成果表 ( 又は基準点成果情報 ) に関する問題である近年では基準点成果表の項目 ( 内容 ) に関する問題よりは平面直角座標系に絡めた問題が出題されているため平面直角座標系の特徴も併せて覚える方か良いここでは水準点を除くものを基準点として記述する基準点について ( : 最重要事項 : 重要事項 : 知っておくと良い )

第2章 位置の計算

第2章位置の計算