• R/O
  • HTTP
  • SSH
  • HTTPS

yubeshi: Commit

source codes and tests


Commit MetaInfo

Revision7011ee5616402d4bb7e10bdc7df066ebc7056e21 (tree)
Time2011-02-28 22:53:46
Authorkimikage <kimikage_ceo@hotm...>
Commiterkimikage

Log Message

球体ベースのENU座標とは別に,日本の平面直角座標を追加

Change Summary

Incremental Difference

--- /dev/null
+++ b/Yubeshi/Japan/Jpc.cs
@@ -0,0 +1,229 @@
1+/*
2+ * Yubeshi GPS Parser
3+ *
4+ * This software is distributed under a zlib-style license.
5+ * See license.txt for more information.
6+ */
7+
8+using System;
9+using System.Collections.Generic;
10+using System.Text;
11+using Llh = Yubeshi.GeodeticCoordinate;
12+
13+namespace Yubeshi.Japan
14+{
15+
16+ public class Jpc
17+ {
18+ // cf.
19+ // http://vldb.gsi.go.jp/sokuchi/surveycalc/algorithm/bl2xy/bl2xy.htm
20+ #region type definition
21+ private static class MeridianArcConstants
22+ {
23+ public static readonly double[] K = new double[]{
24+ 1.0050525018051135299706,
25+ -5.063108622326803921017e-3,
26+ 1.062759032760013166595e-5,
27+ -2.082040715851126522938e-8,
28+ 3.933230470035652928709e-11,
29+ -7.265236169124544217757e-14,
30+ 1.321653368923020906791e-16,
31+ -2.376757353699231655954e-19,
32+ 4.109419740268720338702e-22
33+ };
34+ /*
35+ private static long[] GenerateSquaredBinomial(int length)
36+ {
37+ // cf. A002894 Binomial(2n,n)^2
38+ long k = 1;
39+ long[] b = new long[length];
40+ for(int i = 1; i < length; ++i)
41+ {
42+ b[i - 1] = k * k;
43+ long p = (i + i) * (i + i - 1);
44+ long q = i * i;
45+ k *= p;
46+ k /= q;
47+ }
48+ b[length - 1] = k * k;
49+ return b;
50+ }
51+ */
52+
53+ }
54+ public enum System
55+ {
56+ Unknown = 0,
57+ I = 1,
58+ II = 2,
59+ III = 3,
60+ IV = 4,
61+ V = 5,
62+ VI = 6,
63+ VII = 7,
64+ VIII = 8,
65+ IX = 9,
66+ X = 10,
67+ XI = 11,
68+ XII = 12,
69+ XIII = 13,
70+ XIV = 14,
71+ XV = 15,
72+ XVI = 16,
73+ XVII = 17,
74+ XVIII = 18,
75+ XIX = 19,
76+ }
77+ #endregion
78+
79+ #region fields
80+
81+ public static readonly GeodeticCoordinate[] Origins;
82+ private const double e = 0.00669438002290069;
83+ private const double e2 = 6.7394967754788593e-3;
84+ private const double m0 = 0.9999;
85+ #endregion
86+
87+ #region constructors
88+ public Jpc(GeodeticCoordinate llh, System number)
89+ : this(llh.Latitude, llh.Longitude, number)
90+ {
91+ }
92+
93+ public Jpc(Degree latitude, Degree longitude, System number)
94+ {
95+ GeodeticCoordinate origin = Origins[(int)number];
96+ double n = GetRadiusOfPrimeVertical(latitude);
97+ double cos = Math.Cos(latitude.Radian);
98+ double cos2 = cos * cos;
99+ double c = cos2;
100+
101+ double delta = longitude.Radian - origin.Longitude.Radian;
102+ double delta2 = delta * delta;
103+ double d = delta2;
104+
105+ double t = Math.Tan(latitude.Radian);
106+ double t2 = t * t;
107+ double t4 = t2 * t2;
108+
109+ double eta2 = e2 * cos2;
110+ double eta4 = eta2 * eta2;
111+
112+ double deltaS = GetLengthOfMeridianArc(latitude) -
113+ GetLengthOfMeridianArc(origin.Latitude);
114+
115+ X = cos * t * d / 2.0;
116+
117+ d *= delta2;
118+ c *= cos2;
119+ X += c * t * d / 24.0 *
120+ (5.0 - t2 + 9.0 * eta2 + 4.0 * eta4);
121+
122+ d *= delta2;
123+ c *= cos2;
124+ X -= c * t * d / 720.0 *
125+ (-61.0 + 58.0 * t2 -t4 -270.0 * eta2 + 330.0 * t2 * eta2);
126+
127+ d *= delta2;
128+ c *= cos2;
129+ X -= c * t * d / 40320.0 *
130+ (-1385.0 + 3111.0 * t2 - 543.0 * t4 + t4 * t2);
131+
132+ X *= n;
133+ X = (X + deltaS) * m0;
134+
135+
136+ c = cos;
137+ d = delta;
138+
139+ Y = c * d;
140+
141+ d *= delta2;
142+ c *= cos2;
143+ Y -= c * d / 6.0 * (-1.0 + t2 - eta2);
144+
145+ d *= delta2;
146+ c *= cos2;
147+ Y -= c * d / 120.0 *
148+ (-5.0 + 18.0 * t2 - t4 - 14.0 * eta2 + 58.0 * t2 * eta2);
149+
150+ d *= delta2;
151+ c *= cos2;
152+ Y -= c * d / 5040.0 *
153+ (-61.0 + 479.0 * t2 - 179.0 * t4 + t4 * t2);
154+ Y *= n * m0;
155+ }
156+
157+ public Jpc(double x, double y)
158+ {
159+ X = x;
160+ Y = y;
161+ }
162+
163+ static Jpc()
164+ {
165+ Origins = new GeodeticCoordinate[20]{
166+ null,
167+ new Llh(new Degree(33, 0, 0), new Degree(129, 30, 0)),
168+ new Llh(new Degree(33, 0, 0), new Degree(131, 0, 0)),
169+ new Llh(new Degree(36, 0, 0), new Degree(132, 10, 0)),
170+ new Llh(new Degree(33, 0, 0), new Degree(133, 30, 0)),
171+ new Llh(new Degree(36, 0, 0), new Degree(134, 20, 0)),
172+ new Llh(new Degree(36, 0, 0), new Degree(136, 0, 0)),
173+ new Llh(new Degree(36, 0, 0), new Degree(137, 10, 0)),
174+ new Llh(new Degree(36, 0, 0), new Degree(138, 30, 0)),
175+ new Llh(new Degree(36, 0, 0), new Degree(139, 50, 0)),
176+ new Llh(new Degree(40, 0, 0), new Degree(140, 50, 0)),
177+ new Llh(new Degree(44, 0, 0), new Degree(140, 15, 0)),
178+ new Llh(new Degree(44, 0, 0), new Degree(142, 15, 0)),
179+ new Llh(new Degree(44, 0, 0), new Degree(144, 15, 0)),
180+ new Llh(new Degree(26, 0, 0), new Degree(142, 0, 0)),
181+ new Llh(new Degree(26, 0, 0), new Degree(127, 30, 0)),
182+ new Llh(new Degree(26, 0, 0), new Degree(124, 0, 0)),
183+ new Llh(new Degree(26, 0, 0), new Degree(131, 0, 0)),
184+ new Llh(new Degree(20, 0, 0), new Degree(136, 0, 0)),
185+ new Llh(new Degree(26, 0, 0), new Degree(154, 0, 0)),
186+ };
187+ }
188+
189+ #endregion
190+
191+ #region properties
192+ public double X
193+ {
194+ get;
195+ private set;
196+ }
197+
198+ public double Y
199+ {
200+ get;
201+ private set;
202+ }
203+
204+
205+ #endregion
206+
207+ #region public methods
208+
209+ public static double GetLengthOfMeridianArc(Degree latitude)
210+ {
211+ double s = 0;
212+ double phi = latitude.Radian;
213+ for (int i = MeridianArcConstants.K.Length -1; i > 0; --i)
214+ {
215+ int j = i + i;
216+ s += Math.Sin(j * phi) * MeridianArcConstants.K[i] / j;
217+ }
218+ s += phi * MeridianArcConstants.K[0];
219+ return s * Constants.SemiMajorAxisA * (1 - e);
220+ }
221+
222+ public static double GetRadiusOfPrimeVertical(Degree latitude)
223+ {
224+ double sin = Math.Sin(latitude.Radian);
225+ return Constants.SemiMajorAxisA / Math.Sqrt(1 - e * sin * sin);
226+ }
227+ #endregion
228+ }
229+}
--- a/Yubeshi/Yubeshi.csproj
+++ b/Yubeshi/Yubeshi.csproj
@@ -48,6 +48,7 @@
4848 <Compile Include="Gps\Word.cs" />
4949 <Compile Include="Gps\Subframe.cs" />
5050 <Compile Include="Gps\HandoverWord.cs" />
51+ <Compile Include="Japan\Jpc.cs" />
5152 <Compile Include="Nmea\Parser.cs" />
5253 <Compile Include="OctetString.cs" />
5354 <Compile Include="Parser.cs" />
--- a/YubeshiTest/EnuCoordinateTest.cs
+++ b/YubeshiTest/EnuCoordinateTest.cs
@@ -10,18 +10,17 @@ using System.Collections.Generic;
1010 using System.Text;
1111 using NUnit.Framework;
1212 using Yubeshi;
13+using C = YubeshiTest.SampleCoordinates;
1314
1415 namespace YubeshiTest
1516 {
1617 class EnuCoordinateTest
1718 {
18-
19+ private static double metreError = 1e-9; // 1nm
1920 [Test]
2021 public void ZeroDistance()
2122 {
22- GeodeticCoordinate c =
23- new GeodeticCoordinate(35.710058, 139.810719, 634.0);
24- EnuCoordinate enu = new EnuCoordinate(c, c);
23+ EnuCoordinate enu = new EnuCoordinate(C.SkyTreeTop, C.SkyTreeTop);
2524 Assert.AreEqual(0.0, enu.E);
2625 Assert.AreEqual(0.0, enu.N);
2726 Assert.AreEqual(0.0, enu.U);
@@ -30,24 +29,19 @@ namespace YubeshiTest
3029 [Test]
3130 public void UpOnly()
3231 {
33- GeodeticCoordinate c1 =
34- new GeodeticCoordinate(35.710058, 139.810719, 634.0);
35- GeodeticCoordinate c2 =
36- new GeodeticCoordinate(35.710058, 139.810719, 0.0);
37- EnuCoordinate enu = new EnuCoordinate(c1, c2);
32+ EnuCoordinate enu = new EnuCoordinate(
33+ C.SkyTreeTop, C.SkyTreeBottom);
3834
39- Assert.AreEqual(0.0, enu.E, 1e-9);
40- Assert.AreEqual(0.0, enu.N, 1e-9);
41- Assert.AreEqual(634.0, enu.U, 1e-9);
35+ Assert.AreEqual(0.0, enu.E, metreError);
36+ Assert.AreEqual(0.0, enu.N, metreError);
37+ Assert.AreEqual(634.0, enu.U, metreError);
4238 }
4339
4440 [Test]
4541 public void TwoPoints()
4642 {
47- GeodeticCoordinate c1 =
48- new GeodeticCoordinate(35.710058, 139.810719, 634.0);
49- GeodeticCoordinate c2 =
50- new GeodeticCoordinate(35.658611, 139.745556, 351.0);
43+ GeodeticCoordinate c1 = C.SkyTreeTop;
44+ GeodeticCoordinate c2 = C.TokyoTowerTop;
5145 GeodeticCoordinate c0 = new GeodeticCoordinate(
5246 0.5 * (c1.Latitude + c2.Latitude),
5347 0.5 * (c1.Longitude + c2.Longitude),
@@ -70,5 +64,18 @@ namespace YubeshiTest
7064 Assert.AreEqual(45.946, direction, 1.0 / 60.0);
7165 Assert.AreEqual(283.0, du, 1e-3);
7266 }
67+
68+ [Test]
69+ public void ToEcefCoordinate()
70+ {
71+ EnuCoordinate enu = new EnuCoordinate(
72+ C.SkyTreeTop, C.TokyoTowerTop);
73+ EcefCoordinate ecef1 = enu.ToEcefCoordinate();
74+ EcefCoordinate ecef2 = C.SkyTreeTop.ToEcefCoordinate();
75+
76+ Assert.AreEqual(ecef2.X, ecef1.X, metreError);
77+ Assert.AreEqual(ecef2.Y, ecef1.Y, metreError);
78+ Assert.AreEqual(ecef2.Z, ecef1.Z, metreError);
79+ }
7380 }
7481 }
--- a/YubeshiTest/GpsTimeTest.cs
+++ b/YubeshiTest/GpsTimeTest.cs
@@ -32,7 +32,7 @@ namespace YubeshiTest
3232 }
3333
3434 [Test]
35- public void FromUTC()
35+ public void FromUtc()
3636 {
3737 GpsTime gps = new GpsTime(utc);
3838 Assert.AreEqual(false, gps.Rounded);
--- /dev/null
+++ b/YubeshiTest/JapanTest/JpcTest.cs
@@ -0,0 +1,37 @@
1+/*
2+ * Yubeshi GPS Parser
3+ *
4+ * This software is distributed under a zlib-style license.
5+ * See license.txt for more information.
6+ */
7+
8+using System;
9+using System.Collections.Generic;
10+using System.Text;
11+using NUnit.Framework;
12+using Yubeshi.Japan;
13+using C = YubeshiTest.SampleCoordinates;
14+
15+namespace YubeshiTest.JapanTest
16+{
17+ class JpcTest
18+ {
19+
20+ [Test]
21+ public void LlhToJpc9()
22+ {
23+ Jpc jpc = new Jpc(C.SkyTreeTop, Jpc.System.IX);
24+ Assert.AreEqual(-32167.406, jpc.X, 0.1);
25+ Assert.AreEqual(-2046.184, jpc.Y, 0.1);
26+ }
27+
28+ [Test]
29+ public void GetLengthOfMeridianArc()
30+ {
31+ Assert.AreEqual(0.0, Jpc.GetLengthOfMeridianArc(0.0));
32+ Assert.AreEqual(10001965.7, Jpc.GetLengthOfMeridianArc(90.0), 0.1);
33+ Assert.AreEqual(4984944.4, Jpc.GetLengthOfMeridianArc(45.0), 0.1);
34+ }
35+
36+ }
37+}
--- /dev/null
+++ b/YubeshiTest/SampleCoordinates.cs
@@ -0,0 +1,26 @@
1+/*
2+ * Yubeshi GPS Parser
3+ *
4+ * This software is distributed under a zlib-style license.
5+ * See license.txt for more information.
6+ */
7+
8+using System;
9+using System.Collections.Generic;
10+using System.Text;
11+using Yubeshi;
12+
13+namespace YubeshiTest
14+{
15+ public class SampleCoordinates
16+ {
17+ public static readonly GeodeticCoordinate SkyTreeTop =
18+ new GeodeticCoordinate(35.710058, 139.810719, 634.0);
19+
20+ public static readonly GeodeticCoordinate SkyTreeBottom =
21+ new GeodeticCoordinate(35.710058, 139.810719, 0.0);
22+
23+ public static readonly GeodeticCoordinate TokyoTowerTop =
24+ new GeodeticCoordinate(35.658611, 139.745556, 351.0);
25+ }
26+}
--- a/YubeshiTest/YubeshiTest.csproj
+++ b/YubeshiTest/YubeshiTest.csproj
@@ -48,11 +48,13 @@
4848 <Compile Include="GeodeticCoordinateTest.cs" />
4949 <Compile Include="GpsTest\WordTest.cs" />
5050 <Compile Include="GpsTimeTest.cs" />
51+ <Compile Include="JapanTest\JpcTest.cs" />
5152 <Compile Include="NmeaTest\PacketTest.cs" />
5253 <Compile Include="NmeaTest\ParserTest.cs" />
5354 <Compile Include="NmeaTest\SamplePackets.cs" />
5455 <Compile Include="ParserTest.cs" />
5556 <Compile Include="Properties\AssemblyInfo.cs" />
57+ <Compile Include="SampleCoordinates.cs" />
5658 <Compile Include="UbxTest\MyPacket.cs" />
5759 <Compile Include="UbxTest\NavPacketTest.cs" />
5860 <Compile Include="UbxTest\ParserTest.cs" />
Show on old repository browser