尧图网站建设 
简介面向数字信号处理学习者与MATLAB开发者这份4KB的小型代码包聚焦多速率信号处理中的采样率转换与抗混叠滤波器设计。资源共2个文件包含可运行的main.m脚本与配套README.md说明文档。脚本基于MATLAB信号处理工具箱完成信号插值、抽取及抗混叠低通滤波的仿真配合文档可直观理解resample、upfirdn等核心函数的触发条件、参数设置与结果验证方法并掌握不同采样率转换场景下滤波器阶数与截止频率的选取思路。包体虽小但结构紧凑适合作为课程设计、毕业设计或工程入门实践的参考模板帮助读者在短时间内搭建完整的采样率转换仿真流程并延伸到数字通信、语音处理等实际应用。目前已有65人学习下载可作为信号处理方向快速上手的轻量示例。1. 多速率信号处理与采样率转换为什么设计系统先设计抗混叠滤波器通信基站的ADC输出采样率是65.536 MHz符号速率却是2 Msps音频设备要在48 kHz与44.1 kHz之间来回转换雷达数字中频之后还要做十倍以上的抽取。这些场景指向同一个问题采样率转换。多数人的第一反应是抽点或插值但工程上做完之后的结论通常是反过来的——先算清楚频谱折叠发生在哪个频点再决定抗混叠滤波器长什么样否则所有指标都是空谈。多速率信号处理的核心不是“多采样率”本身而是速率变化前后频谱如何搬移、叠加、镜像。抽取会把高频折叠进低频插值会在通带外侧制造镜像副本这两种伪影都必须靠预滤波器压到目标底噪以下。本文从抽取与插值的频率规则讲起用MATLAB把几种常用采样率转换路径跑通再落到抗混叠滤波器的参数设计与验证上。内容适合刚接触数字信号处理的工程师也适合那些常年用resample但说不清边界条件的老手。2. 采样率转换的频谱折叠原理抽取、插值在MATLAB里如何验证2.1 抽取的频谱扩张折叠进带内的信号来自哪一段抽取就是每M点保留一点x_down[n] x[Mn]。频域对应关系为X_d(e^{jω}) 1/M Σ X(e^{j(ω−2πk)/M})含义是原频谱被拉伸到M倍然后以2π为周期做叠加。当原始信号中最高频率成分超过fs/(2M)时超出部分就会折叠回低频区和有效信号叠在一起。这个折叠不是噪声是确定性的频谱搬移所以混叠具体出现在哪个频点是可计算的。用MATLAB验证这一段规则最直接。生成1 kHz与2.5 kHz两个正弦16 kHz采样M4直接抽取。抽取后采样率变成4 kHz奈奎斯特边界只有2 kHz2.5 kHz会折叠到4 kHz − 2.5 kHz 1.5 kHz而1 kHz的位置不变。fs0 16000; t 0:1/fs0:0.2; x sin(2*pi*1000*t) sin(2*pi*2500*t); M 4; xd x(1:M:end); % 直接抽取没有抗混叠滤波 Nfft 2048; f0 (0:Nfft/2-1)/Nfft*fs0; fd (0:Nfft/2-1)/Nfft*fs0/M; plot(f0, 20*log10(abs(fft(x, Nfft)))); hold on; plot(fd, 20*log10(abs(fft(xd, Nfft))), r); legend(原信号 16 kHz, 直接抽取 4 kHz); xlabel(Hz); ylabel(dB);这段代码里fft长度取了2048直接用正弦做频谱会有泄漏但峰值位置足够说明问题。运行后能看到原始2.5 kHz的谱峰在抽取后移动到1.5 kHz而1 kHz谱峰不变这就是没有抗混叠滤波时混叠发生的典型形态。反过来想如果有用信号恰好落在1.5 kHz那么来自2.5 kHz的混叠分量会直接叠加在它上面后续任何处理都无法再分开。2.2 插值形成的镜像频谱低通截止频率怎么定插值是在每两个原始样本之间插入L−1个零数学形式是x_up[n] x[n/L]当n是L的倍数时频域表现为频谱压缩到原来的1/L并且在原采样率fs的整数倍位置上复制出镜像副本。以16 kHz采样、1 kHz信号为例L4插值后采样率变为64 kHz镜像出现在1 kHz、15 kHz、17 kHz、31 kHz等位置其中15 kHz是−1 kHz的副本。L 4; xu upsample(x, L); % 中间插零不滤波 f_up (0:Nfft/2-1)/Nfft*fs0*L; plot(f_up, 20*log10(abs(fft(xu, Nfft)))); xlabel(Hz); ylabel(dB);插值后的目标是把原始带宽之外的镜像滤掉。低通截止频率取fs/(2L)2 kHz就能保留1 kHz成分并把15 kHz、17 kHz等副本全部压掉。这里的重点是插值本身不会损坏信号损坏发生在你忘记滤波的时候。真正决定重建质量的是滤波器通带平坦度和对第一镜像的抑制量。2.3 抗混叠滤波器的位置先滤波后抽取先插值后滤波把两节结论放在一起抽取之前必须低通把超过目标奈奎斯特边界的分量提前去掉插值之后必须低通把零插入产生的镜像清掉。滤波器放在这两处不是约定俗成而是频谱规则决定的。后面章节里出现的resample、upfirdn、firdecim等函数本质上都是在做“滤波抽取”或“插值滤波”的组合。理解了这个顺序才能看懂MATLAB函数内部为什么要先设计一个多相抗混叠滤波器。提示抽取的滤波和插值的滤波虽然都叫抗混叠但一个针对折叠一个针对镜像。两者对阻带起始频率的要求不同设计滤波器时不要套用同一组参数。3. 用MATLAB实现采样率转换resample、upfirdn与整数倍抽取/插值函数3.1 最快路径用resample完成有理数采样率转换MATLAB里最常被调用的多速率函数是resample一句话写法是y resample(x, p, q)把输入按p/q倍的速率重采样。p和q是互质整数目标采样率与源采样率之比等于p/q。例如48 kHz转44.1 kHz比值是0.91875用rat函数得到分子分母fs_in 48000; fs_out 44100; [p, q] rat(fs_out / fs_in); % p147, q160 y resample(x, p, q);resample内部采用Kaiser窗设计一个多相抗混叠滤波器然后依次做p倍插值、滤波、q倍抽取。对大多数原型验证场景这一条命令就够用不需要自己设计滤波器。resample有两个可选参数值得注意L控制滤波器长度L越大过渡带越窄beta控制Kaiser窗形状间接决定阻带衰减。默认值在多数场合表现尚可但如果你对混叠底噪有明确指标就应该自己设计h并改用3.2节的方式。注意resample输入输出长度会做自动对齐输出长度大致为ceil(length(x)*p/q)但首尾会有瞬态做符号同步或时延敏感处理时不能忽略。3.2 自定义滤波器用upfirdn把插值、滤波、抽取拆开当你需要精确控制抗混叠滤波器指标时resample的黑盒就不够用了。upfirdn(x, h, p, q)是更底层的原语先对x做p倍零插值与h卷积再做q倍抽取。三个操作在一条调用里完成但h需要你预先设计好。resample内部做的事情和它基本一致。p 147; q 160; fs_mid fs_in * p; % 插值后的中间采样率 7.056 MHz fp 20000; % 通带边界 fst 22500; % 阻带起始需覆盖输出Nyquist df [0 fp fst fs_mid/2] / (fs_mid/2); % 以中间采样率Nyquist做归一化 a [1 1 0 0]; % 通带1阻带0 dev [1e-4 1e-5]; % 通带纹波与阻带衰减线性幅度 [n, fo, ao, w] firpmord(df, a, dev); h firpm(n, fo, ao, w); y upfirdn(x, h, p, q);这里最容易被忽略的是h的频率轴。upfirdn中h工作在插值后的采样率fs_mid上滤波器设计必须换算到fs_mid的奈奎斯特频率。如果直接拿输入采样率的Nyquist做归一化实际滤波器通带会宽p倍镜像一个都滤不掉。这段代码里的fst取22500 Hz略高于输出奈奎斯特22050 Hz给过渡带留了少量余量工程上更稳妥。3.3 整数倍抽取与插值firdecim、firinterp和intfilt当转换比为整数时firdecim和firinterp是两个专用入口。firdecim(h, M, x)等价于先filter再downsample但内部对h做多相分解能显著减少计算量。firinterp同理做L倍插值后滤波。它们的典型调用如下M 4; h fir1(256, 1/M); % 截止在fs/(2M)的简单低通 yd firdecim(h, M, x); % 抽取输出采样率fs/M L 4; h2 fir1(256, 1/L); yu firinterp(h2, L, x); % 插值输出采样率fs*Lintfilt是更老的函数用于生成整数倍插值滤波器现在大部分新代码已经被firinterp或designMultirateFIR取代。新版本MATLAB里designMultirateFIR更值得关注它把抗混叠滤波器设计与抽取/插值参数统一成一个接口R2022b之后我一般在正式项目里优先用它。3.4 函数选型对比这四种方式别再混着用函数完成操作关键参数典型场景resample插值滤波抽取有理数转换p, q, L, beta快速原型、算法验证upfirdn插值滤波抽取可自定义hh, p, q需要指定滤波器指标firdecim滤波M倍抽取h, M固定整数倍抽取firinterpL倍插值滤波h, L固定整数倍插值designMultirateFIR设计抽取/插值滤波器filter, decim, dev, w新工程推荐很多人犯的错是用downsample直接做抽取x(1:M:end)得到的结果没有经过任何滤波高频混叠全部落在带内。对比一下误用和正确路径的底噪结论非常明显y_err x(1:4:end); % 直接抽取无滤波 y_ok resample(x, 1, 4); % 完整抗混叠链路 snr_err snr(y_err); % 通常只有20~30 dB snr_ok snr(y_ok); % 能到60 dB以上snr函数在Signal Processing Toolbox里只做指标对比可以接受。真正做音频或通信系统时还要查看混叠落入带内的具体频点而不仅仅是总SNR。4. 抗混叠滤波器设计窗函数、等波纹、多级级联与48 kHz到44.1 kHz实例4.1 滤波器指标怎么从系统需求推导出来设计抗混叠滤波器之前先定三个数通带边界fp、阻带起始频率fst、阻带衰减As。fp通常是信号本身的有效带宽fst取决于目标采样率下的奈奎斯特边界或镜像/混叠区最低频率As由系统底噪决定。以48 kHz到44.1 kHz音频转换为例信号内容只保留到20 kHz输出奈奎斯特是22.05 kHz所以fp20 kHz、fst22.05 kHz阻带衰减通常要求90 dB以上否则可听底噪里会混入镜像。设计方法优点缺点适用场景窗函数法 fir1简单、稳定、线性相位通带纹波与阻带衰减互相牵扯快速滤波、指标不苛刻等波纹 firpm相同阶数阻带抑制更均衡阶数估计不准时需要迭代抽取/插值混叠抑制最小二乘 firls通带纹波总体最低阻带可能出现个别超差对总失真敏感的音频链路4.2 等波纹设计firpmord配合firpm的完整流程等波纹法的优势是阻带内抑制均匀不会出现窗函数法那种靠近截止频率处凹陷更深、远处回升的情况。对混叠抑制来说所需的是整个阻带都低于某个阈值等波纹法更贴合这个需求。用firpmord先估阶数再交给firpm设计是标准流程。% 中间采样率 7.056 MHz归一化以插值后Nyquist为1 fs_mid 48000 * 147; % 7.056 MHz f_norm [0 20000 22050 fs_mid/2] / (fs_mid/2); a [1 1 0 0]; dev [0.0005 1e-5]; % 通带纹波约0.0087 dB阻带100 dB [n, fo, ao, w] firpmord(f_norm, a, dev); h firpm(n, fo, ao, w);这段代码会算出非常高的阶数因为过渡带只有2050 Hz而中间采样率高达7.056 MHz相对过渡带宽度不到0.15%。这是多速率设计里的经典两难转换比越大镜像频点越密集对滤波器的要求越高。实际工程中不会让过渡带卡这么紧一般会把fst放宽到24 kHz甚至更高让滤波器阶数从数千降到几百代价是损失一点高频带宽。4.3 多级级联什么时候拆成两级更划算当抽取因子D超过10并且过渡带相对窄时单级滤波器阶数会爆炸。把D拆成D1×D2先按D1做第一次抽取再用D2做第二次每一级的过渡带都可以放宽总滤波器长度通常只有单级的几分之一。经验法则每级抽取因子尽量控制在10以内过渡带宽度至少占该级奈奎斯特带宽的20%到50%否则级联的优势不明显。% 单级抽取 D12 h12 designMultirateFIR(1, 12, 64, 90); y1 upfirdn(x, h12, 1, 12); % 两级抽取 3×4 h3 designMultirateFIR(1, 3, 64, 90); h4 designMultirateFIR(1, 4, 64, 90); y2 upfirdn(upfirdn(x, h3, 1, 3), h4, 1, 4);4.4 完整示例48 kHz音频转44.1 kHz的自定义滤波器链路把前面的内容串成一个可运行脚本。信号源用20 kHz以内的多音合成方便在频谱上直接看镜像是否被压住fs_in 48000; t 0:1/fs_in:1; x sin(2*pi*1000*t) 0.5*sin(2*pi*10000*t) 0.2*sin(2*pi*19900*t); p 147; q 160; fs_mid fs_in * p; f_norm [0 20000 22050 fs_mid/2] / (fs_mid/2); a [1 1 0 0]; dev [0.0005 1e-5]; [n, fo, ao, w] firpmord(f_norm, a, dev); n n rem(n, 2); % 保证线性相位为偶数阶 h firpm(n, fo, ao, w); y upfirdn(x, h, p, q); y y(1:round(length(x)*p/q)); % 截断到预期输出长度 fs_out fs_in * p / q;脚本中把滤波器阶数向上取偶是为了让线性相位群延迟正好是整数样本点后面做时延补偿时会省事。运行后输出采样率是44100 Hz频谱图在20 kHz以上应该看不到超过−90 dB的镜像残留。如果看到明显的20 kHz以上谱线先查f_norm是不是按fs_mid归一化这是最容易出错的点。5. 多速率系统的多相分解与验证群延迟补偿和混叠残留定位5.1 把滤波器拆开多相分解到底省在哪里抗混叠滤波器通常动辄上百阶直接在整个采样率上做卷积计算量很高。多相分解的思路是把长度N的FIR滤波器按抽取因子M拆成M个子滤波器每个子滤波器只在低速率的某一路并行工作。数学形式是H(z) Σ z^−k E_k(z^M)其中E_k(z^M)对应第k个子多相分量。MATLAB的firdecim和resample内部都是这种结构。手动验证多相结构的等价性可以用一段很短的MATLAB代码M 4; h fir1(96, 1/M); hp reshape([h(:); zeros(M - mod(length(h), M), 1)], length(h)M-mod(length(h),M), []); % hp的每一列就是一个多相子滤波器hp的每一列就是一路低速子滤波器输入序列按相位分配到各路卷积后再合并。子滤波器长度大致是主滤波器的1/M每路乘法次数大幅减少。这个结构在多速率系统里几乎是标配理解它才能理解为什么resample处理长序列时比“先filter后downsample”快那么多。5.2 群延迟怎么补线性相位滤波器的延迟是N/2线性相位FIR滤波器的群延迟是(N−1)/2个样本以滤波器工作采样率计。使用upfirdn时这个延迟折算到输出采样率要除以q。不补偿的话重采样后的信号整体偏移几十个样本做符号同步或者数据对齐时会直接错位。delay_mid (length(h) - 1) / 2; % 中间采样率下的延迟 delay_out delay_mid / q; % 换算到输出采样率 y upfirdn(x, h, p, q); y y(round(delay_out) 1 : end); % 去掉瞬态并补偿群延迟这段代码放在设计流程最后比用xcorr估计延迟更精确因为群延迟是已知的不需要用信号相关性去猜。复数调制信号或最小相位滤波器不适用这条规律最小相位滤波器的群延迟需要用grpdelay逐点计算。5.3 用pwelch和往返测试定位混叠残留验证抗混叠效果不能只看时域波形把转换前后的功率谱放在同一张图里对比是标准做法。pwelch返回的功率谱密度可以直接看到残留镜像落在哪个频段[pxx, f] pwelch(y, hann(4096), 2048, 8192, fs_out); [pxo, fo] pwelch(x, hann(4096), 2048, 8192, fs_in); semilogy(fo, sqrt(pxo), k); hold on; semilogy(f, sqrt(pxx), r); xline(22050/44100*2, --, Nyquist);更严格的做法是往返测试把信号从48 kHz转到44.1 kHz再转回48 kHz与原始信号比对。往返转换会同时暴露出混叠、镜像和滤波器非理想过渡带三类问题。若往返后的SNR低于系统指标把两个方向的滤波器阻带衰减各提高10 dB再测一次通常能找到问题出在哪一级。多速率系统的交付验证多数时候就靠这两步一张频谱图一个往返SNR。本文还有配套的精品资源点击获取 多速率信号处理与采样率转换：抗混叠滤波器设计实战指南
简介面向数字信号处理学习者与MATLAB开发者这份4KB的小型代码包聚焦多速率信号处理中的采样率转换与抗混叠滤波器设计。资源共2个文件包含可运行的main.m脚本与配套README.md说明文档。脚本基于MATLAB信号处理工具箱完成信号插值、抽取及抗混叠低通滤波的仿真配合文档可直观理解resample、upfirdn等核心函数的触发条件、参数设置与结果验证方法并掌握不同采样率转换场景下滤波器阶数与截止频率的选取思路。包体虽小但结构紧凑适合作为课程设计、毕业设计或工程入门实践的参考模板帮助读者在短时间内搭建完整的采样率转换仿真流程并延伸到数字通信、语音处理等实际应用。目前已有65人学习下载可作为信号处理方向快速上手的轻量示例。1. 多速率信号处理与采样率转换为什么设计系统先设计抗混叠滤波器通信基站的ADC输出采样率是65.536 MHz符号速率却是2 Msps音频设备要在48 kHz与44.1 kHz之间来回转换雷达数字中频之后还要做十倍以上的抽取。这些场景指向同一个问题采样率转换。多数人的第一反应是抽点或插值但工程上做完之后的结论通常是反过来的——先算清楚频谱折叠发生在哪个频点再决定抗混叠滤波器长什么样否则所有指标都是空谈。多速率信号处理的核心不是“多采样率”本身而是速率变化前后频谱如何搬移、叠加、镜像。抽取会把高频折叠进低频插值会在通带外侧制造镜像副本这两种伪影都必须靠预滤波器压到目标底噪以下。本文从抽取与插值的频率规则讲起用MATLAB把几种常用采样率转换路径跑通再落到抗混叠滤波器的参数设计与验证上。内容适合刚接触数字信号处理的工程师也适合那些常年用resample但说不清边界条件的老手。2. 采样率转换的频谱折叠原理抽取、插值在MATLAB里如何验证2.1 抽取的频谱扩张折叠进带内的信号来自哪一段抽取就是每M点保留一点x_down[n] x[Mn]。频域对应关系为X_d(e^{jω}) 1/M Σ X(e^{j(ω−2πk)/M})含义是原频谱被拉伸到M倍然后以2π为周期做叠加。当原始信号中最高频率成分超过fs/(2M)时超出部分就会折叠回低频区和有效信号叠在一起。这个折叠不是噪声是确定性的频谱搬移所以混叠具体出现在哪个频点是可计算的。用MATLAB验证这一段规则最直接。生成1 kHz与2.5 kHz两个正弦16 kHz采样M4直接抽取。抽取后采样率变成4 kHz奈奎斯特边界只有2 kHz2.5 kHz会折叠到4 kHz − 2.5 kHz 1.5 kHz而1 kHz的位置不变。fs0 16000; t 0:1/fs0:0.2; x sin(2*pi*1000*t) sin(2*pi*2500*t); M 4; xd x(1:M:end); % 直接抽取没有抗混叠滤波 Nfft 2048; f0 (0:Nfft/2-1)/Nfft*fs0; fd (0:Nfft/2-1)/Nfft*fs0/M; plot(f0, 20*log10(abs(fft(x, Nfft)))); hold on; plot(fd, 20*log10(abs(fft(xd, Nfft))), r); legend(原信号 16 kHz, 直接抽取 4 kHz); xlabel(Hz); ylabel(dB);这段代码里fft长度取了2048直接用正弦做频谱会有泄漏但峰值位置足够说明问题。运行后能看到原始2.5 kHz的谱峰在抽取后移动到1.5 kHz而1 kHz谱峰不变这就是没有抗混叠滤波时混叠发生的典型形态。反过来想如果有用信号恰好落在1.5 kHz那么来自2.5 kHz的混叠分量会直接叠加在它上面后续任何处理都无法再分开。2.2 插值形成的镜像频谱低通截止频率怎么定插值是在每两个原始样本之间插入L−1个零数学形式是x_up[n] x[n/L]当n是L的倍数时频域表现为频谱压缩到原来的1/L并且在原采样率fs的整数倍位置上复制出镜像副本。以16 kHz采样、1 kHz信号为例L4插值后采样率变为64 kHz镜像出现在1 kHz、15 kHz、17 kHz、31 kHz等位置其中15 kHz是−1 kHz的副本。L 4; xu upsample(x, L); % 中间插零不滤波 f_up (0:Nfft/2-1)/Nfft*fs0*L; plot(f_up, 20*log10(abs(fft(xu, Nfft)))); xlabel(Hz); ylabel(dB);插值后的目标是把原始带宽之外的镜像滤掉。低通截止频率取fs/(2L)2 kHz就能保留1 kHz成分并把15 kHz、17 kHz等副本全部压掉。这里的重点是插值本身不会损坏信号损坏发生在你忘记滤波的时候。真正决定重建质量的是滤波器通带平坦度和对第一镜像的抑制量。2.3 抗混叠滤波器的位置先滤波后抽取先插值后滤波把两节结论放在一起抽取之前必须低通把超过目标奈奎斯特边界的分量提前去掉插值之后必须低通把零插入产生的镜像清掉。滤波器放在这两处不是约定俗成而是频谱规则决定的。后面章节里出现的resample、upfirdn、firdecim等函数本质上都是在做“滤波抽取”或“插值滤波”的组合。理解了这个顺序才能看懂MATLAB函数内部为什么要先设计一个多相抗混叠滤波器。提示抽取的滤波和插值的滤波虽然都叫抗混叠但一个针对折叠一个针对镜像。两者对阻带起始频率的要求不同设计滤波器时不要套用同一组参数。3. 用MATLAB实现采样率转换resample、upfirdn与整数倍抽取/插值函数3.1 最快路径用resample完成有理数采样率转换MATLAB里最常被调用的多速率函数是resample一句话写法是y resample(x, p, q)把输入按p/q倍的速率重采样。p和q是互质整数目标采样率与源采样率之比等于p/q。例如48 kHz转44.1 kHz比值是0.91875用rat函数得到分子分母fs_in 48000; fs_out 44100; [p, q] rat(fs_out / fs_in); % p147, q160 y resample(x, p, q);resample内部采用Kaiser窗设计一个多相抗混叠滤波器然后依次做p倍插值、滤波、q倍抽取。对大多数原型验证场景这一条命令就够用不需要自己设计滤波器。resample有两个可选参数值得注意L控制滤波器长度L越大过渡带越窄beta控制Kaiser窗形状间接决定阻带衰减。默认值在多数场合表现尚可但如果你对混叠底噪有明确指标就应该自己设计h并改用3.2节的方式。注意resample输入输出长度会做自动对齐输出长度大致为ceil(length(x)*p/q)但首尾会有瞬态做符号同步或时延敏感处理时不能忽略。3.2 自定义滤波器用upfirdn把插值、滤波、抽取拆开当你需要精确控制抗混叠滤波器指标时resample的黑盒就不够用了。upfirdn(x, h, p, q)是更底层的原语先对x做p倍零插值与h卷积再做q倍抽取。三个操作在一条调用里完成但h需要你预先设计好。resample内部做的事情和它基本一致。p 147; q 160; fs_mid fs_in * p; % 插值后的中间采样率 7.056 MHz fp 20000; % 通带边界 fst 22500; % 阻带起始需覆盖输出Nyquist df [0 fp fst fs_mid/2] / (fs_mid/2); % 以中间采样率Nyquist做归一化 a [1 1 0 0]; % 通带1阻带0 dev [1e-4 1e-5]; % 通带纹波与阻带衰减线性幅度 [n, fo, ao, w] firpmord(df, a, dev); h firpm(n, fo, ao, w); y upfirdn(x, h, p, q);这里最容易被忽略的是h的频率轴。upfirdn中h工作在插值后的采样率fs_mid上滤波器设计必须换算到fs_mid的奈奎斯特频率。如果直接拿输入采样率的Nyquist做归一化实际滤波器通带会宽p倍镜像一个都滤不掉。这段代码里的fst取22500 Hz略高于输出奈奎斯特22050 Hz给过渡带留了少量余量工程上更稳妥。3.3 整数倍抽取与插值firdecim、firinterp和intfilt当转换比为整数时firdecim和firinterp是两个专用入口。firdecim(h, M, x)等价于先filter再downsample但内部对h做多相分解能显著减少计算量。firinterp同理做L倍插值后滤波。它们的典型调用如下M 4; h fir1(256, 1/M); % 截止在fs/(2M)的简单低通 yd firdecim(h, M, x); % 抽取输出采样率fs/M L 4; h2 fir1(256, 1/L); yu firinterp(h2, L, x); % 插值输出采样率fs*Lintfilt是更老的函数用于生成整数倍插值滤波器现在大部分新代码已经被firinterp或designMultirateFIR取代。新版本MATLAB里designMultirateFIR更值得关注它把抗混叠滤波器设计与抽取/插值参数统一成一个接口R2022b之后我一般在正式项目里优先用它。3.4 函数选型对比这四种方式别再混着用函数完成操作关键参数典型场景resample插值滤波抽取有理数转换p, q, L, beta快速原型、算法验证upfirdn插值滤波抽取可自定义hh, p, q需要指定滤波器指标firdecim滤波M倍抽取h, M固定整数倍抽取firinterpL倍插值滤波h, L固定整数倍插值designMultirateFIR设计抽取/插值滤波器filter, decim, dev, w新工程推荐很多人犯的错是用downsample直接做抽取x(1:M:end)得到的结果没有经过任何滤波高频混叠全部落在带内。对比一下误用和正确路径的底噪结论非常明显y_err x(1:4:end); % 直接抽取无滤波 y_ok resample(x, 1, 4); % 完整抗混叠链路 snr_err snr(y_err); % 通常只有20~30 dB snr_ok snr(y_ok); % 能到60 dB以上snr函数在Signal Processing Toolbox里只做指标对比可以接受。真正做音频或通信系统时还要查看混叠落入带内的具体频点而不仅仅是总SNR。4. 抗混叠滤波器设计窗函数、等波纹、多级级联与48 kHz到44.1 kHz实例4.1 滤波器指标怎么从系统需求推导出来设计抗混叠滤波器之前先定三个数通带边界fp、阻带起始频率fst、阻带衰减As。fp通常是信号本身的有效带宽fst取决于目标采样率下的奈奎斯特边界或镜像/混叠区最低频率As由系统底噪决定。以48 kHz到44.1 kHz音频转换为例信号内容只保留到20 kHz输出奈奎斯特是22.05 kHz所以fp20 kHz、fst22.05 kHz阻带衰减通常要求90 dB以上否则可听底噪里会混入镜像。设计方法优点缺点适用场景窗函数法 fir1简单、稳定、线性相位通带纹波与阻带衰减互相牵扯快速滤波、指标不苛刻等波纹 firpm相同阶数阻带抑制更均衡阶数估计不准时需要迭代抽取/插值混叠抑制最小二乘 firls通带纹波总体最低阻带可能出现个别超差对总失真敏感的音频链路4.2 等波纹设计firpmord配合firpm的完整流程等波纹法的优势是阻带内抑制均匀不会出现窗函数法那种靠近截止频率处凹陷更深、远处回升的情况。对混叠抑制来说所需的是整个阻带都低于某个阈值等波纹法更贴合这个需求。用firpmord先估阶数再交给firpm设计是标准流程。% 中间采样率 7.056 MHz归一化以插值后Nyquist为1 fs_mid 48000 * 147; % 7.056 MHz f_norm [0 20000 22050 fs_mid/2] / (fs_mid/2); a [1 1 0 0]; dev [0.0005 1e-5]; % 通带纹波约0.0087 dB阻带100 dB [n, fo, ao, w] firpmord(f_norm, a, dev); h firpm(n, fo, ao, w);这段代码会算出非常高的阶数因为过渡带只有2050 Hz而中间采样率高达7.056 MHz相对过渡带宽度不到0.15%。这是多速率设计里的经典两难转换比越大镜像频点越密集对滤波器的要求越高。实际工程中不会让过渡带卡这么紧一般会把fst放宽到24 kHz甚至更高让滤波器阶数从数千降到几百代价是损失一点高频带宽。4.3 多级级联什么时候拆成两级更划算当抽取因子D超过10并且过渡带相对窄时单级滤波器阶数会爆炸。把D拆成D1×D2先按D1做第一次抽取再用D2做第二次每一级的过渡带都可以放宽总滤波器长度通常只有单级的几分之一。经验法则每级抽取因子尽量控制在10以内过渡带宽度至少占该级奈奎斯特带宽的20%到50%否则级联的优势不明显。% 单级抽取 D12 h12 designMultirateFIR(1, 12, 64, 90); y1 upfirdn(x, h12, 1, 12); % 两级抽取 3×4 h3 designMultirateFIR(1, 3, 64, 90); h4 designMultirateFIR(1, 4, 64, 90); y2 upfirdn(upfirdn(x, h3, 1, 3), h4, 1, 4);4.4 完整示例48 kHz音频转44.1 kHz的自定义滤波器链路把前面的内容串成一个可运行脚本。信号源用20 kHz以内的多音合成方便在频谱上直接看镜像是否被压住fs_in 48000; t 0:1/fs_in:1; x sin(2*pi*1000*t) 0.5*sin(2*pi*10000*t) 0.2*sin(2*pi*19900*t); p 147; q 160; fs_mid fs_in * p; f_norm [0 20000 22050 fs_mid/2] / (fs_mid/2); a [1 1 0 0]; dev [0.0005 1e-5]; [n, fo, ao, w] firpmord(f_norm, a, dev); n n rem(n, 2); % 保证线性相位为偶数阶 h firpm(n, fo, ao, w); y upfirdn(x, h, p, q); y y(1:round(length(x)*p/q)); % 截断到预期输出长度 fs_out fs_in * p / q;脚本中把滤波器阶数向上取偶是为了让线性相位群延迟正好是整数样本点后面做时延补偿时会省事。运行后输出采样率是44100 Hz频谱图在20 kHz以上应该看不到超过−90 dB的镜像残留。如果看到明显的20 kHz以上谱线先查f_norm是不是按fs_mid归一化这是最容易出错的点。5. 多速率系统的多相分解与验证群延迟补偿和混叠残留定位5.1 把滤波器拆开多相分解到底省在哪里抗混叠滤波器通常动辄上百阶直接在整个采样率上做卷积计算量很高。多相分解的思路是把长度N的FIR滤波器按抽取因子M拆成M个子滤波器每个子滤波器只在低速率的某一路并行工作。数学形式是H(z) Σ z^−k E_k(z^M)其中E_k(z^M)对应第k个子多相分量。MATLAB的firdecim和resample内部都是这种结构。手动验证多相结构的等价性可以用一段很短的MATLAB代码M 4; h fir1(96, 1/M); hp reshape([h(:); zeros(M - mod(length(h), M), 1)], length(h)M-mod(length(h),M), []); % hp的每一列就是一个多相子滤波器hp的每一列就是一路低速子滤波器输入序列按相位分配到各路卷积后再合并。子滤波器长度大致是主滤波器的1/M每路乘法次数大幅减少。这个结构在多速率系统里几乎是标配理解它才能理解为什么resample处理长序列时比“先filter后downsample”快那么多。5.2 群延迟怎么补线性相位滤波器的延迟是N/2线性相位FIR滤波器的群延迟是(N−1)/2个样本以滤波器工作采样率计。使用upfirdn时这个延迟折算到输出采样率要除以q。不补偿的话重采样后的信号整体偏移几十个样本做符号同步或者数据对齐时会直接错位。delay_mid (length(h) - 1) / 2; % 中间采样率下的延迟 delay_out delay_mid / q; % 换算到输出采样率 y upfirdn(x, h, p, q); y y(round(delay_out) 1 : end); % 去掉瞬态并补偿群延迟这段代码放在设计流程最后比用xcorr估计延迟更精确因为群延迟是已知的不需要用信号相关性去猜。复数调制信号或最小相位滤波器不适用这条规律最小相位滤波器的群延迟需要用grpdelay逐点计算。5.3 用pwelch和往返测试定位混叠残留验证抗混叠效果不能只看时域波形把转换前后的功率谱放在同一张图里对比是标准做法。pwelch返回的功率谱密度可以直接看到残留镜像落在哪个频段[pxx, f] pwelch(y, hann(4096), 2048, 8192, fs_out); [pxo, fo] pwelch(x, hann(4096), 2048, 8192, fs_in); semilogy(fo, sqrt(pxo), k); hold on; semilogy(f, sqrt(pxx), r); xline(22050/44100*2, --, Nyquist);更严格的做法是往返测试把信号从48 kHz转到44.1 kHz再转回48 kHz与原始信号比对。往返转换会同时暴露出混叠、镜像和滤波器非理想过渡带三类问题。若往返后的SNR低于系统指标把两个方向的滤波器阻带衰减各提高10 dB再测一次通常能找到问题出在哪一级。多速率系统的交付验证多数时候就靠这两步一张频谱图一个往返SNR。本文还有配套的精品资源点击获取