1. 項目緣起當“固定”的濾波器遇上“變化”的世界在數字信號處理的世界里我們常常面臨一個經典的矛盾設計的靈活性與實現的效率。就拿一個最常見的需求來說——你需要一個截止頻率可變的低通濾波器。最直接的想法是什么沒錯每次頻率參數一變我就重新計算一遍濾波器的系數然后更新到我的濾波器結構中。這在離線處理或者對實時性要求不高的場合或許可行但一旦放到FPGA、DSP或者高性能嵌入式系統中頻繁地重新計算并加載一組可能多達幾十甚至上百個的濾波器系數比如一個高階FIR濾波器帶來的計算開銷和延遲往往是不可接受的。尤其是在軟件無線電、雷達信號處理、音頻效果器這類對實時性要求極高的領域參數需要連續、平滑地調整這種“硬切換”系數的方式會引入信號的不連續產生可聞的咔嗒聲或影響系統性能。那么有沒有一種方法能讓濾波器的某個關鍵參數比如截止頻率、分數延遲像擰旋鈕一樣連續可調而無需在每次調整時都大動干戈地重新計算所有系數呢這就是Farrow結構濾波器要解決的核心問題。它本質上是一種實現可變分數延遲濾波器或可變參數濾波器的高效結構。我第一次在項目中接觸它是在為一個多速率采樣系統設計一個采樣率轉換模塊時傳統的多項式插值方法在精度和復雜度上難以平衡而Farrow結構提供了一種優雅的折中方案。今天我們就來徹底拆解這個在工程上極具魅力的結構從它要解決什么問題到它的核心原理再到如何一步步設計并實現它最后分享一些在硬件實現中容易踩的坑。2. Farrow結構的核心思想將參數變化“固化”到結構里要理解Farrow結構我們得先忘掉那些復雜的公式從一個更直觀的視角來看。想象一下一個濾波器的輸出是其系數與輸入信號卷積的結果。如果這個濾波器的某個特性如群延遲需要連續變化傳統做法是讓系數本身成為這個變化參數的函數。也就是說每一個濾波器系數h[k]不再是固定值而是一個關于可變參數d通常代表分數延遲范圍在0到1之間的函數h[k] f_k(d)。Farrow結構的巧妙之處在于它對這個函數f_k(d)做了一個關鍵的假設和簡化每個系數關于可變參數d的變化可以用一個低階多項式來近似。也就是說h[k] ≈ c_{k,0} c_{k,1} * d c_{k,2} * d^2 ... c_{k,L} * d^L這里L是多項式的階數c_{k,l}就是我們需要預先計算并固定下來的“子濾波器”系數。看到了嗎變化的部分d及其冪次被抽離出來了而需要存儲和參與實時卷積運算的變成了固定不變的系數c_{k,l}。這個思想帶來了巨大的優勢實時性當需要改變延遲d時我們不再需要重新計算或加載任何濾波器系數c_{k,l}。我們只需要根據新的d值實時計算出一組“權重”(1, d, d^2, ..., d^L)然后用這組權重去組合那些固定的子濾波器的輸出。結構化與模塊化整個濾波器可以被分解為(L1)個并行的、系數固定的子濾波器每個對應多項式的一項以及一個多項式求值即權重組合模塊。這種結構非常規整特別適合用硬件如FPGA進行并行流水線實現。設計靈活性我們可以通過設計多項式階數L和子濾波器系數c_{k,l}來權衡逼近精度、計算復雜度和濾波器性能。那么這個結構具體長什么樣呢一個典型的、用于實現分數延遲的Farrow結構框圖如下所示我們以三次多項式即L3為例輸入 x[n] | |---- 固定子濾波器 H0(z) (系數為 c_{k,0}) ---- 乘 1 (即 d^0) ---\ |---- 固定子濾波器 H1(z) (系數為 c_{k,1}) ---- 乘 d ----- 加法器 ---- 輸出 y[n] |---- 固定子濾波器 H2(z) (系數為 c_{k,2}) ---- 乘 d^2 -----/ ---- 固定子濾波器 H3(z) (系數為 c_{k,3}) ---- 乘 d^3 -----/關鍵點H0(z),H1(z), ...,H3(z)都是普通的、系數固定的FIR濾波器。它們的輸出分別乘以d^0,d^1,d^2,d^3然后求和就得到了最終具有分數延遲d的輸出y[n]。d可以在每個采樣時刻動態改變而所有子濾波器的系數是焊死在硬件或代碼里的。3. 從零開始Farrow濾波器系數設計方法論理解了思想下一步就是如何得到那些固定的子濾波器系數c_{k,l}。這是Farrow濾波器設計的核心。設計目標通常是讓整個濾波器在感興趣的頻帶內其頻率響應盡可能地逼近一個理想的、延遲為d的分數延遲器。理想分數延遲器的頻率響應是H_ideal(e^{jω}) e^{-jωd}它具有線性相位群延遲恒為d。設計方法主要有兩大類基于最小二乘誤差準則和基于拉格朗日插值。這里我重點講工程上最常用、也相對直觀的最小二乘設計法。3.1 設計問題建模假設我們設計一個長度為N階數為N-1的FIR濾波器用于近似分數延遲d。其傳遞函數為H(z, d) Σ_{k0}^{N-1} h[k] * z^{-k}而根據Farrow結構我們假設h[k]是d的L階多項式h[k] Σ_{l0}^{L} c_{k,l} * d^l我們的目標是在d ∈ [0, 1]和ω ∈ [0, απ]α是帶寬因子比如0.8表示利用80%的奈奎斯特帶寬的范圍內讓實際響應H(e^{jω}, d)逼近理想響應e^{-jωd}。定義一個誤差函數E(ω, d) H(e^{jω}, d) - e^{-jωd}設計目標就是找到一組系數c_{k,l}使得在定義的(ω, d)區域上誤差E的某種范數如平方誤差的積分最小。這是一個標準的優化問題。3.2 具體設計步驟與MATLAB示例雖然推導過程涉及積分和矩陣運算但幸運的是我們可以借助MATLAB等工具來高效完成。下面是一個手把手的步驟步驟一確定設計參數N: 主FIR濾波器長度抽頭數。越大逼近精度越高計算量也越大。L: 多項式階數。越高對參數d變化的擬合能力越強但結構也更復雜。通常L3三次是一個很好的折中。alpha: 歸一化帶寬0到1之間。我們只關心這個帶寬內的性能。設計網格在d和ω的二維平面上需要劃分密集的網格點來進行離散化優化。例如d從0到1步進0.01ω從0到alpha*pi步進pi/500。步驟二構建最小二乘問題對于每一個網格點(ω_i, d_j)我們可以寫出H(e^{jω_i}, d_j) Σ_{k0}^{N-1} (Σ_{l0}^{L} c_{k,l} * d_j^l) * e^{-jω_i k}令C是一個將所有系數c_{k,l}按特定順序例如先按k再按l排列成的列向量。那么對于所有網格點我們可以建立一個巨大的線性方程組A * C ≈ B其中A矩陣的每一行對應一個(ω_i, d_j)點其元素由d_j^l * e^{-jω_i k}構成B向量對應每個點的理想響應值e^{-jω_i d_j}。步驟三求解系數這是一個超定線性方程組我們用最小二乘法求解C (A^H * A) \ (A^H * B)這里^H表示共軛轉置\是MATLAB中的左除運算符用于求解最小二乘解。步驟四驗證與評估求解出系數C后將其重新排列成N x (L1)的矩陣C_mat其中第k行第l列就是c_{k,l}。 然后我們可以固定一個d值如0.5用C_mat生成對應的FIR系數h C_mat * [1; d; d^2; ... d^L]繪制其頻率響應并與理想延遲e^{-jωd}比較。掃描d從0到1觀察濾波器幅頻響應應接近1和群延遲應接近d的變化是否平滑。實操心得在MATLAB中構建矩陣A時一定要注意索引和維度的對應關系這是最容易出錯的地方。一個技巧是先用循環寫一個清晰但低效的版本確保邏輯正確后再嘗試用向量化操作meshgrid,kron等進行加速。另外帶寬因子alpha不要設得太大如0.95以上因為逼近理想延遲器在頻帶邊緣非常困難強求會導致通帶內紋波增大。通常0.8~0.9是穩健的選擇。4. 超越分數延遲Farrow結構的變體與應用擴展雖然Farrow結構最初是為可變分數延遲而生但其“用固定子濾波器組合實現參數可變”的思想可以被推廣到更廣泛的可變參數濾波器設計中。這才是它真正強大和有趣的地方。4.1 可變截止頻率濾波器這是另一個非常實用的場景。假設我們需要一個截止頻率fc可變的低通濾波器。理想情況下濾波器的系數應該是fc的函數。我們可以借鑒Farrow思想將每個系數h[k]用關于歸一化截止頻率ff fc / fsfs為采樣率的多項式來近似h[k](f) ≈ Σ_{l0}^{L} c_{k,l} * f^l設計過程與分數延遲器類似但目標函數變了。此時我們需要讓設計出的濾波器在通帶[0, f]內響應接近1在阻帶[fΔ, 0.5]內響應接近0Δ是過渡帶。這同樣可以轉化為一個在(ω, f)二維區域上的最小二乘優化問題只是誤差權重函數需要精心設計例如在通帶和阻帶賦予高權重在過渡帶賦予低權重。實現結構和經典Farrow結構一模一樣只是輸入參數從延遲d變成了截止頻率f子濾波器系數c_{k,l}是針對截止頻率變化而優化的。4.2 應用于采樣率轉換多項式插值這是Farrow結構最早、也是最成功的應用之一。在異步采樣率轉換中我們需要計算輸入序列在非整數采樣點上的值這本質上就是一個分數延遲問題。例如常用的三次拉格朗日插值器其系數就是關于分數間隔μ相當于d的三次多項式。你可以驗證這些多項式系數正好可以排列成一個4x4的C_mat矩陣從而完美地用Farrow結構實現一個高效的三次插值器。更一般地任何基于多項式的插值核如分段拋物線、樣條插值都可以用Farrow結構來實現。這使得它在數字上下變頻、軟件無線電接收機中成為了標準模塊。4.3 其他可變參數理論上只要濾波器性能是某個參數p的平滑函數并且可以用多項式較好地近似就可以嘗試Farrow結構。例如可變帶寬的帶通濾波器中心頻率固定帶寬可變。可變Q值的諧振器諧振頻率固定品質因數Q可變。可變滾降因子的升余弦濾波器。設計經驗在將這些擴展時最大的挑戰在于多項式近似的有效性。如果濾波器系數隨參數p的變化非常劇烈或非線性那么可能需要很高的多項式階數L才能較好近似這會急劇增加子濾波器的數量L1個和計算量。因此在決定采用Farrow結構前一定要先分析系數隨參數變化的曲線是否相對平滑。通常在參數變化范圍較小如d在0~1之間時Farrow結構表現最佳。5. 硬件實現考量與實戰中的“坑”將Farrow濾波器從算法模型搬到FPGA或ASIC上是另一個充滿細節的戰場。這里分享幾個我趟過的雷區。5.1 計算復雜度與資源優化一個N抽頭、L階多項式的Farrow濾波器需要(L1)個并行的N抽頭FIR子濾波器。直接實現的乘法器數量是N*(L1)這看起來很大。但我們可以利用結構特點進行優化子濾波器合并計算觀察結構每個輸入采樣x[n]需要與所有子濾波器的第一抽頭系數c_{0,0}, c_{0,1}, ..., c_{0,L}相乘然后分別延遲、再與下一組系數相乘。我們可以將同一抽頭位置k上的、屬于不同子濾波器的系數c_{k,0} ... c_{k,L}視為一個向量。當計算該抽頭的貢獻時我們實際上是在計算這個系數向量與權重向量[1, d, d^2, ..., d^L]^T的點積。這可以轉化為一個先乘累加、再乘以輸入信號x[n-k]的過程有時能減少乘法器數量。多項式求值優化計算權重d^l可以使用霍納法則將Σ c_{k,l} * d^l的計算轉化為y_k c_{k,0} d*(c_{k,1} d*(c_{k,2} ... d*c_{k,L})...)這樣只需要L次乘法和L次加法而不是L次冪運算和L次乘法。系數對稱性利用如果目標頻率響應是對稱的如線性相位那么設計出的子濾波器系數c_{k,l}也可能呈現某種對稱性。例如對于可變分數延遲濾波器其沖激響應關于中心點近似對稱。利用這種對稱性可以將乘法器數量幾乎減半。5.2 動態參數更新的時序問題參數d或f是動態更新的。在硬件中這需要仔細處理更新速率d更新的速度不能超過數據處理流水線的“吞吐量”。如果d每個時鐘周期都變那么權重d^l需要每個周期重新計算并且必須確保在用到新權重的時刻數據流水線中對應的是新的輸入數據。通常d的更新速率遠低于采樣率。同步新的d值必須與輸入數據流正確同步。一個常見的做法是使用一個“參數有效”信號當d更新時該信號拉高一個周期標志著從此之后進入濾波器的數據將使用新的d值。這需要在數據路徑上插入相應的控制邏輯。5.3 有限字長效應與精度管理這是硬件實現中最容易出問題的地方尤其是在定點設計中。系數量化設計得到的c_{k,l}通常是高精度浮點數。我們需要將其量化為定點數如Q格式。量化會引入誤差可能導致頻率響應偏離設計目標甚至不穩定。必須進行充分的仿真掃描不同的量化位寬如12位、16位、18位觀察通帶紋波、阻帶衰減等關鍵指標的變化找到滿足性能要求的最小位寬。中間結果位寬擴展在計算d^l以及子濾波器乘累加的過程中中間結果的動態范圍會擴大。例如計算d^2假設d是16位小數可能需要32位來保證精度不損失。在加法樹中位寬擴展更明顯。必須為每條數據路徑仔細規劃位寬防止溢出和精度過度損失。一個安全的方法是先做全精度仿真記錄中間結果的最大最小值再據此確定定點位寬。權重計算的非線性d^l的計算尤其是高次冪對d的量化誤差非常敏感。當d接近0或1時d^l的值可能非常小量化誤差會占據主導。可以考慮采用查找表LUT來存儲d^l的值或者使用分段線性近似等方法來計算高次冪以平衡精度和資源消耗。踩坑實錄在一次通信接收機項目中我們使用Farrow結構做定時恢復。仿真時性能完美但上板后誤碼率總是差一點。用邏輯分析儀抓取中間數據發現當分數延遲d在0.1以下時輸出信號有明顯的失真。最終定位到問題計算d^3的模塊為了節省乘法器我們采用了(d*d)*d的順序計算。由于d很小d*d的結果在定點量化后直接變成了0導致三次項完全失效。教訓是對于可能接近零的小數運算必須保證中間結果的精度或者改變計算順序如先計算高精度浮點再量化或者采用查找表。后來我們改用一個小型的、針對d的LUT來存儲d^2和d^3的值問題得以解決。6. 性能評估與設計實例一個可變延遲線讓我們用一個具體的設計實例把前面所有的點串起來。目標設計一個用于音頻處理的可變分數延遲線要求延遲d在0到1個采樣周期內連續可調音頻帶寬為20kHz采樣率fs48kHz因此alpha ≈ 0.833。設計參數選擇主濾波器長度N32。這是一個折中能提供較好的阻帶抑制。多項式階數L3。三次多項式足以在0-1區間平滑近似。帶寬因子alpha0.85略高于需求留有余量。采用最小二乘準則在d[0:0.02:1]和ω[0:π/1000:0.85π]的網格上設計。設計結果 使用MATLAB的farrow設計函數或按前述步驟自編代碼得到32x4的系數矩陣C_mat。性能評估固定延遲響應取d0.5生成系數h C_mat * [1; 0.5; 0.25; 0.125]。繪制其頻率響應。實測在0-20kHz通帶內幅度起伏 0.01 dB群延遲波動 0.001個樣本非常接近理想的0.5樣本延遲。動態掃描讓d從0線性變化到1觀察濾波器幅頻響應。在整個過程中通帶內幅度保持平坦阻帶20kHz衰減始終大于80dB。群延遲值緊密跟隨d值變化誤差在±0.005個樣本以內。聽感測試用一段正弦波掃頻信號通過該延遲線同時用低頻三角波調制d在0-1之間變化。輸出信號聽起來應該是純凈的、音高無變化因為只有延遲無頻移但帶有周期性相位調制的效果。如果設計不良可能會引入可聞的諧波失真或振幅調制。本例中聽感干凈。硬件實現概要FPGA系數存儲32*4128個系數每個量化為18位有符號定點數存儲在ROM中。數據處理路徑輸入音頻數據x[n]為24位。設計4個并行的32抽頭FIR濾波器H0-H3。利用系數的近似對稱性每個FIR采用轉置結構乘法器數量減半約需(32/2)*4 64個乘法器。參數d為16位無符號小數0x0000~0xFFFF 對應 0~1。權重計算模塊使用一個時鐘周期計算d2 d*d下一個周期計算d3 d2*d均為32位中間結果。同時使用一個三級流水線的霍納法則單元計算每個抽頭k的加權和sum_k c_k0 d*(c_k1 d*(c_k2 d*c_k3))。最終輸出y[n]為24位由所有sum_k * x[n-k]的累加和得到。資源與性能在Xilinx Artix-7 FPGA上占用約70個DSP slice最高運行時鐘頻率可達150MHz遠高于48kHz的音頻采樣率有充足資源進行多通道處理。通過這個實例可以看到一個精心設計的Farrow濾波器能夠在參數連續變化時保持卓越且穩定的性能同時滿足實時硬件處理的要求。它不再是教科書上一個抽象的數學結構而是一個能解決實際工程難題的強力工具。