添加gdal测试
This commit is contained in:
user
parent
fcae99253e
commit
0dc3298069
5
pom.xml
5
pom.xml
|
|
@ -93,6 +93,11 @@
|
|||
<artifactId>gt-opengis</artifactId>
|
||||
<version>29.5</version>
|
||||
</dependency>
|
||||
<dependency>
|
||||
<groupId>org.geotools</groupId>
|
||||
<artifactId>gt-geotiff</artifactId>
|
||||
<version>${geotools.version}</version>
|
||||
</dependency>
|
||||
<!-- https://mvnrepository.com/artifact/org.gdal/gdal -->
|
||||
<dependency>
|
||||
<groupId>org.gdal</groupId>
|
||||
|
|
|
|||
|
|
@ -0,0 +1,288 @@
|
|||
package org.example.test;
|
||||
|
||||
import it.geosolutions.jaiext.range.Range;
|
||||
import org.gdal.gdal.Band;
|
||||
import org.gdal.gdal.Dataset;
|
||||
import org.gdal.gdal.gdal;
|
||||
import org.gdal.gdalconst.gdalconst;
|
||||
import org.gdal.ogr.Layer;
|
||||
import org.gdal.osr.SpatialReference;
|
||||
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.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.data.DefaultTransaction;
|
||||
import org.geotools.data.collection.ListFeatureCollection;
|
||||
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.GeoTiffReader;
|
||||
import org.geotools.geometry.Position2D;
|
||||
import org.geotools.geometry.jts.ReferencedEnvelope;
|
||||
import org.geotools.referencing.CRS;
|
||||
import org.geotools.referencing.crs.DefaultGeographicCRS;
|
||||
import org.geotools.util.factory.Hints;
|
||||
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.opengis.geometry.Envelope;
|
||||
|
||||
import java.io.File;
|
||||
import java.io.IOException;
|
||||
import java.io.Serializable;
|
||||
import java.nio.charset.Charset;
|
||||
import java.util.*;
|
||||
|
||||
public class CreateContour {
|
||||
public static void main(String[] args) throws Exception{
|
||||
loadDem();
|
||||
// 假设我们使用WGS84坐标系统
|
||||
CoordinateReferenceSystem crs = CRS.decode("EPSG:4326");
|
||||
GeometryFactory geometryFactory = new GeometryFactory();
|
||||
// 定义矩形的左下角和右上角坐标（WGS84）
|
||||
List<double[]> 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 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};
|
||||
// PrimitiveIterator.OfDouble iterator = Arrays.stream(stepSet).iterator();
|
||||
// double nextHeight=targetHeight,step = iterator.next();
|
||||
// while (Math.abs(allSum)-100>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);
|
||||
|
||||
}
|
||||
public static double calcBlanceHeight(Coordinate[] coordinates,GeometryFactory geometryFactory,double height) throws Exception{
|
||||
double bufferLen = 100, stepLen = 2;
|
||||
|
||||
List<SimpleFeature> featuresList = new ArrayList<>();
|
||||
SimpleFeatureType TYPE = createFeatureType();
|
||||
SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE);
|
||||
|
||||
Polygon rectangle = geometryFactory.createPolygon(coordinates);
|
||||
LineString boundary = (LineString)rectangle.getBoundary();
|
||||
featureBuilder.add(geometryFactory.createLineString(boundary.getCoordinates()));
|
||||
featureBuilder.add("目标平面的第一条等高线");
|
||||
featureBuilder.add(height);
|
||||
//buildFeature(id)创建一个要素,ID可以是null,当ID为null时会由builder内部生成
|
||||
SimpleFeature feature = featureBuilder.buildFeature("1");
|
||||
featuresList.add(feature);
|
||||
|
||||
// 计算放坡的等高线
|
||||
for(int i=1;i<=bufferLen/stepLen;i++){
|
||||
BufferParameters bufferParameters = new BufferParameters();
|
||||
bufferParameters.setEndCapStyle(BufferParameters.CAP_FLAT);
|
||||
bufferParameters.setJoinStyle(BufferParameters.JOIN_MITRE);
|
||||
BufferOp bufOp = new BufferOp(rectangle,bufferParameters);
|
||||
|
||||
double degree = meter2degree(stepLen*i);
|
||||
Polygon bg = (Polygon)bufOp.getResultGeometry(degree);
|
||||
LineString lineString = (LineString)bg.getBoundary();
|
||||
featureBuilder.add(lineString);
|
||||
featureBuilder.add(("扩展平面的第"+i+"条等高线"));
|
||||
featureBuilder.add(height-i*stepLen/2);
|
||||
//buildFeature(id)创建一个要素,ID可以是null,当ID为null时会由builder内部生成
|
||||
SimpleFeature simpleFeature = featureBuilder.buildFeature("ext"+i);
|
||||
featuresList.add(simpleFeature);
|
||||
}
|
||||
SimpleFeatureCollection collection = new ListFeatureCollection(TYPE, featuresList);
|
||||
ReferencedEnvelope bounds = collection.getBounds();
|
||||
|
||||
ShapefileDataStore newDataStore = getShapefileDataStore(collection);
|
||||
newDataStore.createSchema(TYPE);
|
||||
newDataStore.setCharset(Charset.forName("UTF-8"));
|
||||
Transaction transaction = new DefaultTransaction("create");
|
||||
|
||||
String typeName = newDataStore.getTypeNames()[0];
|
||||
SimpleFeatureSource featureSource = newDataStore.getFeatureSource(typeName);
|
||||
SimpleFeatureType SHAPE_TYPE = featureSource.getSchema();
|
||||
/*
|
||||
* The Shapefile format has a couple limitations:
|
||||
* - "the_geom" is always first, and used for the geometry attribute name
|
||||
* - "the_geom" must be of type Point, MultiPoint, MuiltiLineString, MultiPolygon
|
||||
* - Attribute names are limited in length
|
||||
* - Not all data types are supported (example Timestamp represented as Date)
|
||||
*
|
||||
* Each data store has different limitations so check the resulting SimpleFeatureType.
|
||||
*/
|
||||
System.out.println("SHAPE:" + SHAPE_TYPE);
|
||||
|
||||
if (featureSource instanceof SimpleFeatureStore) {
|
||||
SimpleFeatureStore featureStore = (SimpleFeatureStore) featureSource;
|
||||
/*
|
||||
* SimpleFeatureStore has a method to add features from a
|
||||
* SimpleFeatureCollection object, so we use the ListFeatureCollection
|
||||
* class to wrap our list of features.
|
||||
*/
|
||||
|
||||
featureStore.setTransaction(transaction);
|
||||
try {
|
||||
featureStore.addFeatures(collection);
|
||||
transaction.commit();
|
||||
} catch (Exception problem) {
|
||||
problem.printStackTrace();
|
||||
transaction.rollback();
|
||||
} finally {
|
||||
transaction.close();
|
||||
}
|
||||
// System.exit(0); // success!
|
||||
} else {
|
||||
System.out.println(typeName + " does not support read/write access");
|
||||
System.exit(1);
|
||||
}
|
||||
|
||||
// 注册所有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("矢量数据没有图层！");
|
||||
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<Object> 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;
|
||||
}
|
||||
|
||||
private static ShapefileDataStore getShapefileDataStore(SimpleFeatureCollection collection) throws IOException {
|
||||
// 创建输出路径
|
||||
File newFile = new File("E:/code/CutFillDemo/target/classes/tmp.shp");
|
||||
// 创建shapefile类型的数据存储工厂实例
|
||||
ShapefileDataStoreFactory dataStoreFactory = new ShapefileDataStoreFactory();
|
||||
Map<String, Serializable> params = new HashMap<>();
|
||||
params.put("url", newFile.toURI().toURL());
|
||||
params.put("create spatial index", Boolean.TRUE);
|
||||
// 创建一个空的数据存储
|
||||
ShapefileDataStore newDataStore = (ShapefileDataStore) dataStoreFactory.createNewDataStore(params);
|
||||
return newDataStore;
|
||||
}
|
||||
|
||||
// public static
|
||||
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
|
||||
CreateGrid.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 = CreateGrid.coverage.getGridGeometry();
|
||||
GridCoordinates2D gridPoint = gridGeometry.worldToGrid(position);
|
||||
double[] elevation = new double[1];
|
||||
CreateGrid.coverage.evaluate(gridPoint,elevation);
|
||||
double height = elevation[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;
|
||||
}
|
||||
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", LineString.class);
|
||||
builder.length(50).add("name", String.class); // <- 15 chars width for name field
|
||||
builder.add("ELEV", Double.class);
|
||||
|
||||
// build the type
|
||||
final SimpleFeatureType LOCATION = builder.buildFeatureType();
|
||||
|
||||
return LOCATION;
|
||||
}
|
||||
}
|
||||
|
|
@ -70,7 +70,7 @@ public class CreateGrid {
|
|||
loadDem();
|
||||
// 假设我们使用WGS84坐标系统
|
||||
CoordinateReferenceSystem crs = CRS.decode("EPSG:4326");
|
||||
|
||||
GeometryFactory geometryFactory = new GeometryFactory();
|
||||
// 定义矩形的左下角和右上角坐标（WGS84）
|
||||
List<double[]> coords = new ArrayList(){{
|
||||
add(new double[]{113.52896835361318, 23.654018962759377});
|
||||
|
|
@ -88,8 +88,41 @@ public class CreateGrid {
|
|||
// 闭合多边形，复制第一个坐标作为最后一个坐标
|
||||
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 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};
|
||||
PrimitiveIterator.OfDouble iterator = Arrays.stream(stepSet).iterator();
|
||||
double nextHeight=targetHeight,step = iterator.next();
|
||||
while (Math.abs(allSum)-100>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);
|
||||
}
|
||||
public static double calcBlanceHeight(Coordinate[] coordinates,GeometryFactory geometryFactory,double targetHeight) throws Exception{
|
||||
// 创建GeometryFactory
|
||||
GeometryFactory geometryFactory = new GeometryFactory();
|
||||
|
||||
List<Map<String,Object>> featuresProperties = new ArrayList<>();// 面要素属性集合
|
||||
// List<Map<String,Object>> pointFeaturesProperties = new ArrayList<>();// 扩展点属性集合
|
||||
List<Map<String,Object>> tPointFeaturesProperties = new ArrayList<>();// 目标平面的点集属性集合
|
||||
|
|
@ -98,6 +131,11 @@ public class CreateGrid {
|
|||
List<Map<String,Object>> p3PointFeaturesProperties = new ArrayList<>();// p3面内的扩展点属性集合
|
||||
List<Map<String,Object>> p4PointFeaturesProperties = new ArrayList<>();// p4面内的扩展点属性集合
|
||||
|
||||
List<Map<String,Object>> p1TopPointFeaturesProperties = new ArrayList<>();// p1面内的扩展点属性集合
|
||||
List<Map<String,Object>> p2TopPointFeaturesProperties = new ArrayList<>();// p2面内的扩展点属性集合
|
||||
List<Map<String,Object>> p3TopPointFeaturesProperties = new ArrayList<>();// p3面内的扩展点属性集合
|
||||
List<Map<String,Object>> p4TopPointFeaturesProperties = new ArrayList<>();// p4面内的扩展点属性集合
|
||||
|
||||
// 创建多边形（矩形）
|
||||
Map<String,Object> innerPolygon = new HashMap<>();
|
||||
Polygon rectangle = geometryFactory.createPolygon(coordinates);
|
||||
|
|
@ -106,15 +144,6 @@ public class CreateGrid {
|
|||
innerPolygon.put("id",1);
|
||||
featuresProperties.add(innerPolygon);
|
||||
|
||||
// 计算设计标高，取四个角的高程平均值
|
||||
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);
|
||||
|
||||
// 缓冲区建立，外扩100m
|
||||
double meter = 100;
|
||||
|
|
@ -178,85 +207,105 @@ public class CreateGrid {
|
|||
|
||||
// 打点
|
||||
double outBottomHeight = targetHeight - meter/2;
|
||||
System.out.println("边界高程："+outBottomHeight);
|
||||
double outTopHeight = targetHeight + meter/1.5;
|
||||
|
||||
// System.out.println("边界高程："+outBottomHeight);
|
||||
Coordinate[] coordArr1 = new Coordinate[]{bgcoords[0],bgcoords[1]};
|
||||
List gridsP1 = addPartPoint(coordArr1,geometryFactory,side,p1,p1PointFeaturesProperties,outBottomHeight,"-1");
|
||||
Coordinate[] coordArr2 = new Coordinate[]{bgcoords[1],bgcoords[2]};
|
||||
List gridsP2 = addPartPoint(coordArr2,geometryFactory,side,p2,p2PointFeaturesProperties,outBottomHeight,"-1");
|
||||
Coordinate[] coordArr3 = new Coordinate[]{bgcoords[2],bgcoords[3]};
|
||||
List gridsP3 = addPartPoint(coordArr3,geometryFactory,side,p3,p3PointFeaturesProperties,outBottomHeight,"-1");
|
||||
Coordinate[] coordArr4 = new Coordinate[]{bgcoords[3],bgcoords[0]};
|
||||
|
||||
List gridsP1 = addPartPoint(coordArr1,geometryFactory,side,p1,p1PointFeaturesProperties,outBottomHeight,"-1");
|
||||
List gridsP2 = addPartPoint(coordArr2,geometryFactory,side,p2,p2PointFeaturesProperties,outBottomHeight,"-1");
|
||||
List gridsP3 = addPartPoint(coordArr3,geometryFactory,side,p3,p3PointFeaturesProperties,outBottomHeight,"-1");
|
||||
List gridsP4 = addPartPoint(coordArr4,geometryFactory,side,p4,p4PointFeaturesProperties,outBottomHeight,"-1");
|
||||
|
||||
List gridsTopP1 = addPartPoint(coordArr1,geometryFactory,side,p1,p1TopPointFeaturesProperties,outTopHeight,"1");
|
||||
List gridsTopP2 = addPartPoint(coordArr2,geometryFactory,side,p2,p2TopPointFeaturesProperties,outTopHeight,"1");
|
||||
List gridsTopP3 = addPartPoint(coordArr3,geometryFactory,side,p3,p3TopPointFeaturesProperties,outTopHeight,"1");
|
||||
List gridsTopP4 = addPartPoint(coordArr4,geometryFactory,side,p4,p4TopPointFeaturesProperties,outTopHeight,"1");
|
||||
|
||||
|
||||
// 计算，整数为填，负数为挖
|
||||
targetPanelCalc(tPointFeaturesProperties,side);// 计算矩形的填挖量
|
||||
double p1Result = extPanelCalc(gridsP1,side);
|
||||
double p2Result = extPanelCalc(gridsP2,side);
|
||||
double p3Result = extPanelCalc(gridsP3,side);
|
||||
double p4Result = extPanelCalc(gridsP4,side);
|
||||
double p1Result = extPanelCalc(gridsP1,side,"-1");
|
||||
double p2Result = extPanelCalc(gridsP2,side,"-1");
|
||||
double p3Result = extPanelCalc(gridsP3,side,"-1");
|
||||
double p4Result = extPanelCalc(gridsP4,side,"-1");
|
||||
|
||||
double p1TopResult = extPanelCalc(gridsTopP1,side,"1");
|
||||
double p2TopResult = extPanelCalc(gridsTopP2,side,"1");
|
||||
double p3TopResult = extPanelCalc(gridsTopP3,side,"1");
|
||||
double p4TopResult = extPanelCalc(gridsTopP4,side,"1");
|
||||
|
||||
|
||||
|
||||
|
||||
// // 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<List<Map<String,Object>>> pointList = new ArrayList(){{
|
||||
//// add(pointFeaturesProperties);
|
||||
// add(tPointFeaturesProperties);
|
||||
// add(p1TopPointFeaturesProperties);
|
||||
// add(p2PointFeaturesProperties);
|
||||
// add(p3PointFeaturesProperties);
|
||||
// add(p4PointFeaturesProperties);
|
||||
// }};
|
||||
// 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);
|
||||
// Iterator<SimpleFeature> iterator = featureCollection.stream().iterator();
|
||||
// iterator.next();
|
||||
// collection.add(iterator.next());
|
||||
// SimpleFeatureSource source = new CollectionFeatureSource(collection);
|
||||
// IntersectionBuilder intersectionBuilder = new IntersectionBuilder(TYPE, source);
|
||||
// SimpleFeatureSource grid = Grids.createSquareGrid(collection.getBounds(), side, -1);
|
||||
// Style styleBox = SLD.createSimpleStyle(grid.getSchema());
|
||||
// Layer layerBox = new FeatureLayer(grid, styleBox);
|
||||
//// map.addLayer(layerBox);
|
||||
//
|
||||
// // Now display the map
|
||||
// JMapFrame.showMap(map);
|
||||
|
||||
double sum = (p1Result+p2Result+p3Result+p4Result);
|
||||
System.out.println("allSum:"+sum);
|
||||
// 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<List<Map<String,Object>>> pointList = new ArrayList(){{
|
||||
// add(pointFeaturesProperties);
|
||||
add(tPointFeaturesProperties);
|
||||
add(p1PointFeaturesProperties);
|
||||
add(p2PointFeaturesProperties);
|
||||
add(p3PointFeaturesProperties);
|
||||
add(p4PointFeaturesProperties);
|
||||
}};
|
||||
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);
|
||||
Iterator<SimpleFeature> iterator = featureCollection.stream().iterator();
|
||||
iterator.next();
|
||||
collection.add(iterator.next());
|
||||
SimpleFeatureSource source = new CollectionFeatureSource(collection);
|
||||
IntersectionBuilder intersectionBuilder = new IntersectionBuilder(TYPE, source);
|
||||
SimpleFeatureSource grid = Grids.createSquareGrid(collection.getBounds(), side, -1);
|
||||
Style styleBox = SLD.createSimpleStyle(grid.getSchema());
|
||||
Layer layerBox = new FeatureLayer(grid, styleBox);
|
||||
// map.addLayer(layerBox);
|
||||
|
||||
// Now display the map
|
||||
JMapFrame.showMap(map);
|
||||
|
||||
double sumTop = (p1TopResult+p2TopResult+p3TopResult+p4TopResult);
|
||||
// System.out.println("allSum:"+(sum+sumTop));
|
||||
return sum+sumTop;
|
||||
}
|
||||
|
||||
public static void targetPanelCalc (List<Map<String,Object>> pointFeaturesProperties,double side){
|
||||
|
|
@ -272,10 +321,9 @@ public class CreateGrid {
|
|||
throw new RuntimeException(e);
|
||||
}
|
||||
});
|
||||
System.out.println("sum:"+sum);
|
||||
}
|
||||
|
||||
public static double extPanelCalc(List<List<Map<String,Object>>> gridPolygon,double side){
|
||||
public static double extPanelCalc(List<List<Map<String,Object>>> gridPolygon,double side,String type){
|
||||
double sum = 0;// 土方量总和，整数填
|
||||
for(int i=0;i<gridPolygon.size();i++){
|
||||
List<Map<String,Object>> col = gridPolygon.get(i);
|
||||
|
|
@ -286,7 +334,7 @@ public class CreateGrid {
|
|||
Point p = (Point)pfp.get("geom");
|
||||
try {
|
||||
double elevation = getElevation(p);
|
||||
if(elevation>h){
|
||||
if(type.equals("-1")&&elevation>h||type.equals("1")&&elevation<h){
|
||||
break;
|
||||
}
|
||||
// sum.updateAndGet(v -> new Double((double) (v + (h - elevation) * Math.pow(degree2meter(side), 2))));
|
||||
|
|
@ -296,7 +344,7 @@ public class CreateGrid {
|
|||
}
|
||||
}
|
||||
}
|
||||
System.out.println("gridPolygon sum:"+sum);
|
||||
// System.out.println("gridPolygon sum:"+sum);
|
||||
return sum;
|
||||
}
|
||||
// 添加一个面的点集
|
||||
|
|
@ -388,7 +436,7 @@ public class CreateGrid {
|
|||
if(type.equals("-1")){//下挖的高程计算
|
||||
ph+=(degree2meter(side))*i/2;
|
||||
}else if(type.equals("1")){//
|
||||
// ph-=(degree2meter(side))*i/2;
|
||||
ph-=(degree2meter(side))*i/1.5;
|
||||
}
|
||||
tmpPoint.put("height",ph);
|
||||
pointFeaturesProperties.add(tmpPoint);
|
||||
|
|
|
|||
|
|
@ -0,0 +1,158 @@
|
|||
package org.example.test;
|
||||
|
||||
import it.geosolutions.jaiext.range.Range;
|
||||
import org.gdal.gdal.gdal;
|
||||
import org.geotools.api.coverage.grid.GridCoverageWriter;
|
||||
import org.geotools.api.geometry.Bounds;
|
||||
import org.geotools.api.parameter.GeneralParameterValue;
|
||||
import org.geotools.api.parameter.ParameterValueGroup;
|
||||
import org.geotools.api.referencing.FactoryException;
|
||||
import org.geotools.api.referencing.crs.CoordinateReferenceSystem;
|
||||
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.gce.geotiff.GeoTiffFormat;
|
||||
import org.geotools.gce.geotiff.GeoTiffReader;
|
||||
import org.geotools.gce.geotiff.GeoTiffWriteParams;
|
||||
import org.geotools.geometry.Envelope2DArchived;
|
||||
import org.geotools.geometry.jts.ReferencedEnvelope;
|
||||
import org.geotools.referencing.CRS;
|
||||
import org.geotools.util.factory.Hints;
|
||||
|
||||
import javax.imageio.ImageIO;
|
||||
import javax.imageio.stream.ImageInputStream;
|
||||
import javax.media.jai.RasterFactory;
|
||||
import java.awt.*;
|
||||
import java.awt.image.Raster;
|
||||
import java.awt.image.RenderedImage;
|
||||
import java.awt.image.WritableRaster;
|
||||
import java.io.*;
|
||||
import java.util.Random;
|
||||
import java.util.UUID;
|
||||
|
||||
public class CreateTiff {
|
||||
public static void main(String[] args) throws Exception{
|
||||
tiffWrite();
|
||||
|
||||
// tiffRead();
|
||||
}
|
||||
public static String tiffWrite() throws FactoryException, IOException {
|
||||
int[][] matrix = new int[256][256];
|
||||
Random r = new Random();
|
||||
for (int i = 0; i < matrix.length; i++) {
|
||||
int[] row = matrix[i];
|
||||
for (int j = 0; j < row.length; j++) {
|
||||
matrix[i][j] = Math.abs(r.nextInt() % Byte.MAX_VALUE);
|
||||
//System.out.println(i + "," + j + "," + matrix[i][j]);
|
||||
}
|
||||
}
|
||||
double tile = 6378137 * 2 * Math.PI / (1 << 16);
|
||||
Bounds envelope = new Envelope2DArchived( CRS.decode("EPSG:3857"),12923543.888468, 4852595.145082 - tile, tile, tile);
|
||||
|
||||
WritableRaster raster = createWritableRaster(1);
|
||||
addBound(matrix, raster, 0);
|
||||
// addBound(matrix, raster, 1);
|
||||
// addBound(matrix, raster, 2);
|
||||
// addBound(matrix, raster, 3);
|
||||
GridCoverageFactory f = new GridCoverageFactory();
|
||||
GridCoverage2D coverage = f.create("test1", raster, envelope);
|
||||
GeoTiffFormat format = new GeoTiffFormat();
|
||||
|
||||
/*ByteArrayOutputStream bao = new ByteArrayOutputStream();
|
||||
ImageOutputStream ios = ImageIO.createImageOutputStream(bao);
|
||||
GridCoverageWriter writer = format.getWriter(ios);*/
|
||||
|
||||
File out = new File("E:\\data\\test1.tiff");
|
||||
GridCoverageWriter writer = format.getWriter(out);
|
||||
writer.write(coverage, createWriteParameters());
|
||||
writer.dispose();
|
||||
coverage.dispose(true);
|
||||
|
||||
// int i = gdal.RasterizeLayer();
|
||||
|
||||
|
||||
return "test";
|
||||
}
|
||||
|
||||
public static WritableRaster createWritableRaster(int boundNum) {
|
||||
WritableRaster raster = RasterFactory.createBandedRaster(2, 256, 256, boundNum, (Point) null);
|
||||
return raster;
|
||||
}
|
||||
|
||||
public static void addBound(int[][] matrix, WritableRaster raster, int boundIndex) {
|
||||
for (int i = 0; i < 256; i++) {
|
||||
for (int j = 0; j < 256; j++) {
|
||||
raster.setSample(i, j, boundIndex, matrix[i][j]);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
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;
|
||||
}
|
||||
public static void tiffRead(byte[] tiff) throws IOException, FactoryException {
|
||||
|
||||
ByteArrayInputStream arrayInputStream = new ByteArrayInputStream(tiff);
|
||||
ImageInputStream stream = ImageIO.createImageInputStream(arrayInputStream);
|
||||
File file = new File("src\\main\\resources\\public\\data\\test1.tiff");
|
||||
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
|
||||
// coverage = (GridCoverage2D) Operations.DEFAULT.resample(coverage, targetCRS);
|
||||
|
||||
ReferencedEnvelope envelope = coverage.getEnvelope2D();
|
||||
RenderedImage image = coverage.getRenderedImage();
|
||||
|
||||
double height = envelope.getHeight();
|
||||
double width = envelope.getWidth();
|
||||
double minx = envelope.getMinX();
|
||||
double miny = envelope.getMinY();
|
||||
String crs = CRS.toSRS(envelope.getCoordinateReferenceSystem());
|
||||
int type = image.getSampleModel().getDataType();
|
||||
|
||||
int bandNum = image.getSampleModel().getNumBands();
|
||||
int matrixW = image.getWidth();
|
||||
//int matrixW2 = image.getTileWidth();
|
||||
int matrixH = image.getHeight();
|
||||
//int matrixH2 = image.getTileHeight();
|
||||
|
||||
int[][][] matrix = new int[bandNum][matrixW][matrixH];
|
||||
Raster raster = image.getData();
|
||||
for (int b = 0; b < bandNum; b++) {
|
||||
for (int i = 0; i < matrixW; i++) {
|
||||
for (int j = 0; j < matrixH; j++) {
|
||||
int value = raster.getSample(i, j, b);
|
||||
matrix[b][i][j] = value;
|
||||
if (value > 0){
|
||||
System.out.println(i + "," + j + "," + b + "," + value);
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
//tiffWrite2(crs,minx,miny,width,height,type,bandNum,matrixW,matrixH,matrix);
|
||||
// tiffWrite3(coverage);
|
||||
}
|
||||
|
||||
public static void tiffWrite3(GridCoverage2D coverage) throws IOException {
|
||||
File out2 = new File("/Users/didi/Desktop/data/"+ UUID.randomUUID().toString().substring(0,4)+".tiff");
|
||||
GeoTiffFormat format = new GeoTiffFormat();
|
||||
GridCoverageWriter writer2 = format.getWriter(out2);
|
||||
writer2.write(coverage, createWriteParameters());
|
||||
writer2.dispose();
|
||||
coverage.dispose(true);
|
||||
}
|
||||
}
|
||||
|
|
@ -0,0 +1,45 @@
|
|||
package org.example.test;
|
||||
|
||||
import java.io.BufferedReader;
|
||||
import java.io.InputStreamReader;
|
||||
|
||||
public class ExeMain {
|
||||
public static void main(String[] args) {
|
||||
//需要传入的参数（ps: 多个参数用空格隔开）
|
||||
// String paras = " --task C:\\Users\\Admin\\.cesiumlab\\tasks\\terrain_99254c80b1ef11ea8edfe728504fc83a.json --taskserver tcp://127.0.0.1:9001 --taskname 99254c80b1ef11ea8edfe728504fc83a --log_dir C:\\Users\\Admin\\.cesiumlab\\logs";
|
||||
String paras = "";
|
||||
//调用的exe可执行文件（ps: 调用可执行文件和参数拼接必须要用空格隔开）
|
||||
String cmd = "E:\\code\\contour2dem.gui\\contour2dem.gui.exe" + paras;
|
||||
|
||||
openExe(cmd);
|
||||
|
||||
}
|
||||
public static void openExe(String cmd) {
|
||||
BufferedReader br = null;
|
||||
BufferedReader brError = null;
|
||||
|
||||
try {
|
||||
//执行exe cmd可以为字符串(exe存放路径)也可为数组，调用exe时需要传入参数时，可以传数组调用(参数有顺序要求)
|
||||
Process p = Runtime.getRuntime().exec(cmd);
|
||||
String line = null;
|
||||
//获得子进程的输入流。
|
||||
br = new BufferedReader(new InputStreamReader(p.getInputStream()));
|
||||
//获得子进程的错误流。
|
||||
brError = new BufferedReader(new InputStreamReader(p.getErrorStream()));
|
||||
while ((line = br.readLine()) != null || (line = brError.readLine()) != null) {
|
||||
//输出exe输出的信息以及错误信息
|
||||
System.out.println(line);
|
||||
}
|
||||
} catch (Exception e) {
|
||||
e.printStackTrace();
|
||||
} finally {
|
||||
if (br != null) {
|
||||
try {
|
||||
br.close();
|
||||
} catch (Exception e) {
|
||||
e.printStackTrace();
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
|
@ -50,7 +50,7 @@ public class Lession02_Feature {
|
|||
final SimpleFeatureType TYPE =
|
||||
DataUtilities.createType(
|
||||
"Location",
|
||||
"the_geom:Point:srid=4326,"
|
||||
"the_geom:LineString:srid=4326,"
|
||||
+ // <- the geometry attribute: Point type
|
||||
"name:String,"
|
||||
+ // <- a String attribute
|
||||
|
|
|
|||
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
Binary file not shown.
|
|
@ -0,0 +1 @@
|
|||
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"]]
|
||||
Binary file not shown.
Binary file not shown.
Binary file not shown.
|
|
@ -0,0 +1,11 @@
|
|||
<PAMDataset>
|
||||
<PAMRasterBand band="1">
|
||||
<Metadata>
|
||||
<MDI key="STATISTICS_MAXIMUM">160</MDI>
|
||||
<MDI key="STATISTICS_MEAN">126.01388888889</MDI>
|
||||
<MDI key="STATISTICS_MINIMUM">113</MDI>
|
||||
<MDI key="STATISTICS_STDDEV">13.012273112074</MDI>
|
||||
<MDI key="STATISTICS_VALID_PERCENT">100</MDI>
|
||||
</Metadata>
|
||||
</PAMRasterBand>
|
||||
</PAMDataset>
|
||||
Loading…
Reference in New Issue