降雨等值面生成不正确

Question

问题遇到的现象和发生背景

背景：降雨等值面生成
现象：没有任何等值线生成

问题相关代码，请勿粘贴截图

public class EquiSurface {

private static final double EPSILON = 0.001;
/**
 * 生成等值面
 *
 * @param trainData    训练数据
 * @param dataInterval 数据间隔
 * @param size         大小，宽，高
 * @param boundryFile  四至
 * @param isclip       是否裁剪
 * @return
 */
public String calEquiSurface(double[][] trainData,
                             double[] dataInterval,
                             int[] size,
                             String boundryFile,
                             boolean isclip) {
    String geojsonpogylon = "";
    try {
        double _undefData = -9999.0;
        SimpleFeatureCollection polygonCollection = null;
        List<PolyLine> cPolylineList = new ArrayList<PolyLine>();
        List<Polygon> cPolygonList = new ArrayList<Polygon>();

        int width = size[0],
                height = size[1];
        double[] _X = new double[width];
        double[] _Y = new double[height];

        File file = new File(boundryFile);
        ShapefileDataStore shpDataStore = null;

        shpDataStore = new ShapefileDataStore(file.toURL());
        //设置编码
        Charset charset = Charset.forName("GBK");
        shpDataStore.setCharset(charset);
        String typeName = shpDataStore.getTypeNames()[0];
        SimpleFeatureSource featureSource = shpDataStore.getFeatureSource(typeName);
        SimpleFeatureCollection fc = featureSource.getFeatures();

        double minX = fc.getBounds().getMinX();
        double minY = fc.getBounds().getMinY();
        double maxX = fc.getBounds().getMaxX();
        double maxY = fc.getBounds().getMaxY();

        Interpolate.CreateGridXY_Num(minX, minY, maxX, maxY, _X, _Y);

        double[][] _gridData = new double[width][height];

        int nc = dataInterval.length;

        _gridData = Interpolate.Interpolation_IDW_Neighbor(trainData,
                _X, _Y, dataInterval.length, _undefData);// IDW插值
        int[][] S1 = new int[_gridData.length][_gridData[0].length];
        List<Border> _borders = Contour.tracingBorders(_gridData, _X, _Y,
                S1, _undefData);

        cPolylineList = Contour.tracingContourLines(_gridData, _X, _Y, nc,
                dataInterval, _undefData, _borders, S1);// 生成等值线

        cPolylineList = Contour.smoothLines(cPolylineList);// 平滑

        cPolygonList = Contour.tracingPolygons(_gridData, cPolylineList,
                _borders, dataInterval);

        geojsonpogylon = getPolygonGeoJson(cPolygonList);
        // 根据边界裁剪等值面

        if (isclip) {
            //polygonCollection = readGeoJsonByString(geojsonpogylon);
            //获取等值面要素 SimpleFeatureCollection
            polygonCollection = getSimpleCollectionByGeoJSON(geojsonpogylon);
            SimpleFeatureCollection sm = clipPolygonFeatureCollection(fc, polygonCollection);
            FeatureCollection featureCollection = sm;
            FeatureJSON featureJSON = new FeatureJSON();
            StringWriter writer = new StringWriter();
            featureJSON.writeFeatureCollection(featureCollection, writer);
            geojsonpogylon = writer.toString();
        }
    } catch (Exception e) {
        e.printStackTrace();
    }
    return geojsonpogylon;
}

/**
 * 通过geojson获取要素类
 * @param geojsom
 * @return
 * @throws IOException
 */
private static SimpleFeatureCollection getSimpleCollectionByGeoJSON(String geojsom) throws IOException {
    FeatureJSON featureJSON=new FeatureJSON();
    FeatureCollection featureCollection=featureJSON.readFeatureCollection(geojsom);
    List<SimpleFeature> list=new ArrayList<>();
    FeatureIterator featureIterator=featureCollection.features();
    while (featureIterator.hasNext()){
        SimpleFeature feature=(SimpleFeature) featureIterator.next();
        list.add(feature);
    }
    SimpleFeatureCollection simpleFeatureCollection=new ListFeatureCollection((SimpleFeatureType) featureCollection.getSchema(),list);
    return simpleFeatureCollection;
}


/**
 * 要素类转geojson
 * @param simpleFeatureCollection
 * @return
 * @throws IOException
 */
/*private static String SimpleCollectionToGeoJSON(SimpleFeatureCollection simpleFeatureCollection) throws IOException {
    FeatureCollection featureCollection=(FeatureCollection) simpleFeatureCollection;
    FeatureJSON featureJSON=new FeatureJSON();
    StringWriter writer = new StringWriter();
    featureJSON.writeFeatureCollection(featureCollection,writer);
    return writer.toString();

}*/

/**
 * 裁剪图形
 * @param fc
 * @return
 */
private SimpleFeatureCollection clipPolygonFeatureCollection(FeatureCollection fc,
                                                             SimpleFeatureCollection gs) throws SchemaException, IOException {
    SimpleFeatureCollection simpleFeatureCollection = null;
    SimpleFeatureType TYPE = DataUtilities.createType("polygons",
            "the_geom:MultiPolygon,lvalue:double,hvalue:double");
    List<SimpleFeature> list = new ArrayList<>();
    SimpleFeatureIterator contourFeatureIterator = gs.features();
    FeatureIterator dataFeatureIterator = fc.features();
    while (dataFeatureIterator.hasNext()){
        SimpleFeature dataFeature = (SimpleFeature) dataFeatureIterator.next();
        Geometry dataGeometry = (Geometry) dataFeature.getDefaultGeometry();
        contourFeatureIterator = gs.features();
        while (contourFeatureIterator.hasNext()){
            SimpleFeature contourFeature = contourFeatureIterator.next();
            Geometry contourGeometry = (Geometry) contourFeature.getDefaultGeometry();
            Double lv = (Double) contourFeature.getProperty("lvalue").getValue();
            Double hv = (Double) contourFeature.getProperty("hvalue").getValue();
            if (dataGeometry.getGeometryType() == "MultiPolygon"){
                for (int i = 0; i <dataGeometry.getNumGeometries(); i++){
                    Geometry geom = dataGeometry.getGeometryN(i);
                    if (geom.intersects(contourGeometry)){
                        Geometry geo = geom.intersection(contourGeometry);
                        SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE);
                        featureBuilder.add(geo);
                        featureBuilder.add(lv);
                        featureBuilder.add(hv);
                        SimpleFeature feature = featureBuilder.buildFeature(null);
                        list.add(feature);

                    }
                }

            } else {
                if (dataGeometry.intersects(contourGeometry)){
                    Geometry geo = dataGeometry.intersection(contourGeometry);
                    SimpleFeatureBuilder featureBuilder = new SimpleFeatureBuilder(TYPE);
                    featureBuilder.add(geo);
                    featureBuilder.add(lv);
                    featureBuilder.add(hv);
                    SimpleFeature feature = featureBuilder.buildFeature(null);
                    list.add(feature);

                }

            }

        }
    }

    contourFeatureIterator.close();
    dataFeatureIterator.close();
    simpleFeatureCollection = new ListFeatureCollection(TYPE, list);

    return simpleFeatureCollection;

}


/**
 * Contour面转geojson2
 * @param polygonList
 * @return
 */
/*private static String polygonToGeoJson(List<Polygon> polygonList) {
    String geo;
    String geometry = "{ \"type\":\"Feature\",\"geometry\":";
    String properties = ",\"properties\":{ \"hvalue\":";

    String head = "{\"type\": \"FeatureCollection\"," + "\"features\":[";
    String end = "]}";

    if (polygonList == null || polygonList.isEmpty()) {
        return null;
    }

    StringBuilder geoBuilder = new StringBuilder();

    for (Polygon pPolygon : polygonList) {
        List<Object> ptsTotal = new ArrayList<>();
        List<Object> pts = new ArrayList<>();
        pPolygon.OutLine.PointList.stream().map(point -> {
            List<Double> pt = new ArrayList<>();
            pt.add(point.X);
            pt.add(point.Y);
            return pt;
        }).forEachOrdered(pts::add);
        ptsTotal.add(pts);

        if (pPolygon.HasHoles()) {
            // 有孔洞的多面
            pPolygon.HoleLines.stream().map(cptLine -> {
                List<Object> cpts = new ArrayList<>();

                cptLine.PointList.stream().map(ccptD -> {
                    List<Double> pt = new ArrayList<>();
                    pt.add(ccptD.X);
                    pt.add(ccptD.Y);
                    return pt;
                }).forEachOrdered(cpts::add);
                return cpts;
            }).filter(cpts -> !cpts.isEmpty()).forEachOrdered(ptsTotal::add);
        }

        JSONObject js = new JSONObject();

        js.put("type", "Polygon");
        js.put("coordinates", ptsTotal);

        double hv = pPolygon.HighValue;
        double lv = pPolygon.LowValue;

        if (Math.abs(hv - lv) < EPSILON) {
            if(pPolygon.IsHighCenter){
                hv=hv+0.1;
                lv=lv+0.1;
            }else{
                hv=hv-0.1;
                lv=lv-0.1;
            }
        }else{
            hv=hv-0.1;
            lv=lv+0.1;
        }

        geoBuilder.insert(0, geometry + js.toString() + properties + hv + ", \"lvalue\":" + lv + "}}" + ",");
    }
    geo = geoBuilder.toString();

    if (geo.contains(",")) {
        geo = geo.substring(0, geo.lastIndexOf(","));
    }

    geo = head + geo + end;

    return geo;
}*/

/**
 * 转geojson1
 * @param cPolygonList
 * @return
 */
private String getPolygonGeoJson(List<Polygon> cPolygonList) {
    String geo = null;
    String geometry = " { \"type\":\"Feature\",\"geometry\":";
    String properties = ",\"properties\":{ \"hvalue\":";

    String head = "{\"type\": \"FeatureCollection\"," + "\"features\": [";
    String end = "  ] }";
    if (cPolygonList == null || cPolygonList.size() == 0) {
        return null;
    }
    try {
        for (Polygon pPolygon : cPolygonList) {

            List<Object> ptsTotal = new ArrayList<Object>();
            List<Object> pts = new ArrayList<Object>();

            PolyLine pline = pPolygon.OutLine;

            for (PointD ptD : pline.PointList) {
                List<Double> pt = new ArrayList<Double>();
                pt.add(ptD.X);
                pt.add(ptD.Y);
                pts.add(pt);
            }

            ptsTotal.add(pts);

            if (pPolygon.HasHoles()) {
                for (PolyLine cptLine : pPolygon.HoleLines) {
                    List<Object> cpts = new ArrayList<Object>();
                    for (PointD ccptD : cptLine.PointList) {
                        List<Double> pt = new ArrayList<Double>();
                        pt.add(ccptD.X);
                        pt.add(ccptD.Y);
                        cpts.add(pt);
                    }
                    if (cpts.size() > 0) {
                        ptsTotal.add(cpts);
                    }
                }
            }

            JSONObject js = new JSONObject();
            js.put("type", "Polygon");
            js.put("coordinates", ptsTotal);
            double hv = pPolygon.HighValue;
            double lv = pPolygon.LowValue;

            if (hv == lv) {
                if (pPolygon.IsClockWise) {
                    if (!pPolygon.IsHighCenter) {
                        hv = hv - 0.1;
                        lv = lv - 0.1;
                    }

                } else {
                    if (!pPolygon.IsHighCenter) {
                        hv = hv - 0.1;
                        lv = lv - 0.1;
                    }
                }
            } else {
                if (!pPolygon.IsClockWise) {
                    lv = lv + 0.1;
                } else {
                    if (pPolygon.IsHighCenter) {
                        hv = hv - 0.1;
                    }
                }

            }

            geo = geometry + js.toString() + properties + hv
                    + ", \"lvalue\":" + lv + "} }" + "," + geo;

        }
        if (geo.contains(",")) {
            geo = geo.substring(0, geo.lastIndexOf(","));
        }

        geo = head + geo + end;
    } catch (Exception e) {
        e.printStackTrace();
        return geo;
    }
    return geo;
}




public static SimpleFeatureCollection readGeoJsonByString(String geojsonpogylon) throws IOException {
    FeatureJSON fjson = new FeatureJSON();
    SimpleFeatureCollection featureCollection = (SimpleFeatureCollection) fjson.readFeatureCollection(geojsonpogylon);
    return featureCollection;
}

public void append2File(String strFile, String strJson){
    File f=new File(strFile);//新建一个文件对象，如果不存在则创建一个该文件
    FileWriter fw;
    try {
        fw = new FileWriter(f);
        fw.write(strJson);//将字符串写入到指定的路径下的文件中
        fw.close();
    } catch (IOException e) { e.printStackTrace(); }
}

public static void main(String[] args) {
    long start = System.currentTimeMillis();

    EquiSurface equiSurface = new EquiSurface();

    // 雨量站点经纬度及雨量值
    JSONArray jsonArray = JSONArray.parseArray(arr);
    double[][] trainData = new double[jsonArray.size()][3];
    for (int i = 0; i < jsonArray.size(); i++){
        trainData[i][0] = jsonArray.getJSONObject(i).getDouble("lng");
        trainData[i][1] = jsonArray.getJSONObject(i).getDouble("lat");
        trainData[i][2] = jsonArray.getJSONObject(i).getDouble("rainFall");
    }

    double[] dataInterval = new double[]{0.0001, 2, 4, 6, 9, 10, 20, 30};
    //double[] dataInterval = new double[]{0,5,10,15,20,25,30,35,40,45,50,55,60,65,70,75,80,85,90,95,100,105,110,115,120,125,130,135,140,145,150,155,160};

    String boundryFile = "E:\\开发war包\\桐庐包\\地图文件\\gdmaptonglu.shp";
    //String boundryFile = "E:\\开发war包\\桐庐包\\地图文件\\共享的shp\\dgBoundary.shp";

    int[] size = new int[]{100, 100};

    boolean isclip = true;

    try {
        String strJson = equiSurface.calEquiSurface(trainData, dataInterval, size, boundryFile, isclip);
        String strFile = "E:\\开发war包\\桐庐包\\testdata_json\\jinzhou.json";
        equiSurface.append2File(strFile, strJson);
        System.out.println(strFile + "差值成功, 共耗时" + (System.currentTimeMillis() - start) + "ms");
    } catch (Exception e) {
        e.printStackTrace();
    }
}


public static String arr = "[\n" +
        "        {\n" +
        "            \"siteId\": \"1523\",\n" +
        "            \"siteName\": \"富春江坝上\",\n" +
        "            \"rainFall\": 4.0,\n" +
        "            \"lng\": 119.66022265320088,\n" +
        "            \"lat\": 29.708598670961084,\n" +
        "            \"time\": \"2022-06-21 00:00:00\"\n" +
        "        },\n" +
        "        {\n" +
        "            \"siteId\": \"1525\",\n" +
        "            \"siteName\": \"桐君山\",\n" +
        "            \"rainFall\": 4.0,\n" +
        "            \"lng\": 119.686385,\n" +
        "            \"lat\": 29.809936,\n" +
        "            \"time\": \"2022-06-21 00:00:00\"\n" +
        "        }
        "    ]";

}

运行结果及报错内容

运行无报错，将生成好的json放到了geojson.io平台上渲染出来没有任何的等值线生成

我的解答思路和尝试过的方法

我想要达到的结果

后台生成geojson数据，根据真实雨量数据画出等值线，前端进行颜色填充