我就废话不多说了,大家还是直接看代码吧~
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 | x0 = linspace( 0.1 , 2 , 100 ); % x0,y0验证函数离散点,可以非等间隔 y0 = 1. / x0; h1 = abs(diff(x0)) ; h = [h1 h1(end)]; ht = h; yapp1 = gradient(y0). / ht; % matlab数值近似 yapp2 = del2(y0). / ht; % matlab数值近似 k2 = abs(yapp2). / ( 1 + yapp1.^ 2 ).^( 3 / 2 ); figure plot(k2) title( '曲率曲线' ) [~,maxFlag] = max(k2); % 曲率最大位置 x_max = x0(maxFlag); y_max = y0(maxFlag); % 画出图像 标注曲率最大点 figure plot(x0,y0, '.-' ); hold on; plot(x_max,y_max, 'rp' ) title( '标注最大曲率点' ) xlabel( 'log10((norm(B*Xk-L)))' ) ylabel( 'log10((norm(Xk)))' ) |
补充:MATLAB 插值+计算离散点曲率
思路:点足够密的话直接用 diff、gradient 求曲率,稀疏的话先插值再算曲率。
公式:
点密的情况 输入曲线坐标(1-2)求一、二阶导数(4-9)通过公式求得曲率(10)
1 2 3 4 5 6 7 8 9 10 11 | x = 0 : 0.01 : 7 ; y = cos(x * 0.5 * pi); h1 = abs(diff(x)); h = [h1 h1(end)]; ht = h; y1 = gradient(y). / ht; y2 = gradient(y1). / ht; curv = abs(y2). / sqrt(( 1 + y1.^ 2 ).^ 3 ); plot(x,y, '-' ,x,curv,' - - r); legend( 'Raw Data, ' Curvature ',' Location', "best" ); grid on |
图像与下文理论值图像相同
点稀疏的情况
1、输入散点坐标(1-2)
2、用样条曲线(B-Spline)等方法插值得到拟合曲线(3-4)
3、diff、gradient 函数求拟合曲线的一、二阶导数(6-11)
4、通过公式求得曲率(12)
例:余弦函数取 8 个点,用 B-Spline 插值
1 2 3 4 5 6 7 8 9 10 11 12 13 | x = 0 : 1 : 7 ; y = cos(x * 0.5 * pi); xx = 0 : 0.01 : 7 ; yy = spline(x,y,xx); h1 = abs(diff(xx)); h = [h1 h1(end)]; ht = h; yy1 = gradient(yy). / ht; yy2 = gradient(yy1). / ht; curv = abs(yy2). / sqrt(( 1 + yy1.^ 2 ).^ 3 ); plot(xx,yy, '-' ,xx,curv, '--r' ,x,y, 'o-' ); legend( 'B-Spline' , 'Curvature' , 'Raw Data' , 'Location' , "best" ); grid on |
补充用法
求最大曲率并在图中标出
1 2 3 | [max_val,max_ind] = max(curv); hold on plot(xx(max_ind),yy(max_ind), '*r' ); |
与理论值(余弦函数曲线)对比
曲率对比
几种插值方法对比
列举四种方法,分别为:分段线性插值、三次样条曲线(B-Spline)插值、三次 Hermite 插值(PCHIP)、修正 Akima 分段三次 Hermite 插值(Akima)
Case 1: 三维螺线

三维螺线散点

插值

俯视

侧视
Case 2:二维梯形波

二维梯形波
Case 3:三维不规则折线

三维不规则折线(不等间距)
对比可得:
Case 1:B-Spline>Akima>PCHIP>Linear
Case 2:Linear>PCHIP>Akima>B-Spline
Case 3:Linear≈PCHIP≈Akima>B-Spline
故在插值的时候需要选择适合的计算方法
以上为个人经验,希望能给大家一个参考,也希望大家多多支持自学编程网。如有错误或未考虑完全的地方,望不吝赐教。
- 本文固定链接: https://zxbcw.cn/post/209751/
- 转载请注明:必须在正文中标注并保留原文链接
- QQ群: PHP高手阵营官方总群(344148542)
- QQ群: Yii2.0开发(304864863)