做网站要域名吗,软件开发用什么笔记本,小城镇建设投稿网站,陕西建设机械官方网站文章目录 前言实现思路想法1想法2想法3 补充实现想法1想法2代码 想法3代码 总结 前言
读前须知#xff1a; 首先我们得确保你已经完全知晓相关的基本的数学知识#xff0c;其中包括用最小二乘法拟合曲二次曲面#xff0c;以及曲面的曲率详细求解。若还是没弄清楚#xff0… 文章目录 前言实现思路想法1想法2想法3 补充实现想法1想法2代码 想法3代码 总结 前言
读前须知 首先我们得确保你已经完全知晓相关的基本的数学知识其中包括用最小二乘法拟合曲二次曲面以及曲面的曲率详细求解。若还是没弄清楚则详细请看下面链接。
【点云、图像】学习中 常见的数学知识及其中的关系与python实战[更新中]建议从一个标题上从上往下看比较循序渐进
补充曲率反映曲面在某一点处的弯曲程度它与该点及其邻近点的位置和法向量有关。
以及一些open3d的常见操作 爆肝5万字❤️Open3D 点云数据处理基础Python版
先上结果: 实现思路
不同于常见的缺陷检测如划痕或者斑点这些肉眼可见的缺陷凹凸性缺陷难以肉眼可见甚至得打光照射才能看见凹槽这里我们使用深度摄像机普通相机深度信息来采集深度信息。此时我们把图像称为深度图当然深度图也可以转换为点云。
这里我们仅对点云这种数据进行处理。
要求对异形曲面微细缺陷识别我6月份之前要完成的毕设 这里缺陷主要指凸起和凹槽。
想法1
想法1如果是一个平面上出现凹槽或凸起的话首先确立一个由大部分点拟合的平面然后对不在此平面的点云进行高程分析以确立凹陷或凸起程度。 事实上不会很平于是想法1排除但可以用来做平面来做一个简单示例
想法2
想法2计算点云上每个点的领域曲率来描绘点的弯曲程度。
想法3
想法3计算点云上每个点的高斯曲率和平均曲率来描绘点的弯曲程度。
补充
机器学习——详解KD-Tree原理
实现
想法1
暂无
想法2
想法2 在求取完领域曲率的基础上我们对其曲率的大小进行一个分割并进行可视化。
这里进一步对这个领域曲率的定义进行详解。 首先我们已经在这里的实战1了解到了领域曲率的求法。 基于邻域特征点提取和匹配的点云配准_李新春 可是这个定义我只在其他的博客中和下面这篇论文中找到定义并无在wiki百科中找到相关阐述。
找了半天后依然无法解释其中原因于是在参考了下面两篇博文后 如何理解矩阵特征值 特征值的最大值与最小值 想出了这么一个合理的解释 1、这个曲率定义的优点是它不依赖于法向量的方向而且它的值域是 [0, 1/3]这使得它比较容易进行归一化和可视化。
2、那么为什么用最小的特征值除以特征值和而不是用最大的特征值除以特征值和呢
这是因为最小的特征值对应的特征向量是曲面的法向量而最大的特征值对应的特征向量是曲面的主方向。
如果用最大的特征值除以特征值和那么曲率的值就会与曲面的主方向的弯曲程度成正比而与曲面的法向方向的弯曲程度无关。这样就会忽略掉曲面的凹凸变化导致曲率的计算不准确。
如果用最小的特征值除以特征值和那么曲率的值就会与曲面的法向方向的弯曲程度成正比而与曲面的主方向的弯曲程度无关。这样就可以反映出曲面的凹凸变化提高曲率的计算精度。
因此用最小的特征值除以特征值和而不是用最大的特征值除以特征值和是为了更好地描述曲面的局部形状而不是曲面的整体方向。
于是我们就可以坦然用这个定义来求取了。
代码
怕点云太多算不过来首先对例子中的兔子点云进行了一个下采样。如何获取兔子点云到时候再出教程
import open3d as o3d
import numpy as npdef pca_compute(data, sortTrue): #1、主成分分析average_data np.mean(data, axis0) # 求每一列的平均值即求各个特征的平均值decentration_matrix data - average_data # 去中心化矩阵H np.dot(decentration_matrix.T, decentration_matrix) # 求协方差矩阵 #协方差是衡量两个变量关系的统计量协方差为正表示两个变量正相关为负表示两个变量负相关eigenvectors, eigenvalues, eigenvectors_T np.linalg.svd(H) # 求特征值与特征向量 #H UΣV^T #输出列向量、对角矩阵、行向量if sort:sort eigenvalues.argsort()[::-1] # 从大到小排序 .argsort()是升序排序[::-1]是将数组反转实现降序排序eigenvalues eigenvalues[sort] # 特征值 ## 使用索引来获取排序后的数组return eigenvaluesdef caculate_surface_curvature(radius,pcd):#2、计算点云的表面曲率cloud pcdpoints np.asarray(cloud.points) #点云转换为数组 点云数组形式为[[x1,y1,z1],[x2,y2,z2],...]kdtree o3d.geometry.KDTreeFlann(cloud) #建立KDTreenum_points len(cloud.points) #点云中点的个数curvature [] # 储存表面曲率for i in range(num_points):k, idx, _ kdtree.search_radius_vector_3d(cloud.points[i], radius) #返回邻域点的个数和索引neighbors points[idx] #数组形式为[[x1,y1,z1],[x2,y2,z2],...]w pca_compute(neighbors)#调用第1步 #由降序排序w[2]为最小特征值 #np.zeros_like(w[2])生成与w[2]相同形状的全0数组delt np.divide(w[2], np.sum(w)) #根据公式求取领域曲率curvature.append(delt)curvature np.array(curvature, dtypenp.float64)return curvaturedef curvature_normal():#3、曲率归一化 从0-1/3归到0-1之间curvature caculate_surface_curvature(radius,pcd) #调用第2步c_max max(curvature)c_min min(curvature)cur_normal [(float(i) - c_min) / (c_max - c_min) for i in curvature] return cur_normaldef draw(cur_max,cur_min,pcd):#4、绘图cur_normal curvature_normal()#调用第3步pcd.paint_uniform_color([0.5,0.5,0.5]) #初始化所有颜色为灰色for i in range(len(cur_normal)):if 0 cur_normal[i] cur_min: #归一化后的曲率np.asarray(pcd.colors)[i] [1, 0, 0]#红elif cur_min cur_normal[i] cur_max:np.asarray(pcd.colors)[i] [0, 1, 0]#绿elif cur_max cur_normal[i] 1: np.asarray(pcd.colors)[i] [0, 0, 1]#蓝# 可视化o3d.visualization.draw_geometries([pcd])cur_max 0.7
cur_min 0.3 #曲率分割基准
radius 0.05
voxel_size 0.01 #越小密度越大
pcd o3d.io.read_point_cloud(bunny.pcd)
print(pcd)
pcd pcd.voxel_down_sample(voxel_size) #下采样
draw(cur_max,cur_min,pcd)
结果 这个长度宽度报错可以不管有点子强迫症的可以在可视化改成 o3d.visualization.draw_geometries([pcd],window_name可视化原始点云,width800, height800, left50, top50,mesh_show_back_faceFalse)红色为曲率较低绿色曲率中等蓝色曲率较高。 发现效果上不太行红色一些部分看着曲率也很高甚至还出现了一个初始化时候的灰点可能求邻近点的时候没取到不曾得知。
想法3
首先我们在这里的高斯曲率和平均曲率求解有了一些认识。 附一张图
代码
import open3d as o3d
import numpy as np
from scipy.optimize import curve_fitvoxel_size 0.01 #越小密度越大
radius 0.07
pcd o3d.io.read_point_cloud(bunny.pcd)
pcd pcd.voxel_down_sample(voxel_size) #下采样
cloud pcd
points np.asarray(cloud.points) #点云转换为数组 点云数组形式为[[x1,y1,z1],[x2,y2,z2],...]
kdtree o3d.geometry.KDTreeFlann(cloud) #建立KDTree
num_points len(cloud.points) #点云中点的个数
pcd.paint_uniform_color([1, 0, 0]) # 初始化所有颜色为红色
# 定义非线性函数这里假设是一个二次曲面
def func(x, a, b, c, d, e, f):return a * x[0]**2 b * x[1]**2 c * x[0] * x[1] d * x[0] e * x[1] f
def f(x, y):return popt[0]*x**2 popt[1]*y**2 popt[2]* x*y popt[3]*x popt[4]*y popt[5]
# 定义曲面的梯度函数即一阶偏导数
def gradient(f, x, y):# 使用中心差分法近似求导 #参考https://cloud.tencent.com/developer/article/1685164h 1e-6 # 差分步长可以根据精度要求调整df_dx (f(x h, y) - f(x - h, y)) / (2 * h) # 对x求偏导df_dy (f(x, y h) - f(x, y - h)) / (2 * h) # 对y求偏导return df_dx, df_dy
# 定义曲面的曲率函数即二阶偏导数
def curvature(f, x, y):# 使用中心差分法近似求导h 1e-6 # 差分步长可以根据精度要求调整d2f_dx2 (f(x h, y) - 2 * f(x, y) f(x - h, y)) / (h ** 2) # 对x求二阶偏导d2f_dy2 (f(x, y h) - 2 * f(x, y) f(x, y - h)) / (h ** 2) # 对y求二阶偏导d2f_dxdy (f(x h, y h) - f(x h, y - h) - f(x - h, y h) f(x - h, y - h)) / (4 * h ** 2) # 对xy求混合偏导# 根据公式计算高斯曲率K和平均曲率Hdf_dx, df_dy gradient(f, x, y) # 调用梯度函数求一阶偏导E 1 df_dx ** 2F df_dx * df_dyG 1 df_dy ** 2L d2f_dx2 / np.sqrt(1 df_dx ** 2 df_dy ** 2) #np.sqrt()表示开方M d2f_dxdy / np.sqrt(1 df_dx ** 2 df_dy ** 2)N d2f_dy2 / np.sqrt(1 df_dx ** 2 df_dy ** 2)K (L * N - M ** 2) / (E * G - F ** 2) # 高斯曲率H (E * N G * L - 2 * F * M) / (2 * (E * G - F ** 2)) # 平均曲率return K, Hcurvatures []
for i in range(num_points):k, idx, _ kdtree.search_radius_vector_3d(cloud.points[i], radius) #返回邻域点的个数和索引neighbors points[idx] #数组形式为[[x1,y1,z1],[x2,y2,z2],...]#print(k)Y neighbors[:, 2] # 因变量X neighbors[:, [0,1]] # 自变量 [[2,3],[1,1],[8,9],[11,12],[4,5],[8,9]] 6*2popt, pcov curve_fit(func, xdataX.T,ydata Y)x cloud.points[i][0] # 某一点的x坐标y cloud.points[i][1] # 某一点的y坐标K, H curvature(f, x, y) # 计算该点的曲率curvatures.append([K,H])print(curvatures)
for i in range(len(curvatures)):if -0.05curvatures[i][0] 0.05 and -0.05curvatures[i][1] 0.05: #平坦np.asarray(pcd.colors)[i] [0, 0, 0]#黑elif -0.05curvatures[i][0] 0.05 and curvatures[i][1] 0.05: #凸np.asarray(pcd.colors)[i] [1, 0, 0]#红elif -0.05curvatures[i][0] 0.05 and -0.05curvatures[i][1] 0.05: #凹np.asarray(pcd.colors)[i] [0, 1, 0]#绿elif curvatures[i][0] -0.05 and curvatures[i][1] 0.05: #鞍形脊 大部分凸少部分凹np.asarray(pcd.colors)[i] [0, 0, 1]#蓝elif curvatures[i][0] -0.05 and curvatures[i][1] -0.05: #鞍形谷 大部分凹少部分凸np.asarray(pcd.colors)[i] [0, 1, 1]#青elif curvatures[i][0] 0.05 and curvatures[i][1] 0.05: #峰 np.asarray(pcd.colors)[i] [1, 0, 1]#紫elif curvatures[i][0] 0.05 and curvatures[i][1] -0.05: #阱np.asarray(pcd.colors)[i] [1, 1, 0]#黄#显示点云
o3d.visualization.draw_geometries([pcd])
结果 结果出奇的好每个点都进行了划分比想法1好太多了。这里暂时用兔子点云测试到时候创造一些平面点云再来测试一下。 总结
学习东西都不是一蹴而就的果然还是得一步一步脚踏实地地学才学的明白。chatgpt是个好东西只有你也会点东西时它才会回答的正确不能轻信之。
未完成 ps:1、兔子点云pcd读取 2、创建平面点云 3、返回面积、深度信息
