python - 检测阵列的抛物面形状部分

Question

我正在尝试检测 numpy.ndarray 内的抛物面形状点.

我有一个数据集，使用维纳滤波器对其进行平滑处理(使用 NIFTy4 包):它保存为 .csv 并且可以在 here on Google Drive 中找到。或here on PasteBin .

我正在寻找一种方法来识别数组中描述向上打开的抛物面的部分。

在提供的示例中，我想检测 5 个这样的形状。重点主要在于形状的宽度(开始和结束)而不是实际的最小位置。

提前致谢

Answer 1

使用numpy的解决方案

对于这个解决方案，我将使用numpy:

import numpy as np

数据集

创建数据集

OP Phteven我很友善地提供了要使用的数据集，但由于与数据集的链接往往会消失，我还创建了一个函数来生成类似曲线。

def polyval_points(n):
    """Return random points generated from polyval."""
    x = np.linspace(-20, 20, num=n)
    coef = np.random.normal(-3, 3, size=(5))
    y = np.polyval(coef, x) * np.sin(x)
    return x, y

加载数据集

def dataset_points():
    """Return points loaded from dataset."""
    y = np.loadtxt("data.csv", delimiter=',')
    x = np.linspace(0, 1, num=len(y))
    return x, y

阐述要点

斜率卷积

由于点是离散的，我们必须表示斜率。其中一种方法是通过统一内核。

def convolute_slopes(y, k=3):
    """Return slopes convoluted with an uniform kernel of size k."""
    d2y = np.gradient(np.gradient(y))
    return np.convolve(d2y, np.ones(k)/k, mode="same")

获取抛物面

现在我们可以计算卷积斜率，确定它在哪里切换方向，并选择平均斜率大于绝对平均值乘以描述抛物面“斜率”程度的系数的间隔。

def get_paraboloids(x, y, c=0.2):
    """Return list of points (x,y) that are part of paraboloids in given set.
        x: np.ndarray of floats
        y: np.ndarray of floats
        c: slopyness coefficient
    """
    slopes = convolute_slopes(y)
    mean = np.mean(np.abs(slopes))
    w = np.where(np.diff(slopes > 0) > 0)[0] + 1
    w = np.insert(w, [0, len(w)], [0, len(x)])
    return [(x[lower:upper], y[lower:upper])
            for lower, upper in zip(w[:-1], w[1:])
            if np.mean(slopes[lower:upper]) > mean * c]

如何使用它和可视化

首先我们加载数据集并生成更多数据:

datasets = [dataset_points(), polyval_points(10000), polyval_points(10000)]

然后，迭代每个数据集:

from matplotlib import pyplot as plt

plt.figure(figsize=(10 * len(datasets), 10))
for i, points in enumerate(datasets):
    x, y = points
    plt.subplot(1, len(datasets), i + 1)
    plt.plot(x, y, linewidth=1)
    for gx, gy in get_paraboloids(x, y):
        plt.plot(gx, gy, linewidth=3)
plt.show()

结果