Query regarding calculating frequency and amplitude
이전 댓글 표시
Hello everyone, I am loading a text file which has has a big data (comprising of 2 columns and multiple rows). Data represents (Time and Electrical activity). So far, I have used this syntax to plot the text file :
[fid,msg] = fopen('sample.txt','rt');
assert(fid>=3,msg)
C = textscan(fid, '%f%f', 'CommentStyle','#', 'CollectOutput',true);
fclose(fid);
M = C{1};
plot (M)
I am adding a screenshot of the plot as well as the text file. My next goal is to calculate the frequency and amplitude for a period of 5 seconds before t= 376.6 second (time has a comment in the text file.) Can someone please help me out with this. #Thanks in advance.
답변 (1개)
Eduard Reitmann
2018년 8월 3일
Hope this helps.
%%Create sample data with 5 Hz sinosoid (Do not include this section in
% your code)
n = 1001;
t = linspace(360,380,n)';
M = [t sin(5*2*pi*t)];
%%Reallocate data
t = M(:,1);
X = M(:,2);
% Trim data (5 seconds before 376.6s)
i5 = t >= 376.6-5 & t <= 376.6;
t = t(i5);
X = X(i5);
Fs = 1/(t(2)-t(1));
%%Calculate Fast Fourier Transform
n = numel(X);
Fn=Fs/2; % Nyquist frequency
Y=fft(X,n); % Calculate fft
fftUnique=Y(1:floor(n/2)); % Find unique values (Symetry)
% Magnitude
fftScale=fftUnique/n; % Scale fft
fftComp=fftScale*2; % Compensate for unique values
mag=abs(fftComp); % Calculate magnitude
% Phase (optional)
phaUnique=angle(fftUnique); % Calculate phase
pha=unwrap(phaUnique); % Correct phase angles to produce
% smoother phase plots
% Frequency
f=(Fn*linspace(0,1,numel(mag)))'; % Output frequency
%%Plot
figure;
subplot(2,1,1)
plot(t,X)
subplot(2,1,2)
plot(f,mag)
댓글 수: 19
Akshat Shrivastava
2018년 8월 6일
Eduard Reitmann
2018년 8월 6일
Put this section after you have imported Emma:
t = Emma.VarName1;
X = Emma.VarName2;
indexOfComments = find(~isnan(Emma.VarName3),1,'first');
t_end = t(indexOfComments);
indexTrim = t >= t_end-5 & t <= t_end;
t = t(indexTrim);
X = X(indexTrim);
Fs = 1/(t(2)-t(1));
Akshat Shrivastava
2018년 8월 6일
Eduard Reitmann
2018년 8월 6일
Can you kindly send my a screenshot of?
t = Emma.VarName1;
X = Emma.VarName2;
plot(t,X);
I suspect the data is incorrectly imported. Also, please send me your script from top to bottom. Thanks.
Akshat Shrivastava
2018년 8월 6일
Eduard Reitmann
2018년 8월 6일
The
indexTrim = t >= t_end-5 & t <= t_end;
line should be changed to:
indexTrim = t >= t_end - 5 & t <= t_end + 5;
It basically, finds the indices of all the values larger than t_end-5 and smaller than t_end+5. Thus, a 10 s band around t_end.
With regards to the plots that you sent, the top graph of "plot(f,mag).PNG" does not make sense because the x-axis is from 0 to 5. What is your t_end value? Also, what are the dimensions ("size(t)" and "size(X)") of your t and X variables? There are multiple lines in the top graphs of "plot(f,mag).PNG" which also doesn't make sense.
The "plot(t,X).PNG" graph seems believable.
Can you attach the data file?
Akshat Shrivastava
2018년 8월 7일
Eduard Reitmann
2018년 8월 7일
One problem with the input data was that the time vector restarted at zero a couple of times.

Also, I did not realize that there were multiple events. The updated script (attached) should do the trick.
Akshat Shrivastava
2018년 8월 7일
Akshat Shrivastava
2018년 8월 7일
Akshat Shrivastava
2018년 8월 7일
Akshat Shrivastava
2018년 8월 8일
Eduard Reitmann
2018년 8월 8일
The new 'sample.txt' file to 111.225 still restarted at some point.
Use the attached script with your original data file. The time vector is not important if you are only interested in frequency and the samples are evenly spaced. You only need the sample frequency (Fs). This updated script takes +-5s (10s) band around all the comment positions (assuming it is one continuous data file) and averages all the frequency responses.
Eduard Reitmann
2018년 8월 8일
FYI the positions of the comments might need to be scrutinized. They overlap sometimes and are also at positions where the file seems to have erroneous data points.
Akshat Shrivastava
2018년 8월 8일
Eduard Reitmann
2018년 8월 8일
I would not recommend doing this automatically - it can be done but the effort outweighs the reward. Rather manually define all the sections using the attached script. Just populate "tpos" with all the sections you want, according to the newly defined time vector. Use first plot as reference.
Akshat Shrivastava
2018년 8월 8일
Akshat Shrivastava
2018년 8월 9일
Akshat Shrivastava
2018년 8월 10일
카테고리
도움말 센터 및 File Exchange에서 Bartlett에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!