diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..c41cc9e --- /dev/null +++ b/.gitignore @@ -0,0 +1 @@ +/target \ No newline at end of file diff --git a/pom.xml b/pom.xml index a4efda8..119df35 100644 --- a/pom.xml +++ b/pom.xml @@ -103,6 +103,21 @@ gt-geotiff ${geotools.version} + + org.geotools + gt-geojson + ${geotools.version} + + + org.geotools + gt-geometry + 24.0 + + + com.alibaba + fastjson + 1.2.47 + org.gdal @@ -111,6 +126,14 @@ system ${project.basedir}/${gdal.binddir}/gdal-3.8.4.jar + + opencv + opencv-490 + opencv-490 + system + true + ${project.basedir}/src/main/resources/opencv/opencv-490.jar + junit junit diff --git a/src/main/java/org/example/alpha/AlphaShape.java b/src/main/java/org/example/alpha/AlphaShape.java new file mode 100644 index 0000000..a612f7f --- /dev/null +++ b/src/main/java/org/example/alpha/AlphaShape.java @@ -0,0 +1,134 @@ +package org.example.alpha; + +import org.locationtech.jts.geom.Coordinate; +import org.locationtech.jts.geom.GeometryFactory; +import org.locationtech.jts.geom.LineString; +import org.locationtech.jts.geom.Polygon; +import org.locationtech.jts.operation.polygonize.Polygonizer; + +import java.util.ArrayList; +import java.util.Collection; +import java.util.Collections; +import java.util.List; + +/** + * Alpha Shapes提取离散点凹边轮廓线（java实现） + */ +public class AlphaShape { + + private double radius; + + private int sortCount = 0; + //轮廓探寻的圆半径 + public AlphaShape(double radius) { + this.radius = radius; + } + + public List> boundaryPoints(List points) { + + if (points == null || points.size() == 0) { + return null; + } + //编号 + int no = 1; + for (Vector2D point : points) { + point.setId(no); + no++; + } + + List boundary = new ArrayList<>(); + List lines = new ArrayList<>(); + int edgeNo = 1; + for (int i = 0; i < points.size(); i++) { + + // k从i+1开始，减少重复计算 + for (int k = i + 1; k < points.size(); k++) { + // 跳过距离大于直径的情况 + if (points.get(i).distanceToPoint(points.get(k)) > 2 * radius) { + continue; + } + // 两个圆心 + Vector2D c1, c2; + // 线段中点 + Vector2D center = points.get(i).add(points.get(k)).center(); + // 方向向量 d = (x,y) + Vector2D dir = points.get(i).subtract(points.get(k)); + // 垂直向量 n = (a,b) a*dir.x+b*dir.y = 0; a = -(b*dir.y/dir.x) + Vector2D normal = new Vector2D(); + // 因为未知数有两个，随便给y附一个值5。 + normal.setY(5); + + if (dir.getX() != 0) { + normal.setX(-(normal.getY() * dir.getY()) / dir.getX()); + } else { + // 如果方向平行于y轴 + normal.setX(1); + normal.setY(0); + } + // 法向量单位化 + normal.normalize(); + + double len = Math.sqrt(radius * radius - (0.25 * dir.length() * dir.length())); + c1 = center.add(normal.multiply(len)); + c2 = center.subtract(normal.multiply(len)); + + // b1、b2记录是否在圆C1、C2中找到其他点。 + boolean b1 = false, b2 = false; + for (int m = 0; m < points.size(); m++) { + if (m == i || m == k) { + continue; + } + if (b1 != true && points.get(m).distanceToPoint(c1) < radius) { + b1 = true; + } + if (b2 != true && points.get(m).distanceToPoint(c2) < radius) { + b2 = true; + } + // 如果都有内部点，不必再继续检查了 + if (b1 == true && b2 == true) { + break; + } + } + + if (b1 != true || b2 != true) { + Edge edge = new Edge(); + edge.setId(edgeNo); + edge.setA(points.get(i)); + edge.setB(points.get(k)); + edgeNo++; + boundary.add(points.get(i)); + boundary.add(points.get(k)); + lines.add(edge); + } + } + } + + if (lines.size() <= 1) { + return Collections.emptyList(); + } + + List> polygonPoints = new ArrayList<>(); + Polygonizer polygonizer = new Polygonizer(); + List lineStrings = new ArrayList<>(lines.size()); + for (Edge line : lines) { + Coordinate A = new Coordinate(line.getA().getX(), line.getA().getY()); + Coordinate B = new Coordinate(line.getB().getX(), line.getB().getY()); + LineString lineString = new GeometryFactory().createLineString(new Coordinate[]{A,B}); + lineStrings.add(lineString); + } + //线条多边形化 + polygonizer.add(lineStrings); + Collection polygonList=polygonizer.getPolygons(); + for (Object polObj : polygonList) { + Polygon polygon = (Polygon) polObj; + Coordinate[] coordinates=polygon.getCoordinates(); + List vector2DS = new ArrayList<>(); + for (Coordinate coordinate : coordinates) { + Vector2D vector2D = new Vector2D(coordinate.x, coordinate.y); + vector2DS.add(vector2D); + } + polygonPoints.add(vector2DS); + } + return polygonPoints; + } +} \ No newline at end of file diff --git a/src/main/java/org/example/alpha/Edge.java b/src/main/java/org/example/alpha/Edge.java new file mode 100644 index 0000000..7d8c20e --- /dev/null +++ b/src/main/java/org/example/alpha/Edge.java @@ -0,0 +1,83 @@ +package org.example.alpha; + +public class Edge { + + private int id; + + private Vector2D a; + + private Vector2D b; + + private Edge next; + + private boolean hasReversal=false; + + private boolean hasNext=false; + + public Edge() { + + } + + public Edge(Vector2D a, Vector2D b) { + this.a = a; + this.b = b; + } + + public Edge(int id,Vector2D a, Vector2D b) { + this.id = id; + this.a = a; + this.b = b; + } + + public void reversal() { + Vector2D temp = this.a; + this.a = this.b; + this.b = temp; + hasReversal = true; + } + + public Vector2D getA() { + return a; + } + + public void setA(Vector2D a) { + this.a = a; + } + + public Vector2D getB() { + return b; + } + + public void setB(Vector2D b) { + this.b = b; + } + + public Edge getNext() { + return next; + } + + public void setNext(Edge next) { + if (next != null) { + this.next = next; + hasNext = true; + } + + } + + + public boolean hasReversal() { + return this.hasReversal; + } + + public boolean hasNext() { + return this.hasNext; + } + + public int getId() { + return id; + } + + public void setId(int id) { + this.id = id; + } +} \ No newline at end of file diff --git a/src/main/java/org/example/alpha/Vector2D.java b/src/main/java/org/example/alpha/Vector2D.java new file mode 100644 index 0000000..a547aed --- /dev/null +++ b/src/main/java/org/example/alpha/Vector2D.java @@ -0,0 +1,226 @@ +package org.example.alpha; + +import java.util.Objects; + +public class Vector2D { + + private int id; + private double x; + private double y; + private double distance; + private Vector2D next; + + private double bottomHeight; + private double topHeight; + private String index; + + public Vector2D() { + x = 0; + y = 0; + } + + public Vector2D(double x, double y) { + this.x = x; + this.y = y; + } + + //获取弧度 + public double radian() { + return Math.atan2(y, x); + } + + //获取角度 + public double angle() { + return radian() / Math.PI * 180; + } + + public Vector2D clone() { + return new Vector2D(x, y); + } + + public double length() { + return Math.sqrt(getLengthSQ()); + } + + public double getLengthSQ() { + return x * x + y * y; + } + + //向量置零 + public Vector2D Zero() { + x = 0; + y = 0; + return this; + } + + public boolean isZero() { + return x == 0 && y == 0; + } + + //向量的长度设置为我们期待的value + public void setLength(double value) { + double _angle = angle(); + x = Math.cos(_angle) * value; + y = Math.sin(_angle) * value; + } + + //向量的标准化（方向不变，长度为1） + public Vector2D normalize() { + double length = length(); + x = x / length; + y = y / length; + return this; + } + + //是否已经标准化 + public boolean isNormalized() { + return length() == 1.0; + } + + //向量的方向翻转 + public Vector2D reverse() { + x = -x; + y = -y; + return this; + } + + //2个向量的数量积(点积) + public double dotProduct(Vector2D v) { + return x * v.x + y * v.y; + } + + //2个向量的向量积(叉积) + public double crossProduct(Vector2D v) { + return x * v.y - y * v.x; + } + + //计算2个向量的夹角弧度 + //参考点积公式:v1 * v2 = cos * |v1| *|v2| + public static double radianBetween(Vector2D v1, Vector2D v2) { + if (!v1.isNormalized()) v1 = v1.clone().normalize(); // |v1| = 1 + if (!v2.isNormalized()) v2 = v2.clone().normalize(); // |v2| = 1 + return Math.acos(v1.dotProduct(v2)); + } + + //弧度 = 角度乘以PI后再除以180、 推理可得弧度换算角度的公式 + //弧度转角度 + public static double radian2Angle(double radian) { + return radian / Math.PI * 180; + } + + //向量加 + public Vector2D add(Vector2D v) { + return new Vector2D(x + v.x, y + v.y); + } + + //向量减 + public Vector2D subtract(Vector2D v) { + return new Vector2D(x - v.x, y - v.y); + } + + //向量乘 + public Vector2D multiply(double value) { + return new Vector2D(x * value, y * value); + } + + //向量除 + public Vector2D divide(double value) { + return new Vector2D(x / value, y / value); + } + + + public double distanceToPoint(Vector2D point) { + double dx = this.x - point.getX(); + double dy = this.y - point.getY(); + double distance = Math.sqrt(dx * dx + dy * dy); + return Math.abs(distance); + } + + public double getX() { + return x; + } + + public double getY() { + return y; + } + + public void setX(double x) { + this.x = x; + } + + public void setY(double y) { + this.y = y; + } + + public int getId() { + return id; + } + + public void setId(int id) { + this.id = id; + } + + public double getDistance() { + return distance; + } + + public void setDistance(double distance) { + this.distance = distance; + } + + public Vector2D next() { + return next; + } + + public void setNext(Vector2D next) { + this.next = next; + } + + public double getBottomHeight() { + return bottomHeight; + } + + public double getTopHeight() { + return topHeight; + } + + public void setBottomHeight(double bottomHeight) { + this.bottomHeight = bottomHeight; + } + + public void setTopHeight(double topHeight) { + this.topHeight = topHeight; + } + + public String getIndex() { + return index; + } + + public void setIndex(String index) { + this.index = index; + } + + /** + * 线段中点 + * + * @return + */ + public Vector2D center() { + return new Vector2D(getX() * 0.5, getY() * 0.5); + } + + @Override + public boolean equals(Object o) { + if (this == o) return true; + if (o == null || getClass() != o.getClass()) return false; + Vector2D vector2D = (Vector2D) o; + return id == vector2D.id && + Double.compare(vector2D.x, x) == 0 && + Double.compare(vector2D.y, y) == 0; + } + + @Override + public int hashCode() { + return Objects.hash(x, y); + } +} \ No newline at end of file diff --git a/src/main/java/org/example/controller/RingGridController.java b/src/main/java/org/example/controller/RingGridController.java new file mode 100644 index 0000000..39fd43a --- /dev/null +++ b/src/main/java/org/example/controller/RingGridController.java @@ -0,0 +1,76 @@ +package org.example.controller; + +import com.alibaba.fastjson.JSONArray; +import org.example.dto.FillCutDto; +import org.example.service.RingGridService; +import org.example.test.CreateRingGridFix; +import org.example.utils.CommonResult; +import org.geotools.api.geometry.Position; +import org.geotools.api.referencing.crs.CoordinateReferenceSystem; +import org.geotools.coverage.grid.GridCoordinates2D; +import org.geotools.coverage.grid.GridCoverage2D; +import org.geotools.coverage.grid.GridGeometry2D; +import org.geotools.coverage.processing.Operations; +import org.geotools.gce.geotiff.GeoTiffReader; +import org.geotools.geometry.Position2D; +import org.geotools.referencing.CRS; +import org.geotools.util.factory.Hints; +import org.locationtech.jts.geom.Coordinate; +import org.locationtech.jts.geom.GeometryFactory; +import org.locationtech.jts.geom.Point; +import org.springframework.beans.factory.annotation.Autowired; +import org.springframework.web.bind.annotation.Mapping; +import org.springframework.web.bind.annotation.RequestBody; +import org.springframework.web.bind.annotation.RequestMapping; +import org.springframework.web.bind.annotation.RestController; + +import javax.servlet.http.HttpServletRequest; +import java.io.File; +import java.util.ArrayList; +import java.util.List; +import java.util.Map; + +@RestController +public class RingGridController { + @Autowired + public RingGridService ringGridService; + + @RequestMapping("/fillcut") + public CommonResult> calcFillCut(FillCutDto fcd) throws Exception { + + Coordinate[] coordinates = new Coordinate[5]; + double range = 10; // 容差 + double paramRange = fcd.getRange(); + String paramCoordinate = fcd.getCoordinates(); + if(paramRange != 0) range = paramRange; + if(!paramCoordinate.isEmpty()){ + try{ + // 解析字符串 + JSONArray coordinateArray = (JSONArray)JSONArray.parse(paramCoordinate); + coordinates = new Coordinate[coordinateArray.size()+1]; + for(int i=0;i(false,"解析坐标失败！"); + } + }else{ + return new CommonResult<>(false,"缺少必要参数！"); + } + + + return ringGridService.calcFillCut(coordinates,range); + } + +} diff --git a/src/main/java/org/example/dto/FillCutDto.java b/src/main/java/org/example/dto/FillCutDto.java new file mode 100644 index 0000000..3770a22 --- /dev/null +++ b/src/main/java/org/example/dto/FillCutDto.java @@ -0,0 +1,22 @@ +package org.example.dto; + +public class FillCutDto { + String coordinates; + double range; + + public String getCoordinates() { + return coordinates; + } + + public void setCoordinates(String coordinates) { + this.coordinates = coordinates; + } + + public double getRange() { + return range; + } + + public void setRange(double range) { + this.range = range; + } +} diff --git a/src/main/java/org/example/service/RingGridService.java b/src/main/java/org/example/service/RingGridService.java new file mode 100644 index 0000000..eb10605 --- /dev/null +++ b/src/main/java/org/example/service/RingGridService.java @@ -0,0 +1,10 @@ +package org.example.service; + +import org.example.utils.CommonResult; +import org.locationtech.jts.geom.Coordinate; +import org.springframework.stereotype.Service; + + +public interface RingGridService { + CommonResult calcFillCut(Coordinate[] coordinates, double range) throws Exception; +} diff --git a/src/main/java/org/example/service/impl/RingGridServiceImpl.java b/src/main/java/org/example/service/impl/RingGridServiceImpl.java new file mode 100644 index 0000000..d5d02fc --- /dev/null +++ b/src/main/java/org/example/service/impl/RingGridServiceImpl.java @@ -0,0 +1,38 @@ +package org.example.service.impl; + +import org.example.service.RingGridService; +import org.example.utils.CommonResult; + +import org.locationtech.jts.geom.Coordinate; +import org.locationtech.jts.geom.GeometryFactory; +import org.springframework.stereotype.Service; + +import java.util.Map; + +import static org.example.utils.CalcRingGridFix.loadDem; +import static org.example.utils.CalcRingGridFix.getFillCutBlance; +@Service +public class RingGridServiceImpl implements RingGridService { + /** + * 计算填挖方量 + * @param coordinates 四个点坐标 + * @param range 容差 + * @return + * @throws Exception + */ + public CommonResult calcFillCut(Coordinate[] coordinates, double range) throws Exception{ + CommonResult result = new CommonResult(false,"计算出错！"); + try{ + loadDem(); + Map fcbMap = getFillCutBlance(coordinates,range); + result.setSuccess(true); + result.setResult(fcbMap); + result.setMessage(""); + }catch (Exception e){ + System.out.println("计算出错"); + e.printStackTrace(); + } + return result; + } + +} diff --git a/src/main/java/org/example/test/CreateContour.java b/src/main/java/org/example/test/CreateContour.java index 435869f..6b9ef5a 100644 --- a/src/main/java/org/example/test/CreateContour.java +++ b/src/main/java/org/example/test/CreateContour.java @@ -7,16 +7,22 @@ import org.gdal.gdal.gdal; import org.gdal.gdalconst.gdalconst; import org.gdal.ogr.Layer; import org.gdal.osr.SpatialReference; +import org.geotools.api.coverage.grid.GridCoverageWriter; import org.geotools.api.data.SimpleFeatureSource; import org.geotools.api.data.SimpleFeatureStore; import org.geotools.api.data.Transaction; import org.geotools.api.feature.simple.SimpleFeature; import org.geotools.api.feature.simple.SimpleFeatureType; +import org.geotools.api.geometry.Bounds; import org.geotools.api.geometry.Position; +import org.geotools.api.parameter.GeneralParameterValue; +import org.geotools.api.parameter.ParameterValueGroup; import org.geotools.api.referencing.crs.CoordinateReferenceSystem; import org.geotools.coverage.grid.GridCoordinates2D; import org.geotools.coverage.grid.GridCoverage2D; +import org.geotools.coverage.grid.GridCoverageFactory; import org.geotools.coverage.grid.GridGeometry2D; +import org.geotools.coverage.grid.io.AbstractGridFormat; import org.geotools.coverage.processing.Operations; import org.geotools.data.DefaultTransaction; import org.geotools.data.collection.ListFeatureCollection; @@ -24,7 +30,10 @@ import org.geotools.data.shapefile.ShapefileDataStore; import org.geotools.data.simple.SimpleFeatureCollection; import org.geotools.feature.simple.SimpleFeatureBuilder; import org.geotools.feature.simple.SimpleFeatureTypeBuilder; +import org.geotools.gce.geotiff.GeoTiffFormat; import org.geotools.gce.geotiff.GeoTiffReader; +import org.geotools.gce.geotiff.GeoTiffWriteParams; +import org.geotools.geometry.Envelope2DArchived; import org.geotools.geometry.Position2D; import org.geotools.geometry.jts.ReferencedEnvelope; import org.geotools.referencing.CRS; @@ -34,6 +43,10 @@ import org.locationtech.jts.geom.*; import org.geotools.data.shapefile.ShapefileDataStoreFactory; import org.locationtech.jts.operation.buffer.BufferOp; import org.locationtech.jts.operation.buffer.BufferParameters; +import org.opencv.core.Mat; +import org.opencv.core.Size; +import org.opencv.imgcodecs.Imgcodecs; +import org.opencv.imgproc.Imgproc; import org.opengis.geometry.Envelope; import javax.media.jai.RasterFactory; @@ -56,11 +69,16 @@ public class CreateContour { GeometryFactory geometryFactory = new GeometryFactory(); // 定义矩形的左下角和右上角坐标（WGS84） List coords = new ArrayList(){{ - add(new double[]{113.52896835361318, 23.654018962759377}); - add(new double[]{113.529270627066524, 23.654027766840539 }); - add(new double[]{ 113.529294104616298, 23.653772448486752}); - add(new double[]{ 113.528977157694356, 23.653763644405586 }); - add(new double[]{113.52896835361318, 23.654018962759377}); +// add(new double[]{113.52896835361318, 23.654018962759377}); +// add(new double[]{113.529270627066524, 23.654027766840539 }); +// add(new double[]{ 113.529294104616298, 23.653772448486752}); +// add(new double[]{ 113.528977157694356, 23.653763644405586 }); +// add(new double[]{113.52896835361318, 23.654018962759377}); + add(new double[]{113.529117, 23.65403}); + add(new double[]{113.529297, 23.653906 }); + add(new double[]{ 113.529189, 23.653765}); + add(new double[]{ 113.528963, 23.65389 }); + add(new double[]{113.529117, 23.65403}); }}; // 创建矩形的四个角点坐标 Coordinate[] coordinates = new Coordinate[5]; @@ -138,7 +156,7 @@ public class CreateContour { featuresList.add(simpleFeature); } // 输出等高线shp - String outShpPath = "E:\\code\\CutFillDemo\\target\\classes\\tmp.shp"; + String outShpPath = "E:\\code\\CutFillDemo\\target\\classes\\contour_tmp.shp"; SimpleFeatureCollection collection = new ListFeatureCollection(TYPE, featuresList); ShapefileDataStore newDataStore = getShapefileDataStore(collection,outShpPath); newDataStore.createSchema(TYPE); @@ -165,147 +183,147 @@ public class CreateContour { System.exit(1); } // 调用exe生成 dem文件（ps: 调用可执行文件和参数拼接必须要用空格隔开） - String outDemPath = "E:\\code\\CutFillDemo\\target\\classes\\tmp_dem.tiff"; + String outDemPath = "E:\\code\\CutFillDemo\\target\\classes\\dem_tmp.tiff"; String paras = " -i "+outShpPath+" -o "+outDemPath+" -f ELEV -s 0.0000625"; String cmd = "E:\\code\\contour2dem\\bin\\Debug\\contour2dem.exe" + paras; ExeMain.openExe(cmd); // 读取生成的dem File file = new File(outDemPath); - // Hints.FORCE_LONGITUDE_FIRST_AXIS_ORDER：设置经度为第一轴顺序 GeoTiffReader reader = new GeoTiffReader(file); GridCoverage2D targetCoverage = reader.read(null); // 逐个像素比较 RenderedImage targetImage = targetCoverage.getRenderedImage(); - RenderedImage image = CreateContour.coverage.getRenderedImage(); +// RenderedImage image = CreateContour.coverage.getRenderedImage(); Raster targetRaster = targetImage.getData(); GridGeometry2D targetGridGeometry = targetCoverage.getGridGeometry(); GridGeometry2D gridGeometry = CreateContour.coverage.getGridGeometry(); - Raster data = image.getData(); +// Raster data = image.getData(); + + Bounds outBound = targetCoverage.getEnvelope(); + WritableRaster raster = RasterFactory.createBandedRaster(2, targetRaster.getWidth(), targetRaster.getHeight(), 1, (java.awt.Point) null); + for (int y = 0; y < targetRaster.getHeight(); y++) { for (int x = 0; x < targetRaster.getWidth(); x++) { double elevation1 = targetRaster.getSampleDouble(x, y, 0); // 根据像素坐标获取地理坐标 Position position = targetGridGeometry.gridToWorld(new GridCoordinates2D(x, y)); + // 将地理坐标转换成【数据源】的网格坐标 GridCoordinates2D gridPoint = gridGeometry.worldToGrid(position); double[] elevation = new double[1]; CreateContour.coverage.evaluate(gridPoint,elevation); double diff = elevation1 - elevation[0]; - System.out.println(x + "," + y + "," + 0 + "," + diff); // diffRaster.setSample(x, y, 0, diff); - } - } - - // 寻找零线 - List zeroFeatureList = new ArrayList<>(); - SimpleFeatureType zeroPointType = createPointFeatureType(); - SimpleFeatureBuilder zeroFeatureBuilder = new SimpleFeatureBuilder(zeroPointType); - for (int y = 0; y < targetRaster.getHeight(); y++) { - for (int x = 0; x < targetRaster.getWidth(); x++) { - double elevation1 = targetRaster.getSampleDouble(x, y, 0); - // 根据像素坐标获取地理坐标 - Position2D position = (Position2D)targetGridGeometry.gridToWorld(new GridCoordinates2D(x, y)); - // 将地理坐标转换成【数据源】的网格坐标 - GridCoordinates2D gridPoint = gridGeometry.worldToGrid(position); - double[] elevation = new double[1]; - CreateContour.coverage.evaluate(gridPoint,elevation); - double diff = elevation1 - elevation[0]; - if(diff<0.0001&&diff>=0){ - zeroFeatureBuilder.add(geometryFactory.createPoint(new Coordinate(position.x,position.y))); - zeroFeatureBuilder.add(("零线点")); - zeroFeatureBuilder.add(elevation[0]); - SimpleFeature simpleFeature = zeroFeatureBuilder.buildFeature("zero"+UUID.randomUUID()); - zeroFeatureList.add(simpleFeature); + double[] coord = position.getCoordinate(); + Point point = geometryFactory.createPoint(new Coordinate(coord[0], coord[1])); + if(rectangle.intersects(point)){ + raster.setSample(x,y,0,0); + }else{ + raster.setSample(x,y,0,diff); } - System.out.println(x + "," + y + "," + 0 + "," + diff); } } - String outZeroShpPath = "E:\\code\\CutFillDemo\\target\\classes\\zeroPoint.shp"; - SimpleFeatureCollection zeroCollection = new ListFeatureCollection(zeroPointType, zeroFeatureList); - ShapefileDataStore zeroDataStore = getShapefileDataStore(zeroCollection,outZeroShpPath); - zeroDataStore.createSchema(zeroPointType); - zeroDataStore.setCharset(Charset.forName("UTF-8")); - Transaction zeroTransaction = new DefaultTransaction("create"); - String zeroTypeName = zeroDataStore.getTypeNames()[0]; - SimpleFeatureSource zeroFeatureSource = zeroDataStore.getFeatureSource(zeroTypeName); + // 输出一个tiff + GridCoverageFactory f = new GridCoverageFactory(); + GridCoverage2D coverage = f.create("test1", raster, outBound); + GeoTiffFormat format = new GeoTiffFormat(); + String newDemPath = "E:\\code\\CutFillDemo\\target\\classes\\out_dem.tiff"; + File out = new File(newDemPath); + GridCoverageWriter writer = format.getWriter(out); + writer.write(coverage, createWriteParameters()); + writer.dispose(); + coverage.dispose(true); - if (zeroFeatureSource instanceof SimpleFeatureStore) { - SimpleFeatureStore featureStore = (SimpleFeatureStore) zeroFeatureSource; - featureStore.setTransaction(zeroTransaction); - try { - featureStore.addFeatures(zeroCollection); - zeroTransaction.commit(); - } catch (Exception problem) { - problem.printStackTrace(); - zeroTransaction.rollback(); - } finally { - zeroTransaction.close(); - } - } else { - System.out.println(typeName + " does not support read/write access"); - System.exit(1); - } + // 使用opencv输出边界 + // 加载lib,这个lib的名称 + String path = "D:\\tools\\opencv\\opencv\\build\\java\\x64\\" + System.mapLibraryName("opencv_java490"); + String outPath = "E:\\code\\CutFillDemo\\target\\classes\\out_boundary.png"; + File lib = new File(path); + System.load(lib.getAbsolutePath()); + //读取图片信息 + Mat image = Imgcodecs.imread(new File(newDemPath).getAbsolutePath(),1); + testCanny(image); + Imgcodecs.imwrite(outPath, image); - // 将差异栅格转换为 GridCoverage2D（如果需要） -// GridCoverageFactory gcf = new GridCoverageFactory(); -// GridCoverage2D diffCoverage = gcf.create("Difference", diffRaster, demCoverage1.getEnvelope2D()); +// // 寻找零线 +// List zeroFeatureList = new ArrayList<>(); +// SimpleFeatureType zeroPointType = createPointFeatureType(); +// SimpleFeatureBuilder zeroFeatureBuilder = new SimpleFeatureBuilder(zeroPointType); +// for (int y = 0; y < targetRaster.getHeight(); y++) { +// for (int x = 0; x < targetRaster.getWidth(); x++) { +// double elevation1 = targetRaster.getSampleDouble(x, y, 0); +// // 根据像素坐标获取地理坐标 +// Position2D position = (Position2D)targetGridGeometry.gridToWorld(new GridCoordinates2D(x, y)); +// // 将地理坐标转换成【数据源】的网格坐标 +// GridCoordinates2D gridPoint = gridGeometry.worldToGrid(position); +// double[] elevation = new double[1]; +// CreateContour.coverage.evaluate(gridPoint,elevation); +// double diff = elevation1 - elevation[0]; +// if(diff<0.0001&&diff>=0){ +// zeroFeatureBuilder.add(geometryFactory.createPoint(new Coordinate(position.x,position.y))); +// zeroFeatureBuilder.add(("零线点")); +// zeroFeatureBuilder.add(elevation[0]); +// SimpleFeature simpleFeature = zeroFeatureBuilder.buildFeature("zero"+UUID.randomUUID()); +// zeroFeatureList.add(simpleFeature); +// } +// System.out.println(x + "," + y + "," + 0 + "," + diff); +// } +// } - -// // 注册所有gdal的驱动 -// gdal.AllRegister(); -// Dataset vector = gdal.OpenEx("E:/code/CutFillDemo/target/classes/tmp.shp",gdalconst.OF_VECTOR,null,null,null); -// int layerCount = vector.GetLayerCount(); -// if(layerCount == 0){ -// System.err.println("矢量数据没有图层！"); +// String outZeroShpPath = "E:\\code\\CutFillDemo\\target\\classes\\zeroPoint.shp"; +// SimpleFeatureCollection zeroCollection = new ListFeatureCollection(zeroPointType, zeroFeatureList); +// ShapefileDataStore zeroDataStore = getShapefileDataStore(zeroCollection,outZeroShpPath); +// zeroDataStore.createSchema(zeroPointType); +// zeroDataStore.setCharset(Charset.forName("UTF-8")); +// Transaction zeroTransaction = new DefaultTransaction("create"); +// String zeroTypeName = zeroDataStore.getTypeNames()[0]; +// SimpleFeatureSource zeroFeatureSource = zeroDataStore.getFeatureSource(zeroTypeName); +// +// if (zeroFeatureSource instanceof SimpleFeatureStore) { +// SimpleFeatureStore featureStore = (SimpleFeatureStore) zeroFeatureSource; +// featureStore.setTransaction(zeroTransaction); +// try { +// featureStore.addFeatures(zeroCollection); +// zeroTransaction.commit(); +// } catch (Exception problem) { +// problem.printStackTrace(); +// zeroTransaction.rollback(); +// } finally { +// zeroTransaction.close(); +// } +// } else { +// System.out.println(typeName + " does not support read/write access"); // System.exit(1); // } -// Layer layer = vector.GetLayerByIndex(0); -// SpatialReference crs = layer.GetSpatialRef(); -// double[] boundsArr = layer.GetExtent(); -// double pixelSize = 0.00025; -// int NoData_value = -9999; -// int x_res = (int)((boundsArr[1] - boundsArr[0]) / pixelSize); -// int y_res = (int)((boundsArr[3] - boundsArr[2]) / pixelSize); -// if(x_res == 0){ -// x_res = 2000; -// } -// if(y_res == 0){ -// y_res = 2000; -// } -// // 创建dem输出的数据集 -// String demOutPath = "E:\\data\\demOutFile.tiff"; -// Dataset target_ds = gdal.GetDriverByName("GTiff").Create(demOutPath,x_res,y_res,1, gdalconst.GDT_Byte); -// target_ds.SetGeoTransform(new double[]{ boundsArr[0], pixelSize,0,boundsArr[3],0,-pixelSize }); -// target_ds.SetProjection(crs.ExportToWkt()); -// -// -// Band band = target_ds.GetRasterBand(1); -// band.SetNoDataValue(NoData_value); -//// band.FlushCache(); -// Vector objects = new Vector<>(); -// objects.add("ATTRIBUTE=ELEV"); -// // 将矢量数据栅格化到栅格数据集中 -// double[] burnValues = new double[1]; // 烧录值数组 -// burnValues[0] = 255; // 例如，将矢量要素烧录为255值 -// String[] options = new String[2]; -// options[0] = "ALL_TOUCHED=TRUE"; // 设置选项，例如考虑所有相交的像素 -// options[1] = "BURN_VALUES=" + burnValues[0]; // 设置烧录值 -// -//// gdal.RasterizeLayer(target_ds, new int[]{1}, layer, null, options); -// -// gdal.RasterizeLayer(target_ds,new int[]{1},layer); -// target_ds.FlushCache(); // 确保所有数据都已写入文件 -//// gdal. -// // 释放资源 -// vector.delete(); -// target_ds.delete(); return 0; } - + // canny算子 + public static void testCanny(Mat image) { + //高斯滤波 + Imgproc.GaussianBlur(image, image, new Size(3, 3), 0, 0); + //灰度 + Imgproc.cvtColor(image, image,Imgproc.COLOR_BGR2GRAY); + //一般来说，高阈值maxVal推荐是低阈值minVal的2~3倍 + int lowThresh=45; + //边缘检测 + Imgproc.Canny(image, image,lowThresh, lowThresh*3,3); +// HighGui.imshow("Canny", image); +// HighGui.waitKey(0); + } + public static GeneralParameterValue[] createWriteParameters() { + GeoTiffFormat format = new GeoTiffFormat(); + GeoTiffWriteParams wp = new GeoTiffWriteParams(); + wp.setCompressionMode(GeoTiffWriteParams.MODE_EXPLICIT); + wp.setCompressionType("Deflate"); + ParameterValueGroup params = format.getWriteParameters(); + params.parameter(AbstractGridFormat.GEOTOOLS_WRITE_PARAMS.getName().toString()).setValue(wp); + GeneralParameterValue[] writeParameters = params.values().toArray(new GeneralParameterValue[1]); + return writeParameters; + } private static ShapefileDataStore getShapefileDataStore(SimpleFeatureCollection collection,String outShpPath) throws IOException { // 创建输出路径 File newFile = new File(outShpPath); diff --git a/src/main/java/org/example/test/CreateRingGrid.java b/src/main/java/org/example/test/CreateRingGrid.java index 2b838be..447e485 100644 --- a/src/main/java/org/example/test/CreateRingGrid.java +++ b/src/main/java/org/example/test/CreateRingGrid.java @@ -64,7 +64,8 @@ public class CreateRingGrid { builder.add("the_geom", Point.class); // builder.length(15).add("name", String.class); // <- 15 chars width for name field builder.add("id", String.class); - builder.add("height", double.class); + builder.add("bottomHeight", double.class); + builder.add("topHeight", double.class); return builder.buildFeatureType(); } public static void main(String[] args) throws Exception{ @@ -73,6 +74,7 @@ public class CreateRingGrid { CoordinateReferenceSystem crs = CRS.decode("EPSG:4326"); GeometryFactory geometryFactory = new GeometryFactory(); // 定义矩形的左下角和右上角坐标（WGS84） + // 输入参数：四点坐标，最大挖填方允许差值 List coords = new ArrayList(){{ add(new double[]{113.52896835361318, 23.654018962759377}); add(new double[]{113.529270627066524, 23.654027766840539 }); @@ -100,26 +102,29 @@ public class CreateRingGrid { System.out.println("设计标高："+targetHeight); double allSum = calcBlanceHeight(coordinates,geometryFactory,targetHeight); - double[] stepSet = new double[]{5,3,1,0.5,0.25,0.125,0.0625,0.03125,0.0156125,0.007806125,0.007806125/2,0.007806125/4,0.007806125/8,0.007806125/16,0.007806125/32}; - PrimitiveIterator.OfDouble iterator = Arrays.stream(stepSet).iterator(); - double nextHeight=targetHeight,step = iterator.next(); - while (Math.abs(allSum)-10>0){ - nextHeight = allSum<0?nextHeight-step:nextHeight+step; - System.out.println("nextHeight:"+nextHeight); - double firstSum = calcBlanceHeight(coordinates,geometryFactory,nextHeight); - if(Math.abs(firstSum)-1<0){ - break; - } - nextHeight = firstSum<0?nextHeight-step:nextHeight+step; - System.out.println("nextHeight:"+nextHeight); - System.out.println("firstSum:"+firstSum); - allSum = calcBlanceHeight(coordinates,geometryFactory,nextHeight); - if((firstSum>0&&allSum<0)||(firstSum<0&&allSum>0)){ - step = iterator.next(); - } - System.out.println("allSum:"+allSum); - } - System.out.println("lastHeight:"+nextHeight); +// double[] stepSet = new double[]{5,3,1,0.5,0.25,0.125,0.0625,0.03125,0.0156125,0.007806125,0.007806125/2,0.007806125/4,0.007806125/8,0.007806125/16,0.007806125/32}; +// PrimitiveIterator.OfDouble iterator = Arrays.stream(stepSet).iterator(); +// double nextHeight=targetHeight,step = iterator.next(); +// while (Math.abs(allSum)-10>0){ +// nextHeight = allSum<0?nextHeight-step:nextHeight+step; +// System.out.println("nextHeight:"+nextHeight); +// double firstSum = calcBlanceHeight(coordinates,geometryFactory,nextHeight); +// if(Math.abs(firstSum)-1<0){ +// break; +// } +// nextHeight = firstSum<0?nextHeight-step:nextHeight+step; +// System.out.println("nextHeight:"+nextHeight); +// System.out.println("firstSum:"+firstSum); +// allSum = calcBlanceHeight(coordinates,geometryFactory,nextHeight); +// if((firstSum>0&&allSum<0)||(firstSum<0&&allSum>0)){ +// // fix: 采用线性方式去计算步长 +// step = iterator.next(); +// } +// System.out.println("allSum:"+allSum); +// } +// System.out.println("lastHeight:"+nextHeight); + + // 输出：理想标高、总挖方量、总填方量、填挖方差值、边界数据 } public static double calcBlanceHeight(Coordinate[] coordinates,GeometryFactory geometryFactory,double targetHeight) throws Exception{ // 创建GeometryFactory @@ -184,49 +189,49 @@ public class CreateRingGrid { // System.out.println("targetSum:"+targetSum); double allsum = (double)bottomResult.get("sum")+(double)topResult.get("sum")+targetSum; -// // Create a map content and add our shapefile to it -// MapContent map = new MapContent(); -// map.setTitle("填挖方计算-创建网格"); -// // 构建添加面图层 -// final SimpleFeatureType TYPE =createFeatureType(); -// DefaultFeatureCollection featureCollection = new DefaultFeatureCollection(); -// SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE); -// featuresProperties.forEach(props->{ -// featureBuilder.add(props.get("geom")); -// featureBuilder.add(props.get("name")); -// featureBuilder.add(props.get("id")); -// SimpleFeature feature = featureBuilder.buildFeature((String.valueOf(props.get("id")))); -// featureCollection.add(feature); -// }); -// Style style = SLD.createSimpleStyle(featureCollection.getSchema()); -// Layer layer = new FeatureLayer(featureCollection, style); -// map.addLayer(layer); -// -// // 构建点要素图层 -// List>> pointList = new ArrayList(){{ -//// add(prf); + // Create a map content and add our shapefile to it + MapContent map = new MapContent(); + map.setTitle("填挖方计算-创建网格"); + // 构建添加面图层 + final SimpleFeatureType TYPE =createFeatureType(); + DefaultFeatureCollection featureCollection = new DefaultFeatureCollection(); + SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE); + featuresProperties.forEach(props->{ + featureBuilder.add(props.get("geom")); + featureBuilder.add(props.get("name")); + featureBuilder.add(props.get("id")); + SimpleFeature feature = featureBuilder.buildFeature((String.valueOf(props.get("id")))); + featureCollection.add(feature); + }); + Style style = SLD.createSimpleStyle(featureCollection.getSchema()); + Layer layer = new FeatureLayer(featureCollection, style); + map.addLayer(layer); + + // 构建点要素图层 + List>> pointList = new ArrayList(){{ + add(prf); // add(topResult.get("filterGrid")); -//// add(tPointFeaturesProperties); -// }}; -//// List>> pointList = pointTopRingFeaturesProperties; -// pointList.stream().forEach(list->{ -// if(!list.isEmpty()){ -// final SimpleFeatureType PTYPE =createPointFeatureType(); -// DefaultFeatureCollection pointFeatureCollection = new DefaultFeatureCollection(); -// SimpleFeatureBuilder pointFeatureBuilder = new SimpleFeatureBuilder(PTYPE); -// list.forEach(props->{ -// pointFeatureBuilder.add(props.get("geom")); -// pointFeatureBuilder.add(props.get("id")); -// pointFeatureBuilder.add(props.get("height")); -// SimpleFeature feature = pointFeatureBuilder.buildFeature((String.valueOf(props.get("id")))); -// pointFeatureCollection.add(feature); -// }); -// Style pointStyle = SLD.createSimpleStyle(pointFeatureCollection.getSchema()); -// Layer pointLlayer = new FeatureLayer(pointFeatureCollection, pointStyle); -// map.addLayer(pointLlayer); -// } -// }); -// +// add(tPointFeaturesProperties); + }}; +// List>> pointList = pointTopRingFeaturesProperties; + pointList.stream().forEach(list->{ + if(!list.isEmpty()){ + final SimpleFeatureType PTYPE =createPointFeatureType(); + DefaultFeatureCollection pointFeatureCollection = new DefaultFeatureCollection(); + SimpleFeatureBuilder pointFeatureBuilder = new SimpleFeatureBuilder(PTYPE); + list.forEach(props->{ + pointFeatureBuilder.add(props.get("geom")); + pointFeatureBuilder.add(props.get("id")); + pointFeatureBuilder.add(props.get("height")); + SimpleFeature feature = pointFeatureBuilder.buildFeature((String.valueOf(props.get("id")))); + pointFeatureCollection.add(feature); + }); + Style pointStyle = SLD.createSimpleStyle(pointFeatureCollection.getSchema()); + Layer pointLlayer = new FeatureLayer(pointFeatureCollection, pointStyle); + map.addLayer(pointLlayer); + } + }); + // // 添加网格图层 // // 创建grids点， //// ListFeatureCollection collection = new ListFeatureCollection(TYPE); diff --git a/src/main/java/org/example/test/CreateRingGridFix.java b/src/main/java/org/example/test/CreateRingGridFix.java new file mode 100644 index 0000000..99c3829 --- /dev/null +++ b/src/main/java/org/example/test/CreateRingGridFix.java @@ -0,0 +1,811 @@ +package org.example.test; + +import org.apache.commons.lang3.ArrayUtils; +import org.example.alpha.AlphaShape; +import org.example.alpha.Vector2D; +import org.geotools.api.feature.simple.SimpleFeature; +import org.geotools.api.feature.simple.SimpleFeatureType; +import org.geotools.api.geometry.Position; +import org.geotools.api.referencing.crs.CoordinateReferenceSystem; +import org.geotools.api.style.Style; +import org.geotools.coverage.grid.GridCoordinates2D; +import org.geotools.coverage.grid.GridCoverage2D; +import org.geotools.coverage.grid.GridGeometry2D; +import org.geotools.coverage.processing.Operations; +import org.geotools.feature.DefaultFeatureCollection; +import org.geotools.feature.simple.SimpleFeatureBuilder; +import org.geotools.feature.simple.SimpleFeatureTypeBuilder; +import org.geotools.gce.geotiff.GeoTiffReader; +import org.geotools.geojson.feature.FeatureJSON; +import org.geotools.geojson.geom.GeometryJSON; +import org.geotools.geometry.Position2D; +import org.geotools.map.FeatureLayer; +import org.geotools.map.Layer; +import org.geotools.map.MapContent; +import org.geotools.referencing.CRS; +import org.geotools.referencing.crs.DefaultGeographicCRS; +import org.geotools.styling.SLD; +import org.geotools.swing.JMapFrame; +import org.geotools.util.factory.Hints; +import org.locationtech.jts.geom.*; +import org.locationtech.jts.geom.Point; +import org.locationtech.jts.geom.Polygon; +import org.locationtech.jts.operation.buffer.BufferOp; +import org.locationtech.jts.operation.buffer.BufferParameters; + +import java.awt.*; +import java.io.File; +import java.util.*; +import java.util.List; +import java.util.concurrent.atomic.AtomicReference; + +public class CreateRingGridFix { + public static GridCoverage2D coverage; + public static double extendLen = 100; + public static GeometryFactory geometryFactory = new GeometryFactory(); + public static List> fillPoints; + public static List> cutPoints; + public static double cutAllSum = 0; + public static double fillAllSum = 0; + public static double meter = 2; // 间距,单位为米。 + private static SimpleFeatureType createFeatureType() { + + SimpleFeatureTypeBuilder builder = new SimpleFeatureTypeBuilder(); + builder.setName("Location"); + builder.setCRS(DefaultGeographicCRS.WGS84); // <- Coordinate reference system + + // add attributes in order + builder.add("the_geom", Polygon.class); + builder.length(15).add("name", String.class); // <- 15 chars width for name field + builder.add("id", Integer.class); + + // build the type + final SimpleFeatureType LOCATION = builder.buildFeatureType(); + + return LOCATION; + } + private static SimpleFeatureType createLineStringFeatureType() { + + SimpleFeatureTypeBuilder builder = new SimpleFeatureTypeBuilder(); + builder.setName("Location"); + builder.setCRS(DefaultGeographicCRS.WGS84); // <- Coordinate reference system + + // add attributes in order + builder.add("the_geom", LineString.class); + builder.length(15).add("name", String.class); // <- 15 chars width for name field + builder.add("id", Integer.class); + + // build the type + final SimpleFeatureType LOCATION = builder.buildFeatureType(); + + return LOCATION; + } + private static SimpleFeatureType createPointFeatureType() { + + SimpleFeatureTypeBuilder builder = new SimpleFeatureTypeBuilder(); + builder.setName("Location"); + builder.setCRS(DefaultGeographicCRS.WGS84); // <- Coordinate reference system + + // add attributes in order + builder.add("the_geom", Point.class); +// builder.length(15).add("name", String.class); // <- 15 chars width for name field + builder.add("id", String.class); + builder.add("bottomHeight", double.class); + builder.add("topHeight", double.class); + return builder.buildFeatureType(); + } + public static void main(String[] args) throws Exception{ + loadDem(); + // 假设我们使用WGS84坐标系统 + CoordinateReferenceSystem crs = CRS.decode("EPSG:4326"); + GeometryFactory geometryFactory = new GeometryFactory(); + // 定义矩形的左下角和右上角坐标（WGS84） + // 输入参数：四点坐标，最大挖填方允许差值 + List coords = new ArrayList(){{ + add(new double[]{113.52896835361318, 23.654018962759377}); + add(new double[]{113.529270627066524, 23.654027766840539 }); + add(new double[]{ 113.529294104616298, 23.653772448486752}); + add(new double[]{ 113.528977157694356, 23.653763644405586 }); + add(new double[]{113.52896835361318, 23.654018962759377}); + }}; + // 创建矩形的四个角点坐标 + Coordinate[] coordinates = new Coordinate[5]; + coordinates[0] = new Coordinate(coords.get(0)[0], coords.get(0)[1]); + coordinates[1] = new Coordinate(coords.get(1)[0], coords.get(1)[1]); + coordinates[2] = new Coordinate(coords.get(2)[0], coords.get(2)[1]); + coordinates[3] = new Coordinate(coords.get(3)[0], coords.get(3)[1]); + // 闭合多边形，复制第一个坐标作为最后一个坐标 + coordinates[4] = new Coordinate(coords.get(0)[0], coords.get(0)[1]); + + // 计算设计标高，取四个角的高程平均值 +// double centerHeight = getElevation(rectangle.getCentroid()); + double corner1Height = getElevation( geometryFactory.createPoint(coordinates[0])); + double corner2Height = getElevation( geometryFactory.createPoint(coordinates[1])); + double corner3Height = getElevation( geometryFactory.createPoint(coordinates[2])); + double corner4Height = getElevation( geometryFactory.createPoint(coordinates[3])); + double targetHeight = (corner1Height+corner2Height+corner3Height+corner4Height)/4; + + System.out.println("初始标高："+targetHeight); + double startTime = System.currentTimeMillis(); + double allSum = calcBlanceHeight(coordinates,geometryFactory,targetHeight); + double step = 5;//起始步长 +// double[] stepSet = new double[500]; +// for(int i=0;i<500;i++){ +// stepSet[i] = step/Math.pow(2,i); +// } +// PrimitiveIterator.OfDouble iterator = Arrays.stream(stepSet).iterator(); +// double stepIterator = iterator.next(); +// double range = 10; // 容差 +// double nextHeight=targetHeight; +// while (Math.abs(allSum)-range>0){ +// System.out.println("allSum:"+allSum); +// nextHeight = allSum<0?nextHeight-step:nextHeight+step; +// System.out.println("nextHeight:"+nextHeight); +// double firstSum = calcBlanceHeight(coordinates,geometryFactory,nextHeight); +// System.out.println("firstSum:"+firstSum); +// if(Math.abs(firstSum)-range<0){ +// allSum = firstSum; +// break; +// } +// if((firstSum>0&&allSum<0)||(firstSum<0&&allSum>0)){ +// // fix: 采用线性方式去计算步长 +// step = step*Math.min(Math.abs(allSum),Math.abs(firstSum))/Math.abs(allSum-firstSum); +//// step = iterator.next(); +// } +// double tmpHeight = firstSum<0?nextHeight-step:nextHeight+step; +// System.out.println("nextHeight:"+tmpHeight); +// allSum = calcBlanceHeight(coordinates,geometryFactory,tmpHeight); +// if(Math.abs(tmpHeight-nextHeight)<0.0001){ +// nextHeight = tmpHeight; +// break; +// } +// nextHeight = tmpHeight; +// } +// System.out.println("理想标高:"+nextHeight); +// System.out.println("总挖方量："+cutAllSum); +// System.out.println("总填方量："+fillAllSum); +// System.out.println("填挖方差值："+allSum); +// +// System.out.println("耗时："+(System.currentTimeMillis()-startTime)); +// +// double radius = meter2degree(meter); +// SimpleFeature cutBoundary = searchBourdary(cutPoints,radius,"挖方边界线"); +// SimpleFeature fillBoundary = searchBourdary(fillPoints,radius,"填方边界线"); +// FeatureJSON featureJSON = new FeatureJSON(new GeometryJSON(8)); +// System.out.println("挖方边界线:"+featureJSON.toString(cutBoundary)); +// System.out.println("填方边界线:"+featureJSON.toString(fillBoundary)); + // 构建添加面图层 +// final SimpleFeatureType TYPE =createFeatureType(); +// DefaultFeatureCollection featureCollection = new DefaultFeatureCollection(); +// SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE); +// featuresProperties.forEach(props->{ +// featureBuilder.add(props.get("geom")); +// featureBuilder.add(props.get("name")); +// featureBuilder.add(props.get("id")); +// SimpleFeature feature = featureBuilder.buildFeature((String.valueOf(props.get("id")))); +// featureCollection.add(feature); +// }); + // 输出：理想标高、总挖方量、总填方量、填挖方差值、边界数据 + } + public static double calcBlanceHeight(Coordinate[] coordinates,GeometryFactory geometryFactory,double targetHeight) throws Exception{ + // 创建GeometryFactory + List> featuresProperties = new ArrayList<>();// 面要素属性集合 + List> pointFeaturesProperties = new ArrayList<>();// 扩展点属性集合 + List> tPointFeaturesProperties = new ArrayList<>();// 目标平面的点集属性集合 + + // 清空 + fillPoints = new ArrayList<>(); + cutPoints = new ArrayList<>(); + cutAllSum = 0; + fillAllSum = 0; + + // 间距,单位为米。 +// double meter = 2; + double bottomRatio = 1.0/2; + double topRatio = 1.0/1.5; + // 创建多边形（矩形） + Map innerPolygon = new HashMap<>(); + Polygon rectangle = geometryFactory.createPolygon(coordinates); + innerPolygon.put("geom",rectangle); + innerPolygon.put("name","输入的面"); + innerPolygon.put("id",1); + featuresProperties.add(innerPolygon); + + // 外接矩形打点 + Map outPolygon1 = new HashMap<>(); + Polygon out1 = (Polygon) rectangle.getEnvelope();//geometryFactory.createPolygon(cords4); + outPolygon1.put("geom",out1); + outPolygon1.put("name","输入面的外接矩形"); + outPolygon1.put("id",6); +// featuresProperties.add(outPolygon1); + Coordinate[] outCoords = out1.getCoordinates(); + Coordinate[] outLine = new Coordinate[]{outCoords[0],outCoords[1]}; + addPartPoint(outLine,geometryFactory,meter2degree(meter),rectangle,tPointFeaturesProperties,targetHeight,"0"); + + // 扩展面打点 + List>> pointRingFeaturesProperties = createRingPoint(rectangle,targetHeight,bottomRatio,topRatio,meter,geometryFactory); + + // 中心点 + Point centroid = rectangle.getCentroid(); + // 计算网格 + Map result = calcGrid(pointRingFeaturesProperties,centroid,meter); + List> fillPointList = (List>)result.get("fillPointList"); + List> cutPointList = (List>)result.get("cutPointList"); + double fillSum = (double)result.get("fillSum"); + double cutSum = (double)result.get("cutSum"); + + // 计算设计平面内的填挖方量 + for(int i=0;i pfp = tPointFeaturesProperties.get(i); + double height = (double)pfp.get("height"); + Point p = (Point)pfp.get("geom"); + double elev = getElevation(p); + if(elev>height){ + cutPointList.add(pfp); + cutSum += (elev - height)*Math.pow(meter,2); + }else{ + fillPointList.add(pfp); + fillSum += (elev - height)*Math.pow(meter,2); + } + } + + List> boundary = searchBourdary(cutPointList, meter2degree(meter)); + + // Create a map content and add our shapefile to it + MapContent map = new MapContent(); + map.setTitle("填挖方计算-创建网格"); + // 构建添加面图层 + final SimpleFeatureType TYPE =createFeatureType(); + DefaultFeatureCollection featureCollection = new DefaultFeatureCollection(); + SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE); + featuresProperties.forEach(props->{ + featureBuilder.add(props.get("geom")); + featureBuilder.add(props.get("name")); + featureBuilder.add(props.get("id")); + SimpleFeature feature = featureBuilder.buildFeature((String.valueOf(props.get("id")))); + featureCollection.add(feature); + }); + Style style = SLD.createSimpleStyle(featureCollection.getSchema()); + Layer layer = new FeatureLayer(featureCollection, style); + map.addLayer(layer); + + // 构建点要素图层 +// List>> pointList = new ArrayList(){{ +// add(result.get("fillPointList")); +// add(tPointFeaturesProperties); +// }}; + List>> pointList = pointRingFeaturesProperties; + pointList.stream().forEach(list->{ + if(!list.isEmpty()){ + final SimpleFeatureType PTYPE =createPointFeatureType(); + DefaultFeatureCollection pointFeatureCollection = new DefaultFeatureCollection(); + SimpleFeatureBuilder pointFeatureBuilder = new SimpleFeatureBuilder(PTYPE); + list.forEach(props->{ + pointFeatureBuilder.add(props.get("geom")); + pointFeatureBuilder.add(props.get("id")); + pointFeatureBuilder.add(props.get("bottomHeight")); + pointFeatureBuilder.add(props.get("topHeight")); + SimpleFeature feature = pointFeatureBuilder.buildFeature((String.valueOf(props.get("id")))); + pointFeatureCollection.add(feature); + }); + Style pointStyle = SLD.createSimpleStyle(pointFeatureCollection.getSchema()); + Layer pointLlayer = new FeatureLayer(pointFeatureCollection, pointStyle); + map.addLayer(pointLlayer); + } + }); + + final SimpleFeatureType PTYPE =createPointFeatureType(); + DefaultFeatureCollection pointFeatureCollection = new DefaultFeatureCollection(); + SimpleFeatureBuilder pointFeatureBuilder = new SimpleFeatureBuilder(PTYPE); + List> list = (List>)result.get("cutPointList"); + boundary.forEach(props->{ + pointFeatureBuilder.add(props.get("geom")); + pointFeatureBuilder.add(props.get("id")); + pointFeatureBuilder.add(props.get("bottomHeight")); + pointFeatureBuilder.add(props.get("topHeight")); + SimpleFeature feature = pointFeatureBuilder.buildFeature((String.valueOf(props.get("id")))); + pointFeatureCollection.add(feature); + }); + Style pointStyle = SLD.createSimpleStyle(pointFeatureCollection.getSchema(),Color.red); + Layer pointLlayer = new FeatureLayer(pointFeatureCollection, pointStyle); + map.addLayer(pointLlayer); + + final SimpleFeatureType FillPTYPE =createPointFeatureType(); + DefaultFeatureCollection pointFillFeatureCollection = new DefaultFeatureCollection(); + SimpleFeatureBuilder pointFillFeatureBuilder = new SimpleFeatureBuilder(FillPTYPE); + List> fillList = (List>)result.get("fillPointList"); + fillList.forEach(props->{ + pointFillFeatureBuilder.add(props.get("geom")); + pointFillFeatureBuilder.add(props.get("id")); + pointFillFeatureBuilder.add(props.get("bottomHeight")); + pointFillFeatureBuilder.add(props.get("topHeight")); + SimpleFeature feature = pointFillFeatureBuilder.buildFeature((String.valueOf(props.get("id")))); + pointFillFeatureCollection.add(feature); + }); + Style pointFillStyle = SLD.createSimpleStyle(pointFillFeatureCollection.getSchema(),Color.blue); + Layer pointFillLlayer = new FeatureLayer(pointFillFeatureCollection, pointFillStyle); + map.addLayer(pointFillLlayer); + + JMapFrame.showMap(map); + +// double sum = (p1Result+p2Result+p3Result+p4Result); +// double sumTop = (p1TopResult+p2TopResult+p3TopResult+p4TopResult); +// return sum+sumTop; + cutAllSum = cutSum; + cutPoints = cutPointList; + fillAllSum = fillSum; + fillPoints = fillPointList; + + return fillSum + cutSum; +// return 0; + } + private static List> searchBourdary(List> fillPointList, double radius) { + + AlphaShape alphaShape = new AlphaShape(radius); + List> result = new ArrayList<>(); + // 构建顶点集合对象 + List points = new ArrayList<>(); + fillPointList.forEach(pointMap->{ + Point p = (Point)pointMap.get("geom"); + Vector2D vector2D = new Vector2D(p.getX(), p.getY()); + + points.add(vector2D); + }); + List> lists = alphaShape.boundaryPoints(points); + lists.forEach(pointList->{ + pointList.forEach(point->{ + HashMap tmp = new HashMap<>(); + tmp.put("geom",geometryFactory.createPoint(new Coordinate(point.getX(),point.getY()))); + tmp.put("id",point.getIndex()); + tmp.put("bottomHeight",point.getBottomHeight()); + tmp.put("topHeight",point.getTopHeight()); + result.add(tmp); + }); + }); + + return result; + } + +// /** +// * 返回线要素的 +// * @param fillPointList +// * @param radius +// * @return +// */ +// private static SimpleFeature searchBourdary(List> fillPointList, double radius,String name) { +// +// AlphaShape alphaShape = new AlphaShape(radius); +// // 构建顶点集合对象 +// List points = new ArrayList<>(); +// fillPointList.forEach(pointMap->{ +// Point p = (Point)pointMap.get("geom"); +// Vector2D vector2D = new Vector2D(p.getX(), p.getY()); +//// vector2D.setBottomHeight(Double.valueOf(pointMap.get("bottomHeight").toString())); +//// vector2D.setTopHeight(Double.valueOf(pointMap.get("topHeight").toString())); +//// vector2D.setIndex((String)pointMap.get("id")); +// points.add(vector2D); +// }); +// List> lists = alphaShape.boundaryPoints(points); +// List coordList = new ArrayList(); +// lists.forEach(pointList->{ +// pointList.forEach(point->{ +// coordList.add(new Coordinate(point.getX(),point.getY())); +//// HashMap tmp = new HashMap<>(); +//// tmp.put("geom",geometryFactory.createPoint(new Coordinate(point.getX(),point.getY()))); +//// tmp.put("id",point.getIndex()); +//// tmp.put("bottomHeight",point.getBottomHeight()); +//// tmp.put("topHeight",point.getTopHeight()); +//// result.add(tmp); +// }); +// }); +// Coordinate[] CoordinateArray = coordList.stream().toArray(Coordinate[]::new); +// LineString lineString = geometryFactory.createLineString(CoordinateArray); +// +// // 构建添加线要素层 +// final SimpleFeatureType TYPE =createLineStringFeatureType(); +// SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE); +// String uuid = UUID.randomUUID().toString(); +// featureBuilder.add(lineString); +// featureBuilder.add(name); +// featureBuilder.add(uuid); +// SimpleFeature feature = featureBuilder.buildFeature(uuid); +//// featureCollection.add(feature); +// +// return feature; +// } + + private static Map calcGrid(List>> pointRingFeaturesProperties, Point centroid,double meter) throws Exception { + List> lastRing = pointRingFeaturesProperties.get(pointRingFeaturesProperties.size()-1); + double currDegree = meter2degree(meter/2); + List> fillPointList = new ArrayList<>();// 填方的点集 + List> cutPointList = new ArrayList<>();// 挖方的点集 + double fillSum = 0;// 填方量 + double cutSum = 0;// 挖方量 + + for(Map featuresProperties:lastRing){ +// Map featuresProperties = lastRing.get(300); + Point point = (Point)featuresProperties.get("geom"); + // 创建缓存面 + LineString lineString = geometryFactory.createLineString(new Coordinate[]{new Coordinate(point.getX(),point.getY()),new Coordinate(centroid.getX(),centroid.getY())}); + BufferParameters bufferParameters = new BufferParameters(); + bufferParameters.setEndCapStyle(BufferParameters.CAP_FLAT); + bufferParameters.setJoinStyle(BufferParameters.JOIN_MITRE); + BufferOp bufOp = new BufferOp(lineString,bufferParameters); + Polygon bufferLinePolygon = (Polygon)bufOp.getResultGeometry(currDegree); + + // 计算 + boolean isGoDown = true;// 是否继续往下 + boolean isGoHead = true; + outer: + for(int j=0;j> ringPointMap = pointRingFeaturesProperties.get(j); + for (Map pointMap : ringPointMap) { + Point p = (Point) pointMap.get("geom"); + if (bufferLinePolygon.contains(p)) { + double bottomHeight = (double) pointMap.get("bottomHeight"); + double topHeight = (double) pointMap.get("topHeight"); + double currHeight = getElevation(p); + + if(currHeight>bottomHeight&&currHeighttopHeight&&isGoHead){ + isGoDown =false; + if(!cutPointList.contains(pointMap)){//去除重复点，避免重复计算 + cutPointList.add(pointMap); + cutSum+= (currHeight - topHeight)*Math.pow(meter,2); + } + } + if(currHeightbottomHeight&&isGoDown){ +//// if(fillPointList.isEmpty()){ +//// fillPointList.add(pointMap); +//// } +// isGoDown = false; +// } +// // 实际高度高于顶部高度 +// if(currHeight>=topHeight){ +// isGoDown = false; +// if(!cutPointList.contains(pointMap)){ +// cutPointList.add(pointMap); +// cutSum+= (currHeight - topHeight)*Math.pow(meter,2); +// } +// } + } + } + } + } + + Map result = new HashMap<>(); + result.put("fillPointList",fillPointList); + result.put("cutPointList",cutPointList); + result.put("fillSum",fillSum); + result.put("cutSum",cutSum); +// result.put("bufferLinePolygon",bufferLinePolygon); + return result; + } + + private static Map filterGrid(List>> pointRingFeaturesProperties,GeometryFactory geometryFactory,Point centroid,double meter,String type){ + Map result = new HashMap<>(); + double sum = 0; + List> filterGrid = new ArrayList<>(); + for(int i=0;i> ringPointList = pointRingFeaturesProperties.get(i); + for(int j=0;j pfp = ringPointList.get(j); + double height = (double)pfp.get("height"); + Point p = (Point) pfp.get("geom"); + try { + double elevation = getElevation(p); + // 填 + if((type.equals("-1")&&height - elevation>=0)||(type.equals("1")&&height - elevation<0)){ + List> centerLinePointSet = getCenterLinePoint(geometryFactory,pointRingFeaturesProperties,i,centroid,meter,p); + if(centerLinePointSet.isEmpty()){ + filterGrid.add(pfp); + sum += Math.abs(height - elevation)*Math.pow(meter, 2); + }else{ + int count = 0; + for(int k=0;k lineOnPoint = centerLinePointSet.get(k); + double hp = (double)lineOnPoint.get("height"); + Point geom = (Point) lineOnPoint.get("geom"); + double ev = getElevation(geom); + if(ev>hp){ + count++; + } + } + if(count==0){ + filterGrid.add(pfp); + sum += Math.abs(height - elevation)*Math.pow(meter, 2); + } + } + } + } catch (Exception e) { + throw new RuntimeException(e); + } + } + }; + result.put("sum",Integer.valueOf(type)*sum); + result.put("filterGrid",filterGrid); + return result; + } + private static List> getCenterLinePoint( GeometryFactory geometryFactory,List>> pointFeaturesProperties,int index,Point centroid,double side,Point point){ + // 创建缓存面 + LineString lineString = geometryFactory.createLineString(new Coordinate[]{new Coordinate(point.getX(),point.getY()),new Coordinate(centroid.getX(),centroid.getY())}); + BufferParameters bufferParameters = new BufferParameters(); + bufferParameters.setEndCapStyle(BufferParameters.CAP_FLAT); + bufferParameters.setJoinStyle(BufferParameters.JOIN_MITRE); + BufferOp bufOp = new BufferOp(lineString,bufferParameters); + double currDegree = meter2degree(side/2); + Polygon bufferLine = (Polygon)bufOp.getResultGeometry(currDegree); + + List> result = new ArrayList<>(); + for(int i=0;i> ring = pointFeaturesProperties.get(i); + for(int j=0;j pfp = ring.get(j); + Point p = (Point) pfp.get("geom"); + if(bufferLine.intersects(p)){ + result.add(pfp); + } + } + } + return result; + } + + private static List>> createRingPoint(Polygon rectangle, double targetHeight, double bottomRatio,double topRatio, double side,GeometryFactory geometryFactory) { + + List>> ringPoint = new ArrayList<>(); + // 先在输入面上打点 + List> firstList = new ArrayList<>(); + Coordinate[] coordinates = rectangle.getCoordinates(); + Coordinate[] coordArr1 = new Coordinate[]{coordinates[0],coordinates[1]}; + Coordinate[] coordArr2 = new Coordinate[]{coordinates[1],coordinates[2]}; + Coordinate[] coordArr3 = new Coordinate[]{coordinates[2],coordinates[3]}; + Coordinate[] coordArr4 = new Coordinate[]{coordinates[3],coordinates[0]}; + firstList.addAll(addLinePoint(coordArr1,geometryFactory,side,targetHeight,targetHeight)); + firstList.addAll(addLinePoint(coordArr2,geometryFactory,side,targetHeight,targetHeight)); + firstList.addAll(addLinePoint(coordArr3,geometryFactory,side,targetHeight,targetHeight)); + firstList.addAll(addLinePoint(coordArr4,geometryFactory,side,targetHeight,targetHeight)); + ringPoint.add(firstList); + for(int i=1;i> pointRingFeaturesProperties = new ArrayList<>(); + double bottomHeight = targetHeight - side*bottomRatio*i; + double topHeight = targetHeight + side*topRatio*i; + BufferParameters bufferParameters = new BufferParameters(); + bufferParameters.setEndCapStyle(BufferParameters.CAP_FLAT); + bufferParameters.setJoinStyle(BufferParameters.JOIN_MITRE); + BufferOp bufOp = new BufferOp(rectangle,bufferParameters); + double currDegree = meter2degree(side*i); + Geometry bg = bufOp.getResultGeometry(currDegree); + Coordinate[] boundCoord = bg.getBoundary().getCoordinates(); + Coordinate[] boundArr1 = new Coordinate[]{boundCoord[0],boundCoord[1]}; + Coordinate[] boundArr2 = new Coordinate[]{boundCoord[1],boundCoord[2]}; + Coordinate[] boundArr3 = new Coordinate[]{boundCoord[2],boundCoord[3]}; + Coordinate[] boundArr4 = new Coordinate[]{boundCoord[3],boundCoord[0]}; + pointRingFeaturesProperties.addAll(addLinePoint(boundArr1,geometryFactory,side,bottomHeight,topHeight)); + pointRingFeaturesProperties.addAll(addLinePoint(boundArr2,geometryFactory,side,bottomHeight,topHeight)); + pointRingFeaturesProperties.addAll(addLinePoint(boundArr3,geometryFactory,side,bottomHeight,topHeight)); + pointRingFeaturesProperties.addAll(addLinePoint(boundArr4,geometryFactory,side,bottomHeight,topHeight)); + ringPoint.add(pointRingFeaturesProperties); + } + return ringPoint; + } + // 添加一个面的点集 + public static List>> addPartPoint(Coordinate[] bgcoords, + GeometryFactory geometryFactory, + double side, + Polygon p, + List> pointFeaturesProperties, + double outHeight, + String type){ + List>> gridList = new ArrayList(); + // 计算外边界上的点 + List> headPart = new ArrayList<>(); + double t1 = (bgcoords[1].y - bgcoords[0].y)/(bgcoords[1].x-bgcoords[0].x); + Coordinate startCoord = new Coordinate(bgcoords[0].x,bgcoords[0].y); + Map centerPointMap = new HashMap<>(); + Point center = geometryFactory.createPoint(startCoord); + + centerPointMap.put("geom",center); + centerPointMap.put("id",UUID.randomUUID()); + centerPointMap.put("height",outHeight); + if(!type.equals("0")){ + pointFeaturesProperties.add(centerPointMap); + headPart.add(centerPointMap); + } + List> lastPart = addPointSet(bgcoords,geometryFactory,side,p,pointFeaturesProperties,center,outHeight,type); + headPart.addAll(lastPart); + gridList.add(headPart); + + int len = (int)Math.ceil((Math.sqrt(Math.pow(bgcoords[1].y - bgcoords[0].y,2)+Math.pow(bgcoords[1].x - bgcoords[0].x,2)))/side); + //确定符号 + double xx = bgcoords[1].x-bgcoords[0].x; + double yy = bgcoords[1].y-bgcoords[0].y; + int xSign = xx>=0?1:-1; + int ySign = yy>=0?1:-1; + for(int i=1;i> extHeadPart = new ArrayList<>(); + Coordinate nextCoord = new Coordinate(nextX,nextY); + Map nextPointMap = new HashMap<>(); + Point nextPoint = geometryFactory.createPoint(nextCoord); + nextPointMap.put("geom",nextPoint); + nextPointMap.put("id",UUID.randomUUID()); + nextPointMap.put("height",outHeight); + if(type.equals("0")){ + if(p.intersects(nextPoint)){ + pointFeaturesProperties.add(nextPointMap); + extHeadPart.add(nextPointMap); + } + }else{ + pointFeaturesProperties.add(nextPointMap); + extHeadPart.add(nextPointMap); + } + // 添加垂直于边界线的点阵 + List> extLastPart = addPointSet(bgcoords,geometryFactory,side,p,pointFeaturesProperties,nextPoint,outHeight,type); + extHeadPart.addAll(extLastPart); + gridList.add(extHeadPart); + } + return gridList; + } + public static void targetPanelCalc (List> pointFeaturesProperties,double side){ + + AtomicReference sum = new AtomicReference<>((double) 0);// 土方量总和，负数挖、整数填 + pointFeaturesProperties.stream().forEach(pfp->{ + double h = (double)pfp.get("height"); + Point p = (Point)pfp.get("geom"); + try { + double elevation = getElevation(p); + sum.updateAndGet(v -> new Double((double) (v + (h - elevation) * Math.pow(degree2meter(side), 2)))); + } catch (Exception e) { + throw new RuntimeException(e); + } + }); + } + + public static double extPanelCalc(List>> gridPolygon,double side,String type){ + double sum = 0;// 土方量总和，整数填 + for(int i=0;i> col = gridPolygon.get(i); + // 逆序比较 + for(int j=col.size()-1;j>=0;j--){ + Map pfp = col.get(j); + double h = (double)pfp.get("height"); + Point p = (Point)pfp.get("geom"); + try { + double elevation = getElevation(p); + if(type.equals("-1")&&elevation>h||type.equals("1")&&elevation new Double((double) (v + (h - elevation) * Math.pow(degree2meter(side), 2)))); + sum += (h - elevation) * Math.pow(degree2meter(side), 2); + } catch (Exception e) { + throw new RuntimeException(e); + } + } + } +// System.out.println("gridPolygon sum:"+sum); + return sum; + } + // 添加一个线的点集 + public static List> addLinePoint(Coordinate[] bgcoords, GeometryFactory geometryFactory, double side, double bottomHeight,double topHeight){ + double degree = meter2degree(side); + // 计算外边界上的点 + List> lineList = new ArrayList<>(); + Coordinate startCoord = new Coordinate(bgcoords[0].x,bgcoords[0].y); + Map startPointMap = new HashMap<>(); + Point sp = geometryFactory.createPoint(startCoord); + startPointMap.put("geom",sp); + startPointMap.put("id",UUID.randomUUID().toString()); + startPointMap.put("bottomHeight",bottomHeight); + startPointMap.put("topHeight",topHeight); + lineList.add(startPointMap); + + int len = (int)Math.floor((Math.sqrt(Math.pow(bgcoords[1].y - bgcoords[0].y,2)+Math.pow(bgcoords[1].x - bgcoords[0].x,2)))/degree); + double t1 = (bgcoords[1].y - bgcoords[0].y)/(bgcoords[1].x-bgcoords[0].x); + //确定符号 + double xx = bgcoords[1].x-bgcoords[0].x; + double yy = bgcoords[1].y-bgcoords[0].y; + int xSign = xx>=0?1:-1; + int ySign = yy>=0?1:-1; + for(int i=1;i nextPointMap = new HashMap<>(); + Point nextPoint = geometryFactory.createPoint(nextCoord); + nextPointMap.put("geom",nextPoint); + nextPointMap.put("id",UUID.randomUUID().toString()); + nextPointMap.put("bottomHeight",bottomHeight); + nextPointMap.put("topHeight",topHeight); + lineList.add(nextPointMap); + } + return lineList; + } + public static List> addPointSet(Coordinate[] bgcoords, + GeometryFactory geometryFactory, + double side, + Polygon polygon, + List> pointFeaturesProperties, + Point startPoint, + double boundHeight, + String type){ + List> extendPointProperties = new ArrayList<>(); + + double t1 = (bgcoords[1].y - bgcoords[0].y)/(bgcoords[1].x-bgcoords[0].x); + //确定符号 + double xx = bgcoords[1].x-bgcoords[0].x; + double yy = bgcoords[1].y-bgcoords[0].y; + int xSign = (xx>=0&&yy>0)||(xx<0&&yy>0)?1:-1; + int ySign = (xx<0&&yy<0)||(xx<0&&yy>0)?1:-1; + + for(int i=1;i<300;i++){ + double x1 = xSign * Math.abs(Math.sin(Math.atan(t1)))*side*i + startPoint.getX(); + double y1 = ySign * Math.abs(Math.cos(Math.atan(t1)))*side*i + startPoint.getY(); + Point point = geometryFactory.createPoint(new Coordinate(x1,y1)); + boolean contains = polygon.contains(point); + if(contains){ + Map tmpPoint = new HashMap<>(); + tmpPoint.put("geom",point); + tmpPoint.put("id",UUID.randomUUID()); + double ph = boundHeight;//+(degree2meter(side))*i/2; + if(type.equals("-1")){//下挖的高程计算 + ph+=(degree2meter(side))*i/2; + }else if(type.equals("1")){// + ph-=(degree2meter(side))*i/1.5; + } + tmpPoint.put("height",ph); + pointFeaturesProperties.add(tmpPoint); + extendPointProperties.add(tmpPoint); + }else{ + break; + } + } + return extendPointProperties; + } + public static void loadDem() throws Exception{ + File file = new File("src\\main\\resources\\public\\data\\ASTGTMV003_N23E113_dem.tif"); + // Hints.FORCE_LONGITUDE_FIRST_AXIS_ORDER：设置经度为第一轴顺序 + GeoTiffReader reader = new GeoTiffReader(file, new Hints(Hints.FORCE_LONGITUDE_FIRST_AXIS_ORDER, Boolean.TRUE)); + GridCoverage2D coverage = reader.read(null); + //设置坐标系 + CoordinateReferenceSystem targetCRS = CRS.decode("EPSG:4326"); + // 将 GridCoverage 进行重采样转换为另一个 CRS + CreateRingGridFix.coverage = (GridCoverage2D) Operations.DEFAULT.resample(coverage, targetCRS); + } + public static double getElevation(Point p) throws Exception{ + // 设置经纬度及坐标系等信息 + Position position = new Position2D(p.getY(), p.getX()); + GridGeometry2D gridGeometry = CreateRingGridFix.coverage.getGridGeometry(); + GridCoordinates2D gridPoint = gridGeometry.worldToGrid(position); + double[] elevation = new double[1]; + CreateRingGridFix.coverage.evaluate(gridPoint,elevation); + double height = elevation[0]; // 这里的高度就是查询位置的高程值 +// System.out.println("---------"+height);//---------108.0 + return height; + } + + public static double meter2degree(double meter){ + return meter/(Math.PI*6371004/180); + } + public static double degree2meter(double degree){ + return degree*(Math.PI*6371004)/180; + } +} diff --git a/src/main/java/org/example/test/CreateTiff.java b/src/main/java/org/example/test/CreateTiff.java index 742a2ec..899981e 100644 --- a/src/main/java/org/example/test/CreateTiff.java +++ b/src/main/java/org/example/test/CreateTiff.java @@ -139,12 +139,14 @@ public class CreateTiff { if (value > 0){ System.out.println(i + "," + j + "," + b + "," + value); } + +// coverage. } } } //tiffWrite2(crs,minx,miny,width,height,type,bandNum,matrixW,matrixH,matrix); -// tiffWrite3(coverage); + tiffWrite3(coverage); } public static void tiffWrite3(GridCoverage2D coverage) throws IOException { @@ -155,4 +157,6 @@ public class CreateTiff { writer2.dispose(); coverage.dispose(true); } + + } diff --git a/src/main/java/org/example/test/Demo3.java b/src/main/java/org/example/test/Demo3.java new file mode 100644 index 0000000..0b84a6a --- /dev/null +++ b/src/main/java/org/example/test/Demo3.java @@ -0,0 +1,214 @@ +package org.example.test; + +import java.io.File; + +import org.opencv.core.Core; +import org.opencv.core.CvType; +import org.opencv.core.Mat; +import org.opencv.core.Size; +import org.opencv.highgui.HighGui; +import org.opencv.imgproc.Imgproc; +import org.opencv.imgcodecs.Imgcodecs; +/** + * 创建日期:2018年1月7日 + * 创建时间:上午10:38:06 + * 创建者 :yellowcong + * 机能概要:利用Opencv获取图片的轮廓 + */ +public class Demo3 { + public static void main(String[] args) { + //图片地址 +// String inputImagePath = Demo3.class.getClassLoader().getResource("public/image/shishi.jpg").getFile(); + String inputImagePath = "E:\\data\\test1.tiff"; + + String outPath = "D:/test/demo2.png"; + + //加载lib,这个lib的名称 + String path = "D:\\tools\\opencv\\opencv\\build\\java\\x64\\" + System.mapLibraryName("opencv_java490"); + File lib = new File(path); + System.load(lib.getAbsolutePath()); +// System.loadLibrary(lib.getAbsolutePath()); +// System.loadLibrary(Core.NATIVE_LIBRARY_NAME); + + //读取图片信息 + Mat image = Imgcodecs.imread(new File(inputImagePath).getAbsolutePath(),1); + testCanny(image); +// testSobel(image); +// testScharr(image); +// testLaplacian(image); +// testThreshold(image); +// testAdaptiveThreshold(image); + //将rgb灰化处理 +// Imgproc.cvtColor(image, image,Imgproc.COLOR_BGR2GRAY); +// +// //平滑处理 +// Imgproc.blur(image, image, new Size(2, 2)); +// +// //轮廓 +// //使用Canndy检测边缘 +// double lowThresh =100;//双阀值抑制中的低阀值 +// double heightThresh = 300;//双阀值抑制中的高阀值 +// Imgproc.Canny(image, image,lowThresh, heightThresh); +// +// // 写入到文件 + Imgcodecs.imwrite(outPath, image); + + + } + // canny算子 + public static void testCanny(Mat image) { + //高斯滤波 + Imgproc.GaussianBlur(image, image, new Size(3, 3), 0, 0); + //灰度 + Imgproc.cvtColor(image, image,Imgproc.COLOR_BGR2GRAY); + //一般来说，高阈值maxVal推荐是低阈值minVal的2~3倍 + int lowThresh=45; + //边缘检测 + Imgproc.Canny(image, image,lowThresh, lowThresh*3,3); +// HighGui.imshow("Canny", image); +// HighGui.waitKey(0); + } + + // Sobel算子 + public static void testSobel(Mat image) { + //高斯滤波 + Imgproc.GaussianBlur(image, image, new Size(3, 3), 0, 0); + //灰度 + Imgproc.cvtColor(image, image,Imgproc.COLOR_BGR2GRAY); + //横向 + Mat x = new Mat(); + /*src：源图像 + *dst：检测结果图像 + *ddepth：输出图像的深度 + *dx：x方向上的差分阶数 + *dy：y方向上的差分阶数 + *ksize：sobel核的大小，默认为3 + *scale：缩放因子 + *delta：结果存入输出图像前可选的delta值，默认为0 + *borderType：边界模式，默认BORDER_DEFAULT + * + *其中，输出图像的深度，支持如下src.depth()和ddepth的组合： + *src.depth() = CV_8U ddepth =-1/CV_16S/CV_32F/CV_64F + *src.depth() = CV_16U/CV_16S ddepth =-1/CV_32F/CV_64F + *src.depth() = CV_32F ddepth =-1/CV_32F/CV_64F + *src.depth() = CV_64F ddepth = -1/CV_64F + */ + Imgproc.Sobel(image, x, CvType.CV_16S, 1, 0, 3, 1, 0, Core.BORDER_DEFAULT); + //需要用convertScaleAbs()函数将其转回原来的uint8形式，否则将无法显示图像 + Core.convertScaleAbs(x, x); + HighGui.imshow("x", x); + //竖向 + Mat y = new Mat(); + Imgproc.Sobel(image, y, CvType.CV_16S, 0, 1, 3, 1, 0, Core.BORDER_DEFAULT); + Core.convertScaleAbs(y, y); + HighGui.imshow("y", y); + //横竖向图像融合 + Mat xy = new Mat(); + Core.addWeighted(x, 0.5, y, 0.5, 0, xy); + HighGui.imshow("xy", xy); + HighGui.waitKey(0); + } + // Scharr滤波器 + public static void testScharr(Mat image) { + //高斯滤波 + Imgproc.GaussianBlur(image, image, new Size(3, 3), 0, 0); + //灰度 + Imgproc.cvtColor(image, image,Imgproc.COLOR_BGR2GRAY); + //横向 + Mat x = new Mat(); + /*src：源图像 + *dst：检测结果图像 + *ddepth：输出图像的深度 + *dx：x方向上的差分阶数 + *dy：y方向上的差分阶数 + *scale：缩放因子 + *delta：结果存入输出图像前可选的delta值，默认为0 + *borderType：边界模式，默认BORDER_DEFAULT + * + *其中，输出图像的深度，支持如下src.depth()和ddepth的组合： + *src.depth() = CV_8U ddepth =-1/CV_16S/CV_32F/CV_64F + *src.depth() = CV_16U/CV_16S ddepth =-1/CV_32F/CV_64F + *src.depth() = CV_32F ddepth =-1/CV_32F/CV_64F + *src.depth() = CV_64F ddepth = -1/CV_64F + */ + Imgproc.Scharr(image, x, CvType.CV_16S, 1, 0,1, 0, Core.BORDER_DEFAULT); + //需要用convertScaleAbs()函数将其转回原来的uint8形式，否则将无法显示图像 + Core.convertScaleAbs(x, x); + HighGui.imshow("x", x); + //竖向 + Mat y = new Mat(); + Imgproc.Scharr(image, y, CvType.CV_16S, 0, 1, 1, 0, Core.BORDER_DEFAULT); + Core.convertScaleAbs(y, y); + HighGui.imshow("y", y); + //横竖向图像融合 + Mat xy = new Mat(); + Core.addWeighted(x, 0.5, y, 0.5, 0, xy); + HighGui.imshow("xy", xy); + HighGui.waitKey(0); + } + // Laplace算子 + public static void testLaplacian(Mat image) { + //高斯滤波 + Imgproc.GaussianBlur(image, image, new Size(3, 3), 0, 0); + //灰度 + Imgproc.cvtColor(image, image,Imgproc.COLOR_BGR2GRAY); + //边缘检测 + Mat dst = new Mat(); + /*src：源图像 + *dst：检测结果图像 + *ddepth：输出图像的深度 + *ksize：计算二阶导数的滤波器的孔径尺寸，必须为正奇数 等于1是四邻域算子，大于1改用八邻域算子 + *scale：缩放因子 计算拉普拉斯值时候可选的比例因子 + *delta：结果存入输出图像前可选的delta值，默认为0 + *borderType：边界模式，默认BORDER_DEFAULT + */ + Imgproc.Laplacian(image, dst, CvType.CV_16S, 3, 1, 0, Core.BORDER_DEFAULT); + Core.convertScaleAbs(dst, dst); + HighGui.imshow("Laplacian", dst); + HighGui.waitKey(0); + } + + // Threshold（固定阈值操作） + public static void testThreshold(Mat image) { + /** + *src 是输入的函数图像 + *dst 是输出的函数图像 + *thresh 阈值的具体值 + *maxval 阈值类型的最大值 + *type 阈值类型 + *THRESH_BINARY = 0,过门限的为maxval其他取零 + *THRESH_BINARY_INV = 1,过门限的为取零，其他maxval + *THRESH_TRUNC = 2,过门限的取门限，其他不变 + *THRESH_TOZERO = 3,过门限的不变，其他取零 + *THRESH_TOZERO_INV = 4,过门限的值取零其他不变 + *THRESH_MASK = 7, + *THRESH_OTSU = 8,自动生成阀值，大于阀值的为255 ，小于阀值的为0 + *THRESH_TRIANGLE = 16;和我们设置的阀值大小没有关系，是自动计算出来的 + */ + Imgproc.cvtColor(image, image,Imgproc.COLOR_BGR2GRAY); + HighGui.imshow("gry", image); + Mat m=new Mat(); + Imgproc.threshold(image, m, 120, 255, Imgproc.THRESH_BINARY); + HighGui.imshow("Threshold", m); + HighGui.waitKey(0); + } + + // adaptiveThreshold（自适应阈值操作） + public static void testAdaptiveThreshold(Mat image) { + /** + * src + * dst + * maxValue 给像素赋的满足条件的非零值 + * ADAPTIVE_THRESH_MEAN_C 计算均值时每个像素的权值是相等的 + * ADAPTIVE_THRESH_GAUSSIAN_C 计算均值时每个像素的权值根据其到中心点的距离通过高斯方程得到 + * int blockSize, //计算阈值大小的一个像素的领域尺寸，取值为大于1的奇数 + * double C); //减去加权平均值后的常数值，通常为正数，少数情况下也可为0或负数 + */ + //(Mat src, Mat dst, double maxValue, int adaptiveMethod, int thresholdType, int blockSize, double C) + //CV_8UC1 in function 'cv::adaptiveThreshold' 先灰度处理 + Imgproc.cvtColor(image, image,Imgproc.COLOR_BGR2GRAY); + Imgproc.adaptiveThreshold(image, image, 200, Imgproc.ADAPTIVE_THRESH_MEAN_C, Imgproc.THRESH_BINARY, 5, 9); + HighGui.imshow("adaptiveThreshold", image); + HighGui.waitKey(0); + } +} diff --git a/src/main/java/org/example/test/Test.java b/src/main/java/org/example/test/Test.java new file mode 100644 index 0000000..70bdb0a --- /dev/null +++ b/src/main/java/org/example/test/Test.java @@ -0,0 +1,34 @@ +package org.example.test; + +import java.io.File; + +import org.opencv.core.Mat; +import org.opencv.highgui.HighGui; +import org.opencv.imgcodecs.Imgcodecs; + +/** + * @author Administrator + */ +public class Test { + + static { + String os = System.getProperty("os.name"); + String type = System.getProperty("sun.arch.data.model"); + if (os.toUpperCase().contains("WINDOWS")) { + File lib; + if (type.endsWith("64")) { + lib = new File("D:\\tools\\opencv\\opencv\\build\\java\\x64\\" + System.mapLibraryName("opencv_java490")); + } else { + lib = new File("D:\\tools\\opencv\\opencv\\build\\java\\x86\\" + System.mapLibraryName("opencv_java490")); + } + System.load(lib.getAbsolutePath()); + } + } + + public static void main(String[] args) { + Mat img = Imgcodecs.imread("C:\\Users\\user\\Desktop\\shishi.jpg"); + HighGui.imshow("Test", img); + HighGui.waitKey(0); + } + +} diff --git a/src/main/java/org/example/utils/CalcRingGridFix.java b/src/main/java/org/example/utils/CalcRingGridFix.java new file mode 100644 index 0000000..6f01b0e --- /dev/null +++ b/src/main/java/org/example/utils/CalcRingGridFix.java @@ -0,0 +1,595 @@ +package org.example.utils; + +import org.example.alpha.AlphaShape; +import org.example.alpha.Vector2D; +import org.example.test.CreateRingGridFix; +import org.geotools.api.feature.simple.SimpleFeature; +import org.geotools.api.feature.simple.SimpleFeatureType; +import org.geotools.api.geometry.Position; +import org.geotools.api.referencing.crs.CoordinateReferenceSystem; +import org.geotools.coverage.grid.GridCoordinates2D; +import org.geotools.coverage.grid.GridCoverage2D; +import org.geotools.coverage.grid.GridGeometry2D; +import org.geotools.coverage.processing.Operations; +import org.geotools.feature.simple.SimpleFeatureBuilder; +import org.geotools.feature.simple.SimpleFeatureTypeBuilder; +import org.geotools.gce.geotiff.GeoTiffReader; +import org.geotools.geojson.feature.FeatureJSON; +import org.geotools.geojson.geom.GeometryJSON; +import org.geotools.geometry.Position2D; +import org.geotools.referencing.CRS; +import org.geotools.referencing.crs.DefaultGeographicCRS; +import org.geotools.util.factory.Hints; +import org.locationtech.jts.geom.*; +import org.locationtech.jts.operation.buffer.BufferOp; +import org.locationtech.jts.operation.buffer.BufferParameters; + +import java.io.File; +import java.util.*; +import java.util.concurrent.atomic.AtomicReference; + +public class CalcRingGridFix { + public static GridCoverage2D coverage; + public static String demFilePath = "src\\main\\resources\\public\\data\\ASTGTMV003_N23E113_dem.tif"; + public static double extendLen = 100; + public static GeometryFactory geometryFactory = new GeometryFactory(); + public static List> fillPoints; + public static List> cutPoints; + public static double cutAllSum = 0; + public static double fillAllSum = 0; + public static double meter = 2; // 间距,单位为米。 + public static double range = 10;//容差 + + private static SimpleFeatureType createLineStringFeatureType() { + SimpleFeatureTypeBuilder builder = new SimpleFeatureTypeBuilder(); + builder.setName("Location"); + builder.setCRS(DefaultGeographicCRS.WGS84); // <- Coordinate reference system + + // add attributes in order + builder.add("the_geom", LineString.class); + builder.length(15).add("name", String.class); // <- 15 chars width for name field + builder.add("id", Integer.class); + // build the type + final SimpleFeatureType LOCATION = builder.buildFeatureType(); + return LOCATION; + } + + public static void main(String[] args) throws Exception{ + +// // 假设我们使用WGS84坐标系统 +// CoordinateReferenceSystem crs = CRS.decode("EPSG:4326"); +// GeometryFactory geometryFactory = new GeometryFactory(); + // 定义矩形的左下角和右上角坐标（WGS84） + // 输入参数：四点坐标，最大挖填方允许差值 + List coords = new ArrayList(){{ + add(new double[]{113.52896835361318, 23.654018962759377}); + add(new double[]{113.529270627066524, 23.654027766840539 }); + add(new double[]{ 113.529294104616298, 23.653772448486752}); + add(new double[]{ 113.528977157694356, 23.653763644405586 }); + add(new double[]{113.52896835361318, 23.654018962759377}); + }}; + // 创建矩形的四个角点坐标 + Coordinate[] coordinates = new Coordinate[5]; + coordinates[0] = new Coordinate(coords.get(0)[0], coords.get(0)[1]); + coordinates[1] = new Coordinate(coords.get(1)[0], coords.get(1)[1]); + coordinates[2] = new Coordinate(coords.get(2)[0], coords.get(2)[1]); + coordinates[3] = new Coordinate(coords.get(3)[0], coords.get(3)[1]); + // 闭合多边形，复制第一个坐标作为最后一个坐标 + coordinates[4] = new Coordinate(coords.get(0)[0], coords.get(0)[1]); + + // 计算设计标高，取四个角的高程平均值 +// double centerHeight = getElevation(rectangle.getCentroid()); + double corner1Height = getElevation( geometryFactory.createPoint(coordinates[0])); + double corner2Height = getElevation( geometryFactory.createPoint(coordinates[1])); + double corner3Height = getElevation( geometryFactory.createPoint(coordinates[2])); + double corner4Height = getElevation( geometryFactory.createPoint(coordinates[3])); + double targetHeight = (corner1Height+corner2Height+corner3Height+corner4Height)/4; + + System.out.println("初始标高："+targetHeight); + double startTime = System.currentTimeMillis(); + double allSum = calcBlanceHeight(coordinates,targetHeight); + double step = 5;//起始步长 + double nextHeight=targetHeight; + while (Math.abs(allSum)-range>0){ + System.out.println("allSum:"+allSum); + nextHeight = allSum<0?nextHeight-step:nextHeight+step; + System.out.println("nextHeight:"+nextHeight); + double firstSum = calcBlanceHeight(coordinates,nextHeight); + System.out.println("firstSum:"+firstSum); + if(Math.abs(firstSum)-range<0){ + allSum = firstSum; + break; + } + if((firstSum>0&&allSum<0)||(firstSum<0&&allSum>0)){ + // fix: 采用线性方式去计算步长 + step = step*Math.min(Math.abs(allSum),Math.abs(firstSum))/Math.abs(allSum-firstSum); +// step = iterator.next(); + } + double tmpHeight = firstSum<0?nextHeight-step:nextHeight+step; + System.out.println("nextHeight:"+tmpHeight); + allSum = calcBlanceHeight(coordinates,tmpHeight); + if(Math.abs(tmpHeight-nextHeight)<0.0001){ + nextHeight = tmpHeight; + break; + } + nextHeight = tmpHeight; + } + System.out.println("理想标高:"+nextHeight); + System.out.println("总挖方量："+cutAllSum); + System.out.println("总填方量："+fillAllSum); + System.out.println("填挖方差值："+allSum); + + System.out.println("耗时："+(System.currentTimeMillis()-startTime)); + + double radius = meter2degree(meter); + SimpleFeature cutBoundary = searchBourdary(cutPoints,radius,"挖方边界线"); + SimpleFeature fillBoundary = searchBourdary(fillPoints,radius,"填方边界线"); + FeatureJSON featureJSON = new FeatureJSON(new GeometryJSON(8)); + System.out.println("挖方边界线:"+featureJSON.toString(cutBoundary)); + System.out.println("填方边界线:"+featureJSON.toString(fillBoundary)); + // 输出：理想标高、总挖方量、总填方量、填挖方差值、边界数据 + } + + public static Map getFillCutBlance(Coordinate[] coordinates,double inputRange) throws Exception{ + range = inputRange; + // 计算设计标高，取四个角的高程平均值 + double corner1Height = getElevation( geometryFactory.createPoint(coordinates[0])); + double corner2Height = getElevation( geometryFactory.createPoint(coordinates[1])); + double corner3Height = getElevation( geometryFactory.createPoint(coordinates[2])); + double corner4Height = getElevation( geometryFactory.createPoint(coordinates[3])); + double targetHeight = (corner1Height+corner2Height+corner3Height+corner4Height)/4; + + System.out.println("初始标高："+targetHeight); + double startTime = System.currentTimeMillis(); + double allSum = calcBlanceHeight(coordinates,targetHeight); + double step = 5;//起始步长 + double nextHeight=targetHeight; + while (Math.abs(allSum)-range>0){ + System.out.println("allSum:"+allSum); + nextHeight = allSum<0?nextHeight-step:nextHeight+step; + System.out.println("nextHeight:"+nextHeight); + double firstSum = calcBlanceHeight(coordinates,nextHeight); + System.out.println("firstSum:"+firstSum); + if(Math.abs(firstSum)-range<0){ + allSum = firstSum; + break; + } + if((firstSum>0&&allSum<0)||(firstSum<0&&allSum>0)){ + // fix: 采用线性方式去计算步长 + step = step*Math.min(Math.abs(allSum),Math.abs(firstSum))/Math.abs(allSum-firstSum); +// step = iterator.next(); + } + double tmpHeight = firstSum<0?nextHeight-step:nextHeight+step; + System.out.println("nextHeight:"+tmpHeight); + allSum = calcBlanceHeight(coordinates,tmpHeight); + if(Math.abs(tmpHeight-nextHeight)<0.0001){ + nextHeight = tmpHeight; + break; + } + nextHeight = tmpHeight; + } + double radius = meter2degree(meter); + SimpleFeature cutBoundary = searchBourdary(cutPoints,radius,"挖方边界线"); + SimpleFeature fillBoundary = searchBourdary(fillPoints,radius,"填方边界线"); + FeatureJSON featureJSON = new FeatureJSON(new GeometryJSON(8)); + + Map result = new HashMap<>(); + result.put("理想标高",nextHeight); + result.put("总挖方量",cutAllSum); + result.put("总填方量",fillAllSum); + result.put("填挖方差值",allSum); + result.put("耗时",(System.currentTimeMillis()-startTime)+" 毫秒"); + result.put("挖方边界线",featureJSON.toString(cutBoundary).replace("\"","'")); + result.put("填方边界线",featureJSON.toString(fillBoundary).replace("\"","'")); + return result; + } + + public static double calcBlanceHeight(Coordinate[] coordinates,double targetHeight) throws Exception{ + // 创建GeometryFactory + List> featuresProperties = new ArrayList<>();// 面要素属性集合 + List> pointFeaturesProperties = new ArrayList<>();// 扩展点属性集合 + List> tPointFeaturesProperties = new ArrayList<>();// 目标平面的点集属性集合 + + // 清空 + fillPoints = new ArrayList<>(); + cutPoints = new ArrayList<>(); + cutAllSum = 0; + fillAllSum = 0; + + double bottomRatio = 1.0/2; + double topRatio = 1.0/1.5; + // 创建多边形（矩形） + Map innerPolygon = new HashMap<>(); + Polygon rectangle = geometryFactory.createPolygon(coordinates); + innerPolygon.put("geom",rectangle); + innerPolygon.put("name","输入的面"); + innerPolygon.put("id",1); + featuresProperties.add(innerPolygon); + + // 外接矩形打点 + Map outPolygon1 = new HashMap<>(); + Polygon out1 = (Polygon) rectangle.getEnvelope();//geometryFactory.createPolygon(cords4); + outPolygon1.put("geom",out1); + outPolygon1.put("name","输入面的外接矩形"); + outPolygon1.put("id",6); +// featuresProperties.add(outPolygon1); + Coordinate[] outCoords = out1.getCoordinates(); + Coordinate[] outLine = new Coordinate[]{outCoords[0],outCoords[1]}; + addPartPoint(outLine,geometryFactory,meter2degree(meter),rectangle,tPointFeaturesProperties,targetHeight,"0"); + + // 扩展面打点 + List>> pointRingFeaturesProperties = createRingPoint(rectangle,targetHeight,bottomRatio,topRatio,meter,geometryFactory); + + // 中心点 + Point centroid = rectangle.getCentroid(); + // 计算网格 + Map result = calcGrid(pointRingFeaturesProperties,centroid,meter); + List> fillPointList = (List>)result.get("fillPointList"); + List> cutPointList = (List>)result.get("cutPointList"); + double fillSum = (double)result.get("fillSum"); + double cutSum = (double)result.get("cutSum"); + + // 计算设计平面内的填挖方量 + for(int i=0;i pfp = tPointFeaturesProperties.get(i); + double height = (double)pfp.get("height"); + Point p = (Point)pfp.get("geom"); + double elev = getElevation(p); + if(elev>height){ + cutPointList.add(pfp); + cutSum += (elev - height)*Math.pow(meter,2); + }else{ + fillPointList.add(pfp); + fillSum += (elev - height)*Math.pow(meter,2); + } + } + cutAllSum = cutSum; + cutPoints = cutPointList; + fillAllSum = fillSum; + fillPoints = fillPointList; + return fillSum + cutSum; + } + + /** + * 返回线要素的 + * @param fillPointList + * @param radius + * @return + */ + private static SimpleFeature searchBourdary(List> fillPointList, double radius,String name) { + + AlphaShape alphaShape = new AlphaShape(radius); + // 构建顶点集合对象 + List points = new ArrayList<>(); + fillPointList.forEach(pointMap->{ + Point p = (Point)pointMap.get("geom"); + Vector2D vector2D = new Vector2D(p.getX(), p.getY()); + points.add(vector2D); + }); + List> lists = alphaShape.boundaryPoints(points); + List coordList = new ArrayList(); + lists.forEach(pointList->{ + pointList.forEach(point->{ + coordList.add(new Coordinate(point.getX(),point.getY())); + }); + }); + Coordinate[] CoordinateArray = coordList.stream().toArray(Coordinate[]::new); + LineString lineString = geometryFactory.createLineString(CoordinateArray); + + // 构建添加线要素层 + final SimpleFeatureType TYPE =createLineStringFeatureType(); + SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE); + String uuid = UUID.randomUUID().toString(); + featureBuilder.add(lineString); + featureBuilder.add(name); + featureBuilder.add(uuid); + SimpleFeature feature = featureBuilder.buildFeature(uuid); + return feature; + } + + private static Map calcGrid(List>> pointRingFeaturesProperties, Point centroid,double meter) throws Exception { + List> lastRing = pointRingFeaturesProperties.get(pointRingFeaturesProperties.size()-1); + double currDegree = meter2degree(meter/2); + List> fillPointList = new ArrayList<>();// 填方的点集 + List> cutPointList = new ArrayList<>();// 挖方的点集 + double fillSum = 0;// 填方量 + double cutSum = 0;// 挖方量 + for(Map featuresProperties:lastRing){ + Point point = (Point)featuresProperties.get("geom"); + // 创建缓存面 + LineString lineString = geometryFactory.createLineString(new Coordinate[]{new Coordinate(point.getX(),point.getY()),new Coordinate(centroid.getX(),centroid.getY())}); + BufferParameters bufferParameters = new BufferParameters(); + bufferParameters.setEndCapStyle(BufferParameters.CAP_FLAT); + bufferParameters.setJoinStyle(BufferParameters.JOIN_MITRE); + BufferOp bufOp = new BufferOp(lineString,bufferParameters); + Polygon bufferLinePolygon = (Polygon)bufOp.getResultGeometry(currDegree); + + // 计算 + boolean isGoDown = true;// 是否继续往下 + boolean isGoHead = true; + outer: + for(int j=0;j> ringPointMap = pointRingFeaturesProperties.get(j); + for (Map pointMap : ringPointMap) { + Point p = (Point) pointMap.get("geom"); + if (bufferLinePolygon.contains(p)) { + double bottomHeight = (double) pointMap.get("bottomHeight"); + double topHeight = (double) pointMap.get("topHeight"); + double currHeight = getElevation(p); + + if(currHeight>bottomHeight&&currHeight=topHeight&&isGoHead){ + isGoDown =false; + if(!cutPointList.contains(pointMap)){//去除重复点，避免重复计算 + cutPointList.add(pointMap); + cutSum+= (currHeight - topHeight)*Math.pow(meter,2); + } + } + if(currHeight<=bottomHeight&&isGoDown){ + isGoHead = false; + if(!fillPointList.contains(pointMap)){//去除重复点，避免重复计算 + fillPointList.add(pointMap); + fillSum+= (currHeight - bottomHeight)*Math.pow(meter,2); + } + } + } + } + } + } + + Map result = new HashMap<>(); + result.put("fillPointList",fillPointList); + result.put("cutPointList",cutPointList); + result.put("fillSum",fillSum); + result.put("cutSum",cutSum); + return result; + } + + private static List>> createRingPoint(Polygon rectangle, double targetHeight, double bottomRatio,double topRatio, double side,GeometryFactory geometryFactory) { + + List>> ringPoint = new ArrayList<>(); + // 先在输入面上打点 + List> firstList = new ArrayList<>(); + Coordinate[] coordinates = rectangle.getCoordinates(); + Coordinate[] coordArr1 = new Coordinate[]{coordinates[0],coordinates[1]}; + Coordinate[] coordArr2 = new Coordinate[]{coordinates[1],coordinates[2]}; + Coordinate[] coordArr3 = new Coordinate[]{coordinates[2],coordinates[3]}; + Coordinate[] coordArr4 = new Coordinate[]{coordinates[3],coordinates[0]}; + firstList.addAll(addLinePoint(coordArr1,geometryFactory,side,targetHeight,targetHeight)); + firstList.addAll(addLinePoint(coordArr2,geometryFactory,side,targetHeight,targetHeight)); + firstList.addAll(addLinePoint(coordArr3,geometryFactory,side,targetHeight,targetHeight)); + firstList.addAll(addLinePoint(coordArr4,geometryFactory,side,targetHeight,targetHeight)); + ringPoint.add(firstList); + for(int i=1;i> pointRingFeaturesProperties = new ArrayList<>(); + double bottomHeight = targetHeight - side*bottomRatio*i; + double topHeight = targetHeight + side*topRatio*i; + BufferParameters bufferParameters = new BufferParameters(); + bufferParameters.setEndCapStyle(BufferParameters.CAP_FLAT); + bufferParameters.setJoinStyle(BufferParameters.JOIN_MITRE); + BufferOp bufOp = new BufferOp(rectangle,bufferParameters); + double currDegree = meter2degree(side*i); + Geometry bg = bufOp.getResultGeometry(currDegree); + Coordinate[] boundCoord = bg.getBoundary().getCoordinates(); + Coordinate[] boundArr1 = new Coordinate[]{boundCoord[0],boundCoord[1]}; + Coordinate[] boundArr2 = new Coordinate[]{boundCoord[1],boundCoord[2]}; + Coordinate[] boundArr3 = new Coordinate[]{boundCoord[2],boundCoord[3]}; + Coordinate[] boundArr4 = new Coordinate[]{boundCoord[3],boundCoord[0]}; + pointRingFeaturesProperties.addAll(addLinePoint(boundArr1,geometryFactory,side,bottomHeight,topHeight)); + pointRingFeaturesProperties.addAll(addLinePoint(boundArr2,geometryFactory,side,bottomHeight,topHeight)); + pointRingFeaturesProperties.addAll(addLinePoint(boundArr3,geometryFactory,side,bottomHeight,topHeight)); + pointRingFeaturesProperties.addAll(addLinePoint(boundArr4,geometryFactory,side,bottomHeight,topHeight)); + ringPoint.add(pointRingFeaturesProperties); + } + return ringPoint; + } + // 添加一个面的点集 + public static List>> addPartPoint(Coordinate[] bgcoords, + GeometryFactory geometryFactory, + double side, + Polygon p, + List> pointFeaturesProperties, + double outHeight, + String type){ + List>> gridList = new ArrayList(); + // 计算外边界上的点 + List> headPart = new ArrayList<>(); + double t1 = (bgcoords[1].y - bgcoords[0].y)/(bgcoords[1].x-bgcoords[0].x); + Coordinate startCoord = new Coordinate(bgcoords[0].x,bgcoords[0].y); + Map centerPointMap = new HashMap<>(); + Point center = geometryFactory.createPoint(startCoord); + + centerPointMap.put("geom",center); + centerPointMap.put("id",UUID.randomUUID()); + centerPointMap.put("height",outHeight); + if(!type.equals("0")){ + pointFeaturesProperties.add(centerPointMap); + headPart.add(centerPointMap); + } + List> lastPart = addPointSet(bgcoords,geometryFactory,side,p,pointFeaturesProperties,center,outHeight,type); + headPart.addAll(lastPart); + gridList.add(headPart); + + int len = (int)Math.ceil((Math.sqrt(Math.pow(bgcoords[1].y - bgcoords[0].y,2)+Math.pow(bgcoords[1].x - bgcoords[0].x,2)))/side); + //确定符号 + double xx = bgcoords[1].x-bgcoords[0].x; + double yy = bgcoords[1].y-bgcoords[0].y; + int xSign = xx>=0?1:-1; + int ySign = yy>=0?1:-1; + for(int i=1;i> extHeadPart = new ArrayList<>(); + Coordinate nextCoord = new Coordinate(nextX,nextY); + Map nextPointMap = new HashMap<>(); + Point nextPoint = geometryFactory.createPoint(nextCoord); + nextPointMap.put("geom",nextPoint); + nextPointMap.put("id",UUID.randomUUID()); + nextPointMap.put("height",outHeight); + if(type.equals("0")){ + if(p.intersects(nextPoint)){ + pointFeaturesProperties.add(nextPointMap); + extHeadPart.add(nextPointMap); + } + }else{ + pointFeaturesProperties.add(nextPointMap); + extHeadPart.add(nextPointMap); + } + // 添加垂直于边界线的点阵 + List> extLastPart = addPointSet(bgcoords,geometryFactory,side,p,pointFeaturesProperties,nextPoint,outHeight,type); + extHeadPart.addAll(extLastPart); + gridList.add(extHeadPart); + } + return gridList; + } + public static void targetPanelCalc (List> pointFeaturesProperties,double side){ + + AtomicReference sum = new AtomicReference<>((double) 0);// 土方量总和，负数挖、整数填 + pointFeaturesProperties.stream().forEach(pfp->{ + double h = (double)pfp.get("height"); + Point p = (Point)pfp.get("geom"); + try { + double elevation = getElevation(p); + sum.updateAndGet(v -> new Double((double) (v + (h - elevation) * Math.pow(degree2meter(side), 2)))); + } catch (Exception e) { + throw new RuntimeException(e); + } + }); + } + + public static double extPanelCalc(List>> gridPolygon,double side,String type){ + double sum = 0;// 土方量总和，整数填 + for(int i=0;i> col = gridPolygon.get(i); + // 逆序比较 + for(int j=col.size()-1;j>=0;j--){ + Map pfp = col.get(j); + double h = (double)pfp.get("height"); + Point p = (Point)pfp.get("geom"); + try { + double elevation = getElevation(p); + if(type.equals("-1")&&elevation>h||type.equals("1")&&elevation new Double((double) (v + (h - elevation) * Math.pow(degree2meter(side), 2)))); + sum += (h - elevation) * Math.pow(degree2meter(side), 2); + } catch (Exception e) { + throw new RuntimeException(e); + } + } + } +// System.out.println("gridPolygon sum:"+sum); + return sum; + } + // 添加一个线的点集 + public static List> addLinePoint(Coordinate[] bgcoords, GeometryFactory geometryFactory, double side, double bottomHeight,double topHeight){ + double degree = meter2degree(side); + // 计算外边界上的点 + List> lineList = new ArrayList<>(); + Coordinate startCoord = new Coordinate(bgcoords[0].x,bgcoords[0].y); + Map startPointMap = new HashMap<>(); + Point sp = geometryFactory.createPoint(startCoord); + startPointMap.put("geom",sp); + startPointMap.put("id",UUID.randomUUID().toString()); + startPointMap.put("bottomHeight",bottomHeight); + startPointMap.put("topHeight",topHeight); + lineList.add(startPointMap); + + int len = (int)Math.floor((Math.sqrt(Math.pow(bgcoords[1].y - bgcoords[0].y,2)+Math.pow(bgcoords[1].x - bgcoords[0].x,2)))/degree); + double t1 = (bgcoords[1].y - bgcoords[0].y)/(bgcoords[1].x-bgcoords[0].x); + //确定符号 + double xx = bgcoords[1].x-bgcoords[0].x; + double yy = bgcoords[1].y-bgcoords[0].y; + int xSign = xx>=0?1:-1; + int ySign = yy>=0?1:-1; + for(int i=1;i nextPointMap = new HashMap<>(); + Point nextPoint = geometryFactory.createPoint(nextCoord); + nextPointMap.put("geom",nextPoint); + nextPointMap.put("id",UUID.randomUUID().toString()); + nextPointMap.put("bottomHeight",bottomHeight); + nextPointMap.put("topHeight",topHeight); + lineList.add(nextPointMap); + } + return lineList; + } + public static List> addPointSet(Coordinate[] bgcoords, + GeometryFactory geometryFactory, + double side, + Polygon polygon, + List> pointFeaturesProperties, + Point startPoint, + double boundHeight, + String type){ + List> extendPointProperties = new ArrayList<>(); + + double t1 = (bgcoords[1].y - bgcoords[0].y)/(bgcoords[1].x-bgcoords[0].x); + //确定符号 + double xx = bgcoords[1].x-bgcoords[0].x; + double yy = bgcoords[1].y-bgcoords[0].y; + int xSign = (xx>=0&&yy>0)||(xx<0&&yy>0)?1:-1; + int ySign = (xx<0&&yy<0)||(xx<0&&yy>0)?1:-1; + + for(int i=1;i<300;i++){ + double x1 = xSign * Math.abs(Math.sin(Math.atan(t1)))*side*i + startPoint.getX(); + double y1 = ySign * Math.abs(Math.cos(Math.atan(t1)))*side*i + startPoint.getY(); + Point point = geometryFactory.createPoint(new Coordinate(x1,y1)); + boolean contains = polygon.contains(point); + if(contains){ + Map tmpPoint = new HashMap<>(); + tmpPoint.put("geom",point); + tmpPoint.put("id",UUID.randomUUID()); + double ph = boundHeight;//+(degree2meter(side))*i/2; + if(type.equals("-1")){//下挖的高程计算 + ph+=(degree2meter(side))*i/2; + }else if(type.equals("1")){// + ph-=(degree2meter(side))*i/1.5; + } + tmpPoint.put("height",ph); + pointFeaturesProperties.add(tmpPoint); + extendPointProperties.add(tmpPoint); + }else{ + break; + } + } + return extendPointProperties; + } + public static void loadDem(){ + try { + File file = new File(demFilePath); + // Hints.FORCE_LONGITUDE_FIRST_AXIS_ORDER：设置经度为第一轴顺序 + GeoTiffReader reader = new GeoTiffReader(file, new Hints(Hints.FORCE_LONGITUDE_FIRST_AXIS_ORDER, Boolean.TRUE)); + GridCoverage2D coverage = reader.read(null); + //设置坐标系 + CoordinateReferenceSystem targetCRS = CRS.decode("EPSG:4326"); + // 将 GridCoverage 进行重采样转换为另一个 CRS + CalcRingGridFix.coverage = (GridCoverage2D) Operations.DEFAULT.resample(coverage, targetCRS); + System.out.println("---------dem--------"); + }catch (Exception e){ + throw new RuntimeException("---------加载dem数据失败！---------"); + } + } + public static double getElevation(Point p) throws Exception{ + // 设置经纬度及坐标系等信息 + Position position = new Position2D(p.getY(), p.getX()); + GridGeometry2D gridGeometry = coverage.getGridGeometry(); + GridCoordinates2D gridPoint = gridGeometry.worldToGrid(position); + double[] elevation = new double[1]; + coverage.evaluate(gridPoint,elevation); + double height = elevation[0]; // 这里的高度就是查询位置的高程值 +// System.out.println("---------"+height);//---------108.0 + return height; + } + + public static double meter2degree(double meter){ + return meter/(Math.PI*6371004/180); + } + public static double degree2meter(double degree){ + return degree*(Math.PI*6371004)/180; + } +} diff --git a/src/main/java/org/example/utils/CommonResult.java b/src/main/java/org/example/utils/CommonResult.java new file mode 100644 index 0000000..d46b4aa --- /dev/null +++ b/src/main/java/org/example/utils/CommonResult.java @@ -0,0 +1,39 @@ +package org.example.utils; + +public class CommonResult { + + public boolean success; + public T result; + public String message; + + public CommonResult(boolean success,String message){ + this.success = success; + this.message = message; + } + public CommonResult(){ + + } + public boolean isSuccess() { + return success; + } + + public void setSuccess(boolean success) { + this.success = success; + } + + public T getResult() { + return result; + } + + public void setResult(T result) { + this.result = result; + } + + public String getMessage() { + return message; + } + + public void setMessage(String message) { + this.message = message; + } +} diff --git a/src/main/resources/application.yml b/src/main/resources/application.yml new file mode 100644 index 0000000..0d71399 --- /dev/null +++ b/src/main/resources/application.yml @@ -0,0 +1,4 @@ +server: + port: 8888 + servlet: + context-path: / \ No newline at end of file diff --git a/src/main/resources/opencv/opencv-490.jar b/src/main/resources/opencv/opencv-490.jar new file mode 100644 index 0000000..80b3e7f Binary files /dev/null and b/src/main/resources/opencv/opencv-490.jar differ diff --git a/src/main/resources/public/image/shishi.jpg b/src/main/resources/public/image/shishi.jpg new file mode 100644 index 0000000..e9dd2c2 Binary files /dev/null and b/src/main/resources/public/image/shishi.jpg differ diff --git a/target/classes/org/example/test/CreateContour$1.class b/target/classes/org/example/test/CreateContour$1.class index 05cd405..3d92dad 100644 Binary files a/target/classes/org/example/test/CreateContour$1.class and b/target/classes/org/example/test/CreateContour$1.class differ diff --git a/target/classes/org/example/test/CreateContour.class b/target/classes/org/example/test/CreateContour.class index 82d9411..8f95190 100644 Binary files a/target/classes/org/example/test/CreateContour.class and b/target/classes/org/example/test/CreateContour.class differ diff --git a/target/classes/org/example/test/CreateRingGrid$1.class b/target/classes/org/example/test/CreateRingGrid$1.class index 7c43a1e..9e33187 100644 Binary files a/target/classes/org/example/test/CreateRingGrid$1.class and b/target/classes/org/example/test/CreateRingGrid$1.class differ diff --git a/target/classes/org/example/test/CreateRingGrid.class b/target/classes/org/example/test/CreateRingGrid.class index ab8152a..07c4f4f 100644 Binary files a/target/classes/org/example/test/CreateRingGrid.class and b/target/classes/org/example/test/CreateRingGrid.class differ diff --git a/target/classes/org/example/test/CreateTiff.class b/target/classes/org/example/test/CreateTiff.class index f850852..af6d885 100644 Binary files a/target/classes/org/example/test/CreateTiff.class and b/target/classes/org/example/test/CreateTiff.class differ diff --git a/target/classes/tmp_dem.tiff b/target/classes/tmp_dem.tiff deleted file mode 100644 index 1c234eb..0000000 Binary files a/target/classes/tmp_dem.tiff and /dev/null differ diff --git a/target/classes/zeroPoint.dbf b/target/classes/zeroPoint.dbf deleted file mode 100644 index 52ee84f..0000000 Binary files a/target/classes/zeroPoint.dbf and /dev/null differ diff --git a/target/classes/zeroPoint.fix b/target/classes/zeroPoint.fix deleted file mode 100644 index d758f18..0000000 Binary files a/target/classes/zeroPoint.fix and /dev/null differ diff --git a/target/classes/zeroPoint.prj b/target/classes/zeroPoint.prj deleted file mode 100644 index 5307c8f..0000000 --- a/target/classes/zeroPoint.prj +++ /dev/null @@ -1 +0,0 @@ -GEOGCS["WGS84(DD)", DATUM["WGS84", SPHEROID["WGS84", 6378137.0, 298.257223563]], PRIMEM["Greenwich", 0.0], UNIT["degree", 0.017453292519943295], AXIS["Geodetic longitude", EAST], AXIS["Geodetic latitude", NORTH], AUTHORITY["EPSG","4326"]] \ No newline at end of file diff --git a/target/classes/zeroPoint.shp b/target/classes/zeroPoint.shp deleted file mode 100644 index 3939d7b..0000000 Binary files a/target/classes/zeroPoint.shp and /dev/null differ diff --git a/target/classes/zeroPoint.shx b/target/classes/zeroPoint.shx deleted file mode 100644 index e7fe7e8..0000000 Binary files a/target/classes/zeroPoint.shx and /dev/null differ