使用 geotools 多边形化光栅文件（geotiff）时遇到问题

Question

我有一些大的geotiff文件，大小范围从200MB到10GB。 它们是由无人机按照特定路线生成的，以一定的宽度捕获图像。 这是一个带有马赛克的示例： 示例tiff 文件大小约为200MB。

现在我想获取这个tif文件中的多边形信息，就像图中红线标记的区域一样。 所以我尝试使用 PolygonExtractionProcess 来完成这项工作，但它需要永远，而且我的 cpu 几乎运行到 100%，我不得不放弃。

然后我发现一些“中间”tif 文件，只有白色和黑色，它们的大小非常小（从 50KB 到 10MB），如下所示： 黑白色tiff

在我使用这个小tif作为输入后，多边形化过程持续大约1-2分钟，我得到这样的结果（不完全是我想要的，但非常接近，我只想要黑色部分）： 多边形化结果

代码是这样的（只是非常简单的用法）：

  AbstractGridCoverage2DReader gridCoverage2DReader = gridFormat.getReader(iis, hints);
                GridCoverage2D coverage = gridCoverage2DReader.read(null);
                PolygonExtractionProcess process = new PolygonExtractionProcess();
                SimpleFeatureCollection features = process.execute(coverage, 0, Boolean.TRUE, null, null, null, null);
                SimpleFeatureIterator sfi = features.features();

经过一番搜索，我意识到在多边形化之前，我需要对tiff文件进行灰度或黑白二值化以减小文件大小。但是使用JAVA/GeoTools似乎很难做到这一点，很多解决方案都使用GDAL或GIS软件像 QGIS/ARCGIS。

我的问题是：

1. 我走在正确的道路上吗？如果没有的话有没有标准的解决方案？
2. Java/Geotools 可以高效地完成大型 tif 文件的灰度工作而不占用太多内存吗？
3. 我怎样才能从黑白tiff中只得到黑色部分的多边形？

Answer 1

地理工具

在这个方法中，我使用了JDK 8的并行流处理。虽然按边界划分可以节省内存，但速度仍然很慢。如果有更好的解决方案，我们可以讨论交流。

/**
 * Converts the grid coverage to a SimpleFeatureCollection.
 * This method divides the given grid coverage into multiple small grids, extracts the geometry for each small grid,
 * and then merges these geometries into a SimpleFeatureCollection.
 *
 * @param gridCoverage2D The grid coverage that provides the geographic range and data to be processed.
 * @param nodata         The value used to fill no-data regions.
 * @param handle         A functional interface used for further processing of the merged geometries.
 * @param row            The number of rows to divide.
 * @param cols           The number of columns to divide.
 * @return Returns a SimpleFeatureCollection containing the processed geometries.
 */
public static SimpleFeatureCollection subGrid2SFC(GridCoverage2D gridCoverage2D, double nodata, boolean insideEdges,
                                                  int row, int cols, Function<Geometry, Geometry> handle) {
    // Divide the grid coverage into smaller regions
    List<Envelope2D> envelope2DS = GeoToolsUtils.divideEnvelope(gridCoverage2D.getEnvelope2D(), row, cols);

    PolygonExtractionProcess process = new PolygonExtractionProcess();
    List<Geometry> rets = new CopyOnWriteArrayList<>();
    AtomicInteger atomicInteger = new AtomicInteger();
    // Parallel processing of each small grid
    envelope2DS.parallelStream()
            .forEach(d -> {
                // Log the processing progress
                if (LOG.isDebugEnabled()) {
                    int i = atomicInteger.incrementAndGet();
                    LOG.debug("Current progress of raster to vector conversion: {} out of {}", i, envelope2DS.size());
                }
                // Convert the boundaries of the small grid into polygons
                Polygon polygon = GeoToolsUtils.createPolygon(d);
                // Extract data from the grid coverage using the polygon and create a SimpleFeatureCollection
                SimpleFeatureCollection ret = process.execute(gridCoverage2D, 0, insideEdges, polygon,
                        Collections.singleton(nodata), null, null);
                if (ret.isEmpty()) {
                    return;
                }
                // Convert the SimpleFeatureCollection into a list of geometries
                List<Geometry> geometries = GeoToolsUtils.fcToGeometryList(ret);
                // Merge the list of geometries into a single geometry
                Geometry union = GeoToolsUtils.union(geometries);
                if (Objects.isNull(union) || union.isEmpty()) {
                    return;
                }
                // Apply additional geometry processing logic
                Geometry apply = handle.apply(union);
                rets.add(apply);
            });
    CoordinateReferenceSystem crs = gridCoverage2D.getCoordinateReferenceSystem();
    return GeoToolsUtils.geometries2FeatureCollection(rets, "subGrid2SFC", crs);
}

GDAL

使用 GDAL 速度更快，内存效率更高。

#!/bin/bash

# Check the number of parameters
if [ "$#" -ne 3 ]; then
    echo "Usage: $0 <input_tiff> <band> <output_shapefile>"
    exit 1
fi

# Parameter 1: Path to the input TIFF file
INPUT_TIFF=$1

BAND=$2

# Parameter 2: Path to the output shapefile
OUTPUT_SHP=$3

# Execute the gdal_polygonize.py command
gdal_polygonize.py "$INPUT_TIFF" -b "$BAND" -f 'ESRI Shapefile' "$OUTPUT_SHP" 1 fid

    /**
     * 栅格转矢量
     *
     * @param gridPath        网格路径
     * @param targetShapePath 目标形状路径
     * @param band            波段
     */
    public static void grid2Vector(Path gridPath, Path targetShapePath, int band) {
        String gdalCMD = your_shell_path + "/gdal_polygonize.sh %s %s %s";
        String cmd = String.format(gdalCMD, gridPath.toAbsolutePath(), band, targetShapePath.toAbsolutePath());
        LOG.debug("gdal_polygonize cmd：{}", cmd);
        String rests = RuntimeUtil.execForStr("/bin/bash " + cmd);
        LOG.debug("gdal 栅格转矢量输出 {}", rests);
        if (Files.notExists(targetShapePath)) {
            throw new ServiceException("gdal 栅格转矢量失败！");
        }
    }