import numpy as np
import open3d as o3d
import matplotlib.pyplot as plt
import os
from matplotlib.gridspec import GridSpec
# 设置中文字体为宋体
plt.rcParams["font.family"] = ["SimSun"]
plt.rcParams['axes.unicode_minus'] = False # 解决负号显示问题
class SimplifiedPointCloudProcessor:
def __init__(self, random_seed=42, output_dir="output_models"):
"""初始化点云处理器,使用固定参数计算法向量"""
self.random_seed = random_seed
self.output_dir = output_dir # 输出目录
np.random.seed(random_seed)
# 创建输出目录(如果不存在)
os.makedirs(self.output_dir, exist_ok=True)
# 初始化数据结构
self.data = None
self.points = None
self.labels = None
self.unique_labels = {}
# 可视化几何体存储
self.geometries = []
self.geometries_pcd = []
# 处理参数
self.cluster_index = 0
self.min_points_threshold = 20
self.small_cluster_threshold = 10
self.outlier_std_ratio = 1.2
# 法线计算参数 - 使用固定KNN参数
self.knn_param = 30 # 固定KNN参数
self.consistent_knn = 12 # 法线一致性检查的邻点数
# 可视化参数
self.visualization_pause_time = 0 # 可视化暂停时间,0表示无限等待
self.coord_offset_ratio = 0.3 # 坐标系偏移比例
self.coord_position_choice = 0 # 0:左下, 1:右下, 2:左上, 3:右上
# 三视图参数
self.viewpoint_size = 100 # 三视图点大小
self.save_views = True # 是否保存三视图
def load_data(self, file_path):
"""加载点云数据"""
try:
self.data = np.genfromtxt(
file_path,
delimiter=',',
skip_header=1
)
self.points = self.data[:, :3]
self.labels = self.data[:, 3].astype(int)
for label in np.unique(self.labels):
self.unique_labels[label] = self.points[self.labels == label]
# 计算点云整体边界框,用于放置坐标系
self._calculate_overall_bounding_box()
print(f"成功加载数据,共 {len(self.points)} 个点,{len(self.unique_labels)} 个类簇")
return True
except Exception as e:
print(f"数据加载失败: {str(e)}")
return False
def _calculate_overall_bounding_box(self):
"""计算整个点云的边界框,用于确定坐标系位置"""
if self.points is not None and len(self.points) > 0:
self.overall_min = np.min(self.points, axis=0)
self.overall_max = np.max(self.points, axis=0)
self.overall_center = np.mean(self.points, axis=0)
self.overall_size = self.overall_max - self.overall_min
else:
self.overall_min = np.array([0, 0, 0])
self.overall_max = np.array([100, 100, 100])
self.overall_center = np.array([50, 50, 50])
self.overall_size = np.array([100, 100, 100])
def _filter_outliers(self, pcd):
"""过滤离群点"""
cl, ind = pcd.remove_statistical_outlier(
nb_neighbors=15,
std_ratio=self.outlier_std_ratio
)
return cl, ind
def _orient_normals_externally(self, pcd):
"""将法线方向统一向外"""
points = np.asarray(pcd.points)
normals = np.asarray(pcd.normals)
# 使用包围盒中心方法
bbox = pcd.get_axis_aligned_bounding_box()
center = bbox.get_center()
to_center = points - center
dot_products = np.sum(normals * to_center, axis=1)
flip_mask = dot_products > 0 # 指向中心的需要翻转
normals[flip_mask] *= -1
pcd.normals = o3d.utility.Vector3dVector(normals)
return pcd
def _create_colored_coordinate_frame(self, size=1.0):
"""创建自定义颜色的坐标系:x轴蓝色,y轴绿色,z轴红色"""
frame = o3d.geometry.TriangleMesh.create_coordinate_frame(size=size)
# 设置坐标轴颜色
frame.vertex_colors[0] = [0, 0, 1] # x轴顶点1 - 蓝色
frame.vertex_colors[1] = [0, 0, 1] # x轴顶点2 - 蓝色
frame.vertex_colors[2] = [0, 1, 0] # y轴顶点1 - 绿色
frame.vertex_colors[3] = [0, 1, 0] # y轴顶点2 - 绿色
frame.vertex_colors[4] = [1, 0, 0] # z轴顶点1 - 红色
frame.vertex_colors[5] = [1, 0, 0] # z轴顶点2 - 红色
# 根据选择计算不同的坐标系位置
offset = self.overall_size * self.coord_offset_ratio
if self.coord_position_choice == 0: # 左下
coord_position = self.overall_min - offset
elif self.coord_position_choice == 1: # 右下
coord_position = np.array([self.overall_max[0] + offset[0],
self.overall_min[1] - offset[1],
self.overall_min[2] - offset[2]])
elif self.coord_position_choice == 2: # 左上
coord_position = np.array([self.overall_min[0] - offset[0],
self.overall_max[1] + offset[1],
self.overall_min[2] - offset[2]])
else: # 右上
coord_position = self.overall_max + offset
# 平移坐标系到计算位置
frame.translate(coord_position)
return frame
def _save_3d_model(self, mesh, cluster_index, label):
"""保存三维模型为多种格式"""
if not mesh:
return
# 构建文件名(包含类簇信息)
base_filename = f"Ore_cluster_{label}_sub_{cluster_index}"
# 保存为PLY格式(包含颜色信息)
ply_path = os.path.join(self.output_dir, f"{base_filename}.ply")
o3d.io.write_triangle_mesh(ply_path, mesh)
print(f"已保存PLY模型: {ply_path}")
# 保存为STL格式(广泛用于3D打印)
stl_path = os.path.join(self.output_dir, f"{base_filename}.stl")
o3d.io.write_triangle_mesh(stl_path, mesh)
print(f"已保存STL模型: {stl_path}")
# 保存为OBJ格式(在许多3D软件中兼容)
obj_path = os.path.join(self.output_dir, f"{base_filename}.obj")
o3d.io.write_triangle_mesh(obj_path, mesh)
print(f"已保存OBJ模型: {obj_path}")
def _process_sub_cluster(self, sub_cluster, cluster_index, label):
"""处理单个子簇"""
# 添加微小噪声打破共面性
noise = np.random.normal(0, 1e-6, sub_cluster.shape)
sub_cluster += noise
# 创建点云对象
pcd = o3d.geometry.PointCloud()
pcd.points = o3d.utility.Vector3dVector(sub_cluster)
# 过滤点数太少的子簇
if len(pcd.points) < self.small_cluster_threshold:
print(f"子簇点数不足 {self.small_cluster_threshold},跳过处理")
return None, None
# 不过滤离群点
filtered_pcd = pcd
# 过滤离群点
# filtered_pcd, inliers = self._filter_outliers(pcd)
# if len(filtered_pcd.points) < 0.5 * len(pcd.points):
# print("过滤后点太少,使用原始点云")
# filtered_pcd = pcd
# else:
# print(f"过滤离群点: {len(pcd.points)} -> {len(filtered_pcd.points)}")
print(f"使用KNN参数计算法向量: {self.knn_param}个最近邻")
# 计算法向量 - 使用固定KNN参数
filtered_pcd.estimate_normals(
o3d.geometry.KDTreeSearchParamKNN(knn=self.knn_param),
fast_normal_computation=False # 禁用快速计算,提高精度
)
# 法线方向一致化
filtered_pcd.orient_normals_consistent_tangent_plane(k=self.consistent_knn)
# 设置颜色
color = plt.get_cmap("Accent")(cluster_index % 8)[:3]
filtered_pcd.paint_uniform_color(color)
# 泊松曲面重建
try:
with o3d.utility.VerbosityContextManager(
o3d.utility.VerbosityLevel.Debug
) as cm:
mesh, densities = o3d.geometry.TriangleMesh.create_from_point_cloud_poisson(
filtered_pcd,
depth=8,
width=5,
scale=1.1,
linear_fit=True
)
# 后处理网格
mesh.compute_vertex_normals()
mesh.paint_uniform_color(color)
# 保存三维模型
self._save_3d_model(mesh, cluster_index, label)
return filtered_pcd, mesh
except Exception as e:
print(f"泊松重建失败: {str(e)}")
return filtered_pcd, None
def process_clusters(self):
"""处理所有类簇"""
if not self.unique_labels:
print("没有可处理的类簇数据")
return
# 遍历每个类簇
for label, cluster_points in self.unique_labels.items():
if label == -1: # 跳过噪声点
continue
if len(cluster_points) < self.min_points_threshold:
print(f"类簇 {label} 点数不足 {self.min_points_threshold},跳过处理")
continue
print(f"处理类簇 {label},包含 {len(cluster_points)} 个点")
current_clusters = self._get_sub_clusters(label, cluster_points)
for sub_cluster in current_clusters:
# 传递label参数用于保存文件
pcd, mesh = self._process_sub_cluster(sub_cluster, self.cluster_index, label)
if pcd:
self.geometries_pcd.append(pcd)
if mesh:
self.geometries.append(mesh)
self.cluster_index += 1
print(f"共处理 {self.cluster_index} 个有效类簇")
def _get_sub_clusters(self, label, cluster_points):
"""获取子簇,直接返回原始类簇"""
return [cluster_points]
def visualize_orthographic_views(self):
"""生成并显示三视图(正视图、侧视图、俯视图)"""
if self.points is None or len(self.points) == 0:
print("没有点云数据可生成三视图")
return
# 创建图形和子图
fig = plt.figure(figsize=(15, 10))
gs = GridSpec(2, 2, figure=fig)
# 正视图 (X-Y平面)
ax1 = fig.add_subplot(gs[0, 0])
# 侧视图 (Y-Z平面)
ax2 = fig.add_subplot(gs[0, 1])
# 俯视图 (X-Z平面)
ax3 = fig.add_subplot(gs[1, 0])
# 设置标题
ax1.set_title('正视图 (X-Y平面)')
ax2.set_title('侧视图 (Y-Z平面)')
ax3.set_title('俯视图 (X-Z平面)')
# 设置坐标轴标签
ax1.set_xlabel('X')
ax1.set_ylabel('Y')
ax2.set_xlabel('Y')
ax2.set_ylabel('Z')
ax3.set_xlabel('X')
ax3.set_ylabel('Z')
# 确保各视图比例一致
ax1.set_aspect('equal')
ax2.set_aspect('equal')
ax3.set_aspect('equal')
# 获取颜色映射
cmap = plt.get_cmap("Accent")
# 按类别绘制点
for i, (label, points) in enumerate(self.unique_labels.items()):
if label == -1: # 跳过噪声点
continue
color = cmap(i % 8)[:3]
# 正视图 (X-Y)
ax1.scatter(points[:, 0], points[:, 1], c=[color], s=self.viewpoint_size, alpha=0.6,
label=f'类别 {label}')
# 侧视图 (Y-Z)
ax2.scatter(points[:, 1], points[:, 2], c=[color], s=self.viewpoint_size, alpha=0.6)
# 俯视图 (X-Z)
ax3.scatter(points[:, 0], points[:, 2], c=[color], s=self.viewpoint_size, alpha=0.6)
# 添加图例
ax1.legend(bbox_to_anchor=(1.05, 1), loc='upper left')
# 调整布局
plt.tight_layout()
# 保存三视图(如果启用)
if self.save_views:
views_path = os.path.join(self.output_dir, "point_cloud_orthographic_views.png")
plt.savefig(views_path, dpi=300, bbox_inches='tight')
print(f"已保存三视图: {views_path}")
# 显示三视图
plt.show()
def visualize_results(self):
"""可视化处理结果,包括3D视图和三视图"""
# 先显示三视图
self.visualize_orthographic_views()
if not self.geometries_pcd and not self.geometries:
print("没有可可视化的3D结果")
return
# 创建自定义颜色的坐标系
coord_frame_size = np.max(self.overall_size) * 0.1
colored_coord_frame = self._create_colored_coordinate_frame(size=coord_frame_size)
# 可视化点云(显示法线)
if self.geometries_pcd:
print("显示点云与法线方向 (按Q键关闭窗口)")
vis = o3d.visualization.Visualizer()
vis.create_window(window_name="点云与法线方向", width=1024, height=768)
# 添加自定义坐标系
vis.add_geometry(colored_coord_frame)
# 添加点云
for pcd in self.geometries_pcd:
vis.add_geometry(pcd)
# 设置渲染参数
opt = vis.get_render_option()
opt.point_size = 3
opt.line_width = 1
opt.show_coordinate_frame = False
# 运行可视化
vis.run()
vis.destroy_window()
# 可视化点云和重建网格
if self.geometries:
print("显示点云和泊松重建网格 (按Q键关闭窗口)")
vis = o3d.visualization.Visualizer()
vis.create_window(window_name="点云和泊松重建网格", width=1024, height=768)
# 添加自定义坐标系
vis.add_geometry(colored_coord_frame)
# 添加几何体
for geom in self.geometries + self.geometries_pcd:
vis.add_geometry(geom)
opt = vis.get_render_option()
opt.mesh_show_back_face = True
opt.point_size = 2
opt.show_coordinate_frame = False
vis.run()
vis.destroy_window()
if __name__ == "__main__":
# 创建处理器实例,指定输出目录
processor = SimplifiedPointCloudProcessor(output_dir="3d_models_output")
# 加载数据(请替换为实际文件路径)
file_path = 'midpoints_with_cluster_index_SandStone_filtered_Test_20_2404080.csv'
if processor.load_data(file_path):
# 处理类簇(处理过程中会自动保存模型)
processor.process_clusters()
# 可视化结果(现在包括三视图)
processor.visualize_results()
else:
print("无法继续处理,数据加载失败")
确保保存的体模型为闭合模型
最新发布