插值

离散数据逼近:插 值 插 值 问 题 实 例 1 x 0 1 2 … ┇ ┇ ┇ ┇ ┇ 1.0 0.8413 0.8438 0.8461 … 1.1 0.8643 0.8665 0.8686 … ┇ ┇ ┇ ┇ ┇ 正态分布函数 Φ(x) 求Φ(1.014) Φ(1.014)=0.8438+ (0.8461−0.8438)×0.4=0.8447 插 值 问 题 实 例 2 x y 机翼下 轮廓线 插 值 问 题 的 提 法 已知 n+1个节点 ,,1,0(),( njyx jj L= 其中 jx 互不相同，不妨设 ),10 bxxxa n =<<<= L 求任一点 )(* jxx ≠ 处的值 .*y • • • •• 0x 1x nx 0y 1y 节点可视为由 )(xgy = 产生 g 表达式复杂, 或无封闭形式, 或未知。 ‹ *x *y • 在离散数据的基础上补插连续函数，使 得这条连续曲线通过全部给定的离散数 据点 • 插值是离散函数逼近的重要方法，利用 它可通过函数在有限个点处的取值状 况，估算出函数在其他点处的近似值 • 插值：用来填充图像变换时像素之间的 空隙 • • • •• 0x 1x nx 0y 1y 求解插值问题的基本思路 构造一个(相对简单的)函数 ),(xfy = 通过全部节点, 即 ),1,0()( njyxf jj L== 再用 )(xf 计算插值，即 ).( ** xfy = ‹ *x *y 拉格朗日插值 分段线性插值 三次样条插值 一 维 插 值 一维插值的方法： 1.拉格朗日(Lagrange)多项式插值 1.1 插值多项式 )1()( 01 1 1 axaxaxaxL n n n nn ++++= −− L ⎥⎥ ⎥ ⎦ ⎤ ⎢⎢ ⎢ ⎣ ⎡ = ⎥⎥ ⎥ ⎦ ⎤ ⎢⎢ ⎢ ⎣ ⎡ = ⎥⎥ ⎥ ⎦ ⎤ ⎢⎢ ⎢ ⎣ ⎡ = − − n n n n n n nn y y Y a a A xx xx X MM L L L 0 0 1 1 00 ,, 1 1 )2(,0)det( ∴≠XQ 有唯一解 ),1,0()( njyxL jjn L== )2(YXA =⇒ 1.拉格朗日(Lagrange)多项式插值 1.2 拉格朗日插值多项式 ni xxxxxxxx xxxxxxxxxl niiiiii nii i ,1,0,)())(()( )())(()()( 110 110 LLL LL =−−−− −−−−= +− +− )3()()( 0 xlyxL i n i in ∑ = = jjnji yxLji ji xl =∴⎩⎨ ⎧ ≠ == )( ,0 ,1 )(Q 拉格朗日插值基函数。 又（2）有唯一解，故（3）与（1）相同。 特别地: 两点一次(线性)插值多项式: ( ) 1 01 0 0 10 1 1 yxx xxy xx xxxP − −+− −= 三点二次(抛物)插值多项式: ( ) ( ) ( )( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( ) 21202 10 1 2101 20 0 2010 21 2 yxxxx xxxxy xxxx xxxxy xxxx xxxxxP −⋅− −⋅−+−⋅− −⋅−+−⋅− −⋅−= ( )满足插值条件直接验证可知 xPn 1.拉格朗日(Lagrange)多项式插值 1.3 误差估计 ),(),( )!1( )()()()( 0 )1( baxx n gxLxgxR n j j n nn ∈−+=−= ∏= + ξξ 1 )1( )( + + ≤ nn Mg ξ ∏ = + −+≤⇒ n j j n n xxn MxR 0 1 )!1( )( 1.拉格朗日(Lagrange)多项式插值 1.4 例 将[0,π/2] n等分，用g(x)=cos(x)产生n+1个 节点，作Ln(x)（取n=1,2) ,计算cos(π/6), 估计误差。 解: n=1, (x0,y0)=(0,1), (x1,y1)=(π/2,0), L1(x)=y0l0+y1l1=1-2x/π, cos(π/6)=0.666667（近似） n=2, (x0,y0)=(0,1), (x1,y1)=(π/4,0.7071),(x2,y2)=(π/2,0), L2(x)=y0l0+y1l1+y2l2=8(x-π/4)(x-π/2)/π2-16x(x-π/2)0.7071/π2 cos(π/6)=L2(π/6)=0.8508； 精确值：cos (π/6)=0.866025404 nhhhhxxxxxnhMxR n j jjjnn L324,,2,1:)( 2 0 11 ×<−<<== ∏ = ++ π 1 12 )2)(1(4 32 4)!1( 1)( + + +=+< n n n nn nhhhh n xR πL n 1 2 3 4 )( xR n 0.3 0.04 4.7×10-3 4.7×10-4 1.拉格朗日(Lagrange)多项式插值 1.5 拉格朗日插值多项式的振荡 ?)(?)( ↓⇒⇒↑ xRxLn nn 55, 1 1)( 2 ≤≤−+= xxxg 63.363.3),()(lim ≤≤−=∞→ xxgxLnn Runge现象: -5 0 5 -1 . 5 -1 -0 . 5 0 0 . 5 1 1 . 5 2 y = 1 / (1 + x 2 ) n = 2n = 4 n = 6 n = 8 n = 1 0 2.分段线性插值 • •• •• • xjxj-1 xj+1x0 xn n越大，误差越小. nnn xxxxgxI ≤≤=∞→ 0),()(lim ⎪⎪ ⎪⎪ ⎩ ⎪⎪ ⎪⎪ ⎨ ⎧ ≤≤− − ≤≤− − = = + + + − − − = ∑ 其它,0 , , )( )()( 1 1 1 1 1 1 0 jj jj j jj jj j j n j jjn xxx xx xx xxx xx xx xl xlyxL 66, 1 1)( 2 ≤≤−+= xxxg 用分段线性插值法求插值,并观察插值误差. 思考： 1.在[-6,6]中平均选取5个点作插值 4.在[-6,6]中平均选取41个点作插值 2.在[-6,6]中平均选取11个点作插值 3.在[-6,6]中平均选取21个点作插值 比分段线性插值更光滑。 x y xi-1 xia b 在上，光滑程度的定量描述是：函数(曲 线)的k阶导数存在且连续，则称该曲线具有k阶光 滑性。 光滑性的阶次越高，则越光滑。是否存在较低 次的分段多项式达到较高阶光滑性的方法？三次 样条插值就是一个很好的例子。 3. 三次样条插值 },1],,[),({)( 1 nixxxxsxS iii L=∈= − ],[)()3 ),1,0()()2 ),1()()1 0 2 23 n ii iiiii xxCxS niyxS nidxcxbxaxs ∈ == =+++= L L )1,1()()(),()(),()( 111 −=′′=′′′=′=⇒ +++ nixsxsxsxsxsxs iiiiiiiiiiii 自然边界条件）(0)()()4 0 =′′=′′ nxSxS )(,,,)4)3)2 xSdcba iiii ⇒⇒ )()(lim xgxS n =∞→ 思考： 66, 1 1)( 2 ≤≤−+= xxxg 用三次样条插值选取11个基点计算插值 一维插值方法小结 拉格朗日插值（高次多项式插值）：曲线光滑； 误差估计有表达式；收敛性不能保证（振荡现 象）。常用于理论分析。 1.对于n+1个节点，拉格朗日插值为什么用n次多项式， 次数大于n或小于n会如何？另外，得到的Ln(x)次数 会不会小于n？ 2.若产生n+1个节点的g(x)为m次多项式，Ln(x)与g(x) 的关系如何（mn)？ 分段线性和三次样条插值（低次多项式插值）：曲线不光滑 （三次样条插值已大有改进）；误差估计较难（对三次样条 插值）；收敛性有保证。简单实用，应用广泛。 思考 3.样条插值为什么用3次多项式，而不是2或4次？ 4. 三次样条插值中自然边界条件的几何意义是什么？ 用MATLAB作插值计算 一维插值函数： yi=interp1(x，y，xi，'method') 插值方法被插值点插值节点xi处的插 值结果 ‘linear’ ： 线性插值； ‘spline’ ： 三次样条插值； ‘cubic’ ： 立方插值。 缺省时： 分段线性插值。 注意：所有的插值方法都要求x是单调的，并且xi不 能够超过x的范围。 例：在1-12的11小时内，每隔1小时测量一次温 度，测得的温度依次为：5，8，9，15，25，29， 31，30，22，25，27，24。试估计每隔1/10小时的 温度值。 hours=1:12; temps=[5 8 9 15 25 29 31 30 22 25 27 24]; h=1:0.1:12; t=interp1(hours,temps,h,'spline'); plot(hours,temps,'+',h,t,hours,temps,'r:') xlabel('Hour'),ylabel('Degrees Celsius’) x y 思考：用MATLAB作插值计算 机翼下 轮廓线 X 0 3 5 7 9 11 12 13 14 15 Y 0 1.2 1.7 2.0 2.1 2.0 1.8 1.2 1.0 1.6 已知下轮廓线上数据如下，求x每改变0.1时的y值。 二维插值 • • • • •• • ••• • • • • • x y O 第一种（网格节点）： 已知 m×n个节点 ),2,1;,...,2,1(),,( njmizyx ijji L== 其中 ji yx , 互不相同，不妨设 bxxxa m =<<<= L21 dyyyc n =<<<= L21 构造一个二元函数 ),,( yxfz = 通过全部已知节点,即 ),1,0;,,1,0( ),( njmi zyxf ijji LL == = 再用 ),( yxf 计算插值，即 ).,( *** yxfz = 第二种（散乱节点）： • • • • • • • ••• • • • • • y x0 已知n个节点 ),...,2,1(),,( nizyx iii = 其中 ),( ii yx 互不相同， 构造一个二元函数 ),,( yxfz = 通过全部已知节点,即 ),,1,0(),( nizyxf iii L== 再用 ),( yxf 计算插值，即 ).,( *** yxfz = 网格节点插值法： 最邻近插值 分片线性插值 双线性插值 注意：最邻近插值一般不连续。具有连续性的最简单 的插值是分片线性插值。 最邻近插值 x • • • • •• • ••• • • • • • y (x1, y1) (x1, y2) (x2, y1) (x2, y2) O 二维或高维情形的最邻近插值，与被插值点最邻近的 节点的函数值即为所求。 将四个插值点（矩形的四个顶点）处的函数值依次 简记为： 分片线性插值 x y • • • • •• • ••• • • • • • (xi, yj) (xi, yj+1) (xi+1, yj) (xi+1, yj+1) O f (xi, yj)=f1，f (xi+1, yj)=f2，f (xi+1, yj+1)=f3，f (xi, yj+1)=f4 插值函数为： ji i1i j1j y)xx( xx yy y +−− −≤ + + )yy)(ff()xx)(ff(f)y,x(f j23i121 −−+−−+= 第二片(上三角形区域)：(x, y)满足 ii i1i j1j y)xx( xx yy y +−− −> + + 插值函数为： )xx)(ff()yy)(ff(f)y,x(f i43j141 −−+−−+= 注意：(x, y)当然应该是在插值节点所形成的矩形区 域内。显然，分片线性插值函数是连续的 分两片的函数表达式如下： 第一片(下三角形区域)： (x, y)满足 双线性插值是一片一片的空间二次曲面构成。 双线性插值函数的形式如下： )dcy)(bax()y,x(f ++= 其中有四个待定系数，利用该函数在矩形的四个顶 点（插值节点）的函数值，得到四个代数方程，正 好确定四个系数。 双线性插值 x • • • • •• • ••• • • • • • y (x1, y1) (x1, y2) (x2, y1) (x2, y2) O 要求x0,y0单调；x，y可取为矩阵，或x取 行向量，y取为列向量，x,y的值分别不能超出 x0,y0的范围。 z=interp2(x0,y0,z0,x,y,’method’) 被插值点 插值方法 用MATLAB作网格节点数据的插值 插值 节点 被插值点 的函数值 ‘nearest’ 最邻近插值 ‘linear’ 双线性插值 ‘cubic’ 双三次插值 缺省时, 双线性插值 例：测得平板表面3*5网格点处的温度分别为： 82 81 80 82 84 79 63 61 65 81 84 84 82 85 86 试作出平板表面的温度分布曲面z=f(x,y)的图形。 输入以下命令： x=1:5; y=1:3; temps=[82 81 80 82 84;79 63 61 65 81;84 84 82 85 86]; mesh(x,y,temps) 输入以下命令: xi=1:0.2:3; yi=1:0.2:3; zi=interp2(x,y,temps,xi,yi,'cubic'); 2．以平滑数据,在x、y方向上进行插值 插值函数griddata格式为: cz =griddata（x，y，z，cx，cy，‘method’） 用MATLAB作散点数据的插值计算 要求cx取行向量，cy取为列向量。 被插值点 插值方法 插值 节点 被插值点 的函数值 ‘nearest’ 最邻近插值 ‘linear’ 双线性插值 ‘cubic’ 双三次插值 'v4'- Matlab提供的插值方法 缺省时, 双线性插值 作业: 在某海域测得一些点(x,y)处的水深z由下表 给出，船的吃水深度为5英尺，在矩形区域（75， 200）*（-50，150）里的哪些地方船要避免进入。 x y z 129 140 103.5 88 185.5 195 105 7.5 141.5 23 147 22.5 137.5 85.5 4 8 6 8 6 8 8 x y z 157.5 107.5 77 81 162 162 117.5 -6.5 -81 3 56.5 -66.5 84 -33.5 9 9 8 8 9 4 9 作业 山区地貌： 在某山区测得一些地点的高度如下表。平面区域为 1200<=x<=4000,1200<=y<=3600) 试作出该山区的地貌图和等高线图，并对几种插值方法进行比较。 Y X 1200 1600 2000 2400 2800 3200 3600 4000 1200 1130 1250 1280 1230 1040 900 500 700 1600 1320 1450 1420 1400 1300 700 900 850 2000 1390 1500 1500 1400 900 1100 1060 950 2400 1500 1200 1100 1350 1450 1200 1150 1010 2800 1500 1200 1100 1550 1600 1550 1380 1070 3200 1500 1550 1600 1550 1600 1600 1600 1550 3600 1480 1500 1550 1510 1430 1300 1200 980 通过此例对最近邻点插值、双线性插值方法和双三次插值方法的插 值效果进行比较。