diff --git a/Main.m b/Main.m index 47d5a19..e52c8b3 100644 --- a/Main.m +++ b/Main.m @@ -4,12 +4,55 @@ # #Last modified: 02/03/2023 08:28:32 # -# -# -# ################# pkg load signal; +clc; +close all; +clear all; + + +signal = csvread ('unknownsignal.csv'); + +figure; +plot(signal); +title("raw Signal"); +xlabel("time"); +ylabel("Amplitude"); + +% set signal of interest +SStart = 100; +SEnd = 400; + +Fs = 300; + +lenDomain = 1 + (SEnd - SStart); + +windowed_signal= zeros(lenDomain); + +% using blackman , we get the signal of interest +windowed_signal = signal(SStart:SEnd) .* blackman(lenDomain)'; + +figure; +plot(SStart: SEnd, windowed_signal); +title("Blackman signal windowing"); +xlabel("samples"); +ylabel("Amplitude"); + +frequencySpectrum(windowed_signal, Fs, 0); + +%spectrogram(windowed_signal, Fs, 1/Fs, 1000*length(signal)/Fs); + +% filter using filter and butter + +[val, ind] = max(windowed_signal); + +figure; +[b, a] = butter(6, 10/Fs); +s = filter(b, a, 10*log10(windowed_signal/windowed_signal(ind))); +plot(50:300, s(50:300)); + + +audiowrite("sound.wav", signal, Fs); - x = csvread (filename) \ No newline at end of file diff --git a/frequencySpectrum.m b/frequencySpectrum.m new file mode 100644 index 0000000..94617d5 --- /dev/null +++ b/frequencySpectrum.m @@ -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)') + diff --git a/sound.wav b/sound.wav new file mode 100644 index 0000000..aa44622 Binary files /dev/null and b/sound.wav differ diff --git a/spectrogram.m b/spectrogram.m new file mode 100644 index 0000000..f01e8b5 --- /dev/null +++ b/spectrogram.m @@ -0,0 +1,38 @@ +function spectrogram(signal, samplingFreq, step_size, window_size) +%%%%%%%%%%%%%%%%%%%%%%% +%function spectrogram(signal, samplingFreq, step_size, window_size) +% ex.: spectrogram(signal, samplingFreq, step_size, window_size) +% +% Task: Plot the spectrogram of a given signal +% +% Inputs: +% -signal: temporal signal to analyse +% -samplingFreq: sampling frequency of the temporal signal +% -step_size: how often the power spectrum will be computed in ms +% -window_size: size of the analysing window in ms +% +% Ouput: None +% +% author: Guillaume Gibert (guillaume.gibert@ecam.fr) +% date: 14/03/2023 +%%%%%%%%%%%%%%%%%%%%%%% + +figure; + subplot(2,1,1); +t=0:1/samplingFreq:length(signal)/samplingFreq-1/samplingFreq; +plot(t, signal'); +xlim([0 length(signal)/samplingFreq-1/samplingFreq]); +ylabel('amplitude (norm. unit)'); + subplot(2,1,2); +step = fix(step_size*samplingFreq/1000); % one spectral slice every step_size ms +window = fix(window_size*samplingFreq/1000); % window_size ms data window +fftn = 2^nextpow2(window); % next highest power of 2 +[S, f, t] = specgram(signal, fftn, samplingFreq, window, window-step); +S = abs(S(2:fftn*4000/samplingFreq,:)); % magnitude in range 0