This commit is contained in:
Malo ROUAULT 2023-07-19 16:25:24 +02:00
parent 6d5a12f6aa
commit 65841c8098
4 changed files with 199 additions and 0 deletions

BIN
.DS_Store vendored Normal file

Binary file not shown.

Binary file not shown.

48
frequencySpectrum.m Normal file
View File

@ -0,0 +1,48 @@
function power = frequencySpectrum(signal, fs)
%%%%%%%%%%%%%%%%%%
%function frequencySpectrum(signal, fs)
%
% Task: Display the power spectrum of a given signal
%
% Input:
% - signal: the input signal to process
% - fs: the sampling rate
%
% Output:
% - power: power spectrum of the signal
%
%
% Guillaume Gibert, guillaume.gibert@ecam.fr
% 25/04/2022
%%%%%%%%%%%%%%%%%%
n = length(signal); % number of samples
y = fft(signal, n);% compute DFT of input signal
power = abs(y).^2/n; % power of the DFT
[val, ind] = max(power); % find the mx value of DFT and its index
% plots
figure;
subplot(1,3,1) % time plot
t=0:1/fs:(n-1)/fs; % time range
plot(t, signal)
xticks(0:0.1*fs:n*fs);
xticklabels(0:0.1:n/fs);
xlabel('Time (s)');
ylabel('Amplitude (a.u.)');
subplot(1,3,2) % linear frequency plot
f = (0:n-1)*(fs/n); % frequency range
plot(f,power, 'b*'); hold on;
plot(f,power, 'r');
xlabel('Frequency (Hz)')
ylabel('Power (a.u.)')
subplot(1,3,3) % log frequency plot
plot(f,10*log10(power/power(ind)));
xlabel('Frequency (Hz)')
ylabel('Power (dB)')

151
retake.m Normal file
View File

@ -0,0 +1,151 @@
clear all
close all
clc
signal = csvread('unknownsignal.csv');
samplingFreq = 650;
t = (0:length(signal)-1)/samplingFreq;
signalDuration = 2;
windowDuration = signalDuration/2;
%%%%%%%%%%%%
%RECTANGULAR
%%%%%%%%%%%%
% create rectangular window
rectangularWin = zeros(1, length(t));
for l_sample=1:windowDuration*samplingFreq
rectangularWin(l_sample + signalDuration*samplingFreq/4) = 1;
end
% plot rectangular window
%~ figure;
%~ plot(t, rectangularWin);
% apply the rectangular window
for l_sample=1:length(t)
signal_rect(l_sample) = signal(l_sample) * rectangularWin(l_sample);
end
% plot signal windowed by rectangular window
%~ figure;
%~ plot(t, signal_rect);
% plot the frequency spectrum of this windowed signal
power_rect = frequencySpectrum(signal_rect, samplingFreq);
power = frequencySpectrum(signal, samplingFreq);
%%%%%%%%%%%%
%HAMMING
%%%%%%%%%%%%
hammingWin = zeros(1, length(t));
for l_sample=1:windowDuration*samplingFreq
hammingWin(l_sample+signalDuration*samplingFreq/4) = (0.5 - 0.5*cos(2*pi*(l_sample)/(signalDuration*samplingFreq/2)));
end
% plot Hamming window
%~ figure;
%~ plot(t, hammingWin);
% apply the Hamming window
for l_sample=1:length(t)
signal_hamming(l_sample) = signal(l_sample) * hammingWin(l_sample);
end
% plot signal windowed by rectangular window
%~ figure;
%~ plot(t, signal_hamming);
% plot the frequency spectrum of this windowed signal
power_hamming = frequencySpectrum(signal_hamming, samplingFreq);
%%%%%%%%%%%%
%HANNING
%%%%%%%%%%%%
hanningWin = zeros(1, length(t));
for l_sample=1:windowDuration*samplingFreq
hanningWin(l_sample+signalDuration*samplingFreq/4) = (0.54 - 0.46*cos(2*pi*(l_sample)/(signalDuration*samplingFreq/2)));
end
% plot Hanning window
%~ figure;
%~ plot(t, hanningWin);
% apply the Hanning window
for l_sample=1:length(t)
signal_hanning(l_sample) = signal(l_sample) * hanningWin(l_sample);
end
% plot signal windowed by rectangular window
%~ figure;
%~ plot(t, signal_hanning);
% plot the frequency spectrum of this windowed signal
power_hanning = frequencySpectrum(signal_hanning, samplingFreq);
%%%%%%%%%%%%
%BLACKMAN
%%%%%%%%%%%%
blackmanWin = zeros(1, length(t));
for l_sample=1:windowDuration*samplingFreq
blackmanWin(l_sample+signalDuration*samplingFreq/4) = (0.42 - 0.5 * cos(2*pi*(l_sample)/(signalDuration*samplingFreq/2)) + 0/08*cos(4*pi*(l_sample)/(windowDuration*samplingFreq/2)));
end
% plot Blackman window
%~ figure;
%~ plot(t, blackmanWin);
% apply the Blackman window
for l_sample=1:length(t)
signal_blackman(l_sample) = signal(l_sample) * blackmanWin(l_sample);
end
% plot signal windowed by rectangular window
%~ figure;
%~ plot(t, signal_blackman);
% plot the frequency spectrum of this windowed signal
power_blackman = frequencySpectrum(signal_blackman, samplingFreq);
%%%%%%%%%%%%
% GLOBAL PLOT
%%%%%%%%%%%%
figure;
plot(t, signal_rect, 'r'); hold on;
plot(t, signal_hamming, 'b');
plot(t, signal_hanning, 'g');
plot(t, signal_blackman, 'k');
xlabel('time (s)');
ylabel('amplitude (a.u.)');
legend('Rectangular', 'Hamming', 'Hanning', 'Blackman');
title('Temporal variation of a windowed cosine signal');
figure;
n = length(t);
f = (0:n-1)*(samplingFreq/n); % frequency range
plot(f,10*log10(power_rect/max(power_rect))); hold on;
plot(f,10*log10(power_hamming/max(power_hamming)));
plot(f,10*log10(power_hanning/max(power_hanning)));
plot(f,10*log10(power_blackman/max(power_blackman)));
xlim([0 20]);
ylim([-100 0]);
legend('Rectangular', 'Hamming', 'Hanning', 'Blackman');
xlabel('Frequency (Hz)')
ylabel('Power (dB)')