Wavelet超入門 その7 DaubechiesD20に挑戦


Daubechiesの係数を20にして、waveletを求めた

 1.Daubechiesの係数を20にして、scaling function、mother waveletを描く

    そのためには、まず、Φ関数の、整数における値を求めます。 Φ(n) n: integerts

    n=0から19までの20個の内、 

    n=0 n=19では、  Φ(0)=Φ(19)=0 です。

    残りの n=1〜18までの時の、Φ(n)を計算します。

    これは、1次連立方程式で表現できます。

    

    この行列の内、1の値が入っている行は

      Φ(1)+Φ(2)+Φ(3)+Φ(4)+Φ(5)+Φ(6)+Φ(7)=0  
    
    と言う、規格化のための条件を表現しています。

    この1の値の入った行を取り除くと、規則性が見えてきます。

    
    つまり、 p1-1となるのは、1行1列目、p2-1は、2行2列目です。

    このような規則性を見つけて、Octaveにてプログラムしてみました。

    D20(係数が20個)の場合


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

    こうして、Octaveによるプログラムで、数値が合っている事を確認しました。

    ここの部分の、Octaveでのプログラムは、こんな風です。

##################################################################
#Φ(0:N)の値を求めるプログラム。 Nは整数
##################################################################
close all;

D=input("\n  ==calc Daubechies's phi[integers] values==\n please input D. D is DaubechiesD??\n  D=4,8,20.  D:");


#Daubechies coefficients size D=N*2
p=zeros(1,D);

#今の所、係数の数は、4,8,20のみ。 移すのがメンドイ (^_^;;
switch D 
 case {4}
  p=[0.6830,1.1830,0.3170,-0.1830];
 case {8}
 p=[0.32580343,1.01094572,0.8922014,-0.03957503,-0.26450717,0.0436163,0.0465036,-0.01498699];
 case {20}
 p=[0.03771716,0.26612218,0.74557507,0.97362811,\
    0.39763774,-0.35333620,-0.27710988,0.18012745,\
    0.13160299,-0.10096657,-0.04165925,0.04696981,5.10043697e-3,\
    -0.01517900,1.97332536e-3,2.81768659e-3,-9.69947840e-4,-1.64709006e-4,1.32354367e-4,-1.875841e-5];
 otherwise
 error ("invalid value");
endswitch

a=zeros(D-2,D-2);
sub=ones(1,D-2);

#先ず、連立方程式を行列で表現する。(D-2)*(D-2)の行列を作る
for i=1:D-2
 for k=1:D-2

  if (i==k) retval=p(i+1)-1;
  else
     if(i*2-k+1<=D && (i*2-k+1>0) )
   retval=p(i*2-k+1);
     else retval=0;
     endif
  endif
  a(i,k)=retval;
 endfor
endfor  
# 最上段に1行の追加
 b=[sub;a];

#連立方程式を解く
 x=zeros(D-2,1);
 c=zeros(D-1,1);
 c(1)=1;
# b*x=c, solve
 x= b\c;
#転置行列
 x=x';
#配列の最初と最後に、Φ(0)=Φ(D-1)=0を加える
phi_value=[0,x,0];

#################################################################################################

      pの係数の数が、20にもなると、1行では、収まりません。
 
    で、複数行に渡ってしまう訳ですが、ここで、大きなミスをしていました。

    改行するときは、バックスラッシュを付けなくてはいけない!

    
これを忘れており、3日ほど、苦悶しちりました、涙。


 2.解像度も可変にする

  前節でも、解像度は、可変できましたので、これと、上記のプログラムを合体しました。Octaveを使ってます、感謝。

  全体のmファイルを、ここに置きます、よかったら使ってみて下さい。

 3.DaubechiesD8で、解像度が、1/16 も描いてみる
            
      

   ついでに DaubechiesD4,reso=16も
      

    最後は、最初に出てきたグラフです
        
      


    いけてるようです、バンザイ (@_@)/ 

H20.10.23