信号压缩矢量化，python源代码

自变量 n 控制点的数量，因变量误差 log10(loss) 

压缩和解压效果图，绿线为原始信号，蓝点为压缩后的控制点，红线为解压出来的信号

源代码：

import numpy as np
from matplotlib import pyplot as plt


def init_N(n, k, X, T):
    """
    生成初始化的参数矩阵
    :param n: 控制点数量，int
    :param k: 基函数的阶数，int, 默认为2，建议2,3,4,最小为2，越大越平滑，越大速度越慢
    :param X: 自变量矩阵，取值范围[0,1] (shape=[sample], dtype=float)
    :param T: 控制点系数矩阵，取值范围[0,1] (shape=[n + k + 1], dtype=float)
    :return N:  (shape=[sample, n], dtype=float)
    """ 
    X = X[None]   # shape=[1, sample]
    T = T[:, None]# shape=[n+k+1, 1]
    m = n + k + 1
    N = 1 * (T[:-1] <= X) * (X < T[1:])  # shape=[n+k+1, sample]
    for k in range(1, k + 1):  # 动态规划 k ∈ [1， k]
        n = m - k - 1  # 每次更新的控制点数量
        c1 = T[k:n + k] - T[:n]
        c2 = T[k + 1:n + k + 1] - T[1:n + 1,]
        c1[c1 == 0] = 1
        c2[c2 == 0] = 1
        w1 = (X - T[:n]) / c1
        w2 = (T[k + 1:n + k + 1] - X) / c2
        N = N[:-1] * w1 + N[1:] * w2
    N = N[:n].T # 截取前n个，然后交换行列
    return N  # sample, n


def init_T_X(k, m, sample):
    """
    生成初始化的参数向量
    :param k: 基函数的阶数，int, 默认为2，建议2,3,4,最小为2，越大越平滑，越大速度越慢
    :param m: 系数矩阵的长度，int，m=n+k+1
    :param sample: 采样点数量， int
    :return T: 系数矩阵，取值范围[0,1] (shape[n + k + 1], dtype=float)
    :return X: 自变量矩阵，取值范围[0,1] (shape=[sample], dtype=float)
    """
    # clame 模式
    T = np.concatenate([
        np.zeros(k),
        np.linspace(0, 1, m - 2 * k),
        np.ones(k)]
    )
    gama = 1e-9
    first = T[k]
    last = T[-k - 1] * (1 - gama) + gama * T[-k - 2]
    X = np.linspace(first, last, sample)
    return T, X


def unzip(P, k, sample):
    """
    解压点，生成曲线
    :param P: 控制点矩阵，(shape=[n], dtype=float)
    :param k: 基函数的阶数，int, 默认为2，建议2,3,4,最小为2，越大越平滑，越大速度越慢
    :param sample: 采样点数量， int
    :return data: 解压后的曲线，(shape=[sample], dtype=float)
    """
    n = len(P)
    T, X = init_T_X(k=k, m=k + n + 1, sample=sample)
    N = init_N(n=n, k=k, T=T, X=X)
    numerator = np.sum(P[None] * N, axis=-1)
    denominator = np.sum(N, axis=-1)
    data = numerator / denominator
    return data


def zip(data, n, k):
    """
    压缩曲线，生成控制点矩阵
    :param data: 曲线上的采样点，(shape=[sample], dtype=float)
    :param n: 控制点的数量，int
    :param k: 曲线的阶数，int
    :return P: 控制点矩阵， （shape=[n], dtype=float）
    """
    T, X = init_T_X(k=k, m=k + n + 1, sample=sample)
    N = init_N(n=n, k=k, T=T, X=X)
    # print(N)
    a = N  # r 行 c列， 其中 r 代表采样点数量，c 代表控制点数量
    b = data * np.sum(N, axis=-1)  # r 行
    P = np.linalg.lstsq(a, b, rcond=None)
    # print(P)
    return P[0],P[1]


if __name__ == "__main__":
    sample = 720 # 采样数量
    X = np.linspace(0, 1, sample)
    Y = np.sin(7 * np.pi * X + 6) * (0.2 + np.abs(X - 0.8)) + \
        np.sin(9 * np.pi * X + 9) * (0.2 + np.abs(X - 0.1))  + \
        np.sin(32 * np.pi * X + 1) * (0.2 + np.abs(X - 0.6))         #  np.sin(17 * np.pi * X + 1) * (0.2 + np.abs(X - 0.1)) + \

    k = 2 # 阶数
    
    y_loss = []
    x_loss = list(range(36,128))
    for n in x_loss:
        P,loss = zip(Y, n, k) # 压缩提取控制点矩阵
        y_loss.append( np.log10(loss))
    plt.plot(x_loss, y_loss, 'g-', marker='+', label='') 
    plt.show()
    n = 32*2 + 2
    P,loss = zip(Y, n, k) # 压缩提取控制点矩阵

    unzip_Y = unzip(P, k, sample) # 解压 

    plt.plot(X, Y, 'g-', marker='+', label='origin line')
    plt.plot(np.linspace(0, 1, n), P, 'bo-.', label='zip point')
    plt.plot(X, unzip_Y, 'r-', marker='+', label='unzip line')

    plt.legend()
    plt.show()