How to extract data from txt file and plot spectogram?

16 次查看（过去 30 天）
Rahul
Rahul 2023-6-3
评论： Star Strider 2023-6-15
Hello,
I have a large sample of data (100000 samples) in txt file consisting of time and amplitude of the signal.
I want to plot spectogram. The code is not at all able to extract the .txt file at all along the path.
Kindly advice me as how to extract data from txt file when data is very large as it is in my case and hence to plot the spectogram.
The code is as below:
close all;
freq_sam = 100000; % samples in tt1 data
period_tt1 = 0.005; % time period
freq_tt1 = 1/period_tt1; % frequency
data = importdata('D:\Data_tokamak\data\obp.txt');
x_data = data(:,1);
y_data = data(:,2);
disp(x_data);
disp(y_data);
signal_tt1 = data;
fft_tt1 = fft(signal_tt1); % fourier transform
fft_tt1 = fftshift(fft_tt1); % performing fft shift
f = freq_sam/2*linspace(-1,1,freq_sam);
figure;
plot(f, fft_tt1);

请先登录，再回答此问题。

采纳的回答

Star Strider
Star Strider 2023-6-3
I prefer the pspectrum function for these analyses. See specifically Spectrogram and Reassigned Spectrogram of Chirp for an example of returning the data and plotting it as well.
data = readmatrix('D:\Data_tokamak\data\obp.txt');
x_data = data(:,1);
y_data = data(:,2);
Ts = mean(diff(x_data))
Fs = 1/Ts;
figure
pspectrum(y_data, Fs, 'spectrogram') % Plot 'pspectrum' 'spectrogram'
[sp,fp,tp] = pspectrum(y_data,fs,"spectrogram"); % Return Values, Then Plot
figure
waterfall(fp,tp,sp')
set(gca,XDir="reverse",View=[60 60])
ylabel("Time (s)")
xlabel("Frequency (Hz)")
Be certain that the ‘x_data’ are sampled uniformly. The easiest way to do that is to calculate the standard deviation (std) of ‘diff(x_data)’ and compare it to ‘Ts’. It should be vanishingly small, and preferably less than 1E-6 of the ‘Ts’ value.
It would help to have ‘obp.txt’ to work with.
.
  37 个评论
显示 35更早的评论
Rahul
Rahul 2023-6-15
Hi,
The value of k in my code is never initialized to zero, its going from 1 to numfiles(=12).
I did as you adviced. But none of the plots are getting executed.
Even the spectrogram is also not getting executed.
kindly advice.
with regards.
% Phase (in radians) versus Frequency (in Hz) plot for
% multiple input data from txt files
numfiles = 12;
x_datac{:,1}=[];x_datac{:,2}=[];x_datac{:,3}=[];x_datac{:,4}=[];x_datac{:,5}=[];x_datac{:,6}=[];x_datac{:,7}=[];
x_datac{:,8}=[];x_datac{:,9}=[];x_datac{:,10}=[];x_datac{:,11}=[];x_datac{:,12}=[];
y_datac{:,1}=[];y_datac{:,2}=[];y_datac{:,3}=[];y_datac{:,4}=[];y_datac{:,5}=[];y_datac{:,6}=[];y_datac{:,7}=[];
y_datac{:,8}=[];y_datac{:,9}=[];y_datac{:,10}=[];y_datac{:,11}=[];y_datac{:,12}=[];
mydata = cell(1, numfiles);
f = @(filename)fullfile('D:\Data_tokamak\TT1 data',filename);
for k = 1:numfiles
myfilename = sprintf('OBP%dN.txt', k);
mydata{k} = readmatrix(f(myfilename), 'HeaderLines',8);
mydata{k}(:,1) = mydata{k}(:,1) * 1E-3;
end
for k=1:numfiles
nonfinite_rows{k} = nnz(~any(isfinite(mydata{k}),2)); % Check To See If Any Rows Have Non-Finite Values
if nonfinite_rows{k} > 0
data(~isfinite(mydata{k})) = NaN;
data = fillmissing(mydata{k}, 'linear');
end
x_data(:,k) = mydata{k}(:,1);
y_data(:,k) = mydata{k}(:,2);
Ts = mean(diff(x_data(:,k)));
Tsd = std(diff(x_data(:,k)));
Fs = 1/Ts;
[y_datac{:,k},x_datac{:,k}] = resample(y_data(:,k), x_data(:,k), round(Fs)); % Resample To Constant Sampling Frequency and using cell arrays
end
Fs = round(Fs);
Fn = Fs/2;
figure(1)
for k=1:numfiles
pspectrum(y_datac{:,k}, Fs, 'spectrogram'); % Plot 'pspectrum' 'spectrogram'
colormap(turbo); % Introduced: R2020b
[sp,fp,tp] = pspectrum(y_datac{:,k},Fs,"spectrogram"); % Return Values, Then Plot
end
wp = [39.8 40.2]/Fn; % Passband frequency Normalized
ws = [0.98 1.02].*wp; % Stopband frequency Normalized
Rp = 1; % Passband Ripple
Rs = 60; % Passband Ripple (Attenuation)
[n,wp] = ellipord(wp,ws,Rp,Rs); % Elliptic order calculation
[z,p,k] = ellip(n,Rp,Rs,wp); % Elliptic filter design: Zero-Pole-Gain
[sos,g] = zp2sos(z,p,k); % Second order selection for stability
figure(2)
freqz(sos, 2^18,Fs) ; % Filter Bode Plot
set(subplot(2,1,1), 'XLim',[0 100]) % OPTTIONAL
set(subplot(2,1,2), 'XLim',[0 100]) % OPTTIONAL
[h w] = freqz(sos, 10,Fs) ;
%Implementation of filtfilt function without plotting
k=12;
y_data_filt = filtfilt(sos, g, y_datac{:,k}); % Filter Signal
% Plot for phase (in radians) versus frequency (in hz)
for k=1:numfiles
ax = gca;
ax.YLim=[0,100];
ax.XTick = [0:10:100];
plot(w, pi/180*(angle(h)));
grid on;
xlabel('Frequency (Hz)');
ylabel('Phase (rad)');
hold on
end
Star Strider
Star Strider 2023-6-15
Are you certain this is correct:
myfilename = sprintf('OBP%dN.txt', k);
MATLAB is case-sensitive, and the file names must match exactly, even with respect to case, and the other file names were all lower-case. If that is the problem, the code should throw an error about not being able to find the file. Are there any errors?
So the first thing to do is to check to see what ‘mydata’ contains. Be certain that it is not empty.

请先登录，再进行评论。

更多回答（1 个）

Diwakar Diwakar
Diwakar Diwakar 2023-6-3
Try this code. May be this code will help you.
% Specify the path to your text file
file_path = 'path/to/your/file.txt';
% Read the data from the text file
data = dlmread(file_path);
% Extract the time and amplitude columns
time = data(:, 1);
amplitude = data(:, 2);
% Set the parameters for the spectrogram
window_size = 256; % Size of the window for the spectrogram
overlap = 128; % Overlap between consecutive windows
% Compute the spectrogram
[~, frequency, time_spectrogram, power] = spectrogram(amplitude, window_size, overlap, [], 'onesided', 'yaxis');
% Plot the spectrogram
figure;
imagesc(time_spectrogram, frequency, 10*log10(power));
axis xy; % Flip the y-axis direction
colormap jet; % Choose a colormap
colorbar; % Add a colorbar for reference
xlabel('Time');
ylabel('Frequency');
title('Spectrogram');
  4 个评论
显示 2更早的评论
Rahul
Rahul 2023-6-4
Its not working
it gives error:
Error using pwelch
Expected x to be finite.
Plz advice.
rgds
Star Strider
Star Strider 2023-6-4
@Rahul — We already covered this and I told you how to fix it.
If we are going to solve it, you need to attach your data.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Descriptive Statistics 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by