1.软件版本
MATLAB2013b
2.本算法理论知识
Volterra级数即含记忆的泰勒级数，其数学公式与泰勒级数及其相似，Volterra 级数模型的输出信号，其计算公式是通过输入信号的幂次方来表达的。Volterra级数与泰勒级数不同之处在于Volterra级数具有延迟功能，因此Volterra级数更适合应用在具有记忆效应的功率放大器线性化处理过程中。因此其数学表达式如下：
变量yk(*)的表达式为：
当输入的信号是复数信号的时候，yk公式可以改写为如下表达式：
由于信号的带宽随着非线性阶数的增加而不断增加，当输入的信号通过Volterra模型时，信号的带宽将随着非线性器件阶数的增加而不断增加。以三阶系统为例，其信号带宽将扩大三倍，其最终的输出信号带宽取决于非线性系统的最高阶数。这就导致了一个宽带系统的输出信号在被捕获之前会被滤波器滤除，在这种情况下，单纯的使用Volterra级数模型，系统的输出信号带宽和输入信号带宽将会产生较大的差异性，并最终影响模型的精度。
整个带限Volterra级数的数字预失真模型可以通过如下表达式表示：
输出信号y(n)可以表示为：
3.核心代码
clc;
clear;
close all;
warning off;
addpath 'func\'
%%
%信号产生高旁瓣频谱的信号
load Rs.mat
L = 2^13;%根据电脑速度选择信号长度
fs = 1.8432*10^8;
fc = 7*10^7;%载波70MHZ
ts = 1/fs;
Scal= 8;
Ns = 1;
%original input
Xn = signal(1:L)/max(abs(signal(1:L)))/2;
Xn0 = Xn;
m = length(Xn);
figure;
[psdx,freq] = func_psd(Xn,m,ts,Scal);
plot(freq/1e6,psdx,'g','linewidth',2);
grid on
hold on
axis([-fc/1e6,fc/1e6,-110,0]);
%%
%论文DPD
K = 127;
Wn = [0.01,0.5];%修改0.1的值，获得不同情况下的band limit效果,0.0001~0.1
w = fir1(K,Wn,'bandpass');
K2 = 0;
w2 = [1,1];
%计算C_Lx1,Volterra kernel of the system
G_BL=[ 1.0513+j*0.0904,-0.0542-j*0.2900,-0.9657-j*0.7028,...
-0.0680-j*0.0023, 0.2234+j*0.2317,-0.2451-j*0.3735,...
0.0289-j*0.0054,-0.0621-j*0.0932, 0.1229+j*0.1508];
G_MP=[1.0513+j*0.0904,-0.0542-j*0.2900,-0.9657-j*0.7028,...
-0.0680-j*0.0023, 0.2234+j*0.2317,-0.2451-j*0.3735,...
0.0289-j*0.0054,-0.0621-j*0.0932, 0.1229+j*0.1508];
%计算U_Nx1 = X_NxL * C_Lx1
U_Nx1 = func_volterra0(Xn,w,G_BL,K,Ns);
[psdu,freq] = func_psd(U_Nx1,m,ts,Scal);
plot(freq/1e6,psdu,'b','linewidth',2);
grid on
hold on
axis([-fc/1e6,fc/1e6,-110,0]);
%is the expected inverse output matrix generated from the PA input (the output of the predistorter) u~,
%u~根据公式22计算得到。
%1st:LS
U = func_volterra_Matrix(Xn,w2,G_BL,G_MP,K2,Ns);
%构造Y,LS estimate
Out = func_volterra1(U,w2,G_BL,K2,Ns);
%公式26
ya = (abs(Out.^0)).*Out;
yb = (abs(Out.^2)).*Out;
yc = (abs(Out.^4)).*Out;
Y1 = [ya(3:m);yb(3:m);yc(3:m);ya(2:m-1);yb(2:m-1);yc(2:m-1);ya(1:m-2);yb(1:m-2);yc(1:m-2)];
Y2 = conj(Y1');
Cest0= pinv(Y2)*U(3:m).';
Cest = Cest0;
%LMS算法的参数
step = 0.002;
Xn2 = Xn;
for ii = 1:5
%变步长LMS estimate
U = func_volterra_Matrix(Xn2,w2,G_BL,G_MP,K2,Ns);
Out = func_volterra1(U,w2,G_BL,K2,Ns);
ya = (abs(Out.^0)).*Out;
yb = (abs(Out.^2)).*Out;
yc = (abs(Out.^4)).*Out;
Y1 = [ya(3:m);yb(3:m);yc(3:m);ya(2:m-1);yb(2:m-1);yc(2:m-1);ya(1:m-2);yb(1:m-2);yc(1:m-2)];
Y2 = conj(Y1');
err = (Xn2(3:end) - [Y2*Cest].')/25;
Cest = Cest - [step*err*Y2].';
Xn2 = func_volterra1(U_Nx1,w2,Cest,K2,Ns);
end
Yn2 = func_volterra1(U_Nx1,w2,Cest0,K2,Ns);
%通过功放
U_Nx2 = func_volterra0(Yn2,w,G_BL,K,Ns);
[psdy,freq] = func_psd(U_Nx2,m,ts,Scal);
plot(freq/1e6,psdy,'r','linewidth',2);
grid on
hold on
axis([-fc/1e6,fc/1e6,-110,0]);
legend('标准信号','无预失真时功放频谱','有预失真时功放频谱');
xlabel('Frequency Offset(MHz)');
ylabel('Normalized power Spectral Density (dB)');
save R2.mat freq psdy
4.操作步骤与仿真结论
5.参考文献
[01]P. Singerl, H. Koeppl. A low-rate identification method for digital pre distorters based on Volterra kernel interpolation. The 48th Midwest Symposium on Circuits and Systems. 2005, 2: 1533~1536.
[02]管鲍. 基于FPGA的DPD设计与实现[D]. 武汉: 武汉邮电科学研究院, 2014.
[03]Hao Lin, Taijun Liu, Yan Ye. Optimization Design of FPGA-Based Look-up-Tables or Linearizing RF Power Amplifiers[J]. 2011, 03(21): 8-11.
A01-145
