EllipsoidGeodesic-5b3623dc.js 15 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441442443444445446447448449450451452453454455456457458459460461462463464465466467468469470471472473474475476477478479480481482483484485486487488489490491492493494495496497498499500501502503504505506507508509510511512513514515516517518519520
  1. define(['exports', './Matrix3-41c58dde', './Check-6ede7e26', './defaultValue-fe22d8c0', './Math-0a2ac845'], (function (exports, Matrix3, Check, defaultValue, Math$1) { 'use strict';
  2. function setConstants(ellipsoidGeodesic) {
  3. const uSquared = ellipsoidGeodesic._uSquared;
  4. const a = ellipsoidGeodesic._ellipsoid.maximumRadius;
  5. const b = ellipsoidGeodesic._ellipsoid.minimumRadius;
  6. const f = (a - b) / a;
  7. const cosineHeading = Math.cos(ellipsoidGeodesic._startHeading);
  8. const sineHeading = Math.sin(ellipsoidGeodesic._startHeading);
  9. const tanU = (1 - f) * Math.tan(ellipsoidGeodesic._start.latitude);
  10. const cosineU = 1.0 / Math.sqrt(1.0 + tanU * tanU);
  11. const sineU = cosineU * tanU;
  12. const sigma = Math.atan2(tanU, cosineHeading);
  13. const sineAlpha = cosineU * sineHeading;
  14. const sineSquaredAlpha = sineAlpha * sineAlpha;
  15. const cosineSquaredAlpha = 1.0 - sineSquaredAlpha;
  16. const cosineAlpha = Math.sqrt(cosineSquaredAlpha);
  17. const u2Over4 = uSquared / 4.0;
  18. const u4Over16 = u2Over4 * u2Over4;
  19. const u6Over64 = u4Over16 * u2Over4;
  20. const u8Over256 = u4Over16 * u4Over16;
  21. const a0 =
  22. 1.0 +
  23. u2Over4 -
  24. (3.0 * u4Over16) / 4.0 +
  25. (5.0 * u6Over64) / 4.0 -
  26. (175.0 * u8Over256) / 64.0;
  27. const a1 = 1.0 - u2Over4 + (15.0 * u4Over16) / 8.0 - (35.0 * u6Over64) / 8.0;
  28. const a2 = 1.0 - 3.0 * u2Over4 + (35.0 * u4Over16) / 4.0;
  29. const a3 = 1.0 - 5.0 * u2Over4;
  30. const distanceRatio =
  31. a0 * sigma -
  32. (a1 * Math.sin(2.0 * sigma) * u2Over4) / 2.0 -
  33. (a2 * Math.sin(4.0 * sigma) * u4Over16) / 16.0 -
  34. (a3 * Math.sin(6.0 * sigma) * u6Over64) / 48.0 -
  35. (Math.sin(8.0 * sigma) * 5.0 * u8Over256) / 512;
  36. const constants = ellipsoidGeodesic._constants;
  37. constants.a = a;
  38. constants.b = b;
  39. constants.f = f;
  40. constants.cosineHeading = cosineHeading;
  41. constants.sineHeading = sineHeading;
  42. constants.tanU = tanU;
  43. constants.cosineU = cosineU;
  44. constants.sineU = sineU;
  45. constants.sigma = sigma;
  46. constants.sineAlpha = sineAlpha;
  47. constants.sineSquaredAlpha = sineSquaredAlpha;
  48. constants.cosineSquaredAlpha = cosineSquaredAlpha;
  49. constants.cosineAlpha = cosineAlpha;
  50. constants.u2Over4 = u2Over4;
  51. constants.u4Over16 = u4Over16;
  52. constants.u6Over64 = u6Over64;
  53. constants.u8Over256 = u8Over256;
  54. constants.a0 = a0;
  55. constants.a1 = a1;
  56. constants.a2 = a2;
  57. constants.a3 = a3;
  58. constants.distanceRatio = distanceRatio;
  59. }
  60. function computeC(f, cosineSquaredAlpha) {
  61. return (
  62. (f * cosineSquaredAlpha * (4.0 + f * (4.0 - 3.0 * cosineSquaredAlpha))) /
  63. 16.0
  64. );
  65. }
  66. function computeDeltaLambda(
  67. f,
  68. sineAlpha,
  69. cosineSquaredAlpha,
  70. sigma,
  71. sineSigma,
  72. cosineSigma,
  73. cosineTwiceSigmaMidpoint
  74. ) {
  75. const C = computeC(f, cosineSquaredAlpha);
  76. return (
  77. (1.0 - C) *
  78. f *
  79. sineAlpha *
  80. (sigma +
  81. C *
  82. sineSigma *
  83. (cosineTwiceSigmaMidpoint +
  84. C *
  85. cosineSigma *
  86. (2.0 * cosineTwiceSigmaMidpoint * cosineTwiceSigmaMidpoint - 1.0)))
  87. );
  88. }
  89. function vincentyInverseFormula(
  90. ellipsoidGeodesic,
  91. major,
  92. minor,
  93. firstLongitude,
  94. firstLatitude,
  95. secondLongitude,
  96. secondLatitude
  97. ) {
  98. const eff = (major - minor) / major;
  99. const l = secondLongitude - firstLongitude;
  100. const u1 = Math.atan((1 - eff) * Math.tan(firstLatitude));
  101. const u2 = Math.atan((1 - eff) * Math.tan(secondLatitude));
  102. const cosineU1 = Math.cos(u1);
  103. const sineU1 = Math.sin(u1);
  104. const cosineU2 = Math.cos(u2);
  105. const sineU2 = Math.sin(u2);
  106. const cc = cosineU1 * cosineU2;
  107. const cs = cosineU1 * sineU2;
  108. const ss = sineU1 * sineU2;
  109. const sc = sineU1 * cosineU2;
  110. let lambda = l;
  111. let lambdaDot = Math$1.CesiumMath.TWO_PI;
  112. let cosineLambda = Math.cos(lambda);
  113. let sineLambda = Math.sin(lambda);
  114. let sigma;
  115. let cosineSigma;
  116. let sineSigma;
  117. let cosineSquaredAlpha;
  118. let cosineTwiceSigmaMidpoint;
  119. do {
  120. cosineLambda = Math.cos(lambda);
  121. sineLambda = Math.sin(lambda);
  122. const temp = cs - sc * cosineLambda;
  123. sineSigma = Math.sqrt(
  124. cosineU2 * cosineU2 * sineLambda * sineLambda + temp * temp
  125. );
  126. cosineSigma = ss + cc * cosineLambda;
  127. sigma = Math.atan2(sineSigma, cosineSigma);
  128. let sineAlpha;
  129. if (sineSigma === 0.0) {
  130. sineAlpha = 0.0;
  131. cosineSquaredAlpha = 1.0;
  132. } else {
  133. sineAlpha = (cc * sineLambda) / sineSigma;
  134. cosineSquaredAlpha = 1.0 - sineAlpha * sineAlpha;
  135. }
  136. lambdaDot = lambda;
  137. cosineTwiceSigmaMidpoint = cosineSigma - (2.0 * ss) / cosineSquaredAlpha;
  138. if (!isFinite(cosineTwiceSigmaMidpoint)) {
  139. cosineTwiceSigmaMidpoint = 0.0;
  140. }
  141. lambda =
  142. l +
  143. computeDeltaLambda(
  144. eff,
  145. sineAlpha,
  146. cosineSquaredAlpha,
  147. sigma,
  148. sineSigma,
  149. cosineSigma,
  150. cosineTwiceSigmaMidpoint
  151. );
  152. } while (Math.abs(lambda - lambdaDot) > Math$1.CesiumMath.EPSILON12);
  153. const uSquared =
  154. (cosineSquaredAlpha * (major * major - minor * minor)) / (minor * minor);
  155. const A =
  156. 1.0 +
  157. (uSquared *
  158. (4096.0 + uSquared * (uSquared * (320.0 - 175.0 * uSquared) - 768.0))) /
  159. 16384.0;
  160. const B =
  161. (uSquared *
  162. (256.0 + uSquared * (uSquared * (74.0 - 47.0 * uSquared) - 128.0))) /
  163. 1024.0;
  164. const cosineSquaredTwiceSigmaMidpoint =
  165. cosineTwiceSigmaMidpoint * cosineTwiceSigmaMidpoint;
  166. const deltaSigma =
  167. B *
  168. sineSigma *
  169. (cosineTwiceSigmaMidpoint +
  170. (B *
  171. (cosineSigma * (2.0 * cosineSquaredTwiceSigmaMidpoint - 1.0) -
  172. (B *
  173. cosineTwiceSigmaMidpoint *
  174. (4.0 * sineSigma * sineSigma - 3.0) *
  175. (4.0 * cosineSquaredTwiceSigmaMidpoint - 3.0)) /
  176. 6.0)) /
  177. 4.0);
  178. const distance = minor * A * (sigma - deltaSigma);
  179. const startHeading = Math.atan2(
  180. cosineU2 * sineLambda,
  181. cs - sc * cosineLambda
  182. );
  183. const endHeading = Math.atan2(cosineU1 * sineLambda, cs * cosineLambda - sc);
  184. ellipsoidGeodesic._distance = distance;
  185. ellipsoidGeodesic._startHeading = startHeading;
  186. ellipsoidGeodesic._endHeading = endHeading;
  187. ellipsoidGeodesic._uSquared = uSquared;
  188. }
  189. const scratchCart1 = new Matrix3.Cartesian3();
  190. const scratchCart2 = new Matrix3.Cartesian3();
  191. function computeProperties(ellipsoidGeodesic, start, end, ellipsoid) {
  192. const firstCartesian = Matrix3.Cartesian3.normalize(
  193. ellipsoid.cartographicToCartesian(start, scratchCart2),
  194. scratchCart1
  195. );
  196. const lastCartesian = Matrix3.Cartesian3.normalize(
  197. ellipsoid.cartographicToCartesian(end, scratchCart2),
  198. scratchCart2
  199. );
  200. //>>includeStart('debug', pragmas.debug);
  201. Check.Check.typeOf.number.greaterThanOrEquals(
  202. "value",
  203. Math.abs(
  204. Math.abs(Matrix3.Cartesian3.angleBetween(firstCartesian, lastCartesian)) - Math.PI
  205. ),
  206. 0.0125
  207. );
  208. //>>includeEnd('debug');
  209. vincentyInverseFormula(
  210. ellipsoidGeodesic,
  211. ellipsoid.maximumRadius,
  212. ellipsoid.minimumRadius,
  213. start.longitude,
  214. start.latitude,
  215. end.longitude,
  216. end.latitude
  217. );
  218. ellipsoidGeodesic._start = Matrix3.Cartographic.clone(
  219. start,
  220. ellipsoidGeodesic._start
  221. );
  222. ellipsoidGeodesic._end = Matrix3.Cartographic.clone(end, ellipsoidGeodesic._end);
  223. ellipsoidGeodesic._start.height = 0;
  224. ellipsoidGeodesic._end.height = 0;
  225. setConstants(ellipsoidGeodesic);
  226. }
  227. /**
  228. * Initializes a geodesic on the ellipsoid connecting the two provided planetodetic points.
  229. *
  230. * @alias EllipsoidGeodesic
  231. * @constructor
  232. *
  233. * @param {Cartographic} [start] The initial planetodetic point on the path.
  234. * @param {Cartographic} [end] The final planetodetic point on the path.
  235. * @param {Ellipsoid} [ellipsoid=Ellipsoid.WGS84] The ellipsoid on which the geodesic lies.
  236. */
  237. function EllipsoidGeodesic(start, end, ellipsoid) {
  238. const e = defaultValue.defaultValue(ellipsoid, Matrix3.Ellipsoid.WGS84);
  239. this._ellipsoid = e;
  240. this._start = new Matrix3.Cartographic();
  241. this._end = new Matrix3.Cartographic();
  242. this._constants = {};
  243. this._startHeading = undefined;
  244. this._endHeading = undefined;
  245. this._distance = undefined;
  246. this._uSquared = undefined;
  247. if (defaultValue.defined(start) && defaultValue.defined(end)) {
  248. computeProperties(this, start, end, e);
  249. }
  250. }
  251. Object.defineProperties(EllipsoidGeodesic.prototype, {
  252. /**
  253. * Gets the ellipsoid.
  254. * @memberof EllipsoidGeodesic.prototype
  255. * @type {Ellipsoid}
  256. * @readonly
  257. */
  258. ellipsoid: {
  259. get: function () {
  260. return this._ellipsoid;
  261. },
  262. },
  263. /**
  264. * Gets the surface distance between the start and end point
  265. * @memberof EllipsoidGeodesic.prototype
  266. * @type {number}
  267. * @readonly
  268. */
  269. surfaceDistance: {
  270. get: function () {
  271. //>>includeStart('debug', pragmas.debug);
  272. Check.Check.defined("distance", this._distance);
  273. //>>includeEnd('debug');
  274. return this._distance;
  275. },
  276. },
  277. /**
  278. * Gets the initial planetodetic point on the path.
  279. * @memberof EllipsoidGeodesic.prototype
  280. * @type {Cartographic}
  281. * @readonly
  282. */
  283. start: {
  284. get: function () {
  285. return this._start;
  286. },
  287. },
  288. /**
  289. * Gets the final planetodetic point on the path.
  290. * @memberof EllipsoidGeodesic.prototype
  291. * @type {Cartographic}
  292. * @readonly
  293. */
  294. end: {
  295. get: function () {
  296. return this._end;
  297. },
  298. },
  299. /**
  300. * Gets the heading at the initial point.
  301. * @memberof EllipsoidGeodesic.prototype
  302. * @type {number}
  303. * @readonly
  304. */
  305. startHeading: {
  306. get: function () {
  307. //>>includeStart('debug', pragmas.debug);
  308. Check.Check.defined("distance", this._distance);
  309. //>>includeEnd('debug');
  310. return this._startHeading;
  311. },
  312. },
  313. /**
  314. * Gets the heading at the final point.
  315. * @memberof EllipsoidGeodesic.prototype
  316. * @type {number}
  317. * @readonly
  318. */
  319. endHeading: {
  320. get: function () {
  321. //>>includeStart('debug', pragmas.debug);
  322. Check.Check.defined("distance", this._distance);
  323. //>>includeEnd('debug');
  324. return this._endHeading;
  325. },
  326. },
  327. });
  328. /**
  329. * Sets the start and end points of the geodesic
  330. *
  331. * @param {Cartographic} start The initial planetodetic point on the path.
  332. * @param {Cartographic} end The final planetodetic point on the path.
  333. */
  334. EllipsoidGeodesic.prototype.setEndPoints = function (start, end) {
  335. //>>includeStart('debug', pragmas.debug);
  336. Check.Check.defined("start", start);
  337. Check.Check.defined("end", end);
  338. //>>includeEnd('debug');
  339. computeProperties(this, start, end, this._ellipsoid);
  340. };
  341. /**
  342. * Provides the location of a point at the indicated portion along the geodesic.
  343. *
  344. * @param {number} fraction The portion of the distance between the initial and final points.
  345. * @param {Cartographic} [result] The object in which to store the result.
  346. * @returns {Cartographic} The location of the point along the geodesic.
  347. */
  348. EllipsoidGeodesic.prototype.interpolateUsingFraction = function (
  349. fraction,
  350. result
  351. ) {
  352. return this.interpolateUsingSurfaceDistance(
  353. this._distance * fraction,
  354. result
  355. );
  356. };
  357. /**
  358. * Provides the location of a point at the indicated distance along the geodesic.
  359. *
  360. * @param {number} distance The distance from the inital point to the point of interest along the geodesic
  361. * @param {Cartographic} [result] The object in which to store the result.
  362. * @returns {Cartographic} The location of the point along the geodesic.
  363. *
  364. * @exception {DeveloperError} start and end must be set before calling function interpolateUsingSurfaceDistance
  365. */
  366. EllipsoidGeodesic.prototype.interpolateUsingSurfaceDistance = function (
  367. distance,
  368. result
  369. ) {
  370. //>>includeStart('debug', pragmas.debug);
  371. Check.Check.defined("distance", this._distance);
  372. //>>includeEnd('debug');
  373. const constants = this._constants;
  374. const s = constants.distanceRatio + distance / constants.b;
  375. const cosine2S = Math.cos(2.0 * s);
  376. const cosine4S = Math.cos(4.0 * s);
  377. const cosine6S = Math.cos(6.0 * s);
  378. const sine2S = Math.sin(2.0 * s);
  379. const sine4S = Math.sin(4.0 * s);
  380. const sine6S = Math.sin(6.0 * s);
  381. const sine8S = Math.sin(8.0 * s);
  382. const s2 = s * s;
  383. const s3 = s * s2;
  384. const u8Over256 = constants.u8Over256;
  385. const u2Over4 = constants.u2Over4;
  386. const u6Over64 = constants.u6Over64;
  387. const u4Over16 = constants.u4Over16;
  388. let sigma =
  389. (2.0 * s3 * u8Over256 * cosine2S) / 3.0 +
  390. s *
  391. (1.0 -
  392. u2Over4 +
  393. (7.0 * u4Over16) / 4.0 -
  394. (15.0 * u6Over64) / 4.0 +
  395. (579.0 * u8Over256) / 64.0 -
  396. (u4Over16 - (15.0 * u6Over64) / 4.0 + (187.0 * u8Over256) / 16.0) *
  397. cosine2S -
  398. ((5.0 * u6Over64) / 4.0 - (115.0 * u8Over256) / 16.0) * cosine4S -
  399. (29.0 * u8Over256 * cosine6S) / 16.0) +
  400. (u2Over4 / 2.0 -
  401. u4Over16 +
  402. (71.0 * u6Over64) / 32.0 -
  403. (85.0 * u8Over256) / 16.0) *
  404. sine2S +
  405. ((5.0 * u4Over16) / 16.0 -
  406. (5.0 * u6Over64) / 4.0 +
  407. (383.0 * u8Over256) / 96.0) *
  408. sine4S -
  409. s2 *
  410. ((u6Over64 - (11.0 * u8Over256) / 2.0) * sine2S +
  411. (5.0 * u8Over256 * sine4S) / 2.0) +
  412. ((29.0 * u6Over64) / 96.0 - (29.0 * u8Over256) / 16.0) * sine6S +
  413. (539.0 * u8Over256 * sine8S) / 1536.0;
  414. const theta = Math.asin(Math.sin(sigma) * constants.cosineAlpha);
  415. const latitude = Math.atan((constants.a / constants.b) * Math.tan(theta));
  416. // Redefine in terms of relative argument of latitude.
  417. sigma = sigma - constants.sigma;
  418. const cosineTwiceSigmaMidpoint = Math.cos(2.0 * constants.sigma + sigma);
  419. const sineSigma = Math.sin(sigma);
  420. const cosineSigma = Math.cos(sigma);
  421. const cc = constants.cosineU * cosineSigma;
  422. const ss = constants.sineU * sineSigma;
  423. const lambda = Math.atan2(
  424. sineSigma * constants.sineHeading,
  425. cc - ss * constants.cosineHeading
  426. );
  427. const l =
  428. lambda -
  429. computeDeltaLambda(
  430. constants.f,
  431. constants.sineAlpha,
  432. constants.cosineSquaredAlpha,
  433. sigma,
  434. sineSigma,
  435. cosineSigma,
  436. cosineTwiceSigmaMidpoint
  437. );
  438. if (defaultValue.defined(result)) {
  439. result.longitude = this._start.longitude + l;
  440. result.latitude = latitude;
  441. result.height = 0.0;
  442. return result;
  443. }
  444. return new Matrix3.Cartographic(this._start.longitude + l, latitude, 0.0);
  445. };
  446. exports.EllipsoidGeodesic = EllipsoidGeodesic;
  447. }));