Time response from T parameter matrix
3 次查看(过去 30 天)
显示 更早的评论
Hello,I need to pass a signal(CSV file having time and voltage columns) through a T parameter matrix of 5000 rows of frequency.Can somebody let me know which function to use to process a signal through a T parameter matrix?
10 个评论
Mathieu NOE
2024-2-28
so your input signal looks like a low pass filtered square wave (frequency = 0.77 GHz + harmonics)
the signal is a bit short so I could only do a single buffer fft (no averaging) to have finest resolution (still limited to 0.2 GHz) , so if you pick two peaks from the time plot you can see the time delta (signal period) is around 1.3e-9 s , therefore my frequency = 0.77 GHz.
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1630391/image.png)
the first peak from the fft graph is matching this value (at 0.8 GHz +/- 0.2 GHz resolution)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1630396/image.png)
now if we want to keep the general shape of your signal , that means we need to take care of the upper harmonics as well.
We can see from the fft graph that above 10 GHz the spectrum amplitude drops rapidely , so we can forget about what is above this limit .
so back to our FRF models , wee can say that we don't look for perfect match above 10 GHz
FYI, this is my code to generate the time and fft plots
data=readmatrix('probe.xlsx');
time=data(:,1);
signal=data(:,2);
dt = mean(diff(time));
Fs = 1/dt;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% FFT parameters
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
NFFT = length(signal); % for maximal resolution take NFFT = signal length (but then no averaging possible)
OVERLAP = 0.75; % overlap will be used only if NFFT is shorter than signal length
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% display 1 : time domain plot
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
figure(1),
plot(time,signal,'b');grid on
title(['Time plot / Fs = ' num2str(Fs) ' Hz ']);
xlabel('Time (s)');ylabel('Amplitude');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% display 2 : averaged FFT spectrum
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[freq, spectrum_raw] = myfft_peak(signal,Fs,NFFT,OVERLAP);
ind = freq>0;
freq = freq(ind);
spectrum_raw = spectrum_raw(ind);
figure(2),semilogx(freq,20*log10(spectrum_raw),'b');grid on
df = freq(2)-freq(1); % frequency resolution
title(['Averaged FFT Spectrum / Fs = ' num2str(Fs) ' Hz / Delta f = ' num2str(df,3) ' Hz ']);
xlabel('Frequency (Hz)');ylabel('Amplitude (dB)');
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
Mathieu NOE
2024-2-28
编辑:Mathieu NOE
2024-2-28
now let's try to put a low order model on your FRF
here I opted for a very low order (NA = NB = 1) and IMHO it should suffice
now , as we have the B and A coefficients of the IIR filters , we will use them latter on to apply them on your input signal
NB that the invfreqz operation is done with the sampling rate of your input signal , in order afterwards to make the filtering (with filter) : ( Fs = 1.2800e+12 )
% model for FRF11 / FRF22
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1630636/image.png)
% model for FRF12 / FRF21
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1630641/image.png)
%% code if you have the following toolboxes
% Antenna Toolbox
% RF Toolbox
% % Load the .s2p file using the sparameters function
% s2pObject = sparameters('2Port.s2p');
%
% % Extract frequency and S21 parameter
% frequency = s2pObject.Frequencies;
% S21 = s2pObject.Parameters(:,2,1); % S21 parameter
% % Find the magnitude of S21
% S21_magnitude = abs(S21);
%% my home made solution
% data file is organized s follows :
% !! freq mag(S11) phase(S11) mag(S12) phase(S12)
% ! mag(S21) phase(S21) mag(S22) phase(S22)
% !
% # GHz S MA R 50
data = readmatrix('2Port.s2p','FileType','text');
%freq/mag(S11)/phase(S11)/mag(S12)/phase(S12)/mag(S21)/phase(S21)/mag(S22)/phase(S22)
freq = data(:,1)*1e9; % Hz
FRF11 = data(:,2).*(cosd(data(:,3))+1i*sind(data(:,3))); % complex FRF
FRF12 = data(:,4).*(cosd(data(:,5))+1i*sind(data(:,5))); % complex FRF
FRF21 = data(:,6).*(cosd(data(:,7))+1i*sind(data(:,7))); % complex FRF
FRF22 = data(:,8).*(cosd(data(:,9))+1i*sind(data(:,9))); % complex FRF
%% create model for FRF11 / FRF22
% IIR model
iter = 100;
Fs = 1.2800e+12; % sampling rate of your time signal probe
W = 2*pi*freq/Fs;
ind = freq <20e9; % frequency weighting ;
WT = zeros(size(W));
WT(ind) = 1;
NB = 1;
NA = 1;
[B,A] = invfreqz(FRF11,W,NB,NA,WT,iter);
% my trick to force both values of B numerator to have same value but opposite
% sign, in order to have a true differentiator at the numerator
Bmean = mean(abs(B));
B = Bmean*[1 -1];
% check bode plots
H = freqz(B,A,2*pi*freq/Fs);
figure(1)
subplot(2,1,1),semilogx(freq,20*log10(abs(FRF11)),freq,20*log10(abs(FRF22)),freq,20*log10(abs(H)),'k--');
xlabel('Frequency(Hz)');
ylabel('Modulus (dB)');
legend('FRF11','FRF22','identified model');
subplot(2,1,2),semilogx(freq,180/pi*(angle(FRF11)),freq,180/pi*(angle(FRF22)),freq,180/pi*(angle(H)),'k--');
xlabel('Frequency(Hz)');
ylabel('angle (deg)');
legend('FRF11','FRF22','identified model');
%% create model for FRF12 / FRF21
% IIR model
iter = 100;
Fs = 1.2800e+12; % sampling rate of your time signal probe
W = 2*pi*freq/Fs;
ind = freq <3e9; % frequency weighting ;
WT = zeros(size(W));
WT(ind) = 1;
NB = 1;
NA = 1;
[B,A] = invfreqz(FRF12,W,NB,NA,WT,iter);
% check bode plots
H = freqz(B,A,2*pi*freq/Fs);
figure(2)
subplot(2,1,1),semilogx(freq,20*log10(abs(FRF12)),freq,20*log10(abs(FRF21)),freq,20*log10(abs(H)),'k--');
xlabel('Frequency(Hz)');
ylabel('Modulus (dB)');
legend('FRF12','FRF21','identified model');
subplot(2,1,2),semilogx(freq,180/pi*(angle(FRF12)),freq,180/pi*(angle(FRF21)),freq,180/pi*(angle(H)),'k--');
xlabel('Frequency(Hz)');
ylabel('angle (deg)');
legend('FRF12','FRF21','identified model');
回答(0 个)
另请参阅
类别
在 Help Center 和 File Exchange 中查找有关 Spectral Analysis 的更多信息
产品
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!