FloatFFT_1D
Computes 1D Discrete Fourier Transform (DFT) of complex and real, single
PI
TWO_PI
bk1
bk1l
bk2
bk2l
factors
ip
ipl
n
nBluestein
nBluesteinl
nc
ncl
nl
nw
nwl
plan
useLargeArrays
w
wl
wtable
wtable_r
wtable_rl
wtablel
bluestein_complex()
bluestein_complex(pl.edu.icm.jlargearrays.FloatLargeArray, long, int)
bluestein_complex(float[], int, int)
bluestein_real_forward()
bluestein_real_forward(pl.edu.icm.jlargearrays.FloatLargeArray, long)
bluestein_real_forward(float[], int)
bluestein_real_full()
bluestein_real_full(pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
bluestein_real_full(float[], int, int)
bluestein_real_inverse()
bluestein_real_inverse(pl.edu.icm.jlargearrays.FloatLargeArray, long)
bluestein_real_inverse(float[], int)
bluestein_real_inverse2()
bluestein_real_inverse2(pl.edu.icm.jlargearrays.FloatLargeArray, long)
bluestein_real_inverse2(float[], int)
bluesteini()
bluesteini()
bluesteinil()
bluesteinil()
cfftf()
cfftf(pl.edu.icm.jlargearrays.FloatLargeArray, long, int)
cfftf(float[], int, int)
cffti()
cffti(int, int)
cffti()
cfftil()
cfftil()
complexForward()
complexForward(pl.edu.icm.jlargearrays.FloatLargeArray)
complexForward(float[], int)
complexForward(pl.edu.icm.jlargearrays.FloatLargeArray, long)
complexForward(float[])
complexInverse()
complexInverse(float[], boolean)
complexInverse(pl.edu.icm.jlargearrays.FloatLargeArray, long, boolean)
complexInverse(float[], int, boolean)
complexInverse(pl.edu.icm.jlargearrays.FloatLargeArray, boolean)
passf2()
passf2(int, int, float[], int, float[], int, int, int)
passf2(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, long)
passf3()
passf3(int, int, float[], int, float[], int, int, int)
passf3(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, long)
passf4()
passf4(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, int)
passf4(int, int, float[], int, float[], int, int, int)
passf5()
passf5(int, int, float[], int, float[], int, int, int)
passf5(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, long)
passfg()
passfg(int[], int, int, int, int, float[], int, float[], int, int, int)
passfg(int[], long, long, long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, long)
radb2()
radb2(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radb2(int, int, float[], int, float[], int, int)
radb3()
radb3(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radb3(int, int, float[], int, float[], int, int)
radb4()
radb4(int, int, float[], int, float[], int, int)
radb4(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radb5()
radb5(int, int, float[], int, float[], int, int)
radb5(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radbg()
radbg(int, int, int, int, float[], int, float[], int, int)
radbg(long, long, long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf2()
radf2(int, int, float[], int, float[], int, int)
radf2(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf3()
radf3(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf3(int, int, float[], int, float[], int, int)
radf4()
radf4(int, int, float[], int, float[], int, int)
radf4(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf5()
radf5(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf5(int, int, float[], int, float[], int, int)
radfg()
radfg(int, int, int, int, float[], int, float[], int, int)
radfg(long, long, long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
realForward()
realForward(pl.edu.icm.jlargearrays.FloatLargeArray, long)
realForward(float[])
realForward(pl.edu.icm.jlargearrays.FloatLargeArray)
realForward(float[], int)
realForwardFull()
realForwardFull(pl.edu.icm.jlargearrays.FloatLargeArray, long)
realForwardFull(float[])
realForwardFull(float[], int)
realForwardFull(pl.edu.icm.jlargearrays.FloatLargeArray)
realInverse()
realInverse(float[], boolean)
realInverse(float[], int, boolean)
realInverse(pl.edu.icm.jlargearrays.FloatLargeArray, boolean)
realInverse(pl.edu.icm.jlargearrays.FloatLargeArray, long, boolean)
realInverse2()
realInverse2(pl.edu.icm.jlargearrays.FloatLargeArray, long, boolean)
realInverse2(float[], int, boolean)
realInverseFull()
realInverseFull(pl.edu.icm.jlargearrays.FloatLargeArray, boolean)
realInverseFull(float[], boolean)
realInverseFull(pl.edu.icm.jlargearrays.FloatLargeArray, long, boolean)
realInverseFull(float[], int, boolean)
rfftb()
rfftb(float[], int)
rfftb(pl.edu.icm.jlargearrays.FloatLargeArray, long)
rfftf()
rfftf(float[], int)
rfftf(pl.edu.icm.jlargearrays.FloatLargeArray, long)
rffti()
rffti()
rfftil()
rfftil()
[brief]
[full]
org.jtransforms.fft
Computes 3D Discrete Fourier Transform (DFT) of complex and real, single
BenchmarkDoubleFFT
Benchmark of single precision FFT's
BenchmarkFloatFFT
Benchmark of single precision FFT's
DoubleFFT_1D
Computes 1D Discrete Fourier Transform (DFT) of complex and real, double
DoubleFFT_1DTest
This is a series of JUnit tests for the {@link DoubleFFT_1D}. First, {@link DoubleFFT_1D#complexForward(double[])} is tested by comparison with {@link DoubleFFT_1D#complexForward(double[])} as a reference.
DoubleFFT_2D
Computes 2D Discrete Fourier Transform (DFT) of complex and real, double
DoubleFFT_2DTest
This is a test of the class {@link DoubleFFT_2D}. In this test, a very crude 2d FFT method is implemented (see {@link #complexForward(double[][])}), assuming that {@link DoubleFFT_1D} has been fully tested and validated. This
DoubleFFT_3D
Computes 3D Discrete Fourier Transform (DFT) of complex and real, double
DoubleFFT_3DTest
This is a test of the class {@link DoubleFFT_3D}. In this test, a very crude 3d FFT method is implemented (see {@link #complexForward(double[][][])}), assuming that {@link DoubleFFT_1D} and {@link DoubleFFT_2D} have been fully
FloatFFT_1D
Computes 1D Discrete Fourier Transform (DFT) of complex and real, single
FloatFFT_1DTest
This is a series of JUnit tests for the {@link FloatFFT_1D}. First, {@link FloatFFT_1D#complexForward(float[])} is tested by comparison with {@link FloatFFT_1D#complexForward(float[])} as a reference.
FloatFFT_2D
Computes 2D Discrete Fourier Transform (DFT) of complex and real, single
FloatFFT_2DTest
This is a test of the class {@link FloatFFT_2D}. In this test, a very crude 2d FFT method is implemented (see {@link #complexForward(float[][])}), assuming that {@link FloatFFT_1D} has been fully tested and validated. This
FloatFFT_3D
Computes 3D Discrete Fourier Transform (DFT) of complex and real, single
FloatFFT_3DTest
This is a test of the class {@link FloatFFT_3D}. In this test, a very crude 3d FFT method is implemented (see {@link #complexForward(float[][][])}), assuming that {@link FloatFFT_1D} and {@link FloatFFT_2D} have been fully
RealFFTUtils_2D
@formatter:on
RealFFTUtils_2DTest
Test of the utility class {@link RealFFTUtils_2D}.
RealFFTUtils_3D
@formatter:on
RealFFTUtils_3DTest
Test of the utility class {@link RealFFTUtils_3D}.
[save]
LICENSE
README.md
FloatFFT_1D
Computes 1D Discrete Fourier Transform (DFT) of complex and real, single
PI
TWO_PI
bk1
bk1l
bk2
bk2l
factors
ip
ipl
n
nBluestein
nBluesteinl
nc
ncl
nl
nw
nwl
plan
useLargeArrays
w
wl
wtable
wtable_r
wtable_rl
wtablel
bluestein_complex()
bluestein_complex(pl.edu.icm.jlargearrays.FloatLargeArray, long, int)
bluestein_complex(float[], int, int)
bluestein_real_forward()
bluestein_real_forward(pl.edu.icm.jlargearrays.FloatLargeArray, long)
bluestein_real_forward(float[], int)
bluestein_real_full()
bluestein_real_full(pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
bluestein_real_full(float[], int, int)
bluestein_real_inverse()
bluestein_real_inverse(pl.edu.icm.jlargearrays.FloatLargeArray, long)
bluestein_real_inverse(float[], int)
bluestein_real_inverse2()
bluestein_real_inverse2(pl.edu.icm.jlargearrays.FloatLargeArray, long)
bluestein_real_inverse2(float[], int)
bluesteini()
bluesteini()
bluesteinil()
bluesteinil()
cfftf()
cfftf(pl.edu.icm.jlargearrays.FloatLargeArray, long, int)
cfftf(float[], int, int)
cffti()
cffti(int, int)
cffti()
cfftil()
cfftil()
complexForward()
complexForward(pl.edu.icm.jlargearrays.FloatLargeArray)
complexForward(float[], int)
complexForward(pl.edu.icm.jlargearrays.FloatLargeArray, long)
complexForward(float[])
complexInverse()
complexInverse(float[], boolean)
complexInverse(pl.edu.icm.jlargearrays.FloatLargeArray, long, boolean)
complexInverse(float[], int, boolean)
complexInverse(pl.edu.icm.jlargearrays.FloatLargeArray, boolean)
passf2()
passf2(int, int, float[], int, float[], int, int, int)
passf2(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, long)
passf3()
passf3(int, int, float[], int, float[], int, int, int)
passf3(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, long)
passf4()
passf4(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, int)
passf4(int, int, float[], int, float[], int, int, int)
passf5()
passf5(int, int, float[], int, float[], int, int, int)
passf5(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, long)
passfg()
passfg(int[], int, int, int, int, float[], int, float[], int, int, int)
passfg(int[], long, long, long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long, long)
radb2()
radb2(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radb2(int, int, float[], int, float[], int, int)
radb3()
radb3(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radb3(int, int, float[], int, float[], int, int)
radb4()
radb4(int, int, float[], int, float[], int, int)
radb4(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radb5()
radb5(int, int, float[], int, float[], int, int)
radb5(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radbg()
radbg(int, int, int, int, float[], int, float[], int, int)
radbg(long, long, long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf2()
radf2(int, int, float[], int, float[], int, int)
radf2(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf3()
radf3(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf3(int, int, float[], int, float[], int, int)
radf4()
radf4(int, int, float[], int, float[], int, int)
radf4(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf5()
radf5(long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
radf5(int, int, float[], int, float[], int, int)
radfg()
radfg(int, int, int, int, float[], int, float[], int, int)
radfg(long, long, long, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, pl.edu.icm.jlargearrays.FloatLargeArray, long, long)
realForward()
realForward(pl.edu.icm.jlargearrays.FloatLargeArray, long)
realForward(float[])
realForward(pl.edu.icm.jlargearrays.FloatLargeArray)
realForward(float[], int)
realForwardFull()
realForwardFull(pl.edu.icm.jlargearrays.FloatLargeArray, long)
realForwardFull(float[])
realForwardFull(float[], int)
realForwardFull(pl.edu.icm.jlargearrays.FloatLargeArray)
realInverse()
realInverse(float[], boolean)
realInverse(float[], int, boolean)
realInverse(pl.edu.icm.jlargearrays.FloatLargeArray, boolean)
realInverse(pl.edu.icm.jlargearrays.FloatLargeArray, long, boolean)
realInverse2()
realInverse2(pl.edu.icm.jlargearrays.FloatLargeArray, long, boolean)
realInverse2(float[], int, boolean)
realInverseFull()
realInverseFull(pl.edu.icm.jlargearrays.FloatLargeArray, boolean)
realInverseFull(float[], boolean)
realInverseFull(pl.edu.icm.jlargearrays.FloatLargeArray, long, boolean)
realInverseFull(float[], int, boolean)
rfftb()
rfftb(float[], int)
rfftb(pl.edu.icm.jlargearrays.FloatLargeArray, long)
rfftf()
rfftf(float[], int)
rfftf(pl.edu.icm.jlargearrays.FloatLargeArray, long)
rffti()
rffti()
rfftil()
rfftil()
org.jtransforms.fft
Computes 3D Discrete Fourier Transform (DFT) of complex and real, single
BenchmarkDoubleFFT
Benchmark of single precision FFT's
BenchmarkFloatFFT
Benchmark of single precision FFT's
DoubleFFT_1D
Computes 1D Discrete Fourier Transform (DFT) of complex and real, double
DoubleFFT_1DTest
This is a series of JUnit tests for the {@link DoubleFFT_1D}. First, {@link DoubleFFT_1D#complexForward(double[])} is tested by comparison with {@link DoubleFFT_1D#complexForward(double[])} as a reference.
DoubleFFT_2D
Computes 2D Discrete Fourier Transform (DFT) of complex and real, double
DoubleFFT_2DTest
This is a test of the class {@link DoubleFFT_2D}. In this test, a very crude 2d FFT method is implemented (see {@link #complexForward(double[][])}), assuming that {@link DoubleFFT_1D} has been fully tested and validated. This
DoubleFFT_3D
Computes 3D Discrete Fourier Transform (DFT) of complex and real, double
DoubleFFT_3DTest
This is a test of the class {@link DoubleFFT_3D}. In this test, a very crude 3d FFT method is implemented (see {@link #complexForward(double[][][])}), assuming that {@link DoubleFFT_1D} and {@link DoubleFFT_2D} have been fully
FloatFFT_1D
Computes 1D Discrete Fourier Transform (DFT) of complex and real, single
FloatFFT_1DTest
This is a series of JUnit tests for the {@link FloatFFT_1D}. First, {@link FloatFFT_1D#complexForward(float[])} is tested by comparison with {@link FloatFFT_1D#complexForward(float[])} as a reference.
FloatFFT_2D
Computes 2D Discrete Fourier Transform (DFT) of complex and real, single
FloatFFT_2DTest
This is a test of the class {@link FloatFFT_2D}. In this test, a very crude 2d FFT method is implemented (see {@link #complexForward(float[][])}), assuming that {@link FloatFFT_1D} has been fully tested and validated. This
FloatFFT_3D
Computes 3D Discrete Fourier Transform (DFT) of complex and real, single
FloatFFT_3DTest
This is a test of the class {@link FloatFFT_3D}. In this test, a very crude 3d FFT method is implemented (see {@link #complexForward(float[][][])}), assuming that {@link FloatFFT_1D} and {@link FloatFFT_2D} have been fully
RealFFTUtils_2D
@formatter:on
RealFFTUtils_2DTest
Test of the utility class {@link RealFFTUtils_2D}.
RealFFTUtils_3D
@formatter:on
RealFFTUtils_3DTest
Test of the utility class {@link RealFFTUtils_3D}.
LICENSE
README.md
/* ***** BEGIN LICENSE BLOCK ***** * JTransforms * Copyright (c) 2007 onward, Piotr Wendykier * All rights reserved. * * Redistribution and use in source and binary forms, with or without * modification, are permitted provided that the following conditions are met: * * 1. Redistributions of source code must retain the above copyright notice, this * list of conditions and the following disclaimer. * 2. Redistributions in binary form must reproduce the above copyright notice, * this list of conditions and the following disclaimer in the documentation * and/or other materials provided with the distribution. * * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR * ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. * * ***** END LICENSE BLOCK ***** */package org .jtransforms .fft ;import java .util .concurrent .Future ;import org .jtransforms .utils .CommonUtils ;import java .util .concurrent .ExecutionException ;import java .util .logging .Level ;import java .util .logging .Logger ;import pl.edu.icm.jlargearrays.ConcurrencyUtils ;import pl.edu.icm.jlargearrays.FloatLargeArray ;import pl.edu.icm.jlargearrays.LongLargeArray ;import pl.edu.icm.jlargearrays.LargeArrayUtils ;import static org .apache.commons.math3.util.FastMath. *;import pl.edu.icm.jlargearrays.LargeArray ;/** * Computes 1D Discrete Fourier Transform (DFT) of complex and real, single * precision data. The size of the data can be an arbitrary number. This is a * parallel implementation of split-radix and mixed-radix algorithms optimized * for SMP systems. <br> * <br> * This code is derived from General Purpose FFT Package written by Takuya Ooura * ( http://www.kurims.kyoto-u.ac.jp/~ooura/fft.html) and from JFFTPack written * by Baoshe Zhang ( http://jfftpack.sourceforge.net/) * * @author Piotr Wendykier (piotr.wendykier@gmail.com) */ public final class FloatFFT_1D { private static enum Plans { SPLIT_RADIX , MIXED_RADIX , BLUESTEIN } private int n ; private long nl ; private int nBluestein ; private long nBluesteinl ; private int [] ip ; private LongLargeArray ipl ; private float [] w ; private FloatLargeArray wl ; private int nw ; private long nwl ; private int nc ; private long ncl ; private float [] wtable ; private FloatLargeArray wtablel ; private float [] wtable_r ; private FloatLargeArray wtable_rl ; private float [] bk1 ; private FloatLargeArray bk1l ; private float [] bk2 ; private FloatLargeArray bk2l ; private Plans plan ; private boolean useLargeArrays ; private static final int [] factors = {4, 2, 3, 5}; private static final float PI = 3.14159265358979311599796346854418516f; private static final float TWO_PI = 6.28318530717958623199592693708837032f; /** * Creates new instance of FloatFFT_1D. * * @param n size of data */ public FloatFFT_1D (long n ) { if (n < 1) { throw new IllegalArgumentException ("n must be greater than 0"); } this .useLargeArrays = (CommonUtils.isUseLargeArrays () || 2 * n > LargeArray.getMaxSizeOf32bitArray ()); this .n = (int ) n ; this .nl = n ; if (this .useLargeArrays == false ) { if (!CommonUtils.isPowerOf2 (n )) { if (CommonUtils.getReminder (n , factors ) >= 211) { plan = Plans.BLUESTEIN ; nBluestein = CommonUtils.nextPow2 (this .n * 2 - 1); bk1 = new float [2 * nBluestein ]; bk2 = new float [2 * nBluestein ]; this .ip = new int [2 + (int ) ceil (2 + (1 << (int ) (log (nBluestein + 0.5f) / log (2)) / 2))]; this .w = new float [nBluestein ]; int twon = 2 * nBluestein ; nw = twon >> 2; CommonUtils.makewt (nw , ip , w ); nc = nBluestein >> 2; CommonUtils.makect (nc , w , nw , ip ); bluesteini (); } else { plan = Plans.MIXED_RADIX ; wtable = new float [4 * this .n + 15]; wtable_r = new float [2 * this .n + 15]; cffti (); rffti (); } } else { plan = Plans.SPLIT_RADIX ; this .ip = new int [2 + (int ) ceil (2 + (1 << (int ) (log (n + 0.5f) / log (2)) / 2))]; this .w = new float [this .n ]; int twon = 2 * this .n ; nw = twon >> 2; CommonUtils.makewt (nw , ip , w ); nc = this .n >> 2; CommonUtils.makect (nc , w , nw , ip ); } } else if (!CommonUtils.isPowerOf2 (nl )) { if (CommonUtils.getReminder (nl , factors ) >= 211) { plan = Plans.BLUESTEIN ; nBluesteinl = CommonUtils.nextPow2 (nl * 2 - 1); bk1l = new FloatLargeArray (2l * nBluesteinl ); bk2l = new FloatLargeArray (2l * nBluesteinl ); this .ipl = new LongLargeArray (2l + (long ) ceil (2l + (1l << (long ) (log (nBluesteinl + 0.5f) / log (2.)) / 2))); this .wl = new FloatLargeArray (nBluesteinl ); long twon = 2 * nBluesteinl ; nwl = twon >> 2l; CommonUtils.makewt (nwl , ipl , wl ); ncl = nBluesteinl >> 2l; CommonUtils.makect (ncl , wl , nwl , ipl ); bluesteinil (); } else { plan = Plans.MIXED_RADIX ; wtablel = new FloatLargeArray (4 * nl + 15); wtable_rl = new FloatLargeArray (2 * nl + 15); cfftil (); rfftil (); } } else { plan = Plans.SPLIT_RADIX ; this .ipl = new LongLargeArray (2l + (long ) ceil (2 + (1l << (long ) (log (nl + 0.5f) / log (2)) / 2))); this .wl = new FloatLargeArray (nl ); long twon = 2 * nl ; nwl = twon >> 2l; CommonUtils.makewt (nwl , ipl , wl ); ncl = nl >> 2l; CommonUtils.makect (ncl , wl , nwl , ipl ); } } /** * Computes 1D forward DFT of complex data leaving the result in * <code>a</code>. Complex number is stored as two float values in * sequence: the real and imaginary part, i.e. the size of the input array * must be greater or equal 2*n. The physical layout of the input data has * to be as follows:<br> * * <pre> * a[2*k] = Re[k], * a[2*k+1] = Im[k], 0<=k<n * </pre> * * @param a data to transform */ public void complexForward (float [] a ) { complexForward (a , 0); } /** * Computes 1D forward DFT of complex data leaving the result in * <code>a</code>. Complex number is stored as two float values in * sequence: the real and imaginary part, i.e. the size of the input array * must be greater or equal 2*n. The physical layout of the input data has * to be as follows:<br> * * <pre> * a[2*k] = Re[k], * a[2*k+1] = Im[k], 0<=k<n * </pre> * * @param a data to transform */ public void complexForward (FloatLargeArray a ) { complexForward (a , 0); } /** * Computes 1D forward DFT of complex data leaving the result in * <code>a</code>. Complex number is stored as two float values in * sequence: the real and imaginary part, i.e. the size of the input array * must be greater or equal 2*n. The physical layout of the input data has * to be as follows:<br> * * <pre> * a[offa+2*k] = Re[k], * a[offa+2*k+1] = Im[k], 0<=k<n * </pre> * * @param a data to transform * @param offa index of the first element in array <code>a</code> */ public void complexForward (float [] a , int offa ) { if (useLargeArrays ) { complexForward (new FloatLargeArray (a ), offa ); } else { if (n == 1) { return ; } switch (plan ) { case SPLIT_RADIX : CommonUtils.cftbsub (2 * n , a , offa , ip , nw , w ); break ; case MIXED_RADIX : cfftf (a , offa , -1); break ; case BLUESTEIN : bluestein_complex (a , offa , -1); break ; } } } /** * Computes 1D forward DFT of complex data leaving the result in * <code>a</code>. Complex number is stored as two float values in * sequence: the real and imaginary part, i.e. the size of the input array * must be greater or equal 2*n. The physical layout of the input data has * to be as follows:<br> * * <pre> * a[offa+2*k] = Re[k], * a[offa+2*k+1] = Im[k], 0<=k<n * </pre> * * @param a data to transform * @param offa index of the first element in array <code>a</code> */ public void complexForward (FloatLargeArray a , long offa ) { if (!useLargeArrays ) { if (!a.isLarge () && !a.isConstant () && offa < Integer.MAX_VALUE ) { complexForward (a.getData (), (int ) offa ); } else { throw new IllegalArgumentException ("The data array is too big."); } } else { if (nl == 1) { return ; } switch (plan ) { case SPLIT_RADIX : CommonUtils.cftbsub (2 * nl , a , offa , ipl , nwl , wl ); break ; case MIXED_RADIX : cfftf (a , offa , -1); break ; case BLUESTEIN : bluestein_complex (a , offa , -1); break ; } } } /** * Computes 1D inverse DFT of complex data leaving the result in * <code>a</code>. Complex number is stored as two float values in * sequence: the real and imaginary part, i.e. the size of the input array * must be greater or equal 2*n. The physical layout of the input data has * to be as follows:<br> * * <pre> * a[2*k] = Re[k], * a[2*k+1] = Im[k], 0<=k<n * </pre> * * @param adata to transform * @param scale if true then scaling is performed */ public void complexInverse (float [] a , boolean scale ) { complexInverse (a , 0, scale ); } /** * Computes 1D inverse DFT of complex data leaving the result in * <code>a</code>. Complex number is stored as two float values in * sequence: the real and imaginary part, i.e. the size of the input array * must be greater or equal 2*n. The physical layout of the input data has * to be as follows:<br> * * <pre> * a[2*k] = Re[k], * a[2*k+1] = Im[k], 0<=k<n * </pre> * * @param adata to transform * @param scale if true then scaling is performed */ public void complexInverse (FloatLargeArray a , boolean scale ) { complexInverse (a , 0, scale ); } /** * Computes 1D inverse DFT of complex data leaving the result in * <code>a</code>. Complex number is stored as two float values in * sequence: the real and imaginary part, i.e. the size of the input array * must be greater or equal 2*n. The physical layout of the input data has * to be as follows:<br> * * <pre> * a[offa+2*k] = Re[k], * a[offa+2*k+1] = Im[k], 0<=k<n * </pre> * * @param adata to transform * @param offa index of the first element in array <code>a</code> * @param scale if true then scaling is performed */ public void complexInverse (float [] a , int offa , boolean scale ) { if (useLargeArrays ) { complexInverse (new FloatLargeArray (a ), offa , scale ); } else { if (n == 1) { return ; } switch (plan ) { case SPLIT_RADIX : CommonUtils.cftfsub (2 * n , a , offa , ip , nw , w ); break ; case MIXED_RADIX : cfftf (a , offa , +1); break ; case BLUESTEIN : bluestein_complex (a , offa , 1); break ; } if (scale ) { CommonUtils.scale (n , 1.0f / (float ) n , a , offa , true ); } } } /** * Computes 1D inverse DFT of complex data leaving the result in * <code>a</code>. Complex number is stored as two float values in * sequence: the real and imaginary part, i.e. the size of the input array * must be greater or equal 2*n. The physical layout of the input data has * to be as follows:<br> * * <pre> * a[offa+2*k] = Re[k], * a[offa+2*k+1] = Im[k], 0<=k<n * </pre> * * @param adata to transform * @param offa index of the first element in array <code>a</code> * @param scale if true then scaling is performed */ public void complexInverse (FloatLargeArray a , long offa , boolean scale ) { if (!useLargeArrays ) { if (!a.isLarge () && !a.isConstant () && offa < Integer.MAX_VALUE ) { complexInverse (a.getData (), (int ) offa , scale ); } else { throw new IllegalArgumentException ("The data array is too big."); } } else { if (nl == 1) { return ; } switch (plan ) { case SPLIT_RADIX : CommonUtils.cftfsub (2 * nl , a , offa , ipl , nwl , wl ); break ; case MIXED_RADIX : cfftf (a , offa , +1); break ; case BLUESTEIN : bluestein_complex (a , offa , 1); break ; } if (scale ) { CommonUtils.scale (nl , 1.0f / (float ) nl , a , offa , true ); } } } /** * Computes 1D forward DFT of real data leaving the result in <code>a</code> * . The physical layout of the output data is as follows:<br> * * if n is even then * * <pre> * a[2*k] = Re[k], 0<=k<n/2 * a[2*k+1] = Im[k], 0<k<n/2 * a[1] = Re[n/2] * </pre> * * if n is odd then * * <pre> * a[2*k] = Re[k], 0<=k<(n+1)/2 * a[2*k+1] = Im[k], 0<k<(n-1)/2 * a[1] = Im[(n-1)/2] * </pre> * * This method computes only half of the elements of the real transform. The * other half satisfies the symmetry condition. If you want the full real * forward transform, use <code>realForwardFull</code>. To get back the * original data, use <code>realInverse</code> on the output of this method. * * @param a data to transform */ public void realForward (float [] a ) { realForward (a , 0); } /** * Computes 1D forward DFT of real data leaving the result in <code>a</code> * . The physical layout of the output data is as follows:<br> * * if n is even then * * <pre> * a[2*k] = Re[k], 0<=k<n/2 * a[2*k+1] = Im[k], 0<k<n/2 a[1] = Re[n/2] * </pre> * * if n is odd then * * <pre> * a[2*k] = Re[k], 0<=k<(n+1)/2 * a[2*k+1] = Im[k], 0<k<(n-1)/2 * a[1] = Im[(n-1)/2] * </pre> * * This method computes only half of the elements of the real transform. The * other half satisfies the symmetry condition. If you want the full real * forward transform, use <code>realForwardFull</code>. To get back the * original data, use <code>realInverse</code> on the output of this method. * * @param a data to transform */ public void realForward (FloatLargeArray a ) { realForward (a , 0); } /** * Computes 1D forward DFT of real data leaving the result in <code>a</code> * . The physical layout of the output data is as follows:<br> * * if n is even then * * <pre> * a[offa+2*k] = Re[k], 0<=k<n/2 * a[offa+2*k+1] = Im[k], 0<k<n/2 * a[offa+1] = Re[n/2] * </pre> * * if n is odd then * * <pre> * a[offa+2*k] = Re[k], 0<=k<(n+1)/2 * a[offa+2*k+1] = Im[k], 0<k<(n-1)/2 * a[offa+1] = Im[(n-1)/2] * </pre> * * This method computes only half of the elements of the real transform. The * other half satisfies the symmetry condition. If you want the full real * forward transform, use <code>realForwardFull</code>. To get back the * original data, use <code>realInverse</code> on the output of this method. * * @param a data to transform * @param offa index of the first element in array <code>a</code> */ public void realForward (float [] a , int offa ) { if (useLargeArrays ) { realForward (new FloatLargeArray (a ), offa ); } else { if (n == 1) { return ; } switch (plan ) { case SPLIT_RADIX : float xi ; if (n > 4) { CommonUtils.cftfsub (n , a , offa , ip , nw , w ); CommonUtils.rftfsub (n , a , offa , nc , w , nw ); } else if (n == 4) { CommonUtils.cftx020 (a , offa ); } xi = a [offa ] - a [offa + 1]; a [offa ] += a [offa + 1]; a [offa + 1] = xi ; break ; case MIXED_RADIX : rfftf (a , offa ); for (int k = n - 1; k >= 2; k --) { int idx = offa + k ; float tmp = a [idx ]; a [idx ] = a [idx - 1]; a [idx - 1] = tmp ; } break ; case BLUESTEIN : bluestein_real_forward (a , offa ); break ; } } } /** * Computes 1D forward DFT of real data leaving the result in <code>a</code> * . The physical layout of the output data is as follows:<br> * * if n is even then * * <pre> * a[offa+2*k] = Re[k], 0<=k<n/2 * a[offa+2*k+1] = Im[k], 0<k<n/2 * a[offa+1] = Re[n/2] * </pre> * * if n is odd then * * <pre> * a[offa+2*k] = Re[k], 0<=k<(n+1)/2 * a[offa+2*k+1] = Im[k], 0<k<(n-1)/2 * a[offa+1] = Im[(n-1)/2] * </pre> * * This method computes only half of the elements of the real transform. The * other half satisfies the symmetry condition. If you want the full real * forward transform, use <code>realForwardFull</code>. To get back the * original data, use <code>realInverse</code> on the output of this method. * * @param a data to transform * @param offa index of the first element in array <code>a</code> */ public void realForward (FloatLargeArray a , long offa ) { if (!useLargeArrays ) { if (!a.isLarge () && !a.isConstant () && offa < Integer.MAX_VALUE ) { realForward (a.getData (), (int ) offa ); } else { throw new IllegalArgumentException ("The data array is too big."); } } else { if (nl == 1) { return ; } switch (plan ) { case SPLIT_RADIX : float xi ; if (nl > 4) { CommonUtils.cftfsub (nl , a , offa , ipl , nwl , wl ); CommonUtils.rftfsub (nl , a , offa , ncl , wl , nwl ); } else if (nl == 4) { CommonUtils.cftx020 (a , offa ); } xi = a.getFloat (offa ) - a.getFloat (offa + 1); a.setFloat (offa , a.getFloat (offa ) + a.getFloat (offa + 1)); a.setFloat (offa + 1, xi ); break ; case MIXED_RADIX : rfftf (a , offa ); for (long k = nl - 1; k >= 2; k --) { long idx = offa + k ; float tmp = a.getFloat (idx ); a.setFloat (idx , a.getFloat (idx - 1)); a.setFloat (idx - 1, tmp ); } break ; case BLUESTEIN : bluestein_real_forward (a , offa ); break ; } } } /** * Computes 1D forward DFT of real data leaving the result in <code>a</code> * . This method computes the full real forward transform, i.e. you will get * the same result as from <code>complexForward</code> called with all * imaginary parts equal 0. Because the result is stored in <code>a</code>, * the size of the input array must greater or equal 2*n, with only the * first n elements filled with real data. To get back the original data, * use <code>complexInverse</code> on the output of this method. * * @param a data to transform */ public void realForwardFull (float [] a ) { realForwardFull (a , 0); } /** * Computes 1D forward DFT of real data leaving the result in <code>a</code> * . This method computes the full real forward transform, i.e. you will get * the same result as from <code>complexForward</code> called with all * imaginary parts equal 0. Because the result is stored in <code>a</code>, * the size of the input array must greater or equal 2*n, with only the * first n elements filled with real data. To get back the original data, * use <code>complexInverse</code> on the output of this method. * * @param a data to transform */ public void realForwardFull (FloatLargeArray a ) { realForwardFull (a , 0); } /** * Computes 1D forward DFT of real data leaving the result in <code>a</code> * . This method computes the full real forward transform, i.e. you will get * the same result as from <code>complexForward</code> called with all * imaginary part equal 0. Because the result is stored in <code>a</code>, * the size of the input array must greater or equal 2*n, with only the * first n elements filled with real data. To get back the original data, * use <code>complexInverse</code> on the output of this method. * * @param a data to transform * @param offa index of the first element in array <code>a</code> */ public void realForwardFull (final float [] a , final int offa ) { if (useLargeArrays ) { realForwardFull (new FloatLargeArray (a ), offa ); } else { final int twon = 2 * n ; switch (plan ) { case SPLIT_RADIX : realForward (a , offa ); int nthreads = ConcurrencyUtils.getNumberOfThreads (); if ((nthreads > 1) && (n / 2 > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { Future <?> [] futures = new Future [nthreads ]; int k = n / 2 / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n / 2 : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { int idx1 , idx2 ; for (int k = firstIdx ; k < lastIdx ; k ++) { idx1 = 2 * k ; idx2 = offa + ((twon - idx1 ) % twon ); a [idx2 ] = a [offa + idx1 ]; a [idx2 + 1] = -a [offa + idx1 + 1]; } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { int idx1 , idx2 ; for (int k = 0; k < n / 2; k ++) { idx1 = 2 * k ; idx2 = offa + ((twon - idx1 ) % twon ); a [idx2 ] = a [offa + idx1 ]; a [idx2 + 1] = -a [offa + idx1 + 1]; } } a [offa + n ] = -a [offa + 1]; a [offa + 1] = 0; break ; case MIXED_RADIX : rfftf (a , offa ); int m ; if (n % 2 == 0) { m = n / 2; } else { m = (n + 1) / 2; } for (int k = 1; k < m ; k ++) { int idx1 = offa + twon - 2 * k ; int idx2 = offa + 2 * k ; a [idx1 + 1] = -a [idx2 ]; a [idx1 ] = a [idx2 - 1]; } for (int k = 1; k < n ; k ++) { int idx = offa + n - k ; float tmp = a [idx + 1]; a [idx + 1] = a [idx ]; a [idx ] = tmp ; } a [offa + 1] = 0; break ; case BLUESTEIN : bluestein_real_full (a , offa , -1); break ; } } } /** * Computes 1D forward DFT of real data leaving the result in <code>a</code> * . This method computes the full real forward transform, i.e. you will get * the same result as from <code>complexForward</code> called with all * imaginary part equal 0. Because the result is stored in <code>a</code>, * the size of the input array must greater or equal 2*n, with only the * first n elements filled with real data. To get back the original data, * use <code>complexInverse</code> on the output of this method. * * @param a data to transform * @param offa index of the first element in array <code>a</code> */ public void realForwardFull (final FloatLargeArray a , final long offa ) { if (!useLargeArrays ) { if (!a.isLarge () && !a.isConstant () && offa < Integer.MAX_VALUE ) { realForwardFull (a.getData (), (int ) offa ); } else { throw new IllegalArgumentException ("The data array is too big."); } } else { final long twon = 2 * nl ; switch (plan ) { case SPLIT_RADIX : realForward (a , offa ); int nthreads = ConcurrencyUtils.getNumberOfThreads (); if ((nthreads > 1) && (nl / 2 > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { Future <?> [] futures = new Future [nthreads ]; long k = nl / 2 / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl / 2 : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { long idx1 , idx2 ; for (long k = firstIdx ; k < lastIdx ; k ++) { idx1 = 2 * k ; idx2 = offa + ((twon - idx1 ) % twon ); a.setFloat (idx2 , a.getFloat (offa + idx1 )); a.setFloat (idx2 + 1, -a.getFloat (offa + idx1 + 1)); } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { long idx1 , idx2 ; for (long k = 0; k < nl / 2; k ++) { idx1 = 2 * k ; idx2 = offa + ((twon - idx1 ) % twon ); a.setFloat (idx2 , a.getFloat (offa + idx1 )); a.setFloat (idx2 + 1, -a.getFloat (offa + idx1 + 1)); } } a.setFloat (offa + nl , -a.getFloat (offa + 1)); a.setFloat (offa + 1, 0); break ; case MIXED_RADIX : rfftf (a , offa ); long m ; if (nl % 2 == 0) { m = nl / 2; } else { m = (nl + 1) / 2; } for (long k = 1; k < m ; k ++) { long idx1 = offa + twon - 2 * k ; long idx2 = offa + 2 * k ; a.setFloat (idx1 + 1, -a.getFloat (idx2 )); a.setFloat (idx1 , a.getFloat (idx2 - 1)); } for (long k = 1; k < nl ; k ++) { long idx = offa + nl - k ; float tmp = a.getFloat (idx + 1); a.setFloat (idx + 1, a.getFloat (idx )); a.setFloat (idx , tmp ); } a.setFloat (offa + 1, 0); break ; case BLUESTEIN : bluestein_real_full (a , offa , -1); break ; } } } /** * Computes 1D inverse DFT of real data leaving the result in <code>a</code> * . The physical layout of the input data has to be as follows:<br> * * if n is even then * * <pre> * a[2*k] = Re[k], 0<=k<n/2 * a[2*k+1] = Im[k], 0<k<n/2 * a[1] = Re[n/2] * </pre> * * if n is odd then * * <pre> * a[2*k] = Re[k], 0<=k<(n+1)/2 * a[2*k+1] = Im[k], 0<k<(n-1)/2 * a[1] = Im[(n-1)/2] * </pre> * * This method computes only half of the elements of the real transform. The * other half satisfies the symmetry condition. If you want the full real * inverse transform, use <code>realInverseFull</code>. * * @param adata to transform * * @param scale if true then scaling is performed * */ public void realInverse (float [] a , boolean scale ) { realInverse (a , 0, scale ); } /** * Computes 1D inverse DFT of real data leaving the result in <code>a</code> * . The physical layout of the input data has to be as follows:<br> * * if n is even then * * <pre> * a[2*k] = Re[k], 0<=k<n/2 * a[2*k+1] = Im[k], 0<k<n/2 * a[1] = Re[n/2] * </pre> * * if n is odd then * * <pre> * a[2*k] = Re[k], 0<=k<(n+1)/2 * a[2*k+1] = Im[k], 0<k<(n-1)/2 * a[1] = Im[(n-1)/2] * </pre> * * This method computes only half of the elements of the real transform. The * other half satisfies the symmetry condition. If you want the full real * inverse transform, use <code>realInverseFull</code>. * * @param adata to transform * * @param scale if true then scaling is performed * */ public void realInverse (FloatLargeArray a , boolean scale ) { realInverse (a , 0, scale ); } /** * Computes 1D inverse DFT of real data leaving the result in <code>a</code> * . The physical layout of the input data has to be as follows:<br> * * if n is even then * * <pre> * a[offa+2*k] = Re[k], 0<=k<n/2 * a[offa+2*k+1] = Im[k], 0<k<n/2 * a[offa+1] = Re[n/2] * </pre> * * if n is odd then * * <pre> * a[offa+2*k] = Re[k], 0<=k<(n+1)/2 * a[offa+2*k+1] = Im[k], 0<k<(n-1)/2 * a[offa+1] = Im[(n-1)/2] * </pre> * * This method computes only half of the elements of the real transform. The * other half satisfies the symmetry condition. If you want the full real * inverse transform, use <code>realInverseFull</code>. * * @param adata to transform * @param offa index of the first element in array <code>a</code> * @param scale if true then scaling is performed * */ public void realInverse (float [] a , int offa , boolean scale ) { if (useLargeArrays ) { realInverse (new FloatLargeArray (a ), offa , scale ); } else { if (n == 1) { return ; } switch (plan ) { case SPLIT_RADIX : a [offa + 1] = 0.5f * (a [offa ] - a [offa + 1]); a [offa ] -= a [offa + 1]; if (n > 4) { CommonUtils.rftfsub (n , a , offa , nc , w , nw ); CommonUtils.cftbsub (n , a , offa , ip , nw , w ); } else if (n == 4) { CommonUtils.cftxc020 (a , offa ); } if (scale ) { CommonUtils.scale (n , 1.0f / (n / 2.0f), a , offa , false ); } break ; case MIXED_RADIX : for (int k = 2; k < n ; k ++) { int idx = offa + k ; float tmp = a [idx - 1]; a [idx - 1] = a [idx ]; a [idx ] = tmp ; } rfftb (a , offa ); if (scale ) { CommonUtils.scale (n , 1.0f / n , a , offa , false ); } break ; case BLUESTEIN : bluestein_real_inverse (a , offa ); if (scale ) { CommonUtils.scale (n , 1.0f / n , a , offa , false ); } break ; } } } /** * Computes 1D inverse DFT of real data leaving the result in <code>a</code> * . The physical layout of the input data has to be as follows:<br> * * if n is even then * * <pre> * a[offa+2*k] = Re[k], 0<=k<n/2 * a[offa+2*k+1] = Im[k], 0<k<n/2 * a[offa+1] = Re[n/2] * </pre> * * if n is odd then * * <pre> * a[offa+2*k] = Re[k], 0<=k<(n+1)/2 * a[offa+2*k+1] = Im[k], 0<k<(n-1)/2 * a[offa+1] = Im[(n-1)/2] * </pre> * * This method computes only half of the elements of the real transform. The * other half satisfies the symmetry condition. If you want the full real * inverse transform, use <code>realInverseFull</code>. * * @param adata to transform * @param offa index of the first element in array <code>a</code> * @param scale if true then scaling is performed * */ public void realInverse (FloatLargeArray a , long offa , boolean scale ) { if (!useLargeArrays ) { if (!a.isLarge () && !a.isConstant () && offa < Integer.MAX_VALUE ) { realInverse (a.getData (), (int ) offa , scale ); } else { throw new IllegalArgumentException ("The data array is too big."); } } else { if (nl == 1) { return ; } switch (plan ) { case SPLIT_RADIX : a.setFloat (offa + 1, 0.5f * (a.getFloat (offa ) - a.getFloat (offa + 1))); a.setFloat (offa , a.getFloat (offa ) - a.getFloat (offa + 1)); if (nl > 4) { CommonUtils.rftfsub (nl , a , offa , ncl , wl , nwl ); CommonUtils.cftbsub (nl , a , offa , ipl , nwl , wl ); } else if (nl == 4) { CommonUtils.cftxc020 (a , offa ); } if (scale ) { CommonUtils.scale (nl , 1.0f / (nl / 2.0f), a , offa , false ); } break ; case MIXED_RADIX : for (long k = 2; k < nl ; k ++) { long idx = offa + k ; float tmp = a.getFloat (idx - 1); a.setFloat (idx - 1, a.getFloat (idx )); a.setFloat (idx , tmp ); } rfftb (a , offa ); if (scale ) { CommonUtils.scale (nl , 1.0f / nl , a , offa , false ); } break ; case BLUESTEIN : bluestein_real_inverse (a , offa ); if (scale ) { CommonUtils.scale (nl , 1.0f / nl , a , offa , false ); } break ; } } } /** * Computes 1D inverse DFT of real data leaving the result in <code>a</code> * . This method computes the full real inverse transform, i.e. you will get * the same result as from <code>complexInverse</code> called with all * imaginary part equal 0. Because the result is stored in <code>a</code>, * the size of the input array must greater or equal 2*n, with only the * first n elements filled with real data. * * @param adata to transform * @param scale if true then scaling is performed */ public void realInverseFull (float [] a , boolean scale ) { realInverseFull (a , 0, scale ); } /** * Computes 1D inverse DFT of real data leaving the result in <code>a</code> * . This method computes the full real inverse transform, i.e. you will get * the same result as from <code>complexInverse</code> called with all * imaginary part equal 0. Because the result is stored in <code>a</code>, * the size of the input array must greater or equal 2*n, with only the * first n elements filled with real data. * * @param adata to transform * @param scale if true then scaling is performed */ public void realInverseFull (FloatLargeArray a , boolean scale ) { realInverseFull (a , 0, scale ); } /** * Computes 1D inverse DFT of real data leaving the result in <code>a</code> * . This method computes the full real inverse transform, i.e. you will get * the same result as from <code>complexInverse</code> called with all * imaginary part equal 0. Because the result is stored in <code>a</code>, * the size of the input array must greater or equal 2*n, with only the * first n elements filled with real data. * * @param adata to transform * @param offa index of the first element in array <code>a</code> * @param scale if true then scaling is performed */ public void realInverseFull (final float [] a , final int offa , boolean scale ) { if (useLargeArrays ) { realInverseFull (new FloatLargeArray (a ), offa , scale ); } else { final int twon = 2 * n ; switch (plan ) { case SPLIT_RADIX : realInverse2 (a , offa , scale ); int nthreads = ConcurrencyUtils.getNumberOfThreads (); if ((nthreads > 1) && (n / 2 > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { Future <?> [] futures = new Future [nthreads ]; int k = n / 2 / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n / 2 : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { int idx1 , idx2 ; for (int k = firstIdx ; k < lastIdx ; k ++) { idx1 = 2 * k ; idx2 = offa + ((twon - idx1 ) % twon ); a [idx2 ] = a [offa + idx1 ]; a [idx2 + 1] = -a [offa + idx1 + 1]; } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { int idx1 , idx2 ; for (int k = 0; k < n / 2; k ++) { idx1 = 2 * k ; idx2 = offa + ((twon - idx1 ) % twon ); a [idx2 ] = a [offa + idx1 ]; a [idx2 + 1] = -a [offa + idx1 + 1]; } } a [offa + n ] = -a [offa + 1]; a [offa + 1] = 0; break ; case MIXED_RADIX : rfftf (a , offa ); if (scale ) { CommonUtils.scale (n , 1.0f / n , a , offa , false ); } int m ; if (n % 2 == 0) { m = n / 2; } else { m = (n + 1) / 2; } for (int k = 1; k < m ; k ++) { int idx1 = offa + 2 * k ; int idx2 = offa + twon - 2 * k ; a [idx1 ] = -a [idx1 ]; a [idx2 + 1] = -a [idx1 ]; a [idx2 ] = a [idx1 - 1]; } for (int k = 1; k < n ; k ++) { int idx = offa + n - k ; float tmp = a [idx + 1]; a [idx + 1] = a [idx ]; a [idx ] = tmp ; } a [offa + 1] = 0; break ; case BLUESTEIN : bluestein_real_full (a , offa , 1); if (scale ) { CommonUtils.scale (n , 1.0f / n , a , offa , true ); } break ; } } } /** * Computes 1D inverse DFT of real data leaving the result in <code>a</code> * . This method computes the full real inverse transform, i.e. you will get * the same result as from <code>complexInverse</code> called with all * imaginary part equal 0. Because the result is stored in <code>a</code>, * the size of the input array must greater or equal 2*n, with only the * first n elements filled with real data. * * @param adata to transform * @param offa index of the first element in array <code>a</code> * @param scale if true then scaling is performed */ public void realInverseFull (final FloatLargeArray a , final long offa , boolean scale ) { if (!useLargeArrays ) { if (!a.isLarge () && !a.isConstant () && offa < Integer.MAX_VALUE ) { realInverseFull (a.getData (), (int ) offa , scale ); } else { throw new IllegalArgumentException ("The data array is too big."); } } else { final long twon = 2 * nl ; switch (plan ) { case SPLIT_RADIX : realInverse2 (a , offa , scale ); int nthreads = ConcurrencyUtils.getNumberOfThreads (); if ((nthreads > 1) && (nl / 2 > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { Future <?> [] futures = new Future [nthreads ]; long k = nl / 2 / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl / 2 : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { long idx1 , idx2 ; for (long k = firstIdx ; k < lastIdx ; k ++) { idx1 = 2 * k ; idx2 = offa + ((twon - idx1 ) % twon ); a.setFloat (idx2 , a.getFloat (offa + idx1 )); a.setFloat (idx2 + 1, -a.getFloat (offa + idx1 + 1)); } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { long idx1 , idx2 ; for (long k = 0; k < nl / 2; k ++) { idx1 = 2 * k ; idx2 = offa + ((twon - idx1 ) % twon ); a.setFloat (idx2 , a.getFloat (offa + idx1 )); a.setFloat (idx2 + 1, -a.getFloat (offa + idx1 + 1)); } } a.setFloat (offa + nl , -a.getFloat (offa + 1)); a.setFloat (offa + 1, 0); break ; case MIXED_RADIX : rfftf (a , offa ); if (scale ) { CommonUtils.scale (nl , 1.0f / nl , a , offa , false ); } long m ; if (nl % 2 == 0) { m = nl / 2; } else { m = (nl + 1) / 2; } for (long k = 1; k < m ; k ++) { long idx1 = offa + 2 * k ; long idx2 = offa + twon - 2 * k ; a.setFloat (idx1 , -a.getFloat (idx1 )); a.setFloat (idx2 + 1, -a.getFloat (idx1 )); a.setFloat (idx2 , a.getFloat (idx1 - 1)); } for (long k = 1; k < nl ; k ++) { long idx = offa + nl - k ; float tmp = a.getFloat (idx + 1); a.setFloat (idx + 1, a.getFloat (idx )); a.setFloat (idx , tmp ); } a.setFloat (offa + 1, 0); break ; case BLUESTEIN : bluestein_real_full (a , offa , 1); if (scale ) { CommonUtils.scale (nl , 1.0f / nl , a , offa , true ); } break ; } } } protected void realInverse2 (float [] a , int offa , boolean scale ) { if (useLargeArrays ) { realInverse2 (new FloatLargeArray (a ), offa , scale ); } else { if (n == 1) { return ; } switch (plan ) { case SPLIT_RADIX : float xi ; if (n > 4) { CommonUtils.cftfsub (n , a , offa , ip , nw , w ); CommonUtils.rftbsub (n , a , offa , nc , w , nw ); } else if (n == 4) { CommonUtils.cftbsub (n , a , offa , ip , nw , w ); } xi = a [offa ] - a [offa + 1]; a [offa ] += a [offa + 1]; a [offa + 1] = xi ; if (scale ) { CommonUtils.scale (n , 1.0f / n , a , offa , false ); } break ; case MIXED_RADIX : rfftf (a , offa ); for (int k = n - 1; k >= 2; k --) { int idx = offa + k ; float tmp = a [idx ]; a [idx ] = a [idx - 1]; a [idx - 1] = tmp ; } if (scale ) { CommonUtils.scale (n , 1.0f / n , a , offa , false ); } int m ; if (n % 2 == 0) { m = n / 2; for (int i = 1; i < m ; i ++) { int idx = offa + 2 * i + 1; a [idx ] = -a [idx ]; } } else { m = (n - 1) / 2; for (int i = 0; i < m ; i ++) { int idx = offa + 2 * i + 1; a [idx ] = -a [idx ]; } } break ; case BLUESTEIN : bluestein_real_inverse2 (a , offa ); if (scale ) { CommonUtils.scale (n , 1.0f / n , a , offa , false ); } break ; } } } protected void realInverse2 (FloatLargeArray a , long offa , boolean scale ) { if (!useLargeArrays ) { if (!a.isLarge () && !a.isConstant () && offa < Integer.MAX_VALUE ) { realInverse2 (a.getData (), (int ) offa , scale ); } else { throw new IllegalArgumentException ("The data array is too big."); } } else { if (nl == 1) { return ; } switch (plan ) { case SPLIT_RADIX : float xi ; if (nl > 4) { CommonUtils.cftfsub (nl , a , offa , ipl , nwl , wl ); CommonUtils.rftbsub (nl , a , offa , ncl , wl , nwl ); } else if (nl == 4) { CommonUtils.cftbsub (nl , a , offa , ipl , nwl , wl ); } xi = a.getFloat (offa ) - a.getFloat (offa + 1); a.setFloat (offa , a.getFloat (offa ) + a.getFloat (offa + 1)); a.setFloat (offa + 1, xi ); if (scale ) { CommonUtils.scale (nl , 1.0f / nl , a , offa , false ); } break ; case MIXED_RADIX : rfftf (a , offa ); for (long k = nl - 1; k >= 2; k --) { long idx = offa + k ; float tmp = a.getFloat (idx ); a.setFloat (idx , a.getFloat (idx - 1)); a.setFloat (idx - 1, tmp ); } if (scale ) { CommonUtils.scale (nl , 1.0f / nl , a , offa , false ); } long m ; if (nl % 2 == 0) { m = nl / 2; for (long i = 1; i < m ; i ++) { long idx = offa + 2 * i + 1; a.setFloat (idx , -a.getFloat (idx )); } } else { m = (nl - 1) / 2; for (long i = 0; i < m ; i ++) { long idx = offa + 2 * i + 1; a.setFloat (idx , -a.getFloat (idx )); } } break ; case BLUESTEIN : bluestein_real_inverse2 (a , offa ); if (scale ) { CommonUtils.scale (nl , 1.0f / nl , a , offa , false ); } break ; } } } /* -------- initializing routines -------- */ /*--------------------------------------------------------- cffti: initialization of Complex FFT --------------------------------------------------------*/ void cffti (int n , int offw ) { if (n == 1) { return ; } final int twon = 2 * n ; final int fourn = 4 * n ; float argh ; int idot , ntry = 0, i , j ; float argld ; int i1 , k1 , l1 , l2 , ib ; float fi ; int ld , ii , nf , ipll , nll , nq , nr ; float arg ; int ido , ipm ; nll = n ; nf = 0; j = 0; factorize_loop : while (true ) { j ++; if (j <= 4) { ntry = factors [j - 1]; } else { ntry += 2; } do { nq = nll / ntry ; nr = nll - ntry * nq ; if (nr != 0) { continue factorize_loop ; } nf ++; wtable [offw + nf + 1 + fourn ] = ntry ; nll = nq ; if (ntry == 2 && nf != 1) { for (i = 2; i <= nf ; i ++) { ib = nf - i + 2; int idx = ib + fourn ; wtable [offw + idx + 1] = wtable [offw + idx ]; } wtable [offw + 2 + fourn ] = 2; } } while (nll != 1); break ; } wtable [offw + fourn ] = n ; wtable [offw + 1 + fourn ] = nf ; argh = TWO_PI / (float ) n ; i = 1; l1 = 1; for (k1 = 1; k1 <= nf ; k1 ++) { ipll = (int ) wtable [offw + k1 + 1 + fourn ]; ld = 0; l2 = l1 * ipll ; ido = n / l2 ; idot = ido + ido + 2; ipm = ipll - 1; for (j = 1; j <= ipm ; j ++) { i1 = i ; wtable [offw + i - 1 + twon ] = 1; wtable [offw + i + twon ] = 0; ld += l1 ; fi = 0; argld = ld * argh ; for (ii = 4; ii <= idot ; ii += 2) { i += 2; fi += 1; arg = fi * argld ; int idx = i + twon ; wtable [offw + idx - 1] = (float ) cos (arg ); wtable [offw + idx ] = (float ) sin (arg ); } if (ipll > 5) { int idx1 = i1 + twon ; int idx2 = i + twon ; wtable [offw + idx1 - 1] = wtable [offw + idx2 - 1]; wtable [offw + idx1 ] = wtable [offw + idx2 ]; } } l1 = l2 ; } } final void cffti () { if (n == 1) { return ; } final int twon = 2 * n ; final int fourn = 4 * n ; float argh ; int idot , ntry = 0, i , j ; float argld ; int i1 , k1 , l1 , l2 , ib ; float fi ; int ld , ii , nf , ipll , nll , nq , nr ; float arg ; int ido , ipm ; nll = n ; nf = 0; j = 0; factorize_loop : while (true ) { j ++; if (j <= 4) { ntry = factors [j - 1]; } else { ntry += 2; } do { nq = nll / ntry ; nr = nll - ntry * nq ; if (nr != 0) { continue factorize_loop ; } nf ++; wtable [nf + 1 + fourn ] = ntry ; nll = nq ; if (ntry == 2 && nf != 1) { for (i = 2; i <= nf ; i ++) { ib = nf - i + 2; int idx = ib + fourn ; wtable [idx + 1] = wtable [idx ]; } wtable [2 + fourn ] = 2; } } while (nll != 1); break ; } wtable [fourn ] = n ; wtable [1 + fourn ] = nf ; argh = TWO_PI / (float ) n ; i = 1; l1 = 1; for (k1 = 1; k1 <= nf ; k1 ++) { ipll = (int ) wtable [k1 + 1 + fourn ]; ld = 0; l2 = l1 * ipll ; ido = n / l2 ; idot = ido + ido + 2; ipm = ipll - 1; for (j = 1; j <= ipm ; j ++) { i1 = i ; wtable [i - 1 + twon ] = 1; wtable [i + twon ] = 0; ld += l1 ; fi = 0; argld = ld * argh ; for (ii = 4; ii <= idot ; ii += 2) { i += 2; fi += 1; arg = fi * argld ; int idx = i + twon ; wtable [idx - 1] = (float ) cos (arg ); wtable [idx ] = (float ) sin (arg ); } if (ipll > 5) { int idx1 = i1 + twon ; int idx2 = i + twon ; wtable [idx1 - 1] = wtable [idx2 - 1]; wtable [idx1 ] = wtable [idx2 ]; } } l1 = l2 ; } } final void cfftil () { if (nl == 1) { return ; } final long twon = 2 * nl ; final long fourn = 4 * nl ; float argh ; long idot , ntry = 0, i , j ; float argld ; long i1 , k1 , l1 , l2 , ib ; float fi ; long ld , ii , nf , ipll , nl2 , nq , nr ; float arg ; long ido , ipm ; nl2 = nl ; nf = 0; j = 0; factorize_loop : while (true ) { j ++; if (j <= 4) { ntry = factors [(int ) (j - 1)]; } else { ntry += 2; } do { nq = nl2 / ntry ; nr = nl2 - ntry * nq ; if (nr != 0) { continue factorize_loop ; } nf ++; wtablel.setFloat (nf + 1 + fourn , ntry ); nl2 = nq ; if (ntry == 2 && nf != 1) { for (i = 2; i <= nf ; i ++) { ib = nf - i + 2; long idx = ib + fourn ; wtablel.setFloat (idx + 1, wtablel.getFloat (idx )); } wtablel.setFloat (2 + fourn , 2); } } while (nl2 != 1); break ; } wtablel.setFloat (fourn , nl ); wtablel.setFloat (1 + fourn , nf ); argh = TWO_PI / (float ) nl ; i = 1; l1 = 1; for (k1 = 1; k1 <= nf ; k1 ++) { ipll = (long ) wtablel.getFloat (k1 + 1 + fourn ); ld = 0; l2 = l1 * ipll ; ido = nl / l2 ; idot = ido + ido + 2; ipm = ipll - 1; for (j = 1; j <= ipm ; j ++) { i1 = i ; wtablel.setFloat (i - 1 + twon , 1); wtablel.setFloat (i + twon , 0); ld += l1 ; fi = 0; argld = ld * argh ; for (ii = 4; ii <= idot ; ii += 2) { i += 2; fi += 1; arg = fi * argld ; long idx = i + twon ; wtablel.setFloat (idx - 1, (float ) cos (arg )); wtablel.setFloat (idx , (float ) sin (arg )); } if (ipll > 5) { long idx1 = i1 + twon ; long idx2 = i + twon ; wtablel.setFloat (idx1 - 1, wtablel.getFloat (idx2 - 1)); wtablel.setFloat (idx1 , wtablel.getFloat (idx2 )); } } l1 = l2 ; } } void rffti () { if (n == 1) { return ; } final int twon = 2 * n ; float argh ; int ntry = 0, i , j ; float argld ; int k1 , l1 , l2 , ib ; float fi ; int ld , ii , nf , ipll , nll , is , nq , nr ; float arg ; int ido , ipm ; int nfm1 ; nll = n ; nf = 0; j = 0; factorize_loop : while (true ) { ++j ; if (j <= 4) { ntry = factors [j - 1]; } else { ntry += 2; } do { nq = nll / ntry ; nr = nll - ntry * nq ; if (nr != 0) { continue factorize_loop ; } ++nf ; wtable_r [nf + 1 + twon ] = ntry ; nll = nq ; if (ntry == 2 && nf != 1) { for (i = 2; i <= nf ; i ++) { ib = nf - i + 2; int idx = ib + twon ; wtable_r [idx + 1] = wtable_r [idx ]; } wtable_r [2 + twon ] = 2; } } while (nll != 1); break ; } wtable_r [twon ] = n ; wtable_r [1 + twon ] = nf ; argh = TWO_PI / (float ) (n ); is = 0; nfm1 = nf - 1; l1 = 1; if (nfm1 == 0) { return ; } for (k1 = 1; k1 <= nfm1 ; k1 ++) { ipll = (int ) wtable_r [k1 + 1 + twon ]; ld = 0; l2 = l1 * ipll ; ido = n / l2 ; ipm = ipll - 1; for (j = 1; j <= ipm ; ++j ) { ld += l1 ; i = is ; argld = (float ) ld * argh ; fi = 0; for (ii = 3; ii <= ido ; ii += 2) { i += 2; fi += 1; arg = fi * argld ; int idx = i + n ; wtable_r [idx - 2] = (float ) cos (arg ); wtable_r [idx - 1] = (float ) sin (arg ); } is += ido ; } l1 = l2 ; } } void rfftil () { if (nl == 1) { return ; } final long twon = 2 * nl ; float argh ; long ntry = 0, i , j ; float argld ; long k1 , l1 , l2 , ib ; float fi ; long ld , ii , nf , ipll , nl2 , is , nq , nr ; float arg ; long ido , ipm ; long nfm1 ; nl2 = nl ; nf = 0; j = 0; factorize_loop : while (true ) { ++j ; if (j <= 4) { ntry = factors [(int ) (j - 1)]; } else { ntry += 2; } do { nq = nl2 / ntry ; nr = nl2 - ntry * nq ; if (nr != 0) { continue factorize_loop ; } ++nf ; wtable_rl.setFloat (nf + 1 + twon , ntry ); nl2 = nq ; if (ntry == 2 && nf != 1) { for (i = 2; i <= nf ; i ++) { ib = nf - i + 2; long idx = ib + twon ; wtable_rl.setFloat (idx + 1, wtable_rl.getFloat (idx )); } wtable_rl.setFloat (2 + twon , 2); } } while (nl2 != 1); break ; } wtable_rl.setFloat (twon , nl ); wtable_rl.setFloat (1 + twon , nf ); argh = TWO_PI / (float ) (nl ); is = 0; nfm1 = nf - 1; l1 = 1; if (nfm1 == 0) { return ; } for (k1 = 1; k1 <= nfm1 ; k1 ++) { ipll = (long ) wtable_rl.getFloat (k1 + 1 + twon ); ld = 0; l2 = l1 * ipll ; ido = nl / l2 ; ipm = ipll - 1; for (j = 1; j <= ipm ; ++j ) { ld += l1 ; i = is ; argld = (float ) ld * argh ; fi = 0; for (ii = 3; ii <= ido ; ii += 2) { i += 2; fi += 1; arg = fi * argld ; long idx = i + nl ; wtable_rl.setFloat (idx - 2, (float ) cos (arg )); wtable_rl.setFloat (idx - 1, (float ) sin (arg )); } is += ido ; } l1 = l2 ; } } private void bluesteini () { int k = 0; float arg ; float pi_n = PI / n ; bk1 [0] = 1; bk1 [1] = 0; for (int i = 1; i < n ; i ++) { k += 2 * i - 1; if (k >= 2 * n ) { k -= 2 * n ; } arg = pi_n * k ; bk1 [2 * i ] = (float ) cos (arg ); bk1 [2 * i + 1] = (float ) sin (arg ); } float scale = 1.0f / nBluestein ; bk2 [0] = bk1 [0] * scale ; bk2 [1] = bk1 [1] * scale ; for (int i = 2; i < 2 * n ; i += 2) { bk2 [i ] = bk1 [i ] * scale ; bk2 [i + 1] = bk1 [i + 1] * scale ; bk2 [2 * nBluestein - i ] = bk2 [i ]; bk2 [2 * nBluestein - i + 1] = bk2 [i + 1]; } CommonUtils.cftbsub (2 * nBluestein , bk2 , 0, ip , nw , w ); } private void bluesteinil () { long k = 0; float arg ; float pi_n = PI / nl ; bk1l.setFloat (0, 1); bk1l.setFloat (1, 0); for (int i = 1; i < nl ; i ++) { k += 2 * i - 1; if (k >= 2 * nl ) { k -= 2 * nl ; } arg = pi_n * k ; bk1l.setFloat (2 * i , (float ) cos (arg )); bk1l.setFloat (2 * i + 1, (float ) sin (arg )); } float scale = 1.0f / nBluesteinl ; bk2l.setFloat (0, bk1l.getFloat (0) * scale ); bk2l.setFloat (1, bk1l.getFloat (1) * scale ); for (int i = 2; i < 2 * nl ; i += 2) { bk2l.setFloat (i , bk1l.getFloat (i ) * scale ); bk2l.setFloat (i + 1, bk1l.getFloat (i + 1) * scale ); bk2l.setFloat (2 * nBluesteinl - i , bk2l.getFloat (i )); bk2l.setFloat (2 * nBluesteinl - i + 1, bk2l.getFloat (i + 1)); } CommonUtils.cftbsub (2 * nBluesteinl , bk2l , 0, ipl , nwl , wl ); } private void bluestein_complex (final float [] a , final int offa , final int isign ) { final float [] ak = new float [2 * nBluestein ]; int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; int k = n / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ] - a [idx4 ] * bk1 [idx2 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ] + a [idx4 ] * bk1 [idx1 ]; } } else { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ] + a [idx4 ] * bk1 [idx2 ]; ak [idx2 ] = -a [idx3 ] * bk1 [idx2 ] + a [idx4 ] * bk1 [idx1 ]; } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); k = nBluestein / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? nBluestein : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = -ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] + ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } else { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] - ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftfsub (2 * nBluestein , ak , 0, ip , nw , w ); k = n / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; a [idx3 ] = bk1 [idx1 ] * ak [idx1 ] - bk1 [idx2 ] * ak [idx2 ]; a [idx4 ] = bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } else { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; a [idx3 ] = bk1 [idx1 ] * ak [idx1 ] + bk1 [idx2 ] * ak [idx2 ]; a [idx4 ] = -bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { if (isign > 0) { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ] - a [idx4 ] * bk1 [idx2 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ] + a [idx4 ] * bk1 [idx1 ]; } } else { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ] + a [idx4 ] * bk1 [idx2 ]; ak [idx2 ] = -a [idx3 ] * bk1 [idx2 ] + a [idx4 ] * bk1 [idx1 ]; } } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); if (isign > 0) { for (int i = 0; i < nBluestein ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = -ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] + ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } else { for (int i = 0; i < nBluestein ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] - ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } CommonUtils.cftfsub (2 * nBluestein , ak , 0, ip , nw , w ); if (isign > 0) { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; a [idx3 ] = bk1 [idx1 ] * ak [idx1 ] - bk1 [idx2 ] * ak [idx2 ]; a [idx4 ] = bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } else { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; a [idx3 ] = bk1 [idx1 ] * ak [idx1 ] + bk1 [idx2 ] * ak [idx2 ]; a [idx4 ] = -bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } } } private void bluestein_complex (final FloatLargeArray a , final long offa , final int isign ) { final FloatLargeArray ak = new FloatLargeArray (2 * nBluesteinl ); int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; long k = nl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 ) - a.getFloat (idx4 ) * bk1l.getFloat (idx2 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx1 )); } } else { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx2 )); ak.setFloat (idx2 , -a.getFloat (idx3 ) * bk1l.getFloat (idx2 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx1 )); } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); k = nBluesteinl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nBluesteinl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = -ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } else { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) - ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftfsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); k = nl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; a.setFloat (idx3 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) - bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (idx4 , bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } else { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; a.setFloat (idx3 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (idx4 , -bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { if (isign > 0) { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 ) - a.getFloat (idx4 ) * bk1l.getFloat (idx2 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx1 )); } } else { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx2 )); ak.setFloat (idx2 , -a.getFloat (idx3 ) * bk1l.getFloat (idx2 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx1 )); } } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); if (isign > 0) { for (long i = 0; i < nBluesteinl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = -ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } else { for (long i = 0; i < nBluesteinl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) - ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } CommonUtils.cftfsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); if (isign > 0) { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; a.setFloat (idx3 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) - bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (idx4 , bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } else { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; a.setFloat (idx3 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (idx4 , -bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } } } private void bluestein_real_full (final float [] a , final int offa , final int isign ) { final float [] ak = new float [2 * nBluestein ]; int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; int k = n / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + i ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ]; } } else { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + i ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ]; ak [idx2 ] = -a [idx3 ] * bk1 [idx2 ]; } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); k = nBluestein / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? nBluestein : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = -ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] + ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } else { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] - ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftfsub (2 * nBluestein , ak , 0, ip , nw , w ); k = n / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + idx1 ] = bk1 [idx1 ] * ak [idx1 ] - bk1 [idx2 ] * ak [idx2 ]; a [offa + idx2 ] = bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } else { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + idx1 ] = bk1 [idx1 ] * ak [idx1 ] + bk1 [idx2 ] * ak [idx2 ]; a [offa + idx2 ] = -bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { if (isign > 0) { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + i ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ]; } } else { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + i ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ]; ak [idx2 ] = -a [idx3 ] * bk1 [idx2 ]; } } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); if (isign > 0) { for (int i = 0; i < nBluestein ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = -ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] + ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } else { for (int i = 0; i < nBluestein ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] - ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } CommonUtils.cftfsub (2 * nBluestein , ak , 0, ip , nw , w ); if (isign > 0) { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + idx1 ] = bk1 [idx1 ] * ak [idx1 ] - bk1 [idx2 ] * ak [idx2 ]; a [offa + idx2 ] = bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } else { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + idx1 ] = bk1 [idx1 ] * ak [idx1 ] + bk1 [idx2 ] * ak [idx2 ]; a [offa + idx2 ] = -bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } } } private void bluestein_real_full (final FloatLargeArray a , final long offa , final long isign ) { final FloatLargeArray ak = new FloatLargeArray (2 * nBluesteinl ); int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; long k = nl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + i ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 )); } } else { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + i ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 )); ak.setFloat (idx2 , -a.getFloat (idx3 ) * bk1l.getFloat (idx2 )); } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); k = nBluesteinl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nBluesteinl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = -ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } else { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) - ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftfsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); k = nl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { if (isign > 0) { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + idx1 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) - bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (offa + idx2 , bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } else { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + idx1 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (offa + idx2 , -bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { if (isign > 0) { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + i ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 )); } } else { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + i ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 )); ak.setFloat (idx2 , -a.getFloat (idx3 ) * bk1l.getFloat (idx2 )); } } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); if (isign > 0) { for (long i = 0; i < nBluesteinl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = -ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } else { for (long i = 0; i < nBluesteinl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) - ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } CommonUtils.cftfsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); if (isign > 0) { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + idx1 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) - bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (offa + idx2 , bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } else { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + idx1 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (offa + idx2 , -bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } } } private void bluestein_real_forward (final float [] a , final int offa ) { final float [] ak = new float [2 * nBluestein ]; int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; int k = n / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + i ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ]; ak [idx2 ] = -a [idx3 ] * bk1 [idx2 ]; } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); k = nBluestein / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? nBluestein : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] - ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + i ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ]; ak [idx2 ] = -a [idx3 ] * bk1 [idx2 ]; } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); for (int i = 0; i < nBluestein ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] - ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } CommonUtils.cftfsub (2 * nBluestein , ak , 0, ip , nw , w ); if (n % 2 == 0) { a [offa ] = bk1 [0] * ak [0] + bk1 [1] * ak [1]; a [offa + 1] = bk1 [n ] * ak [n ] + bk1 [n + 1] * ak [n + 1]; for (int i = 1; i < n / 2; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + idx1 ] = bk1 [idx1 ] * ak [idx1 ] + bk1 [idx2 ] * ak [idx2 ]; a [offa + idx2 ] = -bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } else { a [offa ] = bk1 [0] * ak [0] + bk1 [1] * ak [1]; a [offa + 1] = -bk1 [n ] * ak [n - 1] + bk1 [n - 1] * ak [n ]; for (int i = 1; i < (n - 1) / 2; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + idx1 ] = bk1 [idx1 ] * ak [idx1 ] + bk1 [idx2 ] * ak [idx2 ]; a [offa + idx2 ] = -bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } a [offa + n - 1] = bk1 [n - 1] * ak [n - 1] + bk1 [n ] * ak [n ]; } } private void bluestein_real_forward (final FloatLargeArray a , final long offa ) { final FloatLargeArray ak = new FloatLargeArray (2 * nBluesteinl ); int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; long k = nl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + i ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 )); ak.setFloat (idx2 , -a.getFloat (idx3 ) * bk1l.getFloat (idx2 )); } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); k = nBluesteinl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nBluesteinl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) - ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + i ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 )); ak.setFloat (idx2 , -a.getFloat (idx3 ) * bk1l.getFloat (idx2 )); } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); for (long i = 0; i < nBluesteinl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) - ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } CommonUtils.cftfsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); if (nl % 2 == 0) { a.setFloat (offa , bk1l.getFloat (0) * ak.getFloat (0) + bk1l.getFloat (1) * ak.getFloat (1)); a.setFloat (offa + 1, bk1l.getFloat (nl ) * ak.getFloat (nl ) + bk1l.getFloat (nl + 1) * ak.getFloat (nl + 1)); for (long i = 1; i < nl / 2; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + idx1 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (offa + idx2 , -bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } else { a.setFloat (offa , bk1l.getFloat (0) * ak.getFloat (0) + bk1l.getFloat (1) * ak.getFloat (1)); a.setFloat (offa + 1, -bk1l.getFloat (nl ) * ak.getFloat (nl - 1) + bk1l.getFloat (nl - 1) * ak.getFloat (nl )); for (long i = 1; i < (nl - 1) / 2; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + idx1 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (offa + idx2 , -bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } a.setFloat (offa + nl - 1, bk1l.getFloat (nl - 1) * ak.getFloat (nl - 1) + bk1l.getFloat (nl ) * ak.getFloat (nl )); } } private void bluestein_real_inverse (final float [] a , final int offa ) { final float [] ak = new float [2 * nBluestein ]; if (n % 2 == 0) { ak [0] = a [offa ] * bk1 [0]; ak [1] = a [offa ] * bk1 [1]; for (int i = 1; i < n / 2; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ] - a [idx4 ] * bk1 [idx2 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ] + a [idx4 ] * bk1 [idx1 ]; } ak [n ] = a [offa + 1] * bk1 [n ]; ak [n + 1] = a [offa + 1] * bk1 [n + 1]; for (int i = n / 2 + 1; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + 2 * n - idx1 ; int idx4 = idx3 + 1; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ] + a [idx4 ] * bk1 [idx2 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ] - a [idx4 ] * bk1 [idx1 ]; } } else { ak [0] = a [offa ] * bk1 [0]; ak [1] = a [offa ] * bk1 [1]; for (int i = 1; i < (n - 1) / 2; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + idx1 ; int idx4 = offa + idx2 ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ] - a [idx4 ] * bk1 [idx2 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ] + a [idx4 ] * bk1 [idx1 ]; } ak [n - 1] = a [offa + n - 1] * bk1 [n - 1] - a [offa + 1] * bk1 [n ]; ak [n ] = a [offa + n - 1] * bk1 [n ] + a [offa + 1] * bk1 [n - 1]; ak [n + 1] = a [offa + n - 1] * bk1 [n + 1] + a [offa + 1] * bk1 [n + 2]; ak [n + 2] = a [offa + n - 1] * bk1 [n + 2] - a [offa + 1] * bk1 [n + 1]; for (int i = (n - 1) / 2 + 2; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + 2 * n - idx1 ; int idx4 = idx3 + 1; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ] + a [idx4 ] * bk1 [idx2 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ] - a [idx4 ] * bk1 [idx1 ]; } } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; int k = nBluestein / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? nBluestein : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = -ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] + ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftfsub (2 * nBluestein , ak , 0, ip , nw , w ); k = n / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + i ] = bk1 [idx1 ] * ak [idx1 ] - bk1 [idx2 ] * ak [idx2 ]; } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { for (int i = 0; i < nBluestein ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = -ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] + ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } CommonUtils.cftfsub (2 * nBluestein , ak , 0, ip , nw , w ); for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + i ] = bk1 [idx1 ] * ak [idx1 ] - bk1 [idx2 ] * ak [idx2 ]; } } } private void bluestein_real_inverse (final FloatLargeArray a , final long offa ) { final FloatLargeArray ak = new FloatLargeArray (2 * nBluesteinl ); if (nl % 2 == 0) { ak.setFloat (0, a.getFloat (offa ) * bk1l.getFloat (0)); ak.setFloat (1, a.getFloat (offa ) * bk1l.getFloat (1)); for (long i = 1; i < nl / 2; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 ) - a.getFloat (idx4 ) * bk1l.getFloat (idx2 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx1 )); } ak.setFloat (nl , a.getFloat (offa + 1) * bk1l.getFloat (nl )); ak.setFloat (nl + 1, a.getFloat (offa + 1) * bk1l.getFloat (nl + 1)); for (long i = nl / 2 + 1; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + 2 * nl - idx1 ; long idx4 = idx3 + 1; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx2 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 ) - a.getFloat (idx4 ) * bk1l.getFloat (idx1 )); } } else { ak.setFloat (0, a.getFloat (offa ) * bk1l.getFloat (0)); ak.setFloat (1, a.getFloat (offa ) * bk1l.getFloat (1)); for (long i = 1; i < (nl - 1) / 2; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + idx1 ; long idx4 = offa + idx2 ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 ) - a.getFloat (idx4 ) * bk1l.getFloat (idx2 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx1 )); } ak.setFloat (nl - 1, a.getFloat (offa + nl - 1) * bk1l.getFloat (nl - 1) - a.getFloat (offa + 1) * bk1l.getFloat (nl )); ak.setFloat (nl , a.getFloat (offa + nl - 1) * bk1l.getFloat (nl ) + a.getFloat (offa + 1) * bk1l.getFloat (nl - 1)); ak.setFloat (nl + 1, a.getFloat (offa + nl - 1) * bk1l.getFloat (nl + 1) + a.getFloat (offa + 1) * bk1l.getFloat (nl + 2)); ak.setFloat (nl + 2, a.getFloat (offa + nl - 1) * bk1l.getFloat (nl + 2) - a.getFloat (offa + 1) * bk1l.getFloat (nl + 1)); for (long i = (nl - 1) / 2 + 2; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + 2 * nl - idx1 ; long idx4 = idx3 + 1; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 ) + a.getFloat (idx4 ) * bk1l.getFloat (idx2 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 ) - a.getFloat (idx4 ) * bk1l.getFloat (idx1 )); } } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; long k = nBluesteinl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nBluesteinl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = -ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftfsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); k = nl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + i , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) - bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { for (long i = 0; i < nBluesteinl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = -ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } CommonUtils.cftfsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + i , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) - bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); } } } private void bluestein_real_inverse2 (final float [] a , final int offa ) { final float [] ak = new float [2 * nBluestein ]; int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (n >= CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; int k = n / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? n : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + i ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ]; } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); k = nBluestein / nthreads ; for (int i = 0; i < nthreads ; i ++) { final int firstIdx = i * k ; final int lastIdx = (i == (nthreads - 1)) ? nBluestein : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (int i = firstIdx ; i < lastIdx ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = -ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] + ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { for (int i = 0; i < n ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; int idx3 = offa + i ; ak [idx1 ] = a [idx3 ] * bk1 [idx1 ]; ak [idx2 ] = a [idx3 ] * bk1 [idx2 ]; } CommonUtils.cftbsub (2 * nBluestein , ak , 0, ip , nw , w ); for (int i = 0; i < nBluestein ; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; float im = -ak [idx1 ] * bk2 [idx2 ] + ak [idx2 ] * bk2 [idx1 ]; ak [idx1 ] = ak [idx1 ] * bk2 [idx1 ] + ak [idx2 ] * bk2 [idx2 ]; ak [idx2 ] = im ; } } CommonUtils.cftfsub (2 * nBluestein , ak , 0, ip , nw , w ); if (n % 2 == 0) { a [offa ] = bk1 [0] * ak [0] - bk1 [1] * ak [1]; a [offa + 1] = bk1 [n ] * ak [n ] - bk1 [n + 1] * ak [n + 1]; for (int i = 1; i < n / 2; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + idx1 ] = bk1 [idx1 ] * ak [idx1 ] - bk1 [idx2 ] * ak [idx2 ]; a [offa + idx2 ] = bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } } else { a [offa ] = bk1 [0] * ak [0] - bk1 [1] * ak [1]; a [offa + 1] = bk1 [n ] * ak [n - 1] + bk1 [n - 1] * ak [n ]; for (int i = 1; i < (n - 1) / 2; i ++) { int idx1 = 2 * i ; int idx2 = idx1 + 1; a [offa + idx1 ] = bk1 [idx1 ] * ak [idx1 ] - bk1 [idx2 ] * ak [idx2 ]; a [offa + idx2 ] = bk1 [idx2 ] * ak [idx1 ] + bk1 [idx1 ] * ak [idx2 ]; } a [offa + n - 1] = bk1 [n - 1] * ak [n - 1] - bk1 [n ] * ak [n ]; } } private void bluestein_real_inverse2 (final FloatLargeArray a , final long offa ) { final FloatLargeArray ak = new FloatLargeArray (2 * nBluesteinl ); int threads = ConcurrencyUtils.getNumberOfThreads (); if ((threads > 1) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_2Threads ())) { int nthreads = 2; if ((threads >= 4) && (nl > CommonUtils.getThreadsBeginN_1D_FFT_4Threads ())) { nthreads = 4; } Future <?> [] futures = new Future [nthreads ]; long k = nl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + i ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 )); } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); k = nBluesteinl / nthreads ; for (int i = 0; i < nthreads ; i ++) { final long firstIdx = i * k ; final long lastIdx = (i == (nthreads - 1)) ? nBluesteinl : firstIdx + k ; futures [i ] = ConcurrencyUtils.submit (new Runnable () { public void run () { for (long i = firstIdx ; i < lastIdx ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = -ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } }); } try { ConcurrencyUtils.waitForCompletion (futures ); } catch (InterruptedException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } catch (ExecutionException ex ) { Logger.getLogger (FloatFFT_1D .class .getName ()).log (Level.SEVERE , null , ex ); } } else { for (long i = 0; i < nl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; long idx3 = offa + i ; ak.setFloat (idx1 , a.getFloat (idx3 ) * bk1l.getFloat (idx1 )); ak.setFloat (idx2 , a.getFloat (idx3 ) * bk1l.getFloat (idx2 )); } CommonUtils.cftbsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); for (long i = 0; i < nBluesteinl ; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; float im = -ak.getFloat (idx1 ) * bk2l.getFloat (idx2 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx1 ); ak.setFloat (idx1 , ak.getFloat (idx1 ) * bk2l.getFloat (idx1 ) + ak.getFloat (idx2 ) * bk2l.getFloat (idx2 )); ak.setFloat (idx2 , im ); } } CommonUtils.cftfsub (2 * nBluesteinl , ak , 0, ipl , nwl , wl ); if (nl % 2 == 0) { a.setFloat (offa , bk1l.getFloat (0) * ak.getFloat (0) - bk1l.getFloat (1) * ak.getFloat (1)); a.setFloat (offa + 1, bk1l.getFloat (nl ) * ak.getFloat (nl ) - bk1l.getFloat (nl + 1) * ak.getFloat (nl + 1)); for (long i = 1; i < nl / 2; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + idx1 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) - bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (offa + idx2 , bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } } else { a.setFloat (offa , bk1l.getFloat (0) * ak.getFloat (0) - bk1l.getFloat (1) * ak.getFloat (1)); a.setFloat (offa + 1, bk1l.getFloat (nl ) * ak.getFloat (nl - 1) + bk1l.getFloat (nl - 1) * ak.getFloat (nl )); for (long i = 1; i < (nl - 1) / 2; i ++) { long idx1 = 2 * i ; long idx2 = idx1 + 1; a.setFloat (offa + idx1 , bk1l.getFloat (idx1 ) * ak.getFloat (idx1 ) - bk1l.getFloat (idx2 ) * ak.getFloat (idx2 )); a.setFloat (offa + idx2 , bk1l.getFloat (idx2 ) * ak.getFloat (idx1 ) + bk1l.getFloat (idx1 ) * ak.getFloat (idx2 )); } a.setFloat (offa + nl - 1, bk1l.getFloat (nl - 1) * ak.getFloat (nl - 1) - bk1l.getFloat (nl ) * ak.getFloat (nl )); } } /*--------------------------------------------------------- rfftf1: further processing of Real forward FFT --------------------------------------------------------*/ void rfftf (final float a [], final int offa ) { if (n == 1) { return ; } int l1 , l2 , na , kh , nf , ipll , iw , ido , idl1 ; final float [] ch = new float [n ]; final int twon = 2 * n ; nf = (int ) wtable_r [1 + twon ]; na = 1; l2 = n ; iw = twon - 1; for (int k1 = 1; k1 <= nf ; ++k1 ) { kh = nf - k1 ; ipll = (int ) wtable_r [kh + 2 + twon ]; l1 = l2 / ipll ; ido = n / l2 ; idl1 = ido * l1 ; iw -= (ipll - 1) * ido ; na = 1 - na ; switch (ipll ) { case 2: if (na == 0) { radf2 (ido , l1 , a , offa , ch , 0, iw ); } else { radf2 (ido , l1 , ch , 0, a , offa , iw ); } break ; case 3: if (na == 0) { radf3 (ido , l1 , a , offa , ch , 0, iw ); } else { radf3 (ido , l1 , ch , 0, a , offa , iw ); } break ; case 4: if (na == 0) { radf4 (ido , l1 , a , offa , ch , 0, iw ); } else { radf4 (ido , l1 , ch , 0, a , offa , iw ); } break ; case 5: if (na == 0) { radf5 (ido , l1 , a , offa , ch , 0, iw ); } else { radf5 (ido , l1 , ch , 0, a , offa , iw ); } break ; default : if (ido == 1) { na = 1 - na ; } if (na == 0) { radfg (ido , ipll , l1 , idl1 , a , offa , ch , 0, iw ); na = 1; } else { radfg (ido , ipll , l1 , idl1 , ch , 0, a , offa , iw ); na = 0; } break ; } l2 = l1 ; } if (na == 1) { return ; } System.arraycopy (ch , 0, a , offa , n ); } /*--------------------------------------------------------- rfftf1: further processing of Real forward FFT --------------------------------------------------------*/ void rfftf (final FloatLargeArray a , final long offa ) { if (nl == 1) { return ; } long l1 , l2 , na , kh , nf , iw , ido , idl1 ; int ipll ; final FloatLargeArray ch = new FloatLargeArray (nl ); final long twon = 2 * nl ; nf = (long ) wtable_rl.getFloat (1 + twon ); na = 1; l2 = nl ; iw = twon - 1; for (long k1 = 1; k1 <= nf ; ++k1 ) { kh = nf - k1 ; ipll = (int ) wtable_rl.getFloat (kh + 2 + twon ); l1 = l2 / ipll ; ido = nl / l2 ; idl1 = ido * l1 ; iw -= (ipll - 1) * ido ; na = 1 - na ; switch (ipll ) { case 2: if (na == 0) { radf2 (ido , l1 , a , offa , ch , 0, iw ); } else { radf2 (ido , l1 , ch , 0, a , offa , iw ); } break ; case 3: if (na == 0) { radf3 (ido , l1 , a , offa , ch , 0, iw ); } else { radf3 (ido , l1 , ch , 0, a , offa , iw ); } break ; case 4: if (na == 0) { radf4 (ido , l1 , a , offa , ch , 0, iw ); } else { radf4 (ido , l1 , ch , 0, a , offa , iw ); } break ; case 5: if (na == 0) { radf5 (ido , l1 , a , offa , ch , 0, iw ); } else { radf5 (ido , l1 , ch , 0, a , offa , iw ); } break ; default : if (ido == 1) { na = 1 - na ; } if (na == 0) { radfg (ido , ipll , l1 , idl1 , a , offa , ch , 0, iw ); na = 1; } else { radfg (ido , ipll , l1 , idl1 , ch , 0, a , offa , iw ); na = 0; } break ; } l2 = l1 ; } if (na == 1) { return ; } LargeArrayUtils.arraycopy (ch , 0, a , offa , nl ); } /*--------------------------------------------------------- rfftb1: further processing of Real backward FFT --------------------------------------------------------*/ void rfftb (final float a [], final int offa ) { if (n == 1) { return ; } int l1 , l2 , na , nf , ipll , iw , ido , idl1 ; float [] ch = new float [n ]; final int twon = 2 * n ; nf = (int ) wtable_r [1 + twon ]; na = 0; l1 = 1; iw = n ; for (int k1 = 1; k1 <= nf ; k1 ++) { ipll = (int ) wtable_r [k1 + 1 + twon ]; l2 = ipll * l1 ; ido = n / l2 ; idl1 = ido * l1 ; switch (ipll ) { case 2: if (na == 0) { radb2 (ido , l1 , a , offa , ch , 0, iw ); } else { radb2 (ido , l1 , ch , 0, a , offa , iw ); } na = 1 - na ; break ; case 3: if (na == 0) { radb3 (ido , l1 , a , offa , ch , 0, iw ); } else { radb3 (ido , l1 , ch , 0, a , offa , iw ); } na = 1 - na ; break ; case 4: if (na == 0) { radb4 (ido , l1 , a , offa , ch , 0, iw ); } else { radb4 (ido , l1 , ch , 0, a , offa , iw ); } na = 1 - na ; break ; case 5: if (na == 0) { radb5 (ido , l1 , a , offa , ch , 0, iw ); } else { radb5 (ido , l1 , ch , 0, a , offa , iw ); } na = 1 - na ; break ; default : if (na == 0) { radbg (ido , ipll , l1 , idl1 , a , offa , ch , 0, iw ); } else { radbg (ido , ipll , l1 , idl1 , ch , 0, a , offa , iw ); } if (ido == 1) { na = 1 - na ; } break ; } l1 = l2 ; iw += (ipll - 1) * ido ; } if (na == 0) { return ; } System.arraycopy (ch , 0, a , offa , n ); } /*--------------------------------------------------------- rfftb1: further processing of Real backward FFT --------------------------------------------------------*/ void rfftb (final FloatLargeArray a , final long offa ) { if (nl == 1) { return ; } long l1 , l2 , na , nf , iw , ido , idl1 ; int ipll ; FloatLargeArray ch = new FloatLargeArray (nl ); final long twon = 2 * nl ; nf = (long ) wtable_rl.getFloat (1 + twon ); na = 0; l1 = 1; iw = nl ; for (long k1 = 1; k1 <= nf ; k1 ++) { ipll = (int ) wtable_rl.getFloat (k1 + 1 + twon ); l2 = ipll * l1 ; ido = nl / l2 ; idl1 = ido * l1 ; switch (ipll ) { case 2: if (na == 0) { radb2 (ido , l1 , a , offa , ch , 0, iw ); } else { radb2 (ido , l1 , ch , 0, a , offa , iw ); } na = 1 - na ; break ; case 3: if (na == 0) { radb3 (ido , l1 , a , offa , ch , 0, iw ); } else { radb3 (ido , l1 , ch , 0, a , offa , iw ); } na = 1 - na ; break ; case 4: if (na == 0) { radb4 (ido , l1 , a , offa , ch , 0, iw ); } else { radb4 (ido , l1 , ch , 0, a , offa , iw ); } na = 1 - na ; break ; case 5: if (na == 0) { radb5 (ido , l1 , a , offa , ch , 0, iw ); } else { radb5 (ido , l1 , ch , 0, a , offa , iw ); } na = 1 - na ; break ; default : if (na == 0) { radbg (ido , ipll , l1 , idl1 , a , offa , ch , 0, iw ); } else { radbg (ido , ipll , l1 , idl1 , ch , 0, a , offa , iw ); } if (ido == 1) { na = 1 - na ; } break ; } l1 = l2 ; iw += (ipll - 1) * ido ; } if (na == 0) { return ; } LargeArrayUtils.arraycopy (ch , 0, a , offa , nl ); } /*------------------------------------------------- radf2: Real FFT's forward processing of factor 2 -------------------------------------------------*/ void radf2 (final int ido , final int l1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { int i , ic , idx0 , idx1 , idx2 , idx3 , idx4 ; float t1i , t1r , w1r , w1i ; int iw1 ; iw1 = offset ; idx0 = l1 * ido ; idx1 = 2 * ido ; for (int k = 0; k < l1 ; k ++) { int oidx1 = out_off + k * idx1 ; int oidx2 = oidx1 + idx1 - 1; int iidx1 = in_off + k * ido ; int iidx2 = iidx1 + idx0 ; float i1r = in [iidx1 ]; float i2r = in [iidx2 ]; out [oidx1 ] = i1r + i2r ; out [oidx2 ] = i1r - i2r ; } if (ido < 2) { return ; } if (ido != 2) { for (int k = 0; k < l1 ; k ++) { idx1 = k * ido ; idx2 = 2 * idx1 ; idx3 = idx2 + ido ; idx4 = idx1 + idx0 ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; int widx1 = i - 1 + iw1 ; int oidx1 = out_off + i + idx2 ; int oidx2 = out_off + ic + idx3 ; int iidx1 = in_off + i + idx1 ; int iidx2 = in_off + i + idx4 ; float a1i = in [iidx1 - 1]; float a1r = in [iidx1 ]; float a2i = in [iidx2 - 1]; float a2r = in [iidx2 ]; w1r = wtable_r [widx1 - 1]; w1i = wtable_r [widx1 ]; t1r = w1r * a2i + w1i * a2r ; t1i = w1r * a2r - w1i * a2i ; out [oidx1 ] = a1r + t1i ; out [oidx1 - 1] = a1i + t1r ; out [oidx2 ] = t1i - a1r ; out [oidx2 - 1] = a1i - t1r ; } } if (ido % 2 == 1) { return ; } } idx2 = 2 * idx1 ; for (int k = 0; k < l1 ; k ++) { idx1 = k * ido ; int oidx1 = out_off + idx2 + ido ; int iidx1 = in_off + ido - 1 + idx1 ; out [oidx1 ] = -in [iidx1 + idx0 ]; out [oidx1 - 1] = in [iidx1 ]; } } /*------------------------------------------------- radf2: Real FFT's forward processing of factor 2 -------------------------------------------------*/ void radf2 (final long ido , final long l1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { long i , ic , idx0 , idx1 , idx2 , idx3 , idx4 ; float t1i , t1r , w1r , w1i ; long iw1 ; iw1 = offset ; idx0 = l1 * ido ; idx1 = 2 * ido ; for (long k = 0; k < l1 ; k ++) { long oidx1 = out_off + k * idx1 ; long oidx2 = oidx1 + idx1 - 1; long iidx1 = in_off + k * ido ; long iidx2 = iidx1 + idx0 ; float i1r = in.getFloat (iidx1 ); float i2r = in.getFloat (iidx2 ); out.setFloat (oidx1 , i1r + i2r ); out.setFloat (oidx2 , i1r - i2r ); } if (ido < 2) { return ; } if (ido != 2) { for (long k = 0; k < l1 ; k ++) { idx1 = k * ido ; idx2 = 2 * idx1 ; idx3 = idx2 + ido ; idx4 = idx1 + idx0 ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; long widx1 = i - 1 + iw1 ; long oidx1 = out_off + i + idx2 ; long oidx2 = out_off + ic + idx3 ; long iidx1 = in_off + i + idx1 ; long iidx2 = in_off + i + idx4 ; float a1i = in.getFloat (iidx1 - 1); float a1r = in.getFloat (iidx1 ); float a2i = in.getFloat (iidx2 - 1); float a2r = in.getFloat (iidx2 ); w1r = wtable_rl.getFloat (widx1 - 1); w1i = wtable_rl.getFloat (widx1 ); t1r = w1r * a2i + w1i * a2r ; t1i = w1r * a2r - w1i * a2i ; out.setFloat (oidx1 , a1r + t1i ); out.setFloat (oidx1 - 1, a1i + t1r ); out.setFloat (oidx2 , t1i - a1r ); out.setFloat (oidx2 - 1, a1i - t1r ); } } if (ido % 2 == 1) { return ; } } idx2 = 2 * idx1 ; for (long k = 0; k < l1 ; k ++) { idx1 = k * ido ; long oidx1 = out_off + idx2 + ido ; long iidx1 = in_off + ido - 1 + idx1 ; out.setFloat (oidx1 , -in.getFloat (iidx1 + idx0 )); out.setFloat (oidx1 - 1, in.getFloat (iidx1 )); } } /*------------------------------------------------- radb2: Real FFT's backward processing of factor 2 -------------------------------------------------*/ void radb2 (final int ido , final int l1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { int i , ic ; float t1i , t1r , w1r , w1i ; int iw1 = offset ; int idx0 = l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 2 * idx1 ; int idx3 = idx2 + ido ; int oidx1 = out_off + idx1 ; int iidx1 = in_off + idx2 ; int iidx2 = in_off + ido - 1 + idx3 ; float i1r = in [iidx1 ]; float i2r = in [iidx2 ]; out [oidx1 ] = i1r + i2r ; out [oidx1 + idx0 ] = i1r - i2r ; } if (ido < 2) { return ; } if (ido != 2) { for (int k = 0; k < l1 ; ++k ) { int idx1 = k * ido ; int idx2 = 2 * idx1 ; int idx3 = idx2 + ido ; int idx4 = idx1 + idx0 ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; int idx5 = i - 1 + iw1 ; int idx6 = out_off + i ; int idx7 = in_off + i ; int idx8 = in_off + ic ; w1r = wtable_r [idx5 - 1]; w1i = wtable_r [idx5 ]; int iidx1 = idx7 + idx2 ; int iidx2 = idx8 + idx3 ; int oidx1 = idx6 + idx1 ; int oidx2 = idx6 + idx4 ; t1r = in [iidx1 - 1] - in [iidx2 - 1]; t1i = in [iidx1 ] + in [iidx2 ]; float i1i = in [iidx1 ]; float i1r = in [iidx1 - 1]; float i2i = in [iidx2 ]; float i2r = in [iidx2 - 1]; out [oidx1 - 1] = i1r + i2r ; out [oidx1 ] = i1i - i2i ; out [oidx2 - 1] = w1r * t1r - w1i * t1i ; out [oidx2 ] = w1r * t1i + w1i * t1r ; } } if (ido % 2 == 1) { return ; } } for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 2 * idx1 ; int oidx1 = out_off + ido - 1 + idx1 ; int iidx1 = in_off + idx2 + ido ; out [oidx1 ] = 2 * in [iidx1 - 1]; out [oidx1 + idx0 ] = -2 * in [iidx1 ]; } } /*------------------------------------------------- radb2: Real FFT's backward processing of factor 2 -------------------------------------------------*/ void radb2 (final long ido , final long l1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { long i , ic ; float t1i , t1r , w1r , w1i ; long iw1 = offset ; long idx0 = l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 2 * idx1 ; long idx3 = idx2 + ido ; long oidx1 = out_off + idx1 ; long iidx1 = in_off + idx2 ; long iidx2 = in_off + ido - 1 + idx3 ; float i1r = in.getFloat (iidx1 ); float i2r = in.getFloat (iidx2 ); out.setFloat (oidx1 , i1r + i2r ); out.setFloat (oidx1 + idx0 , i1r - i2r ); } if (ido < 2) { return ; } if (ido != 2) { for (long k = 0; k < l1 ; ++k ) { long idx1 = k * ido ; long idx2 = 2 * idx1 ; long idx3 = idx2 + ido ; long idx4 = idx1 + idx0 ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; long idx5 = i - 1 + iw1 ; long idx6 = out_off + i ; long idx7 = in_off + i ; long idx8 = in_off + ic ; w1r = wtable_rl.getFloat (idx5 - 1); w1i = wtable_rl.getFloat (idx5 ); long iidx1 = idx7 + idx2 ; long iidx2 = idx8 + idx3 ; long oidx1 = idx6 + idx1 ; long oidx2 = idx6 + idx4 ; t1r = in.getFloat (iidx1 - 1) - in.getFloat (iidx2 - 1); t1i = in.getFloat (iidx1 ) + in.getFloat (iidx2 ); float i1i = in.getFloat (iidx1 ); float i1r = in.getFloat (iidx1 - 1); float i2i = in.getFloat (iidx2 ); float i2r = in.getFloat (iidx2 - 1); out.setFloat (oidx1 - 1, i1r + i2r ); out.setFloat (oidx1 , i1i - i2i ); out.setFloat (oidx2 - 1, w1r * t1r - w1i * t1i ); out.setFloat (oidx2 , w1r * t1i + w1i * t1r ); } } if (ido % 2 == 1) { return ; } } for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 2 * idx1 ; long oidx1 = out_off + ido - 1 + idx1 ; long iidx1 = in_off + idx2 + ido ; out.setFloat (oidx1 , 2 * in.getFloat (iidx1 - 1)); out.setFloat (oidx1 + idx0 , -2 * in.getFloat (iidx1 )); } } /*------------------------------------------------- radf3: Real FFT's forward processing of factor 3 -------------------------------------------------*/ void radf3 (final int ido , final int l1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { final float taur = -0.5f; final float taui = 0.866025403784438707610604524234076962f; int i , ic ; float ci2 , di2 , di3 , cr2 , dr2 , dr3 , ti2 , ti3 , tr2 , tr3 , w1r , w2r , w1i , w2i ; int iw1 , iw2 ; iw1 = offset ; iw2 = iw1 + ido ; int idx0 = l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx3 = 2 * idx0 ; int idx4 = (3 * k + 1) * ido ; int iidx1 = in_off + idx1 ; int iidx2 = iidx1 + idx0 ; int iidx3 = iidx1 + idx3 ; float i1r = in [iidx1 ]; float i2r = in [iidx2 ]; float i3r = in [iidx3 ]; cr2 = i2r + i3r ; out [out_off + 3 * idx1 ] = i1r + cr2 ; out [out_off + idx4 + ido ] = taui * (i3r - i2r ); out [out_off + ido - 1 + idx4 ] = i1r + taur * cr2 ; } if (ido == 1) { return ; } for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido ; int idx4 = 3 * idx3 ; int idx5 = idx3 + idx0 ; int idx6 = idx5 + idx0 ; int idx7 = idx4 + ido ; int idx8 = idx7 + ido ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; int widx1 = i - 1 + iw1 ; int widx2 = i - 1 + iw2 ; w1r = wtable_r [widx1 - 1]; w1i = wtable_r [widx1 ]; w2r = wtable_r [widx2 - 1]; w2i = wtable_r [widx2 ]; int idx9 = in_off + i ; int idx10 = out_off + i ; int idx11 = out_off + ic ; int iidx1 = idx9 + idx3 ; int iidx2 = idx9 + idx5 ; int iidx3 = idx9 + idx6 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; float i2i = in [iidx2 - 1]; float i2r = in [iidx2 ]; float i3i = in [iidx3 - 1]; float i3r = in [iidx3 ]; dr2 = w1r * i2i + w1i * i2r ; di2 = w1r * i2r - w1i * i2i ; dr3 = w2r * i3i + w2i * i3r ; di3 = w2r * i3r - w2i * i3i ; cr2 = dr2 + dr3 ; ci2 = di2 + di3 ; tr2 = i1i + taur * cr2 ; ti2 = i1r + taur * ci2 ; tr3 = taui * (di2 - di3 ); ti3 = taui * (dr3 - dr2 ); int oidx1 = idx10 + idx4 ; int oidx2 = idx11 + idx7 ; int oidx3 = idx10 + idx8 ; out [oidx1 - 1] = i1i + cr2 ; out [oidx1 ] = i1r + ci2 ; out [oidx2 - 1] = tr2 - tr3 ; out [oidx2 ] = ti3 - ti2 ; out [oidx3 - 1] = tr2 + tr3 ; out [oidx3 ] = ti2 + ti3 ; } } } /*------------------------------------------------- radf3: Real FFT's forward processing of factor 3 -------------------------------------------------*/ void radf3 (final long ido , final long l1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { final float taur = -0.5f; final float taui = 0.866025403784438707610604524234076962f; long i , ic ; float ci2 , di2 , di3 , cr2 , dr2 , dr3 , ti2 , ti3 , tr2 , tr3 , w1r , w2r , w1i , w2i ; long iw1 , iw2 ; iw1 = offset ; iw2 = iw1 + ido ; long idx0 = l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx3 = 2 * idx0 ; long idx4 = (3 * k + 1) * ido ; long iidx1 = in_off + idx1 ; long iidx2 = iidx1 + idx0 ; long iidx3 = iidx1 + idx3 ; float i1r = in.getFloat (iidx1 ); float i2r = in.getFloat (iidx2 ); float i3r = in.getFloat (iidx3 ); cr2 = i2r + i3r ; out.setFloat (out_off + 3 * idx1 , i1r + cr2 ); out.setFloat (out_off + idx4 + ido , taui * (i3r - i2r )); out.setFloat (out_off + ido - 1 + idx4 , i1r + taur * cr2 ); } if (ido == 1) { return ; } for (long k = 0; k < l1 ; k ++) { long idx3 = k * ido ; long idx4 = 3 * idx3 ; long idx5 = idx3 + idx0 ; long idx6 = idx5 + idx0 ; long idx7 = idx4 + ido ; long idx8 = idx7 + ido ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; long widx1 = i - 1 + iw1 ; long widx2 = i - 1 + iw2 ; w1r = wtable_rl.getFloat (widx1 - 1); w1i = wtable_rl.getFloat (widx1 ); w2r = wtable_rl.getFloat (widx2 - 1); w2i = wtable_rl.getFloat (widx2 ); long idx9 = in_off + i ; long idx10 = out_off + i ; long idx11 = out_off + ic ; long iidx1 = idx9 + idx3 ; long iidx2 = idx9 + idx5 ; long iidx3 = idx9 + idx6 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); float i2i = in.getFloat (iidx2 - 1); float i2r = in.getFloat (iidx2 ); float i3i = in.getFloat (iidx3 - 1); float i3r = in.getFloat (iidx3 ); dr2 = w1r * i2i + w1i * i2r ; di2 = w1r * i2r - w1i * i2i ; dr3 = w2r * i3i + w2i * i3r ; di3 = w2r * i3r - w2i * i3i ; cr2 = dr2 + dr3 ; ci2 = di2 + di3 ; tr2 = i1i + taur * cr2 ; ti2 = i1r + taur * ci2 ; tr3 = taui * (di2 - di3 ); ti3 = taui * (dr3 - dr2 ); long oidx1 = idx10 + idx4 ; long oidx2 = idx11 + idx7 ; long oidx3 = idx10 + idx8 ; out.setFloat (oidx1 - 1, i1i + cr2 ); out.setFloat (oidx1 , i1r + ci2 ); out.setFloat (oidx2 - 1, tr2 - tr3 ); out.setFloat (oidx2 , ti3 - ti2 ); out.setFloat (oidx3 - 1, tr2 + tr3 ); out.setFloat (oidx3 , ti2 + ti3 ); } } } /*------------------------------------------------- radb3: Real FFT's backward processing of factor 3 -------------------------------------------------*/ void radb3 (final int ido , final int l1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { final float taur = -0.5f; final float taui = 0.866025403784438707610604524234076962f; int i , ic ; float ci2 , ci3 , di2 , di3 , cr2 , cr3 , dr2 , dr3 , ti2 , tr2 , w1r , w2r , w1i , w2i ; int iw1 , iw2 ; iw1 = offset ; iw2 = iw1 + ido ; for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int iidx1 = in_off + 3 * idx1 ; int iidx2 = iidx1 + 2 * ido ; float i1i = in [iidx1 ]; tr2 = 2 * in [iidx2 - 1]; cr2 = i1i + taur * tr2 ; ci3 = 2 * taui * in [iidx2 ]; out [out_off + idx1 ] = i1i + tr2 ; out [out_off + (k + l1 ) * ido ] = cr2 - ci3 ; out [out_off + (k + 2 * l1 ) * ido ] = cr2 + ci3 ; } if (ido == 1) { return ; } int idx0 = l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 3 * idx1 ; int idx3 = idx2 + ido ; int idx4 = idx3 + ido ; int idx5 = idx1 + idx0 ; int idx6 = idx5 + idx0 ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; int idx7 = in_off + i ; int idx8 = in_off + ic ; int idx9 = out_off + i ; int iidx1 = idx7 + idx2 ; int iidx2 = idx7 + idx4 ; int iidx3 = idx8 + idx3 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; float i2i = in [iidx2 - 1]; float i2r = in [iidx2 ]; float i3i = in [iidx3 - 1]; float i3r = in [iidx3 ]; tr2 = i2i + i3i ; cr2 = i1i + taur * tr2 ; ti2 = i2r - i3r ; ci2 = i1r + taur * ti2 ; cr3 = taui * (i2i - i3i ); ci3 = taui * (i2r + i3r ); dr2 = cr2 - ci3 ; dr3 = cr2 + ci3 ; di2 = ci2 + cr3 ; di3 = ci2 - cr3 ; int widx1 = i - 1 + iw1 ; int widx2 = i - 1 + iw2 ; w1r = wtable_r [widx1 - 1]; w1i = wtable_r [widx1 ]; w2r = wtable_r [widx2 - 1]; w2i = wtable_r [widx2 ]; int oidx1 = idx9 + idx1 ; int oidx2 = idx9 + idx5 ; int oidx3 = idx9 + idx6 ; out [oidx1 - 1] = i1i + tr2 ; out [oidx1 ] = i1r + ti2 ; out [oidx2 - 1] = w1r * dr2 - w1i * di2 ; out [oidx2 ] = w1r * di2 + w1i * dr2 ; out [oidx3 - 1] = w2r * dr3 - w2i * di3 ; out [oidx3 ] = w2r * di3 + w2i * dr3 ; } } } /*------------------------------------------------- radb3: Real FFT's backward processing of factor 3 -------------------------------------------------*/ void radb3 (final long ido , final long l1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { final float taur = -0.5f; final float taui = 0.866025403784438707610604524234076962f; long i , ic ; float ci2 , ci3 , di2 , di3 , cr2 , cr3 , dr2 , dr3 , ti2 , tr2 , w1r , w2r , w1i , w2i ; long iw1 , iw2 ; iw1 = offset ; iw2 = iw1 + ido ; for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long iidx1 = in_off + 3 * idx1 ; long iidx2 = iidx1 + 2 * ido ; float i1i = in.getFloat (iidx1 ); tr2 = 2 * in.getFloat (iidx2 - 1); cr2 = i1i + taur * tr2 ; ci3 = 2 * taui * in.getFloat (iidx2 ); out.setFloat (out_off + idx1 , i1i + tr2 ); out.setFloat (out_off + (k + l1 ) * ido , cr2 - ci3 ); out.setFloat (out_off + (k + 2 * l1 ) * ido , cr2 + ci3 ); } if (ido == 1) { return ; } long idx0 = l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 3 * idx1 ; long idx3 = idx2 + ido ; long idx4 = idx3 + ido ; long idx5 = idx1 + idx0 ; long idx6 = idx5 + idx0 ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; long idx7 = in_off + i ; long idx8 = in_off + ic ; long idx9 = out_off + i ; long iidx1 = idx7 + idx2 ; long iidx2 = idx7 + idx4 ; long iidx3 = idx8 + idx3 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); float i2i = in.getFloat (iidx2 - 1); float i2r = in.getFloat (iidx2 ); float i3i = in.getFloat (iidx3 - 1); float i3r = in.getFloat (iidx3 ); tr2 = i2i + i3i ; cr2 = i1i + taur * tr2 ; ti2 = i2r - i3r ; ci2 = i1r + taur * ti2 ; cr3 = taui * (i2i - i3i ); ci3 = taui * (i2r + i3r ); dr2 = cr2 - ci3 ; dr3 = cr2 + ci3 ; di2 = ci2 + cr3 ; di3 = ci2 - cr3 ; long widx1 = i - 1 + iw1 ; long widx2 = i - 1 + iw2 ; w1r = wtable_rl.getFloat (widx1 - 1); w1i = wtable_rl.getFloat (widx1 ); w2r = wtable_rl.getFloat (widx2 - 1); w2i = wtable_rl.getFloat (widx2 ); long oidx1 = idx9 + idx1 ; long oidx2 = idx9 + idx5 ; long oidx3 = idx9 + idx6 ; out.setFloat (oidx1 - 1, i1i + tr2 ); out.setFloat (oidx1 , i1r + ti2 ); out.setFloat (oidx2 - 1, w1r * dr2 - w1i * di2 ); out.setFloat (oidx2 , w1r * di2 + w1i * dr2 ); out.setFloat (oidx3 - 1, w2r * dr3 - w2i * di3 ); out.setFloat (oidx3 , w2r * di3 + w2i * dr3 ); } } } /*------------------------------------------------- radf4: Real FFT's forward processing of factor 4 -------------------------------------------------*/ void radf4 (final int ido , final int l1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { final float hsqt2 = 0.707106781186547572737310929369414225f; int i , ic ; float ci2 , ci3 , ci4 , cr2 , cr3 , cr4 , ti1 , ti2 , ti3 , ti4 , tr1 , tr2 , tr3 , tr4 , w1r , w1i , w2r , w2i , w3r , w3i ; int iw1 , iw2 , iw3 ; iw1 = offset ; iw2 = offset + ido ; iw3 = iw2 + ido ; int idx0 = l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 4 * idx1 ; int idx3 = idx1 + idx0 ; int idx4 = idx3 + idx0 ; int idx5 = idx4 + idx0 ; int idx6 = idx2 + ido ; float i1r = in [in_off + idx1 ]; float i2r = in [in_off + idx3 ]; float i3r = in [in_off + idx4 ]; float i4r = in [in_off + idx5 ]; tr1 = i2r + i4r ; tr2 = i1r + i3r ; int oidx1 = out_off + idx2 ; int oidx2 = out_off + idx6 + ido ; out [oidx1 ] = tr1 + tr2 ; out [oidx2 - 1 + ido + ido ] = tr2 - tr1 ; out [oidx2 - 1] = i1r - i3r ; out [oidx2 ] = i4r - i2r ; } if (ido < 2) { return ; } if (ido != 2) { for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = idx1 + idx0 ; int idx3 = idx2 + idx0 ; int idx4 = idx3 + idx0 ; int idx5 = 4 * idx1 ; int idx6 = idx5 + ido ; int idx7 = idx6 + ido ; int idx8 = idx7 + ido ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; int widx1 = i - 1 + iw1 ; int widx2 = i - 1 + iw2 ; int widx3 = i - 1 + iw3 ; w1r = wtable_r [widx1 - 1]; w1i = wtable_r [widx1 ]; w2r = wtable_r [widx2 - 1]; w2i = wtable_r [widx2 ]; w3r = wtable_r [widx3 - 1]; w3i = wtable_r [widx3 ]; int idx9 = in_off + i ; int idx10 = out_off + i ; int idx11 = out_off + ic ; int iidx1 = idx9 + idx1 ; int iidx2 = idx9 + idx2 ; int iidx3 = idx9 + idx3 ; int iidx4 = idx9 + idx4 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; float i2i = in [iidx2 - 1]; float i2r = in [iidx2 ]; float i3i = in [iidx3 - 1]; float i3r = in [iidx3 ]; float i4i = in [iidx4 - 1]; float i4r = in [iidx4 ]; cr2 = w1r * i2i + w1i * i2r ; ci2 = w1r * i2r - w1i * i2i ; cr3 = w2r * i3i + w2i * i3r ; ci3 = w2r * i3r - w2i * i3i ; cr4 = w3r * i4i + w3i * i4r ; ci4 = w3r * i4r - w3i * i4i ; tr1 = cr2 + cr4 ; tr4 = cr4 - cr2 ; ti1 = ci2 + ci4 ; ti4 = ci2 - ci4 ; ti2 = i1r + ci3 ; ti3 = i1r - ci3 ; tr2 = i1i + cr3 ; tr3 = i1i - cr3 ; int oidx1 = idx10 + idx5 ; int oidx2 = idx11 + idx6 ; int oidx3 = idx10 + idx7 ; int oidx4 = idx11 + idx8 ; out [oidx1 - 1] = tr1 + tr2 ; out [oidx4 - 1] = tr2 - tr1 ; out [oidx1 ] = ti1 + ti2 ; out [oidx4 ] = ti1 - ti2 ; out [oidx3 - 1] = ti4 + tr3 ; out [oidx2 - 1] = tr3 - ti4 ; out [oidx3 ] = tr4 + ti3 ; out [oidx2 ] = tr4 - ti3 ; } } if (ido % 2 == 1) { return ; } } for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 4 * idx1 ; int idx3 = idx1 + idx0 ; int idx4 = idx3 + idx0 ; int idx5 = idx4 + idx0 ; int idx6 = idx2 + ido ; int idx7 = idx6 + ido ; int idx8 = idx7 + ido ; int idx9 = in_off + ido ; int idx10 = out_off + ido ; float i1i = in [idx9 - 1 + idx1 ]; float i2i = in [idx9 - 1 + idx3 ]; float i3i = in [idx9 - 1 + idx4 ]; float i4i = in [idx9 - 1 + idx5 ]; ti1 = -hsqt2 * (i2i + i4i ); tr1 = hsqt2 * (i2i - i4i ); out [idx10 - 1 + idx2 ] = tr1 + i1i ; out [idx10 - 1 + idx7 ] = i1i - tr1 ; out [out_off + idx6 ] = ti1 - i3i ; out [out_off + idx8 ] = ti1 + i3i ; } } /*------------------------------------------------- radf4: Real FFT's forward processing of factor 4 -------------------------------------------------*/ void radf4 (final long ido , final long l1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { final float hsqt2 = 0.707106781186547572737310929369414225f; long i , ic ; float ci2 , ci3 , ci4 , cr2 , cr3 , cr4 , ti1 , ti2 , ti3 , ti4 , tr1 , tr2 , tr3 , tr4 , w1r , w1i , w2r , w2i , w3r , w3i ; long iw1 , iw2 , iw3 ; iw1 = offset ; iw2 = offset + ido ; iw3 = iw2 + ido ; long idx0 = l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 4 * idx1 ; long idx3 = idx1 + idx0 ; long idx4 = idx3 + idx0 ; long idx5 = idx4 + idx0 ; long idx6 = idx2 + ido ; float i1r = in.getFloat (in_off + idx1 ); float i2r = in.getFloat (in_off + idx3 ); float i3r = in.getFloat (in_off + idx4 ); float i4r = in.getFloat (in_off + idx5 ); tr1 = i2r + i4r ; tr2 = i1r + i3r ; long oidx1 = out_off + idx2 ; long oidx2 = out_off + idx6 + ido ; out.setFloat (oidx1 , tr1 + tr2 ); out.setFloat (oidx2 - 1 + ido + ido , tr2 - tr1 ); out.setFloat (oidx2 - 1, i1r - i3r ); out.setFloat (oidx2 , i4r - i2r ); } if (ido < 2) { return ; } if (ido != 2) { for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = idx1 + idx0 ; long idx3 = idx2 + idx0 ; long idx4 = idx3 + idx0 ; long idx5 = 4 * idx1 ; long idx6 = idx5 + ido ; long idx7 = idx6 + ido ; long idx8 = idx7 + ido ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; long widx1 = i - 1 + iw1 ; long widx2 = i - 1 + iw2 ; long widx3 = i - 1 + iw3 ; w1r = wtable_rl.getFloat (widx1 - 1); w1i = wtable_rl.getFloat (widx1 ); w2r = wtable_rl.getFloat (widx2 - 1); w2i = wtable_rl.getFloat (widx2 ); w3r = wtable_rl.getFloat (widx3 - 1); w3i = wtable_rl.getFloat (widx3 ); long idx9 = in_off + i ; long idx10 = out_off + i ; long idx11 = out_off + ic ; long iidx1 = idx9 + idx1 ; long iidx2 = idx9 + idx2 ; long iidx3 = idx9 + idx3 ; long iidx4 = idx9 + idx4 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); float i2i = in.getFloat (iidx2 - 1); float i2r = in.getFloat (iidx2 ); float i3i = in.getFloat (iidx3 - 1); float i3r = in.getFloat (iidx3 ); float i4i = in.getFloat (iidx4 - 1); float i4r = in.getFloat (iidx4 ); cr2 = w1r * i2i + w1i * i2r ; ci2 = w1r * i2r - w1i * i2i ; cr3 = w2r * i3i + w2i * i3r ; ci3 = w2r * i3r - w2i * i3i ; cr4 = w3r * i4i + w3i * i4r ; ci4 = w3r * i4r - w3i * i4i ; tr1 = cr2 + cr4 ; tr4 = cr4 - cr2 ; ti1 = ci2 + ci4 ; ti4 = ci2 - ci4 ; ti2 = i1r + ci3 ; ti3 = i1r - ci3 ; tr2 = i1i + cr3 ; tr3 = i1i - cr3 ; long oidx1 = idx10 + idx5 ; long oidx2 = idx11 + idx6 ; long oidx3 = idx10 + idx7 ; long oidx4 = idx11 + idx8 ; out.setFloat (oidx1 - 1, tr1 + tr2 ); out.setFloat (oidx4 - 1, tr2 - tr1 ); out.setFloat (oidx1 , ti1 + ti2 ); out.setFloat (oidx4 , ti1 - ti2 ); out.setFloat (oidx3 - 1, ti4 + tr3 ); out.setFloat (oidx2 - 1, tr3 - ti4 ); out.setFloat (oidx3 , tr4 + ti3 ); out.setFloat (oidx2 , tr4 - ti3 ); } } if (ido % 2 == 1) { return ; } } for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 4 * idx1 ; long idx3 = idx1 + idx0 ; long idx4 = idx3 + idx0 ; long idx5 = idx4 + idx0 ; long idx6 = idx2 + ido ; long idx7 = idx6 + ido ; long idx8 = idx7 + ido ; long idx9 = in_off + ido ; long idx10 = out_off + ido ; float i1i = in.getFloat (idx9 - 1 + idx1 ); float i2i = in.getFloat (idx9 - 1 + idx3 ); float i3i = in.getFloat (idx9 - 1 + idx4 ); float i4i = in.getFloat (idx9 - 1 + idx5 ); ti1 = -hsqt2 * (i2i + i4i ); tr1 = hsqt2 * (i2i - i4i ); out.setFloat (idx10 - 1 + idx2 , tr1 + i1i ); out.setFloat (idx10 - 1 + idx7 , i1i - tr1 ); out.setFloat (out_off + idx6 , ti1 - i3i ); out.setFloat (out_off + idx8 , ti1 + i3i ); } } /*------------------------------------------------- radb4: Real FFT's backward processing of factor 4 -------------------------------------------------*/ void radb4 (final int ido , final int l1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { final float sqrt2 = 1.41421356237309514547462185873882845f; int i , ic ; float ci2 , ci3 , ci4 , cr2 , cr3 , cr4 ; float ti1 , ti2 , ti3 , ti4 , tr1 , tr2 , tr3 , tr4 , w1r , w1i , w2r , w2i , w3r , w3i ; int iw1 , iw2 , iw3 ; iw1 = offset ; iw2 = iw1 + ido ; iw3 = iw2 + ido ; int idx0 = l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 4 * idx1 ; int idx3 = idx1 + idx0 ; int idx4 = idx3 + idx0 ; int idx5 = idx4 + idx0 ; int idx6 = idx2 + ido ; int idx7 = idx6 + ido ; int idx8 = idx7 + ido ; float i1r = in [in_off + idx2 ]; float i2r = in [in_off + idx7 ]; float i3r = in [in_off + ido - 1 + idx8 ]; float i4r = in [in_off + ido - 1 + idx6 ]; tr1 = i1r - i3r ; tr2 = i1r + i3r ; tr3 = i4r + i4r ; tr4 = i2r + i2r ; out [out_off + idx1 ] = tr2 + tr3 ; out [out_off + idx3 ] = tr1 - tr4 ; out [out_off + idx4 ] = tr2 - tr3 ; out [out_off + idx5 ] = tr1 + tr4 ; } if (ido < 2) { return ; } if (ido != 2) { for (int k = 0; k < l1 ; ++k ) { int idx1 = k * ido ; int idx2 = idx1 + idx0 ; int idx3 = idx2 + idx0 ; int idx4 = idx3 + idx0 ; int idx5 = 4 * idx1 ; int idx6 = idx5 + ido ; int idx7 = idx6 + ido ; int idx8 = idx7 + ido ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; int widx1 = i - 1 + iw1 ; int widx2 = i - 1 + iw2 ; int widx3 = i - 1 + iw3 ; w1r = wtable_r [widx1 - 1]; w1i = wtable_r [widx1 ]; w2r = wtable_r [widx2 - 1]; w2i = wtable_r [widx2 ]; w3r = wtable_r [widx3 - 1]; w3i = wtable_r [widx3 ]; int idx12 = in_off + i ; int idx13 = in_off + ic ; int idx14 = out_off + i ; int iidx1 = idx12 + idx5 ; int iidx2 = idx13 + idx6 ; int iidx3 = idx12 + idx7 ; int iidx4 = idx13 + idx8 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; float i2i = in [iidx2 - 1]; float i2r = in [iidx2 ]; float i3i = in [iidx3 - 1]; float i3r = in [iidx3 ]; float i4i = in [iidx4 - 1]; float i4r = in [iidx4 ]; ti1 = i1r + i4r ; ti2 = i1r - i4r ; ti3 = i3r - i2r ; tr4 = i3r + i2r ; tr1 = i1i - i4i ; tr2 = i1i + i4i ; ti4 = i3i - i2i ; tr3 = i3i + i2i ; cr3 = tr2 - tr3 ; ci3 = ti2 - ti3 ; cr2 = tr1 - tr4 ; cr4 = tr1 + tr4 ; ci2 = ti1 + ti4 ; ci4 = ti1 - ti4 ; int oidx1 = idx14 + idx1 ; int oidx2 = idx14 + idx2 ; int oidx3 = idx14 + idx3 ; int oidx4 = idx14 + idx4 ; out [oidx1 - 1] = tr2 + tr3 ; out [oidx1 ] = ti2 + ti3 ; out [oidx2 - 1] = w1r * cr2 - w1i * ci2 ; out [oidx2 ] = w1r * ci2 + w1i * cr2 ; out [oidx3 - 1] = w2r * cr3 - w2i * ci3 ; out [oidx3 ] = w2r * ci3 + w2i * cr3 ; out [oidx4 - 1] = w3r * cr4 - w3i * ci4 ; out [oidx4 ] = w3r * ci4 + w3i * cr4 ; } } if (ido % 2 == 1) { return ; } } for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 4 * idx1 ; int idx3 = idx1 + idx0 ; int idx4 = idx3 + idx0 ; int idx5 = idx4 + idx0 ; int idx6 = idx2 + ido ; int idx7 = idx6 + ido ; int idx8 = idx7 + ido ; int idx9 = in_off + ido ; int idx10 = out_off + ido ; float i1r = in [idx9 - 1 + idx2 ]; float i2r = in [idx9 - 1 + idx7 ]; float i3r = in [in_off + idx6 ]; float i4r = in [in_off + idx8 ]; ti1 = i3r + i4r ; ti2 = i4r - i3r ; tr1 = i1r - i2r ; tr2 = i1r + i2r ; out [idx10 - 1 + idx1 ] = tr2 + tr2 ; out [idx10 - 1 + idx3 ] = sqrt2 * (tr1 - ti1 ); out [idx10 - 1 + idx4 ] = ti2 + ti2 ; out [idx10 - 1 + idx5 ] = -sqrt2 * (tr1 + ti1 ); } } /*------------------------------------------------- radb4: Real FFT's backward processing of factor 4 -------------------------------------------------*/ void radb4 (final long ido , final long l1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { final float sqrt2 = 1.41421356237309514547462185873882845f; long i , ic ; float ci2 , ci3 , ci4 , cr2 , cr3 , cr4 ; float ti1 , ti2 , ti3 , ti4 , tr1 , tr2 , tr3 , tr4 , w1r , w1i , w2r , w2i , w3r , w3i ; long iw1 , iw2 , iw3 ; iw1 = offset ; iw2 = iw1 + ido ; iw3 = iw2 + ido ; long idx0 = l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 4 * idx1 ; long idx3 = idx1 + idx0 ; long idx4 = idx3 + idx0 ; long idx5 = idx4 + idx0 ; long idx6 = idx2 + ido ; long idx7 = idx6 + ido ; long idx8 = idx7 + ido ; float i1r = in.getFloat (in_off + idx2 ); float i2r = in.getFloat (in_off + idx7 ); float i3r = in.getFloat (in_off + ido - 1 + idx8 ); float i4r = in.getFloat (in_off + ido - 1 + idx6 ); tr1 = i1r - i3r ; tr2 = i1r + i3r ; tr3 = i4r + i4r ; tr4 = i2r + i2r ; out.setFloat (out_off + idx1 , tr2 + tr3 ); out.setFloat (out_off + idx3 , tr1 - tr4 ); out.setFloat (out_off + idx4 , tr2 - tr3 ); out.setFloat (out_off + idx5 , tr1 + tr4 ); } if (ido < 2) { return ; } if (ido != 2) { for (long k = 0; k < l1 ; ++k ) { long idx1 = k * ido ; long idx2 = idx1 + idx0 ; long idx3 = idx2 + idx0 ; long idx4 = idx3 + idx0 ; long idx5 = 4 * idx1 ; long idx6 = idx5 + ido ; long idx7 = idx6 + ido ; long idx8 = idx7 + ido ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; long widx1 = i - 1 + iw1 ; long widx2 = i - 1 + iw2 ; long widx3 = i - 1 + iw3 ; w1r = wtable_rl.getFloat (widx1 - 1); w1i = wtable_rl.getFloat (widx1 ); w2r = wtable_rl.getFloat (widx2 - 1); w2i = wtable_rl.getFloat (widx2 ); w3r = wtable_rl.getFloat (widx3 - 1); w3i = wtable_rl.getFloat (widx3 ); long idx12 = in_off + i ; long idx13 = in_off + ic ; long idx14 = out_off + i ; long iidx1 = idx12 + idx5 ; long iidx2 = idx13 + idx6 ; long iidx3 = idx12 + idx7 ; long iidx4 = idx13 + idx8 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); float i2i = in.getFloat (iidx2 - 1); float i2r = in.getFloat (iidx2 ); float i3i = in.getFloat (iidx3 - 1); float i3r = in.getFloat (iidx3 ); float i4i = in.getFloat (iidx4 - 1); float i4r = in.getFloat (iidx4 ); ti1 = i1r + i4r ; ti2 = i1r - i4r ; ti3 = i3r - i2r ; tr4 = i3r + i2r ; tr1 = i1i - i4i ; tr2 = i1i + i4i ; ti4 = i3i - i2i ; tr3 = i3i + i2i ; cr3 = tr2 - tr3 ; ci3 = ti2 - ti3 ; cr2 = tr1 - tr4 ; cr4 = tr1 + tr4 ; ci2 = ti1 + ti4 ; ci4 = ti1 - ti4 ; long oidx1 = idx14 + idx1 ; long oidx2 = idx14 + idx2 ; long oidx3 = idx14 + idx3 ; long oidx4 = idx14 + idx4 ; out.setFloat (oidx1 - 1, tr2 + tr3 ); out.setFloat (oidx1 , ti2 + ti3 ); out.setFloat (oidx2 - 1, w1r * cr2 - w1i * ci2 ); out.setFloat (oidx2 , w1r * ci2 + w1i * cr2 ); out.setFloat (oidx3 - 1, w2r * cr3 - w2i * ci3 ); out.setFloat (oidx3 , w2r * ci3 + w2i * cr3 ); out.setFloat (oidx4 - 1, w3r * cr4 - w3i * ci4 ); out.setFloat (oidx4 , w3r * ci4 + w3i * cr4 ); } } if (ido % 2 == 1) { return ; } } for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 4 * idx1 ; long idx3 = idx1 + idx0 ; long idx4 = idx3 + idx0 ; long idx5 = idx4 + idx0 ; long idx6 = idx2 + ido ; long idx7 = idx6 + ido ; long idx8 = idx7 + ido ; long idx9 = in_off + ido ; long idx10 = out_off + ido ; float i1r = in.getFloat (idx9 - 1 + idx2 ); float i2r = in.getFloat (idx9 - 1 + idx7 ); float i3r = in.getFloat (in_off + idx6 ); float i4r = in.getFloat (in_off + idx8 ); ti1 = i3r + i4r ; ti2 = i4r - i3r ; tr1 = i1r - i2r ; tr2 = i1r + i2r ; out.setFloat (idx10 - 1 + idx1 , tr2 + tr2 ); out.setFloat (idx10 - 1 + idx3 , sqrt2 * (tr1 - ti1 )); out.setFloat (idx10 - 1 + idx4 , ti2 + ti2 ); out.setFloat (idx10 - 1 + idx5 , -sqrt2 * (tr1 + ti1 )); } } /*------------------------------------------------- radf5: Real FFT's forward processing of factor 5 -------------------------------------------------*/ void radf5 (final int ido , final int l1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { final float tr11 = 0.309016994374947451262869435595348477f; final float ti11 = 0.951056516295153531181938433292089030f; final float tr12 = -0.809016994374947340240566973079694435f; final float ti12 = 0.587785252292473248125759255344746634f; int i , ic ; float ci2 , di2 , ci4 , ci5 , di3 , di4 , di5 , ci3 , cr2 , cr3 , dr2 , dr3 , dr4 , dr5 , cr5 , cr4 , ti2 , ti3 , ti5 , ti4 , tr2 , tr3 , tr4 , tr5 , w1r , w1i , w2r , w2i , w3r , w3i , w4r , w4i ; int iw1 , iw2 , iw3 , iw4 ; iw1 = offset ; iw2 = iw1 + ido ; iw3 = iw2 + ido ; iw4 = iw3 + ido ; int idx0 = l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 5 * idx1 ; int idx3 = idx2 + ido ; int idx4 = idx3 + ido ; int idx5 = idx4 + ido ; int idx6 = idx5 + ido ; int idx7 = idx1 + idx0 ; int idx8 = idx7 + idx0 ; int idx9 = idx8 + idx0 ; int idx10 = idx9 + idx0 ; int idx11 = out_off + ido - 1; float i1r = in [in_off + idx1 ]; float i2r = in [in_off + idx7 ]; float i3r = in [in_off + idx8 ]; float i4r = in [in_off + idx9 ]; float i5r = in [in_off + idx10 ]; cr2 = i5r + i2r ; ci5 = i5r - i2r ; cr3 = i4r + i3r ; ci4 = i4r - i3r ; out [out_off + idx2 ] = i1r + cr2 + cr3 ; out [idx11 + idx3 ] = i1r + tr11 * cr2 + tr12 * cr3 ; out [out_off + idx4 ] = ti11 * ci5 + ti12 * ci4 ; out [idx11 + idx5 ] = i1r + tr12 * cr2 + tr11 * cr3 ; out [out_off + idx6 ] = ti12 * ci5 - ti11 * ci4 ; } if (ido == 1) { return ; } for (int k = 0; k < l1 ; ++k ) { int idx1 = k * ido ; int idx2 = 5 * idx1 ; int idx3 = idx2 + ido ; int idx4 = idx3 + ido ; int idx5 = idx4 + ido ; int idx6 = idx5 + ido ; int idx7 = idx1 + idx0 ; int idx8 = idx7 + idx0 ; int idx9 = idx8 + idx0 ; int idx10 = idx9 + idx0 ; for (i = 2; i < ido ; i += 2) { int widx1 = i - 1 + iw1 ; int widx2 = i - 1 + iw2 ; int widx3 = i - 1 + iw3 ; int widx4 = i - 1 + iw4 ; w1r = wtable_r [widx1 - 1]; w1i = wtable_r [widx1 ]; w2r = wtable_r [widx2 - 1]; w2i = wtable_r [widx2 ]; w3r = wtable_r [widx3 - 1]; w3i = wtable_r [widx3 ]; w4r = wtable_r [widx4 - 1]; w4i = wtable_r [widx4 ]; ic = ido - i ; int idx15 = in_off + i ; int idx16 = out_off + i ; int idx17 = out_off + ic ; int iidx1 = idx15 + idx1 ; int iidx2 = idx15 + idx7 ; int iidx3 = idx15 + idx8 ; int iidx4 = idx15 + idx9 ; int iidx5 = idx15 + idx10 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; float i2i = in [iidx2 - 1]; float i2r = in [iidx2 ]; float i3i = in [iidx3 - 1]; float i3r = in [iidx3 ]; float i4i = in [iidx4 - 1]; float i4r = in [iidx4 ]; float i5i = in [iidx5 - 1]; float i5r = in [iidx5 ]; dr2 = w1r * i2i + w1i * i2r ; di2 = w1r * i2r - w1i * i2i ; dr3 = w2r * i3i + w2i * i3r ; di3 = w2r * i3r - w2i * i3i ; dr4 = w3r * i4i + w3i * i4r ; di4 = w3r * i4r - w3i * i4i ; dr5 = w4r * i5i + w4i * i5r ; di5 = w4r * i5r - w4i * i5i ; cr2 = dr2 + dr5 ; ci5 = dr5 - dr2 ; cr5 = di2 - di5 ; ci2 = di2 + di5 ; cr3 = dr3 + dr4 ; ci4 = dr4 - dr3 ; cr4 = di3 - di4 ; ci3 = di3 + di4 ; tr2 = i1i + tr11 * cr2 + tr12 * cr3 ; ti2 = i1r + tr11 * ci2 + tr12 * ci3 ; tr3 = i1i + tr12 * cr2 + tr11 * cr3 ; ti3 = i1r + tr12 * ci2 + tr11 * ci3 ; tr5 = ti11 * cr5 + ti12 * cr4 ; ti5 = ti11 * ci5 + ti12 * ci4 ; tr4 = ti12 * cr5 - ti11 * cr4 ; ti4 = ti12 * ci5 - ti11 * ci4 ; int oidx1 = idx16 + idx2 ; int oidx2 = idx17 + idx3 ; int oidx3 = idx16 + idx4 ; int oidx4 = idx17 + idx5 ; int oidx5 = idx16 + idx6 ; out [oidx1 - 1] = i1i + cr2 + cr3 ; out [oidx1 ] = i1r + ci2 + ci3 ; out [oidx3 - 1] = tr2 + tr5 ; out [oidx2 - 1] = tr2 - tr5 ; out [oidx3 ] = ti2 + ti5 ; out [oidx2 ] = ti5 - ti2 ; out [oidx5 - 1] = tr3 + tr4 ; out [oidx4 - 1] = tr3 - tr4 ; out [oidx5 ] = ti3 + ti4 ; out [oidx4 ] = ti4 - ti3 ; } } } /*------------------------------------------------- radf5: Real FFT's forward processing of factor 5 -------------------------------------------------*/ void radf5 (final long ido , final long l1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { final float tr11 = 0.309016994374947451262869435595348477f; final float ti11 = 0.951056516295153531181938433292089030f; final float tr12 = -0.809016994374947340240566973079694435f; final float ti12 = 0.587785252292473248125759255344746634f; long i , ic ; float ci2 , di2 , ci4 , ci5 , di3 , di4 , di5 , ci3 , cr2 , cr3 , dr2 , dr3 , dr4 , dr5 , cr5 , cr4 , ti2 , ti3 , ti5 , ti4 , tr2 , tr3 , tr4 , tr5 , w1r , w1i , w2r , w2i , w3r , w3i , w4r , w4i ; long iw1 , iw2 , iw3 , iw4 ; iw1 = offset ; iw2 = iw1 + ido ; iw3 = iw2 + ido ; iw4 = iw3 + ido ; long idx0 = l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 5 * idx1 ; long idx3 = idx2 + ido ; long idx4 = idx3 + ido ; long idx5 = idx4 + ido ; long idx6 = idx5 + ido ; long idx7 = idx1 + idx0 ; long idx8 = idx7 + idx0 ; long idx9 = idx8 + idx0 ; long idx10 = idx9 + idx0 ; long idx11 = out_off + ido - 1; float i1r = in.getFloat (in_off + idx1 ); float i2r = in.getFloat (in_off + idx7 ); float i3r = in.getFloat (in_off + idx8 ); float i4r = in.getFloat (in_off + idx9 ); float i5r = in.getFloat (in_off + idx10 ); cr2 = i5r + i2r ; ci5 = i5r - i2r ; cr3 = i4r + i3r ; ci4 = i4r - i3r ; out.setFloat (out_off + idx2 , i1r + cr2 + cr3 ); out.setFloat (idx11 + idx3 , i1r + tr11 * cr2 + tr12 * cr3 ); out.setFloat (out_off + idx4 , ti11 * ci5 + ti12 * ci4 ); out.setFloat (idx11 + idx5 , i1r + tr12 * cr2 + tr11 * cr3 ); out.setFloat (out_off + idx6 , ti12 * ci5 - ti11 * ci4 ); } if (ido == 1) { return ; } for (long k = 0; k < l1 ; ++k ) { long idx1 = k * ido ; long idx2 = 5 * idx1 ; long idx3 = idx2 + ido ; long idx4 = idx3 + ido ; long idx5 = idx4 + ido ; long idx6 = idx5 + ido ; long idx7 = idx1 + idx0 ; long idx8 = idx7 + idx0 ; long idx9 = idx8 + idx0 ; long idx10 = idx9 + idx0 ; for (i = 2; i < ido ; i += 2) { long widx1 = i - 1 + iw1 ; long widx2 = i - 1 + iw2 ; long widx3 = i - 1 + iw3 ; long widx4 = i - 1 + iw4 ; w1r = wtable_rl.getFloat (widx1 - 1); w1i = wtable_rl.getFloat (widx1 ); w2r = wtable_rl.getFloat (widx2 - 1); w2i = wtable_rl.getFloat (widx2 ); w3r = wtable_rl.getFloat (widx3 - 1); w3i = wtable_rl.getFloat (widx3 ); w4r = wtable_rl.getFloat (widx4 - 1); w4i = wtable_rl.getFloat (widx4 ); ic = ido - i ; long idx15 = in_off + i ; long idx16 = out_off + i ; long idx17 = out_off + ic ; long iidx1 = idx15 + idx1 ; long iidx2 = idx15 + idx7 ; long iidx3 = idx15 + idx8 ; long iidx4 = idx15 + idx9 ; long iidx5 = idx15 + idx10 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); float i2i = in.getFloat (iidx2 - 1); float i2r = in.getFloat (iidx2 ); float i3i = in.getFloat (iidx3 - 1); float i3r = in.getFloat (iidx3 ); float i4i = in.getFloat (iidx4 - 1); float i4r = in.getFloat (iidx4 ); float i5i = in.getFloat (iidx5 - 1); float i5r = in.getFloat (iidx5 ); dr2 = w1r * i2i + w1i * i2r ; di2 = w1r * i2r - w1i * i2i ; dr3 = w2r * i3i + w2i * i3r ; di3 = w2r * i3r - w2i * i3i ; dr4 = w3r * i4i + w3i * i4r ; di4 = w3r * i4r - w3i * i4i ; dr5 = w4r * i5i + w4i * i5r ; di5 = w4r * i5r - w4i * i5i ; cr2 = dr2 + dr5 ; ci5 = dr5 - dr2 ; cr5 = di2 - di5 ; ci2 = di2 + di5 ; cr3 = dr3 + dr4 ; ci4 = dr4 - dr3 ; cr4 = di3 - di4 ; ci3 = di3 + di4 ; tr2 = i1i + tr11 * cr2 + tr12 * cr3 ; ti2 = i1r + tr11 * ci2 + tr12 * ci3 ; tr3 = i1i + tr12 * cr2 + tr11 * cr3 ; ti3 = i1r + tr12 * ci2 + tr11 * ci3 ; tr5 = ti11 * cr5 + ti12 * cr4 ; ti5 = ti11 * ci5 + ti12 * ci4 ; tr4 = ti12 * cr5 - ti11 * cr4 ; ti4 = ti12 * ci5 - ti11 * ci4 ; long oidx1 = idx16 + idx2 ; long oidx2 = idx17 + idx3 ; long oidx3 = idx16 + idx4 ; long oidx4 = idx17 + idx5 ; long oidx5 = idx16 + idx6 ; out.setFloat (oidx1 - 1, i1i + cr2 + cr3 ); out.setFloat (oidx1 , i1r + ci2 + ci3 ); out.setFloat (oidx3 - 1, tr2 + tr5 ); out.setFloat (oidx2 - 1, tr2 - tr5 ); out.setFloat (oidx3 , ti2 + ti5 ); out.setFloat (oidx2 , ti5 - ti2 ); out.setFloat (oidx5 - 1, tr3 + tr4 ); out.setFloat (oidx4 - 1, tr3 - tr4 ); out.setFloat (oidx5 , ti3 + ti4 ); out.setFloat (oidx4 , ti4 - ti3 ); } } } /*------------------------------------------------- radb5: Real FFT's backward processing of factor 5 -------------------------------------------------*/ void radb5 (final int ido , final int l1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { final float tr11 = 0.309016994374947451262869435595348477f; final float ti11 = 0.951056516295153531181938433292089030f; final float tr12 = -0.809016994374947340240566973079694435f; final float ti12 = 0.587785252292473248125759255344746634f; int i , ic ; float ci2 , ci3 , ci4 , ci5 , di3 , di4 , di5 , di2 , cr2 , cr3 , cr5 , cr4 , ti2 , ti3 , ti4 , ti5 , dr3 , dr4 , dr5 , dr2 , tr2 , tr3 , tr4 , tr5 , w1r , w1i , w2r , w2i , w3r , w3i , w4r , w4i ; int iw1 , iw2 , iw3 , iw4 ; iw1 = offset ; iw2 = iw1 + ido ; iw3 = iw2 + ido ; iw4 = iw3 + ido ; int idx0 = l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = 5 * idx1 ; int idx3 = idx2 + ido ; int idx4 = idx3 + ido ; int idx5 = idx4 + ido ; int idx6 = idx5 + ido ; int idx7 = idx1 + idx0 ; int idx8 = idx7 + idx0 ; int idx9 = idx8 + idx0 ; int idx10 = idx9 + idx0 ; int idx11 = in_off + ido - 1; float i1r = in [in_off + idx2 ]; ti5 = 2 * in [in_off + idx4 ]; ti4 = 2 * in [in_off + idx6 ]; tr2 = 2 * in [idx11 + idx3 ]; tr3 = 2 * in [idx11 + idx5 ]; cr2 = i1r + tr11 * tr2 + tr12 * tr3 ; cr3 = i1r + tr12 * tr2 + tr11 * tr3 ; ci5 = ti11 * ti5 + ti12 * ti4 ; ci4 = ti12 * ti5 - ti11 * ti4 ; out [out_off + idx1 ] = i1r + tr2 + tr3 ; out [out_off + idx7 ] = cr2 - ci5 ; out [out_off + idx8 ] = cr3 - ci4 ; out [out_off + idx9 ] = cr3 + ci4 ; out [out_off + idx10 ] = cr2 + ci5 ; } if (ido == 1) { return ; } for (int k = 0; k < l1 ; ++k ) { int idx1 = k * ido ; int idx2 = 5 * idx1 ; int idx3 = idx2 + ido ; int idx4 = idx3 + ido ; int idx5 = idx4 + ido ; int idx6 = idx5 + ido ; int idx7 = idx1 + idx0 ; int idx8 = idx7 + idx0 ; int idx9 = idx8 + idx0 ; int idx10 = idx9 + idx0 ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; int widx1 = i - 1 + iw1 ; int widx2 = i - 1 + iw2 ; int widx3 = i - 1 + iw3 ; int widx4 = i - 1 + iw4 ; w1r = wtable_r [widx1 - 1]; w1i = wtable_r [widx1 ]; w2r = wtable_r [widx2 - 1]; w2i = wtable_r [widx2 ]; w3r = wtable_r [widx3 - 1]; w3i = wtable_r [widx3 ]; w4r = wtable_r [widx4 - 1]; w4i = wtable_r [widx4 ]; int idx15 = in_off + i ; int idx16 = in_off + ic ; int idx17 = out_off + i ; int iidx1 = idx15 + idx2 ; int iidx2 = idx16 + idx3 ; int iidx3 = idx15 + idx4 ; int iidx4 = idx16 + idx5 ; int iidx5 = idx15 + idx6 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; float i2i = in [iidx2 - 1]; float i2r = in [iidx2 ]; float i3i = in [iidx3 - 1]; float i3r = in [iidx3 ]; float i4i = in [iidx4 - 1]; float i4r = in [iidx4 ]; float i5i = in [iidx5 - 1]; float i5r = in [iidx5 ]; ti5 = i3r + i2r ; ti2 = i3r - i2r ; ti4 = i5r + i4r ; ti3 = i5r - i4r ; tr5 = i3i - i2i ; tr2 = i3i + i2i ; tr4 = i5i - i4i ; tr3 = i5i + i4i ; cr2 = i1i + tr11 * tr2 + tr12 * tr3 ; ci2 = i1r + tr11 * ti2 + tr12 * ti3 ; cr3 = i1i + tr12 * tr2 + tr11 * tr3 ; ci3 = i1r + tr12 * ti2 + tr11 * ti3 ; cr5 = ti11 * tr5 + ti12 * tr4 ; ci5 = ti11 * ti5 + ti12 * ti4 ; cr4 = ti12 * tr5 - ti11 * tr4 ; ci4 = ti12 * ti5 - ti11 * ti4 ; dr3 = cr3 - ci4 ; dr4 = cr3 + ci4 ; di3 = ci3 + cr4 ; di4 = ci3 - cr4 ; dr5 = cr2 + ci5 ; dr2 = cr2 - ci5 ; di5 = ci2 - cr5 ; di2 = ci2 + cr5 ; int oidx1 = idx17 + idx1 ; int oidx2 = idx17 + idx7 ; int oidx3 = idx17 + idx8 ; int oidx4 = idx17 + idx9 ; int oidx5 = idx17 + idx10 ; out [oidx1 - 1] = i1i + tr2 + tr3 ; out [oidx1 ] = i1r + ti2 + ti3 ; out [oidx2 - 1] = w1r * dr2 - w1i * di2 ; out [oidx2 ] = w1r * di2 + w1i * dr2 ; out [oidx3 - 1] = w2r * dr3 - w2i * di3 ; out [oidx3 ] = w2r * di3 + w2i * dr3 ; out [oidx4 - 1] = w3r * dr4 - w3i * di4 ; out [oidx4 ] = w3r * di4 + w3i * dr4 ; out [oidx5 - 1] = w4r * dr5 - w4i * di5 ; out [oidx5 ] = w4r * di5 + w4i * dr5 ; } } } /*------------------------------------------------- radb5: Real FFT's backward processing of factor 5 -------------------------------------------------*/ void radb5 (final long ido , final long l1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { final float tr11 = 0.309016994374947451262869435595348477f; final float ti11 = 0.951056516295153531181938433292089030f; final float tr12 = -0.809016994374947340240566973079694435f; final float ti12 = 0.587785252292473248125759255344746634f; long i , ic ; float ci2 , ci3 , ci4 , ci5 , di3 , di4 , di5 , di2 , cr2 , cr3 , cr5 , cr4 , ti2 , ti3 , ti4 , ti5 , dr3 , dr4 , dr5 , dr2 , tr2 , tr3 , tr4 , tr5 , w1r , w1i , w2r , w2i , w3r , w3i , w4r , w4i ; long iw1 , iw2 , iw3 , iw4 ; iw1 = offset ; iw2 = iw1 + ido ; iw3 = iw2 + ido ; iw4 = iw3 + ido ; long idx0 = l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = 5 * idx1 ; long idx3 = idx2 + ido ; long idx4 = idx3 + ido ; long idx5 = idx4 + ido ; long idx6 = idx5 + ido ; long idx7 = idx1 + idx0 ; long idx8 = idx7 + idx0 ; long idx9 = idx8 + idx0 ; long idx10 = idx9 + idx0 ; long idx11 = in_off + ido - 1; float i1r = in.getFloat (in_off + idx2 ); ti5 = 2 * in.getFloat (in_off + idx4 ); ti4 = 2 * in.getFloat (in_off + idx6 ); tr2 = 2 * in.getFloat (idx11 + idx3 ); tr3 = 2 * in.getFloat (idx11 + idx5 ); cr2 = i1r + tr11 * tr2 + tr12 * tr3 ; cr3 = i1r + tr12 * tr2 + tr11 * tr3 ; ci5 = ti11 * ti5 + ti12 * ti4 ; ci4 = ti12 * ti5 - ti11 * ti4 ; out.setFloat (out_off + idx1 , i1r + tr2 + tr3 ); out.setFloat (out_off + idx7 , cr2 - ci5 ); out.setFloat (out_off + idx8 , cr3 - ci4 ); out.setFloat (out_off + idx9 , cr3 + ci4 ); out.setFloat (out_off + idx10 , cr2 + ci5 ); } if (ido == 1) { return ; } for (long k = 0; k < l1 ; ++k ) { long idx1 = k * ido ; long idx2 = 5 * idx1 ; long idx3 = idx2 + ido ; long idx4 = idx3 + ido ; long idx5 = idx4 + ido ; long idx6 = idx5 + ido ; long idx7 = idx1 + idx0 ; long idx8 = idx7 + idx0 ; long idx9 = idx8 + idx0 ; long idx10 = idx9 + idx0 ; for (i = 2; i < ido ; i += 2) { ic = ido - i ; long widx1 = i - 1 + iw1 ; long widx2 = i - 1 + iw2 ; long widx3 = i - 1 + iw3 ; long widx4 = i - 1 + iw4 ; w1r = wtable_rl.getFloat (widx1 - 1); w1i = wtable_rl.getFloat (widx1 ); w2r = wtable_rl.getFloat (widx2 - 1); w2i = wtable_rl.getFloat (widx2 ); w3r = wtable_rl.getFloat (widx3 - 1); w3i = wtable_rl.getFloat (widx3 ); w4r = wtable_rl.getFloat (widx4 - 1); w4i = wtable_rl.getFloat (widx4 ); long idx15 = in_off + i ; long idx16 = in_off + ic ; long idx17 = out_off + i ; long iidx1 = idx15 + idx2 ; long iidx2 = idx16 + idx3 ; long iidx3 = idx15 + idx4 ; long iidx4 = idx16 + idx5 ; long iidx5 = idx15 + idx6 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); float i2i = in.getFloat (iidx2 - 1); float i2r = in.getFloat (iidx2 ); float i3i = in.getFloat (iidx3 - 1); float i3r = in.getFloat (iidx3 ); float i4i = in.getFloat (iidx4 - 1); float i4r = in.getFloat (iidx4 ); float i5i = in.getFloat (iidx5 - 1); float i5r = in.getFloat (iidx5 ); ti5 = i3r + i2r ; ti2 = i3r - i2r ; ti4 = i5r + i4r ; ti3 = i5r - i4r ; tr5 = i3i - i2i ; tr2 = i3i + i2i ; tr4 = i5i - i4i ; tr3 = i5i + i4i ; cr2 = i1i + tr11 * tr2 + tr12 * tr3 ; ci2 = i1r + tr11 * ti2 + tr12 * ti3 ; cr3 = i1i + tr12 * tr2 + tr11 * tr3 ; ci3 = i1r + tr12 * ti2 + tr11 * ti3 ; cr5 = ti11 * tr5 + ti12 * tr4 ; ci5 = ti11 * ti5 + ti12 * ti4 ; cr4 = ti12 * tr5 - ti11 * tr4 ; ci4 = ti12 * ti5 - ti11 * ti4 ; dr3 = cr3 - ci4 ; dr4 = cr3 + ci4 ; di3 = ci3 + cr4 ; di4 = ci3 - cr4 ; dr5 = cr2 + ci5 ; dr2 = cr2 - ci5 ; di5 = ci2 - cr5 ; di2 = ci2 + cr5 ; long oidx1 = idx17 + idx1 ; long oidx2 = idx17 + idx7 ; long oidx3 = idx17 + idx8 ; long oidx4 = idx17 + idx9 ; long oidx5 = idx17 + idx10 ; out.setFloat (oidx1 - 1, i1i + tr2 + tr3 ); out.setFloat (oidx1 , i1r + ti2 + ti3 ); out.setFloat (oidx2 - 1, w1r * dr2 - w1i * di2 ); out.setFloat (oidx2 , w1r * di2 + w1i * dr2 ); out.setFloat (oidx3 - 1, w2r * dr3 - w2i * di3 ); out.setFloat (oidx3 , w2r * di3 + w2i * dr3 ); out.setFloat (oidx4 - 1, w3r * dr4 - w3i * di4 ); out.setFloat (oidx4 , w3r * di4 + w3i * dr4 ); out.setFloat (oidx5 - 1, w4r * dr5 - w4i * di5 ); out.setFloat (oidx5 , w4r * di5 + w4i * dr5 ); } } } /*--------------------------------------------------------- radfg: Real FFT's forward processing of general factor --------------------------------------------------------*/ void radfg (final int ido , final int ip , final int l1 , final int idl1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { int idij , ipph , j2 , ic , jc , lc , is , nbd ; float dc2 , ai1 , ai2 , ar1 , ar2 , ds2 , dcp , arg , dsp , ar1h , ar2h , w1r , w1i ; int iw1 = offset ; arg = TWO_PI / (float ) ip ; dcp = (float ) cos (arg ); dsp = (float ) sin (arg ); ipph = (ip + 1) / 2; nbd = (ido - 1) / 2; if (ido != 1) { for (int ik = 0; ik < idl1 ; ik ++) { out [out_off + ik ] = in [in_off + ik ]; } for (int j = 1; j < ip ; j ++) { int idx1 = j * l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx2 = k * ido + idx1 ; out [out_off + idx2 ] = in [in_off + idx2 ]; } } if (nbd <= l1 ) { is = -ido ; for (int j = 1; j < ip ; j ++) { is += ido ; idij = is - 1; int idx1 = j * l1 * ido ; for (int i = 2; i < ido ; i += 2) { idij += 2; int idx2 = idij + iw1 ; int idx4 = in_off + i ; int idx5 = out_off + i ; w1r = wtable_r [idx2 - 1]; w1i = wtable_r [idx2 ]; for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido + idx1 ; int oidx1 = idx5 + idx3 ; int iidx1 = idx4 + idx3 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; out [oidx1 - 1] = w1r * i1i + w1i * i1r ; out [oidx1 ] = w1r * i1r - w1i * i1i ; } } } } else { is = -ido ; for (int j = 1; j < ip ; j ++) { is += ido ; int idx1 = j * l1 * ido ; for (int k = 0; k < l1 ; k ++) { idij = is - 1; int idx3 = k * ido + idx1 ; for (int i = 2; i < ido ; i += 2) { idij += 2; int idx2 = idij + iw1 ; w1r = wtable_r [idx2 - 1]; w1i = wtable_r [idx2 ]; int oidx1 = out_off + i + idx3 ; int iidx1 = in_off + i + idx3 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; out [oidx1 - 1] = w1r * i1i + w1i * i1r ; out [oidx1 ] = w1r * i1r - w1i * i1i ; } } } } if (nbd >= l1 ) { for (int j = 1; j < ipph ; j ++) { jc = ip - j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido + idx1 ; int idx4 = k * ido + idx2 ; for (int i = 2; i < ido ; i += 2) { int idx5 = in_off + i ; int idx6 = out_off + i ; int iidx1 = idx5 + idx3 ; int iidx2 = idx5 + idx4 ; int oidx1 = idx6 + idx3 ; int oidx2 = idx6 + idx4 ; float o1i = out [oidx1 - 1]; float o1r = out [oidx1 ]; float o2i = out [oidx2 - 1]; float o2r = out [oidx2 ]; in [iidx1 - 1] = o1i + o2i ; in [iidx1 ] = o1r + o2r ; in [iidx2 - 1] = o1r - o2r ; in [iidx2 ] = o2i - o1i ; } } } } else { for (int j = 1; j < ipph ; j ++) { jc = ip - j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; for (int i = 2; i < ido ; i += 2) { int idx5 = in_off + i ; int idx6 = out_off + i ; for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido + idx1 ; int idx4 = k * ido + idx2 ; int iidx1 = idx5 + idx3 ; int iidx2 = idx5 + idx4 ; int oidx1 = idx6 + idx3 ; int oidx2 = idx6 + idx4 ; float o1i = out [oidx1 - 1]; float o1r = out [oidx1 ]; float o2i = out [oidx2 - 1]; float o2r = out [oidx2 ]; in [iidx1 - 1] = o1i + o2i ; in [iidx1 ] = o1r + o2r ; in [iidx2 - 1] = o1r - o2r ; in [iidx2 ] = o2i - o1i ; } } } } } else { System.arraycopy (out , out_off , in , in_off , idl1 ); } for (int j = 1; j < ipph ; j ++) { jc = ip - j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido + idx1 ; int idx4 = k * ido + idx2 ; int oidx1 = out_off + idx3 ; int oidx2 = out_off + idx4 ; float o1r = out [oidx1 ]; float o2r = out [oidx2 ]; in [in_off + idx3 ] = o1r + o2r ; in [in_off + idx4 ] = o2r - o1r ; } } ar1 = 1; ai1 = 0; int idx0 = (ip - 1) * idl1 ; for (int l = 1; l < ipph ; l ++) { lc = ip - l ; ar1h = dcp * ar1 - dsp * ai1 ; ai1 = dcp * ai1 + dsp * ar1 ; ar1 = ar1h ; int idx1 = l * idl1 ; int idx2 = lc * idl1 ; for (int ik = 0; ik < idl1 ; ik ++) { int idx3 = out_off + ik ; int idx4 = in_off + ik ; out [idx3 + idx1 ] = in [idx4 ] + ar1 * in [idx4 + idl1 ]; out [idx3 + idx2 ] = ai1 * in [idx4 + idx0 ]; } dc2 = ar1 ; ds2 = ai1 ; ar2 = ar1 ; ai2 = ai1 ; for (int j = 2; j < ipph ; j ++) { jc = ip - j ; ar2h = dc2 * ar2 - ds2 * ai2 ; ai2 = dc2 * ai2 + ds2 * ar2 ; ar2 = ar2h ; int idx3 = j * idl1 ; int idx4 = jc * idl1 ; for (int ik = 0; ik < idl1 ; ik ++) { int idx5 = out_off + ik ; int idx6 = in_off + ik ; out [idx5 + idx1 ] += ar2 * in [idx6 + idx3 ]; out [idx5 + idx2 ] += ai2 * in [idx6 + idx4 ]; } } } for (int j = 1; j < ipph ; j ++) { int idx1 = j * idl1 ; for (int ik = 0; ik < idl1 ; ik ++) { out [out_off + ik ] += in [in_off + ik + idx1 ]; } } if (ido >= l1 ) { for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = idx1 * ip ; for (int i = 0; i < ido ; i ++) { in [in_off + i + idx2 ] = out [out_off + i + idx1 ]; } } } else { for (int i = 0; i < ido ; i ++) { for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; in [in_off + i + idx1 * ip ] = out [out_off + i + idx1 ]; } } } int idx01 = ip * ido ; for (int j = 1; j < ipph ; j ++) { jc = ip - j ; j2 = 2 * j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; int idx3 = j2 * ido ; for (int k = 0; k < l1 ; k ++) { int idx4 = k * ido ; int idx5 = idx4 + idx1 ; int idx6 = idx4 + idx2 ; int idx7 = k * idx01 ; in [in_off + ido - 1 + idx3 - ido + idx7 ] = out [out_off + idx5 ]; in [in_off + idx3 + idx7 ] = out [out_off + idx6 ]; } } if (ido == 1) { return ; } if (nbd >= l1 ) { for (int j = 1; j < ipph ; j ++) { jc = ip - j ; j2 = 2 * j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; int idx3 = j2 * ido ; for (int k = 0; k < l1 ; k ++) { int idx4 = k * idx01 ; int idx5 = k * ido ; for (int i = 2; i < ido ; i += 2) { ic = ido - i ; int idx6 = in_off + i ; int idx7 = in_off + ic ; int idx8 = out_off + i ; int iidx1 = idx6 + idx3 + idx4 ; int iidx2 = idx7 + idx3 - ido + idx4 ; int oidx1 = idx8 + idx5 + idx1 ; int oidx2 = idx8 + idx5 + idx2 ; float o1i = out [oidx1 - 1]; float o1r = out [oidx1 ]; float o2i = out [oidx2 - 1]; float o2r = out [oidx2 ]; in [iidx1 - 1] = o1i + o2i ; in [iidx2 - 1] = o1i - o2i ; in [iidx1 ] = o1r + o2r ; in [iidx2 ] = o2r - o1r ; } } } } else { for (int j = 1; j < ipph ; j ++) { jc = ip - j ; j2 = 2 * j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; int idx3 = j2 * ido ; for (int i = 2; i < ido ; i += 2) { ic = ido - i ; int idx6 = in_off + i ; int idx7 = in_off + ic ; int idx8 = out_off + i ; for (int k = 0; k < l1 ; k ++) { int idx4 = k * idx01 ; int idx5 = k * ido ; int iidx1 = idx6 + idx3 + idx4 ; int iidx2 = idx7 + idx3 - ido + idx4 ; int oidx1 = idx8 + idx5 + idx1 ; int oidx2 = idx8 + idx5 + idx2 ; float o1i = out [oidx1 - 1]; float o1r = out [oidx1 ]; float o2i = out [oidx2 - 1]; float o2r = out [oidx2 ]; in [iidx1 - 1] = o1i + o2i ; in [iidx2 - 1] = o1i - o2i ; in [iidx1 ] = o1r + o2r ; in [iidx2 ] = o2r - o1r ; } } } } } /*--------------------------------------------------------- radfg: Real FFT's forward processing of general factor --------------------------------------------------------*/ void radfg (final long ido , final long ip , final long l1 , final long idl1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { long idij , ipph , j2 , ic , jc , lc , is , nbd ; float dc2 , ai1 , ai2 , ar1 , ar2 , ds2 , dcp , arg , dsp , ar1h , ar2h , w1r , w1i ; long iw1 = offset ; arg = TWO_PI / (float ) ip ; dcp = (float ) cos (arg ); dsp = (float ) sin (arg ); ipph = (ip + 1) / 2; nbd = (ido - 1) / 2; if (ido != 1) { for (long ik = 0; ik < idl1 ; ik ++) { out.setFloat (out_off + ik , in.getFloat (in_off + ik )); } for (long j = 1; j < ip ; j ++) { long idx1 = j * l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx2 = k * ido + idx1 ; out.setFloat (out_off + idx2 , in.getFloat (in_off + idx2 )); } } if (nbd <= l1 ) { is = -ido ; for (long j = 1; j < ip ; j ++) { is += ido ; idij = is - 1; long idx1 = j * l1 * ido ; for (long i = 2; i < ido ; i += 2) { idij += 2; long idx2 = idij + iw1 ; long idx4 = in_off + i ; long idx5 = out_off + i ; w1r = wtable_rl.getFloat (idx2 - 1); w1i = wtable_rl.getFloat (idx2 ); for (long k = 0; k < l1 ; k ++) { long idx3 = k * ido + idx1 ; long oidx1 = idx5 + idx3 ; long iidx1 = idx4 + idx3 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); out.setFloat (oidx1 - 1, w1r * i1i + w1i * i1r ); out.setFloat (oidx1 , w1r * i1r - w1i * i1i ); } } } } else { is = -ido ; for (long j = 1; j < ip ; j ++) { is += ido ; long idx1 = j * l1 * ido ; for (long k = 0; k < l1 ; k ++) { idij = is - 1; long idx3 = k * ido + idx1 ; for (long i = 2; i < ido ; i += 2) { idij += 2; long idx2 = idij + iw1 ; w1r = wtable_rl.getFloat (idx2 - 1); w1i = wtable_rl.getFloat (idx2 ); long oidx1 = out_off + i + idx3 ; long iidx1 = in_off + i + idx3 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); out.setFloat (oidx1 - 1, w1r * i1i + w1i * i1r ); out.setFloat (oidx1 , w1r * i1r - w1i * i1i ); } } } } if (nbd >= l1 ) { for (long j = 1; j < ipph ; j ++) { jc = ip - j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx3 = k * ido + idx1 ; long idx4 = k * ido + idx2 ; for (long i = 2; i < ido ; i += 2) { long idx5 = in_off + i ; long idx6 = out_off + i ; long iidx1 = idx5 + idx3 ; long iidx2 = idx5 + idx4 ; long oidx1 = idx6 + idx3 ; long oidx2 = idx6 + idx4 ; float o1i = out.getFloat (oidx1 - 1); float o1r = out.getFloat (oidx1 ); float o2i = out.getFloat (oidx2 - 1); float o2r = out.getFloat (oidx2 ); in.setFloat (iidx1 - 1, o1i + o2i ); in.setFloat (iidx1 , o1r + o2r ); in.setFloat (iidx2 - 1, o1r - o2r ); in.setFloat (iidx2 , o2i - o1i ); } } } } else { for (long j = 1; j < ipph ; j ++) { jc = ip - j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; for (long i = 2; i < ido ; i += 2) { long idx5 = in_off + i ; long idx6 = out_off + i ; for (long k = 0; k < l1 ; k ++) { long idx3 = k * ido + idx1 ; long idx4 = k * ido + idx2 ; long iidx1 = idx5 + idx3 ; long iidx2 = idx5 + idx4 ; long oidx1 = idx6 + idx3 ; long oidx2 = idx6 + idx4 ; float o1i = out.getFloat (oidx1 - 1); float o1r = out.getFloat (oidx1 ); float o2i = out.getFloat (oidx2 - 1); float o2r = out.getFloat (oidx2 ); in.setFloat (iidx1 - 1, o1i + o2i ); in.setFloat (iidx1 , o1r + o2r ); in.setFloat (iidx2 - 1, o1r - o2r ); in.setFloat (iidx2 , o2i - o1i ); } } } } } else { LargeArrayUtils.arraycopy (out , out_off , in , in_off , idl1 ); } for (long j = 1; j < ipph ; j ++) { jc = ip - j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx3 = k * ido + idx1 ; long idx4 = k * ido + idx2 ; long oidx1 = out_off + idx3 ; long oidx2 = out_off + idx4 ; float o1r = out.getFloat (oidx1 ); float o2r = out.getFloat (oidx2 ); in.setFloat (in_off + idx3 , o1r + o2r ); in.setFloat (in_off + idx4 , o2r - o1r ); } } ar1 = 1; ai1 = 0; long idx0 = (ip - 1) * idl1 ; for (long l = 1; l < ipph ; l ++) { lc = ip - l ; ar1h = dcp * ar1 - dsp * ai1 ; ai1 = dcp * ai1 + dsp * ar1 ; ar1 = ar1h ; long idx1 = l * idl1 ; long idx2 = lc * idl1 ; for (long ik = 0; ik < idl1 ; ik ++) { long idx3 = out_off + ik ; long idx4 = in_off + ik ; out.setFloat (idx3 + idx1 , in.getFloat (idx4 ) + ar1 * in.getFloat (idx4 + idl1 )); out.setFloat (idx3 + idx2 , ai1 * in.getFloat (idx4 + idx0 )); } dc2 = ar1 ; ds2 = ai1 ; ar2 = ar1 ; ai2 = ai1 ; for (long j = 2; j < ipph ; j ++) { jc = ip - j ; ar2h = dc2 * ar2 - ds2 * ai2 ; ai2 = dc2 * ai2 + ds2 * ar2 ; ar2 = ar2h ; long idx3 = j * idl1 ; long idx4 = jc * idl1 ; for (long ik = 0; ik < idl1 ; ik ++) { long idx5 = out_off + ik ; long idx6 = in_off + ik ; out.setFloat (idx5 + idx1 , out.getFloat (idx5 + idx1 ) + ar2 * in.getFloat (idx6 + idx3 )); out.setFloat (idx5 + idx2 , out.getFloat (idx5 + idx2 ) + ai2 * in.getFloat (idx6 + idx4 )); } } } for (long j = 1; j < ipph ; j ++) { long idx1 = j * idl1 ; for (long ik = 0; ik < idl1 ; ik ++) { out.setFloat (out_off + ik , out.getFloat (out_off + ik ) + in.getFloat (in_off + ik + idx1 )); } } if (ido >= l1 ) { for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = idx1 * ip ; for (long i = 0; i < ido ; i ++) { in.setFloat (in_off + i + idx2 , out.getFloat (out_off + i + idx1 )); } } } else { for (long i = 0; i < ido ; i ++) { for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; in.setFloat (in_off + i + idx1 * ip , out.getFloat (out_off + i + idx1 )); } } } long idx01 = ip * ido ; for (long j = 1; j < ipph ; j ++) { jc = ip - j ; j2 = 2 * j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; long idx3 = j2 * ido ; for (long k = 0; k < l1 ; k ++) { long idx4 = k * ido ; long idx5 = idx4 + idx1 ; long idx6 = idx4 + idx2 ; long idx7 = k * idx01 ; in.setFloat (in_off + ido - 1 + idx3 - ido + idx7 , out.getFloat (out_off + idx5 )); in.setFloat (in_off + idx3 + idx7 , out.getFloat (out_off + idx6 )); } } if (ido == 1) { return ; } if (nbd >= l1 ) { for (long j = 1; j < ipph ; j ++) { jc = ip - j ; j2 = 2 * j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; long idx3 = j2 * ido ; for (long k = 0; k < l1 ; k ++) { long idx4 = k * idx01 ; long idx5 = k * ido ; for (long i = 2; i < ido ; i += 2) { ic = ido - i ; long idx6 = in_off + i ; long idx7 = in_off + ic ; long idx8 = out_off + i ; long iidx1 = idx6 + idx3 + idx4 ; long iidx2 = idx7 + idx3 - ido + idx4 ; long oidx1 = idx8 + idx5 + idx1 ; long oidx2 = idx8 + idx5 + idx2 ; float o1i = out.getFloat (oidx1 - 1); float o1r = out.getFloat (oidx1 ); float o2i = out.getFloat (oidx2 - 1); float o2r = out.getFloat (oidx2 ); in.setFloat (iidx1 - 1, o1i + o2i ); in.setFloat (iidx2 - 1, o1i - o2i ); in.setFloat (iidx1 , o1r + o2r ); in.setFloat (iidx2 , o2r - o1r ); } } } } else { for (long j = 1; j < ipph ; j ++) { jc = ip - j ; j2 = 2 * j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; long idx3 = j2 * ido ; for (long i = 2; i < ido ; i += 2) { ic = ido - i ; long idx6 = in_off + i ; long idx7 = in_off + ic ; long idx8 = out_off + i ; for (long k = 0; k < l1 ; k ++) { long idx4 = k * idx01 ; long idx5 = k * ido ; long iidx1 = idx6 + idx3 + idx4 ; long iidx2 = idx7 + idx3 - ido + idx4 ; long oidx1 = idx8 + idx5 + idx1 ; long oidx2 = idx8 + idx5 + idx2 ; float o1i = out.getFloat (oidx1 - 1); float o1r = out.getFloat (oidx1 ); float o2i = out.getFloat (oidx2 - 1); float o2r = out.getFloat (oidx2 ); in.setFloat (iidx1 - 1, o1i + o2i ); in.setFloat (iidx2 - 1, o1i - o2i ); in.setFloat (iidx1 , o1r + o2r ); in.setFloat (iidx2 , o2r - o1r ); } } } } } /*--------------------------------------------------------- radbg: Real FFT's backward processing of general factor --------------------------------------------------------*/ void radbg (final int ido , final int ip , final int l1 , final int idl1 , final float in [], final int in_off , final float out [], final int out_off , final int offset ) { int idij , ipph , j2 , ic , jc , lc , is ; float dc2 , ai1 , ai2 , ar1 , ar2 , ds2 , w1r , w1i ; int nbd ; float dcp , arg , dsp , ar1h , ar2h ; int iw1 = offset ; arg = TWO_PI / (float ) ip ; dcp = (float ) cos (arg ); dsp = (float ) sin (arg ); nbd = (ido - 1) / 2; ipph = (ip + 1) / 2; int idx0 = ip * ido ; if (ido >= l1 ) { for (int k = 0; k < l1 ; k ++) { int idx1 = k * ido ; int idx2 = k * idx0 ; for (int i = 0; i < ido ; i ++) { out [out_off + i + idx1 ] = in [in_off + i + idx2 ]; } } } else { for (int i = 0; i < ido ; i ++) { int idx1 = out_off + i ; int idx2 = in_off + i ; for (int k = 0; k < l1 ; k ++) { out [idx1 + k * ido ] = in [idx2 + k * idx0 ]; } } } int iidx0 = in_off + ido - 1; for (int j = 1; j < ipph ; j ++) { jc = ip - j ; j2 = 2 * j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; int idx3 = j2 * ido ; for (int k = 0; k < l1 ; k ++) { int idx4 = k * ido ; int idx5 = idx4 * ip ; int iidx1 = iidx0 + idx3 + idx5 - ido ; int iidx2 = in_off + idx3 + idx5 ; float i1r = in [iidx1 ]; float i2r = in [iidx2 ]; out [out_off + idx4 + idx1 ] = i1r + i1r ; out [out_off + idx4 + idx2 ] = i2r + i2r ; } } if (ido != 1) { if (nbd >= l1 ) { for (int j = 1; j < ipph ; j ++) { jc = ip - j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; int idx3 = 2 * j * ido ; for (int k = 0; k < l1 ; k ++) { int idx4 = k * ido + idx1 ; int idx5 = k * ido + idx2 ; int idx6 = k * ip * ido + idx3 ; for (int i = 2; i < ido ; i += 2) { ic = ido - i ; int idx7 = out_off + i ; int idx8 = in_off + ic ; int idx9 = in_off + i ; int oidx1 = idx7 + idx4 ; int oidx2 = idx7 + idx5 ; int iidx1 = idx9 + idx6 ; int iidx2 = idx8 + idx6 - ido ; float a1i = in [iidx1 - 1]; float a1r = in [iidx1 ]; float a2i = in [iidx2 - 1]; float a2r = in [iidx2 ]; out [oidx1 - 1] = a1i + a2i ; out [oidx2 - 1] = a1i - a2i ; out [oidx1 ] = a1r - a2r ; out [oidx2 ] = a1r + a2r ; } } } } else { for (int j = 1; j < ipph ; j ++) { jc = ip - j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; int idx3 = 2 * j * ido ; for (int i = 2; i < ido ; i += 2) { ic = ido - i ; int idx7 = out_off + i ; int idx8 = in_off + ic ; int idx9 = in_off + i ; for (int k = 0; k < l1 ; k ++) { int idx4 = k * ido + idx1 ; int idx5 = k * ido + idx2 ; int idx6 = k * ip * ido + idx3 ; int oidx1 = idx7 + idx4 ; int oidx2 = idx7 + idx5 ; int iidx1 = idx9 + idx6 ; int iidx2 = idx8 + idx6 - ido ; float a1i = in [iidx1 - 1]; float a1r = in [iidx1 ]; float a2i = in [iidx2 - 1]; float a2r = in [iidx2 ]; out [oidx1 - 1] = a1i + a2i ; out [oidx2 - 1] = a1i - a2i ; out [oidx1 ] = a1r - a2r ; out [oidx2 ] = a1r + a2r ; } } } } } ar1 = 1; ai1 = 0; int idx01 = (ip - 1) * idl1 ; for (int l = 1; l < ipph ; l ++) { lc = ip - l ; ar1h = dcp * ar1 - dsp * ai1 ; ai1 = dcp * ai1 + dsp * ar1 ; ar1 = ar1h ; int idx1 = l * idl1 ; int idx2 = lc * idl1 ; for (int ik = 0; ik < idl1 ; ik ++) { int idx3 = in_off + ik ; int idx4 = out_off + ik ; in [idx3 + idx1 ] = out [idx4 ] + ar1 * out [idx4 + idl1 ]; in [idx3 + idx2 ] = ai1 * out [idx4 + idx01 ]; } dc2 = ar1 ; ds2 = ai1 ; ar2 = ar1 ; ai2 = ai1 ; for (int j = 2; j < ipph ; j ++) { jc = ip - j ; ar2h = dc2 * ar2 - ds2 * ai2 ; ai2 = dc2 * ai2 + ds2 * ar2 ; ar2 = ar2h ; int idx5 = j * idl1 ; int idx6 = jc * idl1 ; for (int ik = 0; ik < idl1 ; ik ++) { int idx7 = in_off + ik ; int idx8 = out_off + ik ; in [idx7 + idx1 ] += ar2 * out [idx8 + idx5 ]; in [idx7 + idx2 ] += ai2 * out [idx8 + idx6 ]; } } } for (int j = 1; j < ipph ; j ++) { int idx1 = j * idl1 ; for (int ik = 0; ik < idl1 ; ik ++) { int idx2 = out_off + ik ; out [idx2 ] += out [idx2 + idx1 ]; } } for (int j = 1; j < ipph ; j ++) { jc = ip - j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido ; int oidx1 = out_off + idx3 ; int iidx1 = in_off + idx3 + idx1 ; int iidx2 = in_off + idx3 + idx2 ; float i1r = in [iidx1 ]; float i2r = in [iidx2 ]; out [oidx1 + idx1 ] = i1r - i2r ; out [oidx1 + idx2 ] = i1r + i2r ; } } if (ido == 1) { return ; } if (nbd >= l1 ) { for (int j = 1; j < ipph ; j ++) { jc = ip - j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido ; for (int i = 2; i < ido ; i += 2) { int idx4 = out_off + i ; int idx5 = in_off + i ; int oidx1 = idx4 + idx3 + idx1 ; int oidx2 = idx4 + idx3 + idx2 ; int iidx1 = idx5 + idx3 + idx1 ; int iidx2 = idx5 + idx3 + idx2 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; float i2i = in [iidx2 - 1]; float i2r = in [iidx2 ]; out [oidx1 - 1] = i1i - i2r ; out [oidx2 - 1] = i1i + i2r ; out [oidx1 ] = i1r + i2i ; out [oidx2 ] = i1r - i2i ; } } } } else { for (int j = 1; j < ipph ; j ++) { jc = ip - j ; int idx1 = j * l1 * ido ; int idx2 = jc * l1 * ido ; for (int i = 2; i < ido ; i += 2) { int idx4 = out_off + i ; int idx5 = in_off + i ; for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido ; int oidx1 = idx4 + idx3 + idx1 ; int oidx2 = idx4 + idx3 + idx2 ; int iidx1 = idx5 + idx3 + idx1 ; int iidx2 = idx5 + idx3 + idx2 ; float i1i = in [iidx1 - 1]; float i1r = in [iidx1 ]; float i2i = in [iidx2 - 1]; float i2r = in [iidx2 ]; out [oidx1 - 1] = i1i - i2r ; out [oidx2 - 1] = i1i + i2r ; out [oidx1 ] = i1r + i2i ; out [oidx2 ] = i1r - i2i ; } } } } System.arraycopy (out , out_off , in , in_off , idl1 ); for (int j = 1; j < ip ; j ++) { int idx1 = j * l1 * ido ; for (int k = 0; k < l1 ; k ++) { int idx2 = k * ido + idx1 ; in [in_off + idx2 ] = out [out_off + idx2 ]; } } if (nbd <= l1 ) { is = -ido ; for (int j = 1; j < ip ; j ++) { is += ido ; idij = is - 1; int idx1 = j * l1 * ido ; for (int i = 2; i < ido ; i += 2) { idij += 2; int idx2 = idij + iw1 ; w1r = wtable_r [idx2 - 1]; w1i = wtable_r [idx2 ]; int idx4 = in_off + i ; int idx5 = out_off + i ; for (int k = 0; k < l1 ; k ++) { int idx3 = k * ido + idx1 ; int iidx1 = idx4 + idx3 ; int oidx1 = idx5 + idx3 ; float o1i = out [oidx1 - 1]; float o1r = out [oidx1 ]; in [iidx1 - 1] = w1r * o1i - w1i * o1r ; in [iidx1 ] = w1r * o1r + w1i * o1i ; } } } } else { is = -ido ; for (int j = 1; j < ip ; j ++) { is += ido ; int idx1 = j * l1 * ido ; for (int k = 0; k < l1 ; k ++) { idij = is - 1; int idx3 = k * ido + idx1 ; for (int i = 2; i < ido ; i += 2) { idij += 2; int idx2 = idij + iw1 ; w1r = wtable_r [idx2 - 1]; w1i = wtable_r [idx2 ]; int idx4 = in_off + i ; int idx5 = out_off + i ; int iidx1 = idx4 + idx3 ; int oidx1 = idx5 + idx3 ; float o1i = out [oidx1 - 1]; float o1r = out [oidx1 ]; in [iidx1 - 1] = w1r * o1i - w1i * o1r ; in [iidx1 ] = w1r * o1r + w1i * o1i ; } } } } } /*--------------------------------------------------------- radbg: Real FFT's backward processing of general factor --------------------------------------------------------*/ void radbg (final long ido , final long ip , final long l1 , final long idl1 , final FloatLargeArray in , final long in_off , final FloatLargeArray out , final long out_off , final long offset ) { long idij , ipph , j2 , ic , jc , lc , is ; float dc2 , ai1 , ai2 , ar1 , ar2 , ds2 , w1r , w1i ; long nbd ; float dcp , arg , dsp , ar1h , ar2h ; long iw1 = offset ; arg = TWO_PI / (float ) ip ; dcp = (float ) cos (arg ); dsp = (float ) sin (arg ); nbd = (ido - 1) / 2; ipph = (ip + 1) / 2; long idx0 = ip * ido ; if (ido >= l1 ) { for (long k = 0; k < l1 ; k ++) { long idx1 = k * ido ; long idx2 = k * idx0 ; for (long i = 0; i < ido ; i ++) { out.setFloat (out_off + i + idx1 , in.getFloat (in_off + i + idx2 )); } } } else { for (long i = 0; i < ido ; i ++) { long idx1 = out_off + i ; long idx2 = in_off + i ; for (long k = 0; k < l1 ; k ++) { out.setFloat (idx1 + k * ido , in.getFloat (idx2 + k * idx0 )); } } } long iidx0 = in_off + ido - 1; for (long j = 1; j < ipph ; j ++) { jc = ip - j ; j2 = 2 * j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; long idx3 = j2 * ido ; for (long k = 0; k < l1 ; k ++) { long idx4 = k * ido ; long idx5 = idx4 * ip ; long iidx1 = iidx0 + idx3 + idx5 - ido ; long iidx2 = in_off + idx3 + idx5 ; float i1r = in.getFloat (iidx1 ); float i2r = in.getFloat (iidx2 ); out.setFloat (out_off + idx4 + idx1 , i1r + i1r ); out.setFloat (out_off + idx4 + idx2 , i2r + i2r ); } } if (ido != 1) { if (nbd >= l1 ) { for (long j = 1; j < ipph ; j ++) { jc = ip - j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; long idx3 = 2 * j * ido ; for (long k = 0; k < l1 ; k ++) { long idx4 = k * ido + idx1 ; long idx5 = k * ido + idx2 ; long idx6 = k * ip * ido + idx3 ; for (long i = 2; i < ido ; i += 2) { ic = ido - i ; long idx7 = out_off + i ; long idx8 = in_off + ic ; long idx9 = in_off + i ; long oidx1 = idx7 + idx4 ; long oidx2 = idx7 + idx5 ; long iidx1 = idx9 + idx6 ; long iidx2 = idx8 + idx6 - ido ; float a1i = in.getFloat (iidx1 - 1); float a1r = in.getFloat (iidx1 ); float a2i = in.getFloat (iidx2 - 1); float a2r = in.getFloat (iidx2 ); out.setFloat (oidx1 - 1, a1i + a2i ); out.setFloat (oidx2 - 1, a1i - a2i ); out.setFloat (oidx1 , a1r - a2r ); out.setFloat (oidx2 , a1r + a2r ); } } } } else { for (long j = 1; j < ipph ; j ++) { jc = ip - j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; long idx3 = 2 * j * ido ; for (long i = 2; i < ido ; i += 2) { ic = ido - i ; long idx7 = out_off + i ; long idx8 = in_off + ic ; long idx9 = in_off + i ; for (long k = 0; k < l1 ; k ++) { long idx4 = k * ido + idx1 ; long idx5 = k * ido + idx2 ; long idx6 = k * ip * ido + idx3 ; long oidx1 = idx7 + idx4 ; long oidx2 = idx7 + idx5 ; long iidx1 = idx9 + idx6 ; long iidx2 = idx8 + idx6 - ido ; float a1i = in.getFloat (iidx1 - 1); float a1r = in.getFloat (iidx1 ); float a2i = in.getFloat (iidx2 - 1); float a2r = in.getFloat (iidx2 ); out.setFloat (oidx1 - 1, a1i + a2i ); out.setFloat (oidx2 - 1, a1i - a2i ); out.setFloat (oidx1 , a1r - a2r ); out.setFloat (oidx2 , a1r + a2r ); } } } } } ar1 = 1; ai1 = 0; long idx01 = (ip - 1) * idl1 ; for (long l = 1; l < ipph ; l ++) { lc = ip - l ; ar1h = dcp * ar1 - dsp * ai1 ; ai1 = dcp * ai1 + dsp * ar1 ; ar1 = ar1h ; long idx1 = l * idl1 ; long idx2 = lc * idl1 ; for (long ik = 0; ik < idl1 ; ik ++) { long idx3 = in_off + ik ; long idx4 = out_off + ik ; in.setFloat (idx3 + idx1 , out.getFloat (idx4 ) + ar1 * out.getFloat (idx4 + idl1 )); in.setFloat (idx3 + idx2 , ai1 * out.getFloat (idx4 + idx01 )); } dc2 = ar1 ; ds2 = ai1 ; ar2 = ar1 ; ai2 = ai1 ; for (long j = 2; j < ipph ; j ++) { jc = ip - j ; ar2h = dc2 * ar2 - ds2 * ai2 ; ai2 = dc2 * ai2 + ds2 * ar2 ; ar2 = ar2h ; long idx5 = j * idl1 ; long idx6 = jc * idl1 ; for (long ik = 0; ik < idl1 ; ik ++) { long idx7 = in_off + ik ; long idx8 = out_off + ik ; in.setFloat (idx7 + idx1 , in.getFloat (idx7 + idx1 ) + ar2 * out.getFloat (idx8 + idx5 )); in.setFloat (idx7 + idx2 , in.getFloat (idx7 + idx2 ) + ai2 * out.getFloat (idx8 + idx6 )); } } } for (long j = 1; j < ipph ; j ++) { long idx1 = j * idl1 ; for (long ik = 0; ik < idl1 ; ik ++) { long idx2 = out_off + ik ; out.setFloat (idx2 , out.getFloat (idx2 ) + out.getFloat (idx2 + idx1 )); } } for (long j = 1; j < ipph ; j ++) { jc = ip - j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx3 = k * ido ; long oidx1 = out_off + idx3 ; long iidx1 = in_off + idx3 + idx1 ; long iidx2 = in_off + idx3 + idx2 ; float i1r = in.getFloat (iidx1 ); float i2r = in.getFloat (iidx2 ); out.setFloat (oidx1 + idx1 , i1r - i2r ); out.setFloat (oidx1 + idx2 , i1r + i2r ); } } if (ido == 1) { return ; } if (nbd >= l1 ) { for (long j = 1; j < ipph ; j ++) { jc = ip - j ; long idx1 = j * l1 * ido ; long idx2 = jc * l1 * ido ; for (long k = 0; k < l1 ; k ++) { long idx3 = k * ido ; for (long i = 2; i < ido ; i += 2) { long idx4 = out_off + i ; long idx5 = in_off + i ; long oidx1 = idx4 + idx3 + idx1 ; long oidx2 = idx4 + idx3 + idx2 ; long iidx1 = idx5 + idx3 + idx1 ; long iidx2 = idx5 + idx3 + idx2 ; float i1i = in.getFloat (iidx1 - 1); float i1r = in.getFloat (iidx1 ); float i2i = in.getFloat (iidx2 - 1); float i2r = in.getFloat (iidx2 ); out.setFloat (oidx1 - 1, i1i - i2r ); out.setFloat (