Using Findpeaks to plot on a scatter plot and return the peak coordinates
3 views (last 30 days)
Show older comments
Taylor Knuth
on 20 Nov 2021
Answered: Star Strider
on 20 Nov 2021
I used a FFT to transform my data from the time domain to the frequency domain. I also filtered the data to create localized frequency spikes, I am trying to identify the amplitude and frequency of the spikes. but when I use findpeaks it will not plot on my existing plot or show the correct peaks. Do I need to call the findpeaks function by a variable name so that I can access a matrix with the peak coordinates?
I have attached my code and data for your views. Thanks so much for your ideas.
%% Universal Constants
Fs=125; %to be checked Frequency Spectrum (62.5)
dt=1/Fs;
NFFT=800; %fft length (1024) changes the filter average
OVERLAP = .75; % between 0 and 1 max ; buffer overlap = OVERLAP*NFFT
%% Load Files
WBR_signal=load('WB0430Right.txt');
%% Samples and Channels for size
[R_samples_WB,R_channels_WB]=size(WBR_signal);
%% Selected Channels
selected_channel_AX=3;
selected_channel_AY=4;
selected_channel_AZ=5;
selected_channel_GX=6;
selected_channel_GY=7;
selected_channel_GZ=8;
selected_channel_MX=9;
selected_channel_MY=10;
selected_channel_MZ=11;
%% Samples and Channels for size
%WorkBoot
%Right
W_signal_RAX=WBR_signal(:,selected_channel_AX);
W_signal_RAY=WBR_signal(:,selected_channel_AY);
W_signal_RAZ=WBR_signal(:,selected_channel_AZ);
W_signal_RGX=WBR_signal(:,selected_channel_GX);
W_signal_RGY=WBR_signal(:,selected_channel_GY);
W_signal_RGZ=WBR_signal(:,selected_channel_GZ);
W_signal_RMX=WBR_signal(:,selected_channel_MX);
W_signal_RMY=WBR_signal(:,selected_channel_MY);
W_signal_RMZ=WBR_signal(:,selected_channel_MZ);
%% Time
%Work Boot
WR_time=(0:R_samples_WB-1)*dt;
%% Myfft specification
%Work Boot
%right
[W_freq_RAX,W_sensor_spectrum_RAX] = myfft_peak(W_signal_RAX,Fs,NFFT,OVERLAP);
[W_freq_RAY,W_sensor_spectrum_RAY] = myfft_peak(W_signal_RAY,Fs,NFFT,OVERLAP);
[W_freq_RAZ,W_sensor_spectrum_RAZ] = myfft_peak(W_signal_RAZ,Fs,NFFT,OVERLAP);
[W_freq_RGX,W_sensor_spectrum_RGX] = myfft_peak(W_signal_RGX,Fs,NFFT,OVERLAP);
[W_freq_RGY,W_sensor_spectrum_RGY] = myfft_peak(W_signal_RGY,Fs,NFFT,OVERLAP);
[W_freq_RGZ,W_sensor_spectrum_RGZ] = myfft_peak(W_signal_LGZ,Fs,NFFT,OVERLAP);
[W_freq_RMX,W_sensor_spectrum_RMX] = myfft_peak(W_signal_RMX,Fs,NFFT,OVERLAP);
[W_freq_RMY,W_sensor_spectrum_RMY] = myfft_peak(W_signal_RMY,Fs,NFFT,OVERLAP);
[W_freq_RMZ,W_sensor_spectrum_RMZ] = myfft_peak(W_signal_RMZ,Fs,NFFT,OVERLAP);
%% Figures
figure
plot(W_freq_RMZ,W_sensor_spectrum_RMZ/12.63)
xlabel('Frequency (Hz)')
ylabel('Amplitude')
ylim([0,1])
%Find Peaks
x=W_freq_RMZ;
y=W_sensor_spectrum_RMZ/12.63;
[peaks,locs]=findpeaks(y,x);
figure
z=plot(x,y,'b');
grid on
xlabel('Frequency (Hz)')
ylabel('Amplitude')
title('Peaks')
axis tight
hold on
scatter(peaks,locs,'^r');
%% Function
function [freq_vector,fft_spectrum] = myfft_peak(signal, Fs, nfft, Overlap)
% FFT peak spectrum of signal (example sinus amplitude 1 = 0 dB after fft).
% Linear averaging
% signal - input signal,
% Fs - Sampling frequency (Hz).
% nfft - FFT window size
% Overlap - buffer percentage of overlap % (between 0 and 0.95)
[samples,channels] = size(signal);
% fill signal with zeros if its length is lower than nfft
if samples<nfft
s_tmp = zeros(nfft,channels);
s_tmp((1:samples),:) = signal;
signal = s_tmp;
samples = nfft;
end
% window : hanning
window = hanning(nfft);
window = window(:);
% compute fft with overlap
offset = fix((1-Overlap)*nfft);
spectnum = 1+ fix((samples-nfft)/offset); % Number of windows
% % for info is equivalent to :
% noverlap = Overlap*nfft;
% spectnum = fix((samples-noverlap)/(nfft-noverlap)); % Number of windows
% main loop
fft_spectrum = 0;
for i=1:spectnum
start = (i-1)*offset;
sw = signal((1+start):(start+nfft),:).*(window*ones(1,channels));
fft_spectrum = fft_spectrum + (abs(fft(sw))*4/nfft); % X=fft(x.*hanning(N))*4/N; % hanning only
end
fft_spectrum = fft_spectrum/spectnum; % to do linear averaging scaling
% one sidded fft spectrum % Select first half
if rem(nfft,2) % nfft odd
select = (1:(nfft+1)/2)';
else
select = (1:nfft/2+1)';
end
fft_spectrum = fft_spectrum(select,:);
freq_vector = (select - 1)*Fs/nfft;
end
0 Comments
Accepted Answer
Star Strider
on 20 Nov 2021
I have not run the code, however I believe the scatter call used to plot the peaks has the arguments reversed.
Try this instead —
scatter(locs,peaks,'^r');
That should work.
.
0 Comments
More Answers (0)
See Also
Categories
Find more on Transforms in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!