FFTを自前で実装して周波数分析(Player10)

2008/08/31

こんばんは。近頃、大雨がすごいですね。夜中に雷を伴う雨が数日続くなんてことは生まれてこのかた、あんまり体験したことないんで驚いています。北京五輪で、雨を降らせないためにロケットを何十発も打ち込んだなんて話をきくと、「まさか、それと関係じなんてないよね?」なんて勘ぐってみたくなってしまう、今日この頃です。

さて、今回は周波数を分析するFFTについてです。Player9からSoundMixer.computeSpectrumっちゅう関数が実装されて、音を波形から周波数成分に分解することができるようになりました。
それで、これはこれでよかったんですが、困ることもありました。それは、たとえばニコニコ動画とかをみながら、このメソッドを実行するswfが埋め込まれたページを開いたとします。すると、埋め込まれたページだけでなく、関係ないニコニコの音まで拾って解析しようとするんです。するとセキュリティサンドボックスエラーみたいなことが起こってうまく実行できなくなってしまうんです。

なので、Adobe標準のものじゃなくて、自前で実装しようと思いました。ただ、FFT(=高速離散フーリエ変換)を自分で一から実装するのは正直難しいです。なので、これのアルゴリズムをどこかから移植しようと思いました。
いろいろとあって、Java版やExcel版の簡単そうなのを書き直してみたりしたんですが、いまいち動きません。私の理解不足なためにそうなってしまうんでしょう。それで、「もうだめだ!」と投げかけてたところに現れたのが、このページです。
→FFT (高速フーリエ・コサイン・サイン変換) の概略と設計法
いや、すごい。本当にすごいです。大浦様、ありがとうございました。ソースはフリーで使ってもよいと書いてあったんでたぶん大丈夫でしょう。
それで、ここから基本となるプログラムを移植して完成したものが、こちらです。
→サンプルページ(要Player10)

今回それなりにフーリエ変換について勉強しました。最後の部分に移植したFFTのプログラムを載せます。これの使い方なんですが、以下のようになります。

・入力データサンプル数は2のべき乗個分用意する。(2, 4, 6, ….1024, 2048)
・出力されたデータは左右対称の形をしている。このうち右半分のデータは不必要なのでカットする。
・カットしたデータが周波数成分表になるが、横軸の最大値がサンプリングレートの半分となる。
(たとえば44.1kHzでサンプリングしているとき、22.05kHzが最大値となる。今回Player10の音生成とからめているので、Palyer10の音生成のサンプリングレートは44.1kHzのため、これにあてはまることになる。)
・縦軸の単位はわからない。が、割合を表していると思う。
(たとえば、棒が2本でていて、それぞれ3,7ぐらいの高さだとすると、その棒の位置の周波数成分が3:7の割合だということ)

とりあえずソースです。今回Player10から登場したVectorクラスを初めて使ってみました。Arrayと比べて、書くのが面倒ですが、なんか計算速いみたいですね!

package
{
  import __AS3__.vec.Vector;

  public class FFT
  {
    public static function getFFT(inputData:Vector.<Number>):Vector.<Number>
    {
      var n:uint = inputData.length; //データの数。必ず2の乗数であること
      var theta:Number = 2 * Math.PI / n;
      var ar:Vector.<Number>, ai:Vector.<Number>;
      var m:int, mh:int, i:int, j:int, k:int, irev:int;
      var wr:Number, wi:Number, xr:Number, xi:Number;

      //データの読み込み。xrは実数部、xiは虚数部
      ar = new Vector.<Number>(n);
      ai = new Vector.<Number>(n);
      for(i = 0; i < n; i++)
      {
        ar[i] = inputData[i];
        ai[i] = 0;
      }

      //scrambler
      i = 0;
      for (j = 1; j < n - 1; j++) {
        for (k = n >> 1; k > (i ^= k); k >>= 1);
        if (j < i) {
          xr = ar[j];
          xi = ai[j];
          ar[j] = ar[i];
          ai[j] = ai[i];
          ar[i] = xr;
          ai[i] = xi;
        }
      }
      for (mh = 1; (m = mh << 1) <= n; mh = m) {
        irev = 0;
        for (i = 0; i < n; i += m) {
          wr = Math.cos(theta * irev);
          wi = Math.sin(theta * irev);
          for (k = n >> 2; k > (irev ^= k); k >>= 1);
          for (j = i; j < mh + i; j++) {
            k = j + mh;
            xr = ar[j] - ar[k];
            xi = ai[j] - ai[k];
            ar[j] += ar[k];
            ai[j] += ai[k];
            ar[k] = wr * xr - wi * xi;
            ai[k] = wr * xi + wi * xr;
          }
        }
      }

      var output:Vector.<Number> = new Vector.<Number>(n);
      for(i = 0; i < n; i++)
      {
         //output[i] = ar[i]; //実数部
         //output[i] = ai[i]; //虚数部
         output[i] = Math.sqrt(Math.pow(ar[i], 2) + Math.pow(ai[i], 2));
      }

      return output;
    }
  }
}

あと、おまけですがクライアント側のプログラムも載せときます。、、、長い。

package
{
  import __AS3__.vec.Vector;

  import fl.controls.Button;
  import fl.controls.ComboBox;
  import fl.controls.Label;
  import fl.controls.Slider;
  import fl.events.SliderEvent;

  import flash.display.Graphics;
  import flash.display.Sprite;
  import flash.display.StageScaleMode;
  import flash.events.Event;
  import flash.events.MouseEvent;
  import flash.events.SampleDataEvent;
  import flash.media.Sound;
  import flash.media.SoundChannel;
  import flash.text.TextField;
  import flash.text.TextFieldAutoSize;
  import flash.text.TextFormat;

  [SWF(width="640", height="350", frameRate="30", backgroundColor="#b6c3be")]
  public class FFTsample extends Sprite
  {
    private var graphBase:Sprite;
    private var graph:Sprite;
    private var asset:FFT_asset;
    private var isPlay:Boolean;
    private var ch1:SoundChannel;
    private var sound:Sound;
    private var waveTypeCombo:ComboBox;
    public static const SINE:String = "sine";
    public static const SAW:String = "saw";
    public static const TRIANGLE:String = "triangle";
    public static const PULSE:String = "pulse";
    public static const NOISE:String = "noise";
    private var waveType:String = SINE;
    private var frequencySlider:Slider;
    private var frequency:uint = 440;

    public function FFTsample()
    {
      super();
      initialize();
    }

    private function initialize():void
    {
      stage.scaleMode = StageScaleMode.NO_SCALE;
      initializeGraph();
      initializeAsset();

      frequency = 440;
      isPlay = true;
      sound = new Sound();
      sound.addEventListener("sampleData", sampleDataComplete);
    }

    private function sampleDataComplete(e:SampleDataEvent):void
    {
      var input:Vector.<Number> = new Vector.<Number>();
      var output:Vector.<Number>;
      var sample:Number;
      //var frequency:uint = 440;
      var f0:Number = frequency / 44100;
      var f1:Number = 44100 / frequency;
      var PAI2:Number = 2 * Math.PI;
      for ( var c:int=0; c < 2048; c++ )
      {
        switch(waveType)
        {
          //サイン波
          case SINE:
          sample = Math.sin(PAI2 * f0 * Number(c + e.position));
          break;

          //矩形波
          case PULSE:
          sample = (Math.sin(PAI2 * f0 * Number(c + e.position)) < 0) ? 1 : -1;
          break;

          //のこぎり波
          case SAW:
          sample = 2 * f0 * ((c + e.position) % f1) - 1;
          break;

          //三角波
          case TRIANGLE:
          sample = (Math.sin(PAI2 * f0 * Number(c + e.position)) < 0) ?
              (4 * f0 * ((c + e.position) % (f1 / 2)) - 1) :
              (-4 * f0 * ((c + e.position) % (f1 / 2)) + 1);
          break;

          //ホワイトノイズ
          case NOISE:
          sample = Math.random() * 2 - 1;
          break;

          default: break;
        }
        input[c] = sample;
        e.data.writeFloat(sample);
        e.data.writeFloat(sample);
      }

      output = FFT.getFFT(input);
      output = output.splice(0, output.length / 2);
      drawFFT(output, graph, 380, 280);
    }

    private function drawFFT(data:Vector.<Number>, target:Sprite, WIDTH:uint, HEIGHT:uint):void
    {
      var g:Graphics = target.graphics;
      var i:uint;
      var len:uint = data.length;
      var rate:Number = WIDTH / len;
      var step:Number = (rate < 1) ? (1 / rate) : 1;
      g.clear();
      g.lineStyle(1, 0xff0000, 1);
      g.moveTo(0, HEIGHT);
      for(i = 0; i < data.length; i += step)
      {
        g.moveTo(i * rate, HEIGHT);
        g.lineTo(i * rate, HEIGHT - HEIGHT * data[i] / 800);
      }
    }

    private function initializeGraph():void
    {
      var g:Graphics;
      var devideNum:uint = 10;
      var i:uint;
      var devideWidth:Number;
      var HEIGHT:uint;
      var tf:TextField;
      var tft:TextFormat;
      var samplingRate:uint = 44.1;
      graphBase = new Sprite();
      addChild(graphBase);
      graph = new Sprite();
      graphBase.addChild(graph);
      g = graphBase.graphics;
      g.beginFill(0x000000, 1);
      g.drawRect(0,0,380,280);
      g.endFill();
      graphBase.x = 15;
      graphBase.y = 10;

      g.lineStyle(1, 0, 1);
      devideWidth = graphBase.width / devideNum;
      HEIGHT = graphBase.height + 5;
      tft = new TextFormat();
      tft.size = 10;
      for(i = 0; i <= devideNum; i++)
      {
        g.moveTo(i * devideWidth, HEIGHT);
        g.lineTo(i * devideWidth, HEIGHT + 5);
        tf = new TextField();
        tf.autoSize = TextFieldAutoSize.CENTER;
        tf.text = (Math.floor(samplingRate / 2 / devideNum * i * 10) / 10).toString();
        addChild(tf);
        tf.setTextFormat(tft);
        tf.x = i * devideWidth;
        tf.y = HEIGHT + 15;
      }

      tf = new TextField();
      tf.autoSize = TextFieldAutoSize.LEFT;
      tf.text = "[kHz]";
      addChild(tf);
      tf.setTextFormat(tft);
      tf.x = graphBase.width + 20;
      tf.y = HEIGHT + 15;
    }

    private function initializeAsset():void
    {
      asset = new FFT_asset();
      addChild(asset);

      var btn:Button = asset.startStop;
      btn.label = "スタート";
      btn.addEventListener(MouseEvent.MOUSE_DOWN, startStopHandler);

      initializeCombobox();
      initializeFrequencySlider();
    }

    private function initializeCombobox():void
    {
      waveTypeCombo = asset.combo;
      waveTypeCombo.addItem( { label: "サイン波", data:SINE } );
      waveTypeCombo.addItem( { label: "矩形波", data:PULSE } );
      waveTypeCombo.addItem( { label: "のこぎり波", data:SAW } );
      waveTypeCombo.addItem( { label: "三角波", data:TRIANGLE } );
      waveTypeCombo.addItem( { label: "ノイズ", data:NOISE } );
      waveTypeCombo.addEventListener(Event.CHANGE, waveTypeComboChangeHandler);
    }

    private function waveTypeComboChangeHandler(e:Event):void
    {
      waveType = waveTypeCombo.value;
    }

    private function initializeFrequencySlider():void
    {
      frequencySlider = asset.freqSlide;
      frequencySlider.value = frequency;
      frequencySlider.width = 200;
      frequencySlider.tickInterval = 10;
      frequencySlider.maximum = 20000;
      frequencySlider.minimum = 20;
      frequencySlider.addEventListener(SliderEvent.THUMB_DRAG, frequencySliderDrag);
      var label:Label = asset.freq_txt;
      label.text = frequency + "Hz";
    }

    private function frequencySliderDrag(e:SliderEvent):void
    {
      frequency = frequencySlider.value;
      var label:Label = asset.freq_txt;
      label.text = frequency + "Hz";
    }

    private function startStopHandler(e:MouseEvent):void
    {
      var btn:Button = asset.startStop;
      if(isPlay)
      {
        btn.label = "ストップ";
      }
      else
      {
        btn.label = "スタート";
      }
      playStopSound();
      isPlay = !isPlay;
    }

    private function playStopSound():void
    {
      if(isPlay)
      {
        ch1 = sound.play();
      }
      else
      {
        ch1.stop();
      }
    }
  }
}

自作iPhoneアプリ 好評発売中!
フォルメモ - シンプルなフォルダつきメモ帳
ジッピー電卓 - 消費税や割引もサクサク計算!

ページトップへ戻る