java使用GeoTools,JTS实现叠加分析、提取分析、近邻分析

最新推荐文章于 2024-08-13 17:33:49 发布

weixin_40538248

最新推荐文章于 2024-08-13 17:33:49 发布

阅读量1.1k

点赞数 6

文章标签： java 开发语言 arcgis

本文链接：https://blog.csdn.net/weixin_40538248/article/details/138892132

版权

java使用GeoTools,JTS实现叠加分析、提取分析、近邻分析

简介

概念：
JTS：JTS Topology Suite (JTS)是一个开源的Java软件库，它提供了平面几何的对象模型和基本的几何函数，符合OGC发布的“Simple Features for SQL”（SFSQL）规范。JTS被设计用作基于矢量地理信息软件的核心组件，还可以用作计算几何的通用算法库。

GeoTools：GeoTools 是一个开源 (LGPL) Java 代码库，它为操作地理空间数据提供符合标准的方法，例如实现地理信息系统。GeoTools 库数据结构基于开放地理空间联盟 (OGC) 规范。
官网地址：
https://www.geotools.org/

你好！这是你第一次使用 Markdown编辑器 所展示的欢迎页。如果你想学习如何使用Markdown编辑器, 可以仔细阅读这篇文章，了解一下Markdown的基本语法知识。

数据操作

**shape：

public static void readShpFile(String shpPath) {
       File shpFile = new File(shpPath);
       try {
           ShapefileDataStore shapefileDataStore = new ShapefileDataStore(shpFile.toURI().toURL());
           // 设置编码,防止属性的中文字符出现乱码
           shapefileDataStore.setCharset(Charset.forName("UTF-8"));
           // 这个typeNamae不传递，默认是文件名称
           FeatureSource featuresource = shapefileDataStore.getFeatureSource(shapefileDataStore.getTypeNames()[0]);
           // 读取bbox
           ReferencedEnvelope bbox  =featuresource.getBounds();
           // 读取投影
           CoordinateReferenceSystem crs = featuresource.getSchema().getCoordinateReferenceSystem();
           // 特征总数
           int count = featuresource.getCount(Query.ALL);
           // 获取当前数据的geometry类型（点、线、面）
           GeometryType geometryType = featuresource.getSchema().getGeometryDescriptor().getType();
           // 读取要素
           SimpleFeatureCollection simpleFeatureCollection = (SimpleFeatureCollection) featuresource.getFeatures();
           // 获取当前矢量数据有哪些属性字段值
           List<AttributeDescriptor> attributes = simpleFeatureCollection.getSchema().getAttributeDescriptors();
           // 拿到迭代器
           SimpleFeatureIterator simpleFeatureIterator = simpleFeatureCollection.features();
           // 遍历每一个要素
           while(simpleFeatureIterator.hasNext()) {
               SimpleFeature simpleFeature = simpleFeatureIterator.next();
               // java8新特性流api
               attributes.stream().forEach((a) -> {
                   // 依次读取这个shape中每一个属性值，当然这个属性值，可以处理其它业务
                   System.out.println(a.getLocalName() + ":" + simpleFeature.getAttribute(a.getLocalName()));
               });
           }
       } catch (IOException e) {
           e.printStackTrace();
       }
       System.out.println("读取完成！");
   }

2. Geoserver服务数据读取

    public static final String getFeatureApi = "{{url}}/geoserver/wfs?service=wfs&version=2.0.0&outputFormat=application/json&request=GetFeature&typeNames={{workspace}}:{{layerName}}";
        /**
     * 获取SimpleFeatureCollection
     */
    public static SimpleFeatureCollection getFeature(
            String url,
            String username,
            String password,
            String workspace,
            String layerName) {
        try {
            String urlStr = buildUrl(getFeatureApi, url, workspace, layerName);
            CredentialsProvider credentialsProvider = new BasicCredentialsProvider();
            credentialsProvider.setCredentials(AuthScope.ANY, new UsernamePasswordCredentials(username, password));
            CloseableHttpClient httpClient = HttpClients.custom().setDefaultCredentialsProvider(credentialsProvider).build();
            HttpGet httpGet = new HttpGet(urlStr);
            HttpResponse response = httpClient.execute(httpGet);
            String jsonStr = EntityUtils.toString(response.getEntity());
            ByteArrayInputStream inputStream = new ByteArrayInputStream(jsonStr.getBytes(StandardCharsets.UTF_8));
            return (SimpleFeatureCollection) new FeatureJSON().readFeatureCollection(inputStream);
        } catch (Exception e) {
            e.printStackTrace();
            return null;
        }
    }

3. wkt数据读取

// ex : MultiPolygon ( ((10 10, 10 20, 20 20, 20 15, 10 10)), ((60 60, 70 70, 80 60, 60 60 )) )
WKTReader wktReader = new WKTReader();
Geometry geometry = wktReader.read(wktStr);

4. geojson读取

//{“type”:“Feature”,“crs”:{“type”:“name”,“properties”:{“name”:“EPSG:2380”}},“geometry”:{“type”:“MultiPolygon”,“coordinates”:[[[[646398.9535,3267941.9664],[649558.7196,3267895.3528],[649674.763,3265683.4124],[646387.8773,3265827.4858],[646398.9535,3267941.9664]]]]},“properties”:{“Id”:0}}
public static void readGeoJson(String jsonPath) throws IOException {
       File file = new File(jsonPath);
       FileInputStream fileInputStream = new FileInputStream(file);
       // 这里可能存在问题，如果是超大文件怎么办，一次性读入会不会报内存
       // 解决方案是不是对大文件进行拆分
       GeoJSONReader geoJSONReader = new GeoJSONReader(fileInputStream);
       SimpleFeatureCollection featureCollection = geoJSONReader.getFeatures();
       SimpleFeatureIterator iterator = featureCollection.features();
       List<AttributeDescriptor> attributes = featureCollection.getSchema().getAttributeDescriptors();
       while (iterator.hasNext()) {
           SimpleFeature simpleFeature = iterator.next();
           System.out.println();
           attributes.stream().forEach((a) -> {
               // 依次读取这个shape中每一个属性值，当然这个属性值，可以处理其它业务
               System.out.println(a.getLocalName() + ":" + simpleFeature.getAttribute(a.getLocalName()));
           });
       }
       fileInputStream.close();
       System.out.println("读取JSON完毕！");
   }

5. 遍历要素及属性

SimpleFeatureCollection simpleFeatureCollection = (SimpleFeatureCollection) featuresource.getFeatures();
// 获取当前矢量数据有哪些属性字段值
List<AttributeDescriptor> attributes = simpleFeatureCollection.getSchema().getAttributeDescriptors();
// 拿到迭代器
SimpleFeatureIterator simpleFeatureIterator = simpleFeatureCollection.features();
// 遍历每一个要素
while(simpleFeatureIterator.hasNext()) {
   SimpleFeature simpleFeature = simpleFeatureIterator.next();
   // java8新特性流api
   attributes.stream().forEach((a) -> {
       // 依次读取这个shape中每一个属性值，当然这个属性值，可以处理其它业务
       //Object jzmcValue = feature.getAttribute("jzmc");
       System.out.println(a.getLocalName() + ":" + simpleFeature.getAttribute(a.getLocalName()));
   });
}

6.新增要素（创建图层）

//1.构造TYPE（定义数据的坐标参考、属性字段）
//Class的取值可以是：几何（Point.class,Polygon.class,MultiPolygon.class），属性（String.class，Integer.class， Double.class）
public static SimpleFeatureType createType(Class<?> c, String layerName) {
     SimpleFeatureTypeBuilder builder = new SimpleFeatureTypeBuilder();
     builder.setCRS(DefaultGeographicCRS.WGS84);
     builder.add("FID",String.class);
     // 需要注意的是几何字段的属性名称是固定的
     builder.add("the_geom", c);
     // 设置了图层的名字
     builder.setName(layerName);
     SimpleFeatureType simpleFeatureType = builder.buildFeatureType();
     return simpleFeatureType;
 }
 
 //2.根据TYPE构建单个要素
SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(simpleFeatureType);
WKTReader wktReader = new WKTReader();
Geometry geometry = wktReader.read(wktStr);
// 这里的添加顺序和上面TYPE的时候保持一致
featureBuilder.add("1");
featureBuilder.add(geometry);
SimpleFeature feature = featureBuilder.buildFeature(null);

//3.创建FeatureCollection
List<SimpleFeature> features = new ArrayList<>();
 // 只需要将上面步骤的单个要素放入循环中创建更多的要素
 features.add(feature);
 SimpleFeatureCollection collection = new ListFeatureCollection(TYPE, features);
 
 //4. 生成shpfile
File shpFile = new File(shpPath);
ShapefileDataStoreFactory dataStoreFactory = new ShapefileDataStoreFactory();
// 创造shpstore需要的参数
Map<String, Serializable> params = new HashMap<>();
params.put("url", shpFile.toURI().toURL());
params.put("create spatial index", Boolean.TRUE);
ShapefileDataStore newDataStore = (ShapefileDataStore) dataStoreFactory.createNewDataStore(params);
newDataStore.createSchema(simpleFeatureType);
Transaction transaction = new DefaultTransaction("create");
String typeName = newDataStore.getTypeNames()[0];
SimpleFeatureSource featureSource = newDataStore.getFeatureSource(typeName);
SimpleFeatureStore featureStore = (SimpleFeatureStore) featureSource;
SimpleFeatureCollection collection = new ListFeatureCollection(simpleFeatureType, features);
featureStore.setTransaction(transaction);
featureStore.addFeatures(collection);
featureStore.setTransaction(transaction);
transaction.commit();
transaction.close();

7.读取TIFF

public static Coverage readTiff(String tiffPath) throws IOException {
     File f = new File(tiffPath);
     ParameterValue<OverviewPolicy> policy = AbstractGridFormat.OVERVIEW_POLICY.createValue();
     policy.setValue(OverviewPolicy.IGNORE);
     ParameterValue<String> gridsize = AbstractGridFormat.SUGGESTED_TILE_SIZE.createValue();
     ParameterValue<Boolean> useJaiRead = AbstractGridFormat.USE_JAI_IMAGEREAD.createValue();
     useJaiRead.setValue(true);
     GridCoverage2D image = new GeoTiffReader(f).read(new GeneralParameterValue[]{policy, gridsize, useJaiRead});
     return image;
 }

生成TIFF

File file = new File(outTiffPath);
GeoTiffWriter geoTiffWriter = new GeoTiffWriter(file);
final GeoTiffFormat format = new GeoTiffFormat();
final GeoTiffWriteParams wp = new GeoTiffWriteParams();
// 设置写出参数
wp.setCompressionMode(GeoTiffWriteParams.MODE_DEFAULT);
wp.setTilingMode(GeoToolsWriteParams.MODE_DEFAULT);
ParameterValueGroup paramWrite = format.getWriteParameters();
paramWrite.parameter(AbstractGridFormat.GEOTOOLS_WRITE_PARAMS.getName().toString()).setValue(wp);
geoTiffWriter.write((GridCoverage) coverage, paramWrite.values().toArray(new GeneralParameterValue[4]) );
geoTiffWriter.dispose();

空间分析

空间关系（Geometry Relationships）：
常见的空间关系（Geometry Relationships）包括：Disjoint、Intersects、Touches、Crosses、Within、Contains、Overlaps、Relates。

空间操作（Geometry Operations）：
常见的空间操作（Geometry Operations）包括：Buffer、Intersection、ConvexHull、Intersection、Union、Difference、SymDifference。

裁剪

    /**
     * 裁剪
     *
     **/
    public static  List<String>  newClip(GeoserverAnalysisClipDTO dto) {
        //获取服务数据
        SimpleFeatureCollection collection = GeoServerAnalysisUtils.getFeature(dto.getGeoserverUrl(), dto.getUsername(), dto.getPassword(), dto.getWorkspace(), dto.getLayerName());

        // 要裁剪的几何形状
        Geometry clipGeometry = SpatialAnalysisUtils.convertToJtsGeometry(dto.getPolygon());

        // 创建一个用于存储裁剪后要素的 SimpleFeatureCollection
        //DefaultFeatureCollection clippedCollection = new DefaultFeatureCollection();
        List<Geometry> list = new ArrayList<>();

        // 遍历每个要素
        try (SimpleFeatureIterator iterator = collection.features()) {
            while (iterator.hasNext()) {
                SimpleFeature feature = iterator.next();
                // 获取要素的几何属性
                Geometry featureGeometry = (Geometry) feature.getDefaultGeometry();
                if (featureGeometry.intersects(clipGeometry)) {                    // 判断要素的几何形状与指定几何形状是否有交集
                    Geometry intersection = featureGeometry.intersection(clipGeometry);
                    if (!intersection.isEmpty()) {
                        // 将交集部分保留下来，作为裁剪后的几何形状
                        //clippedCollection.add(SimpleFeatureBuilder.build(feature.getFeatureType(), new Object[]{intersection}, null));
                        list.add(intersection);
                    }
                }
            }
        }

        return geometryToWKT(list);
    }

筛选分析

 // 筛选分析API接口
    public static final String screeningApi = "{{url}}/geoserver/wfs?service=wfs&version=2.0.0&request=GetFeature&outputFormat=application/json&typeNames={{workspace}}:{{layerName}}";

   public static String screening( String geoserverUrl, String username,String password,String workspace,String layerName, Map<String, Object> filter) {
        try {
            String urlStr = buildUrl(screeningApi, geoserverUrl, workspace, layerName);
            String filterStr = buildCQLFilter(filter);
            urlStr = urlStr + "&" + filterStr;
            CredentialsProvider credentialsProvider = new BasicCredentialsProvider();
            credentialsProvider.setCredentials(AuthScope.ANY, new UsernamePasswordCredentials(username, password));
            CloseableHttpClient httpClient = HttpClients.custom().setDefaultCredentialsProvider(credentialsProvider).build();
            HttpGet httpGet = new HttpGet(urlStr);
            HttpResponse response = httpClient.execute(httpGet);
            HttpEntity entity = response.getEntity();
            return EntityUtils.toString(entity);
        } catch (Exception e) {
            e.printStackTrace();
            return "后台异常";
        }
     }

割裂

    // 割裂API接口
    public static final String splitApi = "{{url}}/geoserver/wfs?service=wfs&version=2.0.0&request=GetFeature&outputFormat=application/json&typeNames={{workspace}}:{{layerName}}&splitField={{splitField}}";
    public static String split(String geoserverUrl, String username,String password,String workspace,String layerName,
            String splitField) {
        try {
            String urlStr = buildUrl(splitApi, geoserverUrl, workspace, layerName, splitField);
            CredentialsProvider credentialsProvider = new BasicCredentialsProvider();
            credentialsProvider.setCredentials(AuthScope.ANY, new UsernamePasswordCredentials(username, password));
            CloseableHttpClient httpClient = HttpClients.custom().setDefaultCredentialsProvider(credentialsProvider).build();
            HttpGet httpGet = new HttpGet(urlStr);
            HttpResponse response = httpClient.execute(httpGet);
            HttpEntity entity = response.getEntity();
            return EntityUtils.toString(entity);
        } catch (Exception e) {
            e.printStackTrace();
            return "后台异常";
        }
    }

交集

intersection

并集

union

差集

difference

缓冲区

    public static String buffer(GeoserverAnalysiBufferDTO dto){
        BufferFeatureCollection bf = new BufferFeatureCollection();
        SimpleFeatureCollection featureCollection = getFeature(dto.getGeoserverUrl(),dto.getUsername(),dto.getPassword(),dto.getWorkspace(),dto.getLayerName());
        //featureCollection.getSchema().getDescriptor();
        String[] attributeHeaders = getAttributeHeaders(featureCollection);
        // 打印属性字段头
        for (String header : attributeHeaders) {
            System.out.println(header);
        }
        SimpleFeatureCollection simpleFeatureCollection = bf.execute(featureCollection, dto.getDistance(), dto.getAttribute());
        return Shp2GeojsonUtils.shape2Geojson(simpleFeatureCollection);
    }
    
   private static String[] getAttributeHeaders(SimpleFeatureCollection featureCollection) {
        // 获取第一个要素
        SimpleFeatureIterator iterator = featureCollection.features();
        Feature feature = iterator.next();
        iterator.close();

        // 获取属性字段
        List<String> attributeHeaders = new ArrayList<>();
        for (Property property : feature.getProperties()) {
            String attributeName = property.getName().getLocalPart();
            attributeHeaders.add(attributeName);
        }

        // 将属性字段头转换为字符串数组并返回
        return attributeHeaders.toArray(new String[0]);
    }

擦除

    public static String newErase(GeoserverAnalysisEraseDTO dto) {
        SimpleFeatureCollection collection = GeoServerAnalysisUtils.getFeature(dto.getGeoserverUrl(), dto.getUsername(), dto.getPassword(), dto.getWorkspace(), dto.getLayerName());
        Geometry oldGeo = SpatialAnalysisUtils.convertToJtsGeometry(dto.getPolygon());
        Geometry resultGeo = null;

        try (SimpleFeatureIterator iterator = collection.features()) {
            while (iterator.hasNext()) {
                SimpleFeature feature = iterator.next();
                Geometry geometry = (Geometry) feature.getDefaultGeometry();
                if (geometry instanceof org.locationtech.jts.geom.Polygon || geometry instanceof org.locationtech.jts.geom.MultiPolygon) {
                    if (geometry.intersects(oldGeo)) {
                        Geometry difference = geometry.difference(oldGeo);
                        if (resultGeo == null) {
                            resultGeo = difference;
                        } else {
                            resultGeo = resultGeo.union(difference);
                        }

                    }
                }
            }
        } catch (Exception e) {
            e.printStackTrace();
        }
        WKTWriter writer = new WKTWriter();
        return  writer.write(resultGeo);
    }

点距离

    /**
     * 两点之间距离计算
     * @Date 2024/4/7 15:19
     **/
    public static double getDistance(double x1, double y1, double x2, double y2) {
        // 84坐标系构造GeodeticCalculator
        GeodeticCalculator geodeticCalculator = new GeodeticCalculator(DefaultGeographicCRS.WGS84);
        // 起点经纬度
        geodeticCalculator.setStartingGeographicPoint(x1, y1);
        // 末点经纬度
        geodeticCalculator.setDestinationGeographicPoint(x2, y2);
        // 计算距离，单位：米
        return geodeticCalculator.getOrthodromicDistance();
    }

近邻分析

    public static GeoServerNeighborVO getNeighbor(GeoServerNeighborDTO dto) {

        GeoServerNeighborVO vo = new GeoServerNeighborVO();
        double miniDist = Double.MAX_VALUE;

        SimpleFeatureCollection collection = GeoServerAnalysisUtils.getFeature(dto.getGeoserverUrl(), dto.getUsername(), dto.getPassword(), dto.getWorkspace(), dto.getLayerName());

        if (collection != null) {
            // 遍历 feature 集合
            try (SimpleFeatureIterator iterator = collection.features()) {
                while (iterator.hasNext()) {
                    SimpleFeature feature = iterator.next();
                    Geometry geometry = (Geometry) feature.getDefaultGeometry();

                    String jzmc = null;

                    // 检查几何类型是否为点
                    if (geometry instanceof org.locationtech.jts.geom.Point) {
                        Coordinate coordinate = geometry.getCoordinate();
                        Point point = new Point(coordinate.getX(), coordinate.getY());
                        // 获取点的坐标
                        double x = point.getX();
                        double y = point.getY();
                        //System.out.println("Point coordinates: (" + x + ", " + y + ")");
                        //计算距离
                        double distance = findShort(dto.getPolygon(), point);
                        if (distance != 0) {
                            if (distance < miniDist) {
                                miniDist = distance;

                                Object jzmcValue = feature.getAttribute("jzmc");
                                if (jzmcValue != null) {
                                    jzmc = (String) jzmcValue;
                                }
                                vo.setPoint(point);
                                vo.setName(jzmc);
                            }
                        }
                    } else {
                        System.out.println("This feature is not a point.");
                    }
                }
            } catch (Exception e) {
                e.printStackTrace();
            }
        }
        vo.setDistance(miniDist * 100);
        System.out.println("近邻分析结果:" + vo.toString());
        return vo;
    }

面邻域

    public static List<String> getSquareNeighbor(GeoServerNeighborDTO dto) {

        SimpleFeatureCollection collection = GeoServerAnalysisUtils.getFeature(dto.getGeoserverUrl(), dto.getUsername(), dto.getPassword(), dto.getWorkspace(), dto.getLayerName());
        Geometry oldGeo = convertToJtsGeometry(dto.getPolygon());

        //GeoServerSquareNeighborVO vo = new GeoServerSquareNeighborVO();
        List<Geometry> list = new ArrayList<>();
        List<String> result = new ArrayList<>();
        if (collection == null) {
            return result ;
        }
        try (SimpleFeatureIterator iterator = collection.features()) {
            while (iterator.hasNext()) {
                SimpleFeature feature = iterator.next();
                Geometry geometry = (Geometry) feature.getDefaultGeometry();
                if (geometry instanceof org.locationtech.jts.geom.Polygon || geometry instanceof org.locationtech.jts.geom.MultiPolygon) {
                    // org.locationtech.jts.geom.MultiPolygon multiPolygon = (org.locationtech.jts.geom.MultiPolygon) geometry;
                    if (geometry.intersects(oldGeo)) {
                        list.add(geometry);
                    }
                }
            }
        } catch (Exception e) {
            e.printStackTrace();
        }

        System.out.println("相交:" +  geometryToWKT(list));

        return geometryToWKT(list);
    }