From 478e77c7c7773af442ca2b1a346a2ec65393f805 Mon Sep 17 00:00:00 2001 From: gwlucastrig Date: Sun, 24 Dec 2023 11:55:40 -0500 Subject: [PATCH] Issue 104: Include shoreline geometry in output shapefiles --- .../java/org/tinfour/contour/Contour.java | 4 + .../tinfour/contour/ContourBuilderForTin.java | 2 +- .../contour/ContourIntegrityCheck.java | 4 +- .../org/tinfour/contour/ContourRegion.java | 1008 +++++++++-------- .../utils/VisvalingamLineSimplification.java | 5 + .../tinfour/gis/shapefile/DbfFieldDouble.java | 10 +- .../gis/shapefile/ShapefileRecord.java | 5 - .../java/org/tinfour/svm/SvmComputation.java | 2 +- .../java/org/tinfour/svm/SvmContourGraph.java | 161 ++- .../tinfour/svm/properties/SvmProperties.java | 13 + .../org/tinfour/svm/SvmTemplate.properties | 7 + 11 files changed, 700 insertions(+), 521 deletions(-) diff --git a/core/src/main/java/org/tinfour/contour/Contour.java b/core/src/main/java/org/tinfour/contour/Contour.java index 1f955f8f..7f566447 100644 --- a/core/src/main/java/org/tinfour/contour/Contour.java +++ b/core/src/main/java/org/tinfour/contour/Contour.java @@ -253,6 +253,10 @@ public int size() { */ public void complete() { if (closedLoop && n > 6) { + // ensure that there is a "closure" vertex included in the contour. + // If there existing endpoints are numerically close, they are adjusted + // slightly to ensure exact matches. Otherwise, an additional + // vertex is added to the contour. double x0 = xy[0]; double y0 = xy[1]; double x1 = xy[n - 2]; diff --git a/core/src/main/java/org/tinfour/contour/ContourBuilderForTin.java b/core/src/main/java/org/tinfour/contour/ContourBuilderForTin.java index 6f7a75f4..4d863167 100644 --- a/core/src/main/java/org/tinfour/contour/ContourBuilderForTin.java +++ b/core/src/main/java/org/tinfour/contour/ContourBuilderForTin.java @@ -1003,7 +1003,7 @@ private void organizeNestedRegions() { double[] xy = rI.getXY(); for (int j = i + 1; j < nRegion; j++) { ContourRegion rJ = regionList.get(j); - if (rJ.contourRegionType == ContourRegionType.Primary) { + if (rJ.contourRegionType == ContourRegionType.Perimeter) { // regions that include perimeter contours are never // enclosed by other regions. continue; diff --git a/core/src/main/java/org/tinfour/contour/ContourIntegrityCheck.java b/core/src/main/java/org/tinfour/contour/ContourIntegrityCheck.java index 2724c4b3..121bbd93 100644 --- a/core/src/main/java/org/tinfour/contour/ContourIntegrityCheck.java +++ b/core/src/main/java/org/tinfour/contour/ContourIntegrityCheck.java @@ -109,7 +109,7 @@ private boolean checkRegionAreas(Listregions){ double adjustedSum = 0; for(ContourRegion r: regions){ ContourRegionType rt = r.getContourRegionType(); - if(rt == ContourRegionType.Primary){ + if(rt == ContourRegionType.Perimeter){ aSum+=r.getArea(); } adjustedSum+=r.getAdjustedArea(); @@ -140,7 +140,7 @@ private boolean checkRegionAreas(Listregions){ private boolean checkContourTraversal(List regions) { for (ContourRegion r : regions) { ContourRegionType rt = r.getContourRegionType(); - if (rt != ContourRegionType.Primary) { + if (rt != ContourRegionType.Perimeter) { continue; } for(ContourRegionMember member: r.memberList){ diff --git a/core/src/main/java/org/tinfour/contour/ContourRegion.java b/core/src/main/java/org/tinfour/contour/ContourRegion.java index 49c7c9b4..3325ac4c 100644 --- a/core/src/main/java/org/tinfour/contour/ContourRegion.java +++ b/core/src/main/java/org/tinfour/contour/ContourRegion.java @@ -1,489 +1,519 @@ -/* -------------------------------------------------------------------- - * Copyright 2019 Gary W. Lucas. - * - * Licensed under the Apache License, Version 2.0 (the "License"); - * you may not use this file except in compliance with the License. - * You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - * --------------------------------------------------------------------- - */ - - /* - * ----------------------------------------------------------------------- - * - * Revision History: - * Date Name Description - * ------ --------- ------------------------------------------------- - * 07/2019 G. Lucas Created - * - * Notes: - * - * ----------------------------------------------------------------------- - */ -package org.tinfour.contour; - -import java.awt.geom.AffineTransform; -import java.awt.geom.Path2D; -import java.awt.geom.Point2D; -import java.util.ArrayList; -import java.util.List; -import org.tinfour.contour.Contour.ContourType; - -/** - * Provides a elements and access methods for a region created through a - * contour-building process. - */ -public class ContourRegion { - - /** - * An enumeration that indicates the type of a contour region. - */ - public enum ContourRegionType { - /** - * All contours lie entirely within the interior of the TIN and do not - * intersect its perimeter. - */ - Interior, - /** - * At least one contour lies on the perimeter of the TIN. Note that - * primary regions are never enclosed by another region. - */ - Primary - } - - final ContourRegionType contourRegionType; - final List memberList = new ArrayList<>(); - final int regionIndex; - final double area; - final double absArea; - final double xTest; - final double yTest; - - final List children = new ArrayList<>(); - ContourRegion parent; - - int applicationIndex; - - ContourRegion(List memberList, int regionIndex) { - if (memberList.isEmpty()) { - throw new IllegalArgumentException( - "An empty specification for a region geometry is not supported"); - } - - this.regionIndex = regionIndex; - this.memberList.addAll(memberList); - double a = 0; - ContourRegionType rType = ContourRegionType.Interior; - for (ContourRegionMember member : memberList) { - if (member.contour.getContourType() == ContourType.Boundary) { - rType = ContourRegionType.Primary; - } - double s = calculateAreaContribution(member.contour); - if (member.forward) { - a += s; - } else { - a -= s; - } - } - - contourRegionType = rType; - - area = a / 2.0; - absArea = Math.abs(area); - - Contour contour = memberList.get(0).contour; - xTest = (contour.xy[0] + contour.xy[2]) / 2.0; - yTest = (contour.xy[1] + contour.xy[3]) / 2.0; - - } - - /** - * Construct a region based on a single, closed-loop contour. - * - * @param contour a valid instance describing a single, closed-loop contour. - */ - ContourRegion(Contour contour) { - if (contour.getContourType() == ContourType.Interior) { - contourRegionType = ContourRegionType.Interior; - } else { - contourRegionType = ContourRegionType.Primary; - } - - assert contour.closedLoop : "Single contour constructor requires closed loop"; - - memberList.add(new ContourRegionMember(contour, true)); - area = calculateAreaContribution(contour) / 2; - absArea = Math.abs(area); - if (area < 0) { - regionIndex = contour.rightIndex; - } else { - regionIndex = contour.leftIndex; - } - - xTest = (contour.xy[0] + contour.xy[2]) / 2.0; - yTest = (contour.xy[1] + contour.xy[3]) / 2.0; - } - - private double calculateAreaContribution(Contour contour) { - double x0 = contour.xy[0]; - double y0 = contour.xy[1]; - int n = contour.n / 2; - double a = 0; - for (int i = 1; i < n; i++) { - double x1 = contour.xy[i * 2]; - double y1 = contour.xy[i * 2 + 1]; - a += x0 * y1 - x1 * y0; - x0 = x1; - y0 = y1; - } - return a; - } - - /** - * Get the XY coordinates for the contour region. Coordinates - * are stored in a one-dimensional array of doubles in the order: - * 
-   * { (x0,y0), (x1,y1), (x2,y2), etc. }.
-   *
- * - * @return a safe copy of the geometry of the contour region. - */ - public double[] getXY() { - Contour contour = memberList.get(0).contour; - if (memberList.size() == 1 && memberList.get(0).contour.isClosed()) { - return contour.getXY(); - } - int n = 0; - for (ContourRegionMember member : memberList) { - n += member.contour.size() - 1; - } - n++; // closure point - double[] xy = new double[n * 2]; - int k = 0; - for (ContourRegionMember member : memberList) { - contour = member.contour; - n = contour.size(); - if (member.forward) { - for (int i = 0; i < n - 1; i++) { - xy[k++] = contour.xy[i * 2]; - xy[k++] = contour.xy[i * 2 + 1]; - } - } else { - for (int i = n - 1; i > 0; i--) { - xy[k++] = contour.xy[i * 2]; - xy[k++] = contour.xy[i * 2 + 1]; - } - } - } - xy[k++] = xy[0]; - xy[k++] = xy[1]; - - // remove any duplicate points - n = 2; - double px = xy[0]; - double py = xy[1]; - boolean removed = false; - for (int i = 2; i < k; i += 2) { - if (xy[i] == px && xy[i + 1] == py) { - removed = true; // it's a duplicte, remove it - } else { - // it does not match the previous point - px = xy[i]; - py = xy[i + 1]; - if (removed) { - xy[n] = px; - xy[n + 1] = py; - } - n += 2; - } - } - if (removed) { - // copy the coordinates to a downsized array - double[] oldXy = xy; - xy = new double[n]; - System.arraycopy(oldXy, 0, xy, 0, n); - } - return xy; - } - - void addChild(ContourRegion region) { - children.add(region); - region.parent = this; - } - - /** - * Sets the reference to the contour region that encloses the region - * represented by this class. A null parent reference indicates that the - * region is not enclosed by another. - * - * @param parent a valid reference; or a null. - */ - void setParent(ContourRegion parent) { - this.parent = parent; - } - - /** - * Gets the parent region for this region, if any. Regions that are - * enclosed in other regions will have a parent. Regions that are - * not enclosed will not have a parent. - * - * @return a valid instance, or a null. - */ - public ContourRegion getParent() { - return parent; - } - - - /** - * Indicates whether the specified point is inside the region - * - * @param x the Cartesian coordinate for the point of interest - * @param y the Cartesian coordinate for the point of interest - * @return true if the point is inside the contour; otherwise, fakse - */ - public boolean isPointInsideRegion(double x, double y) { - - double[] xy = getXY(); - - return isPointInsideRegion(xy, x, y); - } - - /** - * Indicates whether the specified point is inside a closed polygon. - * - * @param xy an array giving the Cartesian coordinates of the closed, simple - * polygon for the region to be tested. - * @param x the Cartesian coordinate for the point of interest - * @param y the Cartesian coordinate for the point of interest - * @return true if the point is inside the contour; otherwise, fakse - */ - public boolean isPointInsideRegion(double[] xy, double x, double y) { - int rCross = 0; - int lCross = 0; - int n = xy.length / 2; - double x0 = xy[0]; - double y0 = xy[1]; - for (int i = 1; i < n; i++) { - double x1 = xy[i * 2]; - double y1 = xy[i * 2 + 1]; - - double yDelta = y0 - y1; - if (y1 > y != y0 > y) { - double xTest = (x1 * y0 - x0 * y1 + y * (x0 - x1)) / yDelta; - if (xTest > x) { - rCross++; - } - } - if (y1 < y != y0 < y) { - double xTest = (x1 * y0 - x0 * y1 + y * (x0 - x1)) / yDelta; - if (xTest < x) { - lCross++; - } - } - x0 = x1; - y0 = y1; - } - - // (rCross%2) != (lCross%2) - if (((rCross ^ lCross) & 0x01) == 1) { - return false; // on border - } else if ((rCross & 0x01) == 1) { - return true; // unambiguously inside - } - return false; // unambiguously outside - } - - /** - * Gets the absolute value of the overall area of the region. - * No adjustment is made for enclosed regions. - * - * @return a positive value. - */ - public double getAbsoluteArea() { - return absArea; - } - - /** - * Get the area for the region excluding that of any enclosed - * regions. The enclosed regions are not, strictly speaking, - * part of this region and, so, are not included in the adjusted area. - * - * @return a positive value. - */ - public double getAdjustedArea() { - double sumArea = absArea; - for (ContourRegion enclosedRegion : children) { - sumArea -= enclosedRegion.getAbsoluteArea(); - } - return sumArea; - } - - /** - * Gets the signed area of the region. If the points that specify the region - * are given in a counter-clockwise order, the region will have a positive - * area. If the points are given in a clockwise order, the region will have a - * negative area. - * - * @return a signed real value. - */ - public double getArea() { - return area; - } - - /** - * Gets the index of the region. The indexing scheme is based on the original - * values of the zContour array used when the contour regions were built. The - * minimum proper region index is zero. - *

- * At this time, regions are not constructed for areas of null data. In future - * implementations, null-data regions will be indicated by a region index of - * -1. - * - * @return a positive integer value, or -1 for null-data regions. - */ - public int getRegionIndex() { - return regionIndex; - } - - /** - * Gets a Path2D suitable for rendering purposes. The path includes only the - * outer polygon for the region and does not include the internal (nested) - * polygons. Because the contour builder class sorts polygons by descending - * area, this method will general attain the same effect as the conventional - * getPath2D() call when rendering the full set of regions. However, this - * method is not suitable in cases where interior regions are to be omitted - * or filled with a semi-transparent color. - * - * @param transform a valid AffineTransform, typically specified to map the - * Cartesian coordinates of the contour to pixel coordinate. - * @return a valid instance - */ - public Path2D getPath2DWithoutNesting(AffineTransform transform) { - AffineTransform af = transform; - if (af == null) { - af = new AffineTransform(); // identity transform - } - double[] xy = getXY(); - Path2D path = new Path2D.Double(); - appendPathForward(af, path, xy); - return path; - } - - /** - * Gets a Path2D suitable for rendering purposes including both the outer - * polygon and any internal (nested child) polygons. In used for fill - * operations, regions that include nested child regions will be rendered with - * "holes" where the child polygons are indicated. - * - * @param transform a valid AffineTransform, typically specified to map the - * Cartesian coordinates of the contour to pixel coordinate. - * @return a valid instance - */ - public Path2D getPath2D(AffineTransform transform) { - AffineTransform af = transform; - if (af == null) { - af = new AffineTransform(); // identity transform - } - double[] xy = getXY(); - Path2D path = new Path2D.Double(); - path.setWindingRule(Path2D.WIND_EVEN_ODD); - appendPathForward(transform, path, xy); - for (ContourRegion child : children) { - xy = child.getXY(); - appendPathForward(af, path, xy); - } - - return path; - } - - private void appendPathForward( - AffineTransform transform, Path2D path, double[] xy) { - int n = xy.length / 2; - if (n < 2) { - return; - } - - double[] c = new double[n * 2]; - transform.transform(xy, 0, c, 0, n); - - path.moveTo(c[0], c[1]); - for (int i = 1; i < n; i++) { - path.lineTo(c[i * 2], c[i * 2 + 1]); - } - path.closePath(); - } - - /** - * Gets a point lying on one of the segments in the region to support testing - * for polygon enclosures. Note that the test point is never one of the - * vertices of the segment. - * - * @return a valid instance of a Point2D object - */ - public Point2D getTestPoint() { - return new Point2D.Double(xTest, yTest); - } - - /** - * Gets a list of regions that are enclosed by this region. - * This list includes only those regions that are immediately - * enclosed by this region, but does not include any that may - * be "nested" within the enclosed regions. - * @return a valid, potentially empty, list. - */ - public List getEnclosedRegions() { - List nList = new ArrayList<>(); - nList.addAll(children); - return nList; - } - - /** Gets an enumerated value indicating what kind of contour - * region this instance represents. - * - * @return a valid enumeration. - */ - public ContourRegionType getContourRegionType(){ - return this.contourRegionType; - } - - - @Override - public String toString() { - String areaString; - if (absArea > 0.1) { - areaString = String.format("%12.3f", area); - } else { - areaString = String.format("%f", area); - } - return String.format("%4d %s %s %3d", - regionIndex, - areaString, - parent == null ? "root " : "child", - children.size()); - } - - /** - * Allows a calling implementation to set an application-defined index value - * to be used to associate supplemental data with the region - * @param applicationIndex an arbitrary, application-specific index value. - */ - public void setApplicationIndex(int applicationIndex){ - this.applicationIndex = applicationIndex; - } - - /** - * Gets an application-defined index value to be used by the calling - * implementation to associate supplemental data with the region. - * @return an arbitaru, application-specific index value. - */ - public int getApplicationIndex(){ - return this.applicationIndex; - } -} +/* -------------------------------------------------------------------- + * Copyright 2019 Gary W. Lucas. + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0 + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * --------------------------------------------------------------------- + */ + + /* + * ----------------------------------------------------------------------- + * + * Revision History: + * Date Name Description + * ------ --------- ------------------------------------------------- + * 07/2019 G. Lucas Created + * + * Notes: + * + * ----------------------------------------------------------------------- + */ +package org.tinfour.contour; + +import java.awt.geom.AffineTransform; +import java.awt.geom.Path2D; +import java.awt.geom.Point2D; +import java.util.ArrayList; +import java.util.List; +import org.tinfour.contour.Contour.ContourType; + +/** + * Provides a elements and access methods for a region created through a + * contour-building process. + */ +public class ContourRegion { + + /** + * An enumeration that indicates the type of a contour region. + */ + public enum ContourRegionType { + /** + * All contours lie entirely within the interior of the TIN and do not + * intersect its perimeter. + */ + Interior, + /** + * At least one contour lies on the perimeter of the TIN. Note that + * primary regions are never enclosed by another region. + */ + Perimeter + } + + final ContourRegionType contourRegionType; + final List memberList = new ArrayList<>(); + final int regionIndex; + final double area; + final double absArea; + final double xTest; + final double yTest; + + final List children = new ArrayList<>(); + ContourRegion parent; + + int applicationIndex; + + ContourRegion(List memberList, int regionIndex) { + if (memberList.isEmpty()) { + throw new IllegalArgumentException( + "An empty specification for a region geometry is not supported"); + } + + this.regionIndex = regionIndex; + this.memberList.addAll(memberList); + double a = 0; + ContourRegionType rType = ContourRegionType.Interior; + for (ContourRegionMember member : memberList) { + if (member.contour.getContourType() == ContourType.Boundary) { + rType = ContourRegionType.Perimeter; + } + double s = calculateAreaContribution(member.contour); + if (member.forward) { + a += s; + } else { + a -= s; + } + } + + contourRegionType = rType; + + area = a / 2.0; + absArea = Math.abs(area); + + Contour contour = memberList.get(0).contour; + xTest = (contour.xy[0] + contour.xy[2]) / 2.0; + yTest = (contour.xy[1] + contour.xy[3]) / 2.0; + + } + + /** + * Construct a region based on a single, closed-loop contour. + * + * @param contour a valid instance describing a single, closed-loop contour. + */ + public ContourRegion(Contour contour) { + if (contour.getContourType() == ContourType.Interior) { + contourRegionType = ContourRegionType.Interior; + } else { + contourRegionType = ContourRegionType.Perimeter; + } + + assert contour.closedLoop : "Single contour constructor requires closed loop"; + + memberList.add(new ContourRegionMember(contour, true)); + area = calculateAreaContribution(contour) / 2; + absArea = Math.abs(area); + if (area < 0) { + regionIndex = contour.rightIndex; + } else { + regionIndex = contour.leftIndex; + } + + xTest = (contour.xy[0] + contour.xy[2]) / 2.0; + yTest = (contour.xy[1] + contour.xy[3]) / 2.0; + } + + private double calculateAreaContribution(Contour contour) { + double x0 = contour.xy[0]; + double y0 = contour.xy[1]; + int n = contour.n / 2; + double a = 0; + for (int i = 1; i < n; i++) { + double x1 = contour.xy[i * 2]; + double y1 = contour.xy[i * 2 + 1]; + a += x0 * y1 - x1 * y0; + x0 = x1; + y0 = y1; + } + return a; + } + + /** + * Get the XY coordinates for the contour region. Coordinates + * are stored in a one-dimensional array of doubles in the order: + * 

+   * { (x0,y0), (x1,y1), (x2,y2), etc. }.
+   *
+ * + * @return a safe copy of the geometry of the contour region. + */ + public double[] getXY() { + Contour contour = memberList.get(0).contour; + if (memberList.size() == 1 && memberList.get(0).contour.isClosed()) { + return contour.getXY(); + } + int n = 0; + for (ContourRegionMember member : memberList) { + n += member.contour.size() - 1; + } + n++; // closure point + double[] xy = new double[n * 2]; + int k = 0; + for (ContourRegionMember member : memberList) { + contour = member.contour; + n = contour.size(); + if (member.forward) { + for (int i = 0; i < n - 1; i++) { + xy[k++] = contour.xy[i * 2]; + xy[k++] = contour.xy[i * 2 + 1]; + } + } else { + for (int i = n - 1; i > 0; i--) { + xy[k++] = contour.xy[i * 2]; + xy[k++] = contour.xy[i * 2 + 1]; + } + } + } + xy[k++] = xy[0]; + xy[k++] = xy[1]; + + // remove any duplicate points + n = 2; + double px = xy[0]; + double py = xy[1]; + boolean removed = false; + for (int i = 2; i < k; i += 2) { + if (xy[i] == px && xy[i + 1] == py) { + removed = true; // it's a duplicte, remove it + } else { + // it does not match the previous point + px = xy[i]; + py = xy[i + 1]; + if (removed) { + xy[n] = px; + xy[n + 1] = py; + } + n += 2; + } + } + if (removed) { + // copy the coordinates to a downsized array + double[] oldXy = xy; + xy = new double[n]; + System.arraycopy(oldXy, 0, xy, 0, n); + } + return xy; + } + + /** + * Adds a child (nested) region to the internal list and sets + * the child-region's parent link. + * @param region a valid reference. + */ + public void addChild(ContourRegion region) { + children.add(region); + region.parent = this; + } + + /** + * Removes all child (nested) regions from the internal list and + * nullifies any of the child-region parent links. + */ + public void removeChildren(){ + for(ContourRegion c: children){ + c.setParent(null); + } + children.clear(); + } + /** + * Sets the reference to the contour region that encloses the region + * represented by this class. A null parent reference indicates that the + * region is not enclosed by another. + * + * @param parent a valid reference; or a null. + */ + public void setParent(ContourRegion parent) { + this.parent = parent; + } + + /** + * Gets the parent region for this region, if any. Regions that are + * enclosed in other regions will have a parent. Regions that are + * not enclosed will not have a parent. + * + * @return a valid instance, or a null. + */ + public ContourRegion getParent() { + return parent; + } + + + /** + * Indicates whether the specified point is inside the region + * + * @param x the Cartesian coordinate for the point of interest + * @param y the Cartesian coordinate for the point of interest + * @return true if the point is inside the contour; otherwise, fakse + */ + public boolean isPointInsideRegion(double x, double y) { + + double[] xy = getXY(); + + return isPointInsideRegion(xy, x, y); + } + + /** + * Indicates whether the specified point is inside a closed polygon. + * + * @param xy an array giving the Cartesian coordinates of the closed, simple + * polygon for the region to be tested. + * @param x the Cartesian coordinate for the point of interest + * @param y the Cartesian coordinate for the point of interest + * @return true if the point is inside the contour; otherwise, fakse + */ + public boolean isPointInsideRegion(double[] xy, double x, double y) { + int rCross = 0; + int lCross = 0; + int n = xy.length / 2; + double x0 = xy[0]; + double y0 = xy[1]; + for (int i = 1; i < n; i++) { + double x1 = xy[i * 2]; + double y1 = xy[i * 2 + 1]; + + double yDelta = y0 - y1; + if (y1 > y != y0 > y) { + double xTest = (x1 * y0 - x0 * y1 + y * (x0 - x1)) / yDelta; + if (xTest > x) { + rCross++; + } + } + if (y1 < y != y0 < y) { + double xTest = (x1 * y0 - x0 * y1 + y * (x0 - x1)) / yDelta; + if (xTest < x) { + lCross++; + } + } + x0 = x1; + y0 = y1; + } + + // (rCross%2) != (lCross%2) + if (((rCross ^ lCross) & 0x01) == 1) { + return false; // on border + } else if ((rCross & 0x01) == 1) { + return true; // unambiguously inside + } + return false; // unambiguously outside + } + + /** + * Gets the absolute value of the overall area of the region. + * No adjustment is made for enclosed regions. + * + * @return a positive value. + */ + public double getAbsoluteArea() { + return absArea; + } + + /** + * Get the area for the region excluding that of any enclosed + * regions. The enclosed regions are not, strictly speaking, + * part of this region and, so, are not included in the adjusted area. + * + * @return a positive value. + */ + public double getAdjustedArea() { + double sumArea = absArea; + for (ContourRegion enclosedRegion : children) { + sumArea -= enclosedRegion.getAbsoluteArea(); + } + return sumArea; + } + + /** + * Gets the signed area of the region. If the points that specify the region + * are given in a counter-clockwise order, the region will have a positive + * area. If the points are given in a clockwise order, the region will have a + * negative area. + * + * @return a signed real value. + */ + public double getArea() { + return area; + } + + /** + * Gets the absolute value of the area for this region. + * @return a positve value greater than zero + */ + public double getAbsArea(){ + return absArea; + } + /** + * Gets the index of the region. The indexing scheme is based on the original + * values of the zContour array used when the contour regions were built. The + * minimum proper region index is zero. + *

+ * At this time, regions are not constructed for areas of null data. In future + * implementations, null-data regions will be indicated by a region index of + * -1. + * + * @return a positive integer value, or -1 for null-data regions. + */ + public int getRegionIndex() { + return regionIndex; + } + + /** + * Gets a Path2D suitable for rendering purposes. The path includes only the + * outer polygon for the region and does not include the internal (nested) + * polygons. Because the contour builder class sorts polygons by descending + * area, this method will general attain the same effect as the conventional + * getPath2D() call when rendering the full set of regions. However, this + * method is not suitable in cases where interior regions are to be omitted + * or filled with a semi-transparent color. + * + * @param transform a valid AffineTransform, typically specified to map the + * Cartesian coordinates of the contour to pixel coordinate. + * @return a valid instance + */ + public Path2D getPath2DWithoutNesting(AffineTransform transform) { + AffineTransform af = transform; + if (af == null) { + af = new AffineTransform(); // identity transform + } + double[] xy = getXY(); + Path2D path = new Path2D.Double(); + appendPathForward(af, path, xy); + return path; + } + + /** + * Gets a Path2D suitable for rendering purposes including both the outer + * polygon and any internal (nested child) polygons. In used for fill + * operations, regions that include nested child regions will be rendered with + * "holes" where the child polygons are indicated. + * + * @param transform a valid AffineTransform, typically specified to map the + * Cartesian coordinates of the contour to pixel coordinate. + * @return a valid instance + */ + public Path2D getPath2D(AffineTransform transform) { + AffineTransform af = transform; + if (af == null) { + af = new AffineTransform(); // identity transform + } + double[] xy = getXY(); + Path2D path = new Path2D.Double(); + path.setWindingRule(Path2D.WIND_EVEN_ODD); + appendPathForward(transform, path, xy); + for (ContourRegion child : children) { + xy = child.getXY(); + appendPathForward(af, path, xy); + } + + return path; + } + + private void appendPathForward( + AffineTransform transform, Path2D path, double[] xy) { + int n = xy.length / 2; + if (n < 2) { + return; + } + + double[] c = new double[n * 2]; + transform.transform(xy, 0, c, 0, n); + + path.moveTo(c[0], c[1]); + for (int i = 1; i < n; i++) { + path.lineTo(c[i * 2], c[i * 2 + 1]); + } + path.closePath(); + } + + /** + * Gets a point lying on one of the segments in the region to support testing + * for polygon enclosures. Note that the test point is never one of the + * vertices of the segment. + * + * @return a valid instance of a Point2D object + */ + public Point2D getTestPoint() { + return new Point2D.Double(xTest, yTest); + } + + /** + * Gets a list of regions that are enclosed by this region. + * This list includes only those regions that are immediately + * enclosed by this region, but does not include any that may + * be "nested" within the enclosed regions. + * @return a valid, potentially empty, list. + */ + public List getEnclosedRegions() { + List nList = new ArrayList<>(); + nList.addAll(children); + return nList; + } + + /** Gets an enumerated value indicating what kind of contour + * region this instance represents. + * + * @return a valid enumeration. + */ + public ContourRegionType getContourRegionType(){ + return this.contourRegionType; + } + + + @Override + public String toString() { + String areaString; + if (absArea > 0.1) { + areaString = String.format("%12.3f", area); + } else { + areaString = String.format("%f", area); + } + return String.format("%4d %s %s %3d", + regionIndex, + areaString, + parent == null ? "root " : "child", + children.size()); + } + + /** + * Allows a calling implementation to set an application-defined index value + * to be used to associate supplemental data with the region + * @param applicationIndex an arbitrary, application-specific index value. + */ + public void setApplicationIndex(int applicationIndex){ + this.applicationIndex = applicationIndex; + } + + /** + * Gets an application-defined index value to be used by the calling + * implementation to associate supplemental data with the region. + * @return an arbitaru, application-specific index value. + */ + public int getApplicationIndex(){ + return this.applicationIndex; + } + + public ContourRegionType getRegionType(){ + return contourRegionType; + } + + public boolean hasChildren(){ + return !children.isEmpty(); + } +} diff --git a/core/src/main/java/org/tinfour/utils/VisvalingamLineSimplification.java b/core/src/main/java/org/tinfour/utils/VisvalingamLineSimplification.java index 45476f77..a3244f8f 100644 --- a/core/src/main/java/org/tinfour/utils/VisvalingamLineSimplification.java +++ b/core/src/main/java/org/tinfour/utils/VisvalingamLineSimplification.java @@ -100,6 +100,11 @@ public int simplify(int nPoints, double[] xy, double areaThreshold) { nOperations++; nProcessed+=nPoints; + if (areaThreshold == 0) { + nRemoved += (nPoints - n); + return nPoints; + } + double[] a = new double[n]; a[0] = Double.POSITIVE_INFINITY; a[nPoints - 1] = Double.POSITIVE_INFINITY; diff --git a/gis/src/main/java/org/tinfour/gis/shapefile/DbfFieldDouble.java b/gis/src/main/java/org/tinfour/gis/shapefile/DbfFieldDouble.java index 1299b25e..19e3f3ec 100644 --- a/gis/src/main/java/org/tinfour/gis/shapefile/DbfFieldDouble.java +++ b/gis/src/main/java/org/tinfour/gis/shapefile/DbfFieldDouble.java @@ -31,6 +31,7 @@ import java.io.IOException; import java.util.Arrays; +import java.util.Locale; import org.tinfour.io.BufferedRandomAccessFile; import org.tinfour.io.BufferedRandomAccessReader; @@ -45,6 +46,7 @@ public class DbfFieldDouble extends DbfField { private String writingFormat; + DbfFieldDouble( String name, char fieldType, @@ -56,9 +58,9 @@ public class DbfFieldDouble extends DbfField { super(name, fieldType, dataAddress, fieldLength, fieldDecimalCount, offset); engineeringNotation = useEngineeringNotation; if (useEngineeringNotation) { - writingFormat = String.format("%%%d.%de", fieldLength, fieldLength-7); + writingFormat = String.format(Locale.US, "%%%d.%de", fieldLength, fieldLength-7); } else { - writingFormat = String.format("%%%d.%df", fieldLength, fieldDecimalCount); + writingFormat = String.format(Locale.US, "%%%d.%df", fieldLength, fieldDecimalCount); } } @@ -175,12 +177,14 @@ void read(BufferedRandomAccessReader brad, long recordFilePos) throws IOExceptio @Override void write(BufferedRandomAccessFile braf) throws IOException { - String s = String.format(writingFormat, value); + String s = String.format(Locale.US, writingFormat, value); + byte[] b = s.getBytes(); if(b.length>fieldLength){ throw new IOException("Formatted output exceeds fieldLength of " +fieldLength+": \""+s+"\""); } + braf.write(b, 0, fieldLength); } diff --git a/gis/src/main/java/org/tinfour/gis/shapefile/ShapefileRecord.java b/gis/src/main/java/org/tinfour/gis/shapefile/ShapefileRecord.java index 3947e010..5f7e67ff 100644 --- a/gis/src/main/java/org/tinfour/gis/shapefile/ShapefileRecord.java +++ b/gis/src/main/java/org/tinfour/gis/shapefile/ShapefileRecord.java @@ -125,11 +125,6 @@ public void addPolygon(int nPointsInput, double[] xyzInput, boolean nested) { } int nCoordinates = this.shapefileType.hasZ() ? 3 : 2; int n = nPointsInput; - int index = (nPointsInput - 1) * nCoordinates; - if (xyzInput[index] == xyzInput[0] && xyzInput[index + 1] == xyzInput[1]) { - // the calling application included a loop-closing point in the specification - n--; - } // beacuse we will be subtracting the anchor coordinates // from all points, the first and last line segments will drop out diff --git a/svm/src/main/java/org/tinfour/svm/SvmComputation.java b/svm/src/main/java/org/tinfour/svm/SvmComputation.java index 16799cc5..d3c29c55 100644 --- a/svm/src/main/java/org/tinfour/svm/SvmComputation.java +++ b/svm/src/main/java/org/tinfour/svm/SvmComputation.java @@ -342,7 +342,7 @@ public void processVolume( // During the flat-fixer loop, the total count of triangles // may bounce around a bit. In the case of coves and similar // features, fixing one layer of flats may expose yet more - // flats to be fixed. Also, not that the counts/areas are + // flats to be fixed. Also, note that the counts/areas are // not the total count/area of flat triangles, but the total // of those triangles that were subject to remediation. Counting // the actual flat area would add too much processing for too diff --git a/svm/src/main/java/org/tinfour/svm/SvmContourGraph.java b/svm/src/main/java/org/tinfour/svm/SvmContourGraph.java index 6169a099..3e6115a8 100644 --- a/svm/src/main/java/org/tinfour/svm/SvmContourGraph.java +++ b/svm/src/main/java/org/tinfour/svm/SvmContourGraph.java @@ -39,12 +39,14 @@ import java.awt.font.TextLayout; import java.awt.geom.AffineTransform; import java.awt.geom.Path2D; +import java.awt.geom.Point2D; import java.awt.geom.Rectangle2D; import java.awt.image.BufferedImage; import java.io.File; import java.io.IOException; import java.io.PrintStream; import java.util.ArrayList; +import java.util.Collections; import java.util.List; import javax.imageio.ImageIO; import org.tinfour.common.IConstraint; @@ -57,6 +59,7 @@ import org.tinfour.contour.ContourBuilderForTin; import org.tinfour.contour.ContourIntegrityCheck; import org.tinfour.contour.ContourRegion; +import org.tinfour.contour.ContourRegion.ContourRegionType; import org.tinfour.gis.shapefile.ShapefileRecord; import org.tinfour.gis.shapefile.ShapefileType; import org.tinfour.gis.shapefile.ShapefileWriter; @@ -251,7 +254,7 @@ static void write( // For different data sets, we will need different contour intervals. // Attempt to create a specification with about 10 contour intervals // using the axis too. Then, take the results from the axis tool - // and create our zContour array. Not that the zContours must be + // and create our zContour array. Note that the zContours must be // greater than the zMin value and less than the shoreReferenceElevation double[] aArray = null; double contourInterval = properties.getContourGraphInterval(); @@ -536,16 +539,59 @@ static void write( + ", " + ioex.getMessage()); } - File shapefileRef = properties.getContourRegionShapefile(); - if (shapefileRef != null) { - removeOldShapefiles(ps, shapefileRef); - ps.println("Writing shapefile " + shapefileRef.getPath()); + File regionShapefileRef = properties.getContourRegionShapefile(); + File contourShapefileRef = properties.getContourLineShapefile(); + boolean appendShoreline = properties.isContourShapefileShorelineEnabled(); + if(regionShapefileRef==null || contourShapefileRef==null){ + // nothing more to do. + return; + } + + List shorelineContours = new ArrayList<>(); + List shorelineRegions = new ArrayList<>(); + if (appendShoreline) { + for (PolygonConstraint pc : boundaryConstraints) { + double pcArea = pc.getArea(); + int leftIndex, rightIndex; + if (pcArea > 0) { + // enclosing polygon + leftIndex = zContour.length; + rightIndex = zContour.length + 1; + } else { + leftIndex = zContour.length; + rightIndex = zContour.length + 1; + } + + Contour contour = new Contour(leftIndex, rightIndex, shoreReferenceElevation, true); + List pcvList = pc.getVertices(); + for (Vertex v : pcvList) { + contour.add(v.getX(), v.getY()); + } + contour.complete(); + ContourRegion cr = new ContourRegion(contour); + shorelineRegions.add(cr); + shorelineContours.add(contour); + contours.add(contour); + } + } + + + + if (regionShapefileRef != null) { + removeOldShapefiles(ps, regionShapefileRef); + + ps.println("Writing shapefile " + regionShapefileRef.getPath()); + if(appendShoreline){ + organizeNestedRegions(regions, shorelineRegions); + regions = shorelineRegions; + } int iRegion = 0; for (ContourRegion region : regions) { - int rIndex = region.getRegionIndex(); - if (rIndex >= iN) { + if(region.getContourRegionType()==ContourRegionType.Perimeter){ + // SVM contouring will produce a perimeter region that is + // constructed based on the convex-hull of the TIN. continue; } iRegion++; @@ -557,18 +603,22 @@ static void write( regionSpec.addIntegerField("feature_id", 8); regionSpec.addIntegerField("parent_id", 8); regionSpec.addIntegerField("band_idx", 4); - regionSpec.addFloatingPointField("band_min", 8, 2, false); - regionSpec.addFloatingPointField("band_max", 8, 2, false); - + regionSpec.addFloatingPointField("band_min", 9, 3, false); + regionSpec.addFloatingPointField("band_max", 9, 3, false); regionSpec.addFloatingPointField("Shape_area", 13, 6, true); regionSpec.setShapefilePrjContent(data.getShapefilePrjContent()); - try (ShapefileWriter regionWriter = new ShapefileWriter(shapefileRef, regionSpec);) { + try (ShapefileWriter regionWriter = new ShapefileWriter(regionShapefileRef, regionSpec);) { for (ContourRegion region : regions) { - int rIndex = region.getRegionIndex(); - if (rIndex >= iN) { + if (region.getContourRegionType() == ContourRegionType.Perimeter) { + // SVM contouring will produce a perimeter region that is + // constructed based on the convex-hull of the TIN. continue; } + int rIndex = region.getRegionIndex(); + if (rIndex > iN) { + continue; // should never happen + } double bandMin = zBandMin[rIndex]; double bandMax = zBandMax[rIndex]; @@ -608,27 +658,32 @@ static void write( } } - shapefileRef = properties.getContourLineShapefile(); - if (shapefileRef != null) { - removeOldShapefiles(ps, shapefileRef); - ps.println("Writing shapefile " + shapefileRef.getPath()); + + if (contourShapefileRef != null) { + removeOldShapefiles(ps, contourShapefileRef); + ps.println("Writing shapefile " + contourShapefileRef.getPath()); ShapefileWriterSpecification contourSpec = new ShapefileWriterSpecification(); contourSpec.setShapefileType(ShapefileType.PolyLine); contourSpec.addIntegerField("feature_id", 8); contourSpec.addIntegerField("cntr_idx", 4); - contourSpec.addFloatingPointField("depth", 8, 2, false); + contourSpec.addFloatingPointField("depth", 9, 3, false); if (produceElevations) { - contourSpec.addFloatingPointField("elevation", 8, 2, false); + contourSpec.addFloatingPointField("elevation", 9, 3, false); } contourSpec.addFloatingPointField("Shape_len", 13, 6, true); + contourSpec.addIntegerField("shore", 1); contourSpec.setShapefilePrjContent(data.getShapefilePrjContent()); int nContour = 0; - try (ShapefileWriter contourWriter = new ShapefileWriter(shapefileRef, contourSpec);) { + try (ShapefileWriter contourWriter = new ShapefileWriter(contourShapefileRef, contourSpec);) { for (Contour contour : contours) { if (contour.isBoundary()) { continue; } + int shoreCode = 0; + if(shorelineContours.contains(contour)){ + shoreCode = 1; + } nContour++; int cIndex = contour.getLeftIndex(); double z = contour.getZ(); @@ -652,6 +707,7 @@ static void write( contourWriter.setDbfFieldValue("elevation", z); } contourWriter.setDbfFieldValue("Shape_len", dSum); + contourWriter.setDbfFieldValue("shore", shoreCode); contourWriter.writeRecord(record); } @@ -732,4 +788,69 @@ private static String matchCase(String source, String target) { return sb.toString(); } + + + static private void organizeNestedRegions( + List regions, List regionList) + { + // add the regions to the shoreline regions + // remove any perimeter regions and nullify the parent reference + // of their child regions (by removing children) + for (ContourRegion region : regions) { + if (region.getContourRegionType() == ContourRegionType.Perimeter) { + region.removeChildren(); + continue; + } else { + regionList.add(region); + } + } + + int nRegion = regionList.size(); + if (nRegion < 2) { + return; + } + + // The nesting concept organizes the regions to identify which + // regions enclose which. The "parent" of a region is the region + // that immediately encloses it. The parent may, in turn, be + // enclosed by its own parent region. Metaphorically, this concept + // resembles the way traditional Russian nesting dolls are configured. + // To establish the nesting structure, we sort the regions into + // descending order of area. Then, we loop through each region + // comparing it to larger regions to see if the larger region encloses + // it. The "parent" reference may be reset multiple times as smaller + // and smaller enclosing regions are discovered. At the end of the + // process, the parent reference points to the smallest region that + // encloses the region of interest. + Collections.sort(regionList, (ContourRegion o1, ContourRegion o2) + -> Double.compare(o2.getAbsArea(), o1.getAbsArea()) // sort largest to smalles + ); + + // renumber the regions based on their sorted order + int iRegion = 0; + for(ContourRegion region: regionList){ + iRegion++; + region.setApplicationIndex(iRegion); + } + + // Establish the nesting relationships (parent-child) + // for any regions that are not currently connected. + // Any existing relationships will be unchanged. + for (int i = 0; i < nRegion - 1; i++) { + ContourRegion rI = regionList.get(i); + if (rI.hasChildren()) { + continue; + } + double[] xy = rI.getXY(); + for (int j = i + 1; j < nRegion; j++) { + ContourRegion rJ = regionList.get(j); + if (rJ.getParent() == null) { + Point2D testPoint = rJ.getTestPoint(); + if (rI.isPointInsideRegion(xy, testPoint.getX(), testPoint.getY())) { + rI.addChild(rJ); + } + } + } + } + } } diff --git a/svm/src/main/java/org/tinfour/svm/properties/SvmProperties.java b/svm/src/main/java/org/tinfour/svm/properties/SvmProperties.java index c0e3dfec..5d74dba3 100644 --- a/svm/src/main/java/org/tinfour/svm/properties/SvmProperties.java +++ b/svm/src/main/java/org/tinfour/svm/properties/SvmProperties.java @@ -79,6 +79,7 @@ public class SvmProperties { private final static String contourRegionShapefileKey = "contourRegionShapefile"; private final static String contourLineShapefileKey = "contourLineShapefile"; + private final static String contourShapefileShorelineEnabledKey = "contourShapefileShorelineEnabled"; private final static String anomalyTableFileKey = "anomalyTableFileName"; @@ -829,6 +830,18 @@ public File getContourLineShapefile() { return null; } + + /** + * Indicates whether the optional shoreline append operation is enabled. + * @return true if the append operation is enabled; otherwise, false. + */ + public boolean isContourShapefileShorelineEnabled() { + if (properties.containsKey(contourShapefileShorelineEnabledKey)) { + String s = properties.getProperty(contourShapefileShorelineEnabledKey); + return s.toLowerCase().startsWith("t"); + } + return false; + } /** * Get the dimensions for the contour graph image file. * diff --git a/svm/src/main/resources/org/tinfour/svm/SvmTemplate.properties b/svm/src/main/resources/org/tinfour/svm/SvmTemplate.properties index 0974d169..68bec5e8 100644 --- a/svm/src/main/resources/org/tinfour/svm/SvmTemplate.properties +++ b/svm/src/main/resources/org/tinfour/svm/SvmTemplate.properties @@ -60,9 +60,16 @@ capacityGraphTitle = SVM Analysis for Alan Henry Reservoir # polygon (_r for "region") and/or line-based features (_l for "line"). # The .prj file from the source data is transcribed to form prj files # for the output. +# The shoreline-enabled option is used to control whether the shoreline +# is included in the output. In cases where the contour lines are very +# closely spaced, there is a chance that the shoreline geometry +# might conflict (intersect) those generated by the contouring logic. +# So this setting allows users to control whether the shoreline +# is included in the output products. By default, it is disabled. #contourRegionShapefile = contour/SvmContourAlanHenry_r.shp #contourLineShapefile = contour/SvmContourAlanHenry_l.shp +#contourShapefileShorelineEnabled = true