IIR+全通濾波器實現相位平衡_matlab仿真

http://blog.csdn.net/shenziheng1/article/details/53409743

1.全通系統實現相位均衡

本文從逆向角度出發，基於搜索最優化策略對全通系統相位平衡問題進行了初步研究。主要思想：在頻域中對它在單位圓上抽樣N個點，按照所要求的頻率特性，找出一個全通濾波器，使其滿足一定的相位特性。最后利用最優化（例如遺傳算法）策略進行迭代，找到滿足要求的最優解。

適用於低通系統相位補償的二階全通濾波系統如下：

2.設計滿足特定相位條件的全通濾波器

2.1 產生實驗信號，組份0.01 0.02（保留） 0.6、0.7（低通濾波濾掉）

2.2 設計低通濾波器，進行濾波，並觀察濾波后的信號是否發生了相位失真

低通濾波器作用下，高頻成分被過濾掉，僅剩下低頻成分0.01、0.02。

我們再看看濾波后，幅值和相位發生了什么變化？

我們發現，幅度圖像有點變形！轉過來看一下相位圖像！！！問題出現了，發生了相位失真！就是說，經過低通濾波時0.01和0.02兩個頻率成分“相對步調”不一致，這就是群延時的意義，說明該濾波器的群延時不是常數。

3.構造特定全通濾波器相位，並利用遺傳算法進行參數搜索設計

經過全通濾波器后，發現相位失真被補償回來了！！！0.01/0.02步調又一致了。

最后，我們分析一下，相位補償過程中，我們做的工作：

這張圖才是相位均衡的核心，黑絲線條是低通巴特沃茲濾波器引進的非線性相位；藍色曲線是我們利用全通濾波器構造的非線性相位，兩者通過系統級聯完成了線性相位的設計！

 

 

function [ delt ] = myfit( X  )
%X(1)  a
%X(2)  b
global AP_Phase
[h, w] = freqz([X(2),X(1),1],[1,X(1),X(2)]);
AP_Phase_temp =unwrap(angle(h))*180/pi;
delt = sum(abs(AP_Phase-AP_Phase_temp));  
end

　　

 

clc;clear all;close all;  
%%  生成模擬模擬信號  
% Fs = ;   
% t = ;   
% %歸一化頻率設計0.01 0.02 0.06  0.07   
% y = ;    
% yreal = ;  
%% 設計一個低通濾波器*********************************************%%  
fp=15/500;%通帶截止頻率 0.03  
fs=40/500;%阻帶截止頻率 0.06  
wp=3;     %通帶最大衰減   
ws=30;    %阻帶最大衰減   
[n,fc]=buttord(fp,fs,wp,ws);   
[b1,a1]=butter(n,fc);    
[h1,w]=freqz(b1,a1);  
LP_phase = unwrap(angle(h1))*180/pi;  %低通濾波器的相位函數  
%% 設計全通濾波器 進行相位補償 ************************************%%  
Ideal_Phase = zeros(512,1);  
Y0 = min(min(LP_phase));  
for i=1:1:512  
    Ideal_Phase(i,1) = (Y0*i)/512;  
end  
%全通濾波器需要補償的相位  
AP_Phase = zeros(512,1);  
for j=1:1:512   
    AP_Phase(j,1) = Ideal_Phase(j,1) - LP_phase(j,1);  
end  
global AP_Phase
%% 利用遺傳算法優化設計二階全通濾波器 ******************************%%  
%濾波器原型為：H(z)=（b+a*z^(-1)+z^(-2)）/(1+a*z^(-1)+b*z^(-2));  
%目標函數：minimum -> delt = sum(abs(AP_Phase-AP_Phase_temp));  
%動態范圍：a->[-100,100]  b->[-100,100]  
options = gaoptimset(); %遺傳算法
options.Generations = 200000; %最大迭代數設為20000  
UP = [50,50]; LOW = [-50,-50];%變量上下限  
[X,FVAL] = ga(@myfit, 2 ,[], [],[],[],LOW,UP,[],options);  
[h3,w]=freqz([X(2),X(1),1],[1,X(1),X(2)]);  
AP_Phase_Lab = unwrap(angle(h3))*180/pi;  %低通濾波器的相位函數  
%% 相位校正后的效果 **********************************************%%  
figure;
subplot(3,1,1)
plot(LP_phase, 'k')
hold on;
plot(AP_Phase);
hold on;
plot(Ideal_Phase, 'r')
subplot(3,1,2)
plot(AP_Phase_Lab,'b')
subplot(3,1,3)
plot(LP_phase+AP_Phase_Lab);