// Copyright 2016, Brooks Mershon and Joy Patel (function (global, factory) { typeof exports === 'object' && typeof module !== 'undefined' ? factory(exports, require('d3-polygon')) : typeof define === 'function' && define.amd ? define(['exports', 'd3-polygon'], factory) : (factory((global.partials = global.partials || {}),global.d3_polygon)); }(this, function (exports,d3Polygon) { 'use strict'; /* Computes the dot product between 2 vectors. */ function dot(a, b) { var res = 0.0; if (a.length != b.length){ throw new Error("Invalid dimensions for dot product."); return null; } for (let i = 0; i < a.length; i++) { res += a[i] * b[i]; } return res; } /* Computes the length of vector u. */ function length(u) { return Math.sqrt(dot(u,u)); } /* Returns angle in radians between AB and AC, where a, b, c are specified in counter-clockwise order. */ function angle(a, b, c) { var u = [b[0] - a[0], b[1] - a[1]], v = [c[0] - a[0], c[1] - a[1]]; if ((length(u) * length(v)) == 0) return 0; return Math.acos(Math.max(Math.min(1, dot(u, v) / (length(u) * length(v))), -1)); } function inDelta(actual, expected, epsilon) { return actual < expected + epsilon && actual > expected - epsilon; } /* Given 2 array inputs a and b, this function computes operation between a[i] and b[i] */ function pairwise(a, b, operation){ var res = []; for(let i = 0; i < a.length; i++){ res[i] = operation(a[i], b[i]); } return res; } /* Given 2 array inputs a and b, this function adds a[i] + b[i] */ function add(a, b) { return pairwise(a, b, function(x, y){ return x + y; }); } /* Given 2 array inputs a and b, this function subtracts a[i] - b[i] */ function sub(a, b) { return pairwise(a, b, function(x, y){ return x - y; }); } function scale(s, a) { return [s * a[0], s * a[1]]; } /* Normalizes vector a. */ function normalize(a) { return scale(1 / length(a), a); } /* Returns point of intersection of AB and CD or null if segments do not intersect. */ function intersect(a, b, c, d) { var u = sub(b, a), v = sub(d, c), q = sub(c, a), denom, s, t, epsilon = 0.0; // Cramer's Rule denom = u[0] * (-v[1]) - (-v[0]) * u[1]; s = (q[0] * (-v[1]) - (-v[0]) * q[1]) / (denom); t = (u[0] * q[1] - q[0] * u[1]) / (denom); return (inDelta(s, 0.5, 0.5 + epsilon) && inDelta(t, 0.5, 0.5 + epsilon)) ? add(a, scale(s, u)) : null; } function project(u, v) { return scale(dot(u, v) / dot(v, v), v); } /* Converts radians to degrees. */ function degree(rad) { return rad * 180 / Math.PI; } function rotation(theta) { var rad = theta * Math.PI / 180.0; var s = Math.sin(rad); var c = Math.cos(rad); return [[c, -1.0 * s, 0], [s, c, 0], [0, 0, 1]]; } /* Returns an n by n identity matrix. */ function identity(n) { var eye = [[0, 0, 0], [0, 0, 0], [0, 0, 0]]; for(let i = 0; i < n; i++){ for(let j = 0; j < n; j++){ eye[i][j] = 0.0; } } for(let i = 0; i < n; i++){ eye[i][i] = 1.0; } return eye; } function translation(delta) { var dx = delta[0]; var dy = delta[1]; var res = identity(3); res[0][2] = dx; res[1][2] = dy; return res; } function row(a, row) { var res = []; for(let i = 0; i < a[0].length; i++){ res[i] = a[row][i]; } return res; } /* Multiplies matrix A by vector b. */ function multiply(A, b) { var res = [], x = b.slice(); // pad for homogenous coordinates for (let i = 0; i < 3 - b.length; i++) { x.push(1); } for(let i = 0; i < A.length; i++){ res[i] = dot(row(A, i), x); } return res; } function triangleArea(a,b,c){ var x = length(sub(b,a)); var y = length(sub(c,a)); var z = length(sub(b,c)); var s = 0.5*(x+y+z); return Math.sqrt(s*(s-a)*(s-b)*(s-c)); } /* Returns column col from matrix a. */ function column(a, col) { var res = []; for(let i = 0; i < a.length; i++){ res[i] = a[i][col]; } return res; } /* Multiplies 2 matrices. */ function Multiply(A, B) { if(A[0].length != B.length){ throw new Error("invalid dimensions"); return null; } var res = [[0, 0, 0], [0, 0, 0], [0, 0, 0]]; for(let i = 0; i < A.length; i++){ for(let j = 0; j < B[0].length; j++){ res[i][j] = dot(row(A, i), column(B, j)); } } return res; } function Polygon(P) { this.transforms = []; this._matrix = identity(3); this._origin = this._target = null; } Polygon.prototype = Object.create(Array.prototype); // subclass Array Polygon.prototype.constructor = Polygon; Polygon.prototype.translate = translate; Polygon.prototype.rotate = rotate; Polygon.prototype.accumulate = accumulate; Polygon.prototype.origin = origin; Polygon.prototype.centroid = centroid; Polygon.prototype.rotation = theta; Polygon.prototype.clone = clone; Polygon.prototype.containsPoint = containsPoint; // checks whether polygon contains this point function containsPoint(point){ var myArea = 0.0; var pointArea = 0.0; var a = point; for(let i = 1; i < this.length; i++){ var b = this[i]; var c = this[(i+1 == this.length) ? 0 : i+1]; myArea += triangleArea(a,b,c); } a = this[0]; for(let i = 1; i < this.length-1; i++){ var b = this[i]; var c = this[i+1]; pointArea += triangleArea(a,b,c); } return Math.abs(myArea - pointArea) < 1e-8; } function translate(T) { for (let i = 0, n = this.length; i < n; i++) { let v = this[i]; this[i] = [v[0] + T[0], v[1] + T[1]]; } return this; } function rotate(theta, pivot) { var R = rotation(theta); if (pivot) { this.translate([-pivot[0], -pivot[1]]); } for (let i = 0, n = this.length; i < n; i++) { this[i] = multiply(R, this[i]); } if (pivot) { this.translate(pivot); } return this; } // Returns a polygon translated into the shape it came from. function origin() { if (this._origin) return this._origin; this._origin = this.accumulate(); return this._origin; } function accumulate() { let P = this.clone(), // do not change THIS polygon's geometry n = this.transforms.length, M = identity(3); for (let i = n - 1; i >= 0; i--) { let transform = this.transforms[i]; if (transform.rotate && transform.pivot) { // pivot required P.rotate(transform.rotate, transform.pivot); M = Multiply(translation(scale(-1, transform.pivot)), M); M = Multiply(rotation(transform.rotate), M); M = Multiply(translation(transform.pivot), M); } else if (transform.translate) { P.translate(transform.translate); M = Multiply(translation(transform.translate), M); } } P._matrix = M; return P; } // This polygon's rigid rotation about centroid relative to original placement. function theta() { var a, b; a = this[0]; b = multiply(this.origin()._matrix, a); b = multiply(translation(sub(this.centroid(), this.origin().centroid())), b); return 180 / Math.PI * angle(this.centroid(), a, b.slice(0, 2)); } function centroid() { return d3Polygon.polygonCentroid(this); } function clone() { return polygon(JSON.parse(JSON.stringify(this.slice()))); } // Create new polygon from array of position tuples. function polygon(positions) { var P = new Polygon(); P.push.apply(P, positions); return P; } /* Takes in array of triangle vertices has properties: { transforms, parent } Returns array of polygons with rotations and translations relative to parent. */ function triangle2Rectangle(P) { var index = 0, A, B, C, D, E, F, G, H, M, BC, BA, MB, MC, ME, BCFD, DGB, FCH, maxAngle = -Infinity, rectangle, polygons; if (d3Polygon.polygonArea(P) == 0) return []; // Find largest angle in triangle T for (let i = 0; i < 3; i++) { let theta = 0; theta = angle(P[i % 3], P[(i + 1) % 3], P[(i + 2) % 3]); if (theta > maxAngle) { maxAngle = theta; index = i; } } A = P[index]; B = P[(index + 1) % 3]; C = P[(index + 2) % 3]; BC = sub(C, B); BA = sub(A, B); M = add(B, project(BA, BC)); E = add(A, scale(0.5, sub(M, A))); MB = sub(B, M); MC = sub(C, M); ME = sub(E, M); G = add(M, add(MB, ME)); H = add(M, add(MC, ME)); D = intersect(E, G, A, B); // pivot F = intersect(E, H, A, C); // pivot BCFD = polygon([B, C, F, D]); DGB = polygon([D, G, B]); FCH = polygon([F, C, H]); // transforms to return polygons to previous positions DGB.transforms.push({rotate: degree(Math.PI), pivot: D}); FCH.transforms.push({rotate: degree(-Math.PI), pivot: F}); polygons = [BCFD, DGB, FCH]; polygons.rectangle = [B, C, H, G]; return polygons; }; /* Returns array of polygons resulting from cutting convex polygon P by segment ab. If the polygon is not cut, an empty array is returned. */ function cutPolygon(a, b, P) { var n = P.length, start, end, _P = [], P_ = [], points = [], e, f = P[n - 1], bounds = [], i = -1, rectangle, polygons = []; // walk through polygon edges, attempting intersections // with exact vertices or line segments while (++i < n) { let x; e = f; f = P[i]; // compare floating-point positions by reference if ((a === e || b === e) && !(points.length && Object.is(points[0], e))) x = e; else if ((a === f || b === f) && !(points.length && !Object.is(points[0], f))) x = f; else x = intersect(a, b, e, f); if (!x) continue; // no intersection if (points.length == 2) break; points.push(x); bounds.push(i); } // stitch together two generated polygons if (points.length == 2 && !Object.is(points[0], points[1])) { // stitch half of polygon _P.push(points[0]); i = bounds[0]; while (i != bounds[1]) { // check point is not an EXACT intersection if (!Object.is(points[0], P[i]) && ! Object.is(points[1], P[i])) _P.push(P[i]); i = (i + 1) % n; } _P.push(points[1]); // stitch other half of polygon P_.push(points[0]); P_.push(points[1]); i = bounds[1]; while (i != bounds[0]) { // check point is not an EXACT intersection if (!Object.is(points[0], P[i]) && ! Object.is(points[1], P[i])) P_.push(P[i]); i = (i + 1) % n; } polygons.push(polygon(_P), polygon(P_)); // transfer parent properties polygons.forEach(function(d) { [].push.apply(d.transforms, P.transforms); }); } return polygons; }; /* Returns array of polygons resulting from cutting convex polygons in the collection by segment ab. Polygons that are cut are removed from the collection; polygons that not cut are returned along with the new polygons. */ function cutCollection(a, b, collection) { var N = collection.length, indices = [], polygons = []; for (let i = 0; i < N; i++) { let P = collection[i], generated; generated = cutPolygon(a, b, P); if (generated.length == 2) { [].push.apply(polygons, generated); indices.push(i); } } // include polygons that were not cut into two new polygons polygons = polygons.concat(collection.filter(function(d, i) { return !indices.includes(i); })); return polygons; }; /* Takes in array of polygons forming a canonical rectangle expects rectangle member on collection with bounding rectangle references */ function rectangle2Square(collection) { var dist = [], A, B, C, D, E, F, G, H, J, K, // follows Tralie's labeling AFGD, KCF, area, s, rectangle = collection.rectangle, // bounding rectangle square, slide, l = Infinity, polygons = []; area = Math.abs(d3Polygon.polygonArea(rectangle)); s = Math.sqrt(area); // square side length // assumes escalator conditions hold B = rectangle[0]; // an invariant throughout E = rectangle[3]; // need reference F = rectangle[1]; // need reference G = rectangle[2]; // the square defined by [A, B, C, D] A = add(B, scale(s, normalize(sub(E, B)))); C = add(B, scale(s, normalize(sub(F, B)))); D = add(A, sub(C, B)); l = length(sub(B, F)); // halving the canonical rectangle for the escalator method while (l > 2 * s) { let a, b, left, halved, slideLeft, slideUp, E_Polygon = [], E_Vertex = [], F_Polygon = [], F_Vertex = []; a = add(A, scale(0.5, sub(F, B))); b = add(B, scale(0.5, sub(F, B))); l = length(sub(b, F)); left = [A, B, b, a]; // "overshoot" cut segment endpoints to ensure intersection halved = cutCollection(a, add(b, sub(C, D)), collection); slideUp = sub(E, B); slideLeft = sub(B, b); // translate polgons resulting from the cut (i.e., "stacking the two halves") halved.forEach(function(d, j) { var centroid, T; centroid = d3Polygon.polygonCentroid(d); // update exact references to vertices used in exact intersections d.forEach(function(V, k) { // scan vertices if (Object.is(E, V)) { E_Polygon.push(j); E_Vertex.push(k); } else if (Object.is(F, V)) { F_Polygon.push(j); F_Vertex.push(k); } }); if (d3Polygon.polygonContains(left, centroid)) { // slide "up" d.translate(slideUp).transforms.push({translate: scale(-1, slideUp)}); // undo translate } else { // slide "left" d.translate(slideLeft).transforms.push({translate: scale(-1, slideLeft)}); // undo translate } }); E = halved[E_Polygon[0]][E_Vertex[0]]; F = halved[F_Polygon[0]][F_Vertex[0]]; // make all vertices at F share the same reference (to ensure intersection) for (let i = 1, n = F_Vertex.length; i < n; i++) { F = halved[F_Polygon[i]][F_Vertex[i]] = halved[F_Polygon[i-1]][F_Vertex[i-1]]; } collection = halved; } polygons = collection; // Following "elevator method" assumes rectangle length < 2 x square height G = add(F, sub(E, B)); J = intersect(A, F, E, G); K = intersect(A, F, D, C) KCF = [K, C, F]; // used to locate elevator pieces AFGD = [A, F, G, D]; polygons = cutCollection(A, F, polygons); polygons = cutCollection(D, add(C, scale(1, sub(C, D))), polygons); // slide new polygons using elevator method polygons.forEach(function(d) { var centroid, T; centroid = d3Polygon.polygonCentroid(d); if (d3Polygon.polygonContains(KCF, centroid)) { T = sub(A, K); d.translate(T).transforms.push({translate: scale(-1, T)}); } else if (d3Polygon.polygonContains(AFGD, centroid)) { T = sub(A, J); d.translate(T).transforms.push({translate: scale(-1, T)}); } }); polygons.square = [A, B, C, D]; return polygons; }; /* Takes in array of polygons forming a canonical rectangle and a side length x for the rectangle to be formed. Expects rectangle member on collection with bounding rectangle references. */ function rectangle2Rectangle(collection, x) { var dist = [], A, B, C, D, E, F, G, H, J, K, // follows Tralie's labeling AFGD, KCF, area, width, height, s, rectangle = collection.rectangle, // bounding rectangle square, slide, l = Infinity, polygons = []; area = Math.abs(d3Polygon.polygonArea(rectangle)); s = area / x; height = Math.min(s, x); // square side length width = Math.max(s, x); console.log(area, width, height); // assumes escalator conditions hold B = rectangle[0]; // an invariant throughout E = rectangle[3]; // need reference F = rectangle[1]; // need reference G = rectangle[2]; // the rectangle (of side length x) defined by [A, B, C, D] A = add(B, scale(height, normalize(sub(E, B)))); C = add(B, scale(width, normalize(sub(F, B)))); D = add(A, sub(C, B)); l = length(sub(B, F)); // halving the canonical rectangle for the escalator method while (l > 2 * height) { let a, b, left, halved, slideLeft, slideUp, E_Polygon = [], E_Vertex = [], F_Polygon = [], F_Vertex = []; a = add(A, scale(0.5, sub(F, B))); b = add(B, scale(0.5, sub(F, B))); l = length(sub(b, F)); left = [A, B, b, a]; // "overshoot" cut segment endpoints to ensure intersection halved = cutCollection(a, add(b, sub(C, D)), collection); slideUp = sub(E, B); slideLeft = sub(B, b); // translate polgons resulting from the cut (i.e., "stacking the two halves") halved.forEach(function(d, j) { var centroid, T; centroid = d3Polygon.polygonCentroid(d); // update exact references to vertices used in exact intersections d.forEach(function(V, k) { // scan vertices if (Object.is(E, V)) { E_Polygon.push(j); E_Vertex.push(k); } else if (Object.is(F, V)) { F_Polygon.push(j); F_Vertex.push(k); } }); if (d3Polygon.polygonContains(left, centroid)) { // slide "up" d.translate(slideUp).transforms.push({translate: scale(-1, slideUp)}); // undo translate } else { // slide "left" d.translate(slideLeft).transforms.push({translate: scale(-1, slideLeft)}); // undo translate } }); E = halved[E_Polygon[0]][E_Vertex[0]]; F = halved[F_Polygon[0]][F_Vertex[0]]; // make all vertices at F share the same reference (to ensure intersection) for (let i = 1, n = F_Vertex.length; i < n; i++) { F = halved[F_Polygon[i]][F_Vertex[i]] = halved[F_Polygon[i-1]][F_Vertex[i-1]]; } collection = halved; } polygons = collection; // Following "elevator method" assumes rectangle length < 2 x square height G = add(F, sub(E, B)); J = intersect(A, F, E, G); K = intersect(A, F, D, C) KCF = [K, C, F]; // used to locate elevator pieces AFGD = [A, F, G, D]; polygons = cutCollection(A, F, polygons); polygons = cutCollection(D, add(C, scale(1, sub(C, D))), polygons); // slide new polygons using elevator method polygons.forEach(function(d) { var centroid, T; centroid = d3Polygon.polygonCentroid(d); if (d3Polygon.polygonContains(KCF, centroid)) { T = sub(A, K); d.translate(T).transforms.push({translate: scale(-1, T)}); } else if (d3Polygon.polygonContains(AFGD, centroid)) { T = sub(A, J); d.translate(T).transforms.push({translate: scale(-1, T)}); } }); polygons.square = [D, add(C, scale(1, sub(C, D)))]; return polygons; }; /* Returns point of intersection of segment ab and line cd or null if segments do not intersect. */ function infiniteIntersection(a, b, c, d) { var u = sub(b, a); var v = sub(d, c); var q = sub(a, c); // Cramer's Rule var denom = v[0] * (-u[1]) - (-u[0]) * v[1]; var t = (q[0] * (-u[1]) - (-u[0]) * q[1]) / (denom); var s = (v[0] * q[1] - q[0] * v[1]) / (denom); if(0.0 <= s && s <= 1){ return add(c,scale(t,v)); } return []; } /* Tests if point p is to the right of edge e. p is a single coordinate. e is an array of 2 coordinates. */ function right(p, e) { var u = sub(p, e[0]); var v = sub(e[1], e[0]); var z = u[0]*v[1] - v[0]*u[1]; return z >= 0; } /* Returns the undo operation for a translation or rotation about a pivot. */ function undo(transform){ var undone; if (transform.rotate && transform.pivot) { // pivot required undone = {rotate: -transform.rotate, pivot: transform.pivot}; } else if (transform.translate) { undone = {translate: scale(-1, transform.translate)}; } return undone; } /* Sutherland-Hodgeman algorithm for subject polygon and clip polygon. */ function clipPolygon(subject, clip){ var clipPolygon, outputList, clipped, transforms; clipPolygon = polygon(clip.reverse()); outputList = subject.slice(0); for (let i = 0; i < clipPolygon.length; i++) { let end, clipEdge, inputList, S; end = (i + 1 == clipPolygon.length) ? 0 : i+1; clipEdge = [clipPolygon[i], clipPolygon[end]]; inputList = outputList.slice(0); outputList = []; S = inputList[inputList.length-1]; for (let j = 0; j < inputList.length; j++) { let E, EContains, SContains, inter; E = inputList[j]; EContains = !right(E, clipEdge); SContains = !right(S, clipEdge); if (EContains){ if (!SContains){ inter = infiniteIntersection(E, S, clipPolygon[i], clipPolygon[end]); if (inter != []){ outputList.push(inter); } } outputList.push(E); } else if (SContains) { inter = infiniteIntersection(E, S, clipPolygon[i], clipPolygon[end]); if (inter != []){ outputList.push(inter); } } S = E; } } // transfer transform history clipped = polygon(outputList); transforms = subject.transforms.slice(); transforms.concat(clip.transforms.slice(0).reverse().map(undo)); clipped.transforms = transforms; return clipped; } /* Refer to: http://gamedev.stackexchange.com/questions/13229/sorting-array-of-points-in-clockwise-order Orients points in counter-clockwise order. */ function CCW(poly){ var c = poly.centroid(); poly.sort(function(a,b){ return Math.atan2(b[1] - c[1],b[0] - c[0]) - Math.atan2(a[1] - c[1],a[0] - c[0]); }); return poly; } /* Takes in arrays of polygons A and B and returns all intersections between the two collections using the Sutherland-Hodgeman algorithm. */ function clipCollection(A, B){ var result = []; for (let i = 0; i < A.length; i++) { for (let j = 0; j < B.length; j++) { var sh = clipPolygon(CCW(A[i]), CCW(B[j])); if (sh != [] && sh!=null && sh.length != 0 && sh[0].length != 0){ result.push(sh); } } } return result; } /* Takes in a source and and a subject triangle; */ function triangle2Triangle(source, subject) { var A, B, squareA, squareB, clipped; // TODO rescale subject triangle about centroid in the case that // its area is not equal to that of the source triangle A = rectangle2Square(triangle2Rectangle(source)); B = rectangle2Square(triangle2Rectangle(subject)); squareA = polygon(A.square); squareB = polygon(B.square); clipped = clipCollection(A, B); return clipped; }; exports.angle = angle; exports.intersect = intersect; exports.triangle2Rectangle = triangle2Rectangle; exports.cutPolygon = cutPolygon; exports.cutCollection = cutCollection; exports.rectangle2Square = rectangle2Square; exports.rectangle2Rectangle = rectangle2Rectangle; exports.triangle2Triangle = triangle2Triangle; exports.clipCollection = clipCollection; exports.polygon = polygon; }));