added files

This commit is contained in:
Guillaume LE CHARTIER 2024-03-18 09:27:00 +01:00
parent 8e22a60da3
commit 14ca4cc845
3 changed files with 100 additions and 0 deletions

61
frequencySpectrum.m Normal file
View File

@ -0,0 +1,61 @@
function power = frequencySpectrum(signal, fs, resolution)
%%%%%%%%%%%%%%%%%%
%function power = frequencySpectrum(signal, fs, pad)
%
% Task: Display the power spectrum (lin and log scale) of a given signal
%
% Input:
% - signal: the input signal to process
% - fs: the sampling rate in Hz
% - resolution: frequency resolution in Hz, signal will be padded with zeros if necessary
%
% Output:
% - power: the power spectrum
%
%
% Guillaume Gibert, guillaume.gibert@ecam.fr
% 15/03/2024
%%%%%%%%%%%%%%%%%%
n = length(signal); % number of samples
current_resolution = fs / n;
if (resolution < current_resolution)
n_original = n;
n = fs / resolution;
signal = [signal zeros(1, n-n_original)];
end
%~ if (pad)
%~ n_original = n;
%~ n = 2^(nextpow2(n));
%~ signal = [signal zeros(1, n-n_original)];
%~ end
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)')

38
spectrogram.m Normal file
View File

@ -0,0 +1,38 @@
function spectrogram(signal, samplingFreq, step_size, window_size)
%%%%%%%%%%%%%%%%%%%%%%%
%function spectrogram(signal, samplingFreq, step_size, window_size)
% ex.: spectrogram(signal, 300, 50, 1000)
%
% 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*samplingFreq/2/samplingFreq,:)); % magnitude in range 0<f<=4000 Hz.
S = S/max(S(:)); % normalize magnitude so that max is 0 dB.
%S = max(S, 10^(-40/10)); % clip below -40 dB.
%S = min(S, 10^(-3/10)); % clip above -3 dB.
imagesc (t, f, log(S)); % display in log scale
set (gca, "ydir", "normal"); % put the 'y' direction in the correct direction
xlabel('time (s)');
ylabel('frequency (Hz)');

1
unknownsignal.csv Normal file

File diff suppressed because one or more lines are too long