+Java Web Start
1/f雑音生成アルゴリズム

細田 隆之
Feb. 1, 2002

1/f雑音は増幅器の低周波域の雑音として有名なピンクノイズとも言われる3dB/オクターブで電力が減衰する雑音である。 近年、1/fのゆらぎが人に心地よいとして、エアコンや扇風機の風量制御や照明などにも応用されている。 そのような用途のために、疑似ホワイトノイズから1/f雑音を生成するアルゴリズムを考案した。
1/fの特性、つまり 3dB/oct. の特性を得るには様々な方法があろうが、ここでは、 fig.1 に示すラグリードフィルタ(lag-lead filter)を縦列接続したラダーフィルタ(ladder filter)に相当する伝達関数により1/f特性を近似した。 このラダーフィルタの伝達関数をデジタルフィルタで実現し、疑似ホワイトノイズにこのフィルターを施すことにより1/f雑音を得ている。 伝達関数 f(s) は次の式で表される。
Fig.1 ラグリードフィルタ伝達関数
lag-lead filter f(s)= \prod_{k=0}^n  \displaystyle\frac{2^{2k}s+1}{2^{2k+1}s+1}
(1) 式の周波数・振幅特性例を Fig.2 に示す。
Fig.2 近似式の振幅・周波数特性例 (k = -4〜4)Fig.3 SPICE シミュレーション例
Download source
Download SPICE circuit netlist

下のCプログラム例は、random() でホワイトノイズを生成した後に pinkfilter() で 3dB/oct のフィルタをかけている。 pinkfilter() は MAX_Z オクターブに渡って 3dB/oct の特性を示すデジタルフィルタをIIR型で構成している。 pinkfilter() を使用する前に、init_pink() を呼び出し、フィルタ係数を初期化する必要がある。

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

#define MAX_Z 16

double  z[MAX_Z];
double  k[MAX_Z];

double pinkfilter(double in) {
    extern double   z[MAX_Z];
    extern double   k[MAX_Z];
    static double   t = 0.0;
    double          q;
    int             i;

    q = in;
    for (i = 0; i < MAX_Z; i++) {
        z[i] = (q * k[i] + z[i] * (1.0 - k[i]));
        q = (q + z[i]) * 0.5;
    }
    return (t = 0.75 * q + 0.25 * t); /* add 1st order LPF */
}

void init_pink() {
    extern double   z[MAX_Z];
    extern double   k[MAX_Z];
    int             i;

    for (i = 0; i < MAX_Z; i++)
        z[i] = 0;
    k[MAX_Z - 1] = 0.5;
    for (i = MAX_Z - 1; i > 0; i--)
        k[i - 1] = k[i] * 0.25;
}

int main() {
    int             i, j;
    double          n;

    init_pink();
    srandomdev();
    for (i = 0; i < 65536; i++) {
        printf("%.8g\n", pinkfilter(random() & 1 ? 1.0 : -1.0));
    }
}
Fig.4, Fig.5 にプログラム出力の波形とFFT結果の例を示す。
Fig.4 出力波形例Fig.5 出力のFFT例 (217ポイント Nuttall 窓*note1)
Download PinkNoise Demo [5kB, Java Web Start]
⚠️ Demo の実行には Java 8 のインストール例外サイトへの追加 が必要です。

*note1: 窓関数
Fig.6 窓関数のDFT特性例
Window functions
Blackman 窓 : 比較的狭帯域でメインローブ近傍のサイドローブが -68dB未満とまぁまぁでよく利用される。
 Blackman_Ex(x) = (7938 - 9240 cos(2πx) + 1430 cos(4πx)) / 18608 | 0 ≤ x ≤ 1

Nuttall(4b) 窓 : メインローブ近傍のサイドローブが約 -90dB と低く減衰率も -18dB/oct. と大きく汎用性が高い。
  Nuttall_4b(x) = 2.355768 - 0.487396 cos(2πx) + 0.144232 cos(4πx) - 0.012604 cos(6πx) | 0 ≤ x ≤ 1

Blackman-Nuttall 窓 : メインローブ近傍からサイドローブが-98dB未満と割りと良く高ダイナミックレンジ向き。
  Blackman-Nuttall(x) = 0.3635819 - 0.4891774 cos(2πx) + 0.1365995 cos(4πx) - 0.0106411 cos(6πx) | 0 ≤ x ≤ 1

Finetune5 窓 : デシメーション用に直流近傍のエイリアス周波数における減衰量を重視して設計したもの。

SEE ALSO

Window functions


このページの内容は、 日本国著作権法および国際条約により保護されています。
www.finetune.co.jp [Mail] © 2000 Takayuki HOSODA.
Powered by
 Finetune