| << Chapter < Page | Chapter >> Page > |
| Listing 21. ForwardRealToComplex01.java. |
|---|
/*File ForwardRealToComplex01.java
Copyright 2004, R.G.BaldwinRev 5/14/04
The static method named transform performs a realto complex Fourier transform.
Does not implement the FFT algorithm. Implementsa straight-forward sampled-data version of the
continuous Fourier transform defined usingintegral calculus. See ForwardRealToComplexFFT01
for an FFT algorithm.Returns real, imag, magnitude, and phase angle in
degrees.Incoming parameters are:
double[]data - incoming real data
double[]realOut - outgoing real data
double[]imagOut - outgoing imaginary data
double[]angleOut - outgoing phase angle in
degreesdouble[] magnitude - outgoing amplitudespectrum
int zero - the index of the incoming datasample that represents zero time
double lowF - Low freq limit as fraction ofsampling frequency
double highF - High freq limit as fraction ofsampling frequency
The frequency increment is the difference betweenhigh and low limits divided by the length of
the magnitude arrayThe magnitude is computed as the square root of
the sum of the squares of the real and imaginaryparts. This value is divided by the incoming
data length, which is given by data.length.Returns a number of points in the frequency
domain equal to the incoming data lengthregardless of the high and low frequency
limits.************************************************/
public class ForwardRealToComplex01{public static void transform(
double[]data,
double[]realOut,
double[]imagOut,
double[]angleOut,
double[]magnitude,
int zero,double lowF,
double highF){double pi = Math.PI;//for convenience
int dataLen = data.length;double delF = (highF-lowF)/data.length;
//Outer loop iterates on frequency// values.
for(int i=0; i<dataLen;i++){
double freq = lowF + i*delF;double real = 0.0;
double imag = 0.0;double ang = 0.0;
//Inner loop iterates on time-// series points.
for(int j=0; j<dataLen; j++){
real += data[j]*Math.cos(
2*pi*freq*(j-zero));imag += data[j]*Math.sin(2*pi*freq*(j-zero));
}//end inner looprealOut[i] = real/dataLen;imagOut[i] = imag/dataLen;magnitude[i] = (Math.sqrt(real*real + imag*imag))/dataLen;
//Calculate and return the phase// angle in degrees.
if(imag == 0.0&&real == 0.0){ang = 0.0;}
else{ang = Math.atan(imag/real)*180.0/pi;}if(real<0.0&&imag == 0.0){ang = 180.0;}
else if(real<0.0&&imag == -0.0){
ang = -180.0;}else if(real<0.0&&imag>0.0){
ang += 180.0;}else if(real<0.0&&imag<0.0){
ang += -180.0;}angleOut[i] = ang;}//end outer loop
}//end transform method}//end class ForwardRealToComplex01 |
| Listing 22. Dsp030.java. |
|---|
/* File Dsp030.java
Copyright 2004, R.G.BaldwinRev 5/14/04
Uses an FFT algorithm to compute and display themagnitude of the spectral content for up to five
sinusoids having different frequencies andamplitudes.
See the program named Dsp028 for a program thatdoes not use an FFT algorithm.
Gets input parameters from a file namedDsp030.txt. If that file doesn't exist in the
current directory, the program uses a set ofdefault parameters.
Each parameter value must be stored as characterson a separate line in the file named Dsp030.txt.
The required parameters are as follows:Data length as type int
Number of spectra as type int. Max value is 5.List of sinusoid frequency values as type double.
List of sinusoid amplitude values as type double.CAUTION: THE DATA LENGTH MUST BE A POWER OF TWO.
OTHERWISE, THIS PROGRAM WILL FAIL TO RUNPROPERLY.
The number of values in each of the lists mustmatch the value for the number of spectra.
Note: All frequency values are specified as adouble representing a fractional part of the
sampling frequency. For example, a value of 0.5specifies a frequency that is half the sampling
frequency.Here is a set of sample parameter values. Don't
allow blank lines at the end of the data in thefile.
128.05
0.10.2
0.30.4
0.4560
7080
90100
The plotting program that is used to plot theoutput data from this program requires that the
program implement GraphIntfc01. For example,the plotting program named Graph03 can be used
to plot the data produced by this program. Whenit is used, the usage information is:
java Graph03 Dsp030A static method named transform belonging to the
class named ForwardRealToComplexFFT01 is used toperform the actual spectral analysis. The method
named transform implements an FFT algorithm. TheFFT algorithm requires that the data length be
a power of two.Tested using SDK 1.4.2 under WinXP.
************************************************/import java.util.*;
import java.io.*;class Dsp030 implements GraphIntfc01{
final double pi = Math.PI;//for simplification//Begin default parameters
int len = 128;//data lengthint numberSpectra = 5;
//Frequencies of the sinusoidsdouble[] freq = {0.1,0.2,0.3,0.4,0.5};//Amplitudes of the sinusoidsdouble[] amp = {60,70,80,90,100};//End default parameters
//Following arrays will contain data that is// input to the spectral analysis process.
double[]data1;
double[]data2;
double[]data3;
double[]data4;
double[]data5;
//Following arrays receive information back// from the spectral analysis that is not used
// in this program.double[] real;double[] imag;double[] angle;//Following arrays receive the magnitude
// spectral information back from the spectral// analysis process.
double[]magnitude1;
double[]magnitude2;
double[]magnitude3;
double[]magnitude4;
double[]magnitude5;
public Dsp030(){//constructor//Get the parameters from a file named
// Dsp030.txt. Use the default parameters// if the file doesn't exist in the current
// directory.if(new File("Dsp030.txt").exists()){
getParameters();}//end if
//Note that this program always processes// five sinusoids, even if fewer than five
// were requested as the input parameter// for numberSpectra. In that case, the
// extras are processed using default values// and simply ignored when the results are
// plotted.//Create the raw data. Note that the
// argument for a sinusoid at half the// sampling frequency would be (2*pi*x*0.5).
// This would represent one half cycle or pi// radians per sample.
//First create empty array objects.double[] data1 = new double[len];
double[]data2 = new double[len];double[] data3 = new double[len];
double[]data4 = new double[len];double[] data5 = new double[len];
//Now populate the array objectsfor(int n = 0;n<len;n++){
data1[n]= amp[0]*Math.cos(2*pi*n*freq[0]);
data2[n]= amp[1]*Math.cos(2*pi*n*freq[1]);
data3[n]= amp[2]*Math.cos(2*pi*n*freq[2]);
data4[n]= amp[3]*Math.cos(2*pi*n*freq[3]);
data5[n]= amp[4]*Math.cos(2*pi*n*freq[4]);
}//end for loop//Compute magnitude spectra of the raw data
// and save it in output arrays. Note that// the real, imag, and angle arrays are not
// used later, so they are discarded each// time a new spectral analysis is performed.
magnitude1 = new double[len];
real = new double[len];
imag = new double[len];
angle = new double[len];
ForwardRealToComplexFFT01.transform(data1,real,imag,angle,magnitude1);
magnitude2 = new double[len];
real = new double[len];
imag = new double[len];
angle = new double[len];
ForwardRealToComplexFFT01.transform(data2,real,imag,angle,magnitude2);
magnitude3 = new double[len];
real = new double[len];
imag = new double[len];
angle = new double[len];
ForwardRealToComplexFFT01.transform(data3,real,imag,angle,magnitude3);
magnitude4 = new double[len];
real = new double[len];
imag = new double[len];
angle = new double[len];
ForwardRealToComplexFFT01.transform(data4,real,imag,angle,magnitude4);
magnitude5 = new double[len];
real = new double[len];
imag = new double[len];
angle = new double[len];
ForwardRealToComplexFFT01.transform(data5,real,imag,angle,magnitude5);
}//end constructor//-------------------------------------------//
//This method gets processing parameters from// a file named Dsp030.txt and stores those
// parameters in instance variables belonging// to the object of type Dsp030.
void getParameters(){int cnt = 0;
//Temporary holding area for strings. Allow// space for a few blank lines at the end
// of the data in the file.String[] data = new String[20];
try{//Open an input stream.
BufferedReader inData =new BufferedReader(new FileReader(
"Dsp030.txt"));//Read and save the strings from each of
// the lines in the file. Be careful to// avoid having blank lines at the end,
// which may cause an ArrayIndexOutOfBounds// exception to be thrown.
while((data[cnt]=
inData.readLine()) != null){cnt++;
}//end whileinData.close();
}catch(IOException e){}//Move the parameter values from the
// temporary holding array into the instance// variables, converting from characters to
// numeric values in the process.cnt = 0;
len = (int)Double.parseDouble(data[cnt++]);
numberSpectra = (int)Double.parseDouble(data[cnt++]);for(int fCnt = 0;fCnt<numberSpectra;
fCnt++){freq[fCnt] = Double.parseDouble(data[cnt++]);}//end for loop
for(int aCnt = 0;aCnt<numberSpectra;
aCnt++){amp[aCnt] = Double.parseDouble(data[cnt++]);}//end for loop
//Print parameter values.System.out.println();
System.out.println("Data length: " + len);System.out.println(
"Number spectra: " + numberSpectra);System.out.println("Frequencies");
for(cnt = 0;cnt<numberSpectra;cnt++){
System.out.println(freq[cnt]);
}//end for loopSystem.out.println("Amplitudes");
for(cnt = 0;cnt<numberSpectra;cnt++){
System.out.println(amp[cnt]);
}//end for loop}//end getParameters
//-------------------------------------------////The following six methods are required by the
// interface named GraphIntfc01. The plotting// program pulls the data values to be plotted
// by calling these methods.public int getNmbr(){
//Return number of functions to// process. Must not exceed 5.
return numberSpectra;}//end getNmbr
//-------------------------------------------//public double f1(double x){
int index = (int)Math.round(x);if(index<0 ||
index>magnitude1.length-1){
return 0;}else{
return magnitude1[index];
}//end else}//end function
//-------------------------------------------//public double f2(double x){
int index = (int)Math.round(x);if(index<0 ||
index>magnitude2.length-1){
return 0;}else{
return magnitude2[index];
}//end else}//end function
//-------------------------------------------//public double f3(double x){
int index = (int)Math.round(x);if(index<0 ||
index>magnitude3.length-1){
return 0;}else{
return magnitude3[index];
}//end else}//end function
//-------------------------------------------//public double f4(double x){
int index = (int)Math.round(x);if(index<0 ||
index>magnitude4.length-1){
return 0;}else{
return magnitude4[index];
}//end else}//end function
//-------------------------------------------//public double f5(double x){
int index = (int)Math.round(x);if(index<0 ||
index>magnitude5.length-1){
return 0;}else{
return magnitude5[index];
}//end else}//end function
//-------------------------------------------//}//end class Dsp030 |
| Listing 23. ForwardRealToComplexFFT01.java. |
|---|
/*File ForwardRealToComplexFFT01.java
Copyright 2004, R.G.BaldwinRev 5/14/04
The static method named transform performs a realto complex Fourier transform using a crippled
complex-to-complex FFT algorithm. It is crippledin the sense that it is not being used to its
full potential as a complex-to-complex forwardor inverse FFT algorithm.See ForwardRealToComplex01 for a slower but more
general approach that does not use an FFTalgorithm.
Returns real, imag, magnitude, and phase angle indegrees.
Incoming parameters are:double[] data - incoming real datadouble[] realOut - outgoing real datadouble[] imagOut - outgoing imaginary datadouble[] angleOut - outgoing phase angle indegrees
double[]magnitude - outgoing amplitude
spectrumThe magnitude is computed as the square root of
the sum of the squares of the real and imaginaryparts. This value is divided by the incoming
data length, which is given by data.length.CAUTION: THE INCOMING DATA LENGTH MUST BE A
POWER OF TWO. OTHERWISE, THIS PROGRAM WILL FAILTO RUN PROPERLY.
Returns a number of points in the frequencydomain equal to the incoming data length. Those
points are uniformly distributed between zero andone less than the sampling frequency.
************************************************/public class ForwardRealToComplexFFT01{
public static void transform(double[] data,double[] realOut,double[] imagOut,double[] angleOut,double[] magnitude){double pi = Math.PI;//for convenience
int dataLen = data.length;//The complexToComplex FFT method does an
// in-place transform causing the output// complex data to be stored in the arrays
// containing the input complex data.// Therefore, it is necessary to copy the
// input data to this method into the real// part of the complex data passed to the
// complexToComplex method.System.arraycopy(data,0,realOut,0,dataLen);
//Perform the spectral analysis. The results// are stored in realOut and imagOut. The +1
// causes it to be a forward transform. A -1// would cause it to be an inverse transform.
complexToComplex(1,dataLen,realOut,imagOut);//Compute the magnitude and the phase angle
// in degrees.for(int cnt = 0;cnt<dataLen;cnt++){
magnitude[cnt]=
(Math.sqrt(realOut[cnt]*realOut[cnt]
+ imagOut[cnt]*imagOut[cnt]))/dataLen;if(imagOut[cnt] == 0.0&&realOut[cnt] == 0.0){angleOut[cnt] = 0.0;}//end if
else{angleOut[cnt] = Math.atan(imagOut[cnt]/realOut[cnt])*180.0/pi;
}//end elseif(realOut[cnt]<0.0&&imagOut[cnt] == 0.0){angleOut[cnt] = 180.0;}else if(realOut[cnt]<0.0&&imagOut[cnt] == -0.0){angleOut[cnt] = -180.0;}else if(realOut[cnt]<0.0&&imagOut[cnt]>0.0){
angleOut[cnt]+= 180.0;
}else if(realOut[cnt]<0.0&&imagOut[cnt]<0.0){
angleOut[cnt]+= -180.0;
}//end else}//end for loop
}//end transform method//-------------------------------------------//
//This method computes a complex-to-complex// FFT. The value of sign must be 1 for a
// forward FFT.public static void complexToComplex(
int sign,int len,
double real[],
double imag[]){
double scale = 1.0;//Reorder the input data into reverse binary
// order.int i,j;
for (i=j=0; i<len; ++i) {
if (j>=i) {
double tempr = real[j]*scale;
double tempi = imag[j]*scale;
real[j]= real[i]*scale;imag[j] = imag[i]*scale;
real[i]= tempr;
imag[i]= tempi;
}//end ifint m = len/2;
while (m>=1&&j>=m) {
j -= m;m /= 2;
}//end while loopj += m;
}//end for loop//Input data has been reordered.
int stage = 0;int maxSpectraForStage,stepSize;
//Loop once for each stage in the spectral// recombination process.
for(maxSpectraForStage = 1,stepSize = 2*maxSpectraForStage;
maxSpectraForStage<len;
maxSpectraForStage = stepSize,stepSize = 2*maxSpectraForStage){
double deltaAngle =sign*Math.PI/maxSpectraForStage;
//Loop once for each individual spectrafor (int spectraCnt = 0;
spectraCnt<maxSpectraForStage;
++spectraCnt){double angle = spectraCnt*deltaAngle;
double realCorrection = Math.cos(angle);double imagCorrection = Math.sin(angle);
int right = 0;for (int left = spectraCnt;
left<len;left += stepSize){
right = left + maxSpectraForStage;double tempReal =
realCorrection*real[right]- imagCorrection*imag[right];double tempImag =
realCorrection*imag[right]+ imagCorrection*real[right];real[right] = real[left]-tempReal;
imag[right]= imag[left]-tempImag;real[left] += tempReal;imag[left] += tempImag;}//end for loop
}//end for loop for individual spectramaxSpectraForStage = stepSize;
}//end for loop for stages}//end complexToComplex method
}//end class ForwardRealToComplexFFT01 |
This section contains a variety of miscellaneous information.
Baldwin explains several different programs used for spectral analysis. He also explains the impact of the sampling frequency and the Nyquist folding frequency on spectral analysis.
Financial : Although the Connexions site makes it possible for you to download a PDF file for thismodule at no charge, and also makes it possible for you to purchase a pre-printed version of the PDF file, you should beaware that some of the HTML elements in this module may not translate well into PDF.
I also want you to know that, I receive no financial compensation from the Connexions website even if you purchase the PDF version of the module.
In the past, unknown individuals have copied my modules from cnx.org, converted them to Kindle books, and placed them for sale on Amazon.com showing me as the author. Ineither receive compensation for those sales nor do I know who does receive compensation. If you purchase such a book, please beaware that it is a copy of a module that is freely available on cnx.org and that it was made and published withoutmy prior knowledge.
Affiliation ::: I am a professor of Computer Information Technology at Austin Community College in Austin, TX.
-end-
Notification Switch
Would you like to follow the 'Digital signal processing - dsp' conversation and receive update notifications?
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|