前言
在森林资源监测中，点云数据已经成为提取森林结构参数（胸径、树高）的重要数据源。无论是来自机载激光雷达（ALS），还是地面激光雷达（TLS），在预处理后，有一个步骤几乎永远排在第一位——地面点分类。那你有没有怀疑过为何非要做地面点分类，聪明的你肯定会想到，因为我要提取的是植被的结构参数，所以当然要先分离地面点得到植被点，以防止干扰后续的数据分析咯。对的，但还要更深层的原因。
今天我们就从这个基础讲起，并带你理解一个非常直观的滤波方法：CSF（Cloth Simulation Filtering）布料模拟算法。
为什么要做地面点滤波？
在森林场景中，地形起伏较大，不同树木的生长位置、高度等信息都与地面密切相关。如果不先提取出地面点，后续的很多关键参数都无从计算：
- 树高：你不能直接用点云中的最高点，而是需要相对于地面的高度。
- 胸径（DBH）：通常在离地1.3米处测量，如果不对点云进行“地面归一化”，这个高度没法统一。
-
数字高程模型（DEM）：地面点就是DEM的基础，没有它就无法还原真实地形。
所以，可以说，地面点滤波是森林点云处理中的基础
地面点分类 --- 布料模拟滤波算法
布料模拟滤波算法CSF（ Cloth Simulation Filtering），是一种用布料模拟方法来提取地面点的算法。听起来是不是很奇怪？布料？
直接来看看它是怎么做的：
- 先把点云倒过来：将点云在Z轴上反转（z → -z），相当于把森林场景点云翻转过来，如下图(a)所示。
- 布料模拟与掉落：在反转后的点云上方初始化布料网格，并进行布料下垂模拟。布料在重力和内部张力的作用下逐步贴合地面点云中凸起部分，形成近似的地面。
-
地面点判断：计算原始点云中各点到布料曲面的距离，根据设定的距离阈值（如0.5m）区分地面点与非地面点。
通俗来说：“CSF 就像先将点云倒置，导致点云的上方盖块布估算地形，再看点离布近不近，近的就是地面点。”
CSF地面点分类实践
更详细内容见链接：https://github.com/jianboqi/CSF，漆老师除了这个CSF算法，还做了一个牛逼的三维辐射传输模型的软件（LESS）https://lessrt.org/，用于三位场景的构建和遥感数据的模拟。
LESS主页介绍：LESS是基于光线追踪的异质地表三维辐射传输模型（计算机模拟模型），可用于复杂森林、城市等地表结构下的多光谱、高光谱、激光雷达、热红外、叶绿素荧光等遥感数据，以及三维光分布、光合作用等数据的模拟。
CSF安装：
pip install cloth-simulation-filter
在实际应用中，我认为CSF 地面点滤波算法主要有以下三个参数需要设置：
1. 高度阈值（class_threshold）
该参数用于判断点是否属于地面点，建议初始设为 0.5 m。如果分类结果中出现过多树干或低矮植被被误判为地面点，可适当调小该阈值（如 0.3 m 或更低），以提高滤波精度。
2. 布料刚性（csf.params.rigidness）
该参数决定布料的柔软程度，从而影响其贴合地形的能力。设置值：
- rigidness = 3：适用于地形相对平缓或无明显坡度的区域；
- rigidness = 2：适用于中等起伏地形，如丘陵、河岸、梯田等；
- rigidness = 1：适用于地形陡峭区域，布料需更柔软以贴合复杂起伏。
3. 坡度后处理（csf.params.bSloopSmooth）
该参数控制是否对地形坡度进行后处理，如果是平地则是false，坡地是true，这个需要自己根据实际情况判断。
更多参数设置可以看这篇文章《An Easy-to-Use Airborne LiDAR Data Filtering Method Based on Cloth Simulation》
import laspy
import numpy as np
import CSF
# 1. 读取原始 las 文件
file_path = r"C:\Users\YSYmc\Desktop\近期内容\lidar_data\ALS_small.las"
als = laspy.read(file_path)
points = np.vstack((als.x, als.y, als.z)).T
# 2. CSF地面点分类
csf = CSF.CSF()
# 参数设置
csf.params.bSloopSmooth = True # 是否进行后处理，平地是TRUE，坡地是false
csf.params.class_threshold = 0.5 # 分类阈值,值越大更多的点被分为地面点
csf.params.cloth_resolution = 0.5 # 网格分辨率，模拟布的网格大小，影响分类精度，默认0.5m即可
csf.params.rigidness = 2 # 值越大布的刚性越大。刚性越高越不容易随地形弯曲，适用于平地
# csf.params.interations = 500 # 模拟迭代次数
# csf.params.time_step = 0.65 # 布料下落时间步长，默认
csf.setPointCloud(points)
ground = CSF.VecInt()
non_ground = CSF.VecInt()
csf.do_filtering(ground, non_ground)
#地面点数组
ground_points=points[np.array(ground)]
# 3. 修改classification属性：地面 = 2，非地面 = 1
als.classification[:] = 0 # 初始化为 0
als.classification[np.array(ground)] = 2
als.classification[np.array(non_ground)] = 1
# 4. 导出las文件
out_path = r"C:\Users\YSYmc\Desktop\近期内容\lidar_data\ALS_small_class.las"
als.write(out_path)
