数字信号处理——多速率信号处理(1)
本文详细介绍了多速率信号处理的核心技术及其应用:
抽取(Decimation)
抽取操作用于减少信号的采样率。通过数学表达式:
y[k] = x\left( \left\lfloor \frac{k}{D} \right\rfloor \right)
实现后减少了数据量并去除了不必要的频率成分。Matlab代码使用downsample函数并配置抗混叠滤波器以确保频谱完整性。
内插(Upsampling)
内插用于增加信号的采样率。通过数学表达式:
y[i] = x\left( i + S\mod D \right)
实现后提高了信号质量并扩展了频谱范围。Matlab代码使用upsample函数并通过零填充技术完成。
采样率分数倍变换
这种变换结合了内插和抽取两种操作:先进行高精度的内插以防止频谱失真,再进行低复杂度的数据缩减以降低硬件开销。通过设计一个等效滤波器确保整体系统的稳定性与效率。
仿真结果验证了每种操作的有效性:提取样本、频谱分析以及相位信息均表明所采用的方法能够有效优化信号处理流程,在保持高质量的同时显著降低了数据流的速度需求。这种方法在现代数字信号处理领域具有重要应用价值,尤其在需要高效资源利用的应用场景中表现突出。

目录
引言
1、抽取
MATLAB仿真代码
仿真结果
抽取滤波器
2、内插
MATLAB仿真代码
仿真结果
内插滤波器
3、采样率分数倍变换
MATLAB仿真源码
仿真结果
参考说明
引言
多速率信号处理在当前数字信号处理领域发挥着关键作用。随着ADC速率的不断提升,在实际应用中能够实现更高的采样率。这种提升使得数字化处理技术能够更接近射频前端的核心环节。从软件无线电架构来看,其核心包含三种基本采样方式:射频直接采样、中频采样和低通采样。其中,中频直接采样方案因其灵活性而得到了广泛应用。通过采用过采样技术,在保持相同数据精度的前提下显著提升了系统的信噪比性能。在ADC完成信号采集后段落中的信噪比计算公式如下所示:

然而,在这种情况下(采样速率过高),对应的数据显示流速度也随之急剧上升。这使得后端信号处理环节面临巨大的压力并造成大量的硬件资源消耗。为此引入抽取滤波器便成为了解决这一矛盾的有效手段——它不仅能够降低数据传输速率还能缓解后端的压力。相反地讲 内插技术则是在降低数据传输速率的基础上通过内插操作来提升采样率从而确保高速DAC能够正常驱动。其中核心组件的主要框图如下所示:

本博文系列基于FPGA的硬件平台进行多速率信号处理的研究与探讨。后续内容将提供详细的仿真设计源代码。
本文主要阐述多速率信号处理的基本理论,并通过MATLAB仿真展示可视化效果。同时提供相应的MATLAB代码。
1、抽取
改写说明

然而,在这 M 个样点中选择哪一个具有一定的自由度,并将其定义为相位偏移参数;当抽取倍数等于 M 时,在这种情况下相位偏移参数的可选值共有 M 种可能的选择方式,并且这些可选值依次为从 0 到 M-1 的所有整数值。
举例说明:
采用LFM信号作为研究对象进行分析。具体参数设置如下:将采样率设置为100MHz,并设定其带宽为2.5兆赫;抽样倍数取值为了10;同时确保原始序列的采样点总数设定为了2000个。
MATLAB仿真代码
> 1. %% ---- ---- ---- ---- ---- Data Down sample Simulation ---- ---- ---- ---- ----
>
> 2. % | Author :Xu Y. B.( 昵称:在路上,正出发)
>
> 3. % | Abstract :
>
> 4. % | -1- Signal Type : Linear Frequency Modulation Signal
>
> 5. % | -2- Signal Bandwidth : 5 MHz (0-5MHz)
>
> 6. % | -3- Sample Rate : 100 MHz
>
> 7. % | -4- Samples : 2000
>
> 8. % | -5- Decimation Multiple: 10
>
> 9. % |
>
> 10. % ---------------------------------------------------------------------------------
>
> 11.
>
> 12. %% Prepare
>
> 13. clc;
>
> 14. clearvars;
>
> 15. close all;
>
> 16. set(0,'defaultfigurecolor','w')
>
> 17.
>
> 18. %% Parameter Setting
>
> 19. B = 2.5e6; % Bandwidth ,unit:Hz
>
> 20. N = 2e3; % Samples
>
> 21. Fs = 100e6; % Sample Rates ,unit:Hz
>
> 22. D = 10; % Decimation Multiple
>
> 23. K = B*Fs/N; % Chirp Rate,unit:Hz/s
>
> 24.
>
> 25. %% Signal Genarating
>
> 26. t = (0:N-1)/Fs; % Time Axis
>
> 27. S = cos(pi*K*t.^2); % Sample Sequence
>
> 28.
>
> 29. % Signal Verifying
>
> 30. figure;
>
> 31. subplot(211)
>
> 32. plot(S,'b')
>
> 33. xlabel('Samples Index')
>
> 34. ylabel('Amplitude')
>
> 35. title('Time Domain Signal Waveform')
>
> 36. subplot(212)
>
> 37. plot((linspace(-Fs/2,Fs/2,N))/1e6,abs(fftshift(fft(S))),'r')
>
> 38. xlabel('Frequence /MHz')
>
> 39. ylabel('Amplitude')
>
> 40. title('The Amplitude Spectrum of LFM Signal')
>
> 41.
>
> 42. %% Down Sampling
>
> 43. % -1- Use <downsample> function
>
> 44. figure;
>
> 45. for i = 0:D-1
>
> 46. subplot(5,2,i+1);
>
> 47. plot(downsample(S,D,i),'b')
>
> 48. axis tight
>
> 49. xlabel('Smaples')
>
> 50. ylabel('Amplitude')
>
> 51. title(['Phase Offset is ',num2str(i)])
>
> 52. end
>
> 53. % -2- Direct Operation
>
> 54. figure;
>
> 55. for i = 1:D
>
> 56. subplot(5,2,i);
>
> 57. plot(S(i:D:end),'r')
>
> 58. axis tight
>
> 59. xlabel('Smaples')
>
> 60. ylabel('Amplitude')
>
> 61. title(['Phase Offset is ',num2str(i-1)])
>
> 62. end
>
> 63. 64. %% Spectrum Analysis After Down Sampling
>
> 65. % Amplitude Spectrum
>
> 66. figure;
>
> 67. for i = 1:D
>
> 68. subplot(5,2,i);
>
> 69. plot(linspace(-Fs/2/D,Fs/2/D,N/D)/1e6,abs(fftshift(fft(S(i:D:end)))),'k')
>
> 70. axis tight
>
> 71. xlabel('Frequence /MHz')
>
> 72. ylabel('Amplitude')
>
> 73. title(['Amplitude Spectrum Analysis After Down Sampling (Phase Offset is ',num2str(i-1),')'])
>
> 74. end
>
> 75. 76. % Phase Spectrum
>
> 77. figure;
>
> 78. for i = 1:D
>
> 79. subplot(5,2,i);
>
> 80. plot(linspace(-Fs/2/D,Fs/2/D,N/D)/1e6,angle(fftshift(fft(S(i:D:end)))),'k')
>
> 81. axis tight
>
> 82. xlabel('Frequence /MHz')
>
> 83. ylabel('Phase')
>
> 84. title(['Phase Spectrum Analysis After Down Sampling (Phase Offset is ',num2str(i-1),')'])
>
> 85. end
>
> 86.
>
> 87.
>
> 88.
>
> 89.
>
> 90.
>
> 91.
>
>
> AI助手
仿真结果
原始采样波形及其频谱验证:

不同相位偏移的下采样结果:


下采样后的幅度谱:

相位谱:

抽取滤波器
提取滤波器又被称作抗混叠滤波器。我们知道,在数字域中,信号的频率是有限制的;无穷大并非实际存在的概念。因此,在采样率fs下工作的采样系统而言,则其允许的最大输入信号频率应设定为fs/2;一旦超越这一临界值,则会导致频谱混叠现象的发生。而抽取操作则是基于原有采样率进行降采样的过程;假设采用抽取倍数D,则经过抽取后的无混叠最大频率应被限制在fs/(2D)范围内;因此,在执行抽取运算之前必须安装一个截止频率为fs/(2D)的理想抗混叠滤波器以避免频谱混叠的发生。至于这一类滤波器的具体实现结构,则将在后续的内容中逐一进行详细阐述。
2、内插
内插即是在已有的采样点序列基础上补充零值点以实现数据扩展的过程。其数学表达式为:

与上述抽取的仿真类似,内插的方针如下:
MATLAB仿真代码
> 1. %% ---- ---- ---- ---- ---- Data Up sample Simulation ---- ---- ---- ---- ----
>
> 2. % | Author :Xu Y. B.( 昵称:在路上,正出发)
>
> 3. % | Abstract :
>
> 4. % | -1- Signal Type : Linear Frequency Modulation Signal
>
> 5. % | -2- Signal Bandwidth : 5 MHz (0-5MHz)
>
> 6. % | -3- Sample Rate : 100 MHz
>
> 7. % | -4- Samples : 2000
>
> 8. % | -5- Interpolation Multiple : 10
>
> 9. % | Matlab Version : R2022a
>
> 10. % |
>
> 11. % ---------------------------------------------------------------------------------
>
> 12.
>
> 13. %% Prepare
>
> 14. clc;
>
> 15. clearvars;
>
> 16. close all;
>
> 17. set(0,'defaultfigurecolor','w')
>
> 18.
>
> 19. %% Parameter Setting
>
> 20. B = 2.5e6; % Bandwidth ,unit:Hz
>
> 21. N = 2e3; % Samples
>
> 22. Fs = 100e6; % Sample Rates ,unit:Hz
>
> 23. I = 10; % Interpolation Multiple
>
> 24. K = B*Fs/N; % Chirp Rate,unit:Hz/s
>
> 25.
>
> 26. %% Signal Genarating
>
> 27. t = (0:N-1)/Fs; % Time Axis
>
> 28. S = cos(pi*K*t.^2); % Sample Sequence
>
> 29.
>
> 30. % Signal Verifying
>
> 31. figure;
>
> 32. subplot(211)
>
> 33. plot(S,'b')
>
> 34. xlabel('Samples Index')
>
> 35. ylabel('Amplitude')
>
> 36. title('Time Domain Signal Waveform')
>
> 37. subplot(212)
>
> 38. plot((linspace(-Fs/2,Fs/2,N))/1e6,abs(fftshift(fft(S))),'r')
>
> 39. xlabel('Frequence /MHz')
>
> 40. ylabel('Amplitude')
>
> 41. title('The Amplitude Spectrum of LFM Signal')
>
> 42.
>
> 43. %% Up Sampling
>
> 44. % -1- Use <upsample> function
>
> 45. figure;
>
> 46. for i = 0:I-1
>
> 47. subplot(5,2,i+1);
>
> 48. plot(upsample(S,I,i),'b')
>
> 49. axis tight
>
> 50. xlabel('Smaples')
>
> 51. ylabel('Amplitude')
>
> 52. title(['Phase Offset is ',num2str(i)])
>
> 53. end
>
> 54. % -2- Direct Operation
>
> 55. figure;
>
> 56. for i = 1:I
>
> 57. subplot(5,2,i);
>
> 58. plot(reshape([zeros(N,i-1),S.',zeros(N,I-i)].',[],1),'r')
>
> 59. axis tight
>
> 60. xlabel('Smaples')
>
> 61. ylabel('Amplitude')
>
> 62. title(['Phase Offset is ',num2str(i-1)])
>
> 63. end
>
> 64. 65. %% Spectrum Analysis After Up Sampling
>
> 66. % Amplitude Spectrum
>
> 67. figure;
>
> 68. for i = 1:I
>
> 69. subplot(5,2,i);
>
> 70. plot(linspace(-Fs/2*I,Fs/2*I,N*I)/1e6,abs(fftshift(fft(upsample(S,I,i-1)))),'k')
>
> 71. axis tight
>
> 72. xlabel('Frequence /MHz')
>
> 73. ylabel('Amplitude')
>
> 74. title(['Amplitude Spectrum Analysis After Up Sampling (Phase Offset is ',num2str(i-1),')'])
>
> 75. end
>
> 76. 77. % Phase Spectrum
>
> 78. figure;
>
> 79. for i = 1:I
>
> 80. subplot(5,2,i);
>
> 81. plot(linspace(-Fs/2*I,Fs/2*I,N*I)/1e6,angle(fftshift(fft(upsample(S,I,i-1)))),'k')
>
> 82. axis tight
>
> 83. xlabel('Frequence /MHz')
>
> 84. ylabel('Phase')
>
> 85. title(['Phase Spectrum Analysis After Up Sampling (Phase Offset is ',num2str(i-1),')'])
>
> 86. end
>
> 87.
>
> 88.
>
> 89.
>
> 90.
>
> 91.
>
> 92.
>
>
> AI助手
仿真结果
原始采样波形和频谱与上述仿真结果一致。在此处不展示原始数据图, 主要展示内插处理后的结果图。
内插后时域波形:


内插后幅度谱:

相谱

内插滤波器
-fs/2,fs/2
-fsI/2,fsI/2
内插倍数设定为 I ,内插运算将引起原本位于【-fs/2, fs/2
3、采样率分数倍变换
小数频率转换器等效于系统基准的小数倍频率转换器。为了防止频谱失真,在实现过程中应当优先实施内插运算,并随后进行抽取运算级联。
其处理结构:

中间的滤波结构本质上是由内插_filter_和抽取_filter_串联而成的一个复合_filter_。该复合_filter_所对应的截止频率等于其两个子filter_截止频率中较低的那个值。具体而言,其数学表达式如下所示;

MATLAB仿真源码
> 1. %% ---- ---- ---- ---- ---- Data Fractional Multiple Sample Rate Simulation ---- ---- ---- ---- ----
>
> 2. % | Author :Xu Y. B.( 昵称:在路上,正出发)
>
> 3. % | Abstract :
>
> 4. % | -1- Signal Type : Linear Frequency Modulation Signal
>
> 5. % | -2- Signal Bandwidth : 5 MHz (0-5MHz)
>
> 6. % | -3- Sample Rate : 100 MHz
>
> 7. % | -4- Samples : 2000
>
> 8. % | -5- Interpolation Multiple : 10
>
> 9. % | -6- Decimation Multiple : 4
>
> 10. % | -7- Phase Offset : 0
>
> 11. % | Matlab Version : R2022a
>
> 12. % |
>
> 13. % ---------------------------------------------------------------------------------
>
> 14.
>
> 15. %% Prepare
>
> 16. clc;
>
> 17. clearvars;
>
> 18. close all;
>
> 19. set(0,'defaultfigurecolor','w')
>
> 20.
>
> 21. %% Parameter Setting
>
> 22. B = 2.5e6; % Bandwidth ,unit:Hz
>
> 23. N = 2e3; % Samples
>
> 24. Fs = 100e6; % Sample Rates ,unit:Hz
>
> 25. I = 10; % Interpolation Multiple
>
> 26. D = 4; % Decimation Multiple
>
> 27. K = B*Fs/N; % Chirp Rate,unit:Hz/s
>
> 28.
>
> 29. %% Signal Genarating
>
> 30. t = (0:N-1)/Fs; % Time Axis
>
> 31. S = cos(pi*K*t.^2); % Sample Sequence
>
> 32.
>
> 33. % Signal Verifying
>
> 34. figure;
>
> 35. subplot(311)
>
> 36. plot(S,'b')
>
> 37. xlabel('Samples Index')
>
> 38. ylabel('Amplitude')
>
> 39. title('Time Domain Signal Waveform')
>
> 40. subplot(312)
>
> 41. plot((linspace(-Fs/2,Fs/2,N))/1e6,abs(fftshift(fft(S))),'r')
>
> 42. xlabel('Frequence /MHz')
>
> 43. ylabel('Amplitude')
>
> 44. title('The Amplitude Spectrum of LFM Signal')
>
> 45. subplot(313)
>
> 46. plot((linspace(-Fs/2,Fs/2,N))/1e6,angle(fftshift(fft(S))),'b')
>
> 47. xlabel('Frequence /MHz')
>
> 48. ylabel('Amplitude')
>
> 49. title('The Phase Spectrum of LFM Signal')
>
> 50.
>
> 51.
>
> 52. %% Up Sampling
>
> 53. S_UP_SAMPLE = upsample(S,I,0);
>
> 54. figure;
>
> 55. subplot(311)
>
> 56. plot(S_UP_SAMPLE,'b')
>
> 57. axis tight
>
> 58. xlabel('Smaples')
>
> 59. ylabel('Amplitude')
>
> 60. title('Time Domain Waveform After Up Sampling')
>
> 61. 62. %% Spectrum Analysis After Up Sampling
>
> 63. % Amplitude Spectrum
>
> 64. subplot(312)
>
> 65. plot(linspace(-Fs/2*I,Fs/2*I,N*I)/1e6,abs(fftshift(fft(S_UP_SAMPLE))),'k')
>
> 66. axis tight
>
> 67. xlabel('Frequence /MHz')
>
> 68. ylabel('Amplitude')
>
> 69. title('Amplitude Spectrum Analysis After Up Sampling ')
>
> 70. 71. % Phase Spectrum
>
> 72. subplot(313)
>
> 73. plot(linspace(-Fs/2*I,Fs/2*I,N*I)/1e6,angle(fftshift(fft(S_UP_SAMPLE))),'r')
>
> 74. axis tight
>
> 75. xlabel('Frequence /MHz')
>
> 76. ylabel('Phase')
>
> 77. title('Phase Spectrum Analysis After Up Sampling ')
>
> 78. 79. 80. %% LP Filtering
>
> 81. Order = 40;
>
> 82. f_stop = min(Fs/2,Fs*I/D/2);
>
> 83. Wn = f_stop/(Fs*I/2);
>
> 84. H = fir1(Order , Wn ,"low");
>
> 85. freqz(H,1);
>
> 86. S_FILTER = filter(H,1,S_UP_SAMPLE);
>
> 87. figure;
>
> 88. subplot(311)
>
> 89. plot(S_FILTER,'b')
>
> 90. axis tight
>
> 91. xlabel('Smaples')
>
> 92. ylabel('Amplitude')
>
> 93. title('Time Domain Waveform After Filtering')
>
> 94. 95. %% Spectrum Analysis After Filtering
>
> 96. % Amplitude Spectrum
>
> 97. subplot(312)
>
> 98. plot(linspace(-Fs/2*I,Fs/2*I,N*I)/1e6,abs(fftshift(fft(S_FILTER))),'k')
>
> 99. axis tight
>
> 100. xlabel('Frequence /MHz')
>
> 101. ylabel('Amplitude')
>
> 102. title('Amplitude Spectrum Analysis After Filtering ')
>
> 103. 104. % Phase Spectrum
>
> 105. subplot(313)
>
> 106. plot(linspace(-Fs/2*I,Fs/2*I,N*I)/1e6,angle(fftshift(fft(S_FILTER))),'r')
>
> 107. axis tight
>
> 108. xlabel('Frequence /MHz')
>
> 109. ylabel('Phase')
>
> 110. title('Phase Spectrum Analysis After Filtering ')
>
> 111. 112. %% Down Sampling
>
> 113. S_DOWN_SAMPLE = downsample(S_FILTER,D,0);
>
> 114. figure;
>
> 115. subplot(311);
>
> 116. plot(S_DOWN_SAMPLE,'b')
>
> 117. axis tight
>
> 118. xlabel('Smaples')
>
> 119. ylabel('Amplitude')
>
> 120. title('Time Domain Waveform After Down Sampling')
>
> 121. 122. 123. %% Spectrum Analysis After Down Sampling
>
> 124. % Amplitude Spectrum
>
> 125. subplot(312);
>
> 126. plot(linspace(-Fs*I/2/D,Fs*I/2/D,N*I/D)/1e6,abs(fftshift(fft(S_DOWN_SAMPLE))),'k')
>
> 127. axis tight
>
> 128. xlabel('Frequence /MHz')
>
> 129. ylabel('Amplitude')
>
> 130. title('Amplitude Spectrum Analysis After Down Sampling ')
>
> 131.
>
> 132. % Phase Spectrum
>
> 133. subplot(313);
>
> 134. plot(linspace(-Fs*I/2/D,Fs*I/2/D,N*I/D)/1e6,angle(fftshift(fft(S_DOWN_SAMPLE))),'r')
>
> 135. axis tight
>
> 136. xlabel('Frequence /MHz')
>
> 137. ylabel('Phase')
>
> 138. title('Phase Spectrum Analysis After Down Sampling')
>
>
>
>
> AI助手
仿真结果
原始采样:

内插:

滤波:

抽取:

参考说明
【1】高亚军.基于FPGA的数字信号处理.
