Java Code Examples for org.apache.commons.math.complex.Complex#add()

The following examples show how to use org.apache.commons.math.complex.Complex#add() . You can vote up the ones you like or vote down the ones you don't like, and go to the original project or source file by following the links above each example. You may check out the related API usage on the sidebar.
Example 1
Source File: LaguerreSolver.java    From cacheonix-core with GNU Lesser General Public License v2.1 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 * 
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws MaxIterationsExceededException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function 
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    MaxIterationsExceededException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw new IllegalArgumentException
            ("Polynomial degree must be positive: degree=" + n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 2
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 *
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws MaxIterationsExceededException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    MaxIterationsExceededException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw MathRuntimeException.createIllegalArgumentException(
              NON_POSITIVE_DEGREE_MESSAGE, n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 3
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 * 
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws ConvergenceException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function 
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    ConvergenceException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw new IllegalArgumentException
            ("Polynomial degree must be positive: degree=" + n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 4
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 * 
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws MaxIterationsExceededException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function 
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    MaxIterationsExceededException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw MathRuntimeException.createIllegalArgumentException(
              "polynomial degree must be positive: degree={0}", n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 5
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 *
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws MaxIterationsExceededException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    MaxIterationsExceededException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw MathRuntimeException.createIllegalArgumentException(
              NON_POSITIVE_DEGREE_MESSAGE, n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 6
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 * 
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws MaxIterationsExceededException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function 
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    MaxIterationsExceededException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw MathRuntimeException.createIllegalArgumentException(
              "polynomial degree must be positive: degree={0}", n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 7
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given
 * coefficients, starting from the given initial value.
 *
 * @param coefficients Polynomial coefficients.
 * @param initial Start value.
 * @return the point at which the function value is zero.
 * @throws org.apache.commons.math.exception.TooManyEvaluationsException
 * if the maximum number of evaluations is exceeded.
 * @throws NullArgumentException if the {@code coefficients} is
 * {@code null}.
 * @throws NoDataException if the {@code coefficients} array is empty.
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) {
    if (coefficients == null) {
        throw new NullArgumentException();
    }
    int n = coefficients.length - 1;
    if (n == 0) {
        throw new NoDataException(LocalizedFormats.POLYNOMIAL);
    }
    // Coefficients for deflated polynomial.
    Complex c[] = new Complex[n + 1];
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // Solve individual roots successively.
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n - i + 1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // Polynomial deflation using synthetic division.
        Complex newc = c[n - i];
        Complex oldc = null;
        for (int j = n - i - 1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
    }

    return root;
}
Example 8
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 *
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws MaxIterationsExceededException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    MaxIterationsExceededException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw MathRuntimeException.createIllegalArgumentException(
              NON_POSITIVE_DEGREE_MESSAGE, n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 9
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 *
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws MaxIterationsExceededException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    MaxIterationsExceededException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw MathRuntimeException.createIllegalArgumentException(
              "polynomial degree must be positive: degree={0}", n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 10
Source File: LaguerreSolver.java    From astor with GNU General Public License v2.0 5 votes vote down vote up 
/**
 * Find all complex roots for the polynomial with the given coefficients,
 * starting from the given initial value.
 *
 * @param coefficients the polynomial coefficients array
 * @param initial the start value to use
 * @return the point at which the function value is zero
 * @throws MaxIterationsExceededException if the maximum iteration count is exceeded
 * or the solver detects convergence problems otherwise
 * @throws FunctionEvaluationException if an error occurs evaluating the
 * function
 * @throws IllegalArgumentException if any parameters are invalid
 */
public Complex[] solveAll(Complex coefficients[], Complex initial) throws
    MaxIterationsExceededException, FunctionEvaluationException {

    int n = coefficients.length - 1;
    int iterationCount = 0;
    if (n < 1) {
        throw MathRuntimeException.createIllegalArgumentException(
              NON_POSITIVE_DEGREE_MESSAGE, n);
    }
    Complex c[] = new Complex[n+1];    // coefficients for deflated polynomial
    for (int i = 0; i <= n; i++) {
        c[i] = coefficients[i];
    }

    // solve individual root successively
    Complex root[] = new Complex[n];
    for (int i = 0; i < n; i++) {
        Complex subarray[] = new Complex[n-i+1];
        System.arraycopy(c, 0, subarray, 0, subarray.length);
        root[i] = solve(subarray, initial);
        // polynomial deflation using synthetic division
        Complex newc = c[n-i];
        Complex oldc = null;
        for (int j = n-i-1; j >= 0; j--) {
            oldc = c[j];
            c[j] = newc;
            newc = oldc.add(newc.multiply(root[i]));
        }
        iterationCount += this.iterationCount;
    }

    resultComputed = true;
    this.iterationCount = iterationCount;
    return root;
}
Example 11
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 12
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 13
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 14
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 15
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 16
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 17
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 18
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 19
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}
Example 20
Source File: FastFourierTransformer.java    From astor with GNU General Public License v2.0 4 votes vote down vote up 
/**
 * Perform the base-4 Cooley-Tukey FFT algorithm (including inverse).
 *
 * @param data the complex data array to be transformed
 * @return the complex transformed array
 * @throws IllegalArgumentException if any parameters are invalid
 */
protected Complex[] fft(Complex data[])
    throws IllegalArgumentException {

    final int n = data.length;
    final Complex f[] = new Complex[n];

    // initial simple cases
    verifyDataSet(data);
    if (n == 1) {
        f[0] = data[0];
        return f;
    }
    if (n == 2) {
        f[0] = data[0].add(data[1]);
        f[1] = data[0].subtract(data[1]);
        return f;
    }

    // permute original data array in bit-reversal order
    int ii = 0;
    for (int i = 0; i < n; i++) {
        f[i] = data[ii];
        int k = n >> 1;
        while (ii >= k && k > 0) {
            ii -= k; k >>= 1;
        }
        ii += k;
    }

    // the bottom base-4 round
    for (int i = 0; i < n; i += 4) {
        final Complex a = f[i].add(f[i+1]);
        final Complex b = f[i+2].add(f[i+3]);
        final Complex c = f[i].subtract(f[i+1]);
        final Complex d = f[i+2].subtract(f[i+3]);
        final Complex e1 = c.add(d.multiply(Complex.I));
        final Complex e2 = c.subtract(d.multiply(Complex.I));
        f[i] = a.add(b);
        f[i+2] = a.subtract(b);
        // omegaCount indicates forward or inverse transform
        f[i+1] = roots.isForward() ? e2 : e1;
        f[i+3] = roots.isForward() ? e1 : e2;
    }

    // iterations from bottom to top take O(N*logN) time
    for (int i = 4; i < n; i <<= 1) {
        final int m = n / (i<<1);
        for (int j = 0; j < n; j += i<<1) {
            for (int k = 0; k < i; k++) {
                //z = f[i+j+k].multiply(roots.getOmega(k*m));
                final int k_times_m = k*m;
                final double omega_k_times_m_real = roots.getOmegaReal(k_times_m);
                final double omega_k_times_m_imaginary = roots.getOmegaImaginary(k_times_m);
                //z = f[i+j+k].multiply(omega[k*m]);
                final Complex z = new Complex(
                    f[i+j+k].getReal() * omega_k_times_m_real -
                    f[i+j+k].getImaginary() * omega_k_times_m_imaginary,
                    f[i+j+k].getReal() * omega_k_times_m_imaginary +
                    f[i+j+k].getImaginary() * omega_k_times_m_real);

                f[i+j+k] = f[j+k].subtract(z);
                f[j+k] = f[j+k].add(z);
            }
        }
    }
    return f;
}