虚拟模拟合成器：带限波形

在这篇文章中，我计划展示 polyBLEPS 的效果并实现一个简短的示例程序，该程序使用此方法生成带限波形，以期减少混叠。此外，还有一些其他可用方法的简短概述。这篇文章是在 Jupyter Notebook 中编写的。Notebook 可作为 .html 或 .ipnyb 使用。在减法合成中，我们希望从具有大量“谐波含量”的波形开始，然后我们可以使用滤波器和效果“雕刻掉”这些波形。这通常意味着锯齿波或方波。相比之下，正弦波的“谐波含量”很少。事实上，它的可能性最小，而且非常基本，可以将其视为任何其他信号的构建块。其中\(x\)，在音频的上下文中，是我们的时间参数，\(L\) 是方波周期，注意\(n\) 被定义为仅是奇数整数（推导）。 #Setup import numpy as np import scipy import scipy.signal import matplotlib.pyplot as plt import ipywidgets as widgets L = 1 / 10 ; # 1 周期 0.1s => 10Hz Fs = 1000 # 采样（稍后解释） t = np . linspace ( 0 , 1 , Fs ) # 1sec with 10000 samples signal = scipy .信号 。 square ( t * ( 2 * np . pi / L )) # 10Hz 方波 ps = np . abs ( np . fft . fft ( 信号 )) ** 2 # 功率谱频率 = np . fft 。 fftfreq (signal .size, d = 1 / Fs) # 绘图频率 #Plot idx = np . argsort ( freqs ) fig , ax = plt .子图 () ax 。轴。 set_ticks ( np . arange ( 1 / L , 100 , 1 / L )) ax . set_xlim ( 0 , 100 ) ax 。 set_title ( "%d Hz 方波的频谱" % ( 1 / L )) ax 。 set_xlabel ("Hz") ax 。 plot ( freqs [ idx ], np . sqrt ( ps [ idx ]), ".-" ) 上图显示了 10Hz 方波的频谱。由此我们可以看到构成方波的各个正弦波。最强的峰值是频率为 10Hz（通常称为基波）的正弦波。下一个峰值可以在 30Hz 处看到，这意味着我们在 30Hz 处有一个正弦波，但幅度较低。下一个峰值为 50Hz，下一个峰值为 70Hz，等等。请注意，峰值在 \(\frac{1}{L}\) (20Hz, 40Hz, ...) 的偶数倍处丢失，并且有在 \(\frac{1}{L}\) (10Hz, 30Hz, ...) 的不均匀倍数处达到峰值，正如 \eqref{fourierSeriesSquare} 中所描述的那样。

理论上，这些峰值一直出现到无限高的频率！但是，在数字世界中处理信号时，必须以特定的时间间隔对信号进行采样，因此在出现混叠之前，信号频率是有限制的。当我们试图捕获的信号的频率高于采样频率的一半时，就会发生混叠。当这种情况发生时，我们不能再完美地表示信号。它被扭曲了。奈奎斯特定理说，要通过采样完美地捕获信号，采样频率必须至少是信号中最高频率的两倍。为了说明混叠，我们可以查看以 40Hz 采样的正弦波。频率为 15Hz 的正弦波将在频谱中的 15Hz 处具有峰值。但是 25Hz 的正弦波高于采样频率一半的阈值，并且会出现混叠。我们可以看到，这意味着信号的频率“折叠”并似乎在频率上下降，落在 15Hz 左右而不是 25Hz，从而错误地表示信号。 def fft_spectrum ( t , signal , Fs ): """函数帮助绘制频率响应""" ps = np . abs ( np . fft . fft ( 信号 )) ** 2 # 功率谱频率 = np . fft 。 fftfreq (signal .size, d = 1 / Fs) # 绘图的频率 idx = np . argsort ( freqs ) return ( freqs [ idx ], ps [ idx ]) def sin_fft ( freq , Fs ): L = 1 / freq ; # 1 周期 0.1s => 10Hz Fs = Fs # 采样数 t = np 。 linspace ( 0 , 1 , Fs ) # 1sec with Fs samples signal = np 。 sin ( t * ( 2 * np . pi / L )) # 1/L Hz 正弦波返回 fft_spectrum ( t , signal , Fs ) #Plot！血小板子图 fig ， ax = plt 。子图 ( 2 , sharex = True ) ax [ 0 ] 。轴。 set_ticks ( np . arange ( 10 , 100 , 5 )) ax [ 0 ] 。 set_xlim (0, 20) ax [1]。 set_xlabel ("Hz") #15 Hz freqs , ps = sin_fft ( 15 , 40 ) ax [ 0 ] 。 plot ( freqs , np . sqrt ( ps )) ax [ 0 ] 。 set_title (“%d Hz 方波的频谱”% (15)) #25 Hz freqs , ps = sin_fft ( 25 , 40 ) ax [ 1 ] 。绘图 (freqs, np.sqrt (ps)) ax [1]。 set_title ( "%d Hz 方波的频率频谱" % ( 25 )) 由此我们可以想象，试图表示一个在无限高频率处具有峰值的方波可能会成为一个挑战。但是为什么原始方波示例中没有混叠？在那种情况下，与方波的周期相比，我们的采样率设置得相当高。较高频率的峰值非常小，以至于混叠的影响几乎可以忽略不计。但是，如果我们希望允许具有更短周期的方波（导致更高的基频），这可能会成为一个问题，并将由于混叠而产生的数字信号引入不需要的内容和噪声。

出于音频目的，通常使用 44.1kHz 的采样率。这意味着我们可以表示频率高达 20kHz 的音频，这超过了人耳所能接收的频率。但是存在混叠的风险，高于 20k kHz 的频率会“折回”并在可听范围内引入新的不需要的频率。如何避免混叠？ 信号需要进行带宽限制。这意味着要表示的信号不需要具有高于采样率一半的频率成分（或至少这种成分具有低幅度）。通过麦克风录制真实世界的音频时，几乎肯定会遇到高于采样率的频率。为了避免混叠，使用低通滤波器来确保较高的频率足够低，从而限制信号频带，从而避免混叠。当以数字方式生成波形时，我们需要确保也是如此。限制无法采样的较高频率。有多种方法可以帮助防止波形混叠。常用方法包括 这 3 种技术被称为准带限，因为它们不能完美地消除混叠。然而，它们仍然有助于显着减少混叠效应。下面是对过采样、加法合成技术和 BLEP 技术的快速概述。 BLEP 技术是我的项目最感兴趣的技术之一，似乎最适合我的需求。我没有机会深入研究其他技术。以下是有关其他 BLITS 和 BLAMPS 的更多信息的一些提示：

在这里，信号只是以足够高的频率采样，不会出现混叠或混叠的影响可以忽略不计。音频通常以 44.1kHz 采样，如果我们希望捕获具有高成分频率的信号，我们必须确保采样率是最高成分频率频率的两倍。但这可能意味着计算机必须完成更多工作，这可能是一个问题，尤其是对于实时应用程序。使用上面方波的傅立叶级数表达式（方程 \eqref{fourierSeriesSquare}），我们可以通过手动执行组成频率的求和并在到达奈奎斯特频率时简单地停止来生成信号的每个样本。这将使我们尽可能接近表示方波，而不会出现由混叠引起的错误探测频率。 def added_square ( N , freq = 10 , Fs = 1000 ): L = 1 / freq ; # 周期 Fs = Fs # # 采样数 t = np . linspace ( 0 , 1 , Fs ) # 1sec with Fs samples signal = np 。数组 ([ 0 ] * len ( t ), float ) for i in range ( 1 , len ( t )): t0 = t [ i ] for n in range ( 1 , N + 1 ): 信号 [ i ] += 4/np。 pi * \ 1 / ( 2 * n - 1 ) * np . sin ( 2 * ( 2 * n - 1 ) * np . pi * t0 / L ) return ( t , signal ) def plot_additive_square ( N , freq , Fs ): L = 1 / freq ;打印（“最高组成频率：%g Hz”\%（（（2*N-1）*np.pi/L）/（np.pi）））t，信号=additive_square（N，freq，Fs）plt .图 ( figsize = ( 14 , 3 )) plt .子图 (1, 2, 1) plt。 plot ( t , 信号 , ".-" ) plt.标题（“%d Hz 处的附加方波 \n \ 具有 %d 项，在 %d Hz 处采样” % ( freq , N , Fs )) plt 。 xlabel ("time [s]") freqs, ps = fft_spectrum (t, signal, Fs) plt.子图 (1, 2, 2) plt。绘图（频率，np。sqrt（ps））plt。标题（“方波频谱在 %d Hz \n \ 与 %d 项，在 %d Hz 处采样” % ( freq , N , Fs )) plt . xlabel ( "频率 [Hz]" ) plt . xlim(0, 120) 小部件。交互 ( plot_additive_square , \ N = widgets . IntSlider ( min = 1 , max = 7 , step = 1 , value = 3 ),\ freq = widgets . fixed ( 18 ), Fs = widgets . fixed ( 200 )) > 最高成分频率：90 Hz 在上图中，我们可以看到在 200Hz 下采样的 N 值不同的附加方波。我们还有这个波的频谱。当我们增加 N 时，峰值会以均匀的间隔出现。但是一旦我们达到采样频率 (100Hz) 的一半，不规则的峰值就会开始出现，这会导致非谐波声音（不太适合制作音乐）。所以我们可以在穿过奈奎斯特频率之前简单地停在 N 的值处，在这种情况下 \(N = 3\)。当以 44.1kHz 采样时，我们希望生成基频为 420kHz 的方波，在这种情况下，N 将为 \(N = 53\)。根据情况，对于实时应用程序中的每个样本可能需要执行过多的迭代。一种常见的方法是使用加法合成预先生成不同频率的所需波形，并将生成的样本存储在波表中，然后可以在运行时查找。但这需要预先计算和内存。该方法利用了从方波的低值到方波的高值的步进在理论上是瞬时的这一事实。因此，无论主方波设置为什么周期或频率，执行跳跃的波部分将始终相同。在数字情况下，这意味着这一步只能在采样率允许的情况下尽可能尖锐并具有尽可能多的频率。

诀窍是让此步骤的版本受带宽限制。然后可以将此“BLEP”替换到具有太陡步的天真生成的波形区域（不考虑混叠），从而避免更高的频率。由于 BLEP 需要以不连续性为中心，我们需要能够提前预测几个样本何时发生阶跃，或者只是稍微延迟信号。可以使用上述加法合成或通过模拟对通过低通滤波器发送的脉冲的影响来生成这些步骤。 BLEP 方法的发展是所谓的 polyBLEP。这里给出了带限阶跃函数的多项式近似。这允许快速计算带限步长并且不需要内存查找。这种近似是通过数学解释将阶跃函数通过低通滤波器的过程得出的。 Valimaki 和 Huovilainen 在他们的论文中描述了这个过程。三角波可以用作低通滤波的函数，这是因为它的频率响应具有 sinc^2 的形状。三角波函数由下式给出，三角脉冲的频率响应可用于使用卷积改变方波。方波与三角脉冲卷积。在此过程中，三角脉冲得到积分。常量 \(\frac{1}{2}\) 项是为了确保连续性。这是多项式限带步骤。要应用于原始生成波形的实际校正由 BLEP 残差给出。也就是说，从 BLEP 中减去输入脉冲。我们得到 $$f_{polyBLEP}(t) = \begin{cases}\frac{t^2}{2} + t + \frac{1}{2} \, , \text{when} \, -1 \leq t \leq 0 \newlinet - \frac{t^2}{2} + \frac{1}{2} \, , \text{when} \, 0 \lt t \leq 1\end{cases} \, .$$

在程序中实现这一点意味着我们需要一些逻辑来找到不连续性，然后应用更正。我们需要提前知道我们的不连续性何时出现才能应用校正，但由于波形是周期性的，这应该不会太困难。下面是在 polyBLEPS 的帮助下生成方波校正的示例。 # PolyBLEP by Tale # http://www.kvraudio.com/forum/viewtopic.php?t=375517 #（由我改编）#double t = 0.; // 0 <= t < 1 #double dt = freq / sample_rate; def poly_blep ( t , dt ): """对天真生成的波形中给出的不连续性进行 polyblep 校正。dt = freq/Fs """ if (t < dt): # 0 <= t < 1 t /= dt ; # 2 * (t - t^2/2 - 0.5) 返回 t + t - t * t - 1 elif ( t > 1. - dt ): # -1 < t < 0 t = ( t - 1. ) / dt ; # 2 * (t^2/2 + t + 0.5) return t * t + t + t + 1 else : # 0 否则返回0 def poly_square ( t , dt , aa ): naive_square = 0 if t <= 0.5 : naive_square = 1 ;如果 t > 0.5 : naive_square = - 1 ; if aa : #with polyBlep 校正。 #一个用于上升，一个用于下降（0.5） return naive_square + poly_blep ( t , dt ) \ - poly_blep (( t + 0.5 ) % 1 , dt ) else : return naive_square t = 0 Fs = 3520 f0 = 240 dt = f0 / Fs out_no_aa = np .数组 ([ 0 ] * Fs , float ) out_aa = np 。 array ([ 0 ] * Fs , float ) # 使锯齿波 f0 次 t = np . linspace ( 0 , f0 , Fs ) for i in range (0, len (t)): t_tmp = t [ i ] % 1 out_no_aa [ i ] = poly_square ( t_tmp , dt , 0 ) out_aa [ i ] = poly_mpt , dt , 1 ) t = np 。 linspace ( 0 , 1 , Fs ) #f0 次在 1 秒内 signal_no_aa = out_no_aa signal_aa = out_aa 我们可以通过查看频率响应来了解这对混叠的影响。 ＃阴谋！血小板图 ( figsize = ( 16 , 20 )) plt .子图 (4, 2, 1) plt。 plot ( t , signal_aa , ".-" ) plt. title ("Square with polyBLEP") plt . xlabel ( "time [s]" ) plt 。 xlim ( 0 , 0.05 );血小板ylim (-2, 2) plt。子图 (4, 2, 2) plt。 plot ( t , signal_no_aa , ".-" ) plt。标题（“没有 polyBLEP 的正方形”） plt 。 xlabel ( "time [s]" ) plt 。 xlim ( 0 , 0.05 );血小板ylim ( - 2 , 2 ) freqs , ps = fft_spectrum ( t , signal_aa , Fs ) plt.子图 (4, 2, 3) plt。绘图（频率，np。sqrt（ps），“.-”）plt。标题（“polyBLEP 的频率响应”） plt 。 xlabel ( "频率 [Hz]" ) plt . xlim (0, 2000) plt。 ylim ( 0 , 400 ) freqs , ps = fft_spectrum ( t , signal_no_aa , Fs ) plt.子图 (4, 2, 4) plt。绘图（频率，np。sqrt（ps），“.-”）plt。标题（“没有 polyBLEP 的频率响应”） plt 。 xlabel ( "频率 [Hz]" ) plt . xlim (0, 2000) plt。 ylim ( 0 , 400 ) 我们有它。在没有 polyBLEPS 的信号中，我们看到了很多不需要的频率。任何高度在 150 或以下的峰都可以被认为是混叠（除了最后一个峰）。使用 polyBLEP 校正调整信号时，我们看到混叠显着减少！下面是两个视频，您可以在其中听到 polyBLEPS 的效果。这里使用了 22050 的采样率。音量警告！

返回顶部