EllipsoidGeodesic-84507801.js 17 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413
  1. /**
  2. * Cesium - https://github.com/AnalyticalGraphicsInc/cesium
  3. *
  4. * Copyright 2011-2017 Cesium Contributors
  5. *
  6. * Licensed under the Apache License, Version 2.0 (the "License");
  7. * you may not use this file except in compliance with the License.
  8. * You may obtain a copy of the License at
  9. *
  10. * http://www.apache.org/licenses/LICENSE-2.0
  11. *
  12. * Unless required by applicable law or agreed to in writing, software
  13. * distributed under the License is distributed on an "AS IS" BASIS,
  14. * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  15. * See the License for the specific language governing permissions and
  16. * limitations under the License.
  17. *
  18. * Columbus View (Pat. Pend.)
  19. *
  20. * Portions licensed separately.
  21. * See https://github.com/AnalyticalGraphicsInc/cesium/blob/master/LICENSE.md for full licensing details.
  22. */
  23. define(['exports', './when-8d13db60', './Check-70bec281', './Math-61ede240', './Cartographic-fe4be337', './Cartesian2-85064f09'], function (exports, when, Check, _Math, Cartographic, Cartesian2) { 'use strict';
  24. function setConstants(ellipsoidGeodesic) {
  25. var uSquared = ellipsoidGeodesic._uSquared;
  26. var a = ellipsoidGeodesic._ellipsoid.maximumRadius;
  27. var b = ellipsoidGeodesic._ellipsoid.minimumRadius;
  28. var f = (a - b) / a;
  29. var cosineHeading = Math.cos(ellipsoidGeodesic._startHeading);
  30. var sineHeading = Math.sin(ellipsoidGeodesic._startHeading);
  31. var tanU = (1 - f) * Math.tan(ellipsoidGeodesic._start.latitude);
  32. var cosineU = 1.0 / Math.sqrt(1.0 + tanU * tanU);
  33. var sineU = cosineU * tanU;
  34. var sigma = Math.atan2(tanU, cosineHeading);
  35. var sineAlpha = cosineU * sineHeading;
  36. var sineSquaredAlpha = sineAlpha * sineAlpha;
  37. var cosineSquaredAlpha = 1.0 - sineSquaredAlpha;
  38. var cosineAlpha = Math.sqrt(cosineSquaredAlpha);
  39. var u2Over4 = uSquared / 4.0;
  40. var u4Over16 = u2Over4 * u2Over4;
  41. var u6Over64 = u4Over16 * u2Over4;
  42. var u8Over256 = u4Over16 * u4Over16;
  43. var a0 = (1.0 + u2Over4 - 3.0 * u4Over16 / 4.0 + 5.0 * u6Over64 / 4.0 - 175.0 * u8Over256 / 64.0);
  44. var a1 = (1.0 - u2Over4 + 15.0 * u4Over16 / 8.0 - 35.0 * u6Over64 / 8.0);
  45. var a2 = (1.0 - 3.0 * u2Over4 + 35.0 * u4Over16 / 4.0);
  46. var a3 = (1.0 - 5.0 * u2Over4);
  47. var distanceRatio = a0 * sigma - a1 * Math.sin(2.0 * sigma) * u2Over4 / 2.0 - a2 * Math.sin(4.0 * sigma) * u4Over16 / 16.0 -
  48. a3 * Math.sin(6.0 * sigma) * u6Over64 / 48.0 - Math.sin(8.0 * sigma) * 5.0 * u8Over256 / 512;
  49. var constants = ellipsoidGeodesic._constants;
  50. constants.a = a;
  51. constants.b = b;
  52. constants.f = f;
  53. constants.cosineHeading = cosineHeading;
  54. constants.sineHeading = sineHeading;
  55. constants.tanU = tanU;
  56. constants.cosineU = cosineU;
  57. constants.sineU = sineU;
  58. constants.sigma = sigma;
  59. constants.sineAlpha = sineAlpha;
  60. constants.sineSquaredAlpha = sineSquaredAlpha;
  61. constants.cosineSquaredAlpha = cosineSquaredAlpha;
  62. constants.cosineAlpha = cosineAlpha;
  63. constants.u2Over4 = u2Over4;
  64. constants.u4Over16 = u4Over16;
  65. constants.u6Over64 = u6Over64;
  66. constants.u8Over256 = u8Over256;
  67. constants.a0 = a0;
  68. constants.a1 = a1;
  69. constants.a2 = a2;
  70. constants.a3 = a3;
  71. constants.distanceRatio = distanceRatio;
  72. }
  73. function computeC(f, cosineSquaredAlpha) {
  74. return f * cosineSquaredAlpha * (4.0 + f * (4.0 - 3.0 * cosineSquaredAlpha)) / 16.0;
  75. }
  76. function computeDeltaLambda(f, sineAlpha, cosineSquaredAlpha, sigma, sineSigma, cosineSigma, cosineTwiceSigmaMidpoint) {
  77. var C = computeC(f, cosineSquaredAlpha);
  78. return (1.0 - C) * f * sineAlpha * (sigma + C * sineSigma * (cosineTwiceSigmaMidpoint +
  79. C * cosineSigma * (2.0 * cosineTwiceSigmaMidpoint * cosineTwiceSigmaMidpoint - 1.0)));
  80. }
  81. function vincentyInverseFormula(ellipsoidGeodesic, major, minor, firstLongitude, firstLatitude, secondLongitude, secondLatitude) {
  82. var eff = (major - minor) / major;
  83. var l = secondLongitude - firstLongitude;
  84. var u1 = Math.atan((1 - eff) * Math.tan(firstLatitude));
  85. var u2 = Math.atan((1 - eff) * Math.tan(secondLatitude));
  86. var cosineU1 = Math.cos(u1);
  87. var sineU1 = Math.sin(u1);
  88. var cosineU2 = Math.cos(u2);
  89. var sineU2 = Math.sin(u2);
  90. var cc = cosineU1 * cosineU2;
  91. var cs = cosineU1 * sineU2;
  92. var ss = sineU1 * sineU2;
  93. var sc = sineU1 * cosineU2;
  94. var lambda = l;
  95. var lambdaDot = _Math.CesiumMath.TWO_PI;
  96. var cosineLambda = Math.cos(lambda);
  97. var sineLambda = Math.sin(lambda);
  98. var sigma;
  99. var cosineSigma;
  100. var sineSigma;
  101. var cosineSquaredAlpha;
  102. var cosineTwiceSigmaMidpoint;
  103. do {
  104. cosineLambda = Math.cos(lambda);
  105. sineLambda = Math.sin(lambda);
  106. var temp = cs - sc * cosineLambda;
  107. sineSigma = Math.sqrt(cosineU2 * cosineU2 * sineLambda * sineLambda + temp * temp);
  108. cosineSigma = ss + cc * cosineLambda;
  109. sigma = Math.atan2(sineSigma, cosineSigma);
  110. var sineAlpha;
  111. if (sineSigma === 0.0) {
  112. sineAlpha = 0.0;
  113. cosineSquaredAlpha = 1.0;
  114. } else {
  115. sineAlpha = cc * sineLambda / sineSigma;
  116. cosineSquaredAlpha = 1.0 - sineAlpha * sineAlpha;
  117. }
  118. lambdaDot = lambda;
  119. cosineTwiceSigmaMidpoint = cosineSigma - 2.0 * ss / cosineSquaredAlpha;
  120. if (isNaN(cosineTwiceSigmaMidpoint)) {
  121. cosineTwiceSigmaMidpoint = 0.0;
  122. }
  123. lambda = l + computeDeltaLambda(eff, sineAlpha, cosineSquaredAlpha,
  124. sigma, sineSigma, cosineSigma, cosineTwiceSigmaMidpoint);
  125. } while (Math.abs(lambda - lambdaDot) > _Math.CesiumMath.EPSILON12);
  126. var uSquared = cosineSquaredAlpha * (major * major - minor * minor) / (minor * minor);
  127. var A = 1.0 + uSquared * (4096.0 + uSquared * (uSquared * (320.0 - 175.0 * uSquared) - 768.0)) / 16384.0;
  128. var B = uSquared * (256.0 + uSquared * (uSquared * (74.0 - 47.0 * uSquared) - 128.0)) / 1024.0;
  129. var cosineSquaredTwiceSigmaMidpoint = cosineTwiceSigmaMidpoint * cosineTwiceSigmaMidpoint;
  130. var deltaSigma = B * sineSigma * (cosineTwiceSigmaMidpoint + B * (cosineSigma *
  131. (2.0 * cosineSquaredTwiceSigmaMidpoint - 1.0) - B * cosineTwiceSigmaMidpoint *
  132. (4.0 * sineSigma * sineSigma - 3.0) * (4.0 * cosineSquaredTwiceSigmaMidpoint - 3.0) / 6.0) / 4.0);
  133. var distance = minor * A * (sigma - deltaSigma);
  134. var startHeading = Math.atan2(cosineU2 * sineLambda, cs - sc * cosineLambda);
  135. var endHeading = Math.atan2(cosineU1 * sineLambda, cs * cosineLambda - sc);
  136. ellipsoidGeodesic._distance = distance;
  137. ellipsoidGeodesic._startHeading = startHeading;
  138. ellipsoidGeodesic._endHeading = endHeading;
  139. ellipsoidGeodesic._uSquared = uSquared;
  140. }
  141. var scratchCart1 = new Cartographic.Cartesian3();
  142. var scratchCart2 = new Cartographic.Cartesian3();
  143. function computeProperties(ellipsoidGeodesic, start, end, ellipsoid) {
  144. var firstCartesian = Cartographic.Cartesian3.normalize(ellipsoid.cartographicToCartesian(start, scratchCart2), scratchCart1);
  145. var lastCartesian = Cartographic.Cartesian3.normalize(ellipsoid.cartographicToCartesian(end, scratchCart2), scratchCart2);
  146. //>>includeStart('debug', pragmas.debug);
  147. Check.Check.typeOf.number.greaterThanOrEquals('value', Math.abs(Math.abs(Cartographic.Cartesian3.angleBetween(firstCartesian, lastCartesian)) - Math.PI), 0.0125);
  148. //>>includeEnd('debug');
  149. vincentyInverseFormula(ellipsoidGeodesic, ellipsoid.maximumRadius, ellipsoid.minimumRadius,
  150. start.longitude, start.latitude, end.longitude, end.latitude);
  151. ellipsoidGeodesic._start = Cartographic.Cartographic.clone(start, ellipsoidGeodesic._start);
  152. ellipsoidGeodesic._end = Cartographic.Cartographic.clone(end, ellipsoidGeodesic._end);
  153. ellipsoidGeodesic._start.height = 0;
  154. ellipsoidGeodesic._end.height = 0;
  155. setConstants(ellipsoidGeodesic);
  156. }
  157. /**
  158. * Initializes a geodesic on the ellipsoid connecting the two provided planetodetic points.
  159. *
  160. * @alias EllipsoidGeodesic
  161. * @constructor
  162. *
  163. * @param {Cartographic} [start] The initial planetodetic point on the path.
  164. * @param {Cartographic} [end] The final planetodetic point on the path.
  165. * @param {Ellipsoid} [ellipsoid=Ellipsoid.WGS84] The ellipsoid on which the geodesic lies.
  166. */
  167. function EllipsoidGeodesic(start, end, ellipsoid) {
  168. var e = when.defaultValue(ellipsoid, Cartesian2.Ellipsoid.WGS84);
  169. this._ellipsoid = e;
  170. this._start = new Cartographic.Cartographic();
  171. this._end = new Cartographic.Cartographic();
  172. this._constants = {};
  173. this._startHeading = undefined;
  174. this._endHeading = undefined;
  175. this._distance = undefined;
  176. this._uSquared = undefined;
  177. if (when.defined(start) && when.defined(end)) {
  178. computeProperties(this, start, end, e);
  179. }
  180. }
  181. Object.defineProperties(EllipsoidGeodesic.prototype, {
  182. /**
  183. * Gets the ellipsoid.
  184. * @memberof EllipsoidGeodesic.prototype
  185. * @type {Ellipsoid}
  186. * @readonly
  187. */
  188. ellipsoid : {
  189. get : function() {
  190. return this._ellipsoid;
  191. }
  192. },
  193. /**
  194. * Gets the surface distance between the start and end point
  195. * @memberof EllipsoidGeodesic.prototype
  196. * @type {Number}
  197. * @readonly
  198. */
  199. surfaceDistance : {
  200. get : function() {
  201. //>>includeStart('debug', pragmas.debug);
  202. Check.Check.defined('distance', this._distance);
  203. //>>includeEnd('debug');
  204. return this._distance;
  205. }
  206. },
  207. /**
  208. * Gets the initial planetodetic point on the path.
  209. * @memberof EllipsoidGeodesic.prototype
  210. * @type {Cartographic}
  211. * @readonly
  212. */
  213. start : {
  214. get : function() {
  215. return this._start;
  216. }
  217. },
  218. /**
  219. * Gets the final planetodetic point on the path.
  220. * @memberof EllipsoidGeodesic.prototype
  221. * @type {Cartographic}
  222. * @readonly
  223. */
  224. end : {
  225. get : function() {
  226. return this._end;
  227. }
  228. },
  229. /**
  230. * Gets the heading at the initial point.
  231. * @memberof EllipsoidGeodesic.prototype
  232. * @type {Number}
  233. * @readonly
  234. */
  235. startHeading : {
  236. get : function() {
  237. //>>includeStart('debug', pragmas.debug);
  238. Check.Check.defined('distance', this._distance);
  239. //>>includeEnd('debug');
  240. return this._startHeading;
  241. }
  242. },
  243. /**
  244. * Gets the heading at the final point.
  245. * @memberof EllipsoidGeodesic.prototype
  246. * @type {Number}
  247. * @readonly
  248. */
  249. endHeading : {
  250. get : function() {
  251. //>>includeStart('debug', pragmas.debug);
  252. Check.Check.defined('distance', this._distance);
  253. //>>includeEnd('debug');
  254. return this._endHeading;
  255. }
  256. }
  257. });
  258. /**
  259. * Sets the start and end points of the geodesic
  260. *
  261. * @param {Cartographic} start The initial planetodetic point on the path.
  262. * @param {Cartographic} end The final planetodetic point on the path.
  263. */
  264. EllipsoidGeodesic.prototype.setEndPoints = function(start, end) {
  265. //>>includeStart('debug', pragmas.debug);
  266. Check.Check.defined('start', start);
  267. Check.Check.defined('end', end);
  268. //>>includeEnd('debug');
  269. computeProperties(this, start, end, this._ellipsoid);
  270. };
  271. /**
  272. * Provides the location of a point at the indicated portion along the geodesic.
  273. *
  274. * @param {Number} fraction The portion of the distance between the initial and final points.
  275. * @param {Cartographic} result The object in which to store the result.
  276. * @returns {Cartographic} The location of the point along the geodesic.
  277. */
  278. EllipsoidGeodesic.prototype.interpolateUsingFraction = function(fraction, result) {
  279. return this.interpolateUsingSurfaceDistance(this._distance * fraction, result);
  280. };
  281. /**
  282. * Provides the location of a point at the indicated distance along the geodesic.
  283. *
  284. * @param {Number} distance The distance from the inital point to the point of interest along the geodesic
  285. * @param {Cartographic} result The object in which to store the result.
  286. * @returns {Cartographic} The location of the point along the geodesic.
  287. *
  288. * @exception {DeveloperError} start and end must be set before calling function interpolateUsingSurfaceDistance
  289. */
  290. EllipsoidGeodesic.prototype.interpolateUsingSurfaceDistance = function(distance, result) {
  291. //>>includeStart('debug', pragmas.debug);
  292. Check.Check.defined('distance', this._distance);
  293. //>>includeEnd('debug');
  294. var constants = this._constants;
  295. var s = constants.distanceRatio + distance / constants.b;
  296. var cosine2S = Math.cos(2.0 * s);
  297. var cosine4S = Math.cos(4.0 * s);
  298. var cosine6S = Math.cos(6.0 * s);
  299. var sine2S = Math.sin(2.0 * s);
  300. var sine4S = Math.sin(4.0 * s);
  301. var sine6S = Math.sin(6.0 * s);
  302. var sine8S = Math.sin(8.0 * s);
  303. var s2 = s * s;
  304. var s3 = s * s2;
  305. var u8Over256 = constants.u8Over256;
  306. var u2Over4 = constants.u2Over4;
  307. var u6Over64 = constants.u6Over64;
  308. var u4Over16 = constants.u4Over16;
  309. var sigma = 2.0 * s3 * u8Over256 * cosine2S / 3.0 +
  310. s * (1.0 - u2Over4 + 7.0 * u4Over16 / 4.0 - 15.0 * u6Over64 / 4.0 + 579.0 * u8Over256 / 64.0 -
  311. (u4Over16 - 15.0 * u6Over64 / 4.0 + 187.0 * u8Over256 / 16.0) * cosine2S -
  312. (5.0 * u6Over64 / 4.0 - 115.0 * u8Over256 / 16.0) * cosine4S -
  313. 29.0 * u8Over256 * cosine6S / 16.0) +
  314. (u2Over4 / 2.0 - u4Over16 + 71.0 * u6Over64 / 32.0 - 85.0 * u8Over256 / 16.0) * sine2S +
  315. (5.0 * u4Over16 / 16.0 - 5.0 * u6Over64 / 4.0 + 383.0 * u8Over256 / 96.0) * sine4S -
  316. s2 * ((u6Over64 - 11.0 * u8Over256 / 2.0) * sine2S + 5.0 * u8Over256 * sine4S / 2.0) +
  317. (29.0 * u6Over64 / 96.0 - 29.0 * u8Over256 / 16.0) * sine6S +
  318. 539.0 * u8Over256 * sine8S / 1536.0;
  319. var theta = Math.asin(Math.sin(sigma) * constants.cosineAlpha);
  320. var latitude = Math.atan(constants.a / constants.b * Math.tan(theta));
  321. // Redefine in terms of relative argument of latitude.
  322. sigma = sigma - constants.sigma;
  323. var cosineTwiceSigmaMidpoint = Math.cos(2.0 * constants.sigma + sigma);
  324. var sineSigma = Math.sin(sigma);
  325. var cosineSigma = Math.cos(sigma);
  326. var cc = constants.cosineU * cosineSigma;
  327. var ss = constants.sineU * sineSigma;
  328. var lambda = Math.atan2(sineSigma * constants.sineHeading, cc - ss * constants.cosineHeading);
  329. var l = lambda - computeDeltaLambda(constants.f, constants.sineAlpha, constants.cosineSquaredAlpha,
  330. sigma, sineSigma, cosineSigma, cosineTwiceSigmaMidpoint);
  331. if (when.defined(result)) {
  332. result.longitude = this._start.longitude + l;
  333. result.latitude = latitude;
  334. result.height = 0.0;
  335. return result;
  336. }
  337. return new Cartographic.Cartographic(this._start.longitude + l, latitude, 0.0);
  338. };
  339. exports.EllipsoidGeodesic = EllipsoidGeodesic;
  340. });