getToneFro​mPSD函数中这一段​代码的含义是什么?

11 views (last 30 days)
dayan
dayan on 3 Sep 2024
Answered: Epsilon on 10 Sep 2024
我在研究MATLAB计算sinad的函数sinad(),其中有一段调用了一个内部函数getToneFromPSD(Pxx,F,rbw,toneFreq),在该函数的最后一段,有“protect against nearby tone invading the window kernel”这一段,请问这一段代码的作用是什么呢?之前有通过psdbandpower函数计算得到一个power值,这里为何要对power比较一下并重新赋值呢?
% protect against nearby tone invading the window kernel
if nargin>2 && power < rbw(1)*colPxx(idxToneScalar)
power = rbw(1)*colPxx(idxToneScalar);
freq = colF(idxToneScalar);
end
下面附上getToneFromPSD函数全部代码:
function [power, freq, idxTone, idxLeft, idxRight] = getToneFromPSD(Pxx, F, rbw, toneFreq)
%GETTONEFROMPSD Retrieve the power and frequency of a windowed sinusoid
%
% This function is for internal use only and may be removed in a future
% release of MATLAB
% Copyright 2013-2019 The MathWorks, Inc.
%#codegen
idxTone = []; idxLeft = []; idxRight = [];
coder.varsize('idxTone','idxLeft','idxRight',[1,1],[1,1]);
% force column vector
colPxx = Pxx(:);
colF = F(:);
if nargin<4
[~, idxTone] = max(colPxx);
elseif colF(1) <= toneFreq(1) && toneFreq(1) <= colF(end)
% find closest bin to specified freq
[~, idxTone] = min(abs(colF-toneFreq(1)));
% look for local peak in vicinity of tone
iLeftBin = max(1,idxTone(1)-1);
iRightBin = min(idxTone(1)+1,numel(colPxx));
[~, idxMax] = max(colPxx(iLeftBin:iRightBin));
idxTone = iLeftBin+idxMax-1;
else
power = NaN('like',Pxx);
freq = NaN('like',Pxx);
idxTone = [];
idxLeft = [];
idxRight = [];
return
end
idxToneScalar = idxTone(1);
% sidelobes treated as noise
idxLeft = idxToneScalar - 1;
idxRight = idxToneScalar + 1;
% roll down slope to left
while idxLeft(1) > 0 && colPxx(idxLeft(1)) <= colPxx(idxLeft(1)+1)
idxLeft = idxLeft - 1;
end
% roll down slope to right
while idxRight(1) <= numel(colPxx) && colPxx(idxRight(1)-1) >= colPxx(idxRight(1))
idxRight = idxRight + 1;
end
% provide indices to the tone border (inclusive)
idxLeft = idxLeft+1;
idxRight = idxRight-1;
idxLeftScalar = idxLeft(1);
idxRightScalar = idxRight(1);
% compute the central moment in the neighborhood of the peak
Ffund = colF(idxLeftScalar:idxRightScalar);
Sfund = colPxx(idxLeftScalar:idxRightScalar);
freq = dot(Ffund, Sfund) ./ sum(Sfund);
% report back the integrated power in this band
if idxLeftScalar<idxRightScalar
% more than one bin
power = bandpower(colPxx(idxLeftScalar:idxRightScalar),colF(idxLeftScalar:idxRightScalar),'psd');
elseif 1 < idxRightScalar && idxRightScalar < numel(colPxx)
% otherwise just use the current bin
power = colPxx(idxRightScalar) * (colF(idxRightScalar+1) - colF(idxRightScalar-1))/2;
else
% otherwise just use the average bin width
power = colPxx(idxRightScalar) * mean(diff(colF));
end
% protect against nearby tone invading the window kernel
if nargin>2 && power < rbw(1)*colPxx(idxToneScalar)
power = rbw(1)*colPxx(idxToneScalar);
freq = colF(idxToneScalar);
end
  2 Comments
Walter Roberson
Walter Roberson on 3 Sep 2024
Approximate translation:
I am studying the MATLAB function sinad() that calculates sinad. There is a section that calls an internal function getToneFromPSD(Pxx,F,rbw,toneFreq). At the end of the function, there is a section called "protect against nearby tone invading the window kernel". What is the purpose of this code? A power value was calculated by the psdbandpower function before. Why do we need to compare the power and reassign it here?
dayan
dayan on 4 Sep 2024
thank you for the translation

Sign in to comment.

Answers (1)

Epsilon
Epsilon on 10 Sep 2024
Hi Dayan,
The internal function getToneFromPSD(Pxx,F,rbw,toneFreq) is used to get the power and frequency of the identified frequency along with its index and the index of the boundaries, from the passed power spectral density estimate and the resolution bandwidth. This would be the then removed from the PSD to isolate the noise and distortion components.
The code in the % protect against nearby tone invading the window kernel” section seems to be there to ensure the integrated power in the bin <power> is not less than the baseline level defined by the product of resolution bandwidth <rbw> and the power at the peak frequency bin <colPxx(idxToneScalar)>. This can happen due to spectral leakages as the tone is not perfectly isolated by the boundaries<idxRight> & <idxLeft>, leading to lower power calculation by the <bandpower> function.
It then also reassigns the frequency to that of the peak frequency bin.
This section of the code is therefore there to maintain the baseline level and act as a safeguard, in the presence of spectral interference from a nearby strong tone.

Products


Release

R2023a

Community Treasure Hunt

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

Start Hunting!