// file: $ISIP_IFC/class/java/SpeechDisplay/Fft.java // version: $Id: Fft.java 10398 2006-01-27 19:40:39Z stanley $// // import necessary java libraries // import java.awt.*; import java.awt.event.*; import java.awt.image.BufferedImage; import java.awt.geom.Line2D; import java.awt.geom.Rectangle2D; import java.net.*; import java.io.*; import java.awt.font.*; import java.awt.geom.*; import javax.swing.*; import java.text.DecimalFormat; import java.lang.*; import java.util.Random; /** * class: Fft.java * * this class has the methods to compute a radix-2 cooley-tukey * decimation-in-frequency algorithm. These methods have beeen * imported from the display tcl/tk code developed at isip. * */ public class Fft { //------------------------------------------------------------------------- // // public data // //------------------------------------------------------------------------- /** * the constant that stores the fft length */ public static int n_d; /** * the array that stores the real weights */ public float rad2_wr_d[]; /** * the array that stores the imaginary weights */ public float rad2_wi_d[]; /** * the array that stores the imaginary values */ public float rad2_temp_d[]; /** * the array that stores the real values */ public float output_a[]; //------------------------------------------------------------------------- // // public methods // //------------------------------------------------------------------------- /** * a dummy constructor. * */ public Fft() { } /** * this function intialises the real and imaginary weight array * needed for the fft calcualtion. * * @param n_a the fft length * * @return returns back a boolean value **/ public boolean rad2_init_c(int n_a) { // declare the local variables // int i; float phi = 0; // calculate the weight // float wn = (float)(2.0*Math.PI/n_a); // the array of real and imaginary weights // rad2_wr_d = new float[n_a]; rad2_wi_d = new float[n_a]; // calcualtes the weights and assigns it to appropriate // indeces on the weight table // for (i = 0; i < n_a; i++) { phi = wn * i; rad2_wi_d[i] = (float)Math.sin(phi); rad2_wr_d[i] = (float)Math.cos(phi); } rad2_temp_d = new float[n_a]; // return a boolean value // return true ; } /** * this method performs a radix-2 cooley-tukey decimation-in-frequency * algorithm. the code follows the exact structure of that in the * book: J.G.Proakis, D.G.Manolakis, "Digital Signal Processing - * Principles, Algorithms, and Applications" 2nd Edition, Macmillan * Publishing Company, New York, pp. 731-732, 1992. This function has * been taken from the tcl/tk code. The original code has been * modified so that the method returns only the first half of the fft * values and the mirror image is not calculated. * * @param input_a a window of data * * @return returns back the fft values. **/ public float [] rad2_real_cc( float [] input_a){ // declare local variables // int i,j,k,l,m,i1,i2; int n1,n2; float real_twid,imag_twid; float tempr,tempi; // get the length of the fft to be calculated // n_d = input_a.length; // initilaise the ouput array // output_a = new float [n_d]; // initialize the lookup tables // for (k = 0; k < n_d; ++k) { output_a[k] = input_a[k]; rad2_temp_d[k] = 0; } n2 = n_d; // looping through the first stage // n1 = n2; n2 = n2 >> 1; i1 = 0; i2 = ((n_d)/n1); // this goes through the first set of butterfly calcualtion // for (j = 0; j < n2; ++j) { // loop through each butterfly // real_twid = rad2_wr_d[i1]; imag_twid = rad2_wi_d[i1]; i1 = i1 + i2; for (k = j; k < n_d; k += n1) { l = k + n2; tempr = output_a[k] - output_a[l]; output_a[k] = output_a[k] + output_a[l]; output_a[l] = (real_twid * tempr); rad2_temp_d[l] = -(imag_twid * tempr); } } // get the number of layers of butterfly calcualtion. For a 8 // point fft, we have three stages of fft calculation // m =(int) Math.round((Math.log(n_d)/(float)Math.log(2))); // after the first stage of butterfly calculation is over, // this for loop loops through the remaining stages to get to // get the final fft values // for (i = 1; i < m ; ++i) { n1 = n2; n2 = n2 >> 1; i1 = 0; i2 = ((n_d)/n1); for (j = 0; j < n2; ++j) { // loop through each butterfly // real_twid = rad2_wr_d[i1]; imag_twid = rad2_wi_d[i1]; i1 = i1 + i2; for (k = j; k < n_d; k += n1) { if (i == m - 1) { l = k + n2; output_a[k] = output_a[k] + output_a[l]; rad2_temp_d[k] = rad2_temp_d[k] + rad2_temp_d[l]; } // perform the butterfly computation; computation // within square braces // in eqn (9.3.34) multiplied by the twiddle factor // else { l = k + n2; tempr = output_a[k] - output_a[l]; output_a[k] = output_a[k] + output_a[l]; tempi = rad2_temp_d[k] - rad2_temp_d[l]; rad2_temp_d[k] = rad2_temp_d[k] + rad2_temp_d[l]; output_a[l] = (real_twid * tempr) + (imag_twid * tempi); rad2_temp_d[l] = (real_twid * tempi) - (imag_twid * tempr); } } } } // procedure for bit reversal // j = 0; n1 = n_d - 1; n2 = n_d >> 1; for (i = 0; i < n1 ; ++i) { if (i < j) { tempr = output_a[j]; output_a[j] = output_a[i]; output_a[i] = tempr; tempr = rad2_temp_d[j]; rad2_temp_d[j] = rad2_temp_d[i]; rad2_temp_d[i] = tempr; } k = n2; while (k <= j) { j -= k; k = k >> 1; } j += k; } // interlace the output and store in the output array // int clone_n_d = n_d ; float[] fft_output = new float [n_d/2]; // this for loop is used to calculate the magnitude by just // taking the square root of sum of the squares of the real // and the imaginary parts // for (i = 0; i < n_d/2; i++) { fft_output[i] = (float)(Math.sqrt(output_a[i]* output_a[i] + rad2_temp_d[i] * rad2_temp_d[i])); } // return the fft values // return fft_output; } //------------------------------------------------------------------------- // // main method for this class // //------------------------------------------------------------------------- /** * test the fft program. * * @param args arguments from the command line. */ public static void main(String[] args) { n_d = 512; Fft fft_test = new Fft(); // declare a input array // float [] input_test = new float [n_d]; // initialise the array // for (int i = 0; i < n_d; i++) { input_test[i] = 1; } // initialise the fft object // fft_test.rad2_init_c(n_d); // compute the fft // float [] output_test = fft_test.rad2_real_cc(input_test); } }