2

C でウェーブレット変換を実装しようとしていますが、これまでに行ったことはありません。私はウェーブレットについていくつか読んで、「成長する部分空間」のアイデアと、マラットの片側フィルターバンクが本質的に同じアイデアであることを理解しています。

しかし、Mallat の高速ウェーブレット変換を実際に実装する方法に行き詰まっています。これは私がこれまでに理解していることです：

  1. ハイパス フィルター h(t) は、詳細係数を提供します。特定のスケール j の場合、これは、マザー ウェーブレット W(t) の反映され、拡張され、正規化されたバージョンです。

  2. g(t) は、差を構成するローパス フィルターです。h(t) の直交ミラーであると想定されます。

  3. 詳細係数、または j 番目のレベルの近似係数を取得するには、信号ブロックをそれぞれ h(t) または g(t) で畳み込み、2^{j} で信号をダウンサンプリングする必要があります (つまり、2^ {j}値)

しかし、これらは私の質問です：

  1. h(t) がわかっている場合、どうすれば g(t) を見つけることができますか?

  2. この変換の逆を計算するにはどうすればよいですか?

参照できる C コードはありますか? （はい、wikiで見つけましたが、役に立ちません）

私がいくつかのコードに言いたいことは次のとおりです。

A. フィルタは次のとおりです。

B. これは変換です (非常に明示的に) C.) これは逆変換です (これもダミーです)

ご辛抱いただきありがとうございますが、Step1 - Step2 - Step3 などの明示的な例を示すガイドはないようです (すべての係数が 1 であり、混乱を招くため、HAAR ではありません)。

4

2 に答える 2

4

fwt のマラット レシピは実にシンプルです。matlab コード (例: Jeffrey Kantor によるスクリプト) を見ると、すべての手順が明らかです。

C ではもう少し手間がかかりますが、それは主に、独自の宣言と割り当てを処理する必要があるためです。

まず、あなたの要約について：

  1. 通常、フィルターhはローパス フィルターであり、スケーリング関数 (父) を表します。
  2. 同様に、gは通常、ウェーブレット (マザー) を表すハイパス フィルターです。
  3. 1 つのフィルタリング + ダウンサンプリング ステップでJレベルの分解を実行することはできません。各レベルで、hでフィルタリングしてダウンサンプリングすることによって近似信号cを作成し、 gでフィルタリングしてダウンサンプリングすることによって詳細信号dを作成し、次のレベルでこれを繰り返します (現在のcを使用) 。

ご質問について:

  1. 直交ウェーブレット基底 [h_1 h_2 .. h_m h_n]のフィルターhの場合、QMF は [h_n -h_m .. h_2 -h_1] です。ここで、 nは偶数で、m == n -1です。
  2. 逆変換は fwt の逆を行います。各レベルで詳細dと近似cをアップサンプリングし、 dをgで、cをhで畳み込み、信号を加算します。対応する matlab scriptを参照してください。

この情報を使用し、タイプxのlenポイントの信号、係数のdoubleスケーリングフィルターhとウェーブレットgフィルター(タイプ)、および分解レベル を指定すると、このコードはマラット fwt を実装します。fdoublelev

double *t=calloc(len+f-1, sizeof(double));
memcpy(t, x, len*sizeof(double));        
for (int i=0; i<lev; i++) {            
    memset(y, 0, len*sizeof(double));
    int len2=len/2;
    for (int j=0; j<len2; j++)           
        for (int k=0; k<f; k++) {          
            y[j]     +=t[2*j+k]*h[k];      
            y[j+len2]+=t[2*j+k]*g[k];      
        }
    len=len2;                            
    memcpy(t, y, len*sizeof(double));    
}
free(t);

次の反復のために近似値(開始する入力信号)tをコピーするための「ワークスペース」という 1 つの追加の配列を使用します。cx

この Cプログラムの例を参照してください。これを使用してコンパイルgcc -std=c99 -fpermissive main.cppおよび実行できます./a.out。

逆もこれらの線に沿ったものでなければなりません。幸運を！

于 2014-09-23T13:00:02.517 に答える
1

唯一欠けているのは、フィルター操作のためのパディングです。台詞

y[j]     +=t[2*j+k]*h[k];      
        y[j+len2]+=t[2*j+k]*g[k];

最初の反復中に t 配列の境界を超え、次の反復中に配列の近似部分を超えます。t 配列の先頭に (f-1) 個の要素を追加する必要があります。

double *t=calloc(len+f-1, sizeof(double));
memcpy(&t[f], x, len*sizeof(double));        
for (int i=0; i<lev; i++) {
memset(t, 0, (f-1)*sizeof(double));        
memset(y, 0, len*sizeof(double));
int len2=len/2;
for (int j=0; j<len2; j++)           
    for (int k=0; k<f; k++) {          
        y[j]     +=t[2*j+k]*h[k];      
        y[j+len2]+=t[2*j+k]*g[k];      
    }
len=len2;                            
memcpy(&t[f], y, len*sizeof(double));    
}
于 2015-06-19T19:19:27.453 に答える