Q Factor is not same

12 次查看(过去 30 天)
HI i have a matlab code but i am unable to find the mistake. i have attched the code.
Basically what i am unable to achieve is based on the Quality factor formula Q=f/f2-f1 we supposed to have the same quality factor in the input values alpha = 2*pi/Q. so, When, I change Amplitude or the frequency the Q factor goes far away from the right Q factor. In, unfortunately, I think I got some errors in the dB scale or something else. I must have Q factor in the output as the same Q factor in the input. If wouldn't mind I submit the Matlab code take look at and fix the problem. I cant figure out where is the error in the code that is messing up with this.
  • Input:
  • Digital Signal:
  • F= 100 Hz
  • Q= 500
  • alpha= (f*pi)/Q
  • A= 10
Result:
Q= fmax/f2-f1
Q= 237.5527
  • Digital Signal:
  • F= 100 Hz
  • Q= 500
  • A= 100
Result:
Q= 99.4541
Any help is appreciated. Thanks in advance.
clear all
clc
close all
fs= 10000
t = (0:1/fs:(1-(1/fs)));
f= 4000
alpha= (f*pi)/500
y = 10*exp(-alpha*t).*sin(2*pi*f*t);
figure
plot(t,y)
N1=10000
S =fft(y,N1);
F1=S(1:N1/2); %half of spectra
PF1=2*F1.*conj(F1)/(fs*N1); %Power spectrum density per 1 kHz
E=sum(PF1)%energy of signal V^2*T
Xdb = 20*log10(S);
LPF1=10*log10(PF1); %Power spectrum in dB scale
w=(1:N1/2);
LP=LPF1(1:N1/2);
w1=fs*w/N1;
figure
plot(w1,LP)
BB= LP
CC= w1
E1 = max(LP)
t1= max(t)
talph1= E / -t1
%alpha1= max(BB*20*log(e))
[max_Z, max_index]=max(BB) %maximum value of impedance
threedb=max_Z*sqrt(2)/2; %the 3db point
[db_Z,index_db] = min(abs(BB.'-threedb)-max_Z)
%db3idx = zci(pwr-hpp);
[pk,loc] = findpeaks(BB,'Npeaks',1,'SortStr','descend');
%db3c= threedb-max_Z
ofst= 12
for k1 = 1:length(pk)
varmtx=[(BB(loc(k1)-ofst:loc(k1)));CC(loc(k1)-ofst:loc(k1));(BB(loc(k1):loc(k1)+ofst));CC(loc(k1):loc(k1)+ofst)]
dBpts(k1,1:2) = interp1(varmtx(1,:), varmtx(2,:), -(db_Z), 'linear','extrap');
dBpts(k1,3:4) = interp1(varmtx(3,:),varmtx(4,:), -(db_Z), 'linear','extrap');
end
figure(3)
plot(CC, BB)
hold on
for k1 = 1:length(pk)
plot(dBpts(k1,:), -(db_Z), 'r+')
end
hold off
grid
F_1 = dBpts(k1,1:2)
F_1 = F_1(1)
F_2 = dBpts(k1,3:4)
F_2 =F_2(1)
BW = (F_2)-(F_1)
Q = abs((loc-1)*(1/BW))

采纳的回答

David Goodmanson
David Goodmanson 2021-2-26
Hello Abdullah
be aware that what you are calling alpha is often called alpha/2, but I don't think that had anything to do with the problem you are having. I stayed with your definition of alpha.
Everything appears to be fine up until the width calculation. I don't know whether you intended that w denote circular frequency, but in this case it is (regular) frequency. To avoid ambiguity I went with frequency and did not use circular frequency. The code below finds the widths and the resulting Q.
clc
close all
fs= 10000;
t = (0:1/fs:(1-(1/fs)));
f= 4000;
alpha= (f*pi)/500;
y = 10*exp(-alpha*t).*sin(2*pi*f*t);
figure(1)
plot(t,y); grid on
Q0 = pi*f/alpha
N1=10000;
S =fft(y,N1);
F1=S(1:N1/2); %half of spectra
PF1=2*F1.*conj(F1)/(fs*N1); %Power spectrum density per 1 kHz
% E=sum(PF1);%energy of signal V^2*T
% Xdb = 20*log10(S);
LPF1=10*log10(PF1); %Power spectrum in dB scale
LP=LPF1(1:N1/2);
fvec=fs*(1:N1/2)/N1; % frequency array
figure(2)
plot(fvec,LP); grid on
xlim([3990 4010])
% find half power points by interpolation of LP
[LPmax ind0] = max(LP);
LPhalf = LPmax - 10*log10(2); % half power points
indx = find(LP>LPhalf);
ind1 = [indx(1)-5:ind0]; % use a few more index points, 5 is somewhat arbitrary
f1 = interp1(LP(ind1),fvec(ind1),LPhalf,'spline');
ind2 = [ind0:indx(end)+5];
f2 = interp1(LP(ind2),fvec(ind2),LPhalf,'spline');
Q = f/(f2-f1)
Q0 = 500.0000
Q = 499.9917

更多回答(0 个)

类别

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

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by