如何在terra中实现类QGIS的多边形范围随机采样点功能？

核心优化思路：跳过全量栅格裁剪，先在多边形矢量范围内生成随机点，再提取点对应位置的栅格像素值。这种方式仅处理采样点对应的极小部分栅格数据，能大幅降低计算量和内存占用。

1. Python + GeoPandas/GDAL 实现 #

适合自定义需求较多的场景，完全开源且内存效率极高：

• 生成多边形内的随机点：

  import geopandas as gpd
  import numpy as np
  from osgeo import gdal
  
  # 读取多边形图层
  poly_gdf = gpd.read_file("polygons.shp")
  # 合并所有不连续多边形为单个几何对象，简化随机点判断
  merged_poly = poly_gdf.unary_union
  # 设置采样点数量
  sample_count = 1000
  
  # 生成指定数量的内部随机点
  random_points = []
  min_x, min_y, max_x, max_y = merged_poly.bounds
  while len(random_points) < sample_count:
      x = np.random.uniform(min_x, max_x)
      y = np.random.uniform(min_y, max_y)
      point = gpd.points_from_xy([x], [y])[0]
      if merged_poly.contains(point):
          random_points.append(point)
  
  # 转为GeoDataFrame统一处理
  points_gdf = gpd.GeoDataFrame(geometry=random_points, crs=poly_gdf.crs)
  

• 提取随机点对应的栅格值：

  # 打开栅格文件（仅读取元数据，不加载全量像素）
  raster_ds = gdal.Open("large_raster.tif")
  raster_band = raster_ds.GetRasterBand(1)
  geo_transform = raster_ds.GetGeoTransform()
  
  # 坐标转栅格像素索引的工具函数
  def coord_to_pixel(x, y, gt):
      px = int((x - gt[0]) / gt[1])
      py = int((y - gt[3]) / gt[5])
      return px, py
  
  # 逐个提取点的栅格值（仅读取对应单个像素）
  points_gdf["raster_value"] = points_gdf.apply(
      lambda row: raster_band.ReadAsArray(
          *coord_to_pixel(row.geometry.x, row.geometry.y, geo_transform),
          1, 1
      )[0][0],
      axis=1
  )
  
  # 保存结果为矢量文件
  points_gdf.to_file("sampled_points_with_values.shp")
  

2. ArcPy 实现 #

适合ArcGIS环境下的快速操作，工具内置优化：

• 先生成多边形约束的随机点：

  arcpy.CreateRandomPoints_management(
      out_path="output.gdb",
      out_name="random_sample_points",
      constraining_feature_class="polygons.shp",
      number_of_points_or_field=1000
  )
  

• 提取点的栅格值（工具仅处理点对应栅格区域，不加载全量栅格）：

  arcpy.sa.ExtractValuesToPoints(
      in_point_features="random_sample_points",
      in_raster="large_raster.tif",
      out_point_features="sampled_points_with_values"
  )
  

3. QGIS 图形界面操作 #

无需代码，适合快速验证和小批量处理：

1. 加载多边形矢量和超大栅格（QGIS会按需加载栅格瓦片，不会一次性读取全量数据）。
2. 打开「矢量」→「研究工具」→「创建随机点」，选择多边形作为约束图层，设置采样数量，生成随机点图层。
3. 打开「栅格」→「提取」→「提取点采样值」，选择生成的随机点和目标栅格，设置输出路径后运行即可。

额外优化建议 #

• 如果多边形数量极多，先合并为单个几何对象能减少随机点生成时的判断次数。
• 若采样点数量超过10万级，建议分批次生成和提取，避免内存过载。
• 确保多边形和栅格的坐标系一致，避免坐标转换带来的额外开销。

内容的提问来源于stack exchange，提问作者tjampolay