xfft_cmpr.m

fft c源程序

%function fft_cmpr;
ax_font_size = 22;
ax_lab_size = 24;
fmt = '-tiff -deps';
print_flag=0;

n = 16;

fname_src = input('FFT input data: ','s');
fname_xfft = input('FFT data: ','s');

sprintf('path = %s\n',fname_src)
sprintf('path = %s\n',fname_xfft)

fp=fopen(fname_src,'r');
x=fscanf(fp,'%f');
fclose(fp);

fp=fopen(fname_xfft,'r');
sim_data=fscanf(fp,'%f');
fclose(fp);

% asemble complex input vector

u = (x(1:2:2*n) + j*x(2:2:2*n))/32768.0;

% compute the reference FFT

ref_fft = fft(u,n)/(n);

% form complex series from simulation results

xsim_data = sim_data(1:2:length(sim_data)) + j*sim_data(2:2:length(sim_data));
xsim_data = xsim_data/32768.0;

figure(1)

b = 0:length(xsim_data)-1;

plot(b,real(xsim_data),b,real(ref_fft),'m')
grid
ax=axis; ax(1)=0; ax(2)=n-1; axis(ax);
title('TRANSFORM OVERLAY - REAL');

figure(2)

plot(b,imag(xsim_data),b,imag(ref_fft),'g')
grid
ax=axis; ax(1)=0; ax(2)=n-1; axis(ax);
title('TRANSFORM OVERLAY - IMAG.');

% RMS error calculations for the real and imaginary parts of
% the result data

er = real(ref_fft)-real(xsim_data);
ei = imag(ref_fft)-imag(xsim_data);

n = length(er);

qr = sqrt(1/n*sum(er.*er));
qi = sqrt(1/n*sum(ei.*ei));

q = 0:length(er)-1;
figure(3)
plot(q,er,'g')
grid
ax=axis; ax(1)=0;; ax(2)=n-1; axis(ax);
title('ERROR SIGNAL - REAL')

figure(4)
plot(q,ei,'b')
grid
ax=axis; ax(1)=0;; ax(2)=n-1;; axis(ax);
title('ERROR SIGNAL - IMAG')

figure(5)
x = (-n/2:n/2-1)/n; % horizontal axis in Hz
mx = max(abs(xsim_data));
H=20*log10(abs(xsim_data)/mx);
plot(x,[H(n/2+1:n); H(1:n/2)],'b');
grid
xlabel('NORMALIZED FREQUENCY')
ylabel('LOG MAGNITUDE RESPONSE (DB)')
title('SIMULATION PSD')

figure(6)
mx = max(abs(ref_fft));
G=20*log10(abs(ref_fft)/mx);
plot(x,[G(n/2+1:n);G(1:n/2)],'r');
grid
xlabel('NORMALIZED FREQUENCY')
ylabel('LOG MAGNITUDE RESPONSE (DB)')
title('REFERENCE PSD')

figure(7)
p0=angle(xsim_data);
plot(x,[p0(n/2+1:n); p0(1:n/2)],'b');
grid
xlabel('NORMALIZED FREQUENCY')
ylabel('ANGLE')
title('SIMULATION PSD');

figure(8)
p1=angle(ref_fft);
plot(x,[p1(n/2+1:n); p1(1:n/2)],'r');
grid
xlabel('NORMALIZED FREQUENCY')
ylabel('ANGLE')
title('REFERENCE PSD');

figure(9)
subplot(2,1,1); plot(b,real(xsim_data)); grid;
ax=axis; ax(1)=0; ax(2)=n-1;-1; axis(ax);
title('TRANSFORM - REAL (xFFT)');
subplot(2,1,2); plot(b,real(ref_fft),'m'); grid
ax=axis; ax(1)=0; ax(2)=n-1;-1; axis(ax);
title('TRANSFORM - REAL (REFERENCE)');

figure(10)
subplot(2,1,1); plot(b,imag(xsim_data)); grid;
ax=axis; ax(1)=0; ax(2)=n-1;-1; axis(ax);
title('TRANSFORM - IMAG. (xFFT)');
subplot(2,1,2); plot(b,imag(ref_fft),'m'); 
ax=axis; ax(1)=0; ax(2)=n-1;-1; axis(ax);
grid
title('TRANSFORM - IMAG. (REFERENCE)');

figure(11)
plot(real(u)); grid; axis([1 length(x) -1 1]);
xlabel('SAMPLE INDEX');
ylabel('MAGNITUDE');
title('INPUT DATA - REAL');
figure(12)
plot(imag(u)); grid; axis([1 length(x) -1 1]);
xlabel('SAMPLE INDEX');
ylabel('MAGNITUDE');
title('INPUT DATA - IMAG');