xarray 反向插值（在坐标上，而不是在数据上）

Question

我有以下 DataArray

arr = xr.DataArray([[0.33, 0.25],[0.55, 0.60],[0.85, 0.71],[0.92,0.85],[1.50,0.96],[2.5,1.1]],[('x',[0.25,0.5,0.75,1.0,1.25,1.5]),('y',[1,2])])

这给出了以下输出

<xarray.DataArray (x: 6, y: 2)>
array([[0.33, 0.25],
       [0.55, 0.6 ],
       [0.85, 0.71],
       [0.92, 0.85],
       [1.5 , 0.96],
       [2.5 , 1.1 ]])
Coordinates:
  * x        (x) float64 0.25 0.5 0.75 1.0 1.25 1.5
  * y        (y) int32 1 2

或为方便起见，将 x 和输出 (z) 并排排列在下方。

x         z (y=1)   z(y=2)
0.25      0.33      0.25
0.50      0.55      0.60
0.75      0.85      0.71
1.00      0.92      0.85
1.25      1.50      0.96
1.50      2.50      1.10

我的数据是几个输入值的结果。其中之一是 x 值。其他输入值还有其他几个维度（例如 y）。我想知道我的输出值 (z) 何时增长大于 1.00，保持其他维度不变并改变 x 值。在上面的二维例子中，我想得到答案 [1.03 1.32]。因为当 y=1 时 x 的值为 1.03 时 z 的值为 1.00，而当 y=2 时 x 的值为 1.32 时 z 的值为 1.00。

编辑：由于输出 z 将随着 x 的增加而增长，因此只有一个点 z 的输出为 1.0。

有没有什么有效的方法可以用 xarray 实现这个？我的实际 table 大得多，有 4 个输入（维度）。

感谢您的帮助！

Answer 1

xarray 有一个非常方便的函数：xr.interp 它将对 xarray 进行分段线性插值。

在您的情况下，您可以使用它来获得 (x, y1) 和 (x, y1) 点的分段插值。完成此操作后，唯一剩下要做的就是将与内插 y1/y2/.. 数组的关闭值相关联的内插 x 数组的值获取到目标数字（在您的示例中为 1.00 ).

这可能是这样的：

y_dims = [0, 1,] 
target_value = 1.0
# create a 'high resolution` version of your data array:
arr_itp = arr.interp(x=np.linspace(arr.x.min(), arr.x.max(), 10000))
for y in y_dims:
    # get the index of closest data
    x_closest = np.abs(arr_itp.isel(y=y) - target_value).argmin()
    print(arr_itp.isel(y=y, x=x_closest))

>>> <xarray.DataArray ()>
>>> array(0.99993199)
>>> Coordinates:
>>>     y        int64 1
>>>     x        float64 1.034
>>> <xarray.DataArray ()>
>>> array(1.00003)
>>> Coordinates:
>>>     y        int64 2
>>>     x        float64 1.321

虽然这可行，但它并不是解决问题的真正有效方法，这里有 2 个不可行的原因：

使用xr.interp 对整个DataArray 进行分段插值。然而，我们只需要最接近目标值的两点之间的插值。
在这里，插值是两点之间的直线。但是，如果我们知道该直线上一点的一个坐标 (y=1.00)，那么我们可以通过求解直线的线性方程来简单地计算出另一个坐标，并通过一些算术运算解决问题。

考虑到这些原因，我们可以针对您的问题制定更有效的解决方案：

# solution of linear function between two points (2. reason)
def lin_itp(p1,p2,tv):
    """Get x coord of point on line

    Determine the x coord. of a point (x, target_value) on the line
    through the points p1, p2.

    Approach:
      - parametrize x, y between p1 and p2: 
          x = p1[0] + t*(p2[0]-p1[0])
          y = p1[1] + t*(p2[1]-p1[1])
      - set y = tv and resolve 2nd eqt for t
          t = (tv - p1[1]) / (p2[1] - p1[1])
      - replace t in 1st eqt with solution for t
          x = p1[0] + (tv - p1[1])*(p2[0] - p1[0])/(p2[1] - p1[1])
    """
    return float(p1[0] + (tv - p1[1])*(p2[0] - p1[0])/(p2[1] - p1[1])) 

# target value:
t_v = 1.0
for y in [0, 1]:
    arr_sd = arr.isel(y=y)
    # get index for the value closest to the target value (but smaller)
    s_udim = int(xr.where(arr_sd - t_v <=0, arr_sd, arr_sd.min()).argmax())
    # I'm explicitly defining the two points here
    ps_itp = arr_sd[s_udim:s_udim+2]
    p1, p2 = (ps_itp.x[0], ps_itp[0]), (ps_itp.x[1], ps_itp[1])
    print(lin_itp(p1,p2,t_v))

>>> 1.0344827586206897
>>> 1.3214285714285714

Answer 2

jojo的答案我遇到的问题是很难在多维度上扩展它并保持xarray结构。因此，我决定进一步研究这个问题。我使用了 jojo 代码中的一些想法来做出以下回答。

我制作了两个数组，一个条件是值小于我寻找的值，另一个条件是它们需要更大。我将第二个在 x 方向上移动负 1。现在我将它们组合成一个正常的线性插值公式。这两个数组仅在条件的 'edge' 处具有重叠的值。如果不移动 -1，则不会有任何值重叠。在最后一行中，我对 x 方向求和，因为所有其他值都是 NaN，我提取正确的值并在此过程中从 DataArray 中删除 x 方向。

def interpolate_dimension_x(arr, target_value, step):
    M0 = arr.where(arr - target_value <= 0)
    M1 = arr.where(arr - target_value > 0).shift(x=-1)

    work_mat = M0.x + step * (target_value - M0) / (M1 - M0)

    return work_mat.sum(dim='x')

interpolate_dimension_x(arr, 1, 0.25)

>>> <xarray.DataArray (y: 2)>
array([1.034483, 1.321429])
Coordinates:
  * y        (y) int32 1 2

我的代码有一些缺点。该代码仅在 M0 和 M1 找到满足条件的值时才有效。否则，该行中的所有值都将设置为 NaN。为了防止 M0 出现问题，我决定让 x 值从 0 开始，因为我的目标值总是大于 0。为了防止 M1 出现问题，我选择足够大的 x 值，这样我就知道我的值在那里.当然，这些都不是理想的解决方案，可能会破坏代码。如果我对 xarray 和 python 有更多的经验，我可能会重写。总之，我有以下要解决的问题：

如何推断 x 范围外的值？我目前只是确保我的 x 范围足够大，以便答案落在其中。
如何使代码对可变步长健壮？
如何编写代码以便动态选择我的维度（现在它仅适用于 'x'）
感谢任何优化。

xarray 反向插值（在坐标上，而不是在数据上）

xarray reverse interpolation (on coordinate, not on data)

python

interpolation

python-xarray