Skip to content

Instantly share code, notes, and snippets.

@gustafsson
Created May 31, 2012 17:00
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save gustafsson/2844750 to your computer and use it in GitHub Desktop.
Save gustafsson/2844750 to your computer and use it in GitHub Desktop.
Illustration of sinusoid in time and frequency domain (Octave script)
function freqtime(f)
% The fourier transform below assumes that y is a periodic signal which it is if 'f' is an integer (i.e we don't do any windowing)
if 0==nargin
f = 2;
end
A = 1;
fs = 100;
T = 20; % Longer signal -> higher frequency resolution
t = linspace(0,T,T*fs+1);
y0 = A*sin(t*f*2*pi);
hz = (0:length(y0)-1)/T;
noise=[0 3*A];
for ni=1:length(noise)
y = y0 - noise(ni).*(rand(size(y0))-0.5);
F = fft(y)*2/numel(t);
subplot(4,2,ni);
time_to_plot=min(2,T);
ti = 1:time_to_plot*fs;
plot(t(ti),y(ti));
title(['Time domain of y_' num2str(ni) ]);
ylabel('y');
xlabel('t [s]');
ylim([-1 1]*max(noise));
subplot(4,2,2+ni);
fi = 1:ceil(length(F)/2); % Plot the entire spectra but skip the redundant part
plot(hz(fi),abs(F(fi)));
title('Frequency domain');
ylabel(['|F\{y_' num2str(ni) '\}|']);
xlabel('f [Hz]');
ylim([0 .5]*max(noise));
subplot(4,2,4+ni);
max_hz = 6; % Just plot the first 6 Hz
fi = 1+ (0:max_hz*T);
plot(hz(fi),abs(F(fi)));
title(['Frequency domain from 0 Hz to ' num2str(max_hz) ' Hz']);
ylabel(['|F\{y_' num2str(ni) '\}|']);
xlabel('f [Hz]');
ylim([0 .5]*max(noise));
subplot(4,2,6+ni);
spectraLines = 40;
spectra = zeros(spectraLines,numel(t));
spectra(1,:) = F;
for li = 2:spectraLines
% Create more signal data similiar to the first part
y = y0 - noise(ni).*(rand(size(y0))-0.5);
spectra(li,:) = fft(y)*2/numel(t);
end
max_hz = 6; % Just plot the first 6 Hz
fi = 1+ (0:max_hz*T);
spectra = abs(spectra(:,fi));
C = zeros(256,3);
C(:,2) = linspace(0,1,256);
colormap( C );
imagesc(hz(fi),((1:spectraLines)*T-0.5*T)/60,(spectra./max(spectra(:))).^0.25);
title('Spectrogram of continued signal');
ylabel('t [min]');
xlabel('f [Hz]');
end
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment