Biphasic Pulse signal with FFT
3 views (last 30 days)
Show older comments
I want to create a biphasic pulse signal as used in Functional electrical stimulation, afterwards I would like to do a FFT to get a better understanding of which frequencies would be good for the stimulation. However my code has some errors. I changed my Input variables and now it says:
"Error using plot
Vectors must be the same length.
Error in BiphasicFFT (line 13)
plot(t,wholeSignal);"
What does it mean and how can I change that? I would like to change my parameters, so that I can analyze different stimulation parameters.
Is my FFT correct?
This is my code:
fs=100000; %sampling freq in Hz
pulsesPerS=60; %frequency in Hz
signalDuration=3; % in seconds
amplitude=22*10^-3; % in mA
pulseWidth=250*10^-6; % in µs
t=0:1/fs:signalDuration-1/fs; %vector over the duration of time and sampled every 1/fs
onePeriod=amplitude*[ones(round(pulseWidth*fs),1); -ones(round(pulseWidth*fs),1); zeros(round(1/pulsesPerS*fs)-2*round(pulseWidth*fs),1)];
wholeSignal=repmat(onePeriod,[pulsesPerS*signalDuration 1]);
figure (1)
plot(t,wholeSignal);
xlabel('Time (seconds)');
ylabel('Amplitude');
title('Biphasic Pulse Train')
fourier=abs(fft(wholeSignal));
figure (2)
plot(fourier);
xlabel('Samples');
ylabel('Magnitude');
title('Frequency-Domain plot - FFT');
N=length(wholeSignal); %length of signal
Freq=((0:N-1)/N)*fs; %frequency axis for FFT
Amp=abs(wholeSignal)/(N/2); %amplitude axis for FFT
figure (3)
plot(Freq,fourier);
%stem(Freq,fourier);
xlabel('Frequency in Hz');
ylabel('Amplitude');
title('FFT')
0 Comments
Accepted Answer
Star Strider
on 22 Feb 2024
The problem was that ‘wholeSignal’ had 60 more elements than ‘t’ so I artificially shortened it:
wholeSignal = wholeSignal(1:numel(t));
so that they both fit.
After that, your code works. I would suggest doing a one-sided fft rather than two-sided, since it is a bit easier (more intuitive) to interpret. See the fft documentation for help with that.
Windowing the fft may also help, and you can do that with:
fourier=abs(fft(wholeSignal.*hann(numel(wholeSignal))));
That produces a slightly cleaner result. I added it (commented-out) so you can experiment with it.
fs=100000; %sampling freq in Hz
pulsesPerS=60; %frequency in Hz
signalDuration=3; % in seconds
amplitude=22*10^-3; % in mA
pulseWidth=250*10^-6; % in µs
t=0:1/fs:signalDuration-1/fs; %vector over the duration of time and sampled every 1/fs
onePeriod=amplitude*[ones(round(pulseWidth*fs),1); -ones(round(pulseWidth*fs),1); zeros(round(1/pulsesPerS*fs)-2*round(pulseWidth*fs),1)];
wholeSignal=repmat(onePeriod,[pulsesPerS*signalDuration 1]);
wholeSignal = wholeSignal(1:numel(t));
figure (1)
plot(t,wholeSignal);
xlabel('Time (seconds)');
ylabel('Amplitude');
title('Biphasic Pulse Train')
fourier=abs(fft(wholeSignal));
% fourier=abs(fft(wholeSignal.*hann(numel(wholeSignal))));
figure (2)
plot(fourier);
xlabel('Samples');
ylabel('Magnitude');
title('Frequency-Domain plot - FFT');
N=length(wholeSignal); %length of signal
Freq=((0:N-1)/N)*fs; %frequency axis for FFT
Amp=abs(wholeSignal)/(N/2); %amplitude axis for FFT
figure (3)
plot(Freq,fourier);
%stem(Freq,fourier);
xlabel('Frequency in Hz');
ylabel('Amplitude');
title('FFT')
Those were the only changes I made to your code.
Good luck with your FES research!
.
0 Comments
More Answers (0)
See Also
Categories
Find more on Fourier Analysis and Filtering 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!