/
Ellipsoid.js
568 lines (497 loc) · 22.3 KB
/
Ellipsoid.js
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
/*global define*/
define([
'./freezeObject',
'./defaultValue',
'./DeveloperError',
'./Math',
'./Cartesian3',
'./Cartographic'
], function(
freezeObject,
defaultValue,
DeveloperError,
CesiumMath,
Cartesian3,
Cartographic) {
"use strict";
/**
* A quadratic surface defined in Cartesian coordinates by the equation
* <code>(x / a)^2 + (y / b)^2 + (z / c)^2 = 1</code>. Primarily used
* by Cesium to represent the shape of planetary bodies.
*
* Rather than constructing this object directly, one of the provided
* constants is normally used.
* @alias Ellipsoid
* @constructor
* @immutable
*
* @param {Number} [x=0] The radius in the x direction.
* @param {Number} [y=0] The radius in the y direction.
* @param {Number} [z=0] The radius in the z direction.
*
* @exception {DeveloperError} All radii components must be greater than or equal to zero.
*
* @see Ellipsoid.fromCartesian3
* @see Ellipsoid.WGS84
* @see Ellipsoid.UNIT_SPHERE
*/
var Ellipsoid = function(x, y, z) {
x = defaultValue(x, 0.0);
y = defaultValue(y, 0.0);
z = defaultValue(z, 0.0);
if (x < 0.0 || y < 0.0 || z < 0.0) {
throw new DeveloperError('All radii components must be greater than or equal to zero.');
}
this._radii = new Cartesian3(x, y, z);
this._radiiSquared = new Cartesian3(x * x,
y * y,
z * z);
this._radiiToTheFourth = new Cartesian3(x * x * x * x,
y * y * y * y,
z * z * z * z);
this._oneOverRadii = new Cartesian3(x === 0.0 ? 0.0 : 1.0 / x,
y === 0.0 ? 0.0 : 1.0 / y,
z === 0.0 ? 0.0 : 1.0 / z);
this._oneOverRadiiSquared = new Cartesian3(x === 0.0 ? 0.0 : 1.0 / (x * x),
y === 0.0 ? 0.0 : 1.0 / (y * y),
z === 0.0 ? 0.0 : 1.0 / (z * z));
this._minimumRadius = Math.min(x, y, z);
this._maximumRadius = Math.max(x, y, z);
this._centerToleranceSquared = CesiumMath.EPSILON1;
};
/**
* Duplicates an Ellipsoid instance.
*
* @memberof Ellipsoid
*
* @param {Ellipsoid} ellipsoid The ellipsoid to duplicate.
* @param {Ellipsoid} [result] The object onto which to store the result, or undefined if a new
* instance should be created.
* @returns {Ellipsoid} The cloned Ellipsoid.
*/
Ellipsoid.clone = function(ellipsoid, result) {
var radii = ellipsoid._radii;
if (typeof result === 'undefined') {
return new Ellipsoid(radii.x, radii.y, radii.z);
}
Cartesian3.clone(radii, result._radii);
Cartesian3.clone(ellipsoid._radiiSquared, result._radiiSquared);
Cartesian3.clone(ellipsoid._radiiToTheFourth, result._radiiToTheFourth);
Cartesian3.clone(ellipsoid._oneOverRadii, result._oneOverRadii);
Cartesian3.clone(ellipsoid._oneOverRadiiSquared, result._oneOverRadiiSquared);
result._minimumRadius = ellipsoid._minimumRadius;
result._maximumRadius = ellipsoid._maximumRadius;
result._centerToleranceSquared = ellipsoid._centerToleranceSquared;
return result;
};
/**
* Computes an Ellipsoid from a Cartesian specifying the radii in x, y, and z directions.
*
* @param {Cartesian3} [radii=Cartesian3.ZERO] The ellipsoid's radius in the x, y, and z directions.
* @return {Ellipsoid} A new Ellipsoid instance.
*
* @exception {DeveloperError} All radii components must be greater than or equal to zero.
*
* @see Ellipsoid.WGS84
* @see Ellipsoid.UNIT_SPHERE
*/
Ellipsoid.fromCartesian3 = function(cartesian) {
if (typeof cartesian === 'undefined') {
return new Ellipsoid();
}
return new Ellipsoid(cartesian.x, cartesian.y, cartesian.z);
};
/**
* An Ellipsoid instance initialized to the WGS84 standard.
* @memberof Ellipsoid
*
* @see czm_getWgs84EllipsoidEC
*/
Ellipsoid.WGS84 = freezeObject(new Ellipsoid(6378137.0, 6378137.0, 6356752.3142451793));
/**
* An Ellipsoid instance initialized to radii of (1.0, 1.0, 1.0).
* @memberof Ellipsoid
*/
Ellipsoid.UNIT_SPHERE = freezeObject(new Ellipsoid(1.0, 1.0, 1.0));
/**
* @memberof Ellipsoid
* @return {Cartesian3} The radii of the ellipsoid.
*/
Ellipsoid.prototype.getRadii = function() {
return this._radii;
};
/**
* @memberof Ellipsoid
* @return {Cartesian3} The squared radii of the ellipsoid.
*/
Ellipsoid.prototype.getRadiiSquared = function() {
return this._radiiSquared;
};
/**
* @memberof Ellipsoid
* @return {Cartesian3} The radii of the ellipsoid raised to the fourth power.
*/
Ellipsoid.prototype.getRadiiToTheFourth = function() {
return this._radiiToTheFourth;
};
/**
* @memberof Ellipsoid
* @return {Cartesian3} One over the radii of the ellipsoid.
*/
Ellipsoid.prototype.getOneOverRadii = function() {
return this._oneOverRadii;
};
/**
* @memberof Ellipsoid
* @return {Cartesian3} One over the squared radii of the ellipsoid.
*/
Ellipsoid.prototype.getOneOverRadiiSquared = function() {
return this._oneOverRadiiSquared;
};
/**
* @memberof Ellipsoid
* @return {Cartesian3} The minimum radius of the ellipsoid.
*/
Ellipsoid.prototype.getMinimumRadius = function() {
return this._minimumRadius;
};
/**
* @memberof Ellipsoid
* @return {Cartesian3} The maximum radius of the ellipsoid.
*/
Ellipsoid.prototype.getMaximumRadius = function() {
return this._maximumRadius;
};
/**
* Computes the unit vector directed from the center of this ellipsoid toward the provided Cartesian position.
* @memberof Ellipsoid
*
* @param {Cartesian3} cartesian The Cartesian for which to to determine the geocentric normal.
* @param {Cartesian3} [result] The object onto which to store the result.
* @return {Cartesian3} The modified result parameter or a new Cartesian3 instance if none was provided.
*
* @exception {DeveloperError} cartesian is required.
*/
Ellipsoid.prototype.geocentricSurfaceNormal = Cartesian3.normalize;
/**
* Computes the normal of the plane tangent to the surface of the ellipsoid at the provided position.
* @memberof Ellipsoid
*
* @param {Cartographic} cartographic The cartographic position for which to to determine the geodetic normal.
* @param {Cartesian3} [result] The object onto which to store the result.
* @return {Cartesian3} The modified result parameter or a new Cartesian3 instance if none was provided.
*
* @exception {DeveloperError} cartographic is required.
*/
Ellipsoid.prototype.geodeticSurfaceNormalCartographic = function(cartographic, result) {
if (typeof cartographic === 'undefined') {
throw new DeveloperError('cartographic is required.');
}
var longitude = cartographic.longitude;
var latitude = cartographic.latitude;
var cosLatitude = Math.cos(latitude);
var x = cosLatitude * Math.cos(longitude);
var y = cosLatitude * Math.sin(longitude);
var z = Math.sin(latitude);
if (typeof result === 'undefined') {
result = new Cartesian3();
}
result.x = x;
result.y = y;
result.z = z;
return Cartesian3.normalize(result, result);
};
/**
* Computes the normal of the plane tangent to the surface of the ellipsoid at the provided position.
* @memberof Ellipsoid
*
* @param {Cartesian3} cartesian The Cartesian position for which to to determine the surface normal.
* @param {Cartesian3} [result] The object onto which to store the result.
* @return {Cartesian3} The modified result parameter or a new Cartesian3 instance if none was provided.
*
* @exception {DeveloperError} cartesian is required.
*/
Ellipsoid.prototype.geodeticSurfaceNormal = function(cartesian, result) {
result = Cartesian3.multiplyComponents(cartesian, this._oneOverRadiiSquared, result);
return Cartesian3.normalize(result, result);
};
var cartographicToCartesianNormal = new Cartesian3();
var cartographicToCartesianK = new Cartesian3();
/**
* Converts the provided cartographic to Cartesian representation.
* @memberof Ellipsoid
*
* @param {Cartographic} cartographic The cartographic position.
* @param {Cartesian3} [result] The object onto which to store the result.
* @return {Cartesian3} The modified result parameter or a new Cartesian3 instance if none was provided.
*
* @exception {DeveloperError} cartographic is required.
*
* @example
* //Create a Cartographic and determine it's Cartesian representation on a WGS84 ellipsoid.
* var position = new Cartographic(Math.toRadians(21), Math.toRadians(78), 5000);
* var cartesianPosition = Ellipsoid.WGS84.cartographicToCartesian(position);
*/
Ellipsoid.prototype.cartographicToCartesian = function(cartographic, result) {
//`cartographic is required` is thrown from geodeticSurfaceNormalCartographic.
var n = cartographicToCartesianNormal;
var k = cartographicToCartesianK;
this.geodeticSurfaceNormalCartographic(cartographic, n);
Cartesian3.multiplyComponents(this._radiiSquared, n, k);
var gamma = Math.sqrt(Cartesian3.dot(n, k));
Cartesian3.divideByScalar(k, gamma, k);
Cartesian3.multiplyByScalar(n, cartographic.height, n);
return Cartesian3.add(k, n, result);
};
/**
* Converts the provided array of cartographics to an array of Cartesians.
* @memberof Ellipsoid
*
* @param {Array} cartographics An array of cartographic positions.
* @param {Array} [result] The object onto which to store the result.
* @return {Array} The modified result parameter or a new Array instance if none was provided.
*
* @exception {DeveloperError} cartographics is required.
*
* @example
* //Convert an array of Cartographics and determine their Cartesian representation on a WGS84 ellipsoid.
* var positions = [new Cartographic(Math.toRadians(21), Math.toRadians(78), 0),
* new Cartographic(Math.toRadians(21.321), Math.toRadians(78.123), 100),
* new Cartographic(Math.toRadians(21.645), Math.toRadians(78.456), 250)
* var cartesianPositions = Ellipsoid.WGS84.cartographicArrayToCartesianArray(positions);
*/
Ellipsoid.prototype.cartographicArrayToCartesianArray = function(cartographics, result) {
if (typeof cartographics === 'undefined') {
throw new DeveloperError('cartographics is required.');
}
var length = cartographics.length;
if (typeof result === 'undefined') {
result = new Array(length);
} else {
result.length = length;
}
for ( var i = 0; i < length; i++) {
result[i] = this.cartographicToCartesian(cartographics[i], result[i]);
}
return result;
};
var cartesianToCartographicN = new Cartesian3();
var cartesianToCartographicP = new Cartesian3();
var cartesianToCartographicH = new Cartesian3();
/**
* Converts the provided cartesian to cartographic representation.
* The cartesian is undefined at the center of the ellipsoid.
* @memberof Ellipsoid
*
* @param {Cartesian3} cartesian The Cartesian position to convert to cartographic representation.
* @param {Cartographic} [result] The object onto which to store the result.
* @return {Cartographic} The modified result parameter, new Cartographic instance if none was provided, or undefined if the cartesian is at the center of the ellipsoid.
*
* @exception {DeveloperError} cartesian is required.
*
* @example
* //Create a Cartesian and determine it's Cartographic representation on a WGS84 ellipsoid.
* var position = new Cartesian(17832.12, 83234.52, 952313.73);
* var cartographicPosition = Ellipsoid.WGS84.cartesianToCartographic(position);
*/
Ellipsoid.prototype.cartesianToCartographic = function(cartesian, result) {
//`cartesian is required.` is thrown from scaleToGeodeticSurface
var p = this.scaleToGeodeticSurface(cartesian, cartesianToCartographicP);
if (typeof p === 'undefined') {
return undefined;
}
var n = this.geodeticSurfaceNormal(p, cartesianToCartographicN);
var h = Cartesian3.subtract(cartesian, p, cartesianToCartographicH);
var longitude = Math.atan2(n.y, n.x);
var latitude = Math.asin(n.z);
var height = CesiumMath.sign(Cartesian3.dot(h, cartesian)) * Cartesian3.magnitude(h);
if (typeof result === 'undefined') {
return new Cartographic(longitude, latitude, height);
}
result.longitude = longitude;
result.latitude = latitude;
result.height = height;
return result;
};
/**
* Converts the provided array of cartesians to an array of cartographics.
* @memberof Ellipsoid
*
* @param {Array} cartesians An array of Cartesian positions.
* @param {Array} [result] The object onto which to store the result.
* @return {Array} The modified result parameter or a new Array instance if none was provided.
*
* @exception {DeveloperError} cartesians is required.
*
* @example
* //Create an array of Cartesians and determine their Cartographic representation on a WGS84 ellipsoid.
* var positions = [new Cartesian(17832.12, 83234.52, 952313.73),
* new Cartesian(17832.13, 83234.53, 952313.73),
* new Cartesian(17832.14, 83234.54, 952313.73)]
* var cartographicPositions = Ellipsoid.WGS84.cartesianArrayToCartographicArray(positions);
*/
Ellipsoid.prototype.cartesianArrayToCartographicArray = function(cartesians, result) {
if (typeof cartesians === 'undefined') {
throw new DeveloperError('cartesians is required.');
}
var length = cartesians.length;
if (typeof result === 'undefined') {
result = new Array(length);
} else {
result.length = length;
}
for ( var i = 0; i < length; ++i) {
result[i] = this.cartesianToCartographic(cartesians[i], result[i]);
}
return result;
};
var scaleToGeodeticSurfaceIntersection;
var scaleToGeodeticSurfaceGradient = new Cartesian3();
/**
* Scales the provided Cartesian position along the geodetic surface normal
* so that it is on the surface of this ellipsoid. If the position is
* at the center of the ellipsoid, this function returns undefined.
* @memberof Ellipsoid
*
* @param {Cartesian3} cartesian The Cartesian position to scale.
* @param {Cartesian3} [result] The object onto which to store the result.
* @return {Cartesian3} The modified result parameter, a new Cartesian3 instance if none was provided, or undefined if the position is at the center.
*
* @exception {DeveloperError} cartesian is required.
*/
Ellipsoid.prototype.scaleToGeodeticSurface = function(cartesian, result) {
if (typeof cartesian === 'undefined') {
throw new DeveloperError('cartesian is required.');
}
var positionX = cartesian.x;
var positionY = cartesian.y;
var positionZ = cartesian.z;
var oneOverRadii = this._oneOverRadii;
var oneOverRadiiX = oneOverRadii.x;
var oneOverRadiiY = oneOverRadii.y;
var oneOverRadiiZ = oneOverRadii.z;
var x2 = positionX * positionX * oneOverRadiiX * oneOverRadiiX;
var y2 = positionY * positionY * oneOverRadiiY * oneOverRadiiY;
var z2 = positionZ * positionZ * oneOverRadiiZ * oneOverRadiiZ;
// Compute the squared ellipsoid norm.
var squaredNorm = x2 + y2 + z2;
var ratio = Math.sqrt(1.0 / squaredNorm);
// As an initial approximation, assume that the radial intersection is the projection point.
var intersection = Cartesian3.multiplyByScalar(cartesian, ratio, scaleToGeodeticSurfaceIntersection);
//* If the position is near the center, the iteration will not converge.
if (squaredNorm < this._centerToleranceSquared) {
return !isFinite(ratio) ? undefined : Cartesian3.clone(intersection, result);
}
var oneOverRadiiSquared = this._oneOverRadiiSquared;
var oneOverRadiiSquaredX = oneOverRadiiSquared.x;
var oneOverRadiiSquaredY = oneOverRadiiSquared.y;
var oneOverRadiiSquaredZ = oneOverRadiiSquared.z;
// Use the gradient at the intersection point in place of the true unit normal.
// The difference in magnitude will be absorbed in the multiplier.
var gradient = scaleToGeodeticSurfaceGradient;
gradient.x = intersection.x * oneOverRadiiSquaredX * 2.0;
gradient.y = intersection.y * oneOverRadiiSquaredY * 2.0;
gradient.z = intersection.z * oneOverRadiiSquaredZ * 2.0;
// Compute the initial guess at the normal vector multiplier, lambda.
var lambda = (1.0 - ratio) * Cartesian3.magnitude(cartesian) / (0.5 * Cartesian3.magnitude(gradient));
var correction = 0.0;
var func;
var denominator;
var xMultiplier;
var yMultiplier;
var zMultiplier;
var xMultiplier2;
var yMultiplier2;
var zMultiplier2;
var xMultiplier3;
var yMultiplier3;
var zMultiplier3;
do {
lambda -= correction;
xMultiplier = 1.0 / (1.0 + lambda * oneOverRadiiSquaredX);
yMultiplier = 1.0 / (1.0 + lambda * oneOverRadiiSquaredY);
zMultiplier = 1.0 / (1.0 + lambda * oneOverRadiiSquaredZ);
xMultiplier2 = xMultiplier * xMultiplier;
yMultiplier2 = yMultiplier * yMultiplier;
zMultiplier2 = zMultiplier * zMultiplier;
xMultiplier3 = xMultiplier2 * xMultiplier;
yMultiplier3 = yMultiplier2 * yMultiplier;
zMultiplier3 = zMultiplier2 * zMultiplier;
func = x2 * xMultiplier2 + y2 * yMultiplier2 + z2 * zMultiplier2 - 1.0;
// "denominator" here refers to the use of this expression in the velocity and acceleration
// computations in the sections to follow.
denominator = x2 * xMultiplier3 * oneOverRadiiSquaredX + y2 * yMultiplier3 * oneOverRadiiSquaredY + z2 * zMultiplier3 * oneOverRadiiSquaredZ;
var derivative = -2.0 * denominator;
correction = func / derivative;
} while (Math.abs(func) > CesiumMath.EPSILON12);
if (typeof result === 'undefined') {
return new Cartesian3(positionX * xMultiplier, positionY * yMultiplier, positionZ * zMultiplier);
}
result.x = positionX * xMultiplier;
result.y = positionY * yMultiplier;
result.z = positionZ * zMultiplier;
return result;
};
/**
* Scales the provided Cartesian position along the geocentric surface normal
* so that it is on the surface of this ellipsoid.
* @memberof Ellipsoid
*
* @param {Cartesian3} cartesian The Cartesian position to scale.
* @param {Cartesian3} [result] The object onto which to store the result.
* @return {Cartesian3} The modified result parameter or a new Cartesian3 instance if none was provided.
*
* @exception {DeveloperError} cartesian is required.
*/
Ellipsoid.prototype.scaleToGeocentricSurface = function(cartesian, result) {
if (typeof cartesian === 'undefined') {
throw new DeveloperError('cartesian is required.');
}
var positionX = cartesian.x;
var positionY = cartesian.y;
var positionZ = cartesian.z;
var oneOverRadiiSquared = this._oneOverRadiiSquared;
var beta = 1.0 / Math.sqrt((positionX * positionX) * oneOverRadiiSquared.x +
(positionY * positionY) * oneOverRadiiSquared.y +
(positionZ * positionZ) * oneOverRadiiSquared.z);
return Cartesian3.multiplyByScalar(cartesian, beta, result);
};
/**
* Transforms a Cartesian X, Y, Z position to the ellipsoid-scaled space by multiplying
* its components by the result of {@link Ellipsoid#getOneOverRadii}.
*
* @memberof Ellipsoid
*
* @param {Cartesian3} position The position to transform.
* @param {Cartesian3} [result] The position to which to copy the result, or undefined to create and
* return a new instance.
* @returns {Cartesian3} The position expressed in the scaled space. The returned instance is the
* one passed as the result parameter if it is not undefined, or a new instance of it is.
*/
Ellipsoid.prototype.transformPositionToScaledSpace = function(position, result) {
return Cartesian3.multiplyComponents(position, this._oneOverRadii, result);
};
/**
* Compares this Ellipsoid against the provided Ellipsoid componentwise and returns
* <code>true</code> if they are equal, <code>false</code> otherwise.
* @memberof Ellipsoid
*
* @param {Ellipsoid} [right] The other Ellipsoid.
* @return {Boolean} <code>true</code> if they are equal, <code>false</code> otherwise.
*/
Ellipsoid.prototype.equals = function(right) {
return (this === right) ||
(typeof right !== 'undefined' &&
Cartesian3.equals(this._radii, right._radii));
};
/**
* Creates a string representing this Ellipsoid in the format '(radii.x, radii.y, radii.z)'.
* @memberof Ellipsoid
*
* @return {String} A string representing this ellipsoid in the format '(radii.x, radii.y, radii.z)'.
*/
Ellipsoid.prototype.toString = function() {
return this._radii.toString();
};
return Ellipsoid;
});