File content as of revision 0:db50a6a3e574:

//---------------------------------------------------------------------------
// 学習同定法による線スペクトル強調器（ALE）
//      係数更新で Leaky LMS 法を使用
//      次数：100 次，遅延器数：5 個
//      作成者，著作権者： 三上直樹 2017/02/21
//---------------------------------------------------------------------------

#include "F446_ADC_Interrupt.hpp"   // AD, DA
#include "AQM1602.hpp"              // LCD 表示器
using namespace Mikami;

const int FS_ = 24000;      // 標本化周波数: 24 kHz

const int ORDER_ = 200;             // 次数
const int DELAY_ = 5;               // 相関除去用の遅延器数
const int N_ALL_ = ORDER_ + DELAY_; // 全体で必要な遅延器数

const float GM_ = 0.996;            // γ
const float S0_ = 0.1;              // 分母が 0 になるのを防止
float alpha_ = 0.2;                 // α

float x_[N_ALL_+1], h_[ORDER_+1];

AdcDual_Intr myAdc_(FS_);   // 参照："F446_ADC_Interrupt.hpp"
DacDual myDac_;             // 参照："F446_DAC.hpp"
DigitalOut ledG_(D10, 1);   // LED 緑色

// ADC 変換終了割り込みに対する割り込みサービス･ルーチン
void AdcIsr()
{
    static float mu, err_mu;
    static float ss = 0.0f;

    float xn1, xn2, yn1, yn2;
    myAdc_.Read(xn1, xn2);          // 入力（xn2 は使わない）

    // 左チャンネルの処理（FIR フィルタの計算）
    x_[0] = xn1;
    yn1 = 0.0;
    for(int k=0; k<=ORDER_; k++) yn1 = yn1 + x_[k+DELAY_]*h_[k];

    // ステップ･サイズ･パラメータの計算
    ss = GM_*ss + xn1*xn1;          // パワに比例する値の推定値
    mu = alpha_/(ss + S0_);         // ステップ･サイズ･パラメータの計算
    // 係数の更新
    err_mu = (xn1 - yn1)*mu;        // 誤差*μ
    for(int k=0; k<=ORDER_; k++)
        h_[k] = GM_*h_[k] + err_mu*x_[k+DELAY_];

    // 遅延器内のデータの移動
    for (int k=N_ALL_; k>0; k--) x_[k] = x_[k-1];

    // 右チャンネル出力は入力をそのまま出力する
    yn2 = xn1;

    myDac_.Write(yn1, yn2);         // 出力
}

int main()
{
    printf("\r\nALE using normalized LMS algorithm\r\n");
    AnalogIn a3In(A3);      // VR からの電圧読み取り用
    Aqm1602 lcd;            // LCD 表示器
    lcd.WriteStringXY("ALE", 0, 0);

    // 出力の LPF の遮断周波数を 10 kHz に設定
    myDac_.ScfClock(10000*100);

    // ALE のバッファのクリア
    for (int k=0; k<=N_ALL_; k++) x_[k] = 0.0f;
    for (int k=0; k<=ORDER_; k++) h_[k] = 0.0f;

    // ADC 変換終了割り込みに対する割り込みサービス･ルーチン割り当て
    myAdc_.SetIntrVec(&AdcIsr);

    float vrOld = -1;
    while(1)
    {
        float vrNow = a3In.read();
        if (fabsf(vrOld - vrNow) > 0.02f)
        {
            alpha_ = powf(10.0f, (vrNow-1)*2);
            vrOld = vrNow;

            printf("a0 = %5.3f\r\n", alpha_);
            lcd.ClearLine(1);
            lcd.WriteStringXY("alpha: ", 0, 1);
            lcd.WriteValue("%6.3f", alpha_);
        }
        wait(0.2f);
    }
}