JP2006064892A - 位置計算装置、位置計算方法、位置計算プログラム、及び位置計算プログラムの記憶媒体 - Google Patents
位置計算装置、位置計算方法、位置計算プログラム、及び位置計算プログラムの記憶媒体 Download PDFInfo
- Publication number
- JP2006064892A JP2006064892A JP2004246064A JP2004246064A JP2006064892A JP 2006064892 A JP2006064892 A JP 2006064892A JP 2004246064 A JP2004246064 A JP 2004246064A JP 2004246064 A JP2004246064 A JP 2004246064A JP 2006064892 A JP2006064892 A JP 2006064892A
- Authority
- JP
- Japan
- Prior art keywords
- latitude
- meridian
- value
- calculation
- scale factor
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Images
Landscapes
- Instructional Devices (AREA)
Abstract
【解決手段】 演算部6が、赤道から当該座標値xまでの子午線弧長Wを算出し、子午線曲率半径を以って緯度φj にて赤道から垂線緯度φj まで積分することにより子午線弧長S(φj )を算出する。演算部6は、緯度φj を与えることによって
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 )を計算して緯度φj+1 を算出する。緯度φj+1 と緯度φj との差の絶対値が基準値以内の場合、
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 式3を計算し、求めた緯度φk+1 を緯度φn とし、演算部6に当該緯度φn における縮尺係数mt+1 を算出させ、当該縮尺係数mt+1 を用いて求めた、子午線弧長Wにて、上記式1の計算を行わせる。
【選択図】 図1
Description
尚、地球楕円体を測量の基準にするためには、楕円体の中心を実際の地球上のどの位置に、またその楕円体の座標軸が実際の地球のどこを通るかということを決める必要があり、この位置と方向が決められた地球楕円体が準拠楕円体とされる。
日本で用いられている平面直角座標は、ガウス・クリューゲルの等角投影法によるもので、座標系原点を通る子午線は等長に、図形は等角の相似形に投影される。しかし、距離については、基準子午線から東西に離れるに従って平面距離が増大していくため、投影距離の誤差を相対的に1/10,000以内に収めるよう座標原点に上記の縮尺係数(0.9999)を与え、かつ、基準子午線より東西130km 以内を適用範囲とした座標系(平面直角座標系は日本全国を19の座標系に区分している。)を設けている(非特許文献1)。
なお、座標系のX軸は、原点において子午線に一致する軸とし、原点から真北に向かう値を正としている。Y軸は原点においてX軸に直交する軸とし、真東に向かう値を正としている。この点、三次元直交座標(X,Y, Z)とは定義が異なる。又、各座標系原点の値はX=Y=0.000 メートルとしている。
このプログラムを利用することにより、web上において、上記座標値x,yの入力により、簡便に、その緯度・経度を算出することができる。
ところが、このプログラムにおいて、座標値x,yが属する座標系の基準子午線から、東西に離れるに従って、その誤差が大きくなる。
この誤差は、国土地理院ホームページの当該「緯度・経度への換算」プログラムで得た緯度・経度を、同じく国土地理院ホームページ「平面直角座標への換算」(非特許文献4)プログラムにて、簡単に確認することができる。
a−1)は当該座標系3の座標原点(0,0)を通る子午線上の点であり、a−2)〜a−5)において、入力する座標値Yの値を基準子午線から西方へ徐々に遠ざかるものとしている。
そして上記「平面直角座標への換算」プログラムにおいて、上記a−1)〜a−5)で入力した座標値と、換算結果が一致するように、試行錯誤で北緯と東経を入力した(北緯と東経の値を種々変えて、出力値が、「平面直角座標への換算」プログラム)における入力値と一致するものを探った)。そして、b−1)〜b−5)に示す結果(平面直角座標上の座標値)を得た。
入力値 X座標 -207504.143 m Y座標 0 m
換算結果 緯度 34°07'45.93167" 経度 132°10'00.00000"
真北方向角+ 0°00'00.00" 縮尺係数 0.99990000
a−2)
入力値 X座標 -207494.938 m Y座標 -59997.338 m
換算結果 緯度 34°07'40.02853" 経度 131°30'58.32911"
真北方向角+ 0°21'53.81" 縮尺係数 0.99994436
a−3)
入力値 X座標 -207483.393 m Y座標 -90082.603 m
換算結果 緯度 34°07'32.62476" 経度 131°11'24.23439"
真北方向角+ 0°32'52.52" 縮尺係数 1.00000000
a−4)
入力値 X座標 -207467.727 m Y座標 -119340.074 m
換算結果 緯度 34°07'22.57910" 経度 130°52'22.58340"
真北方向角+ 0°43'32.98" 縮尺係数 1.00007551
a−5)
入力値 X座標 -207462.646 m Y座標 -127394.933 m
換算結果 緯度 34°07'19.32100" 経度 130°47'08.30605"
真北方向角+ 0°46'29.28" 縮尺係数 1.00010000
入力値 北緯 34°07'45.93167" 東経 132°10'00.00000"
換算結果 X座標 -207504.143 m Y座標 0.0000 m
真北方向角+ 0°00'00.00" 縮尺係数 0.99990000
b−2)
入力値 北緯 34°07'40.02854" 東経 131°30'58.32910"
換算結果 X座標 -207494.938 m Y座標 -59997.338m
真北方向角+ 0°21'53.81" 縮尺係数 0.99994436
b−3)
入力値 北緯 34°07'32.62478" 東経 131°11'24.23440"
換算結果 X座標 -207483.393 m Y座標 -90082.603m
真北方向角+ 0°32'52.52" 縮尺係数 1.00000000
b−4)
入力値 北緯 34°07'22.57909" 東経 130°52'22.58341"
換算結果 X座標 -207467.727 m Y座標 -119340.074 m
真北方向角+ 0°43'32.98" 縮尺係数 1.00007551
b−5)
入力値 北緯 34°07'19.32099" 東経 130°47'08.30603"
換算結果 X座標 -207462.646 m Y座標 -127394.933 m
真北方向角+ 0°46'29.28" 縮尺係数 1.00010000
ところが、a−2)について、上記と同様の作業を行ったb−2)を見ると、a−2)における換算結果と、b−2)における入力値とは一致していない。即ち、前者a−2)が緯度を34°07'40.02853"とし経度を 131°30'58.32911"とするのに対して、後者b−2)が北緯を34°07'40.02854"とし東経を 131°30'58.32910"とし、緯度及び経度について、夫々小数点以下第5位(最下位)において食い違いが生じている。
このような誤差の発生について、a−3〜a−5)と対応するb−3)〜b−5)と見れば、基準子午線から西方へ離れるに従って、特に経度についての当該誤差が大きくなっていることが分かる。
例えば、基準子午線から西方へ約127km離れた地点であるa−5)における換算結果(緯度34°07'19.32100"、経度 130°47'08.30605")と、対応するb−5)における入力値(北緯34°07'19.32099"、東経 130°47'08.30603")とでは、特に、緯度について小数点以下第3位にて既に食い違いが生じている。
上記において、そのような誤差が生じる原因は、プログラムの内容が公開されていないことから、本願発明者において正確なことは分からない。
また、発明者にて可能な範囲で、先行技術文献として、特許公報を調べたが、同種の発明が記載された公報を発見するに至っていない。
そこで、発明者は、この点を熟考し、推測するに、上記国土地理院のプログラムにおいて、緯度について例えば、国土地理院ホームページ「赤道からの子午線弧長を与えて緯度を求める計算」(非特許文献2)などに示される公式を用いて、緯度についての収束計算を行っているのであろうが、基準子午線から離れた場合の前述の縮尺係数の変化が十分に反映するものとなっておらず、その結果、収束(反復計算)により安定状態となった緯度の精度が十分なものでない(ずれた位置に収束してしまっている)と考える。
具体的には、国土地理院の上記プログラムは、非特許文献2に示す国土地理院ホームページの「6.赤道からの子午線弧長を与えて緯度を求める計算」で示されているφn+1 を求める式を用いて反復計算を行うものであると考えられる(この式を用いた結果が上記と一致した)。この式を用いて反復計算を行うと上述の通り、大きな誤差を持った値で数値が安定(収束)する。その理由を検討すると、この式の計算において、縮尺係数は、原点における縮尺係数m0 (=0.9999)に固定されており、基準子午線以外の、求めようとする座標位置における縮尺係数は、原点の縮尺係数m0 (=0.9999)と異なるものであるという事実は、無視されている(尚、国土地理院の式で用いているMは、赤道からの座標系原点までの子午線弧長を示し、この明細書で用いる子午線曲率半径を示すMとは異なる)。
現状では、推測の域を出るものではないが、上記プログラムを用いた結果は、歴然としている。
従って、当該プログラムは、著しく信頼性の欠くものとなっており、厳密性が要求される状況において、その利用が困難となっている。
即ち、この位置計算装置は、座標入力部と、出力部と、演算部と、判定部と、回帰制御部とを備える。座標入力部は、上記座標値x,yの入力を受け付ける。演算部は、直角座標の座標値xが与えられることにより、地球楕円体の赤道から当該座標値x,yが属する座標系原点までの子午線弧長に対し、当該座標値xを仮の縮尺係数mt で除した値を、加算して、赤道から当該座標値xまでの基準子午線上の子午線弧長Wとする。
そして、この子午線弧長Wに対応する基準子午線上の緯度、即ち、上記入力座標から基準子午線上に下ろした垂線の緯度を算出するため、演算部は、地球楕円体の長半径と第1離心率と緯度とから定まる子午線曲率半径について、地球楕円体の長半径aと第1離心率eと共に任意の緯度φj を初期値として与え、別途に子午線曲率半径を以って赤道から上記緯度φj まで積分することにより、赤道からの基準子午線上の子午線弧長S(φj )を算出する。更に、演算部は、基準子午線上の上記緯度φj を与えることによって、式1を計算して基準子午線上の緯度φj+1 を漸近値として算出する。
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 )…式1
判定部は、上記式1にて算出した基準子午線上の上記漸近緯度φj+1 と基準子午線上の緯度φj との差の絶対値と、基準値との比較を行うものである。判定部にて、基準子午線上の、上記漸近緯度φj+1 と緯度φj との差の絶対値が基準値を超えると判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするには十分収束していないと判定した場合、回帰制御部は、この漸近緯度φj+1 を上記基準子午線上の緯度φj とし、演算部に基準子午線上の赤道からの子午線弧長S(φj )を算出させると共に当該子午線弧長S(φj )から上記式1の計算にて再度漸近緯度φj+1 を算出させ、判定部に再度上記比較を行わせる。また、判定部が、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値以内と判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするのに十分収束したと判定した場合、回帰制御部は、当該漸近緯度φj+1 を、入力座標から基準子午線上に下ろした垂線の緯度φk として、仮の縮尺係数mt と、当該垂線緯度φk と、当該垂線緯度φk における子午線曲率半径Mk 及び卯酉線曲率半径Nk と、e’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、式3を計算して緯度φk+1 を演算部に算出させる。
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 …式3
回帰制御部は、上記緯度φk+1 を入力座標に対応する緯度φn として、演算部に、縮尺係数mt+1 を算出させ、前記仮の縮尺係数mt に代え当該縮尺係数mt+1 を用いて求めた、基準子午線上の子午線弧長Wにて、上記式1の計算を行わせる。
ここでいう地球楕円体には、各種の準拠楕円体を含む。
基準子午線とは、当該座標が属する座標系の原点を通る子午線をいう。
また、垂線緯度とは、入力座標から基準子午線に下ろした垂線と、当該基準子午線との交点の緯度をいう。
また、積分には、その近似計算を含む。
上記において、回帰制御部は、演算部に対して直接、計算の指令を下すものに限定するものではなく、間接的に、演算部へ計算を行わせるものも含む。
mt ={1+(y2 /2m02R02)+(y4 /24m04R04)}×0.9999…式2
上記の原点における平均曲率半径R0 は、当該座標が属する座標系の原点の緯度φ0 又はその代替値と共に、地球楕円体の長半径a及び短半径b、或いは当該半径a,bから導かれる第1離心率e、第2離心率e’を用いて定めたものである。
本願第3の発明は、上記本願第2の発明にあって、次の特徴を備えた位置計算装置を提供する。即ち、上記縮尺係数mt+1 は、原点縮尺係数m0 に代え上記縮尺係数mt を用い、且つ、原点における平均曲率半径R0 に代え演算部に式3にて算出させた緯度φn における平均曲率半径Rφn を用いて、演算部が、上記式2により算出するものである。
Δλ=(y/mt )/Nk cos φk
−(y/mt )3 (1+2tk 2 +ηk 2 )/6Nk 3cosφk
+(y/mt )5 (5+28tk 2 +24tk 4 )/120Nk 5cosφk …式4
演算部は、当該座標系が属する原点における経度λ0 に対し当該経度差分Δλを加算することにより、上記座標値x,yにおける経度λを算出するものである。
γ=(y/mt )tk /Nk
−(y/mt )3 tk (1+tk 2 −ηk 2 )/3Nk 3
+(y/mt )5 tk (2+5tk 2 +3tk 4 )/15Nk 5
…式5
即ち、この位置計算プログラムは、座標入力部と、出力部と、演算部と、判定部と、回帰制御部とを備えたコンピュータに導入される。このプログラムが導入されたコンピュータは、座標入力部にて、上記座標値x,yの入力を受け付け、演算部に、直角座標の座標値xを与えて、地球楕円体の赤道から当該座標値x,yが属する座標系原点までの子午線弧長に対し、当該座標値xを仮の縮尺係数mt で除した値を加算して、赤道から当該座標値xまでの基準子午線上の子午線弧長Wとする。
そして、この子午線弧長Wに対応する基準子午線上の緯度、即ち、上記入力座標値から基準子午線上に下ろした垂線の緯度を算出するため、演算部にて、地球楕円体の長半径と第1離心率と緯度とから定まる子午線曲率半径について、地球楕円体の長半径aと第1離心率eと共に任意の緯度φj を初期値として与え、別途に子午線曲率半径を以って赤道から上記緯度φj まで積分することにより、赤道からの基準子午線上の子午線弧長S(φj )を算出し、更に、演算部に、基準子午線上の緯度φj を与えることによって、式1を計算して基準子午線上の漸近緯度φj+1 を漸近値として算出させる。
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 )…式1
そして、判定部に、上記式1にて算出した基準子午線上の漸近緯度φj+1 と基準子午線上の緯度φj との差の絶対値と、基準値との比較を行わせる。判定部にて、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値を超えると判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするには十分収束していないと判定した場合、回帰制御部により、漸近緯度φj+1 を上記基準子午線上の緯度φj とし、演算部に基準子午線上の赤道からの子午線弧長S(φj )を算出させると共に当該子午線弧長S(φj )から上記の式1の計算にて再度漸近緯度φj+1 を算出させ、判定部に再度上記比較を行わせるものである。また、判定部にて、基準子午線上の漸近緯度φj+1 と緯度φj との差の絶対値が基準値以内と判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線とするのに十分収束したと判定した場合、当該漸近緯度φj+1 を入力座標からおろした垂線の緯度φk として、回帰制御部により、仮の縮尺係数mt と、当該垂線緯度φk と、当該垂線緯度φk における子午線曲率半径Mk 及び卯酉線曲率半径Nk と、e’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、式3を計算して緯度φk+1 を演算部に算出させる。
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 …式3
更に、回帰制御部により、上記緯度φk+1 を入力座標に対応する緯度φn として、演算部に、縮尺係数mt+1 を算出させ、前記仮の縮尺係数mt に代え当該縮尺係数mt+1 を用いて求めた、基準子午線上の子午線弧長Wにて、上記式1の計算を行わせるものである。
即ち、この方法は、座標入力部と、出力部と、演算部と、判定部と、回帰制御部とを備えたコンピュータを用いるものである。そして、座標入力部にて、上記座標値x,yの入力を受け付ける。演算部に、直角座標の座標値xを与えて、地球楕円体の赤道から当該座標値x,yが属する座標系原点までの子午線弧長に対し、当該座標値xを仮の縮尺係数mt で除した値を加算して、赤道から当該座標値xまでの基準子午線上の子午線弧長Wとし、演算部にて、地球楕円体の長半径と第1離心率と緯度とから定まる子午線曲率半径について、地球楕円体の長半径aと第1離心率eと共に任意の緯度φj を初期値として与え、別途に子午線曲率半径を以って赤道から上記緯度φj まで積分することにより、赤道からの基準子午線上の子午線弧長S(φj )を算出し、演算部に、基準子午線上の緯度φj を与えることによって、式1を計算して基準子午線上の漸近緯度φj+1 漸近値として算出させる。
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 )…式1
判定部にて、上記式1にて算出した基準子午線上の漸近緯度φj+1 と基準子午線上の緯度φj との差の絶対値と、基準値との比較を行わせるものである。そして、判定部にて、上記基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値を超えると判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線に下ろした垂線の緯度とするには十分収束していないと判定した場合、回帰制御部により、漸近緯度φj+1 を上記基準子午線上の緯度φj とし、演算部に基準子午線上の赤道からの子午線弧長S(φj )を算出させると共に当該子午線弧長S(φj )から上記の式1の計算にて再度漸近緯度φj+1 を算出させ、判定部に再度上記比較を行わせる。判定部にて、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値以内と判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線に下ろした垂線の緯度とするのに十分収束したと判定した場合、回帰制御部は、当該漸近緯度φj+1 を、入力座標から基準子午線上に下ろした垂線の緯度φk として、仮の縮尺係数mt と、上記垂線緯度φk と、当該垂線緯度φk における子午線曲率半径Mk 及び卯酉線曲率半径Nk と、e’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、式3を計算して緯度φk+1 を演算部に算出させる。
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 …式3
そして、回帰制御部により、上記緯度φk+1 を入力座標に対応する緯度φn として、演算部に、縮尺係数mt+1 を算出させ、前記仮の縮尺係数mt に代え当該縮尺係数mt+1 を用いて求めた、基準子午線上の子午線弧長Wにて、上記式1の計算を行わせるものである。
即ち、この発明では、入力座標から基準子午線上に下ろした垂線の緯度(垂線緯度)を算出するに当たり、緯度の収束計算に適した式1をイテレーション式として用いて緯度の計算(反復計算)を行い、得た垂線緯度から式3を用いて当該点(入力座標値に対応する)緯度を算出し直すことにより、より精度の高い緯度を得ることができる。更に、このようにして求めた緯度から縮尺係数を算出して当該縮尺係数の収束を加味した後処理を行うことにより、極めて精度の高い当該点の緯度の算出を可能とした。
また、特に、本願第4の発明では、上記の計算を可能としつつも、基準子午線に近い場合など十分に誤差の小さい場合にまで、後処理として、式3を経て得た上記の緯度で算出した縮尺係数による、式1の計算を強要するものではなく、縮尺係数の収束状況で、後処理を行うか否かを決定し、無駄な反復計算を抑制し、迅速に計算を完了する手段を備えた装置を提供し得た。
更に本願第5の発明では、上記の当該点の緯度と共に、経度λを、また本願第6の発明では、子午線収差を、夫々算出し得る位置計算装置を提供した。
本願発明は、平面直角座標上の座標x,yから、対応する地球楕円体上の緯度、更には、経度、子午線収差を精度よく算出する位置計算方法を提供するものである。即ち、ガウス・クリューゲルの等角投影法により地球楕円体を投影する平面直角座標の座標点の座標値から、当該地球楕円上の対応する点の緯度、経度及び子午線収差を、高精度に算出する方法を提供するものである。
この方法は、垂線緯度の算出のための収束計算において、縮尺係数の収束を加味して、計算結果の精度を向上させたものである。
また、上記後処理工程に先立ち、算出した縮尺係数が収束したかを調べ、その結果当該縮尺係数が収束していた場合、後処理工程を行わずに、計算を終了することもできる。
従って、主工程の処理に続いて後処理工程を行い上記にて算出した縮尺係数を加味した結果を得るものとしてもよく、また、主工程の結果に対して、算出した縮尺係数の収束の判定を行い、当該縮尺係数が収束していない場合にのみ、上記の後処理工程を行うものとしてもよい。
このように、主工程後必ず後処理工程を行う場合も、主工程後の縮尺係数の収束の判定により後処理工程を行う場合も、何れも、垂線緯度の収束計算において、縮尺係数の収束を加味したものとなり、緯度の計算結果について必要な精度を確保することができる。尚、縮尺係数の収束を調べて後処理工程を行うか否かの判断を行う場合、無駄な反復計算(後処理工程)を行う必要がない。例えば、基準子午線付近では、主工程のみにて緯度及び経度の十分な精度を確保している場合がある。このような場合にまで後処理工程を行う必要がないからである。
但し、近年のコンピュータの性能の向上により、主工程後、必ず後処理工程を行うものとしても、実用上、計算に要する時間の増加や、ハードウエアの負荷による問題はない。
また、後処理工程は、1回に限らない。更に、上記と同様、後処理工程を行う毎に縮尺係数の収束の判定を行い、その結果により更に後処理工程を行うか否か判断するものとしても実施可能である。
この位置計算方法は、上記座標値x,yの入力を受け付ける。そして当該座標値xが与えられることにより、地球楕円体の赤道から当該座標値x,yが属する座標系原点までの子午線弧長に対し、当該座標値xを0以外の仮の縮尺係数mt で除した値を、加算して、赤道から当該座標値xまでの基準子午線上の子午線弧長Wとする。また、地球楕円体の長半径と第1離心率と緯度とから定まる子午線曲率半径について、地球楕円体の長半径aと第1離心率eと共に任意の数値(緯度φj )を初期値として与え、別途に子午線曲率半径を以って赤道から上記緯度φj まで積分(近似計算)することにより、赤道からの基準子午線上の子午線弧長S(φj )を算出する。そして、基準子午線上の上記漸近緯度φj から、式1を計算して基準子午線上の緯度φj+1 を漸近値(漸近緯度φj+1 )として算出する。
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 )…式1
上記式1にて算出した基準子午線上の緯度φj+1 と基準子午線上の緯度φj との差の絶対値と、基準値との比較を行い、上記基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値を超えると判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするには十分収束していないと判定した場合、この漸近緯度φj+1 を上記基準子午線上の緯度φj とし、基準子午線上の赤道からの子午線弧長S(φj )を算出すると共に当該子午線弧長S(φj )から上記式1の計算にて再度漸近緯度φj+1 を算出し、再度上記比較を行う。その結果、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値以内と判定した場合、即ち、この漸近緯度φj+1 (以下必要に応じて最近似垂線緯度と呼ぶ。)を、上記入力座標から基準子午線上に下ろした垂線とするのに十分収束したと判定した場合、当該漸近緯度φj+1 を、入力座標から基準子午線上に下ろした垂線の緯度φk として、仮の縮尺係数mt と、当該垂線緯度φk と、当該垂線緯度φk における子午線曲率半径Mk 及び卯酉線曲率半径Nk と、e’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、式3を計算して緯度φk+1 を算出する。
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 …式3
上記にて求めた緯度φk+1 を、当該(入力座標値に対応する)緯度φn として縮尺係数mt+1 を算出し、縮尺係数mt+1 と上記仮の縮尺係数mt との差の絶対値と、所定値との比較を行い(以上主工程)、当該縮尺係数mt+1 と上記仮の縮尺係数mt との差の絶対値が所定値を超えると判定した場合に、当該縮尺係数mt+1 を用いて求めた、基準子午線上の子午線弧長Wにて、上記式1の計算をやり直す(後処理工程)。
その結果得られた緯度を計算結果として出力する。また、緯度の計算の過程にて、経度及び子午線収差を計算する。
上記の式1をイテレーション式(反復式)として行う反復計算は、図3に示す平面座標において、基準子午線(X軸)上の漸近緯度φj+1 を、入力座標点(x,y)から基準子午線(X軸)に下ろした垂線hと、基準子午線(X軸)との交点Gに近似させる(近づける)計算である。基準子午線上において縮尺係数は0.9999である。式1による反復計算により、漸近緯度φj+1 を上記交点Gと取り扱うに十分近接した場合、当該漸近緯度φj+1 を最近似垂線緯度とし、式3の計算により、座標値(x,y)に対応する緯度φk+1 を算出するのである。
また、上記の仮の縮尺係数とは、上記直角座標の座標値yと、原点における縮尺係数m0 又はその代替値と、当該座標が属する座標系の原点における平均曲率半径R0 とから、式2の計算により、算出した縮尺係数mt である。
mt ={1+(y2 /2m02R02)+(y4 /24m04R04)}×0.9999…式2
上記の原点における平均曲率半径R0 は、当該座標が属する座標系の原点の緯度φ0 又はその代替値と共に、地球楕円体の長半径a及び短半径b、或いは当該半径a,bから導かれる第1離心率e、第2離心率e’を用いて定めたものである。
上記の仮の縮尺係数には、式2にて算出するものに限定するものではなく、0以外の任意の数値を採用することが可能である。但し、無駄な後処理工程を減ずる観点から、上記式2を用いて算出したものを採用するのが好ましい。
また上記において、式2を用いる場合も、原点における縮尺係数m0 の代替値として0以外の任意の数値を用いることが可能であるが、上記と同様の観点から、このような代替値ではなく、当該原点における縮尺係数m0 を用いて計算を行うのが好ましい。
また、上記の上記縮尺係数mt+1 は、原点縮尺係数m0 に代え上記縮尺係数mt を用い、且つ、原点における平均曲率半径R0 に代え式3にて算出した緯度φn における平均曲率半径Rφn (後述のRk2)を用いて、上記式2により算出するものである(言い換えると、式2において、mt をm0 に代入し、Rφn をR0 に代入して計算する)。
上記にて算出した縮尺係数mt+1 と、それまで各値を算出するのに用いた縮尺係数mt とについて、|mt+1 −mt |が収束しているか判定する。その結果、収束していれば、計算を終了し、式3、式4(式19)、式5及び式2を経て得た上記の緯度φn 、経度λ、子午線収差角γ、及び縮尺係数mt+1 を、目的の数値として出力する。
上記において|mt+1 −mt |が収束していないと判定した場合、後処理工程を行う。後処理工程において、上記主工程で得た縮尺係数mt+1 を用い、他の値は、(主工程で用いたものと)変更することなく、主工程と同様の処理を行う。当該後処理工程において|mt+1 −mt |が収束していなかった場合、当該後処理工程で得た縮尺係数mt+1 を用いて、更なる後処理工程を行う。
要するに、式2において、上記縮尺係数mt を原点縮尺係数m0 の代わりに用いて、式3にて計算した緯度φn における平均曲率半径Rφn (後述のRk2)を原点における平均曲率半径R0 の代わりに用いることにより、当該点の縮尺係数を算出する。
これを反復させると、当該点の縮尺係数は収束する性質を持ち、完全に収束した場合の縮尺係数は、当該点(入力座標)における縮尺係数である。
図1は、この装置の説明図である。図2は、その処理のフローを示す説明図である。
2−1)位置計算装置の概要
この装置は、オペレータから平面直角座標上の座標値x,yの入力を受付け、自動的に(出力まで人手を要することなく専ら装置のみにて)計算を行い、対応する地球楕円体(準拠楕円体を指す。ここではWGS84楕円体を採用する。但し、ベッセル、GRS80等の他の楕円体を採用しても実施可能である。)上の緯度及び経度を出力するものである。
この位置計算装置は、位置計算プログラムをコンピュータへ導入することにより、構築される。このコンピュータは、入力装置と、演算装置と、制御装置と、記憶装置(内部記憶装置及び外部記憶装置)とを備えた一般的なPC(パーソナルコンピュータ)やWS(ワークステーション)を採用して実施することができる。即ち、キーボード、マウスに代表される入力装置、ディスプレイやプリンタに代表される出力装置、CPU(中央上方処理装置)に代表されるされる演算装置や制御装置(演算制御装置)、RAM(ランダムアクセスメモリ)に代表される内部記憶装置、ハードディスクやフレキシブルディスクその他の磁気ディスク或いは光学ディスクに代表される外部記憶装置を備えた一般的なコンピュータを採用して実施することができる。
上記の位置計算プログラムは、コンピュータの資源を全て占有する専用のソフトウエアであっても実施可能であるが、MS−DOS(登録商標)、Windows(登録商標)やUNIX(登録商標)に代表される周知のOS(オペレーティングシステム)上で動作するものとしても実施可能である。
また、位置計算プログラムは、インタプリンタ型のものであっても、アセンブラやコンパイラを使用して生成されたものであっても何れでもよい。
この実施の形態において、位置計算プログラムは、OSが導入されたコンピュータ上に、導入されて位置計算装置を実現する。以下、OSと、位置計算プログラムとを位置計算ソフトウエアと呼ぶ(必要に応じて単にソフトと呼ぶ)。
この実施の形態に示す位置計算装置は、等角投影法により地球楕円体(準拠楕円体)を投影する平面直角座標X,Y上の座標値x,yを入力することにより、少なくとも当該座標値に対応する地球楕円体上の緯度φn 、経度λ、子午線収差γを出力するものである。
図1に示す通り、この装置は、データ置部10を有する座標入力部1と、出力部2と、定数設定部31と、初期値設定部32と、条件設定部33と、第1仮緯度収容部41と、第2仮緯度収容部42と、第1仮縮尺係数収容部51と、第2仮縮尺係数収容部52と、第1〜第7の7つの演算部61,62,63,64,65,66,67を有する演算部6と、演算制御部7と、主判定部81及び副判定部82を有する判定部8と、回帰制御部9とを備える。
この実施の形態において、座標入力部1は、入力されたデータを一時的に保持するデータ置部10を備える。また、定数設定部31は、算出値置部30を備える。
座標入力部1は、上記ソフトと、コンピュータが有する、上記の入力装置と、演算制御装置とにて構築される。出力部2は、上記ソフトと、コンピュータが有する上記出力装置と、演算制御装置とにて構築される。
定数設定部31、初期値設定部32及び条件設定部33の夫々は、上記ソフトと、コンピュータが有する、上記の演算制御装置、内部記憶装置及び外部記憶装置とにて構築される。
第1仮緯度収容部41、第2仮緯度収容部42、第1仮縮尺係数収容部51及び第2仮縮尺係数収容部52は、夫々、上記ソフトと、コンピュータが有する、上記の演算制御装置と、内部記憶装置とにて構築される。
第1〜第7の7つの演算部61,62,63,64,65,66,67と、演算制御部7と、主判定部81と、副判定部82と、回帰制御部9とは、夫々、上記ソフトと、コンピュータが有する、上記の演算制御装置と、内部記憶装置とにて構築される。
上記の座標入力部1は、オペレータからの上記座標値x,yの入力を受け付け。データ置部10に収容する。
定数設定部31は、上記地球楕円体の、長半径a、短半径b、第1離心率e、第2離心率e’の夫々、又はこれらに準じた値と、上記の座標x,yが属する座標系の原点の、緯度φ0 又はその代替値、経度λ0 の夫々と、当該原点における縮尺係数m0 又はその代替値と、後述する定数A〜Fとが、既知の数値として設定されている。
ここで設定とは、この装置を用いた位置計算の度にオペレータがこれらの値(定数)を上記入力装置を通じて入力するもの、オペレータが予め外部記憶装置にこれらの値(定数)を収容しておくもの、位置計算プログラムの一部としてプログラムが備えプログラムの導入によって自動的に外部記憶装置に収容されるものの何れも含む。
ここでは、必要な定数は、プログラムに含まれるものとして外部記憶装置に導入され、導入後オペレータによってコンピュータが備える上記入力装置を通じて、書き換え可能とする(書き換えられるのは、主として、上記の、原点の緯度φ0 又はその代替値と、経度λ0 と、原点の縮尺係数m0 又はその代替値である)。
定数設定部31は、算出値置部30を備える。算出値置部30は、演算部6が算出した値を収容する。
夫々の演算にて、処理される加減乗除及び、正弦、余弦、正接は、周知の方法によって、コンピュータを用いたソフトウエアによる情報処理(数値計算)にて行われる(近似計算も含む)。
演算1)
c=a2 /b=a/(1−e2 )1/2 =a(1+e'2)1/2 …式6
を、地球楕円体の長半径a、地球楕円体の短半径b、地球楕円体の第1離心率e、地球楕円体の第2離心率e’の各値から計算し、地球楕円体の極の曲率半径cを求めること。尚、上記第1及び第2の離心率は、地球楕円体の長半径aと短半径bとから算定することができる。また、この地球楕円体の極の曲率半径cは、定数として算出によらずに、定数設定部31が保持するものとしてもよい。
演算2)
V=(1+e'2cos2φ)1/2 …式7
上記式7を、緯度φと第2離心率e’の各値から計算し、Vを求めること。
演算3)
R=c/V2 …式8
上記式8を、地球楕円体の極の曲率半径cと、上記Vの値とから計算し、緯度φにおける平均曲率半径を求めること。
演算4)
mn+1 =(1+y2 /2mn 2 R2 +y4 /24mn 4 R4 )×0.9999
…式9
上記式9を、平面直角座標上の座標値x,yのyの値と、仮の縮尺係数mn とから計算し、縮尺係数mn+1 を求めること。
この演算4)は、特許請求の範囲に示す式2を処理する機能である(式9の縮尺係数mn を原点縮尺係数m0 とし、平均曲率半径Rを原点における平均曲率半径R0 としたものが、請求項2に示す式2である)。
演算5)
W=a(1−e2 )・{A・φ−(B/2)・(sin 2φ)
+(C/4)・(sin 4φ)−(D/6)・(sin 6φ)
+(E/8)・(sin 8φ)−(F/10)・(sin 10φ)}
+x/m …式10
上記式10を、第2演算部62は、地球楕円体の長半径aの値と、定数A〜Fの各値と、緯度φの値と、平面直角座標上の座標値x,yのxの値と、0以外の任意の実数mとから計算し、赤道から当該座標値xまでの基準子午線上の子午線弧長Wを求めること。尚、各定数A〜Fは、第1離心率eによって定まるものである(後述のステップ107参照)。
演算6)
S (φ )=a(1−e2 )・{A・φ−(B/2)・(sin 2φ)
+(C/4)・(sin 4φ)−(D/6)・(sin 6φ)
+(E/8)・(sin 8φ)−(F/10)・(sin10 φ)}
…式11
上記の式11を、地球楕円体の長半径aと、第1離心率eと、上記の各定数A〜Fと、緯度φとから計算し、赤道から緯度φまでの子午線弧長S (φ )を求めること。
演算7)
φn+1 =φn
+(W−S (φn ) )(1−e2 sin2φn )3/2 /a(1−e2 )
…式12
上記式12を、仮の垂線緯度φn と、赤道から当該座標値xまでの基準子午線上の子午線弧長Wと、仮の垂線緯度φn における赤道からの子午線弧長S (φn ) と、地球楕円体の長半径aと、第1離心率eの、各値から計算し垂線緯度φn+1 を求めること。
この式12は、特許請求の範囲に示す式1と対応する。
演算8)
t=tan φ …式13
上記式13を、与えられた緯度φ値から、計算し、tを求めること。
演算9)
η=e’cos φ …式14
上記式14を、地球楕円体の第2離心率e’と、緯度φの、夫々の値から、計算し、ηを求めること。
演算10)
M=c/V 3 …式15
上記式15を、地球楕円体の極の曲率半径cの値と、緯度φにて定まるVの値とから計算し、当該緯度φについての子午線曲率半径Mを求めること。
演算11)
N=c/V …式16
上記式16を、地球楕円体の極の曲率半径cの値と、緯度φにて定まるVの値とから計算し、当該緯度φについての卯酉線曲率半径Nを求めること。
演算12)
φn+1 =φn −tn (y/mt )2 /2Mn Nn +tn (y/mt )4 (5+3tn 2 +ηn 2 −9tn 2 ηn 2 −4ηn 4 )/24Mn Nn 3 −tn (y/mt )6 (61+90tn 2 +45tn 4 )/720 Mn Nn 5 …式17
上記式17を、平面直角座標上の座標値x,yのyの値、最近似垂線緯度φn 、縮尺係数mt 、最近似垂線緯度φn についてのtn 、最近似垂線緯度φn についてのηn 、最近似垂線緯度φn についての子午線曲率半径Mn 、及び、最近似垂線緯度φn についての卯酉線曲率半径Nn の、夫々の値から計算し、当該緯度φn+1 を求めること。
この式17は、特許請求の範囲に示す式3に対応する。
演算13)
Δλ=(y/mt )/Nn cos φn
−(y/mt )3 (1+2tn 2 +ηn 2 )/6Nn 3cosφn
+(y/mt )5 (5+28tn 2 +24tn 4 )/120Nn 5cosφn
…式18
上記式18を、平面直角座標上の座標値x,yのyの値、最近似垂線緯度φn 、縮尺係数mt 、最近似垂線緯度φn についてのtn 、最近似垂線緯度φn についてのηn 、及び、最近似垂線緯度φn についての卯酉線曲率半径Nn の、夫々の値から計算し、経度差分Δλを求めること。
この式18は、特許請求の範囲(請求項4)の式4に対応する。
演算14)
λ=λ0 +Δλ …式19
上記式19を、平面直角座標上の座標値x,yの属する原点の経度λ0 と、経度差分の夫々の値とから、計算し、経度λを求めること。
演算15)
γ=(y/mt )/Nn
−(y/mt )3 tn (1+tn 2 −ηn 2 )/3Nn 3
+(y/mt )5 tn (2+5tn 2 +3tn 4 )/15Nn 5
…式20
上記式20を、平面直角座標上の座標値x,yのyの値、縮尺係数mt 、最近似垂線緯度φn についてのtn 、最近似垂線緯度φn についてのηn 、及び、最近似垂線緯度φn についての卯酉線曲率半径Nn の、夫々の値から計算し、子午線収差γを求めること。
この式20は、特許請求の範囲(請求項5)の式5に対応する。
また、演算制御部7は、演算部61〜67の演算の結果得た緯度及び縮尺係数の収束状況について、判定部8に判定を行わせる。
判定部8は、垂線緯度の収束状態の判定を行う主判定部81と、縮尺係数の判定状態の判定を行う副判定部82とを備える。
回帰制御部9は、判定部8の判定を受けて、演算制御部7に、判定に対応した指示を行う。当該指示に従って、演算制御部7は、各演算部61〜67に指令を出す。即ち、回帰演算部9は、判定部8による判定の結果に従って、演算制御部7を通じて、各演算部61〜67を制御する。
条件設定部33は、判定部8の判定(比較)に必要な基準値が設定され、設定された基準値を保持する。
第2仮緯度収容部42は、第1仮緯度収容部41から、第1仮緯度収容部41が収容していた緯度の値を受け取り、当該緯度の値を収容する。第2仮緯度収容部42は、更に第1仮緯度収容部41から緯度の値を受け取ることにより、それまで収容していた緯度の値を受け取った後の緯度の値に書き換えて、当該後の緯度の値を収容する。
第2仮縮尺係数収容部52は、第1仮縮尺係数収容部51から、第1仮縮尺係数収容部51が収容していた縮尺係数の値を受け取り、当該縮尺係数の値を収容する。第2仮緯度収容部42は、更に第1仮縮尺係数収容部51から縮尺係数の値を受け取ることにより、それまで収容していた縮尺係数の値を受け取った後の縮尺係数の値に書き換えて、当該後の縮尺係数の値を収容する。
出力部2は、演算制御部7と回帰演算部9の指令を受けて、演算部6の演算結果を出力する。
以下、この装置を利用した緯度・経度・子午線収差の算出の流れについて、図2を用いて、説明する。
この計算方法は、主工程100を遂行し、更に縮尺係数の収束状況の判定結果により後処理工程を遂行するものである。
また、上記の主工程100の遂行に先立ち、事前に、前処理工程を行っておく。
前処理工程は、オペレータにより、位置計算前に行われる種々の設定工程である。
詳しくは、前処理工程において、オペレータにより、初期値設定部32にて緯度計算の初期値を決定するための数値(初期値決定用助数)が設定される。このような助数として、1”(1秒)をrad (ラジアン)に換算した値である 1/206265が、条件設定部33において設定される。但しこのような数値は変更可能である。
又、前処理工程において、オペレータにより、条件設定部33に、垂線緯度の収束条件IPS(しきい値)が設定される。
この収束条件IPSは、計算の精度を確保する上で重要である。即ち、演算部における各式にて的確な計算がなされても、収束条件IPSが広いと、十分な精度を確保するのが困難になるからである。
ここでは、収束条件IPSは、2”× 10-5/206265rad と設定する。
更に、前処理工程において、オペレータにより、初期値設定部32にて、縮尺係数の収束条件IPS2が設定される。ここでは、IPS2は、座標値の有効数字(最大9桁)を考慮して、4× 10-12 と設定してある。
但し、各収束条件IPS、IPS2の各数値は、上記の値に限定するものではなく、変更可能である。
主工程100以降の工程は、オペレータから座標値x,yを受け付けたこの位置計算装置が、位置計算を実際に行う工程である。
主工程は、次のステップ101〜125を順次遂行する。
オペレータは、座標入力部1により、直角座標の座標値x,y、及び当該座標値が属する座標系Zを入力する。座標系Zの入力については、日本に設定されている19の座標系のうち何れに属するかを特定するものである。
上記の入力作業は、オペレータが行う(これ以降の他の作業は、位置計算装置が行う。尚、この実施の形態において、後述する初期値設定部32の設定、及び定数設定部31内の各代替値の設定は、位置計算前に、オペレータにより行われる)。
座標入力部1は、入力された座標値x,yを、データ置部10に収容し、演算制御部7に座標値の入力があったことを通知する。
演算制御部7は、上記通知を受けて、第1演算部61に、定数設定部31を参照させる。即ち、定数設定部31は、各座標系の原点の緯度及び経度の数値を保持しており、第1演算部61は、入力された(座標値x,yが属する)座標系Zの、原点の緯度φ0 及び経度λ0 を特定する。これにて、以下、各演算部が定数設定部31にて参照する原点の緯度φ0 及び経度λ0 は、確定する。
このステップ102において、演算制御部7は、第1演算部61に、定数設定部31を参照させて、定数設定部31が保持している、地球楕円体の長半径aと短半径bから、或いは、地球楕円体の長半径aと第1離心率eから、或いは、地球楕円体の長半径aと第2離心率e’から、下記の計算1を行わせ、地球楕円体の極の曲率半径cを算出させる。第1演算部61は、算出した極の曲率半径cを、算出値置部30に収容する。
c=a2 /b=a/(1−e2 )1/2 =a(1+e'2)1/2 …計算1
計算1は、第1演算部61の演算1)の式6の計算機能にて実行される。
ステップ103において、演算制御部7は、第1演算部61に、定数設定部31を参照させて、第2離心率e’と当該座標が属する原点の緯度φ0 (以下単に原点緯度φ0 と呼ぶ。)とから、計算2を行わせ、V0 を算出させる。第1演算部61は、算出したVφ0 を、算出値置部30に収容する。
V0 =(1+e'2cos2φ0 )1/2 …計算2
計算2は、第1演算部61の演算2)の式7の計算機能にて実行される。
ステップ104において、演算制御部7は、第1演算部61に、定数設定部31(の算出値置部30)を参照させて、極の曲率半径cとV0 の値とから、計算3を行わせ、当該座標系原点の平均曲率半径R0 を算出させる。第1演算部61は、算出した平均曲率半径R0 を算出値置部30に収容する。
R0 =c/V02…計算3
計算3は、第1演算部61の演算3)の式8の計算機能にて実行される。
ステップ105において、演算制御部7は、定数設定部31の原点における縮尺係数m0 (以下原点縮尺係数m0 と呼ぶ。)を、第1仮縮尺係数収容部51に収容させる(定数設定部31の上記原点縮尺係数m0 を第1仮縮尺係数収容部51へコピーする)。この原点縮尺係数m0 は、0.9999という数値を採る。
そして、演算制御部7は、第1演算部61に、座標入力部1のデータ置部10及び定数設定部31(又は第1仮縮尺係数収容部51)を参照させて、入力された座標値yの値と、原点縮尺係数m0 と、上記平均曲率半径R0 とから、計算4を行わせ、縮尺係数mt を算出させる。
mt =(1+y2 /2m02R02+y4 /24m04R04)×0.9999…計算4
(右辺の0.9999は、数値として原点縮尺係数m0 と同じであるが、これは基準子午線上の縮尺係数を0.9999とするための一定の縮尺係数である。)
計算4は、第1演算部61の演算4)の式9の計算機能にて実行される。
第1演算部61は、算出した縮尺係数mt を第1仮縮尺係数収容部51に収容する。第1仮縮尺係数収容部51は、収容していた原点縮尺係数m0 を第2仮縮尺係数収容部52へ渡し、原点縮尺係数m0 に代え上記の通り上記縮尺係数mt を収容する。
ステップ106において、演算制御部7は、第1演算部61に、定数設定部31と、初期値設定部32とを参照させて、原点緯度φ0 から助数を減じた値、即ち、原点緯度φ0 − 1/206265を計算させる。第1演算部61は、算出した値を緯度φj として、第1仮緯度収容部41に収容する。
このステップ107において、演算制御部7は、第2演算部62に、定数設定部31と、座標入力部1のデータ置部10と、第1仮縮尺係数収容部51とを参照させて、長半径aと、第1離心率eと、原点緯度φ0 と、座標値xと、縮尺係数mt とから、計算5をさせ、赤道から当該座標値xまでの基準子午線上の子午線弧長Wとして算出させる。
W=a(1−e2 )・{A・φ0 −(B/2)・(sin 2φ0 )
+(C/4)・(sin 4φ0 )−(D/6)・(sin 6φ0 )
+(E/8)・(sin 8φ0 )−(F/10)・(sin10 φ0 )}
+x/mt …計算5
計算5の、A〜Fは、次の通りである。
A=1+(3/4)e2 +(45/64)e4 +(175 /256 )e6
+(11025 /16384 )e8 +(43659 /65536 )e10
B= (3/4)e2 +(15/16)e4 +(525 /512 )e6
+(2205/2048)e8 +(72765 /65536 )e10
C= (15/64)e4 +(105 /256 )e6
+(2205/4096)e8 +(10395 /16384 )e10
D= (35/512 )e6
+(315 /2048)e8 +(31185 /131072)e10
E= (315 /16384 )e8 +(3465/65536 )e10
F= (693 /131072)e10
計算5は、第2演算部62の演算5)の式10の計算機能にて実行される。
前述の通り、この定数A〜Fは、定数設定部31に予め収容されており、第2演算部62や他の演算部が参照して利用する。
ここで、計算5の右辺のx/mt 以外の部分が、地球楕円体の赤道から座標値x,yが属する座標系原点までの子午線弧長を示している。また、この縮尺係数mt は、特許請求の範囲でいう、仮の縮尺係数である。
第2演算部62は、算出した、赤道から当該座標値xまでの基準子午線上の子午線弧長Wを定数設定部31の算出値置部30に収容する。
このステップ108において、演算制御部7は、第3演算部63に、定数設定部31と、第1仮緯度収容部41とを参照させ、長半径a、第1離心率e、定数A〜Fと、緯度φj とから、計算6を行わせ、緯度φj に対する赤道からの子午線弧長S (φj ) を算出させる。
S (φj ) =a(1−e2 )・{A・φj −(B/2)・(sin 2φj )
+(C/4)・(sin 4φj )−(D/6)・(sin 6φj )
+(E/8)・(sin 8φj )−(F/10)・(sin10 φj )}
…計算6
即ち、このステップ108は、子午線曲率半径Mを以って赤道から垂線緯度φj まで積分するステップであり、上記の通り計算6は、第3演算部63にて当該積分の近似計算である(ここでいう積分とは、当該近似計算を指している)。
計算6は、第3演算部63の演算6)の式11の計算機能にて実行される。
第3演算部63は、上記にて算出した子午線弧長S (φj ) を、上記算出値置部30に収容する。
ステップ109において、演算制御部7は、第4演算部64に、定数設定部31と、第1仮緯度収容部41とを参照させ、緯度φj 、長半径a、第1離心率e、赤道から当該座標値xまでの基準子午線上の子午線弧長W、子午線弧長S (φj ) とから、計算7を行わせ、垂線緯度φj+1 を算出させる。
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 ) …計算7
計算7は、第4演算部64の演算7)の式12の計算機能にて実行される。
第4演算部64は、算出した緯度φj+1 を、第1仮緯度収容部41へ収容する。第4演算部64から上記の緯度φj+1 を受けた第1仮緯度収容部41は、収容していた緯度φj を第2仮緯度収容部42へ渡し、緯度φj に代え上記緯度φj+1 の値を収容する。
ステップ110において、判定部8の主判定部81は、第1仮緯度収容部41と第2仮緯度収容部42と条件設定部33とを参照し、緯度φj+1 と緯度φj の差の絶対値即ち、|φj+1 −φj |と、緯度の収束条件IPS(しきい値)とを比較する。そして、|φj+1 −φj |>IPSと判定した場合、判定部8(主判定部81)が、その結果を回帰制御部9に通知する。
この通知を受けた回帰制御部9は、第1仮緯度収容部41の緯度φj+1 の値を緯度φj の値として(第1仮緯度収容部41の緯度φj+1 を緯度φj に書き換えて)、装置の処理を再びステップ108に移す。そして、回帰制御部9は、ステップ108からステップ110を、判定部8(主判定部81)が|φj+1 −φj |≦IPSと判定するまで繰り返させる(垂線緯度の反復計算をさせる)。この反復計算において、垂線緯度φj+1 を垂線緯度φj に置き換える(書き換える)のみであり、演算部6(の各演算部)が参照する他の数値を変更しない(当該他の数値については、当初の値を用いて最初から計算される)。
以上により、ステップ110において、判定部8(の主判定部81)は、|φj+1 −φj |≦IPSと判定した場合、その結果を回帰制御部9に通知する。
回帰制御部9は、当該収束の通知を受けて、次のステップ111へ処理を移行させる。即ち、回帰制御部9は、演算制御部7に次のステップ111の処理を行わせる。
ステップ111において、演算制御部7は、第1演算部61に、定数設定部31と、第1仮緯度収容部41とを参照させて、第2離心率e’と垂線緯度φj+1 とから、計算8を行わせ、Vk を算出させる。尚、説明の便宜上、この垂線緯度φj+1 を緯度φk として説明する。
Vk =(1+e'2cos2φk )1/2 …計算8
計算8は、第1演算部61の演算2)の式7の計算機能にて実行される。
第1演算部61は、算出したVk を、算出値置部30に収容する。
ステップ112において、演算制御部7は、第5演算部65に、第1仮緯度収容部41を参照させ、緯度φk の値を用いて(即ち、緯度φj+1 の値を緯度φk の値として用いて)、計算9を行わせ、tk (説明の便宜上第2の値tk と呼ぶ。)を算出させる。
tk =tan φk …計算9
計算9は、第5演算部65の演算8)の式13の計算機能にて実行される。
第5演算部65は、算出した第2の値tk を、算出値置部30へ収容する。
ステップ113において、演算制御部7は、第5演算部65に、定数設定部31と第1仮緯度収容部41を参照させ、第2離心率e’と共に、緯度φk の値を用いて、計算10を行わせ、ηk (説明の便宜上第1の値ηk と呼ぶ。)を算出させる。
ηk =e’cos φk …計算10
計算10は第5演算部65の演算9)の式14の計算機能にて実行される。
第5演算部65は、算出した第1の値ηk を、算出値置部30へ収容する。
ステップ114において、演算制御部7は、第5演算部65に、定数設定部31を参照させ、極の曲率半径cとVk とから、計算11にて、緯度φk における子午線曲率半径Mk を算出させる。
Mk =c/Vk 3 …計算11
計算11は、第5演算部65の演算10)の式15の計算機能にて実行される。
第5演算部65は、算出した子午線曲率半径Mk の値を、算出値置部30へ収容する。
ステップ115において、演算制御部7は、第5演算部65に、定数設定部31を参照させ、極の曲率半径cとVk とから、計算12にて、緯度φk における卯酉線曲率半径Nk を算出させる。
Nk =c/Vk …計算12
計算12は、第5演算部65の演算11)の式16の計算機能にて実行される。
第5演算部65は、算出した卯酉線曲率半径Nk の値を、算出値置部30へ収容する。
ステップ116において、回帰制御部9は、演算部制御部7を通じて、第5演算部65に、緯度φk を与え、計算13を行わせる。具体的には、演算制御部7は、回帰制御部9からの指令を受けて、第5演算部65に、座標入力部1のデータ置部10、第1仮緯度収容部41、第1仮縮尺収容部51及び定数設定部31を参照させて、上記にて、第4演算部64に算出させた最近似垂線緯度φk と、当該点の縮尺係数mt と、最近似垂線緯度φk における子午線曲率半径Mk 及び卯酉線曲率半径Nk と、e’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、計算13の計算にて、緯度φk+1 を計算させる。
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 …計算13
計算13は、第5演算部65の演算12)の式17の計算機能にて実行される。
第5演算部65は、算出した緯度φk+1 を、当該点(入力座標点)の緯度φn として、第1仮緯度収容部41へ、それまで第1仮緯度収容部41に収容されていた値と置き換えて収容する。このとき、第1仮緯度収容部41は、それまで保持していた緯度φk の値を第2仮緯度収容部42へ渡す。
このステップ117は、上記のステップ116と平行として行われる。即ち、ステップ117は、ステップ116で得られた結果(緯度φk+1 )を用いるのではなく、ステップ116の処理に用いたデータを用いて、実行される。
ステップ117において、演算制御部7は、第6演算部66に、座標入力部1のデータ置部10、第1仮縮尺係数収容部51及び算出値置部30を参照させて、最近似垂線緯度φk (ステップ111に至るまでの一連の処理で得られた緯度、即ちステップ111〜116の処理のために用いた緯度φk )と、当該点の縮尺係数mt と、最近似垂線緯度φk における子午線曲率半径Mk (ステップ114の結果)及び卯酉線曲率半径Nk (ステップ115の結果)と、e’cos φk にて定まる第1の値ηk (ステップ113の結果)と、tan φk にて定まる第2の値tk (ステップ112の結果)と、上記座標値yとから、計算14を行わせ、経度差分Δλを算出させる。
Δλ=(y/mt )/Nk cos φk
−(y/mt )3 (1+2tk 2 +ηk 2 )/6Nk 3cosφk
+(y/mt )5 (5+28tk 2 +24tk 4 )/120Nk 5cosφk
…計算14
第6演算部66は、算出した経度差分Δλを、定数設定部31の算出値置部30に、収容する。
このステップ118と当該ステップ118に続くステップ119は、ステップ117と同様、上記のステップ116と平行として行われる。
ステップ118において、演算制御部7は、第6演算部66に、定数設定部31を参照させて、原点経度λ0 と経度差分Δλの和の計算を即ち、計算15をさせて、経度λを算出させる。
第6演算部66は、算出した経度λの値を算出値置部30に収容する。
λ=λ0 +Δλ …計算15
計算15は、第6演算部66の演算14)の式19の計算機能にて実行される。
このステップ119は、上記の通り、ステップ118に続く工程であり、ステップ116と平行して行われる。即ち、ステップ116で得られた結果(緯度φk+1 )を用いるのではなく、ステップ116の処理に用いたデータを用いて、実行される。
ステップ119において、演算制御部7は、第7演算部67に、座標入力部1のデータ置部10、第1仮縮尺係数収容部51及び算出値置部30を参照させて、最近似垂線緯度φk (ステップ111に至るまでの一連の処理で得られた緯度、即ちステップ111〜116の処理のために用いた緯度φk )と、当該点の縮尺係数mt と、最近似垂線緯度φk における卯酉線曲率半径Nk (ステップ115の結果)と、e’cos φk にて定まる第1の値ηk (ステップ113の結果)と、tan φk にて定まる第2の値tk (ステップ112の結果)と、上記座標値yとから、計算16を行わせ、子午線収差γを算出させる。
γ=(y/mt )/Nk
−(y/mt )3 tk (1+tk 2 −ηk 2 )/3Nk 3
+(y/mt )5 tk (2+5tk 2 +3tk 4 )/15Nk 5
…計算16
計算16は、第7演算部67の演算15)の式20の計算機能にて実行される。
第7演算部67は、算出した子午線収差γを算出値置部30に収容する。
ステップ120は、ステップ116に続く工程である。
ステップ120において、演算制御部7は、第1演算部61に、定数設定部31と、第1仮緯度収容部41とを参照させて、第2離心率e’とステップ116で得た当該点の緯度φk+1 とから、計算17を行わせ、Vk2を算出させる。尚、説明の便宜上、この緯度φk+1 を緯度φk2として説明する。
Vk2=(1+e'2cos2φk2)1/2 …計算17
計算17は、第1演算部61の演算2)の式7の計算機能にて実行される。
第1演算部61は、算出したVk2を、算出値置部30に収容する。
ステップ121において、演算制御部7は、第5演算部65に、定数設定部31を参照させ、極の曲率半径cとVk2とから、計算18にて、緯度φk2における子午線曲率半径Mk2を算出させる。
Mk2=c/Vk23 …計算18
計算18は、第5演算部65の演算10)の式15の計算機能にて実行される。
第5演算部65は、算出した子午線曲率半径Mk2の値を、算出値置部30へ収容する。
ステップ122において、演算制御部7は、第5演算部65に、定数設定部31を参照させ、極の曲率半径cとVk2とから、計算19にて、緯度φk2における卯酉線曲率半径Nk2を算出させる。
Nk2=c/Vk2 …計算19
計算19は、第5演算部65の演算11)の式16の計算機能にて実行される。
第5演算部65は、算出した卯酉線曲率半径Nk2の値を、算出値置部30へ収容する。
ステップ123において、演算制御部7は、第1演算部61に、定数設定部31を参照させて、極の曲率半径cとVk2の値とから、計算20を行わせ、当該座標(入力座標)の平均曲率半径Rk2を算出させる。第1演算部61は、算出した平均曲率半径Rk2を算出値置部30に収容する。
Rk2=c/Vk22 …計算20
計算20は、第1演算部61の演算3)の式8の計算機能にて実行される。
ステップ124において、演算制御部7は、第1演算部61に、座標入力部1のデータ置部10及び第1仮縮尺係数収容部51を参照させて、入力された座標値yの値と、縮尺係数mt1(第1仮縮尺係数収容部51に収容されている前記縮尺係数mt を縮尺係数mt1として説明する。以下同じ。)と、上記平均曲率半径Rk2とから、計算21を行わせ、縮尺係数mt2を算出させる。
mt2=(1+y2 /2mt12 Rk22 +y4 /24mt14 Rk24 )×0.9999
…計算21
計算21は、第1演算部61の演算4)の式9の計算機能にて実行される。
第1演算部61は、算出した縮尺係数mt2を第1仮縮尺係数収容部51に収容する。第1仮縮尺係数収容部51は、収容していた原点縮尺係数mt を第2仮縮尺係数収容部52へ渡し、上記の通り縮尺係数mt に代え上記縮尺係数mt2を収容する。
ステップ125において、判定部8の副判定部82は、第1仮縮尺係数収容部51と、第2仮縮尺係数収容部52と、条件設定部33とを参照し、縮尺係数mt1と縮尺係数mt2の差の絶対値即ち|mt2−mt1|と、縮尺係数の収束条件IPS2とを比較する。そして、縮尺係数mt2(と縮尺係数mt1の差の絶対値)が収束したと判定された場合、即ち、副判定部82が|mt2−mt1|≦IPS2と判定した場合、判定部8は、その結果を回帰制御部9に通知する。
上記ステップ125にて、縮尺係数の収束の通知を受けた回帰制御部9は、出力部2により、少なくとも、第1仮緯度収容部41が収容している当該点の緯度φk+1 の値(φn ) と、算出値置部30が収容している経度λの値及び子午線収差γの値とを、出力させる。このとき、必要に応じて、座標入力部1のデータ置部10の入力座標値x,yや、定数設定部31が保持しその算出値置部30が収容している、上記入力座標が属する原点の緯度や経度の数値、地球楕円体の長半径や短半径などの定数、或いは計算工程中に算出した中間の値を、併せて、出力部2に出力させるものとしても、上記の結果がどのような座標値に対応するか、また、どのような経緯を経て算出されたか、オペレータに確認させることができ、便利である。
即ち、上記にて、|mt2−mt1|>IPS2と判定された場合、判定部8は、その結果を回帰制御部9に通知する。この通知を受けた回帰制御部9は、演算制御部7を通じて、上記主工程100の各ステップ(ステップ102〜125)を繰り返させる。この後処理工程において、主工程100と異なるのは、ステップ105の処理において、前述の原点縮尺係数m0 に代え、第1仮縮尺係数収容部51に収容されている縮尺係数mt2を第1演算部61に参照させて、計算4の計算を行わせるということである(原点縮尺係数m0 として、上記主工程100を経て得た縮尺係数mt2を用いる。即ち、後処理工程で原点縮尺係数m0 として用いる値は、上記主工程100を経て得た縮尺係数mt2である)。他の値については、上記主工程100と同様のものを用いる。即ち、後処理工程は、縮尺係数mt2を用いて、主工程100をやり直す工程である。従って、後処理工程において、主工程100の結果として加味されるのは、主工程100で算出された縮尺係数のみであり、主工程100で算出した緯度等の値は用いず、再度主工程100で用いたのと同様の値を用いて計算する(縮尺係数を使わない計算や、縮尺係数を用いる計算であっても、縮尺係数以外の数値については、主工程100で用いたものと同様の数値を用いて計算を行う)。
この後処理工程でのステップ125において、|mt2−mt1|>IPS2と判定された場合、判定部8は、その結果を回帰制御部9に通知する。そして、当該後処理工程を経て得た縮尺係数mt2を原点縮尺係数m0 として、再度後処理工程(主工程100のステップ102からのやり直し)を行う。このように、1回の後処理工程を経ても、縮尺係数mt2(と縮尺係数mt1の差の絶対値)が収束したと判定されない場合、収束するまで、即ち、副判定部82が|mt2−mt1|≦IPS2となるまで、上記の後処理工程を繰り返す。
後処理工程におけるステップ125にて、縮尺係数mt2が収束した即ち|mt2−mt1|≦IPS2と副判定部82に判定された場合、上記の通り、回帰制御部9は、演算制御部7に、ステップ126を遂行させる。
上記の実施の形態において、後処理工程は、判定部8(副判定部82)による縮尺係数の収束不十分、即ち、|mt2−mt1|>IPS2の判定により、行うものとした。
この他、前述の通り、このような縮尺係数の収束の判定を行わずに、設定した回数必ず後処理工程を実行するものとしても実施可能である。例えば、後処理工程を必ず1回(或いは2回以上)行うものとして実施することが可能である。
また、上記の実施の形態において、ステップ107で用いる仮の数値として、ステップ105で算出した縮尺係数mt を用いたが、初回のステップ107(後処理工程でのステップ107は除く。)については、0以外の任意の実数を当該縮尺係数mt の代替値として用いて計算を行うものとしても実施可能である。
但し、ステップ107で用いる仮の縮尺係数として、ステップ105で算出した縮尺係数mt を用いることにより、収束が効率良く行われるので、例えば、ステップ110における判定の結果ステップ108への反復(小ループ)を3回程度に抑えることができ、また、ステップ125における後処理工程の反復(大ループ)を2回以内に抑えることができる。従って、この場合、後処理工程をステップ125の判定無しに2回に固定しても、事実上良好な結果が得られる。
また、ステップ106において、原点緯度φ0 に対して、1秒(−1/206265ラジアン)ずらすものとした。このずらす値については、原点緯度φ0 に対して±5秒を超えない範囲で行うのが好ましい。この範囲であれば、上記の通り、主工程100のステップ110の判定による反復計算(ステップ108〜110)を最大3回とし後処理工程を1回として十分な精度が得られることが発明者によって確認されている。従って、このステップ108〜110の反復計算を3回を限度として実施するのが好ましい。但し、反復計算は、収束条件IPSに依存しているため、このような回数に限定するものではない。
従って、このような逆算によって、本願発明に係る方法で得られた計算結果の精度を検証することができる。
本願発明に係る上記の方法にて得られた結果(出力値)を緯度φ(緯度φn )、経度λとして、下記の計算を行う。
計算22にて、経度差分Δλを逆算する。
Δλ=λ−λ0 …計算22
上記結果の緯度φから、
計算23にて、第1の値ηを求める。
η=e’cos φ…計算23
また、上記の結果の緯度φから、計算24にて、第2の値tを求める。
t=tan φ…計算24
更に、上記の結果の緯度φから、計算25にて、Vの値を求める。
V=(1+e'2cos2φ)1/2 …計算25
求めた当該Vの値から、計算26にて子午線曲率半径Mを、計算27にて、卯酉線曲率半径Nを、更に計算28から平均曲率半径Rを、夫々求める。
M=c/V3 …計算26
N=c/V…計算27
R=c/V2 …計算28
上記結果の緯度φ及び経度λ、更に、上記にて求めた第1の値η、第2の値tとから、計算29にて、当該緯度φ及び経度λにおける縮尺係数mを求める。
m={(1+Δλ2cos2 φ(1+η2 )/2+Δλ4cos4 φ(5−4t2 )/24)}×0.9999…計算29
そして、上記の結果の緯度φと、当該φについての上記子午線曲率半径M、卯酉線曲率半径N、平均曲率半径R、上記の経度差分Δλ、第1の値η、第2の値t、縮尺係数mとから、計算30にて、当該緯度φ及び経度λに対応する平面直角座標上の座標値yを求めることができる(計算29の右辺の0.9999は、数値として原点縮尺係数m0 と同じであるが、前述の通り、これは基準子午線上の縮尺係数を0.9999とするための一定の縮尺係数である)。
y={ΔλNcos φ+Δλ3 Ncos 3 φ(1−t2 +η2 )/6
+Δλ5 Ncos 5 φ(5−18t2 +t4 +14η2 −58η2 t2 )/120 }×m…計算30
S (φ )=a(1−e2 )・{A・φ−(B/2)・(sin 2φ)
+(C/4)・(sin 4φ)−(D/6)・(sin 6φ)
+(E/8)・(sin 8φ)−(F/10)・(sin10 φ)}
…計算31
上記で求めた緯度φに対する赤道からの子午線弧長S (φ )と、原点緯度φ0 に対する赤道からの子午線弧長S (φ0 ) とにより、計算32にて、その差Bを求める。
B=S (φ )−S (φ0 ) …計算32
更に、この差Bと、上記の結果の緯度φと、当該φについての上記子午線曲率半径M、卯酉線曲率半径N、平均曲率半径R、上記の経度差分Δλ、第1の値η、第2の値t、縮尺係数mとから、計算33にて、にて、当該緯度φ及び経度λに対応する平面直角座標上の座標値xを求めることができる。
x={B+Δλ2sinφcos φ/2+Δλ4 Nsin φcos3φ(5−t2 +9η2 +4η4 )/24+Δλ6 Nsin φcos 5 φ(61−58t2 +t4 + 270η2 −330 η2 t2 )/720 }×m…計算33
また、上記の結果の緯度φと、上記の経度差分Δλ、第1の値η、第2の値tとから、計算34にて、当該緯度φ及び経度λにおける子午線収差γを求めることができる。
γ=Δλsin φ+Δλ3 sin φcos2φ(1+3η2 +2η4 )/3
+Δλ5 sin φcos4φ(2−t2 )/15…計算34
更に、具体的な実施例について、以下に掲げる。
結果データ1〜5で得た緯度・経度から、上記の検算法(検算プログラム)にて、座標値についての検証データ1〜5を得た(検証データの番号は、検証した結果データの番号と対応している)。
各データの計算において、WGS84地球楕円体の長半径a、地球楕円体の極の曲率半径c(=a2 /b)、離心率(第1離心率)eの2乗、第2離心率e’の2乗、座標系Z、当該座標系Zの原点緯度φ0 、当該座標系Zの原点経度λ0 、原点縮尺係数m0 について、次の数値を定数として用いた。
a=6378137.000
c=6399593.62576
e2 = 0.006694379990141
e’2 =0.006739496742276
Z=3
φ0 =36- 0- 0.000
λ0 = 132-10- 0.000
m0 =0.9999
また、以下において、
ZAHYOCHI X( 1)は、入力座標値xを示し、単位はメートルである。
ZAHYOCHI Y( 1)は、入力座標値yを示し、単位はメートルである。
PHAIは、最近似垂線緯度φj+1 を示す。
LATITUDEは、当該点の緯度φn (北緯)を示す。例えば 34- 7-45.93167 は、 34 度 7分45.93167秒(34°7'45".93167 )を示している。
LONGITUDE は、当該点の経度λ(東経)を示す。
GAMMA ANGLEは、当該点の子午線収差(真北方向角)γを示す。
MERIDIAN RADIUSは、当該点の子午線曲率半径Mを示す。
ORTHOG. MERIDIAN RADIUS は、当該点の卯酉線曲率半径Nを示す。
MEAN RADIUS は、当該点の平均曲率半径Rを示す。
SCALE FACTOR は、当該点の縮尺係数を示す。
1 ZAHYOCHI X( 1) = -207504.143
ZAHYOCHI Y( 1) = 0.000
PHAI = 34- 7-45.93167
LATITUDE = 34- 7-45.93167
LONGITUDE =132-10- 0.00000
GAMMA ANGLE = 0- 0- 0.00000
MERIDIAN RADIUS = 6355518.63673
ORTHOG. MERIDIAN RADIUS =6384868.10519
MEAN RADIUS = 6370176.46816
SCALE FACTOR = 0.9999000000
2 ZAHYOCHI X( 2) = -207504.143
ZAHYOCHI Y( 2) = -60000.000
PHAI = 34- 7-46.23045
LATITUDE = 34- 7-40.02854
LONGITUDE = 131-30-58.32910
GAMMA ANGLE = - 0-21-53.81142
MERIDIAN RADIUS = 6355516.93662
ORTHOG. MERIDIAN RADIUS =6384867.53587
MEAN RADIUS = 6370175.33215
SCALE FACTOR = 0.9999443586
3 ZAHYOCHI X( 3) = -207504.143
ZAHYOCHI Y( 3) = -90091.612
PHAI = 34- 7-46.60518
LATITUDE = 34-7-32.62478
LONGITUDE = 131-11-24.23440
GAMMA ANGLE = -0-32-52.51513
MERIDIAN RADIUS = 6355514.80440
ORTHODG. MERIDIAN RADIUS = 6384866.82185
MEAN RADIUS = 6370173.90739
SCALE FACTOR = 1.0000000000
4 ZAHYOCHI X( 4) = -207504.143
ZAHYOCHI Y( 4) = -119361.021
PHAI = 34-7-47.111365
LATITUDE = 34-7-22.57909
LONGITUDE = 130-52-22.58341
GAMMA ANGLE = -0-43-32.98171
MERIDIAN RADIUS = 6355511.91142
ORTHOG. MERIDIAN RADIUS = 6384865.85306
MEAN RADIUS = 6370171.97429
SCALE FACTOR = 1.0000755078
5 ZAHYOCHI X( 5) = -207504.143
ZAHYOCHI Y( 5) = -127420.415
PHAI = 34-7-47.27856
LATITUDE = 34-7-19.32099
LONGITUDE = 130-47-8.30603
GAMMA ANGLE = -0-46-29.28380
MERIDIAN RADIUS = 6355510.97317
ORTHOG. MERIDIAN RADIUS = 6384865.53887
MEAN RADIUS = 6370171.34734
SCALE FACTOR = 1.0001000000
1 ZAHYOCHI X( 1) = -207504.143
ZAHYOCHI Y( 1) = 0.000
LATITUDE = 34-7-45.93167
LONGITUDE = 132-10-0.00000
GAMMA ANGLE = 0-0-0.00000
MERIDIAN RADIUS = 6355518.63672
ORTHODG. MERIDIAN RADIUS =6384868.10519
MEAN RADIUS = 6370176.46816
SCALE FACTOR = 0.9999000000
2 ZAHYOCHI X( 2) = -207504.143
ZAHYOCHI Y( 2) = -60000.000
LATITUDE = 34-7-40.02854
LONGITUDE = 131-30-58.32910
GAMMA ANGLE = -0-21-53.81142
MERIDIAN RADIUS = 6355516.93662
ORTHODG. MERIDIAN RADIUS =6384867.53587
MEAN RADIUS = 6370175.33215
SCALE FACTOR = 0.9999443586
3 ZAHYOCHI X( 3) = -207504.143
ZAHYOCHI Y( 3) = -90091.612
LATITUDE = 34-7-32.62478
LONGITUDE = 131-11-24.23440
GAMMA ANGLE = -0-32-52.51514
MERIDIAN RADIUS = 6355514.80440
ORTHODG. MERIDIAN RADIUS = 6384866.82185
MEAN RADIUS = 6370173.90739
SCALE FACTOR = 0.9999999999
4 ZAHYOCHI X( 4) = -207504.143
ZAHYOCHI Y( 4) = -119361.021
LATITUDE = 34-7-22.57909
LONGITUDE = 130-52-22.58341
GAMMA ANGLE = -0-43-32.98172
MERIDIAN RADIUS = 6355511.91142
ORTHOG. MERIDIAN RADIUS = 6384865.85306
MEAN RADIUS = 6370171.97429
SCALE FACTOR = 1.0000755076
5 ZAHYOCHI X( 5) = -207504.143
ZAHYOCHI Y( 5) = -127420.415
LATITUDE = 34-7-19.32099
LONGITUDE = 130-47-8.30603
GAMMA ANGLE = -0-46-29.28381
MERIDIAN RADIUS = 6355510.97317
ORTHOG. MERIDIAN RADIUS= 6384865.53887
MEAN RADIUS = 6370171.34734
SCALE FACTOR = 1.0000999998
また、結果データ5を見れば、縮尺係数を示す SCALE FACTORも、1.0001に、精度よく収束しているのが分かる。
2 主力部
6 演算部
7 演算制御部
8 判定部
9 回帰制御部
Claims (9)
- 等角投影法により地球楕円体を投影する平面直角座標X,Y上の座標値x,yを入力することにより、少なくとも当該座標値に対応する地球楕円体上の緯度φn を出力することが可能な位置計算装置において、
座標入力部と、出力部と、演算部と、判定部と、回帰制御部とを備えるものであり、
座標入力部は、上記座標値x,yの入力を受け付け、
演算部は、直角座標の座標値xが与えられることにより、地球楕円体の赤道から当該座標値x,yが属する座標系原点までの子午線弧長に対し、当該座標値xを仮の縮尺係数mt で除した値を、加算して、赤道から当該座標値xまでの基準子午線上の子午線弧長Wとし、
この子午線弧長Wに対応する基準子午線上の緯度、即ち、上記入力座標から基準子午線上に下ろした垂線の緯度を算出するため、
演算部は、地球楕円体の長半径と第1離心率と緯度とから定まる子午線曲率半径について、地球楕円体の長半径aと第1離心率eと共に基準子午線上の任意の緯度φj を初期値として与え、別途に子午線曲率半径を以って赤道から上記緯度φj まで積分することにより、赤道からの基準子午線上の子午線弧長S(φj )を算出し、
演算部は、基準子午線上の緯度φj を与えることによって、式1を計算して基準子午線上の緯度φj+1 を漸近値として算出し、
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 )…式1
判定部は、上記式1にて算出した基準子午線上の上記漸近緯度φj+1 と基準子午線上の緯度φj との差の絶対値と、基準値との比較を行うものであり、
判定部にて、基準子午線上の、上記漸近緯度φj+1 と緯度φj との差の絶対値が基準値を超えると判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするには十分収束していないと判定した場合、回帰制御部は、漸近緯度φj+1 を上記基準子午線上の緯度φj とし、演算部に基準子午線上の赤道からの子午線弧長S(φj )を算出させると共に当該子午線弧長S(φj )から上記式1の計算にて再度漸近緯度φj+1 を算出させ、判定部に再度上記比較を行わせるものであり、
判定部が、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値を基準値以内と判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするのに十分収束したと判定した場合、回帰制御部は、当該漸近緯度φj+1 を、入力座標から基準子午線上に下ろした垂線の緯度φk として、仮の縮尺係数mt と、当該垂線緯度φk と、当該垂線緯度φk における子午線曲率半径Mk 及び卯酉線曲率半径Nk と、e’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、式3を計算して緯度φk+1 を演算部に算出させ、
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 …式3
回帰制御部は、上記緯度φk+1 を入力座標に対応する緯度φn として、演算部に、縮尺係数mt+1 を算出させ、前記仮の縮尺係数mt に代え当該縮尺係数mt+1 を用いて求めた、基準子午線上の子午線弧長Wにて、上記式1の計算を行わせるものであることを特徴とする位置計算装置。 - 上記の仮の縮尺係数とは、演算部に、直角座標の座標値yを与えることにより、原点における縮尺係数m0 又はその代替値と、当該座標が属する座標系の原点における平均曲率半径R0 とから、式2の計算により、演算部が算出した縮尺係数mt であり、
mt ={1+(y2 /2m02R02)+(y4 /24m04R04)}×0.9999…式2
上記の原点における平均曲率半径R0 は、当該座標が属する座標系の原点の緯度φ0 又はその代替値と共に、地球楕円体の長半径a及び短半径b、或いは当該半径a,bから導かれる第1離心率e、第2離心率e’を用いて定めたものであることを特徴とする請求項1記載の位置計算装置。 - 上記縮尺係数mt+1 は、原点縮尺係数m0 に代え上記縮尺係数mt を用い、且つ、原点における平均曲率半径R0 に代え演算部に式3にて算出させた緯度φn における平均曲率半径Rφn を用いて、演算部が、上記式2により算出するものであることを特徴とする請求項2記載の位置計算装置。
- 上記にあって、判定部が、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値以内と判定した場合に、回帰制御部は、上記緯度φn に対応する縮尺係数mt+1 と上記仮の縮尺係数mt との差の絶対値と、所定値との比較を、更に判定部にさせ、判定部が当該縮尺係数mt+1 と上記仮の縮尺数値mt との差の絶対値が所定値を超えると判定した場合にのみ、上記の通り、前記の仮の縮尺係数mt に代え、当該縮尺係数mt+1 を用いて求めた基準子午線上の子午線弧長Wにて、上記式1の計算を演算部に行わせるものであることを特徴とする請求項1乃至3の何れかに記載の位置計算装置。
- 上記にあって、判定部が、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値以内と判定した場合、回帰制御部は、垂線緯度φk として上記漸近緯度φj+1 を演算部に与えるものであり、演算部は、上記の仮の縮尺係数mt と、上記垂線緯度φk における卯酉線曲率半径Nk と、当該垂線緯度φk と第2離心率とによりe’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、式4にて、経度差分Δλを算出することが可能であり、
Δλ=(y/mt )/Nk cos φk
−(y/mt )3 (1+2tk 2 +ηk 2 )/6Nk 3cosφk
+(y/mt )5 (5+28tk 2 +24tk 4 )/120Nk 5cosφk …式4
演算部は、当該座標系が属する原点における経度λ0 に対し当該経度差分Δλを加算することにより、上記座標値x,yにおける経度λを算出するものであることを特徴とする請求項1乃至4の何れかに記載の位置計算装置。 - 上記にあって、判定部が、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値を基準値以内と判定した場合、回帰制御部は、垂線緯度φk としてこの漸近緯度φj+1 を演算部に与えるものであり、演算部は、上記の縮尺係数mt と、当該垂線緯度φk における、子午線曲率半径Mk 及び卯酉線曲率半径Nk と、当該垂線緯度φk と第2離心率とによりe’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、
γ=(y/mt )tk /Nk
−(y/mt )3 tk (1+tk 2 −ηk 2 )/3Nk 3
+(y/mt )5 tk (2+5tk 2 +3tk 4 )/15Nk 5 …式5
式5にて、座標x,yにおける子午線収差γを算出することが可能なことを特徴とする請求項1乃至5の何れかに記載の位置計算装置。 - 等角投影法により地球楕円体を投影する平面直角座標X,Y上の座標値x,yを入力することにより、少なくとも当該座標値に対応する地球楕円体上の緯度φn を出力することが可能な位置計算プログラムにおいて、
座標入力部と、出力部と、演算部と、判定部と、回帰制御部とを備えたコンピュータに導入されるものであり、
座標入力部にて、上記座標値x,yの入力を受け付け、
演算部に、直角座標の座標値xを与えて、地球楕円体の赤道から当該座標値x,yが属する座標系原点までの子午線弧長に対し、当該座標値xを仮の縮尺係数mt で除した値を加算して、赤道から当該座標値xまでの基準子午線上の子午線弧長Wとし、
この子午線弧長Wに対応する基準子午線上の緯度、即ち、上記入力座標から基準子午線上に下ろした垂線の緯度を算出するため、
演算部にて、地球楕円体の長半径と第1離心率と緯度とから定まる子午線曲率半径について、地球楕円体の長半径aと第1離心率eと共に任意の緯度φj を初期値として与え、別途に子午線曲率半径を以って赤道から上記緯度φj まで積分することにより、赤道からの基準子午線上の子午線弧長S(φj )を算出し、
演算部に、基準子午線上の緯度φj を与えることによって、式1を計算して基準子午線上の緯度φj+1 を漸近値として算出させ、
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 )…式1
判定部に、上記式1にて算出した基準子午線上の緯度φj+1 と基準子午線上の緯度φj との差の絶対値と、基準値との比較を行わせ、
判定部にて、基準子午線上の、上記漸近緯度φj+1 と緯度φj との差の絶対値が基準値を超えると判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするには十分収束していないと判定した場合、回帰制御部により、漸近緯度φj+1 を上記基準子午線上の緯度φj とし、演算部に基準子午線上の赤道からの子午線弧長S(φj )を算出させると共に当該子午線弧長S(φj )から上記の式1の計算にて再度漸近緯度φj+1 を算出させ、判定部に再度上記比較を行わせるものであり、
判定部にて、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値以内と判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするのに十分収束したと判定した場合、
当該漸近緯度φj+1 を、入力座標から基準子午線上に下ろした垂線の緯度φk として、回帰制御部により、仮の縮尺係数mt と、当該垂線緯度φk と、当該垂線緯度φk における子午線曲率半径Mk 及び卯酉線曲率半径Nk と、e’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、式3を計算して緯度φk+1 を演算部に算出させ、
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 …式3
回帰制御部により、上記緯度φk+1 を入力座標に対応する緯度φn として、演算部に、縮尺係数mt+1 を算出させ、前記仮の縮尺係数mt に代え当該縮尺係数mt+1 を用いて求めた、基準子午線上の子午線弧長Wにて、上記式1の計算を行わせるものであることを特徴とする位置計算プログラム。 - 請求項7記載の位置計算プログラムの記録媒体。
- 等角投影法により地球楕円体を投影する平面直角座標X,Y上の座標値x,yを入力することにより、少なくとも当該座標値に対応する地球楕円体上の緯度φn を出力することが可能な位置計算方法において、
座標入力部と、出力部と、演算部と、判定部と、回帰制御部とを備えたコンピュータを用いるものであり、
座標入力部にて、上記座標値x,yの入力を受け付け、
演算部に、直角座標の座標値xを与えて、地球楕円体の赤道から当該座標値x,yが属する座標系原点までの子午線弧長に対し、当該座標値xを仮の縮尺係数mt で除した値を加算して、赤道から当該座標値xまでの基準子午線上の子午線弧長Wとし、
演算部にて、地球楕円体の長半径と第1離心率と緯度とから定まる子午線曲率半径について、地球楕円体の長半径aと第1離心率eと共に任意の緯度φj を初期値として与え、別途に子午線曲率半径を以って赤道から上記緯度φj まで積分することにより、赤道からの基準子午線上の子午線弧長S(φj )を算出し、
演算部に、基準子午線上の緯度φj を与えることによって、式1を計算して基準子午線上の緯度φj+1 を漸近値として算出させ、
φj+1 =φj +(W−S (φj ) )(1−e2 sin2φj )3/2 /a(1−e2 )…式1
判定部にて、上記の式1にて算出した基準子午線上の上記漸近緯度φj+1 と基準子午線上の緯度φj との差の絶対値と、基準値との比較を行わせるものであり、
判定部にて、上記基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値を超えると判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするには十分収束していないと判定した場合、回帰制御部により、漸近緯度φj+1 を上記基準子午線上の緯度φj とし、演算部に基準子午線上の赤道からの子午線弧長S(φj )を算出させると共に当該子午線弧長S(φj )から上記の式1の計算にて再度漸近緯度φj+1 を算出させ、判定部に再度上記比較を行わせるものであり、
判定部にて、基準子午線上の、漸近緯度φj+1 と緯度φj との差の絶対値が基準値以内と判定した場合、即ち、この漸近緯度φj+1 を、上記入力座標から基準子午線上に下ろした垂線の緯度とするのに十分収束したと判定した場合、回帰制御部は、当該漸近緯度φj+1 を、入力座標から基準子午線上に下ろした垂線の緯度φk として、仮の縮尺係数mt と、当該垂線緯度φk と、当該垂線緯度φk における子午線曲率半径Mk 及び卯酉線曲率半径Nk と、e’cos φk にて定まる第1の値ηk と、tan φk にて定まる第2の値tk と、上記座標値yとから、式3を計算して緯度φk+1 を演算部に算出させ、
φk+1 =φk −tk (y/mt )2 /2Mk Nk +tk (y/mt )4 (5+3tk 2 +ηk 2 −9tk 2 ηk 2 −4ηk 4 )/24Mk Nk 3 −tk (y/mt )6 (61+90tk 2 +45tk 4 )/720 Mk Nk 5 …式3
回帰制御部により、上記緯度φk+1 を入力座標に対応する緯度φn として、演算部に、縮尺係数mt+1 を算出させ、前記仮の縮尺係数mt に代え当該縮尺係数mt+1 を用いて求めた、基準子午線上の子午線弧長Wにて、上記式1の計算を行わせるものであることを特徴とする位置計算方法。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2004246064A JP4610968B2 (ja) | 2004-08-26 | 2004-08-26 | 位置計算装置、位置計算方法、位置計算プログラム、及び位置計算プログラムの記憶媒体 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2004246064A JP4610968B2 (ja) | 2004-08-26 | 2004-08-26 | 位置計算装置、位置計算方法、位置計算プログラム、及び位置計算プログラムの記憶媒体 |
Related Child Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2010181971A Division JP5132731B2 (ja) | 2010-08-16 | 2010-08-16 | 位置計算装置、位置計算方法、位置計算プログラム、及び位置計算プログラムの記憶媒体 |
Publications (3)
Publication Number | Publication Date |
---|---|
JP2006064892A true JP2006064892A (ja) | 2006-03-09 |
JP2006064892A5 JP2006064892A5 (ja) | 2006-06-29 |
JP4610968B2 JP4610968B2 (ja) | 2011-01-12 |
Family
ID=36111462
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2004246064A Expired - Fee Related JP4610968B2 (ja) | 2004-08-26 | 2004-08-26 | 位置計算装置、位置計算方法、位置計算プログラム、及び位置計算プログラムの記憶媒体 |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP4610968B2 (ja) |
Cited By (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2014108957A1 (ja) * | 2013-01-08 | 2014-07-17 | 日本電気株式会社 | 座標変換装置、座標変換プログラムが格納された非一時的なコンピュータ可読媒体、座標変換方法 |
WO2018083859A1 (ja) * | 2016-11-01 | 2018-05-11 | 日本電気株式会社 | 地理情報処理装置、地理情報処理方法及びコンピュータ可読媒体 |
CN110017833A (zh) * | 2019-04-19 | 2019-07-16 | 西安应用光学研究所 | 基于像元类地模型的全屏像点地理坐标定位方法 |
CN111028286A (zh) * | 2019-12-17 | 2020-04-17 | 中国人民解放军海军大连舰艇学院 | 精度可控的图斑椭球面积计算方法 |
CN112180412A (zh) * | 2020-09-23 | 2021-01-05 | 中国人民解放军空军工程大学 | 一种基于卫星导航定位***的相对定位定向补偿方法 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH095422A (ja) * | 1995-06-21 | 1997-01-10 | Yokogawa Denshi Kiki Kk | 電子海図システム |
JP2003065787A (ja) * | 2001-07-24 | 2003-03-05 | Navigation Technol Corp | 区分的線形補間を用いる距離計算方法 |
-
2004
- 2004-08-26 JP JP2004246064A patent/JP4610968B2/ja not_active Expired - Fee Related
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH095422A (ja) * | 1995-06-21 | 1997-01-10 | Yokogawa Denshi Kiki Kk | 電子海図システム |
JP2003065787A (ja) * | 2001-07-24 | 2003-03-05 | Navigation Technol Corp | 区分的線形補間を用いる距離計算方法 |
Cited By (8)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2014108957A1 (ja) * | 2013-01-08 | 2014-07-17 | 日本電気株式会社 | 座標変換装置、座標変換プログラムが格納された非一時的なコンピュータ可読媒体、座標変換方法 |
JP2014149807A (ja) * | 2013-01-08 | 2014-08-21 | Nec Corp | 座標変換装置、座標変換プログラム、座標変換方法 |
WO2018083859A1 (ja) * | 2016-11-01 | 2018-05-11 | 日本電気株式会社 | 地理情報処理装置、地理情報処理方法及びコンピュータ可読媒体 |
CN110017833A (zh) * | 2019-04-19 | 2019-07-16 | 西安应用光学研究所 | 基于像元类地模型的全屏像点地理坐标定位方法 |
CN110017833B (zh) * | 2019-04-19 | 2022-11-01 | 西安应用光学研究所 | 基于像元类地模型的全屏像点地理坐标定位方法 |
CN111028286A (zh) * | 2019-12-17 | 2020-04-17 | 中国人民解放军海军大连舰艇学院 | 精度可控的图斑椭球面积计算方法 |
CN112180412A (zh) * | 2020-09-23 | 2021-01-05 | 中国人民解放军空军工程大学 | 一种基于卫星导航定位***的相对定位定向补偿方法 |
CN112180412B (zh) * | 2020-09-23 | 2023-05-02 | 中国人民解放军空军工程大学 | 一种基于卫星导航定位***的相对定位定向补偿方法 |
Also Published As
Publication number | Publication date |
---|---|
JP4610968B2 (ja) | 2011-01-12 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Näher | LEDA, a platform for combinatorial and geometric computing | |
Baselga | Global optimization solution of robust estimation | |
Karney | Algorithms for geodesics | |
US7668700B2 (en) | Adaptive distance field constraint for designing a route for a transport element | |
Wongwathanarat et al. | An axis-free overset grid in spherical polar coordinates for simulating 3D self-gravitating flows | |
Karney | Geodesics on an ellipsoid of revolution | |
Mercan et al. | Solution of the weighted symmetric similarity transformations based on quaternions | |
Alton et al. | An ordered upwind method with precomputed stencil and monotone node acceptance for solving static convex Hamilton-Jacobi equations | |
Schreiner et al. | Direct (re) meshing for efficient surface processing | |
Lian et al. | Stress analysis without meshing: Isogeometric boundary-element method | |
Herholz et al. | Diffusion diagrams: Voronoi cells and centroids from diffusion | |
US20130066473A1 (en) | Generating Thermal Zones | |
Esenbuğa et al. | Comparison of principal geodetic distance calculation methods for automated province assignment in Turkey | |
Yoely et al. | Structural optimization with explicit geometric constraints using a B-spline representation | |
JP2006064892A (ja) | 位置計算装置、位置計算方法、位置計算プログラム、及び位置計算プログラムの記憶媒体 | |
JP5132731B2 (ja) | 位置計算装置、位置計算方法、位置計算プログラム、及び位置計算プログラムの記憶媒体 | |
JP2006064892A5 (ja) | ||
Li | Mesh modification procedures for general 3D non-manifold domains | |
Baldauf | Discontinuous Galerkin solver for the shallow-water equations in covariant form on the sphere and the ellipsoid | |
Fok et al. | A comparative analysis of the performance of iterative and non-iterative solutions to the Cartesian to geodetic coordinate transformation | |
KR20200092715A (ko) | 유전 알고리즘을 이용한 수술용 내비게이션의 표면 정합 방법 및 장치 | |
Habib et al. | Developing low-cost automated tool for integrating maps with GNSS satellite positioning data | |
Mitchell | Deterministic Linear Time for Maximal Poisson‐Disk Sampling using Chocks without Rejection or Approximation | |
US20230401350A1 (en) | Method for visualising a 3d infrastructure design model in a planar coordinate system on an ellipsoid | |
US10529444B1 (en) | System that rapidly generates a solvent-excluded surface |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20060509 |
|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20060509 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20090623 |
|
A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20090824 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20100615 |
|
A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20100816 |
|
TRDD | Decision of grant or rejection written | ||
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20101005 |
|
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 |
|
A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20101013 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20131022 Year of fee payment: 3 |
|
R150 | Certificate of patent or registration of utility model |
Ref document number: 4610968 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20131022 Year of fee payment: 3 |
|
FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20131022 Year of fee payment: 3 |
|
R157 | Certificate of patent or utility model (correction) |
Free format text: JAPANESE INTERMEDIATE CODE: R157 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
LAPS | Cancellation because of no payment of annual fees |