python - 在将数据传递给插值例程之前存储数据

标签 python numpy scipy interpolation

首先我要道歉,因为我对 python 完全陌生(这里是 Fortran 人员),并且一直在即时学习。因此,我的知识中可能存在一些非常明显的漏洞,在阅读了我当前的困境后,这些漏洞可能会很明显。

我有一些数据需要写入文件,然后可以通过插值算法读取。目前,该插值算法可能是 SciPi 中的 RectBivariateSpline。虽然 X 和 Y 的间隔不一样,但它们是有规律的,所以这似乎是理想的。

数据一般是这样的形式

X1,Y1,F(X1,Y1)

X1,Y2,F(X1,Y2)

X1,Y3,F(X1,Y3)

X2,Y1,F(X2,Y1)

X2,Y2,F(X2,Y2)

X2,Y3,F(X2,Y3)

等等...

在这种情况下,F(X,Y) 不是一个显式数学函数,而是 X,Y 处物理量的数据点。

数据是从我基本上无法控制的数据源读取的,使用

 data_array = np.loadtxt(Path/DataFile, dtype = Float, delimiter = ";", usecols = #) 

在这种情况下,可以想象有几个不同的列具有不同的数据,但它们都依赖于 X 和 Y。有一个单独的文件,其中包含有关我读入的 X 和 Y 数量的范围和步长的信息并存储在数组中。至少,我相当确定它是一个数组而不是列表,因为我在某处读到 np.loadtxt 和 np.genfromtxt 都生成 numpy 数组而不是普通的 python 列表。

对于如何使用简单的插值例程以最佳方式存储这些数据以将其传输到另一台机器,我尝试过各种想法,但我可以使用一些建议。我首先想到使用

ArrayExample = np.Empty(xRange,yRange) 
For n in (xRange) 
    For m in (yRange) 
        ArrayExample[n,m] = F(X,Y)

但是,这对于保留与数组中的值关联的实际数量 X 或 Y 没有任何作用,而这些对于插值是必需的,并且对于绘图绝对是必需的。

然后我想到,既然我已经以一种可以很容易读取的形式获得了 X 和 Y 的值,我可以做如下的事情。其中 xvalues 和 yvalues 是保存实际 X 和 Y 值的数组。

ArrayExample = np.Empty(xRange,yRange,1) 
For n in (xRange) 
    For m in (yRange) 
        ArrayExample[n,m,1] = (xvalues[n],yvalues[m],F(X,Y))

之后我想到将 ArrayExample 保存为一个 pickle 文件,以便移动到任何地方,然后在另一端使用 pickle 引入。

但是,一旦我有了它,我真的不知道如何让 RectBivariateSpline 来获取数据。我尝试阅读 scipy 网站上的文档并进行谷歌搜索,但到目前为止我发现的所有内容都毫无帮助。如果有人有任何关于如何使用它的好例子,那将会非常有帮助。

如果您有任何建议、想法或批评,我们将不胜感激。

谢谢!

最佳答案

我之前看到过这个,希望有比我了解更多的人回复。 希望这能为您指明正确的方向

要使用 RectBivariateSpline 类,您需要 x 和 y 作为一维数组,z 值作为二维数组(len(x), len(y))

Numpy 需要任何特定的数组索引为整数,而不是 float 。 x 和 y 坐标需要转换为整数索引,以将 z 数据放入数组中。

import numpy as np
from scipy import interpolate

# Generate x and y arrays
x = np.linspace(45., 70., 70)
y = np.linspace(125/7, 59.876, 29)

# I don't know what your data looks like. Generate a list of dictionary records
data = []

for _x in x:
    for _y in y:
        data.append({'x': _x, 'y': _y, 'z': _x * _x - _y*_y/_x})

data[:100:10]
Out[4]: 
[{'x': 45.0, 'y': 17.857142857142858, 'z': 2017.9138321995465},
 {'x': 45.0, 'y': 32.86387755102041, 'z': 2000.9992344958118},
 {'x': 45.0, 'y': 47.870612244897956, 'z': 1974.075655184414},
 {'x': 45.36231884057971, 'y': 19.357816326530614, 'z': 2049.4792585649284},
 {'x': 45.36231884057971, 'y': 34.364551020408165, 'z': 2031.7068577153198},
 {'x': 45.36231884057971, 'y': 49.37128571428571, 'z': 2004.0054192006007},
 {'x': 45.72463768115942, 'y': 20.858489795918366, 'z': 2081.227345221296},
 {'x': 45.72463768115942, 'y': 35.86522448979592, 'z': 2062.610735570439},
 {'x': 45.72463768115942, 'y': 50.87195918367347, 'z': 2034.1437652565721},
 {'x': 46.08695652173913, 'y': 22.359163265306123, 'z': 2113.159976357177}]

def indexer(v):
    """ Returns a function to index into the array v by value.
        int( result + 0.5 ) to avoid any not quite equal with floats.
    """
    v_lo = v.min()
    v_range = v.max()-v_lo
    n = len(v)-1

    def func( x ):
        """ Returns an index from x. """
        return int(n*(x-v_lo)/v_range+0.5)
    return func

x_index = indexer(x)  # Function to index into the x values 
y_index = indexer(y)  # Function to index into the y values

z = np.empty((len(x), len(y)))

for rec in data:
    ix_x = x_index(rec['x'])  # Map the x value to an index
    ix_y = y_index(rec['y'])  # Map the y value to an index
    z[ix_x, ix_y] = rec['z']  # Place the z value at ix_x, ix_y

inter = interpolate.RectBivariateSpline( x, y, z)

inter(45.3, 18)  # Out[6]: array([[2044.93768211]])

inter(np.arange(45., 70.), 18)
Out[7]: 
   array([[2017.8       ], [2108.95652174], [2202.10638298], [2297.25      ],
          [2394.3877551 ], [2493.52      ], [2594.64705882], [2697.76923077],
          [2802.88679245], [2910.        ], [3019.10909091], [3130.21428571],
          [3243.31578947], [3358.4137931 ], [3475.50847458], [3594.6       ],
          [3715.68852459], [3838.77419355], [3963.85714286], [4090.9375    ],
          [4220.01538462], [4351.09090909], [4484.1641791 ], [4619.23529412],
          [4756.30434783]])

希望至少有一些指示/想法。

关于python - 在将数据传递给插值例程之前存储数据,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/55483827/

相关文章:

python - 在 python 中导入 ngram 时出错

python - 我应该使用 numpy 的随机生成器吗?

Python 积分

python - Scipy:如何将 KD-Tree 距离从查询转换为千米(Python/Pandas)

python - 比较 Pandas 中数据帧的标题

python - 用作外键时如何更改 Django Admin 中的用户表示?

python - 将上三角矩阵转换为对称矩阵的快速方法

python - numpy 数组中所有值的平方根,保留符号

python - 给 numpy 数组赋值

python - 大型二维 numpy 数组中相同元素的高效成对计算