aboutsummaryrefslogtreecommitdiff
path: root/assets/js/d3-geo-voronoi.js
diff options
context:
space:
mode:
Diffstat (limited to 'assets/js/d3-geo-voronoi.js')
-rw-r--r--assets/js/d3-geo-voronoi.js1753
1 files changed, 1753 insertions, 0 deletions
diff --git a/assets/js/d3-geo-voronoi.js b/assets/js/d3-geo-voronoi.js
new file mode 100644
index 0000000..3d0711c
--- /dev/null
+++ b/assets/js/d3-geo-voronoi.js
@@ -0,0 +1,1753 @@
+// https://github.com/Fil/d3-geo-voronoi Version 0.0.5. Copyright 2017 Philippe Riviere.
+(function (global, factory) {
+ typeof exports === 'object' && typeof module !== 'undefined' ? factory(exports, require('d3-array'), require('d3-collection'), require('d3-geo'), require('d3-voronoi')) :
+ typeof define === 'function' && define.amd ? define(['exports', 'd3-array', 'd3-collection', 'd3-geo', 'd3-voronoi'], factory) :
+ (factory((global.d3 = global.d3 || {}),global.d3,global.d3,global.d3,global.d3));
+}(this, (function (exports,d3Array,d3Collection,d3Geo,d3Voronoi) { 'use strict';
+
+//
+// Copyright (c) 2016, by Loren Petrich
+//
+// Distributed under the terms of the MIT License
+//
+// http://lpetrich.org/
+//
+
+// Calculates the spherical Delaunay triangulation of a set of points
+// These points are entered as an array of arrays of coordinates: 0, 1, 2
+// Any extra members are ignored
+
+// FindDelaunayTriangulation(Positions) and
+// FindDelaunayTriangulationIndexed(Positions, Indices)
+// work from an array of points as specified above,
+// the second one working from a set of indices into the array,
+// and return an object with these members:
+
+// positions -- vectors on a unit sphere
+// indices -- of all the vertices
+// triangles -- array of TriangleObject
+// edges -- array of EdgeObject
+// hull -- array of vertex indices -- the convex hull
+// vor_positions -- positions of triangles' circumcircle centers (Voronoi vertices)
+// vor_edges -- pair of indices in vor_positions (empty one: [-1,-1])
+// vor_polygons -- object indexed by vertex index,
+// and containing edges (EdgeObject), triangles (TriangleObject),
+// and boundary (vertices in vor_positions)
+// Open ones have a -1 at each end.
+
+// Warning: ImproveTriangulation() is mysteriously buggy
+// and is effectively disabled for now
+
+function dotprd(x,y)
+{
+ var sum = 0.0;
+ for (var ic=0; ic<3; ic++)
+ sum += x[ic]*y[ic];
+ return sum;
+}
+
+function crsprd(x,y)
+{
+ var prod = new Array(3);
+ for (var ic=0; ic<3; ic++)
+ {
+ var ic1 = ic + 1;
+ if (ic1 >= 3) ic1 -= 3;
+ var ic2 = ic + 2;
+ if (ic2 >= 3) ic2 -= 3;
+ prod[ic] = x[ic1]*y[ic2] - x[ic2]*y[ic1];
+ }
+ return prod;
+}
+
+function triple_prd(x,y,z)
+{
+ return dotprd(crsprd(x,y),z);
+}
+
+// This distance formula has some nice properties:
+// distance and not square of distance;
+// the square roots give better numerical resolution
+// distance of antipode(p) to p' = - (distance of p to p')
+// Range: -2 to +2
+function ptdist(x,y)
+{
+ var dst1 = 0.0;
+ var dst2 = 0.0;
+ for (var ic=0; ic<3; ic++)
+ {
+ var diff1 = y[ic] - x[ic];
+ dst1 += diff1*diff1;
+ var diff2 = y[ic] + x[ic];
+ dst2 += diff2*diff2;
+ }
+ return Math.sqrt(dst1) - Math.sqrt(dst2);
+}
+
+function Normalize(vec)
+{
+ var vecres = new Array(3);
+ var sum = 0.0, nrmult;
+ for (var ic=0; ic<3; ic++)
+ {
+ var val = vec[ic];
+ sum += val*val;
+ }
+ if (sum > 0)
+ nrmult = 1/Math.sqrt(sum);
+ else
+ nrmult = 0;
+ for (var ic=0; ic<3; ic++)
+ {
+ vecres[ic] = nrmult*vec[ic];
+ }
+ return vecres;
+}
+
+function crsprd_ix(Positions,ix,iy)
+{
+ return crsprd(Positions[ix],Positions[iy]);
+}
+
+function triple_prd_ix(Positions,ix,iy,iz)
+{
+ return triple_prd(Positions[ix],Positions[iy],Positions[iz]);
+}
+
+function ptdist_ix(Positions,ix,iy)
+{
+ return ptdist(Positions[ix],Positions[iy]);
+}
+
+// Returns a zero 3-vector
+function zerovec()
+{
+ var vec = new Array(3);
+ for (var ic=0; ic<3; ic++)
+ vec[ic] = 0.0;
+ return vec;
+}
+
+// Implements copying
+function vec_copy(x)
+{
+ var vec = new Array(3);
+ for (var ic=0; ic<3; ic++)
+ vec[ic] = x[ic];
+ return vec;
+}
+
+// Implements x += y
+function vec_add_to(x,y)
+{
+ for (var ic=0; ic<3; ic++)
+ x[ic] += y[ic];
+}
+
+// Implements x *= y
+function vec_mult_scalar_to(x,y)
+{
+ for (var ic=0; ic<3; ic++)
+ x[ic] *= y;
+}
+
+// Implements x - y
+function vec_difference(x,y)
+{
+ var diff = zerovec();
+ for (var ic=0; ic<3; ic++)
+ diff[ic] = x[ic] - y[ic];
+ return diff;
+}
+
+// JavaScript's counterpart of "null" / "None":
+function IsNull(x)
+{
+ return (typeof(x) == 'undefined');
+}
+
+function TrianglesEqual(tr1, tr2)
+{
+ if (IsNull(tr1)) return false;
+ if (IsNull(tr2)) return false;
+
+ for (var iv=0; iv<3; iv++)
+ if (tr1.verts[iv] != tr2.verts[iv])
+ return false;
+
+ return true;
+}
+
+function EdgesEqual(ed1, ed2)
+{
+ if (IsNull(ed1)) return false;
+ if (IsNull(ed2)) return false;
+
+ for (var iv=0; iv<2; iv++)
+ if (ed1.verts[iv] != ed2.verts[iv])
+ return false;
+
+ return true;
+}
+
+function max(x,y)
+{
+ return (y > x) ? y : x;
+}
+
+function TriangleObject(Positions, verts)
+{
+ this.verts = verts;
+ this.edges = new Array(3);
+
+ // Find directions for testing whether a point is inside
+ this.dirs = new Array(3);
+ for (var ic=0; ic<3; ic++)
+ {
+ var ic1 = ic + 1;
+ if (ic1 >= 3) ic1 -= 3;
+ var ic2 = ic + 2;
+ if (ic2 >= 3) ic2 -= 3;
+ this.dirs[ic] = crsprd_ix(Positions,verts[ic1],verts[ic2]);
+ }
+
+ // Tetrahedral volume factor
+ this.vol = triple_prd_ix(Positions,verts[0],verts[1],verts[2]);
+
+ // Adjust to get the signs correct for the point-inside test;
+ // the vertex opposite the edges' vertices ought to give a dot product of 1
+ for (var ic=0; ic<3; ic++)
+ vec_mult_scalar_to(this.dirs[ic],1/this.vol);
+
+ // Circumcircle test
+ var ccdir = zerovec();
+ for (var ic=0; ic<3; ic++)
+ vec_add_to(ccdir,this.dirs[ic]);
+ this.ccdir = Normalize(ccdir);
+
+ var ccdsq = 0;
+ for (var ic=0; ic<3; ic++)
+ ccdsq += ptdist(this.ccdir,Positions[verts[ic]]);
+ ccdsq /= 3;
+ this.ccdsq = ccdsq;
+}
+
+// For copying in vertex info from another triangle
+TriangleObject.prototype.copy_vert_info = function(src)
+{
+ this.verts = src.verts;
+ this.dirs = src.dirs;
+ this.vol = src.vol;
+ this.ccdir = src.ccdir;
+ this.ccdsq = src.ccdsq;
+};
+
+TriangleObject.prototype.IsVertexOrderCorrect = function()
+{
+ return this.vol >= 0;
+};
+
+TriangleObject.prototype.IsPointInside = function(p)
+{
+ for (var ic=0; ic<3; ic++)
+ if (dotprd(p,this.dirs[ic]) < 0) return false;
+
+ return true;
+};
+
+TriangleObject.prototype.IsPointInCircumcircle = function(p)
+{
+ return (ptdist(this.ccdir,p) < this.ccdsq);
+};
+
+TriangleObject.prototype.IsVertex = function(ix)
+{
+ for (var ic=0; ic<3; ic++)
+ if (ix == this.verts[ic]) return true;
+ return false;
+};
+
+TriangleObject.prototype.VertexIndexIn = function(ix)
+{
+ for (var ic=0; ic<3; ic++)
+ if (ix == this.verts[ic]) return ic;
+ return -1;
+};
+
+TriangleObject.prototype.EdgeIndexIn = function(ed)
+{
+ for (var ic=0; ic<3; ic++)
+ if (EdgesEqual(this.edges[ic], ed)) return ic;
+ return -1;
+};
+
+function EdgeObject(verts)
+{
+ this.verts = verts;
+ this.polys = new Array(2);
+}
+
+EdgeObject.prototype.IsVertex = function(ix)
+{
+ for (var ic=0; ic<2; ic++)
+ if (ix == this.verts[ic]) return true;
+ return false;
+};
+
+EdgeObject.prototype.VertexIndexIn = function(ix)
+{
+ for (var ic=0; ic<2; ic++)
+ if (ix == this.verts[ic]) return ic;
+ return -1;
+};
+
+EdgeObject.prototype.PolyIndexIn = function(pl)
+{
+ for (var ic=0; ic<2; ic++)
+ if (TrianglesEqual(this.polys[ic],pl)) return ic;
+ return -1;
+};
+
+function EdgeCheckObject(Positions,verts)
+{
+ this.verts = verts;
+ this.pdst = ptdist_ix(Positions,verts[0],verts[1]);
+ this.direc = Normalize(crsprd_ix(Positions,verts[0],verts[1]));
+ var midpnt = zerovec();
+ vec_add_to(midpnt,Positions[verts[0]]);
+ vec_add_to(midpnt,Positions[verts[1]]);
+ this.midpnt = Normalize(midpnt);
+}
+
+// Check on the possible intersection with another edge-check object
+// return a boolean of whether or not it does
+EdgeCheckObject.prototype.intersects = function(Positions,other)
+{
+ // Assume that sharing a vertex means non-intersecting
+ for (var ic=0; ic<2; ic++)
+ for (var ict=0; ict<2; ict++)
+ if (this.verts[ic] == other.verts[ict]) return false;
+
+ // Find intersection point; will test it and its antipode
+ var itsc = Normalize(crsprd(this.direc, other.direc));
+
+ // Find dot product with midpoints to test if the intersection
+ // is in the near hemispheres of the lines' midpoints.
+ // If it is in both near hemispheres or both far hemispheres, it's OK
+ // In both far hemispheres: antipode is in both near hemispheres
+ var near0 = dotprd(itsc,this.midpnt) > 0;
+ var near1 = dotprd(itsc,other.midpnt) > 0;
+ if (near1 != near0) return false;
+
+ var pd0 = [];
+ for (var ic=0; ic<2; ic++)
+ {
+ var pd = ptdist(itsc, Positions[this.verts[ic]]);
+ pd0.push(pd);
+ }
+ var pd1 = [];
+ for (var ic=0; ic<2; ic++)
+ {
+ var pd = ptdist(itsc, Positions[other.verts[ic]]);
+ pd1.push(pd);
+ }
+ var mxpd0 = max(pd0[0],pd0[1]);
+ var mxpd1 = max(pd1[0],pd1[1]);
+ if ((mxpd0 <= this.pdst) && (mxpd1 <= other.pdst) && near0) return true;
+
+ // Handle its antipode; use antipode-related shortcuts
+ // like reversing the distance value and the hemisphere-presence value
+ vec_mult_scalar_to(itsc, -1);
+ near0 = !near0;
+ for (var ic=0; ic<2; ic++)
+ {
+ pd0[ic] = - pd0[ic];
+ pd1[ic] = - pd1[ic];
+ }
+ mxpd0 = max(pd0[0],pd0[1]);
+ mxpd1 = max(pd1[0],pd1[1]);
+ if ((mxpd0 <= this.pdst) && (mxpd1 <= other.pdst) && near0) return true;
+
+ return false;
+};
+
+// Adds to an array if it was not already present;
+// Must resort to this kludge because JavaScript doesn't have sets
+function AddUnique(arr, x)
+{
+ for (var i=0; i<arr.length; i++)
+ if (x == arr[i]) return;
+ arr.push(x);
+}
+
+// Version for edges, since testing equality of objects
+// doesn't work that well in JavaScript
+function AddUniqueEdge(arr, ed)
+{
+ for (var i=0; i<arr.length; i++)
+ if (EdgesEqual(arr[i],ed)) return;
+ arr.push(ed);
+}
+
+// Find the set intersection
+function FindShared(arr1, arr2)
+{
+ var resarr = [];
+ for (var i1=0; i1<arr1.length; i1++)
+ {
+ var x1 = arr1[i1];
+ for (var i2=0; i2<arr2.length; i2++)
+ {
+ var x2 = arr2[i2];
+ if (x1 == x2)
+ {
+ resarr.push(x1);
+ break;
+ }
+ }
+ }
+
+ return resarr;
+}
+
+// Version for edges
+function FindSharedEdges(arr1, arr2)
+{
+ var resarr = [];
+ for (var i1=0; i1<arr1.length; i1++)
+ {
+ var x1 = arr1[i1];
+ for (var i2=0; i2<arr2.length; i2++)
+ {
+ var x2 = arr2[i2];
+ if (EdgesEqual(x1,x2))
+ {
+ resarr.push(x1);
+ break;
+ }
+ }
+ }
+
+ return resarr;
+}
+
+// Takes all the members of of arr2 out of arr1
+// and ignores the arr2 members not present in arr1
+function FindSetDifference(arr1, arr2)
+{
+ if (arr2.length == 0) return;
+ var diffarr = [];
+ for (var i1=0; i1<arr1.length; i1++)
+ {
+ var x1 = arr1[i1];
+ var AddThisOne = true;
+ for (var i2=0; i2<arr2.length; i2++)
+ {
+ var x2 = arr2[i2];
+ if (x2 == x1)
+ {
+ AddThisOne = false;
+ break;
+ }
+ }
+ if (AddThisOne) diffarr.push(x1);
+ }
+
+ // Clear the array
+ arr1.splice(0,arr1.length);
+
+ for (var i=0; i<diffarr.length; i++)
+ arr1.push(diffarr[i]);
+}
+
+// Version for edges
+function FindSetDifferenceEdges(arr1, arr2)
+{
+ if (arr2.length == 0) return;
+ var diffarr = [];
+ for (var i1=0; i1<arr1.length; i1++)
+ {
+ var x1 = arr1[i1];
+ var AddThisOne = true;
+ for (var i2=0; i2<arr2.length; i2++)
+ {
+ var x2 = arr2[i2];
+ if (EdgesEqual(x1,x2))
+ {
+ AddThisOne = false;
+ break;
+ }
+ }
+ if (AddThisOne) diffarr.push(x1);
+ }
+
+ // Clear the array
+ arr1.splice(0,arr1.length);
+
+ for (var i=0; i<diffarr.length; i++)
+ arr1.push(diffarr[i]);
+}
+
+
+// Specified by index ix; returns whether it was possible to do so
+function AddPointInside(TriSet, ix)
+{
+ var Positions = TriSet.positions;
+ var p = Positions[ix];
+
+ var NumTris = TriSet.triangles.length;
+ for (var j=0; j<NumTris; j++)
+ {
+ var tri = TriSet.triangles[j];
+ if (tri.IsPointInside(p))
+ {
+ // Create three new triangles and their edges
+ var eds = tri.edges;
+ var trixs = [];
+ for (var ic=0; ic<3; ic++)
+ trixs.push(eds[ic].PolyIndexIn(tri));
+
+ var newtris = Array(3);
+ var neweds = Array(3);
+ for (var ic=0; ic<3; ic++)
+ {
+ var ic1 = ic + 1;
+ if (ic1 >= 3) ic1 -= 3;
+ newtris[ic] = new TriangleObject(Positions,[tri.verts[ic],tri.verts[ic1],ix]);
+ neweds[ic] = new EdgeObject([tri.verts[ic],ix]);
+ }
+
+ // Connect those triangles and edges
+ for (var ic=0; ic<3; ic++)
+ {
+ var ic1 = ic + 1;
+ if (ic1 >= 3) ic1 -= 3;
+ newtris[ic].edges[0] = neweds[ic1];
+ newtris[ic].edges[1] = neweds[ic];
+ neweds[ic].polys[0] = newtris[ic];
+ neweds[ic1].polys[1] = newtris[ic];
+ }
+
+ // Find which external edges go with which triangles
+ for (var ic=0; ic<3; ic++)
+ {
+ var ed = eds[ic];
+ var trix = trixs[ic];
+ for (var ict=0; ict<3; ict++)
+ {
+ var newtri = newtris[ict];
+ var numverts = 0;
+ for (var iv=0; iv<2; iv++)
+ {
+ if (newtri.IsVertex(ed.verts[iv]))
+ numverts++;
+ if (numverts == 2)
+ {
+ ed.polys[trix] = newtri;
+ newtri.edges[2] = ed;
+ break;
+ }
+ }
+ }
+ }
+
+ // Insert those triangles and edges into the lists
+ TriSet.triangles[j] = newtris[0];
+ for (var ic=1; ic<3; ic++)
+ TriSet.triangles.push(newtris[ic]);
+ for (var ic=0; ic<3; ic++)
+ TriSet.edges.push(neweds[ic]);
+
+ // All done; indicate that the point was added
+ return true;
+ }
+ }
+
+ // The point was inside no triangle, and thus was not added
+ return false;
+}
+
+function ImproveTriangulation(TriSet)
+{
+ var Positions = TriSet.positions;
+ var quad_verts = new Array(4);
+ for (var itr=0; itr<100; itr++)
+ {
+ var numflips = 0;
+ for (var i=0; i<TriSet.edges.length; i++)
+ {
+ var edge = TriSet.edges[i];
+ var tris = edge.polys;
+
+ // Skip over external edges
+ if (IsNull(tris[0])) continue;
+ if (IsNull(tris[1])) continue;
+
+ // Find the containing quadrangle's vertices
+ for (var ic=0; ic<3; ic++)
+ {
+ var ix = tris[0].verts[ic];
+ if (!edge.IsVertex(ix)) break;
+ }
+ var ic1 = ic + 1;
+ if (ic1 >= 3) ic1 -= 3;
+ var ic2 = ic + 2;
+ if (ic2 >= 3) ic2 -= 3;
+ quad_verts[0] = ix;
+ quad_verts[1] = tris[0].verts[ic1];
+ quad_verts[3] = tris[0].verts[ic2];
+
+ for (var ic=0; ic<3; ic++)
+ {
+ var ix = tris[1].verts[ic];
+ if (!edge.IsVertex(ix)) break;
+ }
+ quad_verts[2] = ix;
+
+ // Are the non-edge points in the other triangles' circumcircles?
+ var incc0 = tris[0].IsPointInCircumcircle(Positions[quad_verts[2]]);
+ var incc1 = tris[1].IsPointInCircumcircle(Positions[quad_verts[0]]);
+ if ((!incc0) && (!incc1)) continue;
+
+ // Are the would-be triangles properly oriented?
+ var newtri0 = new TriangleObject(Positions, [quad_verts[0],quad_verts[1],quad_verts[2]]);
+ if (!newtri0.IsVertexOrderCorrect()) continue;
+ var newtri1 = new TriangleObject(Positions, [quad_verts[0],quad_verts[2],quad_verts[3]]);
+ if (!newtri1.IsVertexOrderCorrect()) continue;
+
+ // If so, then flip
+ numflips++;
+
+ // Adjust the edge and triangle memberships:
+ // 0-3 goes from 0 to 1, 1-2 goes from 1 to 0
+ for (var ic=0; ic<3; ic++)
+ {
+ var ed = tris[0].edges[ic];
+ if (EdgesEqual(ed,edge)) continue;
+ else if (ed.IsVertex(quad_verts[3]))
+ {
+ var ed03 = ed;
+ var edix03 = ic;
+ break;
+ }
+ }
+ for (var ic=0; ic<3; ic++)
+ {
+ var ed = tris[1].edges[ic];
+ if (EdgesEqual(ed,edge)) continue;
+ else if (ed.IsVertex(quad_verts[1]))
+ {
+ var ed12 = ed;
+ var edix12 = ic;
+ break;
+ }
+ }
+
+ var trix0 = ed03.PolyIndexIn(tris[0]);
+ var trix1 = ed12.PolyIndexIn(tris[1]);
+
+ ed03.polys[trix0] = tris[1];
+ ed12.polys[trix1] = tris[0];
+ tris[0].edges[edix03] = ed12;
+ tris[1].edges[edix12] = ed03;
+
+ // Add the vertices
+ tris[0].copy_vert_info(newtri0);
+ tris[1].copy_vert_info(newtri1);
+ edge.verts = [quad_verts[0], quad_verts[2]];
+ }
+ if (numflips == 0) break;
+ }
+}
+
+
+function FindConvexHull(TriSet)
+{
+ // var Positions = TriSet.positions;
+
+ // Find boundary loop -- use as convex hull
+ var NextVertex = new Object;
+ var VtxStart = -1;
+ for (var i=0; i<TriSet.edges.length; i++)
+ {
+ var edge = TriSet.edges[i];
+
+ // Find a boundary one -- look for the triangle that it contains
+ if (IsNull(edge.polys[0]))
+ {
+ if (IsNull(edge.polys[1]))
+ continue;
+ else
+ var tri = edge.polys[1];
+ }
+ else
+ {
+ if (IsNull(edge.polys[1]))
+ var tri = edge.polys[0];
+ else
+ continue;
+ }
+
+ // Ensure that the hull is in the same direction as the triangles
+ var ix0 = edge.verts[0];
+ var ix1 = edge.verts[1];
+ var posdiff = tri.VertexIndexIn(ix1) - tri.VertexIndexIn(ix0);
+ if (posdiff < 0) posdiff += 3;
+ if (posdiff != 1)
+ {
+ var ixs = ix0;
+ ix0 = ix1;
+ ix1 = ixs;
+ }
+
+ NextVertex[ix0] = ix1;
+ VtxStart = ix0;
+ }
+
+ if (VtxStart >= 0)
+ {
+ var ix = VtxStart;
+ var hull = [ix];
+ while(true)
+ {
+ var ixnext = NextVertex[ix];
+ if (ixnext == VtxStart) break;
+ hull.push(ixnext);
+ ix = ixnext;
+ }
+
+ TriSet.hull = hull;
+ }
+}
+
+
+// Finds the dual of the Delaunay triangulation
+// Won't bother to do the sort of connectivity
+// that was necessary for the Delaunay triangulation
+function FindVoronoiDiagram(TriSet)
+{
+ // Special cases: 3 or fewer points
+ if (TriSet.triangles.length == 1)
+ {
+ // A single triangle
+ if (TriSet.hull.length == 3)
+ {
+ var tri = TriSet.triangles[0];
+ TriSet.vor_positions.push(tri.ccdir);
+ for (var k=0; k<3; k++)
+ {
+ var kx = k + 1;
+ if (kx >= 3) kx = 0;
+ var ky = k - 1;
+ if (ky < 0) ky = 2;
+
+ var v1 = TriSet.positions[TriSet.hull[k]];
+ var v2 = TriSet.positions[TriSet.hull[kx]];
+ var posdiff = vec_difference(v2,v1);
+ TriSet.vor_positions.push(Normalize(crsprd(posdiff,tri.ccdir)));
+ TriSet.vor_edges.push([0, k+1, 4]);
+
+ var ix = TriSet.hull[k];
+ TriSet.vor_polygons[ix] = new Object;
+ var vor_poly = TriSet.vor_polygons[ix];
+ var iy = TriSet.hull[ky];
+ for (var l=0; l<3; l++)
+ {
+ var edge = TriSet.edges[l];
+ var shrd = FindShared([iy,ix],edge.verts);
+ if (shrd.length == 2) break;
+ }
+
+ vor_poly.edges = [edge];
+ vor_poly.triangles = [tri];
+ vor_poly.boundary = [0, ky+1, 4, k+1];
+ }
+ var ept = vec_copy(tri.ccdir);
+ vec_mult_scalar_to(ept,-1);
+ TriSet.vor_positions.push(ept);
+ }
+ return;
+ }
+ else if (TriSet.triangles.length == 0)
+ {
+ // A biangle
+ if (TriSet.hull.length == 2)
+ {
+ var v0 = TriSet.positions[TriSet.hull[0]];
+ var v1 = TriSet.positions[TriSet.hull[1]];
+
+ var vt0 = zerovec();
+ vec_add_to(vt0,v0);
+ vec_add_to(vt0,v1);
+ vt0 = Normalize(vt0);
+ TriSet.vor_positions.push(vt0);
+
+ var vt1 = Normalize(crsprd(v0,v1));
+ TriSet.vor_positions.push(vt1);
+
+ var vt2 = vec_copy(vt0);
+ vec_mult_scalar_to(vt2,-1);
+ TriSet.vor_positions.push(vt2);
+
+ var vt3 = vec_copy(vt1);
+ vec_mult_scalar_to(vt3,-1);
+ TriSet.vor_positions.push(vt3);
+
+ TriSet.vor_edges.push([0, 1, 2, 3, 0]);
+
+ edge = TriSet.edges[0];
+ for (var k=0; k<2; k++)
+ {
+ var ix = TriSet.hull[k];
+ TriSet.vor_polygons[ix] = new Object;
+ var vor_poly = TriSet.vor_polygons[ix];
+ vor_poly.edges = [edge];
+ vor_poly.triangles = [0];
+ if (k == 0)
+ vor_poly.boundary = [0, 1, 2, 3];
+ else if (k == 1)
+ vor_poly.boundary = [0, 3, 2, 1];
+ }
+ }
+ return;
+ }
+
+ // Create the array of Voronoi-vertex positions:
+ // Add indices to the triangle objects for convenience
+ for (var i=0; i<TriSet.triangles.length; i++)
+ {
+ var tri = TriSet.triangles[i];
+ tri.index = i;
+ TriSet.vor_positions.push(tri.ccdir);
+ }
+
+ // Voronoi edges: a cinch
+ // Voronoi edges parallel original edges
+ for (var i=0; i<TriSet.edges.length; i++)
+ {
+ var edge = TriSet.edges[i];
+ if (!IsNull(edge.polys[0]) && !IsNull(edge.polys[1]))
+ {
+ var vor_edge = [edge.polys[0].index, edge.polys[1].index];
+ TriSet.vor_edges.push(vor_edge);
+ }
+ }
+
+ // Voronoi polygons: -1 at ends means an open one
+
+ // First, collect the edges and triangles at each vertex
+ // Put them into vor_polygons, because each one
+ // is for each original vertex
+ for (var i=0; i<TriSet.indices.length; i++)
+ {
+ var ix = TriSet.indices[i];
+ TriSet.vor_polygons[ix] = new Object;
+
+ var vor_poly = TriSet.vor_polygons[ix];
+ vor_poly.edges = [];
+ vor_poly.triangles = [];
+ vor_poly.boundary = [];
+ }
+
+ for (var i=0; i<TriSet.edges.length; i++)
+ {
+ var edge = TriSet.edges[i];
+ for (var ic=0; ic<2; ic++)
+ TriSet.vor_polygons[edge.verts[ic]].edges.push(edge);
+ }
+
+ for (var i=0; i<TriSet.triangles.length; i++)
+ {
+ var tri = TriSet.triangles[i];
+ for (var ic=0; ic<3; ic++)
+ TriSet.vor_polygons[tri.verts[ic]].triangles.push(tri);
+ }
+
+ for (var i=0; i<TriSet.indices.length; i++)
+ {
+ var ix = TriSet.indices[i];
+ var vor_poly = TriSet.vor_polygons[ix];
+
+ // First triangle
+ var init_tri = vor_poly.triangles[0];
+ var tri = init_tri;
+ vor_poly.boundary.push(tri.index);
+
+ // First edge
+ for (var ic=0; ic<3; ic++)
+ {
+ var edge = tri.edges[ic];
+ if (edge.IsVertex(ix))
+ break;
+ }
+ var init_edge = edge;
+
+ // The next triangle and edge
+ var IsInside = false;
+ while(true)
+ {
+ var iv = edge.PolyIndexIn(tri);
+ tri = edge.polys[1-iv];
+ if (IsNull(tri)) break;
+ if (TrianglesEqual(tri,init_tri))
+ {
+ IsInside = true;
+ break;
+ }
+
+ vor_poly.boundary.push(tri.index);
+
+ for (var ic=0; ic<3; ic++)
+ {
+ var next_edge = tri.edges[ic];
+ if (EdgesEqual(next_edge,edge)) continue;
+ if (next_edge.IsVertex(ix))
+ {
+ edge = next_edge;
+ break;
+ }
+ }
+ }
+
+ if (!IsInside)
+ {
+ vor_poly.boundary.reverse();
+ tri = init_tri;
+
+ // First edge the other way
+ for (var ic=0; ic<3; ic++)
+ {
+ edge = tri.edges[ic];
+ if (EdgesEqual(edge,init_edge)) continue;
+ if (edge.IsVertex(ix))
+ break;
+ }
+
+ while(true)
+ {
+ var iv = edge.PolyIndexIn(tri);
+ tri = edge.polys[1-iv];
+ if (IsNull(tri)) break;
+
+ vor_poly.boundary.push(tri.index);
+
+ for (var ic=0; ic<3; ic++)
+ {
+ var next_edge = tri.edges[ic];
+ if (EdgesEqual(next_edge,edge)) continue;
+ if (next_edge.IsVertex(ix))
+ {
+ edge = next_edge;
+ break;
+ }
+ }
+ }
+ }
+
+ // Add -1 on ends for open polygon:
+ if (!IsInside)
+ {
+ vor_poly.boundary.reverse();
+ vor_poly.boundary.push(-1);
+ vor_poly.boundary.reverse();
+ vor_poly.boundary.push(-1);
+ }
+ }
+
+ // Handle the area outside of the convex hull
+ if (TriSet.hull.length >= 3)
+ {
+ // Set up the initial boundary lines
+ // The boundary lines contain:
+ // Index of Voronoi vertex / triangle center / intersection (in VorPos)
+ // Indices of original vertices on each side of the line
+ var VorBdLns = new Array();
+ var Positions = TriSet.positions;
+ var hlen = TriSet.hull.length;
+ for (var ic=0; ic<hlen; ic++)
+ {
+ var ix = TriSet.hull[ic];
+ var icx = ic + 1;
+ if (icx >= hlen) icx = 0;
+ var ixa = TriSet.hull[icx];
+ var edset1 = TriSet.vor_polygons[ix].edges;
+ var edset2 = TriSet.vor_polygons[ixa].edges;
+ var edsetshr = FindSharedEdges(edset1,edset2);
+ var edge = edsetshr[0];
+ var tvrt = edge.polys[0].index;
+ var vt0 = Positions[ix];
+ var vt1 = Positions[ixa];
+ var vtdf = vec_difference(vt1,vt0);
+ // Contains: triangle index (Voronoi vertex),
+ // vertex index 1 (Voronoi region), position
+ // vertex index 2 (Voronoi region), position,
+ // great-circle normal
+ var VorBdLn = [tvrt, TriSet.vor_positions[tvrt], ix, vt0, ixa, vt1, vtdf];
+ VorBdLns.push(VorBdLn);
+ }
+
+ // Find intersections
+
+ while (VorBdLns.length > 3)
+ {
+ // Check all combinations of neighbors
+ var n = VorBdLns.length;
+ var itscpts = new Array();
+ var ptitscs = new Array();
+ for (var k=0; k<n; k++)
+ ptitscs.push(new Array());
+ for (var k=0; k<n; k++)
+ {
+ // Find the intersection point; use the convex hull's direction
+ var kx = k + 1;
+ if (kx >= n) kx = 0;
+ var itscpt = Normalize(crsprd(VorBdLns[k][6],VorBdLns[kx][6]));
+ vec_mult_scalar_to(itscpt,-1);
+ ptitscs[k].push(itscpts.length);
+ ptitscs[kx].push(itscpts.length);
+ itscpts.push(itscpt);
+ }
+ // Find the intersection points that are the closest to their parent points
+ for (var k=0; k<n; k++)
+ {
+ var ptitsc = ptitscs[k];
+ if (ptitsc.length >= 2)
+ {
+ var dists = new Array();
+ for (var kp=0; kp<ptitsc.length; kp++)
+ dists.push(ptdist(itscpts[ptitsc[kp]],VorBdLns[k][1]));
+ var dx = 0;
+ var dmin = dists[dx];
+ for (var dxi=0; dxi<dists.length; dxi++)
+ {
+ var dst = dists[dxi];
+ if (dst < dmin)
+ {
+ dx = dxi; dmin = dst;
+ }
+ }
+ var ptitscrd = ptitsc[dx];
+ }
+ else if (ptitsc.length == 1)
+ var ptitscrd = ptitsc[0];
+ else
+ var ptitscrd = -1;
+ ptitscs[k] = ptitscrd;
+ }
+
+ var NewVorBdLns = new Array();
+ for (var k=0; k<n; k++)
+ {
+ // Find all matched intersection points and add them
+ var kx = k + 1;
+ if (kx >= n) kx = 0;
+ var ky = k - 1;
+ if (ky < 0) ky = n - 1;
+ // 0 is lone, 1 is leading, 2 is trailing
+ // vorvtidx is the index of the Voronoi vertex
+ var pstat = 0;
+ var ptitsc = ptitscs[k], ptitsc_next;
+ if (ptitsc != -1)
+ {
+ var ptitsc_prev = ptitscs[ky];
+ if (ptitsc == ptitsc_prev)
+ pstat = 2;
+ else
+ {
+ ptitsc_next = ptitscs[kx];
+ if (ptitsc == ptitsc_next)
+ pstat = 1;
+ }
+ }
+
+ if (pstat == 0)
+ {
+ // Keep the Voronoi line without merging
+ NewVorBdLns.push(VorBdLns[k]);
+ }
+ else if (pstat == 1)
+ {
+ // Merge the Voronoi lines and create a new one
+ var VorBdLn0 = VorBdLns[k];
+ var VorBdLn1 = VorBdLns[kx];
+ var itscpt = itscpts[k];
+ var tvrt0 = VorBdLn0[0];
+ var tvrt1 = VorBdLn1[0];
+ var PointOK = (tvrt1 != tvrt0);
+ if (PointOK)
+ {
+ var nitx = TriSet.vor_positions.length;
+ var ix0 = VorBdLn0[2];
+ var vt0 = VorBdLn0[3];
+ var ix1 = VorBdLn1[4];
+ var vt1 = VorBdLn1[5];
+ var dst_in = undefined;
+ var dst_out = undefined;
+ for (var m=0; m<n; m++)
+ {
+ var ptstm = ptdist(VorBdLns[m][3],itscpt);
+ var mrl = m - k;
+ while (mrl < 0) mrl += n;
+ while (mrl >= n) mrl -= n;
+ if (mrl <= 2)
+ {
+ if (dst_in == undefined)
+ dst_in = ptstm;
+ else if (ptstm < dst_in)
+ dst_in = ptstm;
+ }
+ else
+ {
+ if (dst_out == undefined)
+ dst_out = ptstm;
+ else if (ptstm < dst_out)
+ dst_out = ptstm;
+ }
+ }
+ PointOK = (dst_in < dst_out);
+ }
+ if (PointOK)
+ {
+ var vtdf = vec_difference(vt1,vt0);
+ var VorBdLn = [nitx, itscpt, ix0, vt0, ix1, vt1, vtdf];
+ NewVorBdLns.push(VorBdLn);
+ TriSet.vor_positions.push(itscpt);
+ var ixi = VorBdLn0[4];
+ // Should be equal:
+ // ixi = VorBdLn2[2];
+ TriSet.vor_edges.push([tvrt0, nitx]);
+ TriSet.vor_edges.push([tvrt1, nitx]);
+ // Add the point to the center Voronoi region and close it
+ TriSet.vor_polygons[ixi].boundary.shift();
+ var vpln = TriSet.vor_polygons[ixi].boundary.length;
+ TriSet.vor_polygons[ixi].boundary[vpln-1] = nitx;
+ // Add the point to the left Voronoi region
+ if (TriSet.vor_polygons[ix0].boundary[1] == tvrt0)
+ {
+ TriSet.vor_polygons[ix0].boundary.unshift(-1);
+ TriSet.vor_polygons[ix0].boundary[1] = nitx;
+ }
+ else
+ {
+ vpln = TriSet.vor_polygons[ix0].boundary.length;
+ if (TriSet.vor_polygons[ix0].boundary[vpln-2] == tvrt0)
+ {
+ TriSet.vor_polygons[ix0].boundary.push(-1);
+ vpln = TriSet.vor_polygons[ix0].boundary.length;
+ TriSet.vor_polygons[ix0].boundary[vpln-2] = nitx;
+ }
+ }
+ // Add the point to the right Voronoi region
+ if (TriSet.vor_polygons[ix1].boundary[1] == tvrt1)
+ {
+ TriSet.vor_polygons[ix1].boundary.unshift(-1);
+ TriSet.vor_polygons[ix1].boundary[1] = nitx;
+ }
+ else
+ {
+ vpln = TriSet.vor_polygons[ix1].boundary.length;
+ if (TriSet.vor_polygons[ix1].boundary[vpln-2] == tvrt1)
+ {
+ TriSet.vor_polygons[ix1].boundary.push(-1);
+ vpln = TriSet.vor_polygons[ix1].boundary.length;
+ TriSet.vor_polygons[ix1].boundary[vpln-2] = nitx;
+ }
+ }
+ }
+ else
+ {
+ NewVorBdLns.push(VorBdLn0);
+ NewVorBdLns.push(VorBdLn1);
+ }
+ }
+ /*
+ else if (pstat == 2)
+ {
+ // Do nothing
+ }
+ */
+ }
+
+ if (NewVorBdLns.length == VorBdLns.length) break;
+ VorBdLns = NewVorBdLns;
+ }
+
+ // Special cases: only two or three points left
+ if (VorBdLns.length == 2)
+ {
+ if (VorBdLns[0][0] != VorBdLns[1][0])
+ {
+ var VorLn = [];
+ for (var k=0; k<2; k++)
+ {
+ // Connecting line
+ var kx = VorBdLns[k][0];
+ VorLn.push(kx);
+ // Close the Voronoi region by deleting the end -1's
+ kx = VorBdLns[k][2];
+ TriSet.vor_polygons[kx].boundary.shift();
+ TriSet.vor_polygons[kx].boundary.pop();
+ }
+ TriSet.vor_edges.push(VorLn);
+ }
+ }
+ else if (VorBdLns.length == 3)
+ {
+ var ic0 = VorBdLns[0][0];
+ var ic1 = VorBdLns[1][0];
+ var ic2 = VorBdLns[2][0];
+ if (ic0 != ic1 && ic0 != ic2 && ic1 != ic2)
+ {
+ var nitx = TriSet.vor_positions.length;
+ var v0 = VorBdLns[0][3];
+ var v1 = VorBdLns[1][3];
+ var v2 = VorBdLns[2][3];
+ var itscpt = zerovec();
+ vec_add_to(itscpt,crsprd(v0,v1));
+ vec_add_to(itscpt,crsprd(v1,v2));
+ vec_add_to(itscpt,crsprd(v2,v0));
+ itscpt = Normalize(itscpt);
+ vec_mult_scalar_to(itscpt,-1);
+ TriSet.vor_positions.push(itscpt);
+ for (var k=0; k<3; k++)
+ {
+ // Connecting line
+ var VorBdLn = VorBdLns[k];
+ TriSet.vor_edges.push([VorBdLn[0], nitx]);
+ // Add the point to the Voronoi region and close it
+ var ix = VorBdLn[2];
+ TriSet.vor_polygons[ix].boundary.shift();
+ var vpln = TriSet.vor_polygons[ix].boundary.length;
+ TriSet.vor_polygons[ix].boundary[vpln-1] = nitx;
+ }
+ }
+ }
+ }
+
+ // Adjust the orientations
+ for (var k=0; k<TriSet.vor_polygons.length; k++)
+ {
+ vor_poly = TriSet.vor_polygons[k];
+ if (vor_poly.boundary.length >= 3 && vor_poly.boundary[0] >= 0)
+ {
+ tri = new TriangleObject(TriSet.vor_positions,vor_poly.boundary.slice(0,3));
+ if (!tri.IsVertexOrderCorrect())
+ vor_poly.boundary.reverse();
+ }
+ }
+}
+
+
+function FindDelaunayTriangulationIndexed(Positions, Indices)
+{
+ // Create the triangle-set object
+ var TriSet = new Object;
+ TriSet.positions = Positions;
+ TriSet.indices = Indices;
+ TriSet.triangles = [];
+ TriSet.edges = [];
+ TriSet.hull = [];
+ TriSet.vor_positions = [];
+ TriSet.vor_edges = [];
+ TriSet.vor_polygons = new Object;
+
+ // Create the first triangle, if it is possible to create any
+ if (Indices.length < 3)
+ {
+ if (Indices.length == 2)
+ {
+ TriSet.edges.push(new EdgeObject(Indices));
+ TriSet.hull = Indices;
+ }
+ FindVoronoiDiagram(TriSet);
+ return TriSet;
+ }
+
+ var tri = new TriangleObject(Positions, Indices.slice(0,3));
+ if (!tri.IsVertexOrderCorrect())
+ tri = new TriangleObject(Positions, [Indices[0],Indices[2],Indices[1]]);
+ TriSet.triangles.push(tri);
+
+ var echs = new Array(3);
+
+ for (var ic=0; ic<3; ic++)
+ {
+ var ic1 = ic + 1;
+ if (ic1 >= 3) ic1 -= 3;
+ var ix = Indices[ic];
+ var ix1 = Indices[ic1];
+ var vts = [ix, ix1];
+ var edge = new EdgeObject(vts);
+ var echeck = new EdgeCheckObject(Positions,vts);
+ echeck.edge = edge;
+ echs[ic] = echeck;
+ tri.edges[ic] = edge;
+ edge.polys[0] = tri;
+ TriSet.edges.push(edge);
+ }
+
+ // Place those crossing checkers in a boundary object;
+ // will have to use various kludges since JavaScript doesn't have sets
+ var BoundaryVerts = Indices.slice(0,3);
+ var BoundaryEdges = echs;
+ var Verts = Object;
+ for (var ic=0; ic<3; ic++)
+ {
+ var ic1 = ic + 2;
+ if (ic1 >= 3) ic1 -= 3;
+ var ix = Indices[ic];
+ Verts[ix] = [echs[ic],echs[ic+1]];
+ }
+
+ // Add points until it is no longer possible
+ for (var i=3; i<Indices.length; i++)
+ {
+ var ix = Indices[i];
+
+ // If possible, add the point inside
+ if (AddPointInside(TriSet, ix)) continue;
+
+ // Point was not inside
+ Verts[ix] = [];
+ var NewEdges = [];
+ var VertsAddedTo = [];
+ var EdgesToDelete = [];
+
+ // Find all the non-intersecting edges
+ for (var j=0; j<BoundaryVerts.length; j++)
+ {
+ var ix1 = BoundaryVerts[j];
+ var echk = new EdgeCheckObject(Positions, [ix, ix1]);
+ var DoesIntersect = false;
+ for (var k=0; k<BoundaryEdges.length; k++)
+ {
+ var echk1 = BoundaryEdges[k];
+ DoesIntersect = echk.intersects(Positions,echk1);
+ if (DoesIntersect) break;
+ }
+ if (DoesIntersect) continue;
+
+ var edge = new EdgeObject(echk.verts);
+ echk.edge = edge;
+ AddUniqueEdge(NewEdges,echk);
+ AddUniqueEdge(Verts[ix],echk);
+ AddUnique(VertsAddedTo,ix);
+ AddUniqueEdge(Verts[ix1],echk);
+ AddUnique(VertsAddedTo,ix1);
+ }
+
+ // Add the new vertex itself
+ AddUnique(BoundaryVerts,ix);
+
+ // Find all the triangles
+ for (var j=0; j<BoundaryEdges.length; j++)
+ {
+ var echk = BoundaryEdges[j];
+ var echks = [];
+
+ for (var iv=0; iv<2; iv++)
+ {
+ var vset = FindSharedEdges(Verts[ix],Verts[echk.verts[iv]]);
+ if (vset.length == 0) continue;
+ echks.push(vset[0]);
+ }
+ if (echks.length < 2) continue;
+
+ var empt_indx = -1;
+ for (var iv=0; iv<2; iv++)
+ {
+ if (IsNull(echk.edge.polys[iv]))
+ {
+ empt_indx = iv;
+ break;
+ }
+ }
+ if (empt_indx < 0) continue;
+
+ var oldtri = echk.edge.polys[1-empt_indx];
+ var v0 = echk.verts[0];
+ var i0 = oldtri.VertexIndexIn(v0);
+ var v1 = echk.verts[1];
+ var i1 = oldtri.VertexIndexIn(v1);
+ var i01 = i1 - i0;
+ if (i01 < 0) i01 += 3;
+ if (i01 == 1)
+ {
+ // Order in original: other, v0, v1
+ var NewTriVerts = [ix, v1, v0];
+ }
+ else if (i01 == 2)
+ {
+ // Order in original: other, v1, v0
+ var NewTriVerts = [ix, v0, v1];
+ }
+ var tri = new TriangleObject(Positions, NewTriVerts);
+ if (!tri.IsVertexOrderCorrect()) continue;
+
+ // Add the new triangle
+ // Also, add the new edges,
+ // or remove them from the lists if necessary
+ TriSet.triangles.push(tri);
+ echk.edge.polys[empt_indx] = tri;
+ tri.edges[0] = echk.edge;
+ tri.edges[1] = echks[0].edge;
+ tri.edges[2] = echks[1].edge;
+ AddUniqueEdge(EdgesToDelete, echk);
+
+ for (var iv=0; iv<2; iv++)
+ {
+ var echki = echks[iv];
+ if (IsNull(echki.edge.polys[0]))
+ {
+ echki.edge.polys[0] = tri;
+ TriSet.edges.push(echki.edge);
+ }
+ else
+ {
+ echki.edge.polys[1] = tri;
+ AddUniqueEdge(EdgesToDelete,echki);
+ }
+ }
+ }
+
+ // Add the new edges and remove the edges and vertices
+ // that are now in the interior
+ for (var j=0; j<NewEdges.length; j++)
+ AddUniqueEdge(BoundaryEdges,NewEdges[j]);
+ FindSetDifferenceEdges(BoundaryEdges, EdgesToDelete);
+
+ var BoundaryVertsToRemove = [];
+ for (var j=0; j<VertsAddedTo.length; j++)
+ {
+ var ixa = VertsAddedTo[j];
+ FindSetDifferenceEdges(Verts[ixa],EdgesToDelete);
+ if (Verts[ixa].length == 0)
+ BoundaryVertsToRemove.push(ixa);
+ }
+ FindSetDifference(BoundaryVerts, BoundaryVertsToRemove);
+ }
+
+ // Improve it iteratively
+ ImproveTriangulation(TriSet);
+
+ // Find the boundary line of this region
+ FindConvexHull(TriSet);
+
+ // Find the regions around each point:
+ FindVoronoiDiagram(TriSet);
+
+ return TriSet;
+}
+
+function FindDelaunayTriangulation(Positions)
+{
+ var Indices = new Array(Positions.length);
+ for (var i=0; i<Indices.length; i++)
+ Indices[i] = i;
+
+ return FindDelaunayTriangulationIndexed(Positions, Indices);
+}
+
+//
+// (c) 2016 Philippe Riviere
+//
+// https://github.com/Fil/
+//
+// This software is distributed under the terms of the MIT License
+
+var geoVoronoi = function () {
+ var radians = Math.PI / 180;
+
+ var cartesian = function (spherical) {
+ var lambda = spherical[0] * radians,
+ phi = spherical[1] * radians,
+ cosphi = Math.cos(phi);
+ return [
+ cosphi * Math.cos(lambda),
+ cosphi * Math.sin(lambda),
+ Math.sin(phi)
+ ];
+ };
+
+ var spherical = function (cartesian) {
+ var r = Math.sqrt(cartesian[0] * cartesian[0] + cartesian[1] * cartesian[1]),
+ lat = Math.atan2(cartesian[2], r),
+ lng = Math.atan2(cartesian[1], cartesian[0]);
+ return [lng / radians, lat / radians];
+ };
+
+ var mapline = function (Positions, Verts) {
+ return Verts
+ .map(function (v) {
+ return spherical(Positions[v]);
+ });
+ };
+
+
+ var diagram = d3Voronoi.voronoi()([]);
+
+ var DT = diagram.DT = null,
+ sites = diagram.sites = [],
+ pos = diagram.pos = [],
+ x = function (d) {
+ if (typeof d == 'object' && 'type' in d) {
+ return d3Geo.geoCentroid(d)[0];
+ }
+ if (0 in d) return d[0];
+ },
+ y = function (d) {
+ if (typeof d == 'object' && 'type' in d) {
+ return d3Geo.geoCentroid(d)[1];
+ }
+ if (0 in d) return d[1];
+ };
+
+
+ var voro = function (data) {
+ diagram._hull = diagram._polygons = diagram._links = diagram._triangles = null;
+
+ if (typeof data == 'object' && data.type == 'FeatureCollection') {
+ data = data.features;
+ }
+ sites = data.map(function (site, i) {
+ site.index = i;
+ return site;
+ });
+ pos = data.map(function (site) {
+ return [x(site), y(site)];
+ });
+ DT = FindDelaunayTriangulation(pos.map(cartesian));
+ return diagram;
+ };
+
+ diagram.links = voro.links = function (s) {
+ if (s) voro(s);
+ if (diagram._links) return diagram._links;
+
+ var _index = d3Collection.map();
+
+ var features = DT.edges.map(function (i, n) {
+
+ _index.set(d3Array.extent(i.verts), n);
+
+ var properties = {
+ source: sites[i.verts[0]],
+ target: sites[i.verts[1]],
+ urquhart: true, // will be changed to false later
+ length: d3Geo.geoDistance(pos[i.verts[0]], pos[i.verts[1]])
+ };
+
+ // add left and right sites (?)
+
+ // make GeoJSON
+ return {
+ type: "Feature",
+ geometry: {
+ type: 'LineString',
+ coordinates: [spherical(DT.positions[i.verts[0]]), spherical(DT.positions[i.verts[1]])]
+ },
+ properties: properties
+ };
+ });
+
+ // Urquhart Graph? tag longer link from each triangle
+ DT.triangles.forEach(function (t) {
+ var l = 0,
+ length = 0,
+ remove, v;
+ for (var j = 0; j < 3; j++) {
+ v = d3Array.extent([t.verts[j], t.verts[(j + 1) % 3]]);
+ var n = _index.get(v);
+ length = features[n].properties.length;
+ if (length > l) {
+ l = length;
+ remove = n;
+ }
+ }
+ features[remove].properties.urquhart = false;
+ });
+
+ return diagram._links = {
+ type: "FeatureCollection",
+ features: features
+ };
+ };
+
+ diagram.triangles = voro.triangles = function (s) {
+ if (s) voro(s);
+ if (diagram._triangles) return diagram._triangles;
+
+ var features = DT.triangles
+ .map(function (t) {
+ t.spherical = t.verts.map(function (v) {
+ return DT.positions[v];
+ })
+ .map(spherical);
+
+ // correct winding order
+ if (t.ccdsq < 0) {
+ t.spherical = t.spherical.reverse();
+ t.ccdsq *= -1;
+ }
+
+ return t;
+ })
+
+ // make geojson
+ .map(function (t) {
+ return {
+ type: "Feature",
+ geometry: {
+ type: "Polygon",
+ coordinates: [t.spherical.concat([t.spherical[0]])]
+ },
+ properties: {
+ sites: t.verts.map(function (i) {
+ return sites[i];
+ }),
+ area: t.vol, // steradians
+ circumcenter: spherical(t.ccdir),
+ // ccdsq is *not* the geodesic distance
+ /* circumradius: (2-t.ccdsq) * 53 */
+ }
+ }
+ });
+
+ return diagram._triangles = {
+ type: "FeatureCollection",
+ features: features
+ };
+
+ };
+
+ diagram.polygons = voro.polygons = function (s) {
+ if (s) voro(s);
+ if (diagram._polygons) return diagram._polygons;
+
+ var features = DT.indices.map(function (i, n) {
+ var vor_poly = DT.vor_polygons[DT.indices[i]];
+ var geometry = {};
+
+ if (vor_poly == undefined) {
+ geometry.type = "Sphere";
+ } else {
+ var line = mapline(DT.vor_positions,
+ vor_poly.boundary.concat([vor_poly.boundary[0]])
+ );
+
+ // correct winding order
+ var b = {
+ type: "Polygon",
+ coordinates: [[pos[i], line[0], line[1], pos[i]]]
+ };
+ if (d3Geo.geoArea(b) > 2 * Math.PI + 1e-10) {
+ line = line.reverse();
+ }
+
+ geometry.type = "Polygon";
+ geometry.coordinates = [line];
+ }
+
+ return {
+ type: "Feature",
+ geometry: geometry,
+ properties: {
+ site: sites[i],
+ sitecoordinates: pos[i],
+ neighbours: vor_poly.edges.map(function (e) {
+ return e.verts.filter(function (j) {
+ return j !== i;
+ })[0];
+ })
+ }
+ };
+ });
+
+ return diagram._polygons = {
+ type: "FeatureCollection",
+ features: features
+ };
+ };
+
+ diagram.hull = voro.hull = function (s) {
+ if (s) voro(s);
+ if (diagram._hull) return diagram._hull;
+
+ if (!DT.hull.length) {
+ return null; // What is a null GeoJSON?
+ }
+
+ // seems that DT.hull is always clockwise
+ var hull = DT.hull.reverse();
+
+ // make GeoJSON
+ return diagram._hull = {
+ type: "Feature",
+ geometry: {
+ type: "Polygon",
+ coordinates: [hull.concat([hull[0]]).map(function (i) {
+ return pos[i];
+ })]
+ },
+ properties: {
+ sites: hull.map(function (i) {
+ return sites[i];
+ })
+ }
+ };
+ };
+
+
+ diagram.find = function (x, y, radius) {
+ var features = diagram.polygons().features;
+
+ // optimization: start from most recent result
+ var i, next = diagram.find.found || 0;
+ var cell = features[next] || features[next = 0];
+
+ var dist = d3Geo.geoDistance([x, y], cell.properties.sitecoordinates);
+ do {
+ cell = features[i = next];
+ next = null;
+ cell.properties.neighbours.forEach(function (e) {
+ var ndist = d3Geo.geoDistance([x, y], features[e].properties.sitecoordinates);
+ if (ndist < dist) {
+ dist = ndist;
+ next = e;
+ return;
+ }
+ });
+
+ } while (next !== null);
+ diagram.find.found = i;
+ if (!radius || dist < radius * radius) return cell.properties.site;
+ };
+
+ voro.x = function (f) {
+ if (!f) return x;
+ x = f;
+ return voro;
+ };
+ voro.y = function (f) {
+ if (!f) return y;
+ y = f;
+ return voro;
+ };
+ voro.extent = function (f) {
+ if (!f) return null;
+ return voro;
+ };
+ voro.size = function (f) {
+ if (!f) return null;
+ return voro;
+ };
+
+ return voro;
+
+};
+
+exports.geoVoronoi = geoVoronoi;
+
+Object.defineProperty(exports, '__esModule', { value: true });
+
+})));