globe.yaml

cameras:
    perspective:
        type: perspective
        vanishing_point: [0, 0]
        fov: .00001 # perspective correction
styles:
    projection:
        mix: stereo
    stereo:
        blend_order: 0
        shaders:
            defines:
                EARTH_RADIUS: 6378137.0 //radius of ellipsoid, WGS84
                PI: 3.14159265358979323846
                HALF_PI: 1.570796327
                deg2rad(d): (d)*PI/180.0
                rad2deg(d): (d)*180.0/PI
            blocks:
                global: |
                    // http://wiki.openstreetmap.org/wiki/Mercator
                    float y2lat_m (float y) { return 2.0*atan(exp((y/EARTH_RADIUS)))-HALF_PI; }
                    float x2lon_m (float x) { return x/EARTH_RADIUS; }

                    // rotation matrix transformations
                    mat3 globe_rotateX3D(float phi){
                        return mat3(
                            vec3(1.,0.,0.),
                            vec3(0.,cos(phi),-sin(phi)),
                            vec3(0.,sin(phi),cos(phi))
                            );
                    }
                    mat3 globe_rotateY3D(float theta){
                        return mat3(
                            vec3(cos(theta),0.,-sin(theta)),
                            vec3(0.,1.,0.),
                            vec3(sin(theta),0.,cos(theta)));
                    }
                    mat3 globe_rotateZ3D(float psi){
                        return mat3(
                            vec3(cos(psi),-sin(psi),0.),
                            vec3(sin(psi),cos(psi),0.),
                            vec3(0.,0.,1.)
                            );
                    }

                    // convert the positions from a lat, lon to a position on a sphere.
                    vec3 latLongToVector3(float lat, float lon, float radius) {
                        float phi = lat;
                        float theta = lon;

                        float x = radius * cos(phi) * cos(theta);
                        float y = radius * cos(phi) * sin(theta);
                        float z = radius * sin(phi);
                        // u_map_position.z;
                        return vec3(x,z,y);
                    }

                    float haversine(float lat1, float lon1, float lat2, float lon2) {
                        float phi1 = deg2rad(lat1);
                        float phi2 = deg2rad(lat2);
                        float deltaphi = deg2rad(lat2-lat1);
                        float deltalambda = deg2rad(lon2-lon1);

                        float a = sin(deltaphi/2.0) * sin(deltaphi/2.0) +
                                cos(phi1) * cos(phi2) *
                                sin(deltalambda/2.0) * sin(deltalambda/2.0);
                        float c = 2.0 * atan(sqrt(a), sqrt(1.0-a));
                        //float c = 2.0 * atan(sqrt(1.0-a), sqrt(a));

                        float d = EARTH_RADIUS * c;
                        return d;
                    }

                position: |
                    // mercator position of the current vertex, u_map_position = center of screen,
                    // position.xy = vertex screen position in meters from the center of the screen
                    vec2 mercator = u_map_position.xy + position.xy;
                    float lat = y2lat_m(mercator.y);
                    float lon = x2lon_m(mercator.x);

                    // Latitude_Of_Origin
                    float centerlat = y2lat_m(u_map_position.y);
                    // Central_Meridian
                    float centerlon = x2lon_m(u_map_position.x);

                    // bend map into globe
                    position.xyz = latLongToVector3(lat, lon, 2.) * EARTH_RADIUS;

                    // rotate globe with map navigation
                    // 90 degrees in radians = 1.5708
                    position.xyz *= globe_rotateY3D((-centerlon - 1.5708));
                    position.xyz *= globe_rotateX3D(-(centerlat));

                    // only show the front side of the globe:
                    // if the great circle distance between the center of the view and the current vertex's position is greater than the distance to the horizon, snap the point to the horizon so it doesn't show through the globe
                    // use haversine formula to calculate great circle distance
                    float haversine_distance = haversine(rad2deg(lat), rad2deg(lon), rad2deg(centerlat), rad2deg(centerlon));
                    float quarter_circumference = EARTH_RADIUS * PI / 2.0;
                    // if the point is beyond the horizon
                    if (haversine_distance > quarter_circumference) {
                        // if it's closer to the opposite point on the earth than to the horizon
                        if (haversine_distance > quarter_circumference * 1.75) {
                            // clip it
                            position.xy = vec2(3.402823466e+38);
                        } else {
                            // move point to horizon
                            // position.xy may be thought of as a vector from center of screen to current position, in meters
                            // get difference between distance from center of screen to current position and distance from center of screen to the horizon
                            float factor = EARTH_RADIUS / length(position.xy);
                            // multiply by factor to snap point to horizon
                            position.xy *= vec2(factor * 2.);
                        }
                    }