frequency response of time series data

9 次查看(过去 30 天)
Hello veryone.
i have a time series data of 10 seconds attached and i would like to make Psd vs frequency plot.
below is my code: attached also is my result for data which is also attached a .txt file.
I would like to have it llok better and vsible as shown in the eample plot expected.png so i cant clearly see the natural frequency, and harmonics.
Many thanks in advance
function filecontent = importfile(filename, startRow, endRow)
%IMPORTFILE Import numeric data from a text file as a matrix.
% DCLAIB9CM5CM50MS0090DEG = IMPORTFILE(FILENAME) Reads data from text
% file FILENAME for the default selection.
%
% DCLAIB9CM5CM50MS0090DEG = IMPORTFILE(FILENAME, STARTROW, ENDROW) Reads
% data from rows STARTROW through ENDROW of text file FILENAME.
%
% Example:
% Dclaib9cm5cm50ms0090deg = importfile('2Dclaib_9cm_5cm_50ms_0090deg.txt', 1, 3000);
%
% See also TEXTSCAN.
% Auto-generated by MATLAB on 2019/08/08 17:30:42
%% Initialize variables.
delimiter = '\t';
if nargin<=2
startRow = 1;
endRow = inf;
end
%% Format for each line of text:
% column1: double (%f)
% column2: double (%f)
% column3: double (%f)
% column4: double (%f)
% column5: double (%f)
% column6: double (%f)
% column7: double (%f)
% column8: double (%f)
% column9: double (%f)
% column10: double (%f)
% column11: double (%f)
% For more information, see the TEXTSCAN documentation.
formatSpec = '%f%f%f%f%f%f%f%f%f%f%f%[^\n\r]';
%% Open the text file.
fileID = fopen(filename,'r');
%% Read columns of data according to the format.
% This call is based on the structure of the file used to generate this
% code. If an error occurs for a different file, try regenerating the code
% from the Import Tool.
dataArray = textscan(fileID, formatSpec, endRow(1)-startRow(1)+1, 'Delimiter', delimiter, 'TextType', 'string', 'EmptyValue', NaN, 'HeaderLines', startRow(1)-1, 'ReturnOnError', false, 'EndOfLine', '\r\n');
for block=2:length(startRow)
frewind(fileID);
dataArrayBlock = textscan(fileID, formatSpec, endRow(block)-startRow(block)+1, 'Delimiter', delimiter, 'TextType', 'string', 'EmptyValue', NaN, 'HeaderLines', startRow(block)-1, 'ReturnOnError', false, 'EndOfLine', '\r\n');
for col=1:length(dataArray)
dataArray{col} = [dataArray{col};dataArrayBlock{col}];
end
end
%% Close the text file.
fclose(fileID);
%% Post processing for unimportable data.
% No unimportable data rules were applied during the import, so no post
% processing code is included. To generate code which works for
% unimportable data, select unimportable cells in a file and regenerate the
% script.
%% Create output variable
filecontent = [dataArray{1:end-1}];
timestep = filecontent(:,1);
x=zeros(length(timestep),1);
y=zeros(length(timestep),1);
z=zeros(length(timestep),1);
for j =1:length(timestep)
x(j) =((filecontent(j,5)+filecontent(j,6))-(filecontent(j,4)+filecontent(j,7)))/sum(filecontent(j,4:7));
y(j) =((filecontent(j,6)+filecontent(j,7))-(filecontent(j,4)+filecontent(j,6)))/sum(filecontent(j,4:7));
z(j)= abs(x(j)+y(j)*1i);
end
TimserisData =[x y z];
Data =TimserisData(:,3); % Taking the magnitude of the calculate sigmnal
Fs =1000;%sampling frequecy
T=1/Fs; %sampling period
N =length(TimserisData(:,3)); %length of signal
t=(0:N-1)*T;%Time vector
NFFT=2^nextpow2(N);%length of fft
X=fft(Data,NFFT)/N;%Taking the FFT of the signal
f = Fs/2*linspace(0,1,NFFT/2+1); %frequency vector
%Plotting signal and axes labeling
figure
loglog(f,2*abs(X(1:NFFT/2+1)))
xlabel('Frequency (Hz)')
ylabel('|p(f)|[a.u]')

回答(1 个)

Jyotish Kumar
Jyotish Kumar 2019-12-30
Hi,
To calculate PSD, i would suggest you to try "pwelch" function in MATLAB.
Go through the link below for more detals

类别

Help CenterFile Exchange 中查找有关 Parametric Spectral Estimation 的更多信息

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by