用matlab实现双线性插值旋转图像程序
2008年11月06日星期四 19:13
用matlab实现双线性插值旋转图像程序
2008-04-27 1
2:51
双线性插值：输出像素值是它在输入图像中2*2邻域采样点的平均值，它根据某像素周围4个像素的灰度值在水平和垂直两个方法上对其插值。为了方便理解，先考虑一维情况下的线性插值：对于一个数列c，假设c[a]到c[a+1]之间是线性变化的，那么对于浮点数x(a<=x c(x)=c[a+1]*(x-a)+c[a]*(1+a-x); 把这种插值方式扩展到二维情况：对于一个二维数组c，我们假设对于任意一个浮点数i,c(a,i)到c(a+1,i)之间是线性变化的,c(i,b)到c(i,b+1)之间也是线性变化的(a,b都是整数) ，那么对于浮点数的坐标(x,y)满足(a<=x c(x,b+1): c(x,b) = c[a+1][b]*(x-a)+c[a][b]*(1+a-x); c(x,b+1) = c[a+1][b+1]*(x-a)+c[a][b+1]*(1+a-x); 好，现在已经知道c(x,b)和c(x,b+1)了，而根据假设c(x,b)到c(x,b+1)也是线性变化的，所以: c(x,y) = c(x,b+1)*(y-b)+c(x,b)*(1+b-y) 。 对应以上原理，双线性插值的matlab程序如下： function [imN] = im_R(im,rotation) %im是输入图像，rotation是旋转矩阵close all; im = imread(im); subplot(1,2,1); imshow(im); dim = ndims(im); %确定维数，对灰度和彩色图像分别处理 w = size(im,2); h = size(im,1); if( dim == 3 ) imN = zeros( h, w, 3 ); else imN = zeros(h,w); end im = double(im); for r = 1 : h for c = 1 : w t = rotation * [ r - h / 2; c - w / 2 ] + [ h / 2; w / 2 ];%旋转 if( t(1) > 1 && t(2) > 1 && t(1) < h && t(2) < w )%双线性插值 x = floor( t(1) ); y = floor( t(2) ); v1 = im( x + 1, y, : ) * ( t(1) - x ) + im( x, y, : ) * ( 1 + x - t(1) ); v2 = im( x + 1, y + 1, : ) * ( t(1) - x ) + im( x, y + 1, : ) * ( 1 + x - t(1) ); v = v2 * ( t(2) - y ) + v1 * ( 1 + y - t(2) ); imN( r, c, : ) = round(v); end end end imN = uint8(imN); subplot(1,2,2); imshow(imN); 出处：https://www.doczj.com/doc/a07040241.html,/kr1423/blog/item/89921324ce87862cd507424e.html 三次样条插值法 h0=figure('toolbar','none',... 'position',[200 50 350 450],... 'name','实例86'); h1=axes('parent',h0,... 'position',[0.10 0.45 0.8 0.5],... 'visible','off'); x=0:0.2:2*pi; y=sin(x); plot(x,y) b1=uicontrol('parent',h0,... 'units','points',... 'tag','b1',... 'style','pushbutton',... 'string','三次样条插值',... 'backgroundcolor',[0.75 0.75 0.75],... 'position',[20 60 70 20],... 'callback',[... 'y=0,',... 'sy=0,',... 'strn1=get(e2,''string'');,',... 'n1=str2num(strn1);,',... 'strn2=get(e3,''string'');,',... 'n2=str2num(strn2);,',... 'x=n1:0.2:n2;,',... 'for t=n1:0.2:n2,',... 'y(i)=sin(t);,',... 'sy(i)=san(t,n1,n2);,',... 'i=i+1;,',... 'end,',... 'plot(x,y,''b*'',x,sy,''r-''),',... 'axis([0 7 -1.5 1.5]),',... 'legend(''sin(x)'',''N-Hermite插值'')']); b2=uicontrol('parent',h0,... 'units','points',... 'tag','b2',... 'style','pushbutton',... 'string','误差比较',... 'backgroundcolor',[0.75 0.75 0.75],... 'position',[170 60 70 20],... 'callback',[... 'strdn1=get(e2,''string'');,',... 'n1=str2num(strdn1);,',... 'strdn2=get(e3,''string'');,',... 'n2=str2num(strdn2);,',... 'strdn=get(e1,''string'');,',... 'dn=str2num(strdn);,',... 'dd=abs(sin(dn)-san(dn,n1,n2));,',... 'msgbox([''误差为：'',num2str(dd)],''计算结果'')']); e1=uicontrol('parent',h0,... 'units','points',... 'tag','e1',... 'style','edit',... 'fontsize',12,... 'string','1.20',... 'horizontalalignment','right',... 'backgroundcolor',[1 1 1],... 'position',[200 100 40 20]); t1=uicontrol('parent',h0,... 'units','points',... 'tag','t1',... 'style','text',... 'string','误差点：',... 'fontsize',12,... 'backgroundcolor',[0.75 0.75 0.75],... 'position',[160 100 40 20]); e2=uicontrol('parent',h0,... 'units','points',... 'style','edit',... 'fontsize',12,... 'string','1.00',... 'horizontalalignment','right',... 'backgroundcolor',[1 1 1],... 'position',[20 85 40 20]); t2=uicontrol('parent',h0,... 'units','points',... 'tag','t2',... 'style','text',... 'string','第一节点：',... 'fontsize',12,... 'backgroundcolor',[0.75 0.75 0.75],... 'position',[15 105 50 20]); e3=uicontrol('parent',h0,... 'units','points',... 'tag','e3',... 'style','edit',... 'fontsize',12,... 'string','3.00',... 'horizontalalignment','right',... 'backgroundcolor',[1 1 1],... 'position',[100 85 40 20]); t3=uicontrol('parent',h0,... 'units','points',... 'tag','t3',... 'style','text',... 'string','第二节点：',... 'fontsize',12,... 'backgroundcolor',[0.75 0.75 0.75],... 'position',[95 105 50 20]); b3=uicontrol('parent',h0,... 'units','points',... 'tag','b3',... 'style','pushbutton',... 'string','关闭',... 'backgroundcolor',[0.75 0.75 0.75],... 'position',[100 20 60 20],... Several Typical Interpolation in Matlab Lagrange Interpolation Supposing： If x=175, while y=? Solution： Lagrange Interpolation in Matlab: function y=lagrange(x0,y0,x); n=length(x0);m=length(x); for i=1:m z=x(i); s=0.0; for k=1:n p=1.0; for j=1:n if j~=k p=p*(z-x0(j))/(x0(k)-x0(j)); end end s=p*y0(k)+s; end y(i)=s; end input： x0=[144 169 225] y0=[12 13 15] y=lagrange(x0,y0,175) obtain the answer： x0 = 144 169 225 y0 = 12 13 15 y = 13.2302 Spline Interpolation Solution ： Input x=[ 1 4 9 6]；y=[ 1 4 9 6]；x=[ 1 4 9 6]；pp=spline(x,y) pp = form: 'pp' breaks: [1 4 6 9] coefs: [3x4 double] pieces: 3 order: 4 dim: 1 output : pp.coefs ans = -0.0500 0.5333 -0.8167 1.0000 -0.0500 0.0833 1.0333 2.0000 -0.0500 -0.2167 0.7667 4.0000 It shows the coefficients of cubic spline polynomial , so: S （x ）=, 169,3)9(1484.0)9(0063.0)9(0008.0,94,2)4(2714.0)4(0183.0)4(0008 .0, 41,1)1(4024.0)1(0254.0)1(0008.0232 3 23≥≤+-+---≥≤+-+---≥≤+-+---x x x x x x x x x x x x Newton’s Interpolation Resolve 65 Solution: Newton’s Interpolation in matlab ： function yi=newint(x,y,xi); n=length(x); ny=length(y); if n~=ny error end Y=zeros(n);Y(:,1)=y'; Lagrange插值： x=0:3; y=[-5,-6,-1,16]; n=length(x); syms q; for k=1:n fenmu=1; p=1; for j=1:n if(j~=k) fenmu=fenmu*(x(k)-x(j)) p=conv(p,poly(x(j))) end end c(k,:)=p*y(k)/fenmu end a=zeros(1,n); for i=1:n for j=1:n a(i)=a(i)+c(j,i) end end 输出结果： fenmu = -1 p = 1 -1 fenmu = 2 p = 1 -3 2 fenmu = -6 p = 1 -6 11 -6 c = 0.8333 -5.0000 9.1667 -5.0000 fenmu = 1 p = 1 0 fenmu = -1 p = 1 - 2 0 fenmu = 2 p = 1 -5 6 0 c = 0.8333 -5.0000 9.1667 -5.0000 -3.0000 15.0000 -18.0000 0 fenmu = 2 p = 1 0 fenmu = 2 p = 1 -1 0 fenmu = -2 p = 1 -4 3 0 c = 0.8333 -5.0000 9.1667 -5.0000 -3.0000 15.0000 -18.0000 0 0.5000 -2.0000 1.5000 0 fenmu = 3 p = 1 0 fenmu = 6 p = 1 -1 0 fenmu = 6 p = 1 -3 2 0 c = 0.8333 -5.0000 9.1667 -5.0000 -3.0000 15.0000 -18.0000 0 0.5000 -2.0000 1.5000 0 2.6667 -8.0000 5.3333 0 a = 实验五图像的空间变换 一、实验目的 1、学习图像空间变换，并通过实验体会空间变换的效果，对其作出分析。 2、掌握利用最邻近插值和双线性插值算法(灰度插值)实现图像的缩放。 3、掌握MATLAB编程环境中基本的图像处理函数。 二、实验要求 1.读入图像，对其利用最邻近插值和双线性插值法进行缩放变换，要求先使用IPT函数进行变换，然后自己编写函数实现； 2.对比上述得到的结果。 三、实验原理 图像的空间变换，也称几何变换或几何运算，包括图像的平移、旋转、镜像变换、转置、缩放等。几何运算可改变图像中各物体之间的空间关系，这种运算可以看成是将各物体在图像内移动。 空间变换可如下表示：设(u，v)为源图像上的点，(x，y)为目标图像上的点，则空间变换就是将源图像上(u，v)处的像素值与目标图像上(x，y)处的像素值对应起来,并具有以下关系： x=X(u，v)，y=Y(u，v) （即由(u，v)计算对应(x，y)）(1.1) 或u=U(x，y)，v=V(x，y) （即由(x，y)计算对应(u，v)）(1.2) 其中X(u，v)、Y(u，v)、U(x，y)、V(x，y)均为变换。由(1.1)对应的变换称作向前映射法也叫像素移交法，而由(1.2)对应的变换称作向后映射法也叫像素填充法，向后映射法是向前映射法的逆。 最简单的插值算法是最邻近插值，也称为零阶插值。最邻近插值算法简单，在许多情况 下都能得到令人满意的结果，但是当图像中包含像素之间灰度级有变化的细微结构时，最邻近算法会在图像中产生人为加工的痕迹。双线性插值算法计算量比零阶插值大，但缩放后图像质量高，不会出现像素值不连续的的情况，这样就可以获得一个令人满意的结果。最邻近点插值取插值点的4个邻点中距离最近的邻点灰度值作为该点的灰度值。设插值点（i，j）到周边4个邻点fk（i，j）（k ＝1，2，3，4）的距离为dk（k ＝1，2，3，4），则：g（i，j）＝fk（i，j），dl ＝min{d1，d2，d3，d4}，l＝1，2，3，4 。 双线性插值是利用了需要处理的原始图像像素点周围的四个像素点的相关性，通过双线插值算法计算得出的。对于一个目的坐标，通过后映射法得到其在原始图像的对应的浮点坐标(i+u，j+v)，其中i，j均为非负整数，u，v为[0，l]区间的浮点数，则这个像素的值f(i+u，j+v)可由原图像中坐标为(i，j)、(i+l，j)、(i，j+1)、(i+1，j+1)所对应的周围四个像素的值决定，即：f(i+u，j+v)=(1-u)×(1-v)×f(i，j)+(1-u)×v×f(i，j+1)+u×(1-v)×f(i+l，j)+u×v×f(i+l，j+1)，其中f(i，j)表示源图像(i，j)处的的像素值，以此类推，这就是双线性内插值法。 如下图所示，已知(0，0)、(0，1)、(1，0)、(1，1)四点的的灰度，可以由相邻像素的灰度值f(0，0)和f(1，0)在X方向上线性插值求出(x，0)的灰度f(x，0)，由另外两个相邻像素f(0，1)和f(1，1)在X方向上线性插值可求出(x，1)的灰度f(x，1)，最后由f(x，0)，f(x，1)在Y 方向上进行线性插值就可以得到(x，y)的灰度f(x，y)。 四、实验代码 Matlab中插值函数汇总和使用说明 MATLAB中的插值函数为interp1，其调用格式为： yi= interp1(x,y,xi,'method') 其中x，y为插值点，yi为在被插值点xi处的插值结果；x,y为向量，'method'表示采用的插值方法，MA TLAB提供的插值方法有几种：'method'是最邻近插值，'linear'线性插值；'spline'三次样条插值；'cubic'立方插值．缺省时表示线性插值 注意：所有的插值方法都要求x是单调的，并且xi不能够超过x的范围。 例如：在一天24小时内，从零点开始每间隔2小时测得的环境温度数据分别为12，9，9，10，18 ，24，28，27，25，20，18，15，13， 推测中午12点（即13点）时的温度． x=0:2:24; y=[12 9 9 10 18 24 28 27 25 20 18 15 13]; a=13; y1=interp1(x,y,a,'spline') 结果为：27.8725 若要得到一天24小时的温度曲线，则： xi=0:1/3600:24; yi=interp1(x,y,xi, 'spline'); plot(x,y,'o' ,xi,yi) 命令1 interp1 功能一维数据插值（表格查找）。该命令对数据点之间计算内插值。它找出一元函数f(x)在中间点的数值。其中函数f(x)由所给数据决定。 x：原始数据点 Y：原始数据点 xi：插值点 Yi：插值点 格式 (1)yi = interp1(x,Y,xi) 返回插值向量yi，每一元素对应于参量xi，同时由向量x 与Y 的内插值决定。参量x 指定数据Y 的点。 若Y 为一矩阵，则按Y 的每列计算。yi 是阶数为length(xi)*size(Y,2)的输出矩阵。 (2)yi = interp1(Y,xi) 假定x=1:N，其中N 为向量Y 的长度，或者为矩阵Y 的行数。 (3)yi = interp1(x,Y,xi,method) 用指定的算法计算插值： ’nearest’：最近邻点插值，直接完成计算； ’linear’：线性插值（缺省方式），直接完成计算； ’spline’：三次样条函数插值。对于该方法，命令interp1 调用函数spline、ppval、mkpp、umkpp。这些命令生成一系列用于分段多项式操作的函数。命令spline 用它们执行三次样条函数插值； ’pchip’：分段三次Hermite 插值。对于该方法，命令interp1 调用函数pchip，用于对向量x 与y 执行分段三次内插值。该方法保留单调性与数据的外形； ’cubic’：与’pchip’操作相同； ’v5cubic’：在MATLAB 5.0 中的三次插值。 对于超出x 范围的xi 的分量，使用方法’nearest’、’linear’、’v5cubic’的插值算法，相应地将返回NaN。对其他的方法，interp1 将对超出的分量执行外插值算法。 (4)yi = interp1(x,Y,xi,method,'extrap') 对于超出x 范围的xi 中的分量将执行特殊的外插值法extrap。 (5)yi = interp1(x,Y,xi,method,extrapval) 确定超出x 范围的xi 中的分量的外插值extrapval，其值通常取NaN 或0。 例1 1. 2.>>x = 0:10; y = x.*sin(x); 3.>>xx = 0:.25:10; yy = interp1(x,y,xx); 4.>>plot(x,y,'kd',xx,yy) 复制代码 例2 1. 2.>> year = 1900:10:2010; 3.>> product = [75.995 91.972 105.711 123.203 131.669 150.697 179.323 203.212 226.505 4.249.633 256.344 267.893 ]; 5.>>p1995 = interp1(year,product,1995) 6.>>x = 1900:1:2010; 7.>>y = interp1(year,product,x,'pchip'); 8.>>plot(year,product,'o',x,y) 复制代码 插值结果为： 1. 自己写的Matlab旋转+双线性插值图像函数效果图： 源码： clear all; I = imread('original.jpg'); [Height,Width,RGB] = size(I); II = I;%当角度为0时直接输出 %本程序是以左上角为坐标原点 %angle_j是旋转角度，正值是按顺时针旋转，负值时按逆时针旋转 angle_j = 181; %angle是弧度 angle = 2*pi*angle_j/360; %将angle转成正值 while(angle < 0) angle = 2 * pi + angle; end %约束在0-2π内 while(angle > 2 * pi) angle = angle - 2 * pi; end %tag是判断下面的while循环有没有执行过 tag = 0; while(angle > 0) %超过90度的旋转，都先旋转90度，直到角度在0°-90°之间 %原理是旋转90度整数倍时，信息是不丢失的 if angle >= pi/2 a = pi/2; angle = angle - pi/2; elseif0 < angle < pi/2 a = angle; angle = 0; end if tag == 0 tag = 1; else I = II; [Height,Width,RGB] = size(I);%在旋转后的图像上继续旋转，从而实现大于90° 的旋转 end %正向变换用 sina = sin(a); cosa = cos(a); %逆向变换用_m == _minus sina_m = sin(-a); cosa_m = cos(-a); %旋转后图像的长度和宽度 II_height = round(sina * Width + cosa * Height); II_width = round(sina * Height + cosa * Width); II = ones(II_height,II_width,3); %先转成unit8。或者下面赋值0-1规划一下。否则imshow全是白色。 II = im2uint8(II); %%%%%%%%%%%%%%%%%%%%正向映射%%%%%%%%%%%%%%%%%%%%%%%% 拉格朗日插值的调用函数 function y=lagrange(x0,y0,x) n=length(x0);m=length(x); for i=1:m z=x(i); L=0.0; for j=1:n T=1.0; for k=1:n if k~=j T=T*(z-x0(k))/(x0(j)-x0(k)); end end L=T*y0(j)+L; end y(i)=L; end 四个图在一起： x=[-1:0.05:1]; y=1./(1+25*x.^2); x0=[-1:0.4:1]; y0=1./(1+25*x0.^2); y1=lagrange(x0,y0,x); x0=[-1:0.2:1]; y0=1./(1+25*x0.^2); y2=lagrange(x0,y0,x); x0=[-1:0.1:1]; y0=1./(1+25*x0.^2); y3= lagrange(x0,y0,x); plot(x,y,'-r') hold on plot(x,y1,'-b',x,y2,'-r',x,y3,'-r') l5和fx在一起： x=[-1:0.05:1]; y=1./(1+25*x.^2); x0=[-1:0.4:1]; y0=1./(1+25*x0.^2); y1=lagrange(x0,y0,x); plot(x,y,'-r') hold on plot(x,y1,'-b') l10和fx在一起： x=[-1:0.05:1]; y=1./(1+25*x.^2); x0=[-1:0.2:1]; y0=1./(1+25*x0.^2); y2= lagrange(x0,y0,x); plot(x,y,'-r') hold on plot(x,y2,'-b') l20和fx在一起： x=[-1:0.05:1]; y=1./(1+25*x.^2); x0=[-1:0.1:1]; y0=1./(1+25*x0.^2); y3= lagrange(x0,y0,x); plot(x,y,'-r') hold on plot(x,y3,'-b') 双线性插值法 I_=imread('test.jpg'); I=rgb2gray(I_); A=0.7;B=0.7;%失真像素坐标 [i,j]=size(I); m=round(i*A);n=round(j*B); temp=zeros(m,n);%产生m*n矩阵 G=[A0;0B]; for x=1:m for y=1:n ab=[x,y]/G;%取得x/A,y/B a=ab(1)-floor(ab(1));%权值 b=ab(2)-floor(ab(2)); %防溢出处理 if ab(1)<1 ab(1)=1; end if ab(1)>i ab(1)=i; end if ab(2)<1 ab(2)=1; end if ab(2)>j ab(2)=j; end %定义内插值坐标 ab11=[floor(ab(1))floor(ab(2))]; ab12=[floor(ab(1))ceil(ab(2))]; ab21=[ceil(ab(1))floor(ab(2))]; ab22=[ceil(ab(1))ceil(ab(2))]; temp(x,y)=(1-a)*(1-b)*I(ab11(1),ab11(2))+... a*(1-b)*I(ab12(1),ab12(2))+... (1-a)*b*I(ab21(1),ab21(2))+... a*b*I(ab22(1),ab22(2)); end end imshow(uint8(temp)),title('0.7倍双线性');最近邻法 I_=imread('test.jpg');%读入原始图像 I1=rgb2gray(I_); [i,j]=size(I1); m=round(i*1.5);n=round(j*1.5); m_=round(i*0.7);n_=round(j*0.7); %1.5倍最邻近 TEMP=zeros(m,n);%产生m*n矩阵 for i=1:m for j=1:n TEMP(i,j)=I1(round(i/1.5),round(j/1.5)); end end subplot(1,3,1),imshow(I1),title('原图') TEMP1_5=uint8(TEMP); subplot(1,3,2), imshow(TEMP1_5),title('1.5倍最邻近') 全局预测下的图像分割 I_=imread('test.jpg'); I=rgb2gray(I_); [m,n]=size(I); %统计直方图 zhifangtu=zeros(1,255);% for i=1:1:m for j=1:1:n zhifangtu(I(i,j)+1)= zhifangtu(I(i,j)+1)+1; end end plot(zhifangtu); %阈值处理 final=zeros(m,n); for x=1:1:m for y=1:1:n AA=I(x,y); if AA>120 final(x,y)=255; else final(x,y)=0; end end end imshow(uint8(final)); 数值分析作业 姓名王建忠 学号132080202006 学院能源与动力工程 专业机械电子工程 2013年12月16日 https://www.doczj.com/doc/a07040241.html,grange插值多项式程序 function f=nalagr(x,y,xx) %x为节点向量 %y为节点函数值 %xx是插值点 syms s if(length(x)==length(y)) n=length(x); else disp('x和y的维数不相等!'); return; end f=0.0; for(i=1:n) l=y(i); for(j=1:i-1) l=l*(s-x(j))/(x(i)-x(j)); end; for(j=i+1:n) l=l*(s-x(j))/(x(i)-x(j));%计算拉格朗日基函数end; f=f+l;%计算拉格朗日插值函数 simplify(f); if(i==n) if(nargin==3) f=subs(f,'s');%计算插值点的函数值else f=collect(f);%将插值多项式展开 f=vpa(f,6);%将插值多项式的系数化成6位精度的小数 end end end >>x=[-2,-1,0,1];%已知节点向量y=[3,1,1,6];%节点函数值向量 f=nalagr(x,y) f= 0.5*s^3+ 2.5*s^2+ 2.0*s+ 1.0 >>f=nalagr(x,y,0) f=1 >> 2.牛顿插值多项式程序 function[p2,z]=newTon(x,y,t) %输入参数中x,y为元素个数相等的向量，t为待估计的点，可以为数字或向量。%输出参数中p2为所求得的牛顿插值多项式，z为利用多项式所得的t的函数值。 n=length(x); chaS(1)=y(1); for i=2:n x1=x;y1=y; x1(i+1:n)=[]; y1(i+1:n)=[]; n1=length(x1); s1=0; for j=1:n1 t1=1; for k=1:n1 if k==j continue; else t1=t1*(x1(j)-x1(k)); end end s1=s1+y1(j)/t1; end chaS(i)=s1; end b(1,:)=[zeros(1,n-1)chaS(1)]; cl=cell(1,n-1); for i=2:n u1=1; for j=1:i-1 u1=conv(u1,[1-x(j)]); cl{i-1}=u1; end cl{i-1}=chaS(i)*cl{i-1}; b(i,:)=[zeros(1,n-i),cl{i-1}]; end p2=b(1,:); for j=2:n p2=p2+b(j,:); end if length(t)==1 rm=0; 《工程常用算法》综合实践作业二 完成日期： 2013年 4月 14 日 班级 学号 姓名 主要工作说明 自评成绩 0718 2010071826 崔洪亮 算式与程序的编写 18 0718 2010071815 侯闰上 流程图的编辑，程序的审查 0718 2010071809 赵化川 报告的整理汇总 一.作业题目：三次样条插值与分段插值 已知飞机下轮廓线数据如下： x 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 飞机下轮廓线形状大致如下图所示： 要求分别用拉格朗日插值法、Newton 插值法、分段线性插值法和三次样条插值法计算x 每改变0.5时y 的值，即x 取 0.5, 1, 1.5, … , 14.5 时对应的y 值。比较采用不同方法的计算工作量、计算结果和优缺点。 二.程序流程图及图形 1.拉格朗日插值法 开始 x,y,x0 Length (x)==l Ength (y)? n=length (x) i=1:n,l=1。 j=1:i-1&j=i+1:n l=l.*(x0-x(j)/x(i)-x(j) f=f+l*y(i) 结束 否 是 机翼 下轮廓线 2.牛顿插值法 开始 x,y,xi Length(x)==l ength(y)? n=length(x)Y=zeros （n),Y （：1）=y,f=0 a=1:n-1,b=1:n-a,Y(b,a+1)=（Y （b+1，a)-Y(b,a))/(x （b+a)-x(b)) i=1:n,z=1 结束 j=1:i-1,z=z.*(xi-x(j)) f=f+Y(1,i)*z 否 是 3.分段线性插值法 开始 x ，y ，x0 length （x ）==length(y)? k=1:n-1 x(k)<=x0&x0《=x(k+1) temp=x(k)-x(k+1) f=(x0-x(k+1))/temp*y(k)+(x0-x(k))/(-temp)*y(k+1) 结束 否否 是 是 三.matlab 程序及简要的注释（m 文件） 1.拉格朗日插值法 2.牛顿插值法 function f=newdun(x,y,xi) %x 为已知数据点的x 坐标向量 %y 为已知数据点的y 坐标向量 function f=lang(x,y,x0) %x 为已知数据点的x 坐标向量 %y 为已知数据点的y 坐标向量 双线性插值，又称为双线性内插。在数学上，双线性插值是有两个变量的插值函数的线性插值扩展，其核心思想是在两个方向分别进行一次线性插值。例如已知的红色数据点与待插值得到的绿色点如图1所示： 图1 假如我们想得到未知函数在点的值，假设我们已知函数 在, , , 及四个点的值。首先在 x 方向进行线性插值，得到 然后在 y 方向进行线性插值，得到 这样就得到所要的结果, 如果选择一个坐标系统使得的四个已知点坐标分别为(0, 0)、(0, 1)、(1, 0) 和(1, 1)，那么插值公式就可以化简为 或者用矩阵运算表示为 这就是双线性内插值法。双线性内插值法计算量大，但缩放后图像质量高，不会出现像素值 不连续的的情况。由于双线性插值具有低通滤波器的性质，使高频分量受损，所以可能会使 图像轮廓在一定程度上变得模糊。 双线性插值法的MATLAB源代码为： I=imread('lena.jpg'); %读入原图像 [nrows,ncols,z]=size(I); %读取图像矩阵大小，方便后面操作 K = str2double(inputdlg('please input scale factor (must between 0.2 - 5.0)', 'INPUT scale factor', 1, {'0.5'})); width = K * nrows; height = K * ncols; J = uint8(zeros(width,height,z)); widthScale = nrows/width; heightScale = ncols/height; for x = 5:width - 5 % 5是为了防止矩阵超出边界溢出 for y = 5:height - 5 for z=1:3 xx = x * widthScale; % xx, yy为原坐标，x,y为新坐标 yy = y * heightScale; if((xx/double(uint16(xx))==1.0)&&(yy/double(uint16(yy))==1.0)) J(x,y,z) = I(int16(xx),int16(yy),z); %若xx，yy为整数，直接赋值 else a = double(uint16(xx)); b = double(uint16(yy)); x11 = double(I(a,b,z)); % x11 <- I(a,b) x12 = double(I(a,b+1,z)); % x12 <- I(a,b+1) x21 = double(I(a+1,b,z)); % x21 <- I(a+1,b) x22 = double(I(a+1,b+1,z));% x22 <- I(a+1,b+1) J(x,y,z) = uint8((b+1-yy)*((xx-a)*x21+(a+1-xx)*x11)+(yy-b)* ((xx-a)*x22+(a+1-xx)*x12)); %用双线性插值计算公式计算 end end end end 对图像进行旋转，使用最近邻插值法，双线性插值，三次卷积插值三种方法进行插值。 源码： clc;clear all;close all; Img=imread('test1.bmp'); Img=double(Img); [h w]=size(Img); alpha=pi/6; %逆时针旋转的角度 wnew=w*cos(alpha)+h*sin(alpha); %新图像的宽width hnew=w*sin(alpha)+h*cos(alpha); %新图像的高heighth wnew=ceil(wnew); %取整 hnew=ceil(hnew); u0=w*sin(alpha); %平移量 T=[cos(alpha),sin(alpha);-sin(alpha),cos(alpha)]; %变换矩阵 Imgnew1=zeros(hnew,wnew); Imgnew2=zeros(hnew,wnew); Imgnew3=zeros(hnew,wnew); for u=1:hnew %u和v是新图像坐标，变换到原图像坐标x和y中。 for v=1:wnew tem=T*([u;v]-[u0;0]); x=tem(1); y=tem(2); if x>=1 & x<=h & y>=1 & y<=w %若变换出的x和y在原图像范围内 x_low=floor(x); x_up=ceil(x); y_low=floor(y); y_up=ceil(y); if (x-x_low)<=(x_up-x) %采用最近点法，选取距离最近点的像素赋给新图像x=x_low; else x=x_up; end if (y-y_low)<=(y_up-y) y=y_low; else y=y_up; end p1=Img(x_low,y_low); %双线性插值，p1到p4是（x,y）周围的四个点p2=Img(x_up,y_low); p3=Img(x_low,y_low); p4=Img(x_up,y_up); s=x-x_low; t=y-y_low; Imgnew1(u,v)=Img(x,y); Imgnew2(u,v)=(1-s)*(1-t)*p1+(1-s)*t*p3+(1-t)*s*p2+s*t*p4; end if x>=2 & x<=h-2 & y>=2 & y<=w-2 %若变换出的x和y在原图像范围内x_1=floor(x)-1; x_2=floor(x); x_3=floor(x)+1; x_4=floor(x)+2; y_1=floor(y)-1; y_2=floor(y); y_3=floor(y)+1; y_4=floor(y)+2; A=[sw(1+x-x_2),sw(x-x_2),sw(1-(x-x_2)),sw(2-(x-x_2))]; C=[sw(1+y-y_2),sw(y-y_2),sw(1-(y-y_2)),sw(2-(y-y_2))]; B=[ Img(x_1,y_1),Img(x_1,y_2),Img(x_1,y_3),Img(x_1,y_4); Img(x_2,y_1),Img(x_2,y_2),Img(x_2,y_3),Img(x_2,y_4); Img(x_3,y_1),Img(x_3,y_2),Img(x_3,y_3),Img(x_3,y_4); Img(x_4,y_1),Img(x_4,y_2),Img(x_4,y_3),Img(x_4,y_4)]; Imgnew3(u,v)=A*B*C'; end 用M A T L A B实现拉格朗 日插值和分段线性插值 The Standardization Office was revised on the afternoon of December 13, 2020 用MATLAB实现拉格朗日插值和分段线性插值 1、实验内容： 用MATLAB实现拉格朗日插值和分段线性插值。 2、实验目的： 1）学会使用MATLAB软件； 2）会使用MATLAB软件进行拉格朗日插值算法和分段线性 差值算法； 3、实验原理： 利用拉格朗日插值方法进行多项式插值，并将图形显式出来。 4、实验步骤及运行结果 （1）实现lagrange插值 1）定义函数： f = 1/(x^2+1) 将其保存在文件中，具体程序如 下： function y = f1(x) y = 1./(x.^2+1); 2）定义拉格朗日插值函数：将其保存在文件中，具体实现程序 编程如下： function y = lagrange(x0,y0,x) m = length(x); /区间长度/ n = length(x0); for i = 1:n l(i) = 1; end for i = 1:m for j = 1:n for k = 1:n if j == k continue; end l(j) = ( x(i) -x0(k))/( x0(j) - x0(k) )*l(j); end end end y = 0; for i = 1:n y = y0(i) * l(i) + y; end 3）建立测试程序，保存在文件中，实现画图： x=-5::5; y=(1+x.^2).^-1; p=polyfit(x,y,n); py=vpa(poly2sym(p),10) plot_x=-5::5; f1=polyval(p,plot_x); figure pl ot(x,y,‘r',plot_x,f1) 输入n=6,出现下面的图形： 插值方法 晚上做一个曲线拟合，结果才开始用最小二乘法拟合时，拟合出来的东西太难看了！ 于是尝试用其他方法。 经过一番按图索骥，终于发现做曲线拟合的话，采用插值法是比较理想的方法。尤其是样条插值，插完后线条十分光滑。 方法付后，最关键的问题是求解时要积分，放这里想要的时候就可以直接过来拿，不用死去搜索啦。呵呵 插值方法的Matlab实现 一维数据插值 MATLAB中用函数interp1来拟合一维数据，语法是YI = INTERP1(X,Y,XI,方法) 其中（X，Y）是已给的数据点，XI 是插值点， 其中方法主要有 'linear' -线性插值，默认 'pchip' -逐段三次Hermite插值 'spline' -逐段三次样条函数插值 其中最后一种插值的曲线比较平滑 例： x=0:.12:1; x1=0:.02:1;%（其中x=0:.12:1表示显示的插值点，x1=0:.02:1表示插值的步长） y=(x.^2-3*x+5).*exp(-5*x).*sin(x); plot(x,y,'o'); hold on; y1=interp1(x,y,x1,'spline'); plot(x1,y1,':') 如果要根据样本点求函数的定积分，而函数又是比较光滑的，则可以用样条函数进行插值后再积分，在MATLAB中可以编写如下程序： function y=quadspln(x0,y0,a,b) f=inline('interp1(x0,y0,x,''spline'')','x','x0','y0'); y=quadl(f,a,b,1e-8,[],x0,y0); 现求sin(x)在区间[0，pi]上的定积分，只取5点 x0=[0,0.4,1,2,pi]; y0=sin(x0); I=quadspln(x0,y0,0,pi) 结果得到的值为2.01905，精确值为2 求一段matlab插值程序 悬赏分：20 - 解决时间：2009-12-26 19:57 已知5个数据点：x=[0.25 0.5 0.75 1] y=[0 0.3104 0.6177 0.7886 1] ，求一段matlab插值程序，求过这5个数据点的插值多项式，并在x-y坐标中画出y=f（x）图形，并且求出f （x）与x轴围成图形的面积（积分），不胜感激！ 使用Lagrange 插值多项式的方法： 首先把下面的代码复制到M文件中，保存成lagran function [C,L]=lagran(X,Y) % input - X is a vector that contains a list of abscissas % - Y is a vector that contains a list of ordinates 6.1 插值问题及其误差 6.1.2 与插值有关的MATLAB 函数 (一) POLY2SYM函数 调用格式一：poly2sym (C) 调用格式二：f1=poly2sym(C,'V') 或f2=poly2sym(C, sym ('V') ), (二) POLYVAL函数 调用格式：Y = polyval(P,X) (三) POLY函数 调用格式：Y = poly (V) (四) CONV函数 调用格式：C =conv (A, B) 例 6.1.2求三个一次多项式、和的积.它们的零点分别依次为0.4,0.8,1.2. 解我们可以用两种MATLAB程序求之. 方法1如输入MATLAB程序 >> X1=[0.4,0.8,1.2]; l1=poly(X1), L1=poly2sym (l1) 运行后输出结果为 l1 = 1.0000 - 2.4000 1.7600 -0.3840 L1 = x^3-12/5*x^2+44/25*x-48/125 方法2如输入MATLAB程序 >> P1=poly(0.4);P2=poly(0.8);P3=poly(1.2); C =conv (conv (P1, P2), P3) , L1=poly2sym (C) 运行后输出的结果与方法1相同. (五) DECONV 函数 调用格式：[Q,R] =deconv (B,A) (六) roots(poly(1:n))命令 调用格式：roots(poly(1:n)) (七) det(a*eye(size (A)) - A)命令 调用格式：b=det(a*ey e(size (A)) - A) 6.2 拉格朗日（Lagrange）插值及其MATLAB程序 6.2.1 线性插值及其MATLAB程序 例 6.2.1 已知函数在上具有二阶连续导数，，且满足条件 .求线性插值多项式和函数值，并估计其误差. 解输入程序 >> X=[1,3];Y=[1,2]; l01= poly(X(2))/( X(1)- X(2)), l11= poly(X(1))/( X(2)- X(1)), l0=poly2sym (l01),l1=poly2sym (l11), P = l01* Y(1)+ l11* Y(2), L=poly2sym (P),x=1.5; Y = polyval(P,x) 运行后输出基函数l0和l1及其插值多项式的系数向量P（略）、插值多项式L和插值Y为l0 = l1 = L = Y = -1/2*x+3/2 1/2*x-1/2 1/2*x+1/2 1.2500 输入程序 >> M=5;R1=M*abs((x-X(1))* (x-X(2)))/2 m a t l a b实现插值法和 曲线拟合 插值法和曲线拟合 电子科技大学 摘要：理解拉格朗日多项式插值、分段线性插值、牛顿前插，曲线拟合，用matlab编程求解函数，用插值法和分段线性插值求解同一函数，比较插值余项；用牛顿前插公式计算函数，计算函数值；对于曲线拟 合，用不同曲线拟合数据。 关键字：拉格朗日插值多项式；分段线性插值；牛顿前插；曲线拟合 引言： 在数学物理方程中，当给定数据是不同散点时，无法确定函数表达式，求解函数就需要很大的计算量，我们有多种方法对给定的表格函数进行求解，我们这里，利用插值法和曲线拟合对函数进行求解，进一步了解函数性质，两种方法各有利弊，适合我们进行不同的散点函数求解。 正文： 一、插值法和分段线性插值 1拉格朗日多项式原理 对某个多项式函数，已知有给定的k + 1个取值点： 其中对应着自变量的位置，而对应着函数在这个位置的取值。 假设任意两个不同的x j都互不相同，那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为： 其中每个为拉格朗日基本多项式（或称插值基函数），其表达式为： [3] 拉格朗日基本多项式的特点是在上取值为1，在其它的点 上取值为0。 2分段线性插值原理 给定区间[a,b], 将其分割成a=x 0 实验四用MATLAB实现拉格朗日插值、分段线性插值一、实验目的： 1）学会使用MATLAB软件； 2）会使用MATLAB软件进行拉格朗日插值算法和分段线性差值算法； 二、实验内容： 1用MATLAB实现y = 1./(x.^2+1);（-1<=x<=1）的拉格朗日插值、分段线性 2.选择以下函数，在n个节点上分别用分段线性和三次样条插值的方法，计算m个插值点的函数值，通过数值和图形的输出，将插值结果与精确值进行比较，适当增加n，再作比较，由此作初步分析： (1).y=sinx;( 0≤x≤2π) (2).y=(1-x^2)（-1≤x≤1） 三、实验方法与步骤： 问题一用拉格朗日插值法 1）定义函数：y = 1./(x.^2+1);将其保存在f.m 文件中，程序如下： function y = f1(x) y = 1./(x.^2+1); 2）定义拉格朗日插值函数：将其保存在lagrange.m 文件中，具体实现程序编程如下：function y = lagrange(x0,y0,x) m = length(x); /区间长度/ n = length(x0); for i = 1:n l(i) = 1; end for i = 1:m for j = 1:n for k = 1:n if j == k continue; end l(j) = ( x(i) -x0(k))/( x0(j) - x0(k) )*l(j); end end end y = 0; for i = 1:n y = y0(i) * l(i) + y; end 3）建立测试程序，保存在text.m文件中，实现画图：x=-1:0.001:1; y = 1./(x.^2+1); p=polyfit(x,y,n); py=vpa(poly2sym(p),10) plot_x=-5:0.001:5; f1=polyval(p,plot_x); figure plot(x,y,‘r',plot_x,f1)matlab插值法实例
lagrange插值分段线性插值matlab代码
实验5 双线性插值
插值算法与matlab代码
matlab旋转+双线性插值
拉格朗日插值matlab程序
数字图像处理(Matlab复习代码)
matlab插值法,迭代法程序
matlab计算拉格朗日牛顿及分段线性插值的程序
双线性插值
matlab旋转实现(最近邻值,双线性,三次卷积插值实现插值)
用MATLAB实现拉格朗日插值和分段线性插值
matlab插值计算
函数的插值方法及matlab程序
matlab实现插值法和曲线拟合电子教案
实验四用MATLAB实现拉格朗日插值、分段线性插值
相关主题
文本预览