zoukankan      html  css  js  c++  java
  • 拉格朗日插值编程实现

    拉格朗日插值原理:

    拉格朗日插值的具体介绍网址:https://zh.wikipedia.org/wiki/%E6%8B%89%E6%A0%BC%E6%9C%97%E6%97%A5%E6%8F%92%E5%80%BC%E6%B3%95

    翻译成人话就是,该曲线是由多个n次多项式的和构成的,n是参与插值的点的个数。每个n次多项式的计算方法如上图所示。转化成程序的话,就是要保存每个多项式中(x-xi)中的每一项xi。然后就是系数yi/(xj-x0)(xj-x1)...(xj-xk)

    程序实现:

    #include <iostream>
    #include <opencv2/opencv.hpp>
    
    using namespace cv;
    
    template<class T>
    struct Point
    {
        T x;
        T y;
        
        Point()
        {x=0;y=0;}
        
        Point(T a,T b)
        {x=a;y=b;}
        
        Point operator+(Point p)
        {return Point(x+p.x,y+p.y);}
    };
    /*
     拉格朗日曲线,包含多个拉格朗日多项式和曲线的x坐标最大最小值
     */
    struct LagrangeCurve
    {
        vector<double>factor;
        vector<double>minuend;
        int minX;
        int maxX;
    };
    
    void calLagrangePolynomial(LagrangeCurve&curve,Point<int>*curvePoints,int pointNum)
    {
        int i,j;
        curve.maxX=curvePoints[0].x;
        curve.minX=curvePoints[0].x;
        for(i=0;i<pointNum;++i)
        {
            curve.minuend.push_back(curvePoints[i].x);
            curve.factor.push_back(curvePoints[i].y);
            //获取曲线两端的x值
            if(curvePoints[i].x<curve.minX)
                curve.minX=curvePoints[i].x;
            if(curvePoints[i].x>curve.maxX)
                curve.maxX=curvePoints[i].x;
            //计算(xi-x0)*...(xi-x5),并且把x0...x5保存
            for(j=0;j<pointNum;++j)
            {
                if(j!=i)
                {curve.factor[i]/=(curvePoints[i].x-curvePoints[j].x+1e-8);}
            }
        }
    }
    
    double getInterpolation(int x,LagrangeCurve& lagrangeCurve)
    {
        int i,j;
        double y=0;
        for(i=0;i<lagrangeCurve.minuend.size();++i)
        {
            double temp=lagrangeCurve.factor[i];
            for(j=0;j<lagrangeCurve.minuend.size();++j)
            {
                if(j!=i)
                {temp*=(x-lagrangeCurve.minuend[j]);}
            }
            y+=temp;
        }
        return y;
    }
    
    
    int main(int argc, const char * argv[])
    {
        Mat img(256,256,CV_8UC3);
        img=Scalar::all(255);
        Point<int>*points=new Point<int>[3];
        points[0].x=10;points[0].y=30;
        points[1].x=100;points[1].y=70;;
        points[2].x=250;points[2].y=210;
        
        circle(img, cvPoint(points[0].x,points[0].y),5,Scalar::all(0),-1);
        circle(img, cvPoint(points[1].x,points[1].y),5,Scalar::all(0),-1);
        circle(img, cvPoint(points[2].x,points[2].y),5,Scalar::all(0),-1);
        
        LagrangeCurve lagrangeCurve;
        calLagrangePolynomial(lagrangeCurve,points,3);
        for(int i=lagrangeCurve.minX;i<lagrangeCurve.maxX;++i)
        {
            double y=getInterpolation(i,lagrangeCurve);
            img.at<Vec3b>(static_cast<int>(y+0.5),i)=Vec3b(0,0,0);
        }
        imshow("img", img);
        waitKey(-1);
        return 0;
    }

    上述代码使用到了opencv画图,运行结果图如下:

  • 相关阅读:
    [转]人生以快乐为本
    不用iTunes也能添加音乐到iPod
    设计很有意思的U盘
    PhotoFunia 在线生成趣味图片
    [转]关于项目管理的一点杂感
    MVC视频序列和Demo的下载地址
    视频测试序列的下载地址
    fatal error LNK1104: 无法打开文件“LIBC.lib”错误
    ORACLE数据库性能优化概述
    oracle常用hint
  • 原文地址:https://www.cnblogs.com/hrlnw/p/5016700.html
Copyright © 2011-2022 走看看