narginchk(1,2) if isempty(x) error(message('signal:lpc:Empty')); end [m,n] = size(x); if (n>1) && (m==1) x = x(:); [m,n] = size(x); end if nargin < 2 N = m-1; elseif N < 0 % Check for N positive error(message('signal:lpc:negativeOrder')); end % Check the input data type. Single precision is not supported. try chkinputdatatype(x,N); catch ME throwAsCaller(ME); end if (N > m) error(message('signal:lpc:orderTooLarge', 'X must be a vector with length greater or equal to the prediction order.', 'If X is a matrix, the length of each column must be greater or equal to', 'the prediction order.')); end % Compute autocorrelation vector or matrix X = fft(x,2^nextpow2(2*size(x,1)-1)); R = ifft(abs(X).^2); R = R./m; % Biased autocorrelation estimate [a,e] = levinson(R,N); % Return only real coefficients for the predictor if the input is real for k = 1:n if isreal(x(:,k)) a(k,:) = real(a(k,:)); end end 将以上代码改写为c#

以下是将Matlab代码改写为C#的版本：

void Lpc(double[] x, int N, out double[] a, out double e)
{
    int m = x.Length;
    int n = 1;

    if (m == 0)
        throw new ArgumentException("Input vector cannot be empty.");

    if (N < 0)
        throw new ArgumentException("Prediction order must be a non-negative integer.");

    if (N > m)
        throw new ArgumentException("Prediction order cannot be greater than the length of input vector.");

    if (N >= m)
        throw new ArgumentException("Prediction order must be less than the length of input vector.");

    if (N == 0)
        N = m - 1;

    if (N == 1)
    {
        a = new double[] { -x[1] / x[0] };
        e = x[0] + x[1] * a[0];
        return;
    }

    if (n > 1 && m == 1)
    {
        n = m;
        x = new double[] { x[0] };
    }

    double[] X = new double[2 * m - 1];
    for (int i = 0; i < m; i++)
        X[i] = x[i];

    for (int i = m; i < X.Length; i++)
        X[i] = 0;

    Fourier.Forward(X, FourierOptions.Matlab);

    double[] R = new double[X.Length];

    for (int i = 0; i < X.Length; i++)
        R[i] = X[i] * X[i];

    Fourier.Inverse(R, FourierOptions.Matlab);

    e = R[0] / m;

    a = Levinson(R, N);

    for (int k = 0; k < n; k++)
    {
        if (IsReal(x[k]))
            a[k] = Real(a[k]);
    }
}

double[] Levinson(double[] R, int N)
{
    int p = R.Length - 1;
    double[] a = new double[N + 1];
    double[] k = new double[N + 1];
    double[] alpha = new double[N + 1];
    double[] beta = new double[N + 1];

    a[0] = 1;
    k[1] = -R[1] / R[0];
    alpha[1] = R[0] + k[1] * R[1];
    beta[1] = -k[1];

    for (int i = 1; i < N; i++)
    {
        double sum = 0;
        for (int j = 0; j <= i; j++)
            sum += R[i + 1 - j] * a[j];

        k[i + 1] = -sum / alpha[i];
        a[i + 1] = k[i + 1];

        for (int j = 1; j <= i; j++)
            a[j] += k[i + 1] * a[i + 1 - j];

        alpha[i + 1] = alpha[i] * (1 - k[i + 1] * k[i + 1]);

        for (int j = 1; j <= i; j++)
            beta[j] = a[i + 1 - j];

        for (int j = 1; j <= i; j++)
            a[j] = beta[j] + k[i + 1] * beta[i + 1 - j];

        beta[i + 1] = a[0];
    }

    return a;
}

bool IsReal(double[] x)
{
    for (int i = 0; i < x.Length; i++)
        if (x[i] != Math.Round(x[i]))
            return false;

    return true;
}

double Real(double x)
{
    return Math.Round(x);
}