LMS自适应滤波算法的C++实现

头文件：

/** Copyright (c) 2008-2011 Zhang Ming (M. Zhang), zmjerry@163.com** This program is free software; you can redistribute it and/or modify it* under the terms of the GNU General Public License as published by the* Free Software Foundation, either version 2 or any later version.** Redistribution and use in source and binary forms, with or without* modification, are permitted provided that the following conditions are met:** 1. Redistributions of source code must retain the above copyright notice,*    this list of conditions and the following disclaimer.** 2. Redistributions in binary form must reproduce the above copyright*    notice, this list of conditions and the following disclaimer in the*    documentation and/or other materials provided with the distribution.** This program is distributed in the hope that it will be useful, but WITHOUT* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for* more details. A copy of the GNU General Public License is available at:* http://www.fsf.org/licensing/licenses*//******************************************************************************                                   lms.h** Least Mean Square Adaptive Filter.** Least mean squares (LMS) algorithms are a class of adaptive filter used* to mimic a desired filter by finding the filter coefficients that relate* to producing the least mean squares of the error signal (difference* between the desired and the actual signal). It is a stochastic gradient* descent method in that the filter is only adapted based on the error at* the current time.** This file implement three types of the LMS algorithm: conventional LMS,* algorithm, LMS-Newton algorhm and normalized LMS algorithm.** Zhang Ming, 2010-10, Xi'an Jiaotong University.*****************************************************************************/#ifndef LMS_H#define LMS_H#include #include namespace splab{templateType lms( const Type&, const Type&, Vector&, const Type& );templateType lmsNewton( const Type&, const Type&, Vector&,const Type&, const Type&, const Type& );templateType lmsNormalize( const Type&, const Type&, Vector&,const Type&, const Type& );#include }// namespace splab#endif// LMS_H

实现文件：

/** Copyright (c) 2008-2011 Zhang Ming (M. Zhang), zmjerry@163.com** This program is free software; you can redistribute it and/or modify it* under the terms of the GNU General Public License as published by the* Free Software Foundation, either version 2 or any later version.** Redistribution and use in source and binary forms, with or without* modification, are permitted provided that the following conditions are met:** 1. Redistributions of source code must retain the above copyright notice,*    this list of conditions and the following disclaimer.** 2. Redistributions in binary form must reproduce the above copyright*    notice, this list of conditions and the following disclaimer in the*    documentation and/or other materials provided with the distribution.** This program is distributed in the hope that it will be useful, but WITHOUT* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for* more details. A copy of the GNU General Public License is available at:* http://www.fsf.org/licensing/licenses*//******************************************************************************                                lms-impl.h** Implementation for LMS Adaptive Filter.** Zhang Ming, 2010-10, Xi'an Jiaotong University.*****************************************************************************//*** The conventional LMS algorithm, which is sensitive to the scaling of its* input "xn". The filter order p = wn.size(), where "wn" is the Weight* Vector, and "mu" is the iterative setp size, for stability "mu" should* belong to (0, Rr[Rxx]).*/template Type lms( const Type &xk, const Type &dk, Vector &wn, const Type &mu ){int filterLen = wn.size();static Vector xn(filterLen);// update input signalfor( int i=filterLen; i>1; --i )xn(i) = xn(i-1);xn(1) = xk;// get the outputType yk = dotProd( wn, xn );// update the Weight Vectorwn += 2*mu*(dk-yk) * xn;return yk;}/*** The LMS-Newton is a variant of the LMS algorithm which incorporate* estimates of the second-order statistics of the environment signals is* introduced. The objective of the algorithm is to avoid the slowconvergence* of the LMS algorithm when the input signal is highly correlated. The* improvement is achieved at the expense of an increased computational* complexity.*/template Type lmsNewton( const Type &xk, const Type &dk, Vector &wn,const Type &mu, const Type &alpha, const Type &delta ){assert( 0 < alpha );assert( alpha <= Type(0.1) );int filterLen = wn.size();Type beta = 1-alpha;Vector vP(filterLen);Vector vQ(filterLen);static Vector xn(filterLen);// initialize the Correlation Matrix's inversestatic Matrix invR = eye( filterLen, Type(1.0/delta) );// update input signalfor( int i=filterLen; i>1; --i )xn(i) = xn(i-1);xn(1) = xk;Type yk = dotProd(wn,xn);// update the Correlation Matrix's inversevQ = invR * xn;vP = vQ / (beta/alpha+dotProd(vQ,xn));invR = (invR - multTr(vQ,vP)) / beta;// update the Weight Vectorwn += 2*mu * (dk-yk) * (invR*xn);return yk;}/*** The conventional LMS is very hard to choose a learning rate "mu" that* guarantees stability of the algorithm. The Normalised LMS is a variant* of the LMS that solves this problem by normalising with the power of* the input. For stability, the parameter "rho" should beong to (0,2),* and "gamma" is a small number to prevent  == 0.*/template Type lmsNormalize( const Type &xk, const Type &dk, Vector &wn,const Type &rho, const Type &gamma ){assert( 0 < rho );assert( rho < 2 );int filterLen = wn.size();static Vector sn(filterLen);// update input signalfor( int i=filterLen; i>1; --i )sn(i) = sn(i-1);sn(1) = xk;// get the outputType yk = dotProd( wn, sn );// update the Weight Vectorwn += rho*(dk-yk)/(gamma+dotProd(sn,sn)) * sn;return yk;}

测试代码：

/******************************************************************************                                  lms_test.cpp** LMS adaptive filter testing.** Zhang Ming, 2010-10, Xi'an Jiaotong University.*****************************************************************************/#define BOUNDS_CHECK#include #include #include using namespace std;using namespace splab;typedef double  Type;const   int     N = 50;const   int     order = 1;const   int     dispNumber = 10;int main(){Vector dn(N), xn(N), yn(N), wn(order+1);for( int k=0; k
运行结果：
 
The last 10 iterations of Conventional-LMS:observed        desired         output          adaptive filter-0.2225         -0.9749         -0.8216         -0.5683 1.08100.6235          -0.7818         -0.5949         -0.5912 1.08921.0000          -0.0000         0.0879          -0.6084 1.07840.6235          0.7818          0.6991          -0.5983 1.0947-0.2225         0.9749          0.8156          -0.6053 1.1141-0.9010         0.4339          0.2974          -0.6294 1.1082-0.9010         -0.4339         -0.4314         -0.6289 1.1086-0.2225         -0.9749         -0.8589         -0.6239 1.12910.6235          -0.7818         -0.6402         -0.6412 1.13531.0000          -0.0000         0.0667          -0.6543 1.1271The last 10 iterations of LMS-Newton:observed        desired         output          adaptive filter-0.2225         -0.9749         -0.9739         -0.7958 1.27790.6235          -0.7818         -0.7805         -0.7964 1.27831.0000          -0.0000         0.0006          -0.7966 1.27830.6235          0.7818          0.7816          -0.7966 1.2784-0.2225         0.9749          0.9743          -0.7968 1.2787-0.9010         0.4339          0.4334          -0.7971 1.2787-0.9010         -0.4339         -0.4340         -0.7971 1.2787-0.2225         -0.9749         -0.9747         -0.7971 1.27880.6235          -0.7818         -0.7816         -0.7972 1.27891.0000          -0.0000         0.0001          -0.7973 1.2789The last 10 iterations of Normalized-LMS:observed        desired         output          adaptive filter-0.2225         -0.9749         -0.9749         -0.7975 1.27900.6235          -0.7818         -0.7818         -0.7975 1.27901.0000          -0.0000         -0.0000         -0.7975 1.27900.6235          0.7818          0.7818          -0.7975 1.2790-0.2225         0.9749          0.9749          -0.7975 1.2790-0.9010         0.4339          0.4339          -0.7975 1.2790-0.9010         -0.4339         -0.4339         -0.7975 1.2790-0.2225         -0.9749         -0.9749         -0.7975 1.27900.6235          -0.7818         -0.7818         -0.7975 1.27901.0000          -0.0000         -0.0000         -0.7975 1.2790The theoretical optimal filter is:              -0.7972 1.2788Process returned 0 (0x0)   execution time : 0.109 sPress any key to continue.
 
上篇 (10个月前)：Kalman滤波算法的C++实现 
下篇 (10个月前)：各种类型RLS自适应滤波算法的C++实现