Sun termination line

[hpinkos]

Requested by two users on the forum: https://groups.google.com/forum/?hl=en#!topic/cesium-dev/8-H_Nb1vVNo

The day/night lighting is very nice, but it would also be nice to show where the sunlight is simply by drawing a sun termination line, which is the border between day and night, instead of relying on the shadow to demarcate that border.

[jdurrie]

Anybody know the math for the sunlight/penumbra and penumbra/umbra lines?

Here's a wikipedia article showing the geometry: link

[ghost]

• Terminator (solar)
• Solar Terminator Implementation in Javascript / D3.js

Guess you could replace the sun position code with cesium's sun position?

[jdurrie]

I can't find the actual position of the sun in cesium. Anyone know where that can be found? I could calculate it myself, but would prefer to reuse existing code if possible.

I'm trying to produce sunlight/penumbra and penumbra/umbra lines. The terminator line in those links is similar but not the same.

[ghost]

Here's how Cesium calculates the sun (and moon) positions: https://github.com/AnalyticalGraphicsInc/cesium/blob/bd3ddefcc1f38a8444d3bddfd0d093867a979454/Source/Renderer/UniformState.js#L840

[jdurrie]

Thanks. I'll see if I can figure out how to access that.

[hpinkos]

@jdurrie You might be able to do scene.context.uniformState.sunPositionWC Or you could call Simon1994PlanetaryPositions.computeSunPositionInEarthInertialFrame from whatever method you write to add the termination line

[jdurrie]

Thanks. scene.context.uniformState.sunPositionWC is what I need for sun position.

I still looking for penumbra equations, if anyone knows where I can find them. Otherwise, I'll try figuring it out from the geometry.

[ghost]

Check the image on the left on http://mysite.du.edu/~jcalvert/astro/shadows.htm Another: http://mathscinotes.com/2010/10/solar-eclipse-math/

Both are for solar eclipses... not sure if it helps you.

[jdurrie]

Awesome, I think the image on the left of the first link gives me what I want. I'll crunch some numbers.

[jdurrie]

Here's my math:

As shown on the first figure on the left in thinkpadder1's link above (http://mysite.du.edu/~jcalvert/astro/shadows.htm), the penumbra/umbra line can be calculated as the intersection of the surface of the earth and a cone with it's apex at the center of the earth and projecting away from the sun (base is farther from the sun than the apex). The interior half-angle of the cone is the other angle, I'll call a, in the right triangle containing theta, where a=pi-theta or cos(a)=(R-r)/D.

The sunlight/penumbra line is shown by the crossing lines in that figure. Let b be the angle in that figure between the line connecting the circle centers and one of the crossing lines, and c be the other angle in a right triangle containing b and the circle radius perpendicular to the crossing line (not shown). Then the sunlight/penumbra line can be described by the intersection of the earth's surface and a cone with it's apex at the center of the earth, projecting towards the sun, and an interior angle of c, where sin(b)=(r+R)/D and c=pi-b or cos(c)=(r+R)/D.

[jdurrie]

Here's an implementation for the sandcastle. There are probably performance tweaks that could be made.

var viewer = new Cesium.Viewer('cesiumContainer');

var penumbraSunlightPolygon = viewer.entities.add({
  name: 'sunlight/penumbra polygon',
  polygon: {
    hierarchy: Cesium.Cartesian3.fromDegreesArray([
        0, 0,
        1, 0,
        1, 1,
        0, 1
    ]),
    fill: false,
    outline: true,
    outlineWidth: 2,
    outlineColor: Cesium.Color.YELLOW,
    show: true
  }
});
var penumbraUmbraPolygon = viewer.entities.add({
  name: 'umbra/penumbra polygon',
  polygon: {
    hierarchy: Cesium.Cartesian3.fromDegreesArray([
        0, 0,
        1, 0,
        1, 1,
        0, 1
    ]),
    fill: false,
    outline: true,
    outlineWidth: 2,
    outlineColor: Cesium.Color.BLUE,
    show: true
  }
});
var penumbraVertexCount = 6;
var penumbraVertexRotation = 2*Cesium.Math.PI / penumbraVertexCount;
var simUpdateStep = 1; //update time step in seconds; will be scaled by clock multiplier
var sunRadius = 696.3; //mean radius in 10^6 meters
var earthRadius = 6.372797560856; //mean radius in 10^6 meters
var lastUpdateTime = new Cesium.JulianDate(0,0);
function updatePenumbraLines(clock) {
    //updating too often brings cesium to a halt
    if(Cesium.JulianDate.equalsEpsilon(lastUpdateTime, clock.currentTime, Math.abs(clock.multiplier)*simUpdateStep)) {return;}
    //sun position
    var spwc = viewer.scene.context.uniformState.sunPositionWC;
    var D = Cesium.Cartesian3.magnitude(spwc)/1000000.0; //distance from earth to sun in 10^6 meters
    if(D) {//prevent NaN errors
        lastUpdateTime = clock.currentTime;
        var sv = Cesium.Cartesian3.normalize(spwc, new Cesium.Cartesian3());
        //vertex rotation
        var q = Cesium.Quaternion.fromAxisAngle(sv, penumbraVertexRotation);
        var R = Cesium.Matrix3.fromQuaternion(q);
        //sunlight line
        var angle = Math.acos((sunRadius+earthRadius)/D); //sunlight cone inner half-angle
        var v = new Cesium.Cartesian3(0,0,1);
        Cesium.Cartesian3.cross(sv, v, v); //get vector orthogonal to sv
        var qs = Cesium.Quaternion.fromAxisAngle(v, angle);
        var Rs = Cesium.Matrix3.fromQuaternion(qs);
        Cesium.Matrix3.multiplyByVector(Rs, sv, v);
        var pos = [];
        var w;
        for (var i=0; i<penumbraVertexCount; ++i){
            Cesium.Matrix3.multiplyByVector(R, v, v);
            w = Cesium.Cartesian3.clone(v);
            viewer.scene.globe.ellipsoid.scaleToGeocentricSurface(w, w);
            pos.push(w);
        }
        penumbraSunlightPolygon.polygon.hierarchy = new Cesium.PolygonHierarchy(pos);
        //umbra line
        Cesium.Cartesian3.negate(sv, sv); //umbra cone points away from sun
        angle = Math.acos((sunRadius-earthRadius)/D); //umbra cone inner half-angle
        v = new Cesium.Cartesian3(0,0,1);
        Cesium.Cartesian3.cross(sv, v, v); //get vector orthogonal to sv
        var qu = Cesium.Quaternion.fromAxisAngle(v, angle);
        var Ru = Cesium.Matrix3.fromQuaternion(qu);
        Cesium.Matrix3.multiplyByVector(Ru, sv, v);
        pos = [];
        for (i=0; i<penumbraVertexCount; ++i){
            Cesium.Matrix3.multiplyByVector(R, v, v);
            w = Cesium.Cartesian3.clone(v);
            viewer.scene.globe.ellipsoid.scaleToGeocentricSurface(w, w);
            pos.push(w);
        }
        penumbraUmbraPolygon.polygon.hierarchy = new Cesium.PolygonHierarchy(pos);
    }
}

//add the event handler
viewer.clock.onTick.addEventListener(updatePenumbraLines);

[hpinkos]

Wow, looks great @jdurrie thanks! Using a polyline instead of a polygon outline should improve performance drastically. Even though you're only using the outline, there are a lot of extra steps that happen behind the scenes for the polygon

[jdurrie]

I tried polylines and they aren't closed for some reason, even though I'm using loop: true. Maybe being drawn under the earth surface?

var viewer = new Cesium.Viewer('cesiumContainer');

var penumbraSunlightPolyline = viewer.entities.add({
    name: 'sunlight/penumbra polygon',
    polyline: {
        positions: Cesium.Cartesian3.fromDegreesArray([
            0, 0,
            1, 0
        ]),
        width: 2,
        material: Cesium.Color.YELLOW,
        show: true,
        loop: true
    }
});
var penumbraUmbraPolyline = viewer.entities.add({
    name: 'sunlight/penumbra polygon',
    polyline: {
        positions: Cesium.Cartesian3.fromDegreesArray([
            0, 0,
            1, 0
        ]),
        width: 2,
        material: Cesium.Color.BLUE,
        show: true,
        loop: true
    }
});
var penumbraVertexCount = 6;
var penumbraVertexRotation = 2*Cesium.Math.PI / penumbraVertexCount;
var simUpdateStep = 1; //update time step in seconds; will be scaled by clock multiplier
var sunRadius = 696.3; //mean radius in 10^6 meters
var earthRadius = 6.372797560856; //mean radius in 10^6 meters
var lastUpdateTime = new Cesium.JulianDate(0,0);
function updatePenumbraLines(clock) {
    //updating too often brings cesium to a halt
    if(Cesium.JulianDate.equalsEpsilon(lastUpdateTime, clock.currentTime, Math.abs(clock.multiplier)*simUpdateStep)) {return;}
    //sun position
    var spwc = viewer.scene.context.uniformState.sunPositionWC;
    var D = Cesium.Cartesian3.magnitude(spwc)/1000000.0; //distance from earth to sun in 10^6 meters
    if(D) {//prevent NaN errors
        lastUpdateTime = clock.currentTime;
        var sv = Cesium.Cartesian3.normalize(spwc, new Cesium.Cartesian3());
        //vertex rotation
        var q = Cesium.Quaternion.fromAxisAngle(sv, penumbraVertexRotation);
        var R = Cesium.Matrix3.fromQuaternion(q);
        //sunlight line
        var angle = Math.acos((sunRadius+earthRadius)/D); //sunlight cone inner half-angle
        var v = new Cesium.Cartesian3(0,0,1);
        Cesium.Cartesian3.cross(sv, v, v); //get vector orthogonal to sv
        var qs = Cesium.Quaternion.fromAxisAngle(v, angle);
        var Rs = Cesium.Matrix3.fromQuaternion(qs);
        Cesium.Matrix3.multiplyByVector(Rs, sv, v);
        var pos = [];
        var w;
        for (var i=0; i<penumbraVertexCount; ++i){
            Cesium.Matrix3.multiplyByVector(R, v, v);
            w = Cesium.Cartesian3.clone(v);
            viewer.scene.globe.ellipsoid.scaleToGeocentricSurface(w, w);
            pos.push(w);
        }
        penumbraSunlightPolyline.polyline.positions = pos;
        //umbra line
        Cesium.Cartesian3.negate(sv, sv); //umbra cone points away from sun
        angle = Math.acos((sunRadius-earthRadius)/D); //umbra cone inner half-angle
        v = new Cesium.Cartesian3(0,0,1);
        Cesium.Cartesian3.cross(sv, v, v); //get vector orthogonal to sv
        var qu = Cesium.Quaternion.fromAxisAngle(v, angle);
        var Ru = Cesium.Matrix3.fromQuaternion(qu);
        Cesium.Matrix3.multiplyByVector(Ru, sv, v);
        pos = [];
        for (i=0; i<penumbraVertexCount; ++i){
            Cesium.Matrix3.multiplyByVector(R, v, v);
            w = Cesium.Cartesian3.clone(v);
            viewer.scene.globe.ellipsoid.scaleToGeocentricSurface(w, w);
            pos.push(w);
        }
        penumbraUmbraPolyline.polyline.positions = pos;;
    }
}

//add the event handler
viewer.clock.onTick.addEventListener(updatePenumbraLines);

[jdurrie]

Adding another point to the polylines to force the start and end points to be the same fixes the closure problem. But I'm not seeing any performance differences for the polygons vs polylines.