相対的な誤差楕円も表示する平面の測量網計算プログラムの開発
松井佐久夫、木村洋昭 （RIKEN）
木内 淳(SPring‐8 サービス（株）)

1）測量の網平均計算
2）BASICプログラムの検証と改良、使用実績
3)相対的な誤差楕円4）シミュレーションにも便利
5）Excel VBAへの書換え
6）まとめ1）測量の網平均計算

巻 改訂版 基準点測量 ﹄細野武庸 井内登 日本測量協会﹃測量叢書 第水準網調整計算プログラム

1（1992）

細野武庸 他：『測量叢書 第1巻 改訂版 基準点測量』,
日本測量協会（1992） ｐ２５２〜２５６
標準偏差の計算
正規方程式の組立て観測件の読み込み
正規方程式の解

角の観測方程式 シュライバーの消去法の組立て
サブルーチン等

距離の観測方程式の組立て

2）BASICプログラムの検証と改良、使用• プログラムの改良
・入力部：名前をそのまま使う
・計算結果から推測される角度、距離と実測値の差の打ち出し
・（絶対・相対）誤差楕円部の追加
・逆行列の計算時 桁落ち防止のための軸だし・既知点として一軸のみ固定で計算できるように改良・グラフィック部の追加• アプリケーション N88BASICから FBASICへ・Windows7でも走る

測量会社の使用するプログラムの場合• SPring‐8 建設時

検定番号（日本測量協会）
No.61‐96 S 62.02.27
登録番号（国土地理院）
No.62.A‐12 S 62.06.04
コピー検定番号（日本測量協会） No. 4‐20 H 4.5.1
コピー登録番号（国土地理院） No. 4.A‐39 H 4.5.28

電算プログラム検定記録書

（検証1）（株）山口測量による結果との比較
重み 角度(3.5秒)、距離(1cm)、入力データ同じにして計算
結果：良く一致した

（検証２） 蓄積リング dr 半径方向の距離

USMNによる解析との比較
Spatial Analyzer(SA) （New River Kinematics社）(USNM:Unified Spatial Metrology Network)

dr(mm)

2.0

測定個数 角度 1394
距離 1394

1.00.0-1.0

リング全周-2.006121824

2014 SA-USMN
2014 F-BASIC
30

Cell number

364248

（使用実績1） SPring8蓄積リングの基準点とモニュメント

（使用実績２） XSBT系 モニュメント設置と網平均計算
SACLAからSPring‐8のシンクロトロンへのBT系
高さの差〜９ｍあったが水平面内での計算
傾斜部での測量の様子
シンクロトロン

180ｍ

SACLA光源棟

（使用実績3） XFEL SACLA モニュメントの網平均計算
直線に近くても網平均で計算
最初は機器を据え付けるため、その後は変化をモニターしている

TDA5005 による測量→FBASICによる計算

3) 誤差楕円の計算方法
Ａ=観測方程式の行列
Ｐ=重みの行列
Ｌ=設定値から求まる値と観測値の差の行列
Ｘ=解
最小二乗法の正規方程式
ＡTＰＡＸ＋ＡTＰＬ=0

Ｘ=‐(ＡTＰＡ)‐1(ＡTＰＬ)
相対的な誤差楕円(ＡTＰＡ)‐1 : 分散、共分散行列シミュレーションと誤差楕円の例
SACLA光源棟モニュメント測量の場合
測定値のばらつき＝解析時のエラーの場合2σの誤差楕円の内側に86％が入る

シミュレーション
測定値400セット

4）シミュレーションの効用 その１
4‐1 測量の実力値の推測結果から推測した距離と実測値の差

4）シミュレーションの効用 その２
4‐2 測量網の最適化
相対的誤差楕円の計算より
リングの測量網

絶対的な誤差楕円絶対的と相対的な誤差楕円シミュレーションの条件
測量の精度
角度 0.5秒
距離10μm+0.4ppm

楕円の半径：１σシミュレーション例 １）

絶対的な誤差楕円 測量の角度精度依存
角度精度 １秒

0.7秒

0.5秒

変わる！

σ〜0.5mmあるが周長は２０回のシミュ
レーションでは‐0.03±0.07mm

距離精度 10μm+距離×0.4×10‐6

シミュレーション例 ２）

絶対的な誤差楕円 距離精度依存

距離精度
20μm+距離×5×10‐6楕円の短軸は少し変わるが長軸の大きさはあまり変わらない
30μm+距離×5×10‐6

角度精度 １秒

20μm+距離×1×10‐6シミュレーション例 3）

相対的な誤差楕円 角度精度依存
角度精度 １秒

0.8秒変わる！

距離精度 20μm+距離×0.4×10‐60.5秒

シミュレーション例 3）

相対的な誤差楕円 距離の固定部精度依存
近距離の測定も多いため
距離精度20μm+距離×0.4×10‐6変わる！

1σ=33μm

10μm+距離×0.4×10‐6

1σ=24μm

角度精度 0.5秒

5）Excel VBAへの書換え
Excel シートに解析のエラーを入力
角度 （sec）
距離 （固定部 (mm)）（距離に依存する部分 ppm）
半固定する Ｘ軸、Ｙ軸、なし
繰り返し回数既知点 ＸＹ座標
未知点 ＸＹ座標
角度 データ シミュレーションの場合
データなし
距離 データ解析 繰返し回数 時刻

相対的誤差楕円計算結果

結果から推測されるきょり、角度と実測値との差絶対的誤差楕円VBA 計算の一部（プロシーディングには計算の主要部はすべて載せている）
シミュレーションの場合
一様乱数→正規乱数'‐‐ Calculation of Coefficients of observation eq.
c1 = XJ(j0, iter) ‐ XJ(j1, iter): d1 = YJ(j0, iter) ‐ YJ(j1, iter)
If (c1 = 0) * (d1 = 0) ThenElse: e1 = d1 / (c1 ^ 2 + d1 ^ 2): f1 = ‐c1 / (c1 ^ 2 + d1 ^ 2)
End Ifc2 = XJ(j0, iter) ‐ XJ(j2, iter): d2 = YJ(j0, iter) ‐ YJ(j2, iter)
If (c2 = 0) * (d2 = 0) ThenElse: e2 = d2 / (c2 ^ 2 + d2 ^ 2): F2 = ‐c2 / (c2 ^ 2 + d2 ^ 2)
End Ifx = ‐c1: y = ‐d1: t1 = Atn(y / x)
If x < 0 Then
t1 = t1 + pi
ElseIf y < 0 Then
t1 = t1 + pi * 2
End IfFunction bunp(sigma As Double) As Double' Trans Random Number from flat distribution toNormal one
Dim rn As Double, sum As Double, xran As Double,
dxran As Double, comp As Double, pn As Double
Dim iran As Integerpi = 3.141592653: sqpi = Sqr(pi)rn = Rnd ‐ 0.5: pn = Sgn(rn): If rn < 0 Then rn = ‐rn
sum = 0#: xran = 0#: dxran = 0.01
comp = rn * sqpi / dxran
For iran = 1 To 30000
sum = sum + Exp(‐xran * xran)If sum > comp Then rnn = xran: GoTo 14210
xran = xran + dxran
Next iran14210: bunp = pn * rnn * sigma / 0.707107
End Function

解析結果に続いて
For i = nr + 1 To nellipse
誤差楕円の図の表示

With ActiveSheet.Shapes.AddLine(x1 + xc, y1 + yc, x2 + xc, y2 + yc)
End Withxc = XJ(i, maxiter) * xmulti + shiftx: yc = YJ(i, maxiter) * ymulti + shifty: aa = SMAX_abs(j) * multi: bb = SMIN_abs(j)* multi: tilt = ‐ANGLE(j) / ra * degrad
tilt = Angle_abs(i) * ra
For t = 0 To 360 Step 0.2: tt = t * degrad
x1 = aa * Cos(tt) * Sin(tilt) + bb * Sin(tt) * Cos(tilt): y1 = aa * Cos(tt) * Cos(tilt) ‐ bb * Sin(tt) * Sin(tilt)
x2 = aa * Cos(tt + dt) * Sin(tilt) + bb * Sin(tt + dt) * Cos(tilt): y2 = aa * Cos(tt + dt) * Cos(tilt) ‐ bb * Sin(tt + dt)* Sin(tilt)
With ActiveSheet.Shapes.AddLine(x1 + xc, y1 + yc, x2 + xc, y2 + yc)
End With
Next tj=j+1
Next i

計算時間 （Windows７）

場所

測量点

角度

距離

FBASIC

VBA

光源棟50

111

133

<1秒

<1秒

蓄積リン
グ

486

1300

1400

26分
(36分)

10分
(19分)

()Windows XP

まとめ・教科書『基準点測量』の水平面内の網平均計算の具体的なプログラムを基に改良も加えこれまでのSPring8の多くのところで使用し実用上問題はなかった。
・シミュレーションや相対的な誤差楕円の計算など目的や設定が決まっている場合は測量器に付属しているソフトに比べ使いやすい。
・測量機器の精度、空気などの環境やターゲット、測量など含めた実力値はシミュレーションとの比較からある程度推測可能だと考えられる。
・シミュレーションで周長のエラーが評価できる。現状の測量網ではランダムエラーの影響は全周ではキャンセルされシステマティックエラーのみによることがわかった。
・シミュレーション、特に相対的な誤差楕円の値から測量網の最適化が検討できる。ライカ製
AT402の実力を十分発揮できれば、蓄積リングの半径方向の相対的なエラーは現状の測量網でも20μm(1σ)程度で収まることがわかった。
BASIC系ソフト （Excel VBA含む）
・計算と同時に誤差楕円も描けるので便利である。
・BASICからExcel VBAへの変換は同じ文法が使える等スムーズに実現できた。
・Excel （VBA）で解析でき使い勝手が向上する。
・計算時間はＦＢＡＳＩＣコンパイラーよりExcel VBAの方がかなり速い。測点数500、測定数3000程度で１回の計算にはWindows7で10分程度である。
・できればダウンロードできるところに置き、使用する人でさらに使いやすいものに改良してもらえるといいのではないか、と考えている。