function(x[],y[],m,n,x1[],y1[])
   x[]---已知点横坐标
   y[]---已知点纵坐标
   m-----已知点数
   n-----插值点数
   x1[]--插值后点横坐标
   y1[]--插值后点纵坐标 请问有这样的函数么?

解决方案 »

  1.   

    //文件名:splineinsert.cpp
    /*************************************************************************
    *
    * 函数名称:
    *   splineinsert()
    *
    * 参数:
    * InContourSide        - 存储内轮廓的节点链表
    * 返回值:
    *   BOOL               - 运算成功返回TRUE,否则返回FALSE。
    *
    * 说明:
    * 该函数用于输出提取出的轮廓数据。

    * 要求轮廓图象已经过扫描,建立起轮廓图象,对轮廓数据进行三次样条插值。
    ************************************************************************///通过此函数调用三次样条插值
    BOOL WINAPI splineinsert(ChainHead*  InContourSide,ChainHead* ReplaceSide)
    {
    int i,j,k;
    int direction;
    int LongSize;
    int TotalSize,other;
    ChainHead*  pChain;
    Node*       pNode = NULL; 
    Node*       mNode = NULL;
    Node*       sNode = NULL;
    Node*       BreakNode = NULL; int         MidFlag;
    Node*       MinNode = NULL;
        Node*       MaxNode = NULL;
    //边界条件,即假定边界为切矢无穷大
    double  yp1 = 1.7E30;
    double  ypn = 1.7E30;

    /*
    此处代码省,目的是把轮廓分成上下两条
    因为样条插值要求沿一个方向是增大/减小的
    */ //x,y 是输入的坐标值 n是数据点的个数 
    //yp1,ypn是端始点的二阶导数值 y2为输出的二阶导数值
    spline(Xb, Yb, k, yp1, ypn, y2); LongSize = (int)(MaxNode->x-MinNode->x)+1;  //插值点的个数
    double*    sig = new double[LongSize];
    double*    xcha = new double[LongSize];

    xcha[0] = Xb[0];
    for(i=1; i<LongSize; i++)
    {
    xcha[i] = xcha[i-1]+1;
    }

    for(i=0; i<=LongSize; i++)
    {
    splint(Xb, Yb, y2, k, xcha[i], sig[i]);
    }/*
    此处代码省,目的是把轮廓分成上下两条
    因为样条插值要求沿一个方向是增大/减小的
    */
    //x,y 是输入的坐标值 n是数据点的个数 
    //yp1,ypn是端始点的二阶导数值 y3为输出的二阶导数值
    spline(Xc, Yc, other, yp1, ypn, y3); for(i=0; i<=LongSize; i++)
    {
    splint(Xc, Yc, y3, other, xcha[i], sig[i]);
    }

    for(i=1,j=0; i<LongSize; i++,j++)
    {
    mNode = new Node[1];
    pNode->nextNode = mNode;
    pNode = pNode->nextNode;
    pNode->x = xcha[j];
    pNode->y = sig[j];
    }
    pNode->nextNode = NULL;
    double    ypoint;
    Node*  tempNode=NULL; return true;
    }//计算出二次导数
    void WINAPI spline(double x[], double y[], int n, double yp1,double ypn, double y2[])
    {
        double u[100],aaa,sig,p,bbb,ccc,qn,un;
    int i,k;
        if( yp1 > 9.9e+29 )
    {
            y2[1] = 0;
            u[1] = 0;
    }
        else
    {
            y2[1] = -0.5;
            aaa = (y[2] - y[1]) / (x[2] - x[1]);
            u[1] = (3.0 / (x[2] - x[1])) * (aaa - yp1);
        }

        for (i = 2; i<=n-1; i++)
    {
            sig = (x[i] - x[i - 1]) / (x[i + 1] - x[i - 1]);
            p = sig * y2[i - 1] + 2.0;
            y2[i] = (sig - 1.0) / p;
            aaa = (y[i + 1] - y[i]) / (x[i + 1] - x[i]);
            bbb = (y[i] - y[i - 1]) / (x[i] - x[i - 1]);
            ccc = x[i + 1] - x[i - 1];
            u[i] = (6.0 * (aaa - bbb) / ccc - sig * u[i - 1]) / p;
        }

        if (ypn > 9.9e+29 )
    {
            qn = 0.0;
            un = 0.0;
    }
        else
    {
            qn = 0.5;
            aaa = ypn - (y[n] - y[n - 1]) / (x[n] - x[n - 1]);
            un = (3.0/ (x[n] - x[n - 1])) * aaa;
        }
        y2[n] = (un - qn * u[n - 1]) / (qn * y2[n - 1] + 1.0);

        for (k = n - 1;k>=1;k--)
            y2[k] = y2[k] * y2[k + 1] + u[k];
    }
    //计算出相应的三次样条插值点
    void WINAPI splint(double xa[], double ya[], double y2a[], int n, double& x, double& y)
    {
    int klo,khi,k;
    double h,a,b,aaa,bbb;
        klo = 1;
        khi = n;
    loop:   if (khi - klo > 1 )
    {
    k = (khi + klo) / 2;
    if (xa[k] > x)
    khi = k;
    else
    {
    klo = k;
    }
    goto loop;
    }
    h = xa[khi] - xa[klo];
    ASSERT(h!=0); a = (xa[khi] - x) / h;
    b = (x - xa[klo]) / h;
    aaa = a * ya[klo] + b * ya[khi];
    bbb = (a*a*a - a) * y2a[klo] + (b*b*b - b) * y2a[khi];
    y = aaa + bbb * h*h /6.0;
    }你自己看看吧!
      

  2.   


    struct tag_ZYCurve{
    double v0;
    double v1;
    long tag;//附带数据
    };
    typedef CArray<tag_ZYCurve,tag_ZYCurve&> ZYCurve;
    //用三次样条来平滑曲线
    //inCurve输入曲线
    //outCurve输出曲线
    //step步长
    //标志量 flag=1表示输出的曲线中保留输入曲线的点;
    //  flag=2表示输出的曲线中不保留输入曲线的点;
    BOOL _SmoothEpsl(ZYCurve &inCurve, double step,ZYCurve &outCurve,DWORD flag)
    {
    //小于3个数据不能平滑
    if(inCurve.GetSize()<3) return FALSE;
    outCurve.RemoveAll(); int last_row=inCurve.GetUpperBound(); float data1=0;
    float data2=0; //==================================计算n-2个节点上的一阶导数值==================================
    //求hj 个数为n-1
    CArray<double,double> hj;
    hj.SetSize(inCurve.GetSize()-1);
    for(int row=0;row<inCurve.GetSize()-1;row++)
    {
    double hj_cell;
    hj_cell=inCurve[row+1].v0-inCurve[row].v0; if(fabs(hj_cell)<0.0001) return FALSE;
    hj.SetAt(row,hj_cell);
    }

    //求alfa   个数为n-2
    CArray<double,double> alfa;
    alfa.SetSize(hj.GetSize()-1);
    for(int j=1;j<hj.GetSize();j++)
    {
    double alfa_cell;
    alfa_cell=hj.GetAt(j-1)/(hj.GetAt(j-1)+hj.GetAt(j));
    alfa.SetAt(j-1,alfa_cell);
    } //求beta   个数为n-2
    CArray<double,double> beta;
    beta.SetSize(inCurve.GetSize()-2);
    for(row=1;row<inCurve.GetSize()-1;row++)
    {
    int alfa_j=row-1;
    int hj_j=row; double beta_cell;
    beta_cell=3*((1-alfa.GetAt(alfa_j))*(inCurve[row].v1-inCurve[row-1].v1)/hj.GetAt(hj_j-1)+
    alfa.GetAt(alfa_j)*(inCurve[row+1].v1-inCurve[row].v1)/hj.GetAt(hj_j));
    beta.SetAt(row-1,beta_cell);
    } //求a   个数为n-1
    CArray<double,double> a;
    a.SetSize(alfa.GetSize()+1);
    a.SetAt(0,0);
    for(j=0;j<alfa.GetSize();j++)
    {
    double a_cell;
    a_cell=-alfa.GetAt(j)/(2+(1-alfa.GetAt(j))*a.GetAt(j));
    a.SetAt(j+1,a_cell);
    } //求b   个数为n-1
    CArray<double,double>b;
    b.SetSize(a.GetSize());
    b.SetAt(0,data1);
    for(j=0;j<a.GetSize()-1;j++)
    {
    double b_cell;
    b_cell=(beta.GetAt(j)-(1-alfa.GetAt(j))*b.GetAt(j))/(2+(1-alfa.GetAt(j))*a.GetAt(j));
    b.SetAt(j+1,b_cell);
    } //求yx_1 y对x的一阶导数
    CArray<double,double> yx_1;
    yx_1.SetSize(inCurve.GetSize());
    for(row=0;row<yx_1.GetSize();row++)  yx_1.SetAt(row,0);
    yx_1.SetAt(0,data1);
    yx_1.SetAt(yx_1.GetSize()-1,data2);
    for(row=inCurve.GetSize()-2;row>0;row--)
    {
    double yx_1_cell;
    yx_1_cell=a.GetAt(row)*yx_1.GetAt(row+1)+b.GetAt(row); yx_1.SetAt(row,yx_1_cell);
    } double dBeg=inCurve[0].v0;
    double dEnd=inCurve[inCurve.GetUpperBound()].v0; int outNum=int((dEnd-dBeg)/step);
    int outIndex(0);
    while(dBeg<dEnd){
    double x=dBeg; //内插
    //找到位置index
    if(x<inCurve[0].v0) return FALSE;
    if(x>inCurve[inCurve.GetUpperBound()].v0) return -1;

    int index=0;
    while(x>=inCurve[index].v0) 
    {
    index++;
    if(index==inCurve.GetUpperBound()) break;
    }
    index--; double x0=inCurve[index].v0;
    double x1=inCurve[index+1].v0;
    double y0=inCurve[index].v1;;
    double y1=inCurve[index+1].v1; double fdata1=(3*pow(x1-x,2)/pow(hj.GetAt(index),2)-2*pow(x1-x,3)/pow(hj.GetAt(index),3))*inCurve[index].v1;
    double fdata2=(3*pow(x-x0,2)/pow(hj.GetAt(index),2)-2*pow(x-x0,3)/pow(hj.GetAt(index),3))*inCurve[index+1].v1; double fdata3=((pow(x1-x,2)/pow(hj.GetAt(index),2)-pow(x1-x,3)/pow(hj.GetAt(index),3))*yx_1.GetAt(index));
    double fdata4=((pow(x-x0,2)/pow(hj.GetAt(index),2)-pow(x-x0,3)/pow(hj.GetAt(index),3))*yx_1.GetAt(index+1)); double data=fdata1+fdata2+fdata3-fdata4;
    //在输出曲线中写入数据
    tag_ZYCurve ctag;
    ctag.v0=x;
    ctag.v1=data;
    ctag.tag=-1;//-1表示不是输入曲线的值
    outCurve.Add(ctag); dBeg+=step;
    outIndex++;
    } //如果要输出曲线中保留输入曲线
    if(flag==1){
    //输入曲线===>输出曲线,
    //注意在输出曲线中可能已经包含了输入曲线的某些点
    int iIn(0),iOut(0);//两个序列的索引值 //两个序列都没有遍历到最后才有继续下去的必要
    while(iIn<=inCurve.GetUpperBound()  &&  iOut<=outCurve.GetUpperBound()){
    tag_ZYCurve cin=inCurve[iIn];
    tag_ZYCurve cout=outCurve[iOut]; //(cin.v0-cout.v0)的平方小于0.00001认为两者相等
    if((cin.v0-cout.v0)*(cin.v0-cout.v0)<0.0001){ //相等
    cout.tag=iIn;//表明此点是输入曲线上的点
    outCurve.SetAt(iOut,cout);//更新数据 iIn++;iOut++;//继续往后比较
    }
    else if(cin.v0<cout.v0){   //小于
    cin.tag=iIn;//将输入曲线的序号

    outCurve.InsertAt(iOut,cin);
    //插入后outCurve的长度增加了1,因此且在当前索引前加的,因此索引加1
    iOut++; iIn++;
    }
    else if(cin.v0>cout.v0){
    iOut++;
    }
    }
    //如果iOut已经到最后了,则把iIn剩下的数据直接插入outCurve的最后
    if(iOut>=outCurve.GetSize()){
    for(int i(iIn);i<inCurve.GetSize();i++){
    tag_ZYCurve cin=inCurve[i];
    cin.tag=iIn;
    outCurve.Add(cin);
    }
    }
    } return TRUE;
    }