首页 > 编程知识 正文

python克里金插值代码(python 克里金插值)

时间:2023-12-07 13:39:21 阅读:312841 作者:FVZI

本文目录一览:

  • 1、python 拉格朗日插值 不能超过多少个值
  • 2、01-02怎么用python实现克里金插值
  • 3、克里金插值算法
  • 4、双线性插值法原理 python实现

python 拉格朗日插值 不能超过多少个值

拉格朗日插值Python代码实现

1. 数学原理

对某个多项式函数有已知的k+1个点,假设任意两个不同的都互不相同,那么应用拉格朗日插值公式所得到的拉格朗日插值多项式为:

其中每个lj(x)为拉格朗日基本多项式(或称插值基函数),其表达式为:

2. 轻量级实现

利用

直接编写程序,可以直接插值,并且得到对应的函数值。但是不能得到系数,也不能对其进行各项运算。

123456789101112

def h(x,y,a):    ans=0.0    for i in range(len(y)):        t=y[i]        for j in range(len(y)):            if i !=j:                t*=(a-x[j])/(x[i]-x[j])        ans +=t    return ansx=[1,0]y=[0,2]print(h(x,y,2))

上述代码中,h(x,y,a)就是插值函数,直接调用就行。参数说明如下:

x,y分别是对应点的x值和y值。具体详解下解释。

a为想要取得的函数的值。

事实上,最简单的拉格朗日插值就是两点式得到的一条直线。

例如:

p点(1,0)q点(0,2)

这两个点决定了一条直线,所以当x=2的时候,y应该是-2

该代码就是利用这两个点插值,然后a作为x=2调用函数验证的。

3. 引用库

3.1 库的安装

主要依赖与 scipy。官方网站见:

安装的方法很简单,就是使用pip install scipy 如果失败,则将whl文件下载到本地再利用命令进行安装。

可能如果没有安装numpy

3.2 库的使用

from scipy.interplotate import lagrange

直接调用lagrange(x,y)这个函数即可,返回 一个对象。

参数x,y分别是对应各个点的x值和y值。

例如:(1,2) (3,5) (5,9)这三个点,作为函数输入应该这么写:

x=[1,3,5]

y =[2, 5, 9]

a=lagrange(x,y)

直接输出该对象,就能看到插值的函数。

利用该对象,能得到很多特性。具体参见:

a.order得到阶

a[]得到系数

a()得到对应函数值

此外可以对其进行加减乘除运算

3.3 代码实现

1234567   from scipy.interpolate import lagrangex=[1,2,3,4,7]y=[5,7,10,3,9]a=lagrange(x,y)print(a)print(a(1),a(2),a(3))print(a[0],a[2],a[3])   

结果是:

class 'numpy.lib.polynomial.poly1d' 4

4            3              2

0.5472 x - 7.306 x + 30.65 x - 47.03 x + 28.13

5.0 7.0 10.0

28.1333333333 30.6527777778 -7.30555555556

解释:

class 'numpy.lib.polynomial.poly1d' 4

这一行是输出a的类型,以及最高次幂。

4            3              2

0.5472 x - 7.306 x + 30.65 x - 47.03 x + 28.13

第二行和第三行就是插值的结果,显示出的函数。

第二行的数字是对应下午的x的幂,如果对应不齐,则是排版问题。

5.0 7.0 10.0

第四行是代入的x值,得到的结果。

也就是说,用小括号f(x)的这种形式,可以直接得到计算结果。

28.1333333333 30.6527777778 -7.30555555556

01-02怎么用python实现克里金插值

只为了输出的话倒是还好

1

print '['+', '.join(map(lambda x: '0'+str(x) if len(str(x))==1 else str(x),range(20)))+']'

克里金插值算法

根据项目对数据处理的要求,采用了优化的克里金插值算法,将等值线地化数据插值转换为格网数据,以便实现地化数据的三维显示(王家华等,1999)。其主要实现过程如下:

第一步,计算半变异图,用非线性最小二乘拟合半变异函数系数;

第二步,数据点进行四叉树存储;

第三步,对每一格网点搜索邻近数据点;

第四步,由待预测网格点和邻近数据点计算克里金算法中系数矩阵,及右端常数向量;

第五步,对矩阵进行LU分解,回代求解待预测点的预测值。

克里金插值算法主要包括半变异函数和邻近点搜索的计算,实现方法如下。

(1)半变异函数计算

半变异函数是地质统计学中区域化变量理论的基础。地质统计学主要完成2方面的任务:利用半变异函数生成半变异图来量化研究对象的空间结构;通过插值方法利用半变异图中拟合模型和研究对象周围的实测值来对未知值进行预测。

半变异函数是用来描述区域化变量结构性和随机性并存这一空间特征而提出的。在满足假设的条件下,随机函数z(x)和z(x+h)为某一物理参数测定值的一一对应的2组函数,h为每对数之间的距离。半变异函数γ(h)可用下式来计算:

γ(h)= 1/2E{[z(x)-z(x +h)]2}

4种基本的半变异函数模式(除了这4种基本模式以外,还有很多模式),包括:

1)线形模式(Linear Model)

浙江省农业地质环境GIS设计与实现

2)球面模式(Spherical Model)

浙江省农业地质环境GIS设计与实现

3)指数模式(Exponential Model)

浙江省农业地质环境GIS设计与实现

4)高斯模式(Gaussian Model)

浙江省农业地质环境GIS设计与实现

半变异函数γ(h)会随距离h增大而增大,并逐渐逼近一定值(C0 +C),称为基台值(Sill);而逼近基台值所对应的距离,称为影响范围(Range),表示空间中两位置间的距离小于影响范围时,是空间相关性的。在线性和球面模式中,影响范围等于a;在指数和高斯模式中,影响范围则分别等于3a和

。而模式于半变异函数轴的截距(C0)成为块金系数(Nugget Effect),产生的原因主要是样本测定的误差和最小采样间距内的变异。在应用上,为探讨说明空间变异在不同方向上的差距,也可利用非等向性的变异函数模式。半变异图拟合半变异函数模式的拟合方法可采用非线性最小二乘法拟合。

(2)邻近点搜索算法

由于矩阵LU分解求解方程的算法会随着矩阵维数的增加计算量增大,所以针对大量采样数据点时不能采用全部数据进行估计,必须采用插值点的临近点数据进行计算,即采用局部数据进行克里金算法进行计算。搜索邻近点可采用四叉树结构存储总数据,以提高搜索邻近点的速度。

对于选取邻近点的数目要有所限制,因该值的大小选择会影响插值的计算结果。若太大,则内插结果过于平滑;太小,则无法反映地表的变化;距离预测点较远的实测点可能与待估样点已经不存在自相关关系,也不能参与插值计算。采取以插值点为圆心,以R为半径的圆来确定取样的范围和参加计算的实测样点数目(如果存在各向异性,则可考虑划定一椭圆作为研究区域)。为了避免方向上的偏差,将圆平均地分为4个扇区,每个扇区内实测点数目在2~5之间,这样总共参与每个待估点预测的实测点数目平均达到8个。

区域内临近点的选择,存在着两种策略。

1)以邻近点的个数为基准。通常情况下,邻近点的个数以8~12个为宜,并且个数不能少于2个。此时计算出来的图像较为光滑。

2)以邻近点的半径尺度为基准。通常情况下,选择5~10 倍栅格间距的距离为宜。此时必须定义选择邻近点的最小和最大个数,当在一定半径内查找的邻近点个数小于最小个数时,应扩大搜索半径,使之达到最小查找个数;反之在一定半径内查找的邻近点个数大于最大个数时,应缩小搜索半径,使之小于最大查找个数。通常情况下最大最小个数分别可以定为20和4。

克里金算法的优点在于它基于一些可被验证的统计假设。根据这些假设,克里金算法产生的栅格节点估计量是最佳的,所有的估计量都依赖于可获得的观测值,并且平均误差最小。克里金算法提供了方差误差分析的表达式,可以表明每一个栅格节点的估计精度。

双线性插值法原理 python实现

码字不易,如果此文对你有所帮助,请帮忙点赞,感谢!

一. 双线性插值法原理:

        ① 何为线性插值?

        插值就是在两个数之间插入一个数,线性插值原理图如下:

        ② 各种插值法:

        插值法的第一步都是相同的,计算目标图(dstImage)的坐标点对应原图(srcImage)中哪个坐标点来填充,计算公式为:

        srcX = dstX * (srcWidth/dstWidth)

        srcY = dstY * (srcHeight/dstHeight)

        (dstX,dstY)表示目标图像的某个坐标点,(srcX,srcY)表示与之对应的原图像的坐标点。srcWidth/dstWidth 和 srcHeight/dstHeight 分别表示宽和高的放缩比。

        那么问题来了,通过这个公式算出来的 srcX, scrY 有可能是小数,但是原图像坐标点是不存在小数的,都是整数,得想办法把它转换成整数才行。

        不同插值法的区别就体现在 srcX, scrY 是小数时,怎么将其变成整数去取原图像中的像素值。

        最近邻插值(Nearest-neighborInterpolation):看名字就很直白,四舍五入选取最接近的整数。这样的做法会导致像素变化不连续,在目标图像中产生锯齿边缘。

        双线性插值(Bilinear Interpolation):双线性就是利用与坐标轴平行的两条直线去把小数坐标分解到相邻的四个整数坐标点。权重与距离成反比。

        双三次插值(Bicubic Interpolation):与双线性插值类似,只不过用了相邻的16个点。但是需要注意的是,前面两种方法能保证两个方向的坐标权重和为1,但是双三次插值不能保证这点,所以可能出现像素值越界的情况,需要截断。

        ③ 双线性插值算法原理

        假如我们想得到未知函数 f 在点 P = (x, y) 的值,假设我们已知函数 f 在 Q11 = (x1, y1)、Q12 = (x1, y2), Q21 = (x2, y1) 以及 Q22 = (x2, y2) 四个点的值。最常见的情况,f就是一个像素点的像素值。首先在 x 方向进行线性插值,然后再在 y 方向上进行线性插值,最终得到双线性插值的结果。

    ④ 举例说明

二. python实现灰度图像双线性插值算法:

灰度图像双线性插值放大缩小

import numpy as np

import math

import cv2

def double_linear(input_signal, zoom_multiples):

    '''

    双线性插值

    :param input_signal: 输入图像

    :param zoom_multiples: 放大倍数

    :return: 双线性插值后的图像

    '''

    input_signal_cp = np.copy(input_signal)  # 输入图像的副本

    input_row, input_col = input_signal_cp.shape # 输入图像的尺寸(行、列)

    # 输出图像的尺寸

    output_row = int(input_row * zoom_multiples)

    output_col = int(input_col * zoom_multiples)

    output_signal = np.zeros((output_row, output_col)) # 输出图片

    for i in range(output_row):

        for j in range(output_col):

            # 输出图片中坐标 (i,j)对应至输入图片中的最近的四个点点(x1,y1)(x2, y2),(x3, y3),(x4,y4)的均值

            temp_x = i / output_row * input_row

            temp_y = j / output_col * input_col

            x1 = int(temp_x)

            y1 = int(temp_y)

            x2 = x1

            y2 = y1 + 1

            x3 = x1 + 1

            y3 = y1

            x4 = x1 + 1

            y4 = y1 + 1

            u = temp_x - x1

            v = temp_y - y1

            # 防止越界

            if x4 = input_row:

                x4 = input_row - 1

                x2 = x4

                x1 = x4 - 1

                x3 = x4 - 1

            if y4 = input_col:

                y4 = input_col - 1

                y3 = y4

                y1 = y4 - 1

                y2 = y4 - 1

            # 插值

            output_signal[i, j] = (1-u)*(1-v)*int(input_signal_cp[x1, y1]) + (1-u)*v*int(input_signal_cp[x2, y2]) + u*(1-v)*int(input_signal_cp[x3, y3]) + u*v*int(input_signal_cp[x4, y4])

    return output_signal

# Read image

img = cv2.imread("../paojie_g.jpg",0).astype(np.float)

out = double_linear(img,2).astype(np.uint8)

# Save result

cv2.imshow("result", out)

cv2.imwrite("out.jpg", out)

cv2.waitKey(0)

cv2.destroyAllWindows()

三. 灰度图像双线性插值实验结果:

四. 彩色图像双线性插值python实现

def BiLinear_interpolation(img,dstH,dstW):

    scrH,scrW,_=img.shape

    img=np.pad(img,((0,1),(0,1),(0,0)),'constant')

    retimg=np.zeros((dstH,dstW,3),dtype=np.uint8)

    for i in range(dstH-1):

        for j in range(dstW-1):

            scrx=(i+1)*(scrH/dstH)

            scry=(j+1)*(scrW/dstW)

            x=math.floor(scrx)

            y=math.floor(scry)

            u=scrx-x

            v=scry-y

            retimg[i,j]=(1-u)*(1-v)*img[x,y]+u*(1-v)*img[x+1,y]+(1-u)*v*img[x,y+1]+u*v*img[x+1,y+1]

    return retimg

im_path='../paojie.jpg'

image=np.array(Image.open(im_path))

image2=BiLinear_interpolation(image,image.shape[0]*2,image.shape[1]*2)

image2=Image.fromarray(image2.astype('uint8')).convert('RGB')

image2.save('3.png')

五. 彩色图像双线性插值实验结果:

六. 最近邻插值算法和双三次插值算法可参考:

        ① 最近邻插值算法:

        

        ② 双三次插值算法:

        

七. 参考内容:

        

        

版权声明:该文观点仅代表作者本人。处理文章:请发送邮件至 三1五14八八95#扣扣.com 举报,一经查实,本站将立刻删除。