Skip to content

Commit 6f4fadf

Browse files
authored
Merge pull request #8986 from CesiumGS/ellipsoid-surfacearea
Added Ellipsoid surface area
2 parents 365560f + 4e730ef commit 6f4fadf

File tree

3 files changed

+170
-1
lines changed

3 files changed

+170
-1
lines changed

CHANGES.md

+3-1
Original file line numberDiff line numberDiff line change
@@ -8,11 +8,13 @@
88

99
##### Additions :tada:
1010

11+
- Added `Ellipsoid.surfaceArea` for computing the approximate surface area of a rectangle on the surface of an ellipsoid. [#8986](https://github.com/CesiumGS/cesium/pull/8986)
12+
- Added support for PolylineVolume in CZML. [#8841](https://github.com/CesiumGS/cesium/pull/8841)
1113
- Added `Color.toCssHexString` for getting the CSS hex string equivalent for a color. [#8987](https://github.com/CesiumGS/cesium/pull/8987)
12-
- Added support for PolylineVolume in CZML [#8841](https://github.com/CesiumGS/cesium/pull/8841)
1314

1415
##### Fixes :wrench:
1516

17+
- Fixed a divide-by-zero bug in `Ellipsoid.geodeticSurfaceNormal` when given the origin as input. `undefined` is returned instead. [#8986](https://github.com/CesiumGS/cesium/pull/8986)
1618
- Fixed error with `WallGeoemtry` when there were adjacent positions with very close values [#8952](https://github.com/CesiumGS/cesium/pull/8952)
1719
- Fixed artifact for skinned model when log depth is enabled. [#6447](https://github.com/CesiumGS/cesium/issues/6447)
1820
- Fixed a bug where certain rhumb arc polylines would lead to a crash. [#8787](https://github.com/CesiumGS/cesium/pull/8787)

Source/Core/Ellipsoid.js

+113
Original file line numberDiff line numberDiff line change
@@ -362,6 +362,11 @@ Ellipsoid.prototype.geodeticSurfaceNormalCartographic = function (
362362
* @returns {Cartesian3} The modified result parameter or a new Cartesian3 instance if none was provided.
363363
*/
364364
Ellipsoid.prototype.geodeticSurfaceNormal = function (cartesian, result) {
365+
if (
366+
Cartesian3.equalsEpsilon(cartesian, Cartesian3.ZERO, CesiumMath.EPSILON14)
367+
) {
368+
return undefined;
369+
}
365370
if (!defined(result)) {
366371
result = new Cartesian3();
367372
}
@@ -688,4 +693,112 @@ Ellipsoid.prototype.getSurfaceNormalIntersectionWithZAxis = function (
688693

689694
return result;
690695
};
696+
697+
var abscissas = [
698+
0.14887433898163,
699+
0.43339539412925,
700+
0.67940956829902,
701+
0.86506336668898,
702+
0.97390652851717,
703+
0.0,
704+
];
705+
var weights = [
706+
0.29552422471475,
707+
0.26926671930999,
708+
0.21908636251598,
709+
0.14945134915058,
710+
0.066671344308684,
711+
0.0,
712+
];
713+
714+
/**
715+
* Compute the 10th order Gauss-Legendre Quadrature of the given definite integral.
716+
*
717+
* @param {Number} a The lower bound for the integration.
718+
* @param {Number} b The upper bound for the integration.
719+
* @param {Ellipsoid~RealValuedScalarFunction} func The function to integrate.
720+
* @returns {Number} The value of the integral of the given function over the given domain.
721+
*
722+
* @private
723+
*/
724+
function gaussLegendreQuadrature(a, b, func) {
725+
//>>includeStart('debug', pragmas.debug);
726+
Check.typeOf.number("a", a);
727+
Check.typeOf.number("b", b);
728+
Check.typeOf.func("func", func);
729+
//>>includeEnd('debug');
730+
731+
// The range is half of the normal range since the five weights add to one (ten weights add to two).
732+
// The values of the abscissas are multiplied by two to account for this.
733+
var xMean = 0.5 * (b + a);
734+
var xRange = 0.5 * (b - a);
735+
736+
var sum = 0.0;
737+
for (var i = 0; i < 5; i++) {
738+
var dx = xRange * abscissas[i];
739+
sum += weights[i] * (func(xMean + dx) + func(xMean - dx));
740+
}
741+
742+
// Scale the sum to the range of x.
743+
sum *= xRange;
744+
return sum;
745+
}
746+
747+
/**
748+
* A real valued scalar function.
749+
* @callback Ellipsoid~RealValuedScalarFunction
750+
*
751+
* @param {Number} x The value used to evaluate the function.
752+
* @returns {Number} The value of the function at x.
753+
*
754+
* @private
755+
*/
756+
757+
/**
758+
* Computes an approximation of the surface area of a rectangle on the surface of an ellipsoid using
759+
* Gauss-Legendre 10th order quadrature.
760+
*
761+
* @param {Rectangle} rectangle The rectangle used for computing the surface area.
762+
* @returns {Number} The approximate area of the rectangle on the surface of this ellipsoid.
763+
*/
764+
Ellipsoid.prototype.surfaceArea = function (rectangle) {
765+
//>>includeStart('debug', pragmas.debug);
766+
Check.typeOf.object("rectangle", rectangle);
767+
//>>includeEnd('debug');
768+
var minLongitude = rectangle.west;
769+
var maxLongitude = rectangle.east;
770+
var minLatitude = rectangle.south;
771+
var maxLatitude = rectangle.north;
772+
773+
while (maxLongitude < minLongitude) {
774+
maxLongitude += CesiumMath.TWO_PI;
775+
}
776+
777+
var radiiSquared = this._radiiSquared;
778+
var a2 = radiiSquared.x;
779+
var b2 = radiiSquared.y;
780+
var c2 = radiiSquared.z;
781+
var a2b2 = a2 * b2;
782+
return gaussLegendreQuadrature(minLatitude, maxLatitude, function (lat) {
783+
// phi represents the angle measured from the north pole
784+
// sin(phi) = sin(pi / 2 - lat) = cos(lat), cos(phi) is similar
785+
var sinPhi = Math.cos(lat);
786+
var cosPhi = Math.sin(lat);
787+
return (
788+
Math.cos(lat) *
789+
gaussLegendreQuadrature(minLongitude, maxLongitude, function (lon) {
790+
var cosTheta = Math.cos(lon);
791+
var sinTheta = Math.sin(lon);
792+
return Math.sqrt(
793+
a2b2 * cosPhi * cosPhi +
794+
c2 *
795+
(b2 * cosTheta * cosTheta + a2 * sinTheta * sinTheta) *
796+
sinPhi *
797+
sinPhi
798+
);
799+
})
800+
);
801+
});
802+
};
803+
691804
export default Ellipsoid;

Specs/Core/EllipsoidSpec.js

+54
Original file line numberDiff line numberDiff line change
@@ -2,6 +2,7 @@ import { Cartesian3 } from "../../Source/Cesium.js";
22
import { Cartographic } from "../../Source/Cesium.js";
33
import { Ellipsoid } from "../../Source/Cesium.js";
44
import { Math as CesiumMath } from "../../Source/Cesium.js";
5+
import { Rectangle } from "../../Source/Cesium.js";
56
import createPackableSpecs from "../createPackableSpecs.js";
67

78
describe("Core/Ellipsoid", function () {
@@ -137,6 +138,12 @@ describe("Core/Ellipsoid", function () {
137138
);
138139
});
139140

141+
it("geodeticSurfaceNormal returns undefined when given the origin", function () {
142+
var ellipsoid = Ellipsoid.WGS84;
143+
var returnedResult = ellipsoid.geodeticSurfaceNormal(Cartesian3.ZERO);
144+
expect(returnedResult).toBeUndefined();
145+
});
146+
140147
it("geodeticSurfaceNormal works with a result parameter", function () {
141148
var ellipsoid = Ellipsoid.WGS84;
142149
var result = new Cartesian3();
@@ -708,6 +715,53 @@ describe("Core/Ellipsoid", function () {
708715
expect(ellipsoid._squaredXOverSquaredZ).toEqual(squaredXOverSquaredZ);
709716
});
710717

718+
it("surfaceArea throws without rectangle", function () {
719+
expect(function () {
720+
return Ellipsoid.WGS84.surfaceArea(undefined);
721+
}).toThrowDeveloperError();
722+
});
723+
724+
it("computes surfaceArea", function () {
725+
// area of an oblate spheroid
726+
var ellipsoid = new Ellipsoid(4, 4, 3);
727+
var a2 = ellipsoid.radiiSquared.x;
728+
var c2 = ellipsoid.radiiSquared.z;
729+
var e = Math.sqrt(1.0 - c2 / a2);
730+
var area =
731+
CesiumMath.TWO_PI * a2 +
732+
CesiumMath.PI * (c2 / e) * Math.log((1.0 + e) / (1.0 - e));
733+
expect(
734+
ellipsoid.surfaceArea(
735+
new Rectangle(
736+
-CesiumMath.PI,
737+
-CesiumMath.PI_OVER_TWO,
738+
CesiumMath.PI,
739+
CesiumMath.PI_OVER_TWO
740+
)
741+
)
742+
).toEqualEpsilon(area, CesiumMath.EPSILON3);
743+
744+
// area of a prolate spheroid
745+
ellipsoid = new Ellipsoid(3, 3, 4);
746+
a2 = ellipsoid.radiiSquared.x;
747+
c2 = ellipsoid.radiiSquared.z;
748+
e = Math.sqrt(1.0 - a2 / c2);
749+
var a = ellipsoid.radii.x;
750+
var c = ellipsoid.radii.z;
751+
area =
752+
CesiumMath.TWO_PI * a2 + CesiumMath.TWO_PI * ((a * c) / e) * Math.asin(e);
753+
expect(
754+
ellipsoid.surfaceArea(
755+
new Rectangle(
756+
-CesiumMath.PI,
757+
-CesiumMath.PI_OVER_TWO,
758+
CesiumMath.PI,
759+
CesiumMath.PI_OVER_TWO
760+
)
761+
)
762+
).toEqualEpsilon(area, CesiumMath.EPSILON3);
763+
});
764+
711765
createPackableSpecs(Ellipsoid, Ellipsoid.WGS84, [
712766
Ellipsoid.WGS84.radii.x,
713767
Ellipsoid.WGS84.radii.y,

0 commit comments

Comments
 (0)