【PostGIS】POSTGIS实现聚类统计提取外轮廓

项目需求根据某些条件进行聚类统计，然后返回聚类的外轮廓，这里主要用到POSTGIS的两个算法，一个是聚类统计功能，一个是提取外轮廓的功能。

1. 聚类统计

Postgis主要实现并提供了四种聚类方法，前两个为窗口函数，后两个为聚合函数：

ST_ClusterKMeans – 该函数是窗口函数，主要是用K-means(K均值聚类)算法进行聚类，算法原理比较简单，容易实现，主要适用于点样本数据，如果是多边形则用多边形的中心点计算，只有一个参数K簇，就是事先明确了要把这一堆样本数据聚成K个类，然后去计算一个目标函数达到最优解，因此K值对结果影响比较大，不太好选取，而且该算法而且适合凸的数据集，这也算法容易收敛；

ST_ClusterDBSCAN – 该函数也是窗口函数，主要利用DBSCAN算法对输入的地理要空间素进行聚类，该算法基于密度进行聚类，直觉上更加符合认知，主要是通过地理要素分布的紧密程度决定，同一类别的样本之间是紧密相连的，不同样本是分离的。该算法相对复杂，比较常用，聚类效果较好，。

ST_ClusterIntersecting – 该函数是一个聚合函数，方法比较好理解，顾名思义就是把相交关联关系的地理要素作为一个类簇。具体的函数声明和用法，看官网的例子介绍很容易看懂，也可以拿数据进行测试实验，不是很常用，不做重点；

ST_ClusterWithin – 该函数也是一个聚合函数，也很好理解，顾名思义通过一个距离参数，将距离之内的地理要素聚成一个类别，距离之外的就是非同类。详细可参考官网，也不做重点；

综合自己的需求，采用ST_ClusterDBSCAN函数
函数定义，有三个参数，（窗口函数中的每一个地理要素，搜索半径，最小点数），根据以上分析应该不难理解：

integer ST_ClusterDBSCAN(geometry winset geom, float8 eps, integer minpoints);
例：

SELECT ST_ClusterDBSCAN(geom, eps := 0.000179, minpoints := 2) OVER () AS cluster_id,* FROM 表名

这里的eps需要根据表的srid来调整，我这边表的srid2为4490，是经纬度坐标，所以这里的0.000179也是经纬度的距离，大约20米，返回的cluster_id就是聚类后的类id，可以通过它将同一聚类的记录取出来。

2.外轮廓计算

计算凸包
计算凸包十分简单，使用ST_ConvexHull即可，但同时可能失去很多特征点，生成的任务区比较粗，只能是一个大致的范围。

select ST_ConvexHull ((select ST_Collect(shape) from 表名))

计算凹包
计算凹包十分简单相对复杂一些，使用ST_ConcaveHull，但相对于凸包会多保留一些特征点，生成的任务区保留了一些弯曲变化。

这个函数有三个参数：

geometry geom ，几何体
float target_percent，特征点抛弃的比例，0.1-0.99，越小特征保留的越多，但计算效率越低。当等于1时，和凸包的计算效果是一样的。
boolean allow_holes ，是否允许有洞

select ST_ConcaveHull ((select ST_Collect(shape) from 表名),0.7)

根据需求，我们选择凹包的计算方式

3.实现需求

这里聚类和计算外轮廓分步进行计算

首先进行聚类查询

SELECT ST_ClusterDBSCAN(geom, eps := 0.000179, minpoints := 2) OVER () AS cluster_id,* FROM allrecord

这里allrecord我是用了with子句的查询结果，目的是实现按照项目需求条件进行查询，这里不再展示

在java程序中完成聚类统计

Map<Integer, List<实体类>> collect = 聚类查询结果list.stream().filter(tcql -> tcql.getCluster_id() != null).collect(Collectors.groupingBy(实体类::getCluster_id));

然后分别将上述map中的id分别取出来，后续计算外轮廓是用

计算外轮廓
这里直接将mapper的配置写出来

 <select id="selectConcaveHull" resultType="java.util.Map">
        with concaveResult as(
        select ST_ConcaveHull ((select ST_Collect(geom) from 表名 where id in
        <foreach collection="idList" item="id" separator="," open="(" close=")">
            #{id}
        </foreach>
        ),0.7) geom
        )
        select st_astext(geom) geomtext,st_area(st_transform(geom,4528)) geomarea from concaveResult

    </select>