Skip to content

Commit

Permalink
Merge pull request #407 from AnalyticalGraphicsInc/scaleToGeodeticSur…
Browse files Browse the repository at this point in the history
…face

Fixes for Ellipsoid.scaleToGeodeticSurface
  • Loading branch information
kring committed Jan 3, 2013
2 parents ef0cfc0 + 3f0b8c5 commit d5b082c
Show file tree
Hide file tree
Showing 5 changed files with 131 additions and 72 deletions.
133 changes: 79 additions & 54 deletions Source/Core/Ellipsoid.js
Expand Up @@ -66,6 +66,8 @@ define([
this._minimumRadius = Math.min(x, y, z);

this._maximumRadius = Math.max(x, y, z);

this._centerToleranceSquared = CesiumMath.EPSILON1;
};

/**
Expand Down Expand Up @@ -285,11 +287,12 @@ define([

/**
* 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 or a new Cartographic instance if none was provided.
* @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.
*
Expand All @@ -301,6 +304,11 @@ define([
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);

Expand Down Expand Up @@ -351,14 +359,18 @@ define([
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.
* 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 or a new Cartesian3 instance if none was provided.
* @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.
*/
Expand All @@ -371,74 +383,87 @@ define([
var positionY = cartesian.y;
var positionZ = cartesian.z;

var oneOverRadiiSquared = this._oneOverRadiiSquared;
var oneOverRadiiSquaredX = oneOverRadiiSquared.x;
var oneOverRadiiSquaredY = oneOverRadiiSquared.y;
var oneOverRadiiSquaredZ = oneOverRadiiSquared.z;

var radiiSquared = this._radiiSquared;
var radiiSquaredX = radiiSquared.x;
var radiiSquaredY = radiiSquared.y;
var radiiSquaredZ = radiiSquared.z;

var radiiToTheFourth = this._radiiToTheFourth;
var radiiToTheFourthX = radiiToTheFourth.x;
var radiiToTheFourthY = radiiToTheFourth.y;
var radiiToTheFourthZ = radiiToTheFourth.z;
var oneOverRadii = this._oneOverRadii;
var oneOverRadiiX = oneOverRadii.x;
var oneOverRadiiY = oneOverRadii.y;
var oneOverRadiiZ = oneOverRadii.z;

var beta = 1.0 / Math.sqrt(
(positionX * positionX) * oneOverRadiiSquaredX +
(positionY * positionY) * oneOverRadiiSquaredY +
(positionZ * positionZ) * oneOverRadiiSquaredZ);
var x2 = positionX * positionX * oneOverRadiiX * oneOverRadiiX;
var y2 = positionY * positionY * oneOverRadiiY * oneOverRadiiY;
var z2 = positionZ * positionZ * oneOverRadiiZ * oneOverRadiiZ;

var x = beta * positionX * oneOverRadiiSquaredX;
var y = beta * positionY * oneOverRadiiSquaredY;
var z = beta * positionZ * oneOverRadiiSquaredZ;
// Compute the squared ellipsoid norm.
var squaredNorm = x2 + y2 + z2;
var ratio = Math.sqrt(1.0 / squaredNorm);

var n = Math.sqrt(x * x + y * y + z * z);
var alpha = (1.0 - beta) * (Cartesian3.magnitude(cartesian) / n);
// As an initial approximation, assume that the radial intersection is the projection point.
var intersection = Cartesian3.multiplyByScalar(cartesian, ratio, scaleToGeodeticSurfaceIntersection);

var x2 = positionX * positionX;
var y2 = positionY * positionY;
var z2 = positionZ * positionZ;
//* 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 da = 0.0;
var db = 0.0;
var dc = 0.0;
var oneOverRadiiSquared = this._oneOverRadiiSquared;
var oneOverRadiiSquaredX = oneOverRadiiSquared.x;
var oneOverRadiiSquaredY = oneOverRadiiSquared.y;
var oneOverRadiiSquaredZ = oneOverRadiiSquared.z;

var s = 0.0;
var dSdA = 1.0;
// 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 {
alpha -= (s / dSdA);
lambda -= correction;

xMultiplier = 1.0 / (1.0 + lambda * oneOverRadiiSquaredX);
yMultiplier = 1.0 / (1.0 + lambda * oneOverRadiiSquaredY);
zMultiplier = 1.0 / (1.0 + lambda * oneOverRadiiSquaredZ);

da = 1.0 + (alpha * oneOverRadiiSquaredX);
db = 1.0 + (alpha * oneOverRadiiSquaredY);
dc = 1.0 + (alpha * oneOverRadiiSquaredZ);
xMultiplier2 = xMultiplier * xMultiplier;
yMultiplier2 = yMultiplier * yMultiplier;
zMultiplier2 = zMultiplier * zMultiplier;

var da2 = da * da;
var db2 = db * db;
var dc2 = dc * dc;
xMultiplier3 = xMultiplier2 * xMultiplier;
yMultiplier3 = yMultiplier2 * yMultiplier;
zMultiplier3 = zMultiplier2 * zMultiplier;

var da3 = da * da2;
var db3 = db * db2;
var dc3 = dc * dc2;
func = x2 * xMultiplier2 + y2 * yMultiplier2 + z2 * zMultiplier2 - 1.0;

s = x2 / (radiiSquaredX * da2) + y2 / (radiiSquaredY * db2) + z2 / (radiiSquaredZ * dc2) - 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;

dSdA = -2.0 * (x2 / (radiiToTheFourthX * da3) + y2 / (radiiToTheFourthY * db3) + z2 / (radiiToTheFourthZ * dc3));
} while (Math.abs(s) > CesiumMath.EPSILON10);
var derivative = -2.0 * denominator;

x = positionX / da;
y = positionY / db;
z = positionZ / dc;
correction = func / derivative;
} while (Math.abs(func) > CesiumMath.EPSILON12);

if (typeof result === 'undefined') {
return new Cartesian3(x, y, z);
return new Cartesian3(positionX * xMultiplier, positionY * yMultiplier, positionZ * zMultiplier);
}
result.x = x;
result.y = y;
result.z = z;
result.x = positionX * xMultiplier;
result.y = positionY * yMultiplier;
result.z = positionZ * zMultiplier;
return result;
};

Expand Down
5 changes: 5 additions & 0 deletions Source/Core/EllipsoidTangentPlane.js
Expand Up @@ -24,13 +24,15 @@ define([
/**
* A plane tangent to the provided ellipsoid at the provided origin.
* If origin is not on the surface of the ellipsoid, it's surface projection will be used.
* If origin as at the center of the ellipsoid, an exception will be thrown.
* @alias EllipsoidTangentPlane
* @constructor
*
* @param {Ellipsoid} ellipsoid The ellipsoid to use.
* @param {Cartesian3} origin The point on the surface of the ellipsoid where the tangent plane touches.
*
* @exception {DeveloperError} origin is required.
* @exception {DeveloperError} origin must not be at the center of the ellipsoid.
*/
var EllipsoidTangentPlane = function(origin, ellipsoid) {
if (typeof origin === 'undefined') {
Expand All @@ -40,6 +42,9 @@ define([
ellipsoid = defaultValue(ellipsoid, Ellipsoid.WGS84);

origin = ellipsoid.scaleToGeodeticSurface(origin);
if (typeof origin === 'undefined') {
throw new DeveloperError('origin must not be at the center of the ellipsoid.');
}
var eastNorthUp = Transforms.eastNorthUpToFixedFrame(origin, ellipsoid);
this._ellipsoid = ellipsoid;
this._origin = Cartesian3.clone(origin);
Expand Down
24 changes: 13 additions & 11 deletions Source/Scene/Billboard.js
Expand Up @@ -6,6 +6,7 @@ define([
'../Core/Cartesian2',
'../Core/Cartesian3',
'../Core/Cartesian4',
'../Core/Cartographic',
'../Core/Math',
'./HorizontalOrigin',
'./VerticalOrigin',
Expand All @@ -17,6 +18,7 @@ define([
Cartesian2,
Cartesian3,
Cartesian4,
Cartographic,
CesiumMath,
HorizontalOrigin,
VerticalOrigin,
Expand Down Expand Up @@ -548,36 +550,36 @@ define([
};

var tempCartesian4 = new Cartesian4();
var tempCartographic = new Cartographic();
Billboard._computeActualPosition = function(position, frameState, morphTime, modelMatrix) {
var mode = frameState.mode;

if (mode === SceneMode.SCENE3D) {
return position;
}

var projection = frameState.scene2D.projection;
var cartographic, projectedPosition;

modelMatrix.multiplyByPoint(position, tempCartesian4);

if (mode === SceneMode.MORPHING) {
cartographic = projection.getEllipsoid().cartesianToCartographic(tempCartesian4);
projectedPosition = projection.project(cartographic);
var projection = frameState.scene2D.projection;
var cartographic = projection.getEllipsoid().cartesianToCartographic(tempCartesian4, tempCartographic);
if (typeof cartographic === 'undefined') {
return undefined;
}

var projectedPosition = projection.project(cartographic);
if (mode === SceneMode.MORPHING) {
var x = CesiumMath.lerp(projectedPosition.z, tempCartesian4.x, morphTime);
var y = CesiumMath.lerp(projectedPosition.x, tempCartesian4.y, morphTime);
var z = CesiumMath.lerp(projectedPosition.y, tempCartesian4.z, morphTime);
return new Cartesian3(x, y, z);
}

cartographic = projection.getEllipsoid().cartesianToCartographic(tempCartesian4);
projectedPosition = projection.project(cartographic);

if (mode === SceneMode.SCENE2D) {
return new Cartesian3(0.0, projectedPosition.x, projectedPosition.y);
} else if (mode === SceneMode.COLUMBUS_VIEW) {
}
if (mode === SceneMode.COLUMBUS_VIEW) {
return new Cartesian3(projectedPosition.z, projectedPosition.x, projectedPosition.y);
}
return undefined;
};

Billboard._computeScreenSpacePosition = function(modelMatrix, position, eyeOffset, pixelOffset, clampToPixel, uniformState) {
Expand Down
16 changes: 9 additions & 7 deletions Source/Scene/BillboardCollection.js
Expand Up @@ -848,17 +848,19 @@ define([
}

var length = billboards.length;
var positions = new Array(length);
for (var i = 0; i < length; ++i) {
var positions = [];
for ( var i = 0; i < length; ++i) {
var billboard = billboards[i];
var position = billboard.getPosition();
var actualPosition = Billboard._computeActualPosition(position, frameState, morphTime, modelMatrix);
billboard._setActualPosition(actualPosition);
if (typeof actualPosition !== 'undefined') {
billboard._setActualPosition(actualPosition);

if (recomputeBoundingVolume) {
positions[i] = actualPosition;
} else {
boundingVolume.expand(actualPosition, boundingVolume);
if (recomputeBoundingVolume) {
positions.push(actualPosition);
} else {
boundingVolume.expand(actualPosition, boundingVolume);
}
}
}

Expand Down
25 changes: 25 additions & 0 deletions Specs/Core/EllipsoidSpec.js
Expand Up @@ -150,6 +150,24 @@ defineSuite([
expect(returnedResult).toEqualEpsilon(surfaceCartographic, CesiumMath.EPSILON8);
});

it('cartesianToCartographic works close to center', function() {
var expected = new Cartographic(9.999999999999999e-11, 1.0067394967422763e-20, -6378137.0);
var returnedResult = Ellipsoid.WGS84.cartesianToCartographic(new Cartesian3(1e-50, 1e-60, 1e-70));
expect(returnedResult).toEqual(expected);
});

it('cartesianToCartographic return undefined very close to center', function() {
var ellipsoid = Ellipsoid.WGS84;
var returnedResult = ellipsoid.cartesianToCartographic(new Cartesian3(1e-150, 1e-150, 1e-150));
expect(returnedResult).toBeUndefined();
});

it('cartesianToCartographic return undefined at center', function() {
var ellipsoid = Ellipsoid.WGS84;
var returnedResult = ellipsoid.cartesianToCartographic(Cartesian3.ZERO);
expect(returnedResult).toBeUndefined();
});

it('cartesianArrayToCartographicArray works without a result parameter', function() {
var ellipsoid = Ellipsoid.WGS84;
var returnedResult = ellipsoid.cartesianArrayToCartographicArray([spaceCartesian, surfaceCartesian]);
Expand Down Expand Up @@ -254,6 +272,13 @@ defineSuite([
expect(result).toEqualEpsilon(expected, CesiumMath.EPSILON16);
});

it('scaleToGeodeticSurface returns undefined at center', function() {
var ellipsoid = new Ellipsoid(1.0, 2.0, 3.0);
var cartesian = new Cartesian3(0.0, 0.0, 0.0);
var returnedResult = ellipsoid.scaleToGeodeticSurface(cartesian);
expect(returnedResult).toBeUndefined();
});

it('equals works in all cases', function() {
var ellipsoid = new Ellipsoid(1.0, 0.0, 0.0);
expect(ellipsoid.equals(new Ellipsoid(1.0, 0.0, 0.0))).toEqual(true);
Expand Down

0 comments on commit d5b082c

Please sign in to comment.