//---------------------------------------------------------------------------
//	適応線スペクトル強調器（ALE）のクラス：AleLMS
//		LMS アルゴリズムによる ALE
//
//		作成者：三上直樹，2018/02/16 作成，ⓒ三上直樹 2018
//---------------------------------------------------------------------------

using System;
using System.Linq;

namespace Mikami.Filter
{
	public class AleLMS
	{
		private double[] xm_, xk_, hk_;
		private int order_, delay_;
		private double mu_;

		// コンストラクタ
		public AleLMS(int delay, int order, double mu)
		{
			xm_ = new double[delay];		// 相関除去のための遅延器数
			xk_ = new double[order+1];		// FIR フィルタの遅延器
			hk_ = new double[order+1];		// FIR フィルタの係数
			Array.Clear(xm_, 0, xm_.Length);
			Array.Clear(xk_, 0, xk_.Length);
			Array.Clear(hk_, 0, hk_.Length);
			order_ = order;		// フィルタの次数
			delay_ = delay;		// 遅延器の数
			mu_ = mu;			// μ
		}

		// ALE を実行する
		public float Execute(float xn)
		{
			xm_[0] = xn;					// 遅延器の初段に入力信号を格納
			for (int k=order_-1; k>=0; k--)	// フィルタ内の信号のシフト
				xk_[k+1] = xk_[k];
			xk_[0] = xm_[delay_-1];			// フィルタの初段に遅延器の最終段の信号を格納
			for (int k=delay_-2; k>=0; k--)	// 遅延器内の信号のシフト
				xm_[k+1] = xm_[k];

			// FIRフィルタの計算
			double yn = 0.0;
			yn = xk_.Zip(hk_, (x, h) => x*h).Sum();

			// 係数の更新
			double error = xn - yn;			// 誤差信号の計算
			hk_ = hk_.Zip(xk_, (h, x) => h + mu_*error*x).ToArray();

			return (float)yn;
		}
	}
}
