Преобразование быстрого алгоритма MCLT в Matlab в Java

Question

У меня есть Java-код, который дает мне FFT-вывод от реальных входных данных. Мне нужно выполнить MCLT. В настоящее время у меня есть выход FFT со мной в следующем формате. Я видел какой-то быстрый аллогрит MCLT (https://www.microsoft.com/en-us/research/wp-content/uploads/2016/02/tr-2005-02.pdf),, закодированный в Matlab, но не могу понять это идеально. Может кто-нибудь помочь мне в написании соответствующего кода Java.

Java Code Начальная точка:

int dtLength =  data.length/2;
double[] realPart = new double[dtLength];
double[] imagPart = new double[dtLength];

Код Matlab:

function X = fmclt(x)
% FMCLT - Compute MCLT of a vector via double-length FFT
%
% H. Malvar, September 2001 -- (c) 1998-2001 Microsoft Corp.
%
% Syntax: X = fmclt(x)
%
% Input: x : real-valued input vector of length 2*M
%
% Output: X : complex-valued MCLT coefficients, M subbands
% in Matlab, by default j = sqrt(-1)
% determine # of subbands, M
L = length(x);
M = L/2;
% normalized FFT of input
U = sqrt(1/(2*M)) * fft(x);
% compute modulation function
k = [0:M]';
c = W(8,2*k+1) .* W(4*M,k);
% modulate U into V
V = c .* U(1:M+1);
% compute MCLT coefficients
X = j * V(1:M) + V(2:M+1);
return;
% Local function: complex exponential
function w = W(M,r)
w = exp(-j*2*pi*r/M);
return; 

Answer 1

Несмотря на то, что этот вопрос является своего рода пограничным для SO, статья была довольно интересной, поэтому я решил потратить некоторое время, читая ее и пытаясь преобразовать код Matlab в Java. Вот результат:

import org.apache.commons.math3.complex.Complex;

public class MCLT
{
    public static void main(String args[])
    {
        Complex[] x = new Complex[16];

        for (int i = 1; i <= 16; ++i)
            x[(i - 1)] = new Complex((double)i, 0.0d);

        Complex[] result = fmclt(x);

        for (int i = 0; i < result.length; ++i)
            System.out.println(result[i]);
    }

    public static Complex[] fmclt(Complex[] x)
    {
        int L = x.length;
        int M = L / 2;

        double z = Math.sqrt(1.0d / (2.0d * M));

        Complex[] F = fft(x);
        Complex[] U = new Complex[F.length];

        for (int i = 0; i < F.length; ++i)
            U[i] = F[i].multiply(z);

        double[] k = new double[(M + 1)];

        for (int i = 0; i <= M; ++i)
            k[i] = (double)i;

        Complex[] c = new Complex[(M + 1)];

        for (int i = 0; i <= M; ++i)
            c[i] = W(8.0d, ((2.0d * k[i]) + 1.0d)).multiply(W((4.0d * M), k[i]));

        Complex[][] V = new Complex[(M + 1)][];

        for (int i = 0; i <= M; ++i)
        {
            V[i] = new Complex[(M + 1)];

            for (int j = 0; j <= M; ++j)
                V[i][j] = c[i].multiply(U[j]);
        }

        Complex[] V1 = new Complex[M];

        for (int i = 0; i < M; ++i)
            V1[i] = V[i][0];

        Complex[] V2 = new Complex[M];

        for (int i = 1; i <= M; ++i)
            V2[(i - 1)] = V[i][0];

        Complex b = new Complex(0.0d, 1.0d);
        Complex[] result = new Complex[M];

        for (int i = 0; i < M; ++i)
            result[i] = b.multiply(V1[i]).add(V2[i]); 

        return result;
    }

    public static Complex[] fft(Complex[] x)
    {
        int n = x.length;

        if (n == 1)
            return new Complex[] { x[0] };

        if ((n % 2) != 0)
            throw new IllegalArgumentException("Invalid length.");

        int nh = n / 2;

        Complex[] even = new Complex[nh];

        for (int i = 0; i < nh; ++i)
            even[i] = x[(2 * i)];

        Complex[] q = fft(even);

        Complex[] odd  = even;

        for (int i = 0; i < nh; ++i)
            odd[i] = x[((2 * i) + 1)];

        Complex[] r = fft(odd);

        Complex[] y = new Complex[n];

        for (int i = 0; i < nh; ++i)
        {
            double kth = -2.0d * i * (Math.PI / n);
            Complex wk = new Complex(Math.cos(kth), Math.sin(kth));

            y[i] = q[i].add(wk.multiply(r[i]));
            y[(i + nh)] = q[i].subtract(wk.multiply(r[i]));
        }

        return y;
    }

    public static Complex W(double M, double r)
    {
        Complex j = (new Complex(0.0d, 1.0d)).multiply(-1.0d);
        double z = 2.0d * Math.PI * (r / M);

        return j.multiply(z).exp();
    }
}

Использование отдельных двойных массивов для реальных и мнимых частей, на мой взгляд, не было хорошим выбором для проектирования, поэтому я решил основать свой код на Complex классе Apache Commons вместо этого.

Чтобы рассчитать быстрое преобразование Фурье, я решил использовать некоторый готовый код. Моя fft функция основана на этой реализации , которая кажется очень надежной и использует вышеупомянутый класс Complex.

Используя один и тот же вектор значений, оба кода Matlab и Java возвращают один и тот же вывод. Вы можете протестировать код онлайн, скопировав его на на этом сайте , но вам также необходимо установить библиотеку Apache Commons, прежде чем вы сможете успешно запустить его. Нажмите кнопку Add External Library (from Maven Repo), расположенную внизу, а затем вставьте следующие параметры в форму ввода:

<!-- https://mvnrepository.com/artifact/org.apache.commons/commons-math3 -->
<dependency>
    <groupId>org.apache.commons</groupId>
    <artifactId>commons-math3</artifactId>
    <version>3.6.1</version>
</dependency>