1. 引言为什么需要 SciPy 处理空间数据在数据科学和工程领域空间数据Spatial Data无处不在——从地理信息系统GIS中的经纬度坐标到医学影像中的三维体素再到计算机图形学中的点云和网格。处理这类数据需要专门的数学工具和算法而 SciPy 正是 Python 科学计算生态中处理空间数据的利器。SciPy 的scipy.spatial模块提供了一系列高效的空间数据结构和算法包括距离计算欧氏距离、曼哈顿距离、余弦距离等空间索引KDTree、cKDTree 用于快速最近邻搜索几何结构凸包Convex Hull、Voronoi 图、Delaunay 三角剖分空间变换旋转矩阵、四元数、刚体变换本文将带你深入scipy.spatial的各个功能通过丰富的代码实例展示如何在实际项目中应用这些工具。2. 环境准备与数据生成首先确保已安装 SciPy。如果尚未安装可以使用 pippip install scipy numpy matplotlib让我们生成一些示例空间数据用于后续演示import numpy as np import matplotlib.pyplot as plt from scipy import spatial 生成随机二维点集 np.random.seed(42) points_2d np.random.rand(50, 2) * 10 # 50个点范围[0, 10] 生成随机三维点集 points_3d np.random.rand(30, 3) * 5 # 30个点范围[0, 5] 生成带有噪声的曲线数据 t np.linspace(0, 4*np.pi, 100) curve_points np.column_stack([t, np.sin(t) np.random.normal(0, 0.1, 100)]) print(f二维点集形状: {points_2d.shape}) print(f三维点集形状: {points_3d.shape}) print(f曲线点集形状: {curve_points.shape})3. 距离计算与度量scipy.spatial.distance模块提供了多种距离度量函数。让我们看看几个常用距离的计算方法from scipy.spatial import distance 定义两个向量 v1 np.array([1, 2, 3, 4]) v2 np.array([4, 3, 2, 1]) 计算各种距离 print(欧氏距离:, distance.euclidean(v1, v2)) print(曼哈顿距离:, distance.cityblock(v1, v2)) print(余弦距离:, distance.cosine(v1, v2)) print(切比雪夫距离:, distance.chebyshev(v1, v2)) print(闵可夫斯基距离(p3):, distance.minkowski(v1, v2, p3)) 计算点集之间的成对距离矩阵 points np.array([[0, 0], [1, 1], [2, 2], [3, 3]]) dist_matrix distance.pdist(points, euclidean) print(\n压缩距离矩阵:, dist_matrix) print(方形距离矩阵:\n, distance.squareform(dist_matrix))4. KDTree高效的空间索引与最近邻搜索KDTreek-dimensional tree是一种用于组织 k 维空间中点的空间数据结构能够大幅加速最近邻搜索。SciPy 提供了两种实现纯 Python 的KDTree和 C 实现的cKDTree更快。# 使用 cKDTree 构建空间索引 tree spatial.cKDTree(points_2d) 查询单个点的最近邻 query_point np.array([5, 5]) dist, idx tree.query(query_point, k3) # 查找最近的3个点 print(f查询点 {query_point} 的最近3个点索引: {idx}) print(f对应距离: {dist}) 批量查询多个点 query_points np.array([[2, 2], [8, 8], [1, 9]]) batch_dist, batch_idx tree.query(query_points, k2) print(f\n批量查询结果:) for i, (dists, idxs) in enumerate(zip(batch_dist, batch_idx)): print(f 查询点 {query_points[i]}: 最近点索引 {idxs}, 距离 {dists}) 半径搜索查找指定半径内的所有点 indices tree.query_ball_point([5, 5], r2.0) print(f\n半径2.0内的点索引: {indices}) print(f这些点的坐标:\n{points_2d[indices]}) 两棵树之间的所有点对距离小于阈值的点 tree2 spatial.cKDTree(points_2d 0.5) # 稍微偏移的树 pairs tree.query_ball_tree(tree2, r1.5) print(f\n树间距离小于1.5的点对数量: {sum(len(p) for p in pairs)})5. 凸包计算与应用凸包Convex Hull是包含一组点的最小凸多边形。在计算机图形学、模式识别和路径规划中都有广泛应用。# 计算二维点集的凸包 hull_2d spatial.ConvexHull(points_2d) print(f凸包顶点索引: {hull_2d.vertices}) print(f凸包顶点坐标:\n{points_2d[hull_2d.vertices]}) print(f凸包面积: {hull_2d.volume:.2f}) # 二维时volume表示面积 print(f凸包周长: {hull_2d.area:.2f}) # 二维时area表示周长 可视化凸包 plt.figure(figsize(10, 5)) plt.subplot(1, 2, 1) plt.plot(points_2d[:, 0], points_2d[:, 1], o, labelPoints) for simplex in hull_2d.simplices: plt.plot(points_2d[simplex, 0], points_2d[simplex, 1], r-) plt.title(2D Convex Hull) plt.legend() 计算三维点集的凸包 hull_3d spatial.ConvexHull(points_3d) print(f\n三维凸包顶点数: {len(hull_3d.vertices)}) print(f三维凸包体积: {hull_3d.volume:.2f}) print(f三维凸包表面积: {hull_3d.area:.2f}) plt.subplot(1, 2, 2, projection3d) plt.gca().scatter(points_3d[:, 0], points_3d[:, 1], points_3d[:, 2]) for simplex in hull_3d.simplices: for i in range(3): for j in range(i1, 3): plt.plot([points_3d[simplex[i], 0], points_3d[simplex[j], 0]], [points_3d[simplex[i], 1], points_3d[simplex[j], 1]], [points_3d[simplex[i], 2], points_3d[simplex[j], 2]], r-) plt.title(3D Convex Hull) plt.tight_layout() plt.show()6. Delaunay 三角剖分与 Voronoi 图Delaunay 三角剖分和 Voronoi 图是计算几何中的一对对偶概念在有限元分析、地理信息系统和模式识别中广泛应用。# Delaunay 三角剖分 tri spatial.Delaunay(points_2d[:20]) # 使用前20个点避免过于密集 print(f三角形数量: {tri.simplices.shape[0]}) print(f前5个三角形顶点索引:\n{tri.simplices[:5]}) Voronoi 图 vor spatial.Voronoi(points_2d[:15]) # 使用前15个点 print(f\nVoronoi 区域数量: {len(vor.point_region)}) print(fVoronoi 顶点数量: {len(vor.vertices)}) 可视化 plt.figure(figsize(12, 5)) plt.subplot(1, 2, 1) plt.triplot(points_2d[:20, 0], points_2d[:20, 1], tri.simplices, g-, alpha0.5) plt.plot(points_2d[:20, 0], points_2d[:20, 1], o) plt.title(Delaunay Triangulation) plt.subplot(1, 2, 2) spatial.voronoi_plot_2d(vor, axplt.gca(), show_verticesFalse) plt.plot(points_2d[:15, 0], points_2d[:15, 1], ro) plt.title(Voronoi Diagram) plt.tight_layout() plt.show() 查找包含给定点的三角形 test_point np.array([5, 5]) triangle_index tri.find_simplex(test_point) if triangle_index ! -1: print(f\n点 {test_point} 在三角形 {triangle_index} 内) print(f三角形顶点: {tri.simplices[triangle_index]}) else: print(f点 {test_point} 不在任何三角形内)7. 空间变换旋转与对齐SciPy 也提供了一些空间变换工具特别是在scipy.spatial.transform模块中from scipy.spatial.transform import Rotation 创建旋转对象 从欧拉角创建ZYX顺序单位度 rot1 Rotation.from_euler(zyx, [30, 45, 60], degreesTrue) 从旋转矩阵创建 rot_matrix np.array([[0, -1, 0], [1, 0, 0], [0, 0, 1]]) rot2 Rotation.from_matrix(rot_matrix) 从四元数创建 rot3 Rotation.from_quat([0, 0, np.sin(np.pi/4), np.cos(np.pi/4)]) 应用旋转 vector np.array([1, 0, 0]) rotated_vector rot1.apply(vector) print(f向量 {vector} 旋转后: {rotated_vector}) 旋转组合 rot_combined rot1 * rot2 # 先rot2再rot1 print(f组合旋转的欧拉角: {rot_combined.as_euler(zyx, degreesTrue)}) 旋转插值SLERP rot_start Rotation.from_euler(z, 0, degreesTrue) rot_end Rotation.from_euler(z, 90, degreesTrue) times np.linspace(0, 1, 5) rot_interp Rotation.concatenate([rot_start.slerp(rot_end, t) for t in times]) print(\n球面线性插值结果:) for i, rot in enumerate(rot_interp): print(f t{times[i]:.2f}: 欧拉角{rot.as_euler(z, degreesTrue)[0]:.1f}°)8. 实战案例点云配准与表面重建让我们结合多个 SciPy 空间功能实现一个简单的点云配准和表面重建示例def point_cloud_registration(source_points, target_points): 使用 Procrustes 分析进行点云配准 # 中心化 source_centered source_points - np.mean(source_points, axis0) target_centered target_points - np.mean(target_points, axis0) # 计算旋转矩阵使用SVD H source_centered.T target_centered U, _, Vt np.linalg.svd(H) R Vt.T U.T 确保是右手系旋转det(R) 1 if np.linalg.det(R) 0: Vt[-1, :] * -1 R Vt.T U.T 计算平移 t np.mean(target_points, axis0) - np.mean(source_points, axis0) R return R, t def reconstruct_surface_from_points(points, alpha0.1): 使用 alpha shape 进行表面重建简化示例 计算凸包作为基础表面 hull spatial.ConvexHull(points) 使用 Delaunay 三角剖分 tri spatial.Delaunay(points) 过滤太长的边简化版的 alpha shape triangles [] for simplex in tri.simplices: # 计算三角形边长 edge_lengths [] for i in range(3): for j in range(i1, 3): dist np.linalg.norm(points[simplex[i]] - points[simplex[j]]) edge_lengths.append(dist) # 如果所有边都小于阈值保留该三角形 if max(edge_lengths) lt; alpha * np.max(edge_lengths): triangles.append(simplex) return np.array(triangles) 生成测试点云一个扭曲的球面 phi np.linspace(0, np.pi, 20) theta np.linspace(0, 2*np.pi, 40) phi, theta np.meshgrid(phi, theta) x np.sin(phi) * np.cos(theta) np.random.normal(0, 0.05, phi.shape) y np.sin(phi) * np.sin(theta) np.random.normal(0, 0.05, phi.shape) z np.cos(phi) np.random.normal(0, 0.05, phi.shape) points_cloud np.column_stack([x.ravel(), y.ravel(), z.ravel()]) 表面重建 triangles reconstruct_surface_from_points(points_cloud[:200], alpha0.15) print(f原始点数量: {len(points_cloud)}) print(f重建表面三角形数量: {len(triangles)}) 可视化 fig plt.figure(figsize(10, 8)) ax fig.add_subplot(111, projection3d) ax.scatter(points_cloud[:200, 0], points_cloud[:200, 1], points_cloud[:200, 2], cblue, alpha0.3, s10, labelPoints) 绘制重建的表面 for triangle in triangles[:100]: # 只绘制前100个三角形避免过于密集 for i in range(3): for j in range(i1, 3): ax.plot([points_cloud[triangle[i], 0], points_cloud[triangle[j], 0]], [points_cloud[triangle[i], 1], points_cloud[triangle[j], 1]], [points_cloud[triangle[i], 2], points_cloud[triangle[j], 2]], r-, alpha0.2) ax.set_xlabel(X) ax.set_ylabel(Y) ax.set_zlabel(Z) ax.set_title(Point Cloud Surface Reconstruction) ax.legend() plt.show()9. 性能优化技巧与最佳实践处理大规模空间数据时性能至关重要。以下是一些优化建议import time 1. 使用 cKDTree 而不是 KDTree points_large np.random.rand(100000, 3) * 100 # 10万个3D点 start time.time() tree_slow spatial.KDTree(points_large[:10000]) # 只使用1万个点 print(fKDTree 构建时间: {time.time() - start:.3f}s) start time.time() tree_fast spatial.cKDTree(points_large[:10000]) print(fcKDTree 构建时间: {time.time() - start:.3f}s) 2. 批量查询优于循环单次查询 query_points np.random.rand(1000, 3) * 100 start time.time() for q in query_points: tree_fast.query(q, k1) print(f循环单次查询时间: {time.time() - start:.3f}s) start time.time() tree_fast.query(query_points, k1) print(f批量查询时间: {time.time() - start:.3f}s) 3. 适当选择叶子大小leafsize start time.time() tree_default spatial.cKDTree(points_large[:50000]) print(f默认 leafsize 构建时间: {time.time() - start:.3f}s) SciPy 空间数据处理全攻略：从基础操作到实战案例
1. 引言为什么需要 SciPy 处理空间数据在数据科学和工程领域空间数据Spatial Data无处不在——从地理信息系统GIS中的经纬度坐标到医学影像中的三维体素再到计算机图形学中的点云和网格。处理这类数据需要专门的数学工具和算法而 SciPy 正是 Python 科学计算生态中处理空间数据的利器。SciPy 的scipy.spatial模块提供了一系列高效的空间数据结构和算法包括距离计算欧氏距离、曼哈顿距离、余弦距离等空间索引KDTree、cKDTree 用于快速最近邻搜索几何结构凸包Convex Hull、Voronoi 图、Delaunay 三角剖分空间变换旋转矩阵、四元数、刚体变换本文将带你深入scipy.spatial的各个功能通过丰富的代码实例展示如何在实际项目中应用这些工具。2. 环境准备与数据生成首先确保已安装 SciPy。如果尚未安装可以使用 pippip install scipy numpy matplotlib让我们生成一些示例空间数据用于后续演示import numpy as np import matplotlib.pyplot as plt from scipy import spatial 生成随机二维点集 np.random.seed(42) points_2d np.random.rand(50, 2) * 10 # 50个点范围[0, 10] 生成随机三维点集 points_3d np.random.rand(30, 3) * 5 # 30个点范围[0, 5] 生成带有噪声的曲线数据 t np.linspace(0, 4*np.pi, 100) curve_points np.column_stack([t, np.sin(t) np.random.normal(0, 0.1, 100)]) print(f二维点集形状: {points_2d.shape}) print(f三维点集形状: {points_3d.shape}) print(f曲线点集形状: {curve_points.shape})3. 距离计算与度量scipy.spatial.distance模块提供了多种距离度量函数。让我们看看几个常用距离的计算方法from scipy.spatial import distance 定义两个向量 v1 np.array([1, 2, 3, 4]) v2 np.array([4, 3, 2, 1]) 计算各种距离 print(欧氏距离:, distance.euclidean(v1, v2)) print(曼哈顿距离:, distance.cityblock(v1, v2)) print(余弦距离:, distance.cosine(v1, v2)) print(切比雪夫距离:, distance.chebyshev(v1, v2)) print(闵可夫斯基距离(p3):, distance.minkowski(v1, v2, p3)) 计算点集之间的成对距离矩阵 points np.array([[0, 0], [1, 1], [2, 2], [3, 3]]) dist_matrix distance.pdist(points, euclidean) print(\n压缩距离矩阵:, dist_matrix) print(方形距离矩阵:\n, distance.squareform(dist_matrix))4. KDTree高效的空间索引与最近邻搜索KDTreek-dimensional tree是一种用于组织 k 维空间中点的空间数据结构能够大幅加速最近邻搜索。SciPy 提供了两种实现纯 Python 的KDTree和 C 实现的cKDTree更快。# 使用 cKDTree 构建空间索引 tree spatial.cKDTree(points_2d) 查询单个点的最近邻 query_point np.array([5, 5]) dist, idx tree.query(query_point, k3) # 查找最近的3个点 print(f查询点 {query_point} 的最近3个点索引: {idx}) print(f对应距离: {dist}) 批量查询多个点 query_points np.array([[2, 2], [8, 8], [1, 9]]) batch_dist, batch_idx tree.query(query_points, k2) print(f\n批量查询结果:) for i, (dists, idxs) in enumerate(zip(batch_dist, batch_idx)): print(f 查询点 {query_points[i]}: 最近点索引 {idxs}, 距离 {dists}) 半径搜索查找指定半径内的所有点 indices tree.query_ball_point([5, 5], r2.0) print(f\n半径2.0内的点索引: {indices}) print(f这些点的坐标:\n{points_2d[indices]}) 两棵树之间的所有点对距离小于阈值的点 tree2 spatial.cKDTree(points_2d 0.5) # 稍微偏移的树 pairs tree.query_ball_tree(tree2, r1.5) print(f\n树间距离小于1.5的点对数量: {sum(len(p) for p in pairs)})5. 凸包计算与应用凸包Convex Hull是包含一组点的最小凸多边形。在计算机图形学、模式识别和路径规划中都有广泛应用。# 计算二维点集的凸包 hull_2d spatial.ConvexHull(points_2d) print(f凸包顶点索引: {hull_2d.vertices}) print(f凸包顶点坐标:\n{points_2d[hull_2d.vertices]}) print(f凸包面积: {hull_2d.volume:.2f}) # 二维时volume表示面积 print(f凸包周长: {hull_2d.area:.2f}) # 二维时area表示周长 可视化凸包 plt.figure(figsize(10, 5)) plt.subplot(1, 2, 1) plt.plot(points_2d[:, 0], points_2d[:, 1], o, labelPoints) for simplex in hull_2d.simplices: plt.plot(points_2d[simplex, 0], points_2d[simplex, 1], r-) plt.title(2D Convex Hull) plt.legend() 计算三维点集的凸包 hull_3d spatial.ConvexHull(points_3d) print(f\n三维凸包顶点数: {len(hull_3d.vertices)}) print(f三维凸包体积: {hull_3d.volume:.2f}) print(f三维凸包表面积: {hull_3d.area:.2f}) plt.subplot(1, 2, 2, projection3d) plt.gca().scatter(points_3d[:, 0], points_3d[:, 1], points_3d[:, 2]) for simplex in hull_3d.simplices: for i in range(3): for j in range(i1, 3): plt.plot([points_3d[simplex[i], 0], points_3d[simplex[j], 0]], [points_3d[simplex[i], 1], points_3d[simplex[j], 1]], [points_3d[simplex[i], 2], points_3d[simplex[j], 2]], r-) plt.title(3D Convex Hull) plt.tight_layout() plt.show()6. Delaunay 三角剖分与 Voronoi 图Delaunay 三角剖分和 Voronoi 图是计算几何中的一对对偶概念在有限元分析、地理信息系统和模式识别中广泛应用。# Delaunay 三角剖分 tri spatial.Delaunay(points_2d[:20]) # 使用前20个点避免过于密集 print(f三角形数量: {tri.simplices.shape[0]}) print(f前5个三角形顶点索引:\n{tri.simplices[:5]}) Voronoi 图 vor spatial.Voronoi(points_2d[:15]) # 使用前15个点 print(f\nVoronoi 区域数量: {len(vor.point_region)}) print(fVoronoi 顶点数量: {len(vor.vertices)}) 可视化 plt.figure(figsize(12, 5)) plt.subplot(1, 2, 1) plt.triplot(points_2d[:20, 0], points_2d[:20, 1], tri.simplices, g-, alpha0.5) plt.plot(points_2d[:20, 0], points_2d[:20, 1], o) plt.title(Delaunay Triangulation) plt.subplot(1, 2, 2) spatial.voronoi_plot_2d(vor, axplt.gca(), show_verticesFalse) plt.plot(points_2d[:15, 0], points_2d[:15, 1], ro) plt.title(Voronoi Diagram) plt.tight_layout() plt.show() 查找包含给定点的三角形 test_point np.array([5, 5]) triangle_index tri.find_simplex(test_point) if triangle_index ! -1: print(f\n点 {test_point} 在三角形 {triangle_index} 内) print(f三角形顶点: {tri.simplices[triangle_index]}) else: print(f点 {test_point} 不在任何三角形内)7. 空间变换旋转与对齐SciPy 也提供了一些空间变换工具特别是在scipy.spatial.transform模块中from scipy.spatial.transform import Rotation 创建旋转对象 从欧拉角创建ZYX顺序单位度 rot1 Rotation.from_euler(zyx, [30, 45, 60], degreesTrue) 从旋转矩阵创建 rot_matrix np.array([[0, -1, 0], [1, 0, 0], [0, 0, 1]]) rot2 Rotation.from_matrix(rot_matrix) 从四元数创建 rot3 Rotation.from_quat([0, 0, np.sin(np.pi/4), np.cos(np.pi/4)]) 应用旋转 vector np.array([1, 0, 0]) rotated_vector rot1.apply(vector) print(f向量 {vector} 旋转后: {rotated_vector}) 旋转组合 rot_combined rot1 * rot2 # 先rot2再rot1 print(f组合旋转的欧拉角: {rot_combined.as_euler(zyx, degreesTrue)}) 旋转插值SLERP rot_start Rotation.from_euler(z, 0, degreesTrue) rot_end Rotation.from_euler(z, 90, degreesTrue) times np.linspace(0, 1, 5) rot_interp Rotation.concatenate([rot_start.slerp(rot_end, t) for t in times]) print(\n球面线性插值结果:) for i, rot in enumerate(rot_interp): print(f t{times[i]:.2f}: 欧拉角{rot.as_euler(z, degreesTrue)[0]:.1f}°)8. 实战案例点云配准与表面重建让我们结合多个 SciPy 空间功能实现一个简单的点云配准和表面重建示例def point_cloud_registration(source_points, target_points): 使用 Procrustes 分析进行点云配准 # 中心化 source_centered source_points - np.mean(source_points, axis0) target_centered target_points - np.mean(target_points, axis0) # 计算旋转矩阵使用SVD H source_centered.T target_centered U, _, Vt np.linalg.svd(H) R Vt.T U.T 确保是右手系旋转det(R) 1 if np.linalg.det(R) 0: Vt[-1, :] * -1 R Vt.T U.T 计算平移 t np.mean(target_points, axis0) - np.mean(source_points, axis0) R return R, t def reconstruct_surface_from_points(points, alpha0.1): 使用 alpha shape 进行表面重建简化示例 计算凸包作为基础表面 hull spatial.ConvexHull(points) 使用 Delaunay 三角剖分 tri spatial.Delaunay(points) 过滤太长的边简化版的 alpha shape triangles [] for simplex in tri.simplices: # 计算三角形边长 edge_lengths [] for i in range(3): for j in range(i1, 3): dist np.linalg.norm(points[simplex[i]] - points[simplex[j]]) edge_lengths.append(dist) # 如果所有边都小于阈值保留该三角形 if max(edge_lengths) lt; alpha * np.max(edge_lengths): triangles.append(simplex) return np.array(triangles) 生成测试点云一个扭曲的球面 phi np.linspace(0, np.pi, 20) theta np.linspace(0, 2*np.pi, 40) phi, theta np.meshgrid(phi, theta) x np.sin(phi) * np.cos(theta) np.random.normal(0, 0.05, phi.shape) y np.sin(phi) * np.sin(theta) np.random.normal(0, 0.05, phi.shape) z np.cos(phi) np.random.normal(0, 0.05, phi.shape) points_cloud np.column_stack([x.ravel(), y.ravel(), z.ravel()]) 表面重建 triangles reconstruct_surface_from_points(points_cloud[:200], alpha0.15) print(f原始点数量: {len(points_cloud)}) print(f重建表面三角形数量: {len(triangles)}) 可视化 fig plt.figure(figsize(10, 8)) ax fig.add_subplot(111, projection3d) ax.scatter(points_cloud[:200, 0], points_cloud[:200, 1], points_cloud[:200, 2], cblue, alpha0.3, s10, labelPoints) 绘制重建的表面 for triangle in triangles[:100]: # 只绘制前100个三角形避免过于密集 for i in range(3): for j in range(i1, 3): ax.plot([points_cloud[triangle[i], 0], points_cloud[triangle[j], 0]], [points_cloud[triangle[i], 1], points_cloud[triangle[j], 1]], [points_cloud[triangle[i], 2], points_cloud[triangle[j], 2]], r-, alpha0.2) ax.set_xlabel(X) ax.set_ylabel(Y) ax.set_zlabel(Z) ax.set_title(Point Cloud Surface Reconstruction) ax.legend() plt.show()9. 性能优化技巧与最佳实践处理大规模空间数据时性能至关重要。以下是一些优化建议import time 1. 使用 cKDTree 而不是 KDTree points_large np.random.rand(100000, 3) * 100 # 10万个3D点 start time.time() tree_slow spatial.KDTree(points_large[:10000]) # 只使用1万个点 print(fKDTree 构建时间: {time.time() - start:.3f}s) start time.time() tree_fast spatial.cKDTree(points_large[:10000]) print(fcKDTree 构建时间: {time.time() - start:.3f}s) 2. 批量查询优于循环单次查询 query_points np.random.rand(1000, 3) * 100 start time.time() for q in query_points: tree_fast.query(q, k1) print(f循环单次查询时间: {time.time() - start:.3f}s) start time.time() tree_fast.query(query_points, k1) print(f批量查询时间: {time.time() - start:.3f}s) 3. 适当选择叶子大小leafsize start time.time() tree_default spatial.cKDTree(points_large[:50000]) print(f默认 leafsize 构建时间: {time.time() - start:.3f}s)