Functions for interacting with d3.js globe

versor-1.8.1.js

// Version 1.8.1 Copyright 2020 Rohit Kulkarni
(function(global, factory) {
  typeof exports === 'object' && typeof module !== 'undefined' ? module.exports = factory() :
  typeof define === 'function' && define.amd ? define(factory) :
  (global.versor = factory());
}(this, (function() {'use strict';

var acos = Math.acos,
    asin = Math.asin,
    atan2 = Math.atan2,
    cos = Math.cos,
    max = Math.max,
    min = Math.min,
    PI = Math.PI,
    sin = Math.sin,
    pow = Math.pow,
    sqrt = Math.sqrt,
    radians = PI / 180,
    degrees = 180 / PI;

// This function computes the cross product of two vectors v0 & v1
function cross(v0, v1) {
  return [
    v0[1] * v1[2] - v0[2] * v1[1],
    v0[2] * v1[0] - v0[0] * v1[2],
    v0[0] * v1[1] - v0[1] * v1[0]
  ];
}

// This function computes the dot product of two vectors v0 & v1
function dot(v0, v1) {
  return v0[0] * v1[0] + v0[1] * v1[1] + v0[2] * v1[2];
}

// This functions converts euler angles to quaternion
// [ yaw, pitch, roll ] = [ top, left, front ] = [ lambda, phi, gamma ]
// mapping to bloch sphere representation [φ, ψ, θ]
// Returns the unit quaternion for the given Euler rotation angles [λ, φ, γ]
function versor(e) {
  var l = e[0] / 2 * radians, sl = sin(l), cl = cos(l), // λ / 2
      p = e[1] / 2 * radians, sp = sin(p), cp = cos(p), // φ / 2
      g = e[2] / 2 * radians, sg = sin(g), cg = cos(g); // γ / 2
  return [
    cl * cp * cg + sl * sp * sg,
    sl * cp * cg - cl * sp * sg,
    cl * sp * cg + sl * cp * sg,
    cl * cp * sg - sl * sp * cg
  ];
}

// This function converts a [lon, lat] coordinates into a [x, y, z] coordinate
// the [x, y, z] is cartesian, with origin at [lon, lat] = [0,0] center of the sphere
// Returns cartesian coordinates [x, y, z] given spherical coordinates [λ, φ]
versor.cartesian = function(e) {
  var l = e[0] * radians,
      p = e[1] * radians;
  return [cos(p) * cos(l), cos(p) * sin(l), sin(p)];
};

// This function computes quaternion to euler angles
// Returns the euler rotation angles [λ, φ, γ] for the given quaternion
versor.rotation = function(q) {
  return [
    atan2(2 * (q[0] * q[1] + q[2] * q[3]), 1 - 2 * (q[1] * q[1] + q[2] * q[2])) * degrees,
    asin(max(-1, min(1, 2 * (q[0] * q[2] - q[3] * q[1])))) * degrees,
    atan2(2 * (q[0] * q[3] + q[1] * q[2]), 1 - 2 * (q[2] * q[2] + q[3] * q[3])) * degrees
  ];
};

// This function computes a quaternion representation for the rotation between to vectors
// Returns the quaternion to rotate between two cartesian points on the sphere
versor.delta = function(v0, v1) {
  var w = cross(v0, v1),
      l = sqrt(dot(w, w));
  if (!l)
      return [1, 0, 0, 0];
  var t = acos(max(-1, min(1, dot(v0, v1)))) / 2, // t = θ / 2
      s = sin(t);
  return [cos(t), w[2] / l * s, -w[1] / l * s, w[0] / l * s];
};

// This functions computes a quaternion multiply
// Returns the quaternion that represents q0 * q1
versor.multiply = function(q0, q1) {
  return [
    q0[0] * q1[0] - q0[1] * q1[1] - q0[2] * q1[2] - q0[3] * q1[3],
    q0[0] * q1[1] + q0[1] * q1[0] + q0[2] * q1[3] - q0[3] * q1[2],
    q0[0] * q1[2] - q0[1] * q1[3] + q0[2] * q1[0] + q0[3] * q1[1],
    q0[0] * q1[3] + q0[1] * q1[2] - q0[2] * q1[1] + q0[3] * q1[0]
  ];
};

// This function computes the euler angles when given two vectors and a rotation
// e0 & e1 - starting and ending coordinates in [lon, lat]
// o0 - the projection rotation in euler angles at starting coordinate e0
versor.cartetion = function(e0, e1, o0) {
  // first calculate the quaternion rotation between the two vectors e0 & e1
  // then multiply this rotation onto the original rotation at e0
  // finally convert the resulted quaternion angle back to euler angles for rotation
  var q = versor.multiply(versor(o0), versor.delta(versor.cartesian(e0), versor.cartesian(e1)));
  return versor.rotation(q);
};

// This function computes the metric distance when given two vectors e0 & e1
// Returns the distance using r0 as elevation from center of sphere
versor.distance = function(e0, e1, r0) {
  var l0 = e0[0] * radians, p0 = e0[1] * radians,
      l1 = e1[0] * radians, p1 = e1[1] * radians;
  // calculate using haversine formula
  const a = pow(sin((p1-p0)/2),2) + cos(p1) * cos(p0) * pow(sin((l1-l0)/2),2);
  const c = 2 * atan2(sqrt(a), sqrt(1-a)); //c = 2 * asin(min(1, sqrt(a)));
  return c * r0;
};

return versor;
})));