精华内容
下载资源
问答
  • B样条插值算法

    万次阅读 2017-02-04 21:48:55
    K阶B样条插值应用非常广泛,其中三阶样条插值意义最大。本文简单实现了K阶B样条插值算法

    K阶B样条插值应用非常广,其中函数性质也是对称的,通过矩阵求逆,很容易得到系数矩阵。从而得到任意点的值,以及n阶导数。

    K阶B样条函数Ω(k,x)的递推性质:Ω(k,x)=1/k*(x+(k+1)/2)*Ω(k-1,x+0.5)-1/k*(x-(k+1)/2)*Ω(k-1,x-0.5)

    K阶B样条函数的导数性质:Ω‘(k,x)=Ω(k-1,x+0.5)-Ω(k-1,x-0.5)

    在求解系数矩阵的时候需要额外的K-1个方程条件,一般情况下设置首尾一阶导数值或是二阶导数值。具体的条件可以自己确定,从而得到完整的方阵。

    代码很简单,如下:

    /**

    * K阶B样条函数 Ω(k,x)

    * Ω(k,x)的导数性质:Ω‘(k,x)=Ω(k-1,x+0.5)-Ω(k-1,x-0.5)

    * Ω(k,x)的递推性质:Ω(k,x)=1/k*(x+(k+1)/2)*Ω(k-1,x+0.5)-
    1/k*(x-(k+1)/2)*Ω(k-1,x-0.5)
    * 根据递推性质即可求得函数
    * @param x 
    * @param k B样条阶数
    * @param n n阶导数,n=0时得到y值
    * @return
    */
    private static double getOkx(double x,double k,int n){
    if(n==0){
    if(k==0){
    double abs=Math.abs(x);
    if(abs==0.5)
    return  0.5f;
    else if(abs>0.5)
    return 0;
    else
    return 1;
    }else{

    return 1/k*(x+(k+1)/2)*getOkx( x+0.5, k-1,0)-
    1/k*(x-(k+1)/2)*getOkx( x-0.5, k-1,0)
    ;
    }
    }else{
    return getOkx( x+0.5f, k-1,n-1)-getOkx( x-0.5, k-1,n-1);
    }

    }


    /**
    * 根据 S(x)=∑(C*Ω(k,x))可得:
    * ∑(C*Ω(k,0))=y[0]
    * ∑(C*Ω(k,1))=y[1]
    * ∑(C*Ω(k,2))=y[2]
    * ......
    * ∑(C*Ω(k,n))=y[n]


    * 以及加入条件方程组,得到N+K矩阵

    * A[i]=∑Ω(k,i),V=y

    * A*C=V --> C=A逆*V

    * @return C
    */
    public double[] getCi(){
    double[][]A=new double[N+K][N+K];
    double[][]V=new double[N+K][1];
    for(int i=0;i<=N;i++){
    for(int j=0;j<N+K;j++){
    A[i][j]=getOkx(i-j+(K-1)/2,K,0);
    }
    V[i][0]=y[i];
    }

    //还有K-1个条件(n阶导数为0)
    //自定义K-1个条件方程
    //这里是x0的导数=y0d,xn的导数=ynd;

    int m=0;
    for(int i=N+1;i<N+K;i++){
    if(m==0){
    for(int j=0;j<N+K;j++){
    A[i][j]=getOkx(0-j+(K-1)/2,K,1);//x0的一阶导数系数
    }
    V[i][0]=y0d;
    }else if(m==1){
    for(int j=0;j<N+K;j++){
    A[i][j]=getOkx(N-j+(K-1)/2,K,1);//xn的一阶导数导数系数
    }
    V[i][0]=ynd;
    }else{
    for(int j=0;j<N+K;j++){
    A[i][j]=getOkx(0-j+(K-1)/2,K,m);
    }
    V[i][0]=0;
    }
    m++;
    }

    double[][]A1=inverseMatrix(A);//矩阵求逆
    double [][]t=times(A1, V);//矩阵相乘
    double []C=new double[N+K];
    for(int i=0;i<t.length;i++){
    C[i]=t[i][0];
    }
    return C;//返回系数

    }


    /**
    * 获取曲线上的点S(x)=∑(C*Ω(k,x))
    * @param x (x,y)
    * @param C B样条函数系数数组
    * @param n n阶导数,n=0时得到S(x)值
    * @return S(x), S'(x),S''(x)....
    */
    public static double getY(double x,double[]C,int n){
    double rst=0;
    if(x<0){
    return 0;
    }
    int size=0;
    for(int i=(int) x;i<C.length;i++){
    double temp=getOkx(x-i+(K-1)/2,K,n);
    rst+=C[i]*temp;
    size++;
    if(size==K+1){
    break;
    }
    }
    return rst;
    }


    完整的例子见:http://download.csdn.net/detail/zhq1314zhq/9746788



    展开全文
  • 基于非均匀B样条插值算法的图像放大
  • K次B样条插值算法

    2017-02-04 11:11:32
    该算法实现了K次B样条插值算法java实现,K可配置。自己根据项目需要写的。
  • 三次四阶b样条插值算法(Deboor算法)

    热门讨论 2010-03-05 23:49:17
    三次四阶b样条插值算法(Deboor算法)可以实现B样条曲线的控制点的拟合
  • B样条曲面插值算法

    千次阅读 2010-08-30 13:48:00
    前一阵子在搞曲面插值,网上的代码很少,自己就参考B样条曲面公式写了一个。代码如下:矩阵类:public class MatrixUtil { /// /// 构造函数 构造方阵 /// /// 行列数相等矩阵 public MatrixUtil(int ...

     

    前一阵子在搞曲面插值,网上的代码很少,自己就参考B样条曲面公式写了一个。代码如下:

    矩阵类:

     

     

    计算方法:

     

    注:计算曲面化数据方法CalSurfaceData为计算的主要方法。举例:

    我传入一个61*61的二维数组float[][] data,密度density为3。曲面插值后的数据就应该为183*183。上面代码我只对Z值进行插值,因为是规则的矩形化网格,X坐标。Y坐标只进行了普通插值。下图为插值实现方式:

    即由周围16个红点位置插值出灰色区域内的XYZ值。 所有插值完成之后 还要进行补插周围一圈的2*2矩阵值。CalBoundVal方法即为边界插值。因为我插值的数据是河道数据,可能包含陆地点,所以我在上述代码中加入了一些关于陆地点插值的判断。代码写的比较多,难免有疏漏,就是提供个参考吧。

    展开全文
  • 非均匀三次B样条曲线插值的Jacobi-PIA算法
  • 基于Coons-Gordon造型原理,研究了插值两族相交截面线采样点的B样条曲面双向插值造型算法。参数化各采样点并计算每条截面线的节点矢量,估算每条截面线对应的曲面参数,根据每条截面线的节点分布以及另一族截面线...
  • 通过扫描线算法和不连续B样条进行图像插值
  • B样条曲线插值

    2020-07-05 09:07:21
    B样条曲线反算控制点 1 De Boor算法 设u∈[uj,uj+1)u\in[u_j,u_{j+1})u∈[uj​,uj+1​),Vi,0=ViV_{i,0}=V_iVi,0​=Vi​,对于i=j−p,⋯ ,ji=j-p,\cdots,ji=j−p,⋯,j 令 Vi,k=ui+p+1−k−uui+p+1−k−uiVi−1,k−1...

    B样条曲线反算控制点

    1 De Boor算法

    u[uj,uj+1)u\in[u_j,u_{j+1})Vi,0=ViV_{i,0}=V_i,对于i=jp,,ji=j-p,\cdots,j

    Vi,k=ui+p+1kuui+p+1kuiVi1,k1+uuiui+p+1kuiVi,k1,k=1,,p,i=jp+k,,jV_{i,k}=\dfrac{u_{i+p+1-k}-u}{u_{i+p+1-k}-u_i}V_{i-1,k-1}+\dfrac{u-u_i}{u_{i+p+1-k}-u_i}V_{i,k-1},\quad k=1,\cdots,p,\quad i=j-p+k,\cdots,j

    其中ViV_i为控制点,pp为B样条的幂次,P(u)P(u)为B样条曲线,则
    P(u)=Vj,p. P(u)=V_{j,p}.

    2 B样条曲线上点的计算公式

    u[uj,uj+1)u\in[u_j,u_{j+1}),则

    p=3p=3时,
    P(u)=(uj+1u)3(uj+1uj)(uj+1uj1)(uj+1uj2)Vj3+[(uj+1u)2(uuj2)(uj+1uj)(uj+1uj1)(uj+1uj2)+(uj+1u)(uuj1)(uj+2u)(uj+1uj)(uj+1uj1)(uj+2uj1)+(uuj)(uj+2u)2(uj+1uj)(uj+2uj)(uj+2uj1)]Vj2+[(uj+1u)(uuj1)2(uj+1uj)(uj+1uj1)(uj+2uj1)+(uuj)(uj+2u)(uuj1)(uj+1uj)(uj+2uj)(uj+2uj1)+(uuj)2(uj+3u)(uj+1uj)(uj+2uj)(uj+3uj)]Vj1+(uuj)3(uj+1uj)(uj+2uj)(uj+3uj)Vj. \begin{aligned} P(u)=&\dfrac{(u_{j+1}-u)^3}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})(u_{j+1}-u_{j-2})}V_{j-3} +\left[\dfrac{(u_{j+1}-u)^2(u-u_{j-2})}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})(u_{j+1}-u_{j-2})}\right. \\ & \left. +\dfrac{(u_{j+1}-u)(u-u_{j-1})(u_{j+2}-u)}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})(u_{j+2}-u_{j-1})} +\dfrac{(u-u_j)(u_{j+2}-u)^2}{(u_{j+1}-u_j)(u_{j+2}-u_j)(u_{j+2}-u_{j-1})}\right]V_{j-2} \\ &+\left[\dfrac{(u_{j+1}-u)(u-u_{j-1})^2}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})(u_{j+2}-u_{j-1})} +\dfrac{(u-u_j)(u_{j+2}-u)(u-u_{j-1})}{(u_{j+1}-u_j)(u_{j+2}-u_j)(u_{j+2}-u_{j-1})} \right. \\ &\left. +\dfrac{(u-u_j)^2(u_{j+3}-u)}{(u_{j+1}-u_j)(u_{j+2}-u_j)(u_{j+3}-u_j)}\right]V_{j-1} +\dfrac{(u-u_j)^3}{(u_{j+1}-u_j)(u_{j+2}-u_j)(u_{j+3}-u_j)}V_j. \end{aligned}
    p=2p=2时,
    P(u)=(uj+1u)2(uj+1uj)(uj+1uj1)Vj2+[(uj+1u)(uuj1)(uj+1uj)(uj+1uj1)+(uuj)(uj+2u)(uj+1uj)(uj+2uj)]Vj1+(uuj)2(uj+1uj)(uj+2uj)Vj. \begin{aligned} P(u)=&\dfrac{(u_{j+1}-u)^2}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})}V_{j-2} +\left[\dfrac{(u_{j+1}-u)(u-u_{j-1})}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})} \right. \\ &\left.+\dfrac{(u-u_j)(u_{j+2}-u)}{(u_{j+1}-u_j)(u_{j+2}-u_j)}\right]V_{j-1} +\dfrac{(u-u_j)^2}{(u_{j+1}-u_j)(u_{j+2}-u_j)}V_j. \end{aligned}
    p=1p=1时,
    P(u)=uj+1uuj+1ujVj1+uujuj+1ujVj. P(u)=\dfrac{u_{j+1}-u}{u_{j+1}-u_j}V_{j-1}+\dfrac{u-u_j}{u_{j+1}-u_j}V_j.

    3 3次B样条曲线反求控制点

    3.1 3次B样条曲线

    p=3p=3并且u[uj,uj+1)u\in[u_j,u_{j+1})时,
    P(u)=i=0nNi,3(u)Vi=i=j3jNi,3(u)Vi. P(u)=\sum^n_{i=0}N_{i,3}(u)V_i=\sum^j_{i=j-3}N_{i,3}(u)V_i.
    其中,
    Nj3,3(u)=(uj+1u)3(uj+1uj)(uj+1uj1)(uj+1uj2).Nj2,3(u)=(uj+1u)2(uuj2)(uj+1uj)(uj+1uj1)(uj+1uj2)+(uj+1u)(uuj1)(uj+2u)(uj+1uj)(uj+1uj1)(uj+2uj1)+(uuj)(uj+2u)2(uj+1uj)(uj+2uj)(uj+2uj1).Nj1,3(u)=(uj+1u)(uuj1)2(uj+1uj)(uj+1uj1)(uj+2uj1)+(uuj)(uj+2u)(uuj1)(uj+1uj)(uj+2uj)(uj+2uj1)+(uuj)2(uj+3u)(uj+1uj)(uj+2uj)(uj+3uj).Nj,3(u)=(uuj)3(uj+1uj)(uj+2uj)(uj+3uj). \begin{aligned} N_{j-3,3}(u)=&\dfrac{(u_{j+1}-u)^3}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})(u_{j+1}-u_{j-2})}. \\ N_{j-2,3}(u)=&\dfrac{(u_{j+1}-u)^2(u-u_{j-2})}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})(u_{j+1}-u_{j-2})} +\dfrac{(u_{j+1}-u)(u-u_{j-1})(u_{j+2}-u)}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})(u_{j+2}-u_{j-1})} \\ &+\dfrac{(u-u_j)(u_{j+2}-u)^2}{(u_{j+1}-u_j)(u_{j+2}-u_j)(u_{j+2}-u_{j-1})}. \\ N_{j-1,3}(u)=&\dfrac{(u_{j+1}-u)(u-u_{j-1})^2}{(u_{j+1}-u_j)(u_{j+1}-u_{j-1})(u_{j+2}-u_{j-1})} +\dfrac{(u-u_j)(u_{j+2}-u)(u-u_{j-1})}{(u_{j+1}-u_j)(u_{j+2}-u_j)(u_{j+2}-u_{j-1})} \\ &+\dfrac{(u-u_j)^2(u_{j+3}-u)}{(u_{j+1}-u_j)(u_{j+2}-u_j)(u_{j+3}-u_j)}. \\ N_{j,3}(u)=&\dfrac{(u-u_j)^3}{(u_{j+1}-u_j)(u_{j+2}-u_j)(u_{j+3}-u_j)}. \end{aligned}
    于是,可以得到
    Nj3,3(uj)=(uj+1uj)2(uj+1uj1)(uj+1uj2).Nj2,3(uj)=(uj+1uj)(ujuj2)(uj+1uj1)(uj+1uj2)+(ujuj1)(uj+2uj)(uj+1uj1)(uj+2uj1).Nj1,3(uj)=(ujuj1)2(uj+1uj1)(uj+2uj1).Nj,3(uj)=0. \begin{aligned} N_{j-3,3}(u_j)=&\dfrac{(u_{j+1}-u_j)^2}{(u_{j+1}-u_{j-1})(u_{j+1}-u_{j-2})}. \\ N_{j-2,3}(u_j)=&\dfrac{(u_{j+1}-u_j)(u_j-u_{j-2})}{(u_{j+1}-u_{j-1})(u_{j+1}-u_{j-2})} +\dfrac{(u_j-u_{j-1})(u_{j+2}-u_j)}{(u_{j+1}-u_{j-1})(u_{j+2}-u_{j-1})}. \\ N_{j-1,3}(u_j)=&\dfrac{(u_j-u_{j-1})^2}{(u_{j+1}-u_{j-1})(u_{j+2}-u_{j-1})}. \\ N_{j,3}(u_j)=&0. \end{aligned}
    则有
    P(uj)=i=j3jNi,3(uj)Vi=Nj3,3(uj)Vj3+Nj2,3(uj)Vj2+Nj1,3(uj)Vj1. P(u_j)=\sum^j_{i=j-3}N_{i,3}(u_j)V_i=N_{j-3,3}(u_j)V_{j-3}+N_{j-2,3}(u_j)V_{j-2} +N_{j-1,3}(u_j)V_{j-1}.

    3.2 根据型值点列出方程组

    假设

    节点矢量{0,0,0,u1,u2,,un1,un,1,1,1}\{0,0,0,u_1,u_2,\cdots,u_{n-1},u_n,1,1,1\}

    型值点序列{P1,,Pn}\{P_1,\cdots,P_n\}

    控制点序列{V1,,Vn+2}\{V_1,\cdots,V_{n+2}\}


    u1=0,un=1.V1=P1,Vn+2=Pn. \begin{aligned} &u_1=0,\quad u_n=1. \\ &V_1=P_1,\quad V_{n+2}=P_n. \\ \end{aligned}
    并且
    (N2,3(u2)N3,3(u2)N4,3(u2)Nn1,3(un1)Nn,3(un1)Nn+1,3(un1))(V2Vn+1)=(P2Pn+1) \begin{aligned} \begin{pmatrix} N_{2,3}(u_2) & N_{3,3}(u_2) & N_{4,3}(u_2) & & \\ & \ddots & \ddots & \ddots & \\ & & N_{n-1,3}(u_{n-1}) & N_{n,3}(u_{n-1}) & N_{n+1,3}(u_{n-1}) \\ \end{pmatrix} \begin{pmatrix} V_2 \\ \vdots \\ V_{n+1} \\ \end{pmatrix} =\begin{pmatrix} P_2 \\ \vdots \\ P_{n+1} \\ \end{pmatrix} \end{aligned}
    上式共有n2n-2个方程,但有nn个控制点需要求解,因此需要补充端点条件。

    3.3 端点条件

    如果给定端点切矢,P1=T1P_1'=T_1Pn=TnP_n'=T_n,即
    T1=P1=3u2u1(V2V1).Tn=Pn=3unun1(Vn+2Vn+1). \begin{aligned} &T_1=P_1'=\dfrac{3}{u_2-u_1}(V_2-V_1). \\ &T_n=P_n'=\dfrac{3}{u_n-u_{n-1}}(V_{n+2}-V_{n+1}). \\ \end{aligned}
    如果取自由端点条件,则P1=0P_1''=0Pn=0P_n''=0,即
    0=P1=6(u2u1)2(V1V2)6(u2u1)(u3u1)(V2V3).0=Pn=6(unun1)(unun2)(VnVn+1)6(unun1)2(Vn+1Vn+2). \begin{aligned} &0=P_1''=\dfrac{6}{(u_2-u_1)^2}(V_1-V_2)-\dfrac{6}{(u_2-u_1)(u_3-u_1)}(V_2-V_3). \\ &0=P_n''=\dfrac{6}{(u_n-u_{n-1})(u_n-u_{n-2})}(V_n-V_{n+1})-\dfrac{6}{(u_n-u_{n-1})^2}(V_{n+1}-V_{n+2}). \\ \end{aligned}
    于是得到如下方程组
    (b1c1N2,3(u2)N3,3(u2)N4,3(u2)Nn1,3(un1)Nn,3(un1)Nn+1,3(un1)anbn)(V2V3VnVn+1)=(d1P2Pn1dn) \begin{aligned} \begin{pmatrix} b_1 & c_1 & & & \\ N_{2,3}(u_2) & N_{3,3}(u_2) & N_{4,3}(u_2) & & \\ & \ddots & \ddots & \ddots & \\ & & N_{n-1,3}(u_{n-1}) & N_{n,3}(u_{n-1}) & N_{n+1,3}(u_{n-1}) \\ & & & a_n & b_n \\ \end{pmatrix} \begin{pmatrix} V_2 \\ V_3 \\ \vdots \\ V_n \\ V_{n+1} \\ \end{pmatrix} =\begin{pmatrix} d_1 \\ P_2 \\ \vdots \\ P_{n-1} \\ d_n \\ \end{pmatrix} \end{aligned}
    其中,

    当给定端点切矢时,
    b1=3u2u1.c1=0.d1=P1+3u2u1P1.an=0.bn=3unun1.dn=3unun1PnPn. \begin{aligned} &b_1=\dfrac{3}{u_2-u_1}. \\ &c_1=0. \\ &d_1=P_1'+\dfrac{3}{u_2-u_1}P_1. \\ &a_n=0. \\ &b_n=\dfrac{3}{u_n-u_{n-1}}. \\ &d_n=\dfrac{3}{u_n-u_{n-1}}P_n-P_n'. \end{aligned}
    当取自由端点条件时,
    b1=6(u2u1)2+6(u2u1)(u3u1).c1=6(u2u1)(u3u1).d1=6(u2u1)2P1.an=6(unun1)(unun2).bn=6(unun1)(unun2)+6(unun1)2.dn=6(unun1)2Pn. \begin{aligned} &b_1=\dfrac{6}{(u_2-u_1)^2}+\dfrac{6}{(u_2-u_1)(u_3-u_1)}. \\ &c_1=-\dfrac{6}{(u_2-u_1)(u_3-u_1)}. \\ &d_1=\dfrac{6}{(u_2-u_1)^2}P_1. \\ &a_n=-\dfrac{6}{(u_n-u_{n-1})(u_n-u_{n-2})}. \\ &b_n=\dfrac{6}{(u_n-u_{n-1})(u_n-u_{n-2})}+\dfrac{6}{(u_n-u_{n-1})^2}. \\ &d_n=\dfrac{6}{(u_n-u_{n-1})^2}P_n. \\ \end{aligned}
    两种端点条件可以混合使用。但首、末端点都分别与首、末型值点重合,即
    V1=P1,Vn+2=Pn. V_1=P_1,\quad V_{n+2}=P_n.

    3.4 追赶法求解方程组

    A=(a1c1b2a2c2bn1an1cn1bnan)A=\begin{pmatrix}a_1 & c_1 & & & \\b_2 & a_2 & c_2 & & \\& \ddots & \ddots & \ddots & \\ & & b_{n-1} & a_{n-1} & c_{n-1} \\ & & & b_n & a_n \end{pmatrix},对于一切ij>1|i-j|>1Aij=0A_{ij}=0

    如果AA满足Doolittle分解,即
    LU=A=(1p21pn11pn1)(q1c1q2c2qn1cn1qn) \begin{aligned} LU=A=\begin{pmatrix} 1 & & & & \\ p_2 & 1 & & & \\ & \ddots & \ddots & & \\ & & p_{n-1} & 1 & \\ & & & p_n & 1 \end{pmatrix} \begin{pmatrix} q_1 & c_1 & & & \\ & q_2 & c_2 & & \\ & & \ddots & \ddots & \\ & & & q_{n-1} & c_{n-1} \\ & & & & q_n \\ \end{pmatrix} \end{aligned}
    其中
    {q1=a1,pi=biqi1,i=2,,nqi=aipici1,i=2,,n \begin{aligned} \begin{cases} q_1=a_1, \\ p_i=\dfrac{b_i}{q_{i-1}},\quad i=2,\cdots,n \\ q_i=a_i-p_ic_{i-1},\quad i=2,\cdots,n \\ \end{cases} \end{aligned}
    于是,求Ax=bAx=b等价于求{Ly=bUx=y\begin{cases}Ly=b \\Ux=y\end{cases},其中b=(b1,,bn)Tb=(b_1,\cdots,b_n)^T,故有
    (1p21pn1)(y1y2yn)=(b1b2bn) \begin{aligned} \begin{pmatrix} 1 & & & \\ p_2 & 1 & & \\ & \ddots & \ddots & \\ & & p_n & 1 \\ \end{pmatrix} \begin{pmatrix} y_1 \\ y_2 \\ \vdots \\ y_n \\ \end{pmatrix} =\begin{pmatrix} b_1 \\ b_2 \\ \vdots \\ b_n \\ \end{pmatrix} \end{aligned}
    解得{y1=b1,yi=bipiyi1,i=2,,n.\begin{cases}y_1=b_1, \\y_i=b_i-p_iy_{i-1},\quad i=2,\cdots,n.\end{cases}

    再由
    (q1c1qn1cn1qn)(x1xn1xn)=(y1yn1yn) \begin{aligned} \begin{pmatrix} q_1 & c_1 & & \\ & \ddots & \ddots & \\ & & q_{n-1} & c_{n-1} \\ & & & q_n \\ \end{pmatrix} \begin{pmatrix} x_1 \\ \vdots \\ x_{n-1} \\ x_n \\ \end{pmatrix} =\begin{pmatrix} y_1 \\ \vdots \\ y_{n-1} \\ y_n \\ \end{pmatrix} \end{aligned}
    解得{xn=ynqn,xi=yicixi+1qi,i=n1,,1.\begin{cases}x_n=\dfrac{y_n}{q_n},\\ x_i=\dfrac{y_i-c_ix_{i+1}}{q_i},\quad i=n-1,\cdots,1.\end{cases}

    展开全文
  • 一个简单的 双三次B样条算法 实现的 图像放大 B样条算法图像放大的学习和实现 一、总体设计思路 本次B样条图像放大的实现建立在之前做的图像的显示,双线性插值图像放大的基础上。因为有之前的基础,所以实际上...

    一个简单的 双三次B样条算法 实现的 图像放大

    B样条算法图像放大的学习和实现

    一、总体设计思路

    本次B样条图像放大的实现建立在之前做的图像的显示,双线性插值图像放大的基础上。因为有之前的基础,所以实际上本次我只是要实现基于B样条的插值即可。本文也只是讨论B样条插值的实现。本次实现是基于双三次B样条曲面的,且对每一个像素点都采用B样条插值的方法来确定颜色值。B样条插值有以下几个步骤:

    Step 1:取型值点 16个

    Step2:计算每个型值点的(u,v)

    Step3:根据B样条曲面矩阵形式的公式(由控制点及(u,v)得到曲面上的点),利用已有的型值点及(u,v),反算出16个控制点。

    Step4:计算出欲求像素点的(u,v)。得到对应的像素值,并赋给对应的位置。

    二、原理

    1、理论支持

    B样条算法图像放大的学习和实现

    B样条算法图像放大的学习和实现

    我的整个实现过程都是基于这两个公式的,而经过试验。这两个公式都是好用的。 整个

    2、公式推导

    其实这个实现过程中,最关键的部分就是由该公式由16个型值点构建方程组反求出16个控制点,之后再用这里的公式求点即可。

    对此通过计算易得:

    设[u3,u2,u,1]XM1=pre[4]; M2X[v3,v2,v,1]T=behind[4];

    pre的四个元素为 a b c d bdhind的四个元素为e f g h.易得

    p(u,v)=(e(aP00+bP10+cP20+dP30)+f(aP01+bP11+cP21+dP31)+g(aP02+bP12+cP22+dP32)+h(aP03+bP13+cP23+dP33))/36;

    也即

    一个简单的 双三次B样条算法 实现的 图像放大

    p(u,v)=(aeP00+beP10+ceP20+deP30+afP01+bfP11+cfP21+dfP31+agP02+bgP12+cgP22+dgP32+ahP03+bhP13+chP23+dhP33)/36;

    由该式即可设立一个 16个未知数的方程。

    由16个点就可以得到16个方程并解得16个控制点。

    三、具体实现方案及细节

    1、具体实现方案

    ①型值点的选取

    先将新图中的坐标转换成原图中的坐标.也即x=X/m;y=Y/m;(m为边长需要放大的倍数)取周边的16个点。按得到的坐标的值就可以往左右上下各延伸2个单位。若上或下不足(左或右)那么就朝另一边发展。

    ②(u,v)的确定

    在我的实现中,(u,v)的确定比较简单。(在平面上确定)。因为型值点是取的一个4X4的矩阵。假设型值点在这个矩阵中的下标为(i,j);则u=i/3.0;v=j/3.0; 假设最小型值点的x为preX,y为preY.则新图中像素点(X,Y)的u=(X-preX)/3.0;v=(Y-preY)/3.0;

    ③解方程

    1) Ax=b->x=A^-1B;因为我的(u,v)对于每一个像素点来说是不变的。(参看(u,v)的确定)所以该式中的A是不变的。

    因此给出一个方法就是直接得到A之后对A求逆并保存下来之后每一次就只要做矩阵的乘法即可得到结果。

    ④根据公式加上16个的控制点和对应的(u,v)即可求出像素值。

    2、实现过程中的遇到的问题和解决

    ①不会确定(u,v),心里想着只能靠(x,y)但是不知道到底怎么弄。最后还是老师指点空间可以参数化到平面上,才确定了现在取型值点的方法。

    ②解方程我自己先写了一个高斯消元法,但是效率太低所以换成了方案中的矩阵求逆。然后矩阵运算。

    ③求逆我也写了一个类,解决了行列式的解、求伴随矩阵、求逆序数等问题。最后因为矩阵16X16行列式定义法求解直接卡死了。最后是借助matlab完成了这一步工作。同时这个过程中用代码实现,矩阵的产生和整理 也帮了我很多忙。

    ④解完方程之后,在程序能运行速度大幅提升的时候心里很高兴。但是这股高兴劲马上就被出来的效果给浇灭了。出来的效果总是无缘无故多出一些点。经查,是因为我没有处理求出的像素值超过范围的情况。也即求出的像素值超过了255或小于0。解决方案就简单了,将大于255的改为255,小于0的改为0.



    展开全文
  • 三次B样条曲线拟合算法

    万次阅读 多人点赞 2017-01-17 22:10:28
    三次B样条曲线方程B样条曲线分为近似拟合和插值拟合,所谓近似拟合就是不过特征点,而插值拟合就是通过特征点,但是插值拟合需要经过反算得到控制点再拟合出过特征点的B样条曲线方程。这里会一次介绍两种拟合算法。...
  • 整理中。。。 
  • 对于B样条曲线,您可以从几种结和参数选择算法中进行选择,并且可以自由选择控制点的数量。 即使使用“误差极限设置”也可以找到误差最小的曲线。 另外,您需要为面板提供包含您的点的输入文件。 该文件必须为格式或...
  • 为解决空间点列曲线插值...此算法具有明确的物理意义,数学模型简单,生成的曲线结构紧凑,弥补了B样条插值法的不足,能够生成形态复杂的空间曲线,在计算机造型、反求工程及运动轨迹描述等工程领域有重要的实用价值。
  • 并借助控制多边形在每次加细过程中新旧控制顶点对应的几何位置关系,给出一种新的,z次均匀B样条曲线细分算法,基于该算法构造出带有形状参数的局部插值约束的奇次均匀B样条细分曲线。通过理论和算例说明,该算法...
  • QT绘制B样条曲线 三次样条插值算法 三次样条插值(Cubic Spline Interpolation)及代码实现(C语言)
  • 多级B样条曲线的分散数据插值 该库提供了用C ++ 11实现的[1]中的自适应MBA算法。 这是用于分散N维数据插值和逼近的快速算法。 还提供了Python绑定。 C ++中的2D插值示例: # include int main () { // ...
  • 降阶算法B样条曲线和曲面设计的一个基本算法,它广泛应用于组合曲线,蒙皮或扫描曲面等设计中。Piegl与Tiller曾给出B样条曲线的降阶方法。本文给出了解决更一般的非端点插值B样条曲线降阶的方法。新的方法主要是通过...
  • 对10*10的二维高程矩阵进行插值,得到较平滑的三维山地图,使用A星算法得到起点到终点的最优路径,使用B样条插值法对路径进行平滑处理。
  • 用数学归纳法能够证明B样条插值算法的正确性,但比较繁复,限于篇幅不做讨论。 有关B样条插值算法正确性反而是在Matlab实现中发现的,以下列出目前插值算法实现代码:function [ new_vertices, new_knots ] = ...
  • 介绍图像基于卷积计算三种双立方插值函数,分别是三角线性分布、Bell钟型分布、B样条曲线分布 通过它们实现图像双立方插值放大。
  • B样条插值点反求控制点

    千次阅读 2005-04-27 21:56:00
    三次周期B样条曲线的算法 0 £ u和四个控制点p0,p1,p2和p3. 设P(u)是一个三次周期B样条,满足条件: P(0) = (p0 + 4p1 + p2)/6, P(1) = (p1 + 4p2 + p3)/6, P¢(0) = (p2 – p0)/2, P¢(1) = (p3 – p1)/2. ...
  • 指出了Piegl与Tiller所述的B样条曲线升阶方法中的问题,提出了解决问题的新方法,即一个新的端点插值方法,利用此方法对Piegl与Tiller的升阶方法进行改进,使之能够解决所有均匀及非均匀B样条曲线的升阶问题。
  • 在做相关项目需要解决B样条插值问题,记录如下 感谢以下参考资料的帮助 % ref: 闭合 B 样条曲线控制点的快速求解算法及应用 % http://www.doc88.com/p-5714423317458.html % ...
  • 用VC和OpenGL写的B样条空间曲线。算法包括曲线的绘制与交互编辑。
  • 三种常见的图像处理双三次插值算法双立方插值计算涉及16像素,间(i’, j’)像中的包括小数部分的像素坐标。dx表示X方向的小数坐标。dy表示Y方向的小数坐标。...常见有基于三角取值、Bell分布表达、B样条曲...
  • matlab实现b样条线拼接

    2012-05-31 08:39:58
    jlu的作业必备啊,*chexueyuan做b样条作业,要的,你懂的.用用Matlab实现了3次样条曲线插值算法边界条件取为自然.doc ,就是三次b样条线的算法,下了就知道

空空如也

空空如也

1 2 3 4 5
收藏数 89
精华内容 35
关键字:

b样条插值算法