音声1


覚えること


目次

  1. 音声生成過程
  2. 線形予測分析
  3. ケプストラム
  4. 課題
  5. 参考文献

6-1音声生成過程

音のスペクトル解析を行った所では、振幅スペクトルと位相スペクトルについて
紹介して来ました。横軸が周波数で、縦軸が振幅や位相の成分でした。音のデータ
を特に音声データに限ると特徴的なスペクトルが現れます。特徴的に現れるスペク
トルの「音声のパラメータ」について今回は学びます。

まず音声生成の仕組みを見てみましょう。
音声成過程において、肺(lungs)、声門(glottis)は音源、声道(vocaltracts)は共振
口唇(lips)は放射としての役割を持ちます。 声門と声道(声門から口唇まで)のイメージ図は以下の通りです。声道は鼻につながる
経路がありますが図では省略してあります。


  1. 音源
    音源は声門は緊張して閉じた状態になり呼吸流により断続的な振動を起こします。
    男性の場合1秒間に100から150回、女声の場合250から300回ほど振動します。
    この振動数をピッチ周波数、あるいは単にピッチと呼び、声の高さを表します。
  2. 共振
    共振(共鳴)は物理学的には外部エネルギーによって固有振動が起こる状態を言います。
    音叉(おんさ)をたたいて発する音は小さいものですが、これを箱状の共鳴器に接続
    すると、驚くほど音が大きくなります。これも共振(共鳴)現象を利用したもので、
    音叉の振動数に合わせたサイズの共鳴器によって音圧が著しく高まるのです。

    声道を通るときに「ア」「イ」などの響きがつけられます。これは、口が形を変えて
    共振の特性を変えているから「ア」や「イ」に聴えるのです。この振幅スペクトルの
    ピークである共振周波数をフォルマント周波数と呼び、低い方から順に第1フォルマン
    ト(F1)、第2フォルマント(F2)、 第3フォルマント(F3)、、、と名付けられています
    (下図参照:ただし図は振幅スペクトルの包絡線のみを描いている)。

    「ア」などの母音はフォルマント周波数で特徴つけることが出来ます。次の図は
    F1とF2平面に日本語の5母音を配置したものです。フォルマント周波数は性別、年齢、
    話者によってかなり変動します。図よりF1-F2平面だけでは5母音を区別出来ない
    ことが分かります。単独で発話される母音はF1-F2-F3空間でほぼ完全に区別出来ます。

    実際のスペクトルは以下の図の様になります。これは/o/を発音しているときのスペクトルです。
    青線がFFTによる振幅スペクトル。緑線がケプストラム分析による振幅スペクトルの包絡線です。
    F1が400Hz付近に、F2が1100Hz付近にピークとして現れている様子が分かります。

  3. 周波数スペクトルに見る音源と共振

    上で少し紹介した様に、振幅スペクトル包絡のピーク(共振部分)は母音を区別するための
    特徴でした。子音の場合は共振に加えて、鼻へ通じる声道の分岐によって生じる反共振部分
    (スペクトル包絡の谷)も特徴量となります。この様に振幅スペクトルにおいて、スペクトル
    包絡、つまりスペクトルの空間周波数成分の低周波の部分が、声道の共振を表すと考えられ
    ます。

    一方、音源の特徴はスペクトルの空間周波数成分の高周波の部分に現れます。下の音声ファイルのサウンド
    スペクトログラムを見て下さい。赤い色はスペクトルの振幅値が大きいことを示しています。
    母音でピッチがはっきりと出るところは、周波数上で等間隔に赤い線が並んでいます。
    上の図でしたら、ほぼ等間隔に出現している青線の小刻みな高周波のピークがスペクトロ
    グラムでの赤い色になるわけです。このことから、音源を周波数軸上に等間隔に配置した
    インパルスでモデル化する方法が一般的に使われています。


  4. ではスペクトルの空間周波数成分の低周波(共振情報)や高周波(音源情報)をパラメータで表せないでしょうか。
    音声のパラメータのうち代表的なものとして、LPC(線形予測分析)とケプストラム
    があります。今回はケプストラムを中心に音声分析の具体例を見て行きます。

6-2LPC(線形予測分析)

LPCについては火曜日のゼミの本で詳しくやりますので、今回は説明を省略します。


6-3 ケプストラム

音声波x(n)は音源波g(n)と声道のインパルス応答v(n)との畳み込みとして表現出来る(放射の影響は
ないものとする)と仮定する。すなわち

x(n)=g(n)*v(n)

この信号x(n)に離散フーリエ変換の絶対値の対数演算操作Dを掛けると、対数は掛け算を和の形に
するので以下の様になる。

D{x(n)}=D{g(n)*v(n)}=D{g(n)}+D{v(n)}

(注:)積が和に置き換えられるのは対数演算のためである。

このときD{x(n)}の逆離散フーリエ変換をケプストラムc(n)と呼ぶ。対数演算によりスペクトルの
高周波と低周波が分離され、逆離散フーリエ変換を施すと時間軸の低次のところにスペクトル低周波
の共振情報が、時間軸の高次のところに音源情報が現れる。(逆離散フーリエ変換を行うので、ケプス
トラム変数は時間の次元を持つ。)


課題

次のMATLABプログラムを実行しなさい。またこのプログラムについて以下の質問に答えなさい。

============================================================================

フレーム表示のMOVIEを作るプログラム

  1. 波形表示 (WavFrame.m)
  2. FFT振幅スペクトル表示 (SpecFrame.m)
  3. ケプストラム表示 (CepFrame.m)
  4. 低次ケプストラムによる振幅スペクトル表示 (CepSpecFrame.m)

============================================================================

課題が出来たら、答えを山本(eli@sys.wakayama-u.ac.jp)までメイルで送って下さい。

締め切りは12/7(木)17:00までです。


プログラム1(WavFrame.m)

  1. コード
    %===================================================================
    close all;  % 画面クリア
    clear all;  % データクリア
    
    [x,fs]=aiffread('kousyouAD16.aiff'); %データ読み込み
    
    t=[0:1/fs:length(x)/fs-1/fs]; %時間軸の設定
    wlen=fs*0.032; %窓長
    wshift=fs*0.008; % 窓シフト
    
    wnum=length(x)/wshift-3; %フレーム数の計算
    wnum=floor(wnum); % 整数に直す
    
    f=figure; % 図のID取得
    set(f,'doublebuffer','on');%アニメーション表示設定
    M=moviein(wnum,f); %ムービーバッファ初期化
    
    hw=hanning(wlen);
    
    
    for n=0:wnum-1
    
      S=n*wshift; %フレーム開始点
      xtmp=x(S+1:S+wlen); % フレームデータ取得
      xtmp2=hw.*xtmp; %窓掛け
      Z(1:wlen)=  xtmp2'; % データの転置
       
      subplot(2,1,1);plot(t(1:S+wlen),x(1:S+wlen), ...
                  t(S+wlen+1:end),x(S+wlen+1:end));
      xlabel('time');ylabel('Amplitude');
        
      subplot(2,1,2);plot(Z);axis([1 wlen -10^4 10^4]);
      xlabel('Data Points');ylabel('Amplitude');
            
      M(:,n+1)=getframe(f);%ムービー画像の作成
      disp(['n=' int2str(n)]); %フレーム番号表示
    
    end
    
    movie(f,M);%ムービーの再生
    


    %===================================================================

  2. 画面
  3. 注意事項

プログラム2(SpecFrame.m)

  1. コード
    %===================================================================
    close all;  % 画面クリア
    clear all;  % データクリア
    
    [x,fs]=aiffread('kousyouAD16.aiff'); %データ読み込み
    
    t=[0:1/fs:length(x)/fs-1/fs]; %時間軸の設定
    wlen=fs*0.032; %窓長
    wshift=fs*0.008; % 窓シフト
    
    wnum=length(x)/wshift-3; %フレーム数の計算
    wnum=floor(wnum); % 整数に直す
    
    f=figure; % 図のID取得
    set(f,'doublebuffer','on');%アニメーション表示設定
    M=moviein(wnum,f); %ムービーバッファ初期化
    
    hw=hanning(wlen);
    fq=(1:wlen)*fs/wlen;  fq=fq(1:wlen/2);%周波数軸の設定
    
    
    for n=0:wnum-1
    
      S=n*wshift; %フレーム開始点
      xtmp=x(S+1:S+wlen); % フレームデータ取得
      xtmp2=hw.*xtmp;
      
      y=fft(xtmp2');
      Z= 20*log10(abs(y(1:wlen/2))) ;
         
      subplot(2,1,1);plot(t(1:S+wlen),x(1:S+wlen), ...
                  t(S+wlen+1:end),x(S+wlen+1:end));
      xlabel('time');ylabel('Amplitude');
       
      subplot(2,1,2);plot(fq,Z);axis([1 fs/2 0 120]);
      xlabel('Frequency');ylabel('Magnitude'); 
              
      M(:,n+1)=getframe(f);%ムービー画像の作成
      disp(['n=' int2str(n)]); %フレーム番号表示
    
    end
    
    movie(f,M);%ムービーの再生
    


    %===================================================================


  2. 画面
  3. 注意事項


プログラム3(CepFrame.m)

  1. コード
    %===================================================================
    close all;  % 画面クリア
    clear all;  % データクリア
    
    [x,fs]=aiffread('kousyouAD16.aiff'); %データ読み込み
    
    t=[0:1/fs:length(x)/fs-1/fs]; %時間軸の設定
    wlen=fs*0.032; %窓長
    wshift=fs*0.008; % 窓シフト
    
    wnum=length(x)/wshift-3; %フレーム数の計算
    wnum=floor(wnum); % 整数に直す
    
    f=figure; % 図のID取得
    set(f,'doublebuffer','on');%アニメーション表示設定
    M=moviein(wnum,f); %ムービーバッファ初期化
    
    hw=hanning(wlen);
    
    
    for n=0:wnum-1
    
      S=n*wshift; %フレーム開始点
      xtmp=x(S+1:S+wlen); % フレームデータ取得
      xtmp2=hw.*xtmp;
      
      y=fft(xtmp2');
      Z=real(ifft(log(abs(y))));% ケプストラム
       
      subplot(2,1,1);plot(t(1:S+wlen),x(1:S+wlen), ...
                  t(S+wlen+1:end),x(S+wlen+1:end));
      xlabel('time');ylabel('Amplitude');
        
      subplot(2,1,2);plot(Z);axis([1 wlen -1 1]);
      xlabel('Quefrency');ylabel('Amplitude');
            
      M(:,n+1)=getframe(f);%ムービー画像の作成
      disp(['n=' int2str(n)]); %フレーム番号表示
    
    end
    
    movie(f,M);%ムービーの再生
    


    %===================================================================

  2. 画面
  3. 注意事項

プログラム4(CepSpec.m)

  1. コード
    %===================================================================
    close all;  % 画面クリア
    clear all;  % データクリア
    
    [x,fs]=aiffread('kousyouAD16.aiff'); %データ読み込み
    
    t=[0:1/fs:length(x)/fs-1/fs]; %時間軸の設定
    wlen=fs*0.032; %窓長
    wshift=fs*0.008; % 窓シフト
    
    wnum=length(x)/wshift-3; %フレーム数の計算
    wnum=floor(wnum); % 整数に直す
    
    f=figure; % 図のID取得
    set(f,'doublebuffer','on');%アニメーション表示設定
    M=moviein(wnum,f); %ムービーバッファ初期化
    
    hw=hanning(wlen);
    fq=(1:wlen)*fs/wlen;  fq=fq(1:wlen/2);%周波数軸の設定
    
    for n=0:wnum-1
    
      S=n*wshift; %フレーム開始点
      xtmp=x(S+1:S+wlen); % フレームデータ取得
      xtmp2=hw.*xtmp;
      
      y=fft(xtmp2');
      N=30;
      cep=real(ifft(log(abs(y))));
      lifter=[ones(1,N) zeros(1,length(cep)-2*N) ones(1,N)];%リフタの作成
      lowcep=cep.*lifter;%リフタの重みつけによる低次ケプストラムの抽出
       
      yy=fft(lowcep);
    
      Z1= 20*log10(abs(y(1:wlen/2))) ;
      Z2= 20*log10(exp(abs(yy(1:wlen/2))) );
    
       
      subplot(2,1,1);plot(t(1:S+wlen),x(1:S+wlen), ...
                  t(S+wlen+1:end),x(S+wlen+1:end));
      xlabel('time');ylabel('Amplitude');
        
      subplot(2,1,2);plot(fq,Z1,fq,Z2);axis([1 fs/2 0 120]);
      xlabel('Frequency');ylabel('Magnitude');
            
      M(:,n+1)=getframe(f);%ムービー画像の作成
      disp(['n=' int2str(n)]); %フレーム番号表示
    
    end
    
    movie(f,M);%ムービーの再生
    


    %===================================================================

  2. 画面
  3. 注意事項


参考文献

MATLAB Getting Started with MATLAB
ユーザーズディジタル信号処理