基於K-means 切割多邊形 JAVA實現

時間 2019-12-06

原文原文鏈接

基於K-means 切割多邊形 JAVA實現

依賴

geotools
ekmeans

<!--geotools-->
    <!-- https://mvnrepository.com/artifact/org.locationtech.jts/jts-core -->
    <dependency>
      <groupId>org.locationtech.jts</groupId>
      <artifactId>jts-core</artifactId>
      <version>1.16.0</version>
    </dependency>
    <dependency>
      <groupId>junit</groupId>
      <artifactId>junit</artifactId>
      <version>4.12</version>
      <scope>test</scope>
    </dependency>

    <dependency>
      <groupId>org.geotools</groupId>
      <artifactId>gt-shapefile</artifactId>
      <version>${geotools.version}</version>
    </dependency>

    <!-- https://mvnrepository.com/artifact/org.geotools/gt-main -->
    <dependency>
      <groupId>org.geotools</groupId>
      <artifactId>gt-main</artifactId>
      <version>${geotools.version}</version>
    </dependency>

    <dependency>
      <groupId>org.geotools</groupId>
      <artifactId>gt-swing</artifactId>
      <version>${geotools.version}</version>
    </dependency>
    <dependency>
      <groupId>org.geotools</groupId>
      <artifactId>gt-geojson</artifactId>
      <version>${geotools.version}</version>
    </dependency>
    <dependency>
      <groupId>org.geotools</groupId>
      <artifactId>gt-geometry</artifactId>
      <version>${geotools.version}</version>
    </dependency>
    <!-- https://mvnrepository.com/artifact/com.vividsolutions/jts -->
    <dependency>
      <groupId>com.vividsolutions</groupId>
      <artifactId>jts</artifactId>
      <version>1.13</version>
    </dependency>


    <dependency>
      <groupId>org.geotools</groupId>
      <artifactId>gt-geojson</artifactId>
      <version>${geotools.version}</version>
    </dependency>
<!--kmeans 實現-->
   <dependency>
      <groupId>ca.pjer</groupId>
      <artifactId>ekmeans</artifactId>
      <version>2.0.0</version>
      <scope>compile</scope>
    </dependency>

實現過程

結果類 KmeanPolygonResult

package com.huifer.planar.asEt.entity;

import java.util.ArrayList;
import java.util.List;
import lombok.AllArgsConstructor;
import lombok.Data;
import lombok.NoArgsConstructor;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.Point;
import org.locationtech.jts.geom.Polygon;

/**
 * <p>Title : KmeanPolygonResult </p>
 * <p>Description : Kmean+泰森多邊形平分多邊形結果</p>
 *
 * @author huifer
 * @date 2019-01-16
 */
@Data
@NoArgsConstructor
@AllArgsConstructor
public class KmeanPolygonResult {

    /**
     * 原始平面
     */
    private Polygon polygon;
    /**
     * 隨機點集合
     */
    private ArrayList<Point> pointList;
    /**
     * 分組後組id
     */
    private int[] assignments;
    /**
     * 中心集合
     */
    private double[][] centroids;
    /**
     * xlist
     */
    private ArrayList<Double> xlist;
    /**
     * ylist
     */
    private ArrayList<Double> ylist;
    /**
     * 構造泰森多邊形
     */
    private List<Geometry> voronoi;
}

k-means 聚合類

package com.huifer.planar.asEt.utils;

import ca.pjer.ekmeans.EKmeans;
import lombok.Data;

/**
 * <p>Title : Kmeans </p>
 * <p>Description : Kmeans 聚合</p>
 *
 * @author huifer
 * @date 2019-01-15
 */
@Data
public class Kmeans {

    /***
     * 聚合總量
     */
    private int n;
    /***
     * 數據集合
     */
    private double[][] points;
    /**
     * 簇族數量
     */
    private int k;

    /**
     * 中心集合
     */
    private double[][] centroids;
    /**
     * 分組後組id
     */
    private int[] assignments;

    public Kmeans(double[][] points, int k) {
        this.n = points.length;
        this.points = points;
        this.k = k;
        this.centroids = new double[k][2];
        init();
    }

    /**
     * 求得分組數據
     */
    public void init() {
        EKmeans eKmeans = new EKmeans(centroids, points);
        eKmeans.setIteration(128);
        eKmeans.setEqual(true);
        eKmeans.setDistanceFunction(EKmeans.EUCLIDEAN_DISTANCE_FUNCTION);
        eKmeans.run();
        int[] assignments = eKmeans.getAssignments();
        this.assignments = assignments;
        double[][] centroids = eKmeans.getCentroids();
        this.centroids = centroids;
    }

}

泰森多邊形構造

package com.huifer.planar.asEt.algo.impl;

import com.huifer.planar.asEt.algo.VoronoiInterface;
import java.util.ArrayList;
import java.util.Collection;
import java.util.List;
import org.geotools.geometry.jts.JTSFactoryFinder;
import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Envelope;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.Point;
import org.locationtech.jts.triangulate.VoronoiDiagramBuilder;
import org.locationtech.jts.triangulate.quadedge.QuadEdgeSubdivision;

/**
 * <p>Title : VoronoiInterfaceImpl </p>
 * <p>Description : 泰森多邊形 & 德勞內三角形 </p>
 *
 * @author huifer
 * @date 2019-01-16
 */
public class VoronoiInterfaceImpl implements VoronoiInterface {

    @Override
    public List<Geometry> voronoi(double[][] doubles) {
        return voronoi(doublesToCoordinate(doubles));
    }

    @Override
    public Collection delaunay(double[][] doubles) {
        return delaunay(doublesToCoordinate(doubles));
    }

    @Override
    public List<Geometry> voronoi(ArrayList<Point> points) {
        List<Geometry> voronoi = voronoi(pointToCoordinate(points));
        return voronoi;
    }

    @Override
    public Collection delaunay(ArrayList<Point> points) {
        Collection delaunay = delaunay(pointToCoordinate(points));
        return delaunay;
    }


    private List<Coordinate> doublesToCoordinate(double[][] doubles) {
        List<Coordinate> coords = new ArrayList<Coordinate>();
        for (int i = 0; i < doubles.length; i++) {
            Coordinate coord = new Coordinate(doubles[i][0], doubles[i][1], i);
            coords.add(coord);
        }
        return coords;
    }


    private List<Coordinate> pointToCoordinate(ArrayList<Point> points) {
        List<Coordinate> coords = new ArrayList<Coordinate>();
        for (int i = 0; i < points.size(); i++) {
            Coordinate coord = new Coordinate(points.get(i).getX(), points.get(i).getY(), i);
            coords.add(coord);
        }
        return coords;
    }


    private VoronoiDiagramBuilder getVoronoiDiagramBuilder(List<Coordinate> coords) {
        VoronoiDiagramBuilder voronoiDiagramBuilder = new VoronoiDiagramBuilder();
        Envelope clipEnvelpoe = new Envelope();
        voronoiDiagramBuilder.setSites(coords);
        voronoiDiagramBuilder.setClipEnvelope(clipEnvelpoe);
        return voronoiDiagramBuilder;
    }


    @Override
    public List<Geometry> voronoi(List<Coordinate> coords) {
        VoronoiDiagramBuilder voronoiDiagramBuilder = getVoronoiDiagramBuilder(
                coords);
        Geometry geom = voronoiDiagramBuilder
                .getDiagram(JTSFactoryFinder.getGeometryFactory());

        List<Geometry> voronoi = new ArrayList<>();
        int numGeometries = geom.getNumGeometries();
        for (int i = 0; i < numGeometries; i++) {
            Geometry geometryN = geom.getGeometryN(i);
            voronoi.add(geometryN);
        }
        return voronoi;
    }

    @Override
    public Collection delaunay(List<Coordinate> coords) {
        VoronoiDiagramBuilder voronoiDiagramBuilder = getVoronoiDiagramBuilder(coords);
        QuadEdgeSubdivision subdivision = voronoiDiagramBuilder.getSubdivision();
        Collection delaunay = subdivision.getEdges();
        return delaunay;
    }
}

Kmeans 切割面類

package com.huifer.planar.asEt.algo.impl;

import com.huifer.planar.asEt.algo.KmeanPolygonSplitInterface;
import com.huifer.planar.asEt.algo.VoronoiInterface;
import com.huifer.planar.asEt.entity.KmeanPolygonResult;
import com.huifer.planar.asEt.utils.Kmeans;
import com.huifer.planar.asEt.utils.shptools.overlay.Operation;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;
import org.locationtech.jts.geom.Coordinate;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.Point;
import org.locationtech.jts.geom.Polygon;
import org.locationtech.jts.io.ParseException;

/**
 * <p>Title : KmeanPolygonSplitCore </p>
 * <p>Description : Kmean split polygon core </p>
 *
 * @author huifer
 * @date 2019-01-16
 */
public class KmeanPolygonSplitCore implements KmeanPolygonSplitInterface {

    private static double random(double max, double min) {
        double d = (Math.random() * (max - min) + min);
        return d;
    }

    @Override
    public KmeanPolygonResult splitPolygon(String wkt, int setp, int k) throws ParseException {
        Operation op = new Operation();
        KmeanPolygonResult result = new KmeanPolygonResult();
        Polygon polygon = op.createPolygonByWKT(wkt);
        Coordinate[] coordinates = polygon.getCoordinates();
        ArrayList<Double> xList = new ArrayList<>();
        ArrayList<Double> yList = new ArrayList<>();

        Arrays.stream(coordinates).forEach(
                s -> {
                    xList.add(s.x);
                    yList.add(s.y);
                }
        );
        // xy 最大最小值
        Double xMax = xList.stream().reduce(Double::max).get();
        Double xMin = xList.stream().reduce(Double::min).get();
        Double yMax = yList.stream().reduce(Double::max).get();
        Double yMin = yList.stream().reduce(Double::min).get();

        // 當前點數量
        int pointCount = 0;
        ArrayList<Point> pointArrayList = new ArrayList<>();

        for (int i = 0; i < Integer.MAX_VALUE; i++) {
            // 最大最小值隨機
            if (pointCount <= setp) {
                double rx = random(xMax, xMin);
                double ry = random(yMax, yMin);
                Point nowPoint = op.createPointByWkt("POINT(" + rx + " " + ry + ")");
                boolean contains = polygon.contains(nowPoint);
                if (contains) {
                    pointArrayList.add(nowPoint);
                    pointCount++;
                }
            } else {
                break;
            }
        }
        // k-means 數據 構造
        double[][] kmData = new double[pointArrayList.size()][2];
        for (int i = 0; i < pointArrayList.size(); i++) {
            Point point = pointArrayList.get(i);
            double[] oneData = new double[2];
            oneData[0] = point.getX();
            oneData[1] = point.getY();
            kmData[i] = oneData;
        }

        // k-means 結果

        Kmeans kmeans = new Kmeans(kmData, k);

        // 構造泰森多邊形
        VoronoiInterface vo = new VoronoiInterfaceImpl();
        List<Geometry> voronoi = vo.voronoi(kmeans.getCentroids());

        result.setPolygon(polygon);
        result.setPointList(pointArrayList);
        result.setAssignments(kmeans.getAssignments());
        result.setCentroids(kmeans.getCentroids());
        result.setXlist(xList);
        result.setYlist(yList);
        result.setVoronoi(voronoi);
        return result;
    }

}