# thd
总谐波失真
函数库: TySignalProcessing
# 语法
r, = thd(x)
r, = thd(x,fs,n)
r, = thd(pxx,f,"psd")
r, = thd(pxx,f,n,"psd")
r, = thd(sxx,f,rbw,"power")
r, = thd(sxx,f,rbw,n,"power")
r, = thd(___,harmType = "aliased")
r,harmpow,harmfreq = thd(___)
r,harmpow,harmfreq = thd(___;plotfig = Bool)
# 说明
r, = thd(x) 返回实值正弦信号 x 的总谐波失真（thd），单位为 dBc。总谐波失真使用与输入信号长度相同的修改周期图从基频和前五次谐波确定。修改后的周期图使用 β=38 的 Kaiser窗。 示例
r, = thd(x, fs,n) 指定采样率 fs 和用于 thd 计算的谐波数（包括基波）。 示例
r, = thd(pxx,f,"psd") 将输入 pxx 指定为单边功率谱密度（psd）估计值。f 是与 pxx 中的 PSD 估计相对应的频率向量。
r, = thd(pxx,f,n,"psd") 指定用于 thd 计算的谐波数（包括基波）。 示例
r, = thd(sxx,f,rbw,"power") 将输入指定为单边功率谱。rbw 是对每个功率估计进行积分的分辨率带宽。
r, = thd(sxx,f,rbw,n,"power") 指定在总谐波失真（THD）计算中使用的谐波数量（包括基波）。
r, = thd(___,harmType = "aliased") 报告混叠到奈奎斯特范围内的基波谐波。当输入信号采样不足时，使用此选项。如果未指定此选项，或者将其设置为 "omitaliases"，则函数将忽略超出奈奎斯特范围的基频谐波。
r,harmpow,harmfreq = thd(___) 返回谐波的功率（单位：dB）和频率，包括基波。 示例
r,harmpow,harmfreq = thd(___;plotfig = Bool) 绘制信号的频谱，在当前图形窗口中注释谐波。 示例
# 示例
确定具有两个谐波的信号的总谐波失真
该示例明确显示了如何计算由基波和二次谐波组成的信号的总谐波失真（dBc）。根据 thd 返回的结果检查显式计算。
创建以 1 kHz 采样的信号。信号由振幅为 2 的 100Hz 基波和振幅为 0.01 和 0.005 的 200 和 300Hz 的两个谐波组成。使用 thd 明确获得总谐波失真。
using TySignalProcessing
t = [0:0.001:(1 - 0.001);]
x = @. 2 * cos(2 * pi * 100 * t) +
0.01 * cos(2 * pi * 200 * t) +
0.005 * cos(2 * pi * 300 * t)
tharmdist = 10 * log10((0.01^2 + 0.005^2) / 2^2)
tharmdist = -45.05149978319906
r, = thd(x)
r = -45.051499783200384
指定谐波数
创建以 1 kHz 采样的信号。该信号由振幅为 2 的 100 Hz 基波和振幅为 0.01、0.005 和 0.0025 的 200、300 和 400 Hz 的三个谐波组成。
将谐波数设置为 3。这包括基波。因此，在 thd 计算中使用 100 Hz、200 Hz 和 300 Hz 的功率。
using TySignalProcessing
t = [0:0.001:(1 - 0.001);]
x = @. 2 * cos(2 * pi * 100 * t) +
0.01 * cos(2 * pi * 200 * t) +
0.005 * cos(2 * pi * 300 * t) +
0.0025 * sin(2 * pi * 400 * t)
r, = thd(x, 1000, 3)
r = -45.05149978320038
指定谐波数等于 3 将忽略 thd 计算中 400 Hz 的功率。
指定谐波数（PSD 输入）
创建以 1 kHz 采样的信号。该信号由振幅为 2 的 100 Hz 基波和振幅为 0.01、0.005 和 0.0025 的 200、300 和 400 Hz 的三个谐波组成。
获得信号的周期图 PSD 估计值，并使用 PSD 估计值作为 thd 的输入。将谐波数设置为 3。这包括基波。因此，在 thd 计算中使用 100 Hz、200 Hz 和 300 Hz 的功率。
using TySignalProcessing
t = [0:0.001:(1 - 0.001);]
fs = 1000;
x = @. 2 * cos(2 * pi * 100 * t) +
0.01 * cos(2 * pi * 200 * t) +
0.005 * cos(2 * pi * 300 * t) +
0.0025 * sin(2 * pi * 400 * t)
pxx, f = periodogram(x, rectwin(length(x)), length(x), fs; nargout = 2);
r, = thd(pxx, f, 3, "psd")
r = -45.05149978320068
功率谱的总谐波失真
通过输入用汉明窗口获得的功率谱和窗口的分辨率带宽来确定 thd。
创建以 10 kHz 采样的信号。信号由振幅为 2 的 100 Hz 基波和振幅为 0.01、0.005 和 0.0025 的 300、500 和 700 Hz 的三个奇数谐波组成。指定谐波数为 7。确定 thd。
using TySignalProcessing
fs = 10000
t = [0:(1 / fs):(1 - 1 / fs);]
x = @. 2 * cos(2 * pi * 100 * t) +
0.01 * cos(2 * pi * 300 * t) +
0.005 * cos(2 * pi * 500 * t) +
0.0025 * sin(2 * pi * 700 * t)
sxx, f = periodogram(
x, hamming(length(x)), length(x), fs, "power"; nargout = 2
);
rbw = enbw(hamming(length(x)), fs)
r, = thd(sxx, f, rbw, 7, "power")
r = -44.83960278600343
带和不带混叠谐波的总谐波失真
生成一个信号，该信号类似于具有 2.1 kHz 音调作为输入的弱非线性放大器的输出。信号在 10 kHz 下采样 1 s 。
using TySignalProcessing
using TyPlot
using TyMath
Fs = 10000
f = 2100
rng = MT19937ar(1234)
t = [0:(1 / Fs):1;]
x = tanh.(sin.(2 .* pi .* f .* t) .+ 0.1) .+ 0.001 .* randn(rng, length(t))
periodogram(x,kaiser(length(x),38),[],Fs; plotfig = true)
计算信号的总谐波失真。默认情况下，thd 将混叠谐波视为噪声的一部分。
thd(x, Fs, 7; plotfig=true)
重复计算，但现在将混叠谐波视为信号的一部分。
figure()
thd(x, Fs, 7; plotfig=true, harmType="aliased")
谐波功率和相应频率
创建以 10 kHz 采样的信号。信号由振幅为 2 的 100 Hz 基波和振幅为 0.01、0.005 和 0.0025 的 300、500 和 700 Hz 的三个奇数谐波组成。指定谐波数为 7。确定 thd、谐波功率和相应频率。
using TySignalProcessing
fs = 10000
t = [0:(1/fs):(1-1/fs);]
x = @. 2 * cos(2 * pi * 100 * t) +
0.01 * cos(2 * pi * 300 * t) +
0.005 * cos(2 * pi * 500 * t) +
0.0025 * sin(2 * pi * 700 * t) +
1e-10 * cos(2 * pi * 200 * t) +
1e-10 * cos(2 * pi * 400 * t) +
1e-10 * cos(2 * pi * 600 * t)
r, harmpow, harmfreq = thd(x, 10000, 7)
[harmfreq round.(harmpow, sigdigits=3)]
[harmfreq harmpow] = 7×2 Matrix{Float64}:
100.0 3.01
200.0 -203.0
300.0 -43.0
400.0 -203.0
500.0 -49.0
600.0 -203.0
700.0 -55.1
偶数次谐波的功率约为 −200 dB，对应于
放大信号的总谐波失真
生成以 50 kHz 采样的频率为 2.5 kHz 的正弦曲线。将标准偏差为 0.00005 的高斯白噪声添加到信号中。将结果通过弱非线性放大器。绘制 thd。
using TySignalProcessing
using TyMath
fs = 5e4
f0 = 2.5e3
N = 1024
t = [(0:(N - 1)) / fs;]
ct = @. cos(2 * pi * f0 * t)
rng = MT19937ar(1234);
cd = ct + 0.00005 * randn(rng, size(ct))
amp = [1e-5, 5e-6, -1e-3, 6e-5, 1, 25e-3]
sgn = polyval(amp, cd)
thd(sgn, fs; plotfig=true)
该图显示了用于计算比率的频谱以及被视为噪声的区域。计算中不包括直流电平。标注基波和谐波。
# 输入参数
x - 实正弦双音信号向量
实正弦双音信号，指定为向量。
数据类型： Float64
fs - 采样频率1（默认值）| 正实标量
采样率，指定为正实标量。采样率是每单位时间的采样数。如果时间单位为 s，则采样率的单位为 Hz。
数据类型： Float64
n - 谐波数量正整数
谐波数，指定为正整数。
数据类型： Int
pxx - 单边 PSD 估计向量
单边功率谱密度估计，指定为实值、非负向量。
功率谱密度必须用线性单位表示，而不是分贝。使用 db2pow 将分贝值转换为功率值。
数据类型： Float
f - 循环频率向量
与指定为向量的单侧 PSD 估计值 pxx 相对应的周期频率。f 的第一个元素必须为 0。
数据类型： Float
sxx - 功率谱非负实向量
功率谱，指定为实值非负向量。
功率谱必须用线性单位表示，而不是分贝。使用 db2pow 将分贝值转换为功率值。
数据类型： Float
rbw - 分辨率带宽正标量
分辨率带宽，指定为正标量。分辨率带宽是离散傅里叶变换的频率分辨率和窗口的等效噪声带宽的乘积。
数据类型： Float
# 输出参数
r - 总谐波失真（dBc）实标量
以 dBc 为单位的总谐波失真，作为实标量返回。
数据类型： Float
harmpow - 谐波功率实标量 | 向量
谐波的功率，返回为以 dB 表示的实标量或向量。harmpow 是标量还是向量取决于指定为输入参数 n 的谐波数。
数据类型： Float
harmfreq - 谐波频率非负标量 | 向量
谐波的频率，作为非负标量或或向量返回。harmfreq 是标量还是或向量取决于指定为输入参数 n 的谐波数。
数据类型： Float
# 更多
p - 平均频带功率
函数 thd、thd、thd 和 snr 测量正弦激励的弱非线性系统的响应。
当给定时域输入时，thd 使用具有大旁瓣衰减的 Kaiser 窗口执行周期图。为了找到基频，该算法在周期图中搜索最大的非零频谱分量。然后，它计算从最大值单调减少的所有相邻箱的中心矩。为了可检测，基波至少应在第二频率仓中。高次谐波是基频的整数倍。如果一个谐波位于另一个谐波附近的单调递减区域内，则其功率被认为属于较大的谐波。这个较大的谐波可能是，也可能不是基波。
如果基波不是信号中最高的频谱分量，则 thd 失败。
确保频率分量相距足够远，以适应 Kaiser 窗口的旁瓣宽度。如果这不可行，您可以使用 "power" 标志并使用不同的窗口计算周期图。