diff --git a/LICENCE.md b/LICENCE.md
index b5aae37..9d4e163 100644
--- a/LICENCE.md
+++ b/LICENCE.md
@@ -47,3 +47,21 @@ shp.ts uses the following third-party code:
> LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
> OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
> SOFTWARE.
+
+### [mapbox earcut](https://github.com/mapbox/earcut)
+
+> ISC License
+>
+> Copyright (c) 2016, Mapbox
+>
+> Permission to use, copy, modify, and/or distribute this software for any purpose
+> with or without fee is hereby granted, provided that the above copyright notice
+> and this permission notice appear in all copies.
+>
+> THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL WARRANTIES WITH
+> REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED WARRANTIES OF MERCHANTABILITY AND
+> FITNESS. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT,
+> INDIRECT, OR CONSEQUENTIAL DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS
+> OF USE, DATA OR PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER
+> TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR PERFORMANCE OF
+> THIS SOFTWARE.
diff --git a/README.md b/README.md
index c32ea42..46936e5 100644
--- a/README.md
+++ b/README.md
@@ -1,13 +1,12 @@
# SHP.ts πΊοΈ
-β οΈβοΈ Careful, still hot, very early stages of development, consume with caution
-
-TypeScript package for loading Esri Shapefiles.
+TypeScript package for loading Esri Shapefiles, primary developed for for WebGL applications
- β
returns a geojson-like representation
- β
supports all shape types (including MultiPatch) per [Esri Shapefile specification](https://www.esri.com/content/dam/esrisites/sitecore-archive/Files/Pdfs/library/whitepapers/pdfs/shapefile.pdf)
- β
supports X, Y, Z, and M coordinates
- β
uses vitest π§ͺ for testing
+- β
includes mapbox's earcut triangulation
## Install from [npm](https://www.npmjs.com/package/shpts)
@@ -15,6 +14,11 @@ TypeScript package for loading Esri Shapefiles.
npm install shpts
```
+| Branch | |
+| ------- | ------------------------------------------------------------------------------------------------------------------------------------------------------------------------ |
+| Release | [![SHPts CI Release](https://github.com/vojtatom/shpts/actions/workflows/ci.yaml/badge.svg?branch=release)](https://github.com/vojtatom/shpts/actions/workflows/ci.yaml) |
+| Dev | [![SHPts CI Dev](https://github.com/vojtatom/shpts/actions/workflows/ci.yaml/badge.svg?branch=dev)](https://github.com/vojtatom/shpts/actions/workflows/ci.yaml) |
+
## Usage
```typescript
@@ -64,6 +68,23 @@ const properties = reader.readRecord(index);
console.log(properties);
```
+triangulate a polygon (expects a set of rings, where the first one is the outer ring and the rest are holes):
+
+```typescript
+import { triangulate, Ring } from 'shpts';
+
+const input: Ring[] = [
+ [
+ [0, 0],
+ [1, 0],
+ [1, 1],
+ [0, 1],
+ ],
+];
+const output = triangulate(input, CoordType.XY);
+//Float32Array([1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 1, 1])
+```
+
## Credits
- insipred by https://github.com/oyvindi/ts-shapefile (MIT Licence), uses all of its test data, partially uses its code
diff --git a/package-lock.json b/package-lock.json
index bdfdba7..57c7794 100644
--- a/package-lock.json
+++ b/package-lock.json
@@ -1,13 +1,16 @@
{
"name": "shpts",
- "version": "1.0.2",
+ "version": "1.0.5",
"lockfileVersion": 3,
"requires": true,
"packages": {
"": {
"name": "shpts",
- "version": "1.0.2",
+ "version": "1.0.5",
"license": "MIT",
+ "dependencies": {
+ "gl-matrix": "^3.4.3"
+ },
"devDependencies": {
"typescript": "^4.9.3",
"vite": "^4.1.0",
@@ -1088,6 +1091,11 @@
"node": "*"
}
},
+ "node_modules/gl-matrix": {
+ "version": "3.4.3",
+ "resolved": "https://registry.npmjs.org/gl-matrix/-/gl-matrix-3.4.3.tgz",
+ "integrity": "sha512-wcCp8vu8FT22BnvKVPjXa/ICBWRq/zjFfdofZy1WSpQZpphblv12/bOQLBC1rMM7SGOFS9ltVmKOHil5+Ml7gA=="
+ },
"node_modules/glob-parent": {
"version": "5.1.2",
"resolved": "https://registry.npmjs.org/glob-parent/-/glob-parent-5.1.2.tgz",
diff --git a/package.json b/package.json
index 5927102..d16b9c3 100644
--- a/package.json
+++ b/package.json
@@ -1,7 +1,7 @@
{
"name": "shpts",
"private": false,
- "version": "1.0.4",
+ "version": "1.0.5",
"type": "module",
"repository": {
"type": "git",
@@ -43,5 +43,8 @@
"vite": "^4.1.0",
"vite-plugin-dts": "^2.0.0-beta.3",
"vitest": "^0.30.1"
+ },
+ "dependencies": {
+ "gl-matrix": "^3.4.3"
}
}
diff --git a/shpts/geometry/multipatch.ts b/shpts/geometry/multipatch.ts
index 9e617ad..c649ec6 100644
--- a/shpts/geometry/multipatch.ts
+++ b/shpts/geometry/multipatch.ts
@@ -93,7 +93,9 @@ export class MultiPatchRecord extends BaseRingedRecord {
let offsetFirst = 1;
let offsetSecond = 2;
for (let i = 0; i < coords.length - 2; i++) {
- polygons.push([[coords[i], coords[i + offsetFirst], coords[i + offsetSecond]]]);
+ polygons.push([
+ [coords[i], coords[i + offsetFirst], coords[i + offsetSecond], coords[i]],
+ ]);
offsetFirst = offsetFirst === 1 ? 2 : 1;
offsetSecond = offsetSecond === 1 ? 2 : 1;
}
@@ -103,7 +105,7 @@ export class MultiPatchRecord extends BaseRingedRecord {
private static triangleFanToPolygon(coords: Coord[]) {
const polygons = [];
for (let i = 1; i < coords.length - 1; i++) {
- polygons.push([[coords[0], coords[i], coords[i + 1]]]);
+ polygons.push([[coords[0], coords[i], coords[i + 1], coords[0]]]);
}
return polygons;
}
diff --git a/shpts/reader/featureReader.ts b/shpts/reader/featureReader.ts
index 3668329..49575f4 100644
--- a/shpts/reader/featureReader.ts
+++ b/shpts/reader/featureReader.ts
@@ -56,6 +56,8 @@ export class FeatureReader {
let attrs = [];
if (this.dbfReader != null) attrs = this.dbfReader.readRecord(index);
+ if (geom == null) return null;
+
return new Feature(geom, attrs, this.dbfReader?.fields);
}
@@ -63,7 +65,7 @@ export class FeatureReader {
const collection = new FeatureCollection();
for (let i = 0; i < this.featureCount; i++) {
const feature = this.readFeature(i);
- collection.features.push(feature);
+ if (feature != null) collection.features.push(feature);
}
return collection;
}
diff --git a/shpts/reader/shpReader.ts b/shpts/reader/shpReader.ts
index e7845eb..de7c03a 100644
--- a/shpts/reader/shpReader.ts
+++ b/shpts/reader/shpReader.ts
@@ -92,6 +92,14 @@ export class ShapeReader {
}
readGeom(geomIndex: number) {
+ try {
+ return this.readGeometryData(geomIndex);
+ } catch (e: any) {
+ return null;
+ }
+ }
+
+ private readGeometryData(geomIndex: number) {
const offset = this.getShpIndex(geomIndex);
this.shpStream.seek(offset);
const recHead = this.readGeomHeader();
diff --git a/shpts/shpts.ts b/shpts/shpts.ts
index 677c46c..4eb3b55 100644
--- a/shpts/shpts.ts
+++ b/shpts/shpts.ts
@@ -11,6 +11,7 @@ import { MultiPatchRecord } from './geometry/multipatch';
import { DbfRecord } from './table/record';
import { DbfFieldDescr, DbfFieldType } from './types/dbfTypes';
import { Coord, CoordType } from './types/coordinate';
+import { triangulate } from './utils/triangulation';
export {
DbfReader,
@@ -25,6 +26,7 @@ export {
MultiPointRecord,
MultiPatchRecord,
CoordType,
+ triangulate,
};
export type { DbfFieldType, DbfFieldDescr, Coord };
diff --git a/shpts/types/coordinate.ts b/shpts/types/coordinate.ts
index c268d3e..05fa15d 100644
--- a/shpts/types/coordinate.ts
+++ b/shpts/types/coordinate.ts
@@ -14,5 +14,5 @@ export type PointCoord = Coord;
export type MultiPointCoord = Coord[];
export type PolyLineCoord = Coord[][];
export type Ring = Coord[];
-export type PolygonCoord = Coord[][][];
-export type MultiPatchCoord = Coord[][][];
+export type PolygonCoord = Ring[][];
+export type MultiPatchCoord = Ring[][];
diff --git a/shpts/utils/mapbox/earcut.js b/shpts/utils/mapbox/earcut.js
new file mode 100644
index 0000000..1d758f5
--- /dev/null
+++ b/shpts/utils/mapbox/earcut.js
@@ -0,0 +1,771 @@
+'use strict';
+
+export function earcut(data, holeIndices, dim) {
+ dim = dim || 2;
+
+ var hasHoles = holeIndices && holeIndices.length,
+ outerLen = hasHoles ? holeIndices[0] * dim : data.length,
+ outerNode = linkedList(data, 0, outerLen, dim, true),
+ triangles = [];
+
+ if (!outerNode || outerNode.next === outerNode.prev) return triangles;
+
+ var minX, minY, maxX, maxY, x, y, invSize;
+
+ if (hasHoles) outerNode = eliminateHoles(data, holeIndices, outerNode, dim);
+
+ // if the shape is not too simple, we'll use z-order curve hash later; calculate polygon bbox
+ if (data.length > 80 * dim) {
+ minX = maxX = data[0];
+ minY = maxY = data[1];
+
+ for (var i = dim; i < outerLen; i += dim) {
+ x = data[i];
+ y = data[i + 1];
+ if (x < minX) minX = x;
+ if (y < minY) minY = y;
+ if (x > maxX) maxX = x;
+ if (y > maxY) maxY = y;
+ }
+
+ // minX, minY and invSize are later used to transform coords into integers for z-order calculation
+ invSize = Math.max(maxX - minX, maxY - minY);
+ invSize = invSize !== 0 ? 32767 / invSize : 0;
+ }
+
+ earcutLinked(outerNode, triangles, dim, minX, minY, invSize, 0);
+
+ return triangles;
+}
+
+// create a circular doubly linked list from polygon points in the specified winding order
+function linkedList(data, start, end, dim, clockwise) {
+ var i, last;
+
+ if (clockwise === signedArea(data, start, end, dim) > 0) {
+ for (i = start; i < end; i += dim) last = insertNode(i, data[i], data[i + 1], last);
+ } else {
+ for (i = end - dim; i >= start; i -= dim) last = insertNode(i, data[i], data[i + 1], last);
+ }
+
+ if (last && equals(last, last.next)) {
+ removeNode(last);
+ last = last.next;
+ }
+
+ return last;
+}
+
+// eliminate colinear or duplicate points
+function filterPoints(start, end) {
+ if (!start) return start;
+ if (!end) end = start;
+
+ var p = start,
+ again;
+ do {
+ again = false;
+
+ if (!p.steiner && (equals(p, p.next) || area(p.prev, p, p.next) === 0)) {
+ removeNode(p);
+ p = end = p.prev;
+ if (p === p.next) break;
+ again = true;
+ } else {
+ p = p.next;
+ }
+ } while (again || p !== end);
+
+ return end;
+}
+
+// main ear slicing loop which triangulates a polygon (given as a linked list)
+function earcutLinked(ear, triangles, dim, minX, minY, invSize, pass) {
+ if (!ear) return;
+
+ // interlink polygon nodes in z-order
+ if (!pass && invSize) indexCurve(ear, minX, minY, invSize);
+
+ var stop = ear,
+ prev,
+ next;
+
+ // iterate through ears, slicing them one by one
+ while (ear.prev !== ear.next) {
+ prev = ear.prev;
+ next = ear.next;
+
+ if (invSize ? isEarHashed(ear, minX, minY, invSize) : isEar(ear)) {
+ // cut off the triangle
+ triangles.push((prev.i / dim) | 0);
+ triangles.push((ear.i / dim) | 0);
+ triangles.push((next.i / dim) | 0);
+
+ removeNode(ear);
+
+ // skipping the next vertex leads to less sliver triangles
+ ear = next.next;
+ stop = next.next;
+
+ continue;
+ }
+
+ ear = next;
+
+ // if we looped through the whole remaining polygon and can't find any more ears
+ if (ear === stop) {
+ // try filtering points and slicing again
+ if (!pass) {
+ earcutLinked(filterPoints(ear), triangles, dim, minX, minY, invSize, 1);
+
+ // if this didn't work, try curing all small self-intersections locally
+ } else if (pass === 1) {
+ ear = cureLocalIntersections(filterPoints(ear), triangles, dim);
+ earcutLinked(ear, triangles, dim, minX, minY, invSize, 2);
+
+ // as a last resort, try splitting the remaining polygon into two
+ } else if (pass === 2) {
+ splitEarcut(ear, triangles, dim, minX, minY, invSize);
+ }
+
+ break;
+ }
+ }
+}
+
+// check whether a polygon node forms a valid ear with adjacent nodes
+function isEar(ear) {
+ var a = ear.prev,
+ b = ear,
+ c = ear.next;
+
+ if (area(a, b, c) >= 0) return false; // reflex, can't be an ear
+
+ // now make sure we don't have other points inside the potential ear
+ var ax = a.x,
+ bx = b.x,
+ cx = c.x,
+ ay = a.y,
+ by = b.y,
+ cy = c.y;
+
+ // triangle bbox; min & max are calculated like this for speed
+ var x0 = ax < bx ? (ax < cx ? ax : cx) : bx < cx ? bx : cx,
+ y0 = ay < by ? (ay < cy ? ay : cy) : by < cy ? by : cy,
+ x1 = ax > bx ? (ax > cx ? ax : cx) : bx > cx ? bx : cx,
+ y1 = ay > by ? (ay > cy ? ay : cy) : by > cy ? by : cy;
+
+ var p = c.next;
+ while (p !== a) {
+ if (
+ p.x >= x0 &&
+ p.x <= x1 &&
+ p.y >= y0 &&
+ p.y <= y1 &&
+ pointInTriangle(ax, ay, bx, by, cx, cy, p.x, p.y) &&
+ area(p.prev, p, p.next) >= 0
+ )
+ return false;
+ p = p.next;
+ }
+
+ return true;
+}
+
+function isEarHashed(ear, minX, minY, invSize) {
+ var a = ear.prev,
+ b = ear,
+ c = ear.next;
+
+ if (area(a, b, c) >= 0) return false; // reflex, can't be an ear
+
+ var ax = a.x,
+ bx = b.x,
+ cx = c.x,
+ ay = a.y,
+ by = b.y,
+ cy = c.y;
+
+ // triangle bbox; min & max are calculated like this for speed
+ var x0 = ax < bx ? (ax < cx ? ax : cx) : bx < cx ? bx : cx,
+ y0 = ay < by ? (ay < cy ? ay : cy) : by < cy ? by : cy,
+ x1 = ax > bx ? (ax > cx ? ax : cx) : bx > cx ? bx : cx,
+ y1 = ay > by ? (ay > cy ? ay : cy) : by > cy ? by : cy;
+
+ // z-order range for the current triangle bbox;
+ var minZ = zOrder(x0, y0, minX, minY, invSize),
+ maxZ = zOrder(x1, y1, minX, minY, invSize);
+
+ var p = ear.prevZ,
+ n = ear.nextZ;
+
+ // look for points inside the triangle in both directions
+ while (p && p.z >= minZ && n && n.z <= maxZ) {
+ if (
+ p.x >= x0 &&
+ p.x <= x1 &&
+ p.y >= y0 &&
+ p.y <= y1 &&
+ p !== a &&
+ p !== c &&
+ pointInTriangle(ax, ay, bx, by, cx, cy, p.x, p.y) &&
+ area(p.prev, p, p.next) >= 0
+ )
+ return false;
+ p = p.prevZ;
+
+ if (
+ n.x >= x0 &&
+ n.x <= x1 &&
+ n.y >= y0 &&
+ n.y <= y1 &&
+ n !== a &&
+ n !== c &&
+ pointInTriangle(ax, ay, bx, by, cx, cy, n.x, n.y) &&
+ area(n.prev, n, n.next) >= 0
+ )
+ return false;
+ n = n.nextZ;
+ }
+
+ // look for remaining points in decreasing z-order
+ while (p && p.z >= minZ) {
+ if (
+ p.x >= x0 &&
+ p.x <= x1 &&
+ p.y >= y0 &&
+ p.y <= y1 &&
+ p !== a &&
+ p !== c &&
+ pointInTriangle(ax, ay, bx, by, cx, cy, p.x, p.y) &&
+ area(p.prev, p, p.next) >= 0
+ )
+ return false;
+ p = p.prevZ;
+ }
+
+ // look for remaining points in increasing z-order
+ while (n && n.z <= maxZ) {
+ if (
+ n.x >= x0 &&
+ n.x <= x1 &&
+ n.y >= y0 &&
+ n.y <= y1 &&
+ n !== a &&
+ n !== c &&
+ pointInTriangle(ax, ay, bx, by, cx, cy, n.x, n.y) &&
+ area(n.prev, n, n.next) >= 0
+ )
+ return false;
+ n = n.nextZ;
+ }
+
+ return true;
+}
+
+// go through all polygon nodes and cure small local self-intersections
+function cureLocalIntersections(start, triangles, dim) {
+ var p = start;
+ do {
+ var a = p.prev,
+ b = p.next.next;
+
+ if (
+ !equals(a, b) &&
+ intersects(a, p, p.next, b) &&
+ locallyInside(a, b) &&
+ locallyInside(b, a)
+ ) {
+ triangles.push((a.i / dim) | 0);
+ triangles.push((p.i / dim) | 0);
+ triangles.push((b.i / dim) | 0);
+
+ // remove two nodes involved
+ removeNode(p);
+ removeNode(p.next);
+
+ p = start = b;
+ }
+ p = p.next;
+ } while (p !== start);
+
+ return filterPoints(p);
+}
+
+// try splitting polygon into two and triangulate them independently
+function splitEarcut(start, triangles, dim, minX, minY, invSize) {
+ // look for a valid diagonal that divides the polygon into two
+ var a = start;
+ do {
+ var b = a.next.next;
+ while (b !== a.prev) {
+ if (a.i !== b.i && isValidDiagonal(a, b)) {
+ // split the polygon in two by the diagonal
+ var c = splitPolygon(a, b);
+
+ // filter colinear points around the cuts
+ a = filterPoints(a, a.next);
+ c = filterPoints(c, c.next);
+
+ // run earcut on each half
+ earcutLinked(a, triangles, dim, minX, minY, invSize, 0);
+ earcutLinked(c, triangles, dim, minX, minY, invSize, 0);
+ return;
+ }
+ b = b.next;
+ }
+ a = a.next;
+ } while (a !== start);
+}
+
+// link every hole into the outer loop, producing a single-ring polygon without holes
+function eliminateHoles(data, holeIndices, outerNode, dim) {
+ var queue = [],
+ i,
+ len,
+ start,
+ end,
+ list;
+
+ for (i = 0, len = holeIndices.length; i < len; i++) {
+ start = holeIndices[i] * dim;
+ end = i < len - 1 ? holeIndices[i + 1] * dim : data.length;
+ list = linkedList(data, start, end, dim, false);
+ if (list === list.next) list.steiner = true;
+ queue.push(getLeftmost(list));
+ }
+
+ queue.sort(compareX);
+
+ // process holes from left to right
+ for (i = 0; i < queue.length; i++) {
+ outerNode = eliminateHole(queue[i], outerNode);
+ }
+
+ return outerNode;
+}
+
+function compareX(a, b) {
+ return a.x - b.x;
+}
+
+// find a bridge between vertices that connects hole with an outer ring and and link it
+function eliminateHole(hole, outerNode) {
+ var bridge = findHoleBridge(hole, outerNode);
+ if (!bridge) {
+ return outerNode;
+ }
+
+ var bridgeReverse = splitPolygon(bridge, hole);
+
+ // filter collinear points around the cuts
+ filterPoints(bridgeReverse, bridgeReverse.next);
+ return filterPoints(bridge, bridge.next);
+}
+
+// David Eberly's algorithm for finding a bridge between hole and outer polygon
+function findHoleBridge(hole, outerNode) {
+ var p = outerNode,
+ hx = hole.x,
+ hy = hole.y,
+ qx = -Infinity,
+ m;
+
+ // find a segment intersected by a ray from the hole's leftmost point to the left;
+ // segment's endpoint with lesser x will be potential connection point
+ do {
+ if (hy <= p.y && hy >= p.next.y && p.next.y !== p.y) {
+ var x = p.x + ((hy - p.y) * (p.next.x - p.x)) / (p.next.y - p.y);
+ if (x <= hx && x > qx) {
+ qx = x;
+ m = p.x < p.next.x ? p : p.next;
+ if (x === hx) return m; // hole touches outer segment; pick leftmost endpoint
+ }
+ }
+ p = p.next;
+ } while (p !== outerNode);
+
+ if (!m) return null;
+
+ // look for points inside the triangle of hole point, segment intersection and endpoint;
+ // if there are no points found, we have a valid connection;
+ // otherwise choose the point of the minimum angle with the ray as connection point
+
+ var stop = m,
+ mx = m.x,
+ my = m.y,
+ tanMin = Infinity,
+ tan;
+
+ p = m;
+
+ do {
+ if (
+ hx >= p.x &&
+ p.x >= mx &&
+ hx !== p.x &&
+ pointInTriangle(hy < my ? hx : qx, hy, mx, my, hy < my ? qx : hx, hy, p.x, p.y)
+ ) {
+ tan = Math.abs(hy - p.y) / (hx - p.x); // tangential
+
+ if (
+ locallyInside(p, hole) &&
+ (tan < tanMin ||
+ (tan === tanMin && (p.x > m.x || (p.x === m.x && sectorContainsSector(m, p)))))
+ ) {
+ m = p;
+ tanMin = tan;
+ }
+ }
+
+ p = p.next;
+ } while (p !== stop);
+
+ return m;
+}
+
+// whether sector in vertex m contains sector in vertex p in the same coordinates
+function sectorContainsSector(m, p) {
+ return area(m.prev, m, p.prev) < 0 && area(p.next, m, m.next) < 0;
+}
+
+// interlink polygon nodes in z-order
+function indexCurve(start, minX, minY, invSize) {
+ var p = start;
+ do {
+ if (p.z === 0) p.z = zOrder(p.x, p.y, minX, minY, invSize);
+ p.prevZ = p.prev;
+ p.nextZ = p.next;
+ p = p.next;
+ } while (p !== start);
+
+ p.prevZ.nextZ = null;
+ p.prevZ = null;
+
+ sortLinked(p);
+}
+
+// Simon Tatham's linked list merge sort algorithm
+// http://www.chiark.greenend.org.uk/~sgtatham/algorithms/listsort.html
+function sortLinked(list) {
+ var i,
+ p,
+ q,
+ e,
+ tail,
+ numMerges,
+ pSize,
+ qSize,
+ inSize = 1;
+
+ do {
+ p = list;
+ list = null;
+ tail = null;
+ numMerges = 0;
+
+ while (p) {
+ numMerges++;
+ q = p;
+ pSize = 0;
+ for (i = 0; i < inSize; i++) {
+ pSize++;
+ q = q.nextZ;
+ if (!q) break;
+ }
+ qSize = inSize;
+
+ while (pSize > 0 || (qSize > 0 && q)) {
+ if (pSize !== 0 && (qSize === 0 || !q || p.z <= q.z)) {
+ e = p;
+ p = p.nextZ;
+ pSize--;
+ } else {
+ e = q;
+ q = q.nextZ;
+ qSize--;
+ }
+
+ if (tail) tail.nextZ = e;
+ else list = e;
+
+ e.prevZ = tail;
+ tail = e;
+ }
+
+ p = q;
+ }
+
+ tail.nextZ = null;
+ inSize *= 2;
+ } while (numMerges > 1);
+
+ return list;
+}
+
+// z-order of a point given coords and inverse of the longer side of data bbox
+function zOrder(x, y, minX, minY, invSize) {
+ // coords are transformed into non-negative 15-bit integer range
+ x = ((x - minX) * invSize) | 0;
+ y = ((y - minY) * invSize) | 0;
+
+ x = (x | (x << 8)) & 0x00ff00ff;
+ x = (x | (x << 4)) & 0x0f0f0f0f;
+ x = (x | (x << 2)) & 0x33333333;
+ x = (x | (x << 1)) & 0x55555555;
+
+ y = (y | (y << 8)) & 0x00ff00ff;
+ y = (y | (y << 4)) & 0x0f0f0f0f;
+ y = (y | (y << 2)) & 0x33333333;
+ y = (y | (y << 1)) & 0x55555555;
+
+ return x | (y << 1);
+}
+
+// find the leftmost node of a polygon ring
+function getLeftmost(start) {
+ var p = start,
+ leftmost = start;
+ do {
+ if (p.x < leftmost.x || (p.x === leftmost.x && p.y < leftmost.y)) leftmost = p;
+ p = p.next;
+ } while (p !== start);
+
+ return leftmost;
+}
+
+// check if a point lies within a convex triangle
+function pointInTriangle(ax, ay, bx, by, cx, cy, px, py) {
+ return (
+ (cx - px) * (ay - py) >= (ax - px) * (cy - py) &&
+ (ax - px) * (by - py) >= (bx - px) * (ay - py) &&
+ (bx - px) * (cy - py) >= (cx - px) * (by - py)
+ );
+}
+
+// check if a diagonal between two polygon nodes is valid (lies in polygon interior)
+function isValidDiagonal(a, b) {
+ return (
+ a.next.i !== b.i &&
+ a.prev.i !== b.i &&
+ !intersectsPolygon(a, b) && // dones't intersect other edges
+ ((locallyInside(a, b) &&
+ locallyInside(b, a) &&
+ middleInside(a, b) && // locally visible
+ (area(a.prev, a, b.prev) || area(a, b.prev, b))) || // does not create opposite-facing sectors
+ (equals(a, b) && area(a.prev, a, a.next) > 0 && area(b.prev, b, b.next) > 0))
+ ); // special zero-length case
+}
+
+// signed area of a triangle
+function area(p, q, r) {
+ return (q.y - p.y) * (r.x - q.x) - (q.x - p.x) * (r.y - q.y);
+}
+
+// check if two points are equal
+function equals(p1, p2) {
+ return p1.x === p2.x && p1.y === p2.y;
+}
+
+// check if two segments intersect
+function intersects(p1, q1, p2, q2) {
+ var o1 = sign(area(p1, q1, p2));
+ var o2 = sign(area(p1, q1, q2));
+ var o3 = sign(area(p2, q2, p1));
+ var o4 = sign(area(p2, q2, q1));
+
+ if (o1 !== o2 && o3 !== o4) return true; // general case
+
+ if (o1 === 0 && onSegment(p1, p2, q1)) return true; // p1, q1 and p2 are collinear and p2 lies on p1q1
+ if (o2 === 0 && onSegment(p1, q2, q1)) return true; // p1, q1 and q2 are collinear and q2 lies on p1q1
+ if (o3 === 0 && onSegment(p2, p1, q2)) return true; // p2, q2 and p1 are collinear and p1 lies on p2q2
+ if (o4 === 0 && onSegment(p2, q1, q2)) return true; // p2, q2 and q1 are collinear and q1 lies on p2q2
+
+ return false;
+}
+
+// for collinear points p, q, r, check if point q lies on segment pr
+function onSegment(p, q, r) {
+ return (
+ q.x <= Math.max(p.x, r.x) &&
+ q.x >= Math.min(p.x, r.x) &&
+ q.y <= Math.max(p.y, r.y) &&
+ q.y >= Math.min(p.y, r.y)
+ );
+}
+
+function sign(num) {
+ return num > 0 ? 1 : num < 0 ? -1 : 0;
+}
+
+// check if a polygon diagonal intersects any polygon segments
+function intersectsPolygon(a, b) {
+ var p = a;
+ do {
+ if (
+ p.i !== a.i &&
+ p.next.i !== a.i &&
+ p.i !== b.i &&
+ p.next.i !== b.i &&
+ intersects(p, p.next, a, b)
+ )
+ return true;
+ p = p.next;
+ } while (p !== a);
+
+ return false;
+}
+
+// check if a polygon diagonal is locally inside the polygon
+function locallyInside(a, b) {
+ return area(a.prev, a, a.next) < 0
+ ? area(a, b, a.next) >= 0 && area(a, a.prev, b) >= 0
+ : area(a, b, a.prev) < 0 || area(a, a.next, b) < 0;
+}
+
+// check if the middle point of a polygon diagonal is inside the polygon
+function middleInside(a, b) {
+ var p = a,
+ inside = false,
+ px = (a.x + b.x) / 2,
+ py = (a.y + b.y) / 2;
+ do {
+ if (
+ p.y > py !== p.next.y > py &&
+ p.next.y !== p.y &&
+ px < ((p.next.x - p.x) * (py - p.y)) / (p.next.y - p.y) + p.x
+ )
+ inside = !inside;
+ p = p.next;
+ } while (p !== a);
+
+ return inside;
+}
+
+// link two polygon vertices with a bridge; if the vertices belong to the same ring, it splits polygon into two;
+// if one belongs to the outer ring and another to a hole, it merges it into a single ring
+function splitPolygon(a, b) {
+ var a2 = new Node(a.i, a.x, a.y),
+ b2 = new Node(b.i, b.x, b.y),
+ an = a.next,
+ bp = b.prev;
+
+ a.next = b;
+ b.prev = a;
+
+ a2.next = an;
+ an.prev = a2;
+
+ b2.next = a2;
+ a2.prev = b2;
+
+ bp.next = b2;
+ b2.prev = bp;
+
+ return b2;
+}
+
+// create a node and optionally link it with previous one (in a circular doubly linked list)
+function insertNode(i, x, y, last) {
+ var p = new Node(i, x, y);
+
+ if (!last) {
+ p.prev = p;
+ p.next = p;
+ } else {
+ p.next = last.next;
+ p.prev = last;
+ last.next.prev = p;
+ last.next = p;
+ }
+ return p;
+}
+
+function removeNode(p) {
+ p.next.prev = p.prev;
+ p.prev.next = p.next;
+
+ if (p.prevZ) p.prevZ.nextZ = p.nextZ;
+ if (p.nextZ) p.nextZ.prevZ = p.prevZ;
+}
+
+function Node(i, x, y) {
+ // vertex index in coordinates array
+ this.i = i;
+
+ // vertex coordinates
+ this.x = x;
+ this.y = y;
+
+ // previous and next vertex nodes in a polygon ring
+ this.prev = null;
+ this.next = null;
+
+ // z-order curve value
+ this.z = 0;
+
+ // previous and next nodes in z-order
+ this.prevZ = null;
+ this.nextZ = null;
+
+ // indicates whether this is a steiner point
+ this.steiner = false;
+}
+
+// return a percentage difference between the polygon area and its triangulation area;
+// used to verify correctness of triangulation
+earcut.deviation = function (data, holeIndices, dim, triangles) {
+ var hasHoles = holeIndices && holeIndices.length;
+ var outerLen = hasHoles ? holeIndices[0] * dim : data.length;
+
+ var polygonArea = Math.abs(signedArea(data, 0, outerLen, dim));
+ if (hasHoles) {
+ for (var i = 0, len = holeIndices.length; i < len; i++) {
+ var start = holeIndices[i] * dim;
+ var end = i < len - 1 ? holeIndices[i + 1] * dim : data.length;
+ polygonArea -= Math.abs(signedArea(data, start, end, dim));
+ }
+ }
+
+ var trianglesArea = 0;
+ for (i = 0; i < triangles.length; i += 3) {
+ var a = triangles[i] * dim;
+ var b = triangles[i + 1] * dim;
+ var c = triangles[i + 2] * dim;
+ trianglesArea += Math.abs(
+ (data[a] - data[c]) * (data[b + 1] - data[a + 1]) -
+ (data[a] - data[b]) * (data[c + 1] - data[a + 1])
+ );
+ }
+
+ return polygonArea === 0 && trianglesArea === 0
+ ? 0
+ : Math.abs((trianglesArea - polygonArea) / polygonArea);
+};
+
+function signedArea(data, start, end, dim) {
+ var sum = 0;
+ for (var i = start, j = end - dim; i < end; i += dim) {
+ sum += (data[j] - data[i]) * (data[i + 1] + data[j + 1]);
+ j = i;
+ }
+ return sum;
+}
+
+// turn a polygon in a multi-dimensional array form (e.g. as in GeoJSON) into a form Earcut accepts
+earcut.flatten = function (data) {
+ var dim = data[0][0].length,
+ result = { vertices: [], holes: [], dimensions: dim },
+ holeIndex = 0;
+
+ for (var i = 0; i < data.length; i++) {
+ for (var j = 0; j < data[i].length; j++) {
+ for (var d = 0; d < dim; d++) result.vertices.push(data[i][j][d]);
+ }
+ if (i > 0) {
+ holeIndex += data[i - 1].length;
+ result.holes.push(holeIndex);
+ }
+ }
+ return result;
+};
diff --git a/shpts/utils/orientation.ts b/shpts/utils/orientation.ts
index dfd7802..5880845 100644
--- a/shpts/utils/orientation.ts
+++ b/shpts/utils/orientation.ts
@@ -1,7 +1,5 @@
import { Coord, PolygonCoord, Ring } from '@shpts/types/coordinate';
-//TODO first and last points are the same, is it in the code?
-
export function isClockwise(ring: Ring) {
if (ring.length < 3) return false; //is hole if degenerate ring
let area = ring[ring.length - 1][1] * ring[0][0] - ring[ring.length - 1][0] * ring[0][1];
@@ -56,7 +54,7 @@ export function assemblePolygonsWithHoles(coords: Coord[][]) {
}
if (unusedHoldes.length) {
- console.warn('Some holes are not in any polygon, inserting as individual polygons.');
+ //console.warn('Some holes are not in any polygon, inserting as individual polygons.');
polygons.push(...unusedHoldes.map((ring) => [ring]));
}
diff --git a/shpts/utils/triangulation.ts b/shpts/utils/triangulation.ts
new file mode 100644
index 0000000..612acf3
--- /dev/null
+++ b/shpts/utils/triangulation.ts
@@ -0,0 +1,126 @@
+import { CoordType, Ring } from '@shpts/types/coordinate';
+import { earcut } from './mapbox/earcut';
+import { vec3, vec2 } from 'gl-matrix';
+
+//find a normal Newell's method
+function normal(points: number[], dim: number) {
+ const normal: vec3 = [0, 0, 0];
+ const n = points.length;
+ let p1x, p1y, p1z, p2x, p2y, p2z, p1, p2;
+ for (let i = 0; i < n; i += dim) {
+ p1 = i;
+ p2 = (i + dim) % n;
+ p1x = points[p1];
+ p1y = points[p1 + 1];
+ p1z = points[p1 + 2];
+ p2x = points[p2];
+ p2y = points[p2 + 1];
+ p2z = points[p2 + 2];
+ normal[0] += (p1y - p2y) * (p1z + p2z);
+ normal[1] += (p1z - p2z) * (p1x + p2x);
+ normal[2] += (p1x - p2x) * (p1y + p2y);
+ }
+
+ vec3.normalize(normal, normal);
+ return normal;
+}
+
+function centroid(points: number[], dim: number) {
+ const n = points.length;
+ let x = 0,
+ y = 0,
+ z = 0;
+
+ for (let i = 0; i < n; i += dim) {
+ x += points[i];
+ y += points[i + 1];
+ z += points[i + 2];
+ }
+
+ return [x / (n / dim), y / (n / dim), z / (n / dim)] as vec3;
+}
+
+function find2DBasis(normal: vec3) {
+ let u: vec3;
+ if (Math.abs(normal[0]) <= Math.abs(normal[1]) && Math.abs(normal[0]) <= Math.abs(normal[2])) {
+ u = [0, -normal[2], normal[1]];
+ } else if (
+ Math.abs(normal[1]) <= Math.abs(normal[0]) &&
+ Math.abs(normal[1]) <= Math.abs(normal[2])
+ ) {
+ u = [-normal[2], 0, normal[0]];
+ } else {
+ u = [-normal[1], normal[0], 0];
+ }
+
+ const v = [
+ normal[1] * u[2] - normal[2] * u[1],
+ normal[2] * u[0] - normal[0] * u[2],
+ normal[0] * u[1] - normal[1] * u[0],
+ ];
+
+ return [u, v] as [vec3, vec3];
+}
+
+function projectTo2D(x: number, y: number, z: number, origin: vec3, basis: [vec3, vec3]) {
+ const u = basis[0];
+ const v = basis[1];
+ const w: vec3 = [x - origin[0], y - origin[1], z - origin[2]];
+ return [vec3.dot(u, w), vec3.dot(v, w)] as vec2;
+}
+
+function isDegenerate(normal: vec3) {
+ return isNaN(vec3.length(normal)) || vec3.length(normal) < 1e-10;
+}
+
+function project(points: number[], type: CoordType) {
+ const dim = type === CoordType.XYZM ? 4 : type === CoordType.XY ? 2 : 3;
+
+ if (dim === 2) return points;
+
+ const normalVector = normal(points, dim);
+ if (isDegenerate(normalVector)) return [];
+
+ const centroidVector = centroid(points, dim);
+ const basis = find2DBasis(normalVector);
+
+ const projected: number[] = [];
+ for (let i = 0; i < points.length; i += dim) {
+ const [x, y] = projectTo2D(points[i], points[i + 1], points[i + 2], centroidVector, basis);
+ projected.push(x, y);
+ }
+
+ return projected;
+}
+
+// Triangulate a polygon
+// Expects a polygon as an array of rings, where each ring is an array of points, where each point is an array of numbers
+// Returns an array of vertices organized as triangles
+// The first and last point of each ring has to be the same
+export function triangulate(rings: Ring[], type: CoordType, ignoreM = false) {
+ const coordinates: number[] = [];
+ const holeIndices: number[] = [];
+ const dim = type === CoordType.XYZM ? 4 : type === CoordType.XY ? 2 : 3;
+
+ rings.forEach((ring, index) => {
+ if (index > 0) holeIndices.push(coordinates.length / dim);
+ ring.forEach((coord, index) => {
+ if (index === ring.length - 1) return;
+ coordinates.push(...coord.slice());
+ });
+ });
+
+ const projected = project(coordinates, type);
+ if (projected.length === 0) return new Float32Array([]); //handle degenerate case
+ const indices = earcut(projected, holeIndices, 2);
+
+ let dimLoc = dim;
+ if (ignoreM) dimLoc = Math.min(dim, 3);
+
+ const vertices = new Float32Array(indices.length * dimLoc);
+ for (let i = 0; i < indices.length; i++) {
+ const index = indices[i];
+ for (let j = 0; j < dimLoc; j++) vertices[i * dimLoc + j] = coordinates[index * dim + j];
+ }
+ return vertices;
+}
diff --git a/test/features.test.ts b/test/features.test.ts
index e1b8b68..8672b15 100644
--- a/test/features.test.ts
+++ b/test/features.test.ts
@@ -1,6 +1,6 @@
import { expect, test } from 'vitest';
import { expectPointsEqual, openFileAsArray } from './utils';
-import { FeatureReader, PolyLineRecord } from '@shpts/shpts';
+import { Feature, FeatureReader, PolyLineRecord } from '@shpts/shpts';
test('Read feature OK', async () => {
const shp = openFileAsArray('testdata/featureclass.shp');
@@ -13,7 +13,8 @@ test('Read feature OK', async () => {
expect(reader.fields?.at(0)?.name).toEqual('Id');
expect(reader.fields?.at(1)?.name).toEqual('name');
- const feature = reader.readFeature(1);
+ const feature = reader.readFeature(1) as Feature;
+ expect(feature).not.toBeNull();
expect(feature.geom).not.toBeNull();
expect(feature.properties).not.toBeNull();
diff --git a/test/multipatch.test.ts b/test/multipatch.test.ts
index b3a78e6..1115956 100644
--- a/test/multipatch.test.ts
+++ b/test/multipatch.test.ts
@@ -16,6 +16,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 0, y: 0, z: 0, m: NaN },
{ x: 0, y: 0, z: 3, m: NaN },
{ x: 5, y: 0, z: 0, m: NaN },
+ { x: 0, y: 0, z: 0, m: NaN },
]);
polygon = geom.coords[1];
@@ -23,6 +24,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 0, y: 0, z: 3, m: NaN },
{ x: 5, y: 0, z: 3, m: NaN },
{ x: 5, y: 0, z: 0, m: NaN },
+ { x: 0, y: 0, z: 3, m: NaN },
]);
polygon = geom.coords[2];
@@ -30,6 +32,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 5, y: 0, z: 0, m: NaN },
{ x: 5, y: 0, z: 3, m: NaN },
{ x: 5, y: 5, z: 0, m: NaN },
+ { x: 5, y: 0, z: 0, m: NaN },
]);
polygon = geom.coords[3];
@@ -37,6 +40,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 5, y: 0, z: 3, m: NaN },
{ x: 5, y: 5, z: 3, m: NaN },
{ x: 5, y: 5, z: 0, m: NaN },
+ { x: 5, y: 0, z: 3, m: NaN },
]);
polygon = geom.coords[4];
@@ -44,6 +48,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 5, y: 5, z: 0, m: NaN },
{ x: 5, y: 5, z: 3, m: NaN },
{ x: 0, y: 5, z: 0, m: NaN },
+ { x: 5, y: 5, z: 0, m: NaN },
]);
polygon = geom.coords[5];
@@ -51,6 +56,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 5, y: 5, z: 3, m: NaN },
{ x: 0, y: 5, z: 3, m: NaN },
{ x: 0, y: 5, z: 0, m: NaN },
+ { x: 5, y: 5, z: 3, m: NaN },
]);
polygon = geom.coords[6];
@@ -58,6 +64,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 0, y: 5, z: 0, m: NaN },
{ x: 0, y: 5, z: 3, m: NaN },
{ x: 0, y: 0, z: 0, m: NaN },
+ { x: 0, y: 5, z: 0, m: NaN },
]);
polygon = geom.coords[7];
@@ -65,6 +72,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 0, y: 5, z: 3, m: NaN },
{ x: 0, y: 0, z: 3, m: NaN },
{ x: 0, y: 0, z: 0, m: NaN },
+ { x: 0, y: 5, z: 3, m: NaN },
]);
polygon = geom.coords[8];
@@ -72,6 +80,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 2.5, y: 2.5, z: 5, m: NaN },
{ x: 0, y: 0, z: 3, m: NaN },
{ x: 5, y: 0, z: 3, m: NaN },
+ { x: 2.5, y: 2.5, z: 5, m: NaN },
]);
polygon = geom.coords[9];
@@ -79,6 +88,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 2.5, y: 2.5, z: 5, m: NaN },
{ x: 5, y: 0, z: 3, m: NaN },
{ x: 5, y: 5, z: 3, m: NaN },
+ { x: 2.5, y: 2.5, z: 5, m: NaN },
]);
polygon = geom.coords[10];
@@ -86,6 +96,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 2.5, y: 2.5, z: 5, m: NaN },
{ x: 5, y: 5, z: 3, m: NaN },
{ x: 0, y: 5, z: 3, m: NaN },
+ { x: 2.5, y: 2.5, z: 5, m: NaN },
]);
polygon = geom.coords[11];
@@ -93,6 +104,7 @@ test('Reading MultiPatchRecord with Z', async () => {
{ x: 2.5, y: 2.5, z: 5, m: NaN },
{ x: 0, y: 5, z: 3, m: NaN },
{ x: 0, y: 0, z: 3, m: NaN },
+ { x: 2.5, y: 2.5, z: 5, m: NaN },
]);
});
diff --git a/test/polygon.test.ts b/test/polygon.test.ts
index 7843ddc..b29715a 100644
--- a/test/polygon.test.ts
+++ b/test/polygon.test.ts
@@ -241,3 +241,11 @@ test('Reading PolygonRecord with Z', async () => {
{ x: 73.2019891835879, y: -126.20161418254617, z: 0, m: NaN },
]);
});
+
+test('Reading PolygonZ Terrain Example', async () => {
+ const shpBuffer = openFileAsArray('testdata/terrain/ter.shp');
+ const shxBuffer = openFileAsArray('testdata/terrain/ter.shx');
+
+ const reader = await ShapeReader.fromArrayBuffer(shpBuffer, shxBuffer);
+ expect(reader.recordCount).toBe(42798);
+});
diff --git a/test/triangulate.test.ts b/test/triangulate.test.ts
new file mode 100644
index 0000000..002472f
--- /dev/null
+++ b/test/triangulate.test.ts
@@ -0,0 +1,184 @@
+import { CoordType, triangulate } from '@shpts/shpts';
+import { Ring } from '@shpts/types/coordinate';
+import { expect, test } from 'vitest';
+
+test('triangulate 2D', async () => {
+ const input: Ring[] = [
+ [
+ [0, 0],
+ [1, 0],
+ [1, 1],
+ [0, 1],
+ [0, 0],
+ ],
+ ];
+ const output = triangulate(input, CoordType.XY);
+ expect(output).toEqual(new Float32Array([1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 1, 1]));
+
+ const output2 = triangulate(input, CoordType.XY, true);
+ expect(output2).toEqual(new Float32Array([1, 1, 0, 1, 0, 0, 0, 0, 1, 0, 1, 1]));
+});
+
+test('triangulate 2D with hole', async () => {
+ const input: Ring[] = [
+ [
+ [0, 0],
+ [1, 0],
+ [1, 1],
+ [0, 1],
+ [0, 0],
+ ],
+ [
+ [0.25, 0.25],
+ [0.75, 0.25],
+ [0.75, 0.75],
+ [0.25, 0.75],
+ [0.25, 0.25],
+ ],
+ ];
+
+ const output = triangulate(input, CoordType.XY);
+ expect(output).toEqual(
+ new Float32Array([
+ 0, 0, 0.25, 0.25, 0.25, 0.75, 0.75, 0.25, 0.25, 0.25, 0, 0, 0, 1, 0, 0, 0.25, 0.75,
+ 0.75, 0.25, 0, 0, 1, 0, 1, 1, 0, 1, 0.25, 0.75, 0.75, 0.75, 0.75, 0.25, 1, 0, 1, 1,
+ 0.25, 0.75, 0.75, 0.75, 0.75, 0.75, 1, 0, 1, 1,
+ ])
+ );
+});
+
+test('triangulate 3D', async () => {
+ const input: Ring[] = [
+ [
+ [0, 0, 0, NaN],
+ [0, 1, 0, NaN],
+ [0, 1, 1, NaN],
+ [0, 0, 1, NaN],
+ [0, 0, 0, NaN],
+ ],
+ ];
+ const output = triangulate(input, CoordType.XYZM);
+ expect(output).toEqual(
+ new Float32Array([
+ 0,
+ 1,
+ 1,
+ NaN,
+ 0,
+ 0,
+ 1,
+ NaN,
+ 0,
+ 0,
+ 0,
+ NaN,
+ 0,
+ 0,
+ 0,
+ NaN,
+ 0,
+ 1,
+ 0,
+ NaN,
+ 0,
+ 1,
+ 1,
+ NaN,
+ ])
+ );
+});
+
+test('triangulate 3D, ignore M', async () => {
+ const input: Ring[] = [
+ [
+ [0, 0, 0, NaN],
+ [0, 1, 0, NaN],
+ [0, 1, 1, NaN],
+ [0, 0, 1, NaN],
+ [0, 0, 0, NaN],
+ ],
+ ];
+ const output = triangulate(input, CoordType.XYZM, true);
+ expect(output).toEqual(
+ new Float32Array([0, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 1, 1])
+ );
+});
+
+test('triangulate 3D with hole', async () => {
+ const input: Ring[] = [
+ [
+ [0, 0, 0, NaN],
+ [0, 1, 0, NaN],
+ [0, 1, 1, NaN],
+ [0, 0, 1, NaN],
+ [0, 0, 0, NaN],
+ ],
+ [
+ [0, 0.25, 0.25, NaN],
+ [0, 0.75, 0.25, NaN],
+ [0, 0.75, 0.75, NaN],
+ [0, 0.25, 0.75, NaN],
+ [0, 0.25, 0.25, NaN],
+ ],
+ ];
+
+ const output = triangulate(input, CoordType.XYZM, true);
+ expect(output).toEqual(
+ new Float32Array([
+ 0, 1, 0, 0, 0.75, 0.25, 0, 0.25, 0.25, 0, 0.75, 0.75, 0, 0.75, 0.25, 0, 1, 0, 0, 0, 0,
+ 0, 1, 0, 0, 0.25, 0.25, 0, 0.75, 0.75, 0, 1, 0, 0, 1, 1, 0, 0, 1, 0, 0, 0, 0, 0.25,
+ 0.25, 0, 0.25, 0.75, 0, 0.75, 0.75, 0, 1, 1, 0, 0, 1, 0, 0.25, 0.25, 0, 0.25, 0.75, 0,
+ 0.25, 0.75, 0, 1, 1, 0, 0, 1,
+ ])
+ );
+});
+
+test('degenerate cases', async () => {
+ //overlaping points
+ let input: Ring[] = [
+ [
+ [0, 0, 0, NaN],
+ [0, 1, 0, NaN],
+ [0, 1, 0, NaN],
+ [0, 0, 0, NaN],
+ ],
+ ];
+
+ let output = triangulate(input, CoordType.XYZM);
+ expect(output).toEqual(new Float32Array([]));
+
+ //multiple overlaping points
+ input = [
+ [
+ [0, 0, 0, NaN],
+ [0, 0, 0, NaN],
+ [0, 0, 0, NaN],
+ [0, 0, 0, NaN],
+ [0, 0, 0, NaN],
+ ],
+ ];
+ output = triangulate(input, CoordType.XYZM);
+ expect(output).toEqual(new Float32Array([]));
+
+ //empty input
+ input = [[]];
+ output = triangulate(input, CoordType.XYZM);
+ expect(output).toEqual(new Float32Array([]));
+
+ //single point
+ input = [[[0, 0, 0, NaN]]];
+ output = triangulate(input, CoordType.XYZM);
+ expect(output).toEqual(new Float32Array([]));
+
+ //colinear points
+ input = [
+ [
+ [0, 0, 0, NaN],
+ [0, 2, 0, NaN],
+ [0, 1, 0, NaN],
+ [0, 0, 0, NaN],
+ ],
+ ];
+ output = triangulate(input, CoordType.XYZM);
+ expect(output).toEqual(new Float32Array([]));
+});
diff --git a/test/utils.ts b/test/utils.ts
index 56f051d..ab41ed2 100644
--- a/test/utils.ts
+++ b/test/utils.ts
@@ -45,6 +45,7 @@ export function expectGeometry(
): TRecord {
let geom = reader.readGeom(index);
expect(geom).toBeDefined();
+ if (!geom) throw new Error('geom is undefined');
expect(geom.coordType).toBe(coordType);
expect(geom instanceof type).toBe(true);
return geom as TRecord;
diff --git a/testdata/terrain/ter.CPG b/testdata/terrain/ter.CPG
new file mode 100755
index 0000000..3ad133c
--- /dev/null
+++ b/testdata/terrain/ter.CPG
@@ -0,0 +1 @@
+UTF-8
\ No newline at end of file
diff --git a/testdata/terrain/ter.dbf b/testdata/terrain/ter.dbf
new file mode 100755
index 0000000..8d970c8
Binary files /dev/null and b/testdata/terrain/ter.dbf differ
diff --git a/testdata/terrain/ter.prj b/testdata/terrain/ter.prj
new file mode 100755
index 0000000..98ce6bb
--- /dev/null
+++ b/testdata/terrain/ter.prj
@@ -0,0 +1 @@
+PROJCS["S-JTSK_Krovak_East_North",GEOGCS["GCS_S_JTSK",DATUM["D_S_JTSK",SPHEROID["Bessel_1841",6377397.155,299.1528128]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Krovak"],PARAMETER["False_Easting",0.0],PARAMETER["False_Northing",0.0],PARAMETER["Pseudo_Standard_Parallel_1",78.5],PARAMETER["Scale_Factor",0.9999],PARAMETER["Azimuth",30.28813975277778],PARAMETER["Longitude_Of_Center",24.83333333333333],PARAMETER["Latitude_Of_Center",49.5],PARAMETER["X_Scale",-1.0],PARAMETER["Y_Scale",1.0],PARAMETER["XY_Plane_Rotation",90.0],UNIT["Meter",1.0]]
\ No newline at end of file
diff --git a/testdata/terrain/ter.sbn b/testdata/terrain/ter.sbn
new file mode 100755
index 0000000..3c1a9b2
Binary files /dev/null and b/testdata/terrain/ter.sbn differ
diff --git a/testdata/terrain/ter.sbx b/testdata/terrain/ter.sbx
new file mode 100755
index 0000000..4d91744
Binary files /dev/null and b/testdata/terrain/ter.sbx differ
diff --git a/testdata/terrain/ter.shp b/testdata/terrain/ter.shp
new file mode 100755
index 0000000..a183bf6
Binary files /dev/null and b/testdata/terrain/ter.shp differ
diff --git a/testdata/terrain/ter.shx b/testdata/terrain/ter.shx
new file mode 100755
index 0000000..7daedaf
Binary files /dev/null and b/testdata/terrain/ter.shx differ
diff --git a/tsconfig.json b/tsconfig.json
index 3cb401e..eeb32e1 100644
--- a/tsconfig.json
+++ b/tsconfig.json
@@ -3,7 +3,7 @@
"target": "ESNext",
"useDefineForClassFields": true,
"lib": ["DOM", "DOM.Iterable", "ESNext"],
- "allowJs": false,
+ "allowJs": true,
"skipLibCheck": true,
"esModuleInterop": false,
"allowSyntheticDefaultImports": true,