c++ - ソフトしきい値による 2d Haar ウェーブレット変換を使用した画像のノイズ除去 (C++)

Question

Haar 変換による基本的な 2D ウェーブレット変換の 1 つを実装しようとしています。

これを画像のノイズ除去問題に適用しました。

復元された結果には、いくつかの黒いブロックといくつかの白いブロックがあります。

正規化せずにソフトしきい値の部分にこだわったと思います。

これが私のコードです:

   #include "StdAfx.h"
   #include "WaveletDenoising.h"
   #include <cmath>


   WaveletDenoising::WaveletDenoising(void)
   {
   }

   WaveletDenoising::~WaveletDenoising(void)
   {
   }

   /* Forward Haar wavelet transform: */
   void WaveletDenoising::ForwardHaar1D(double* data, int length)
   {
    const float inv_sqrt2 = 1/sqrt((double)2.0);

    float norm = 1.0f/sqrt((double)length);

    for(int i=0; i < length; i++) {
        data[i] *= norm;
    }

    float *tmp = new float[length];

    while(length > 1) {
        length /= 2;

        for(int i=0; i < length; i++) {
            tmp[i] = (data[2*i] + data[2*i+1]) * inv_sqrt2;
            tmp[length + i] = (data[2*i] - data[2*i+1]) * inv_sqrt2;
        }

        memcpy(data, tmp, length*2*sizeof(float));
    }

    delete [] tmp;
   }

   /* Transpose matrix: */   
   void WaveletDenoising::Transpose(double *data, int width, int height)
   {
    double *B = new double[width*height];

    for(int y=0; y < height; y++) {
        for(int x=0; x < width; x++) {
            B[x*height + y] = data[y*width + x];
        }
    }

    memcpy(data, B, sizeof(double)*width*height);

    delete [] B;
   }

   /* Forward 2d Haar wavelet transform: */
   void WaveletDenoising::ForwardHaar2D(double* data, int width, int height)
   {
    for(int i=0; i < height; i++) 
        ForwardHaar1D(&data[i*width], width);

    Transpose(data, width, height);

    for(int i=0; i < width; i++) 
        ForwardHaar1D(&data[i*height], height);

    Transpose(data, height, width);
   }

   /* Inverse 1d Haar transform */
   void WaveletDenoising::InverseHaar1D(double* data, int length)
   {
    const float inv_sqrt2 = 1/sqrt((double)2.0);
    float inv_norm = sqrt((double)length);

    float *tmp = new float[length];
    int k = 1;

    while(k < length)  {
        for(int i=0; i < k; i++) {
            tmp[2*i] = (data[i] + data[k+i]) * inv_sqrt2;
            tmp[2*i+1] = (data[i] - data[k+i]) * inv_sqrt2;
        }

        memcpy(data, tmp, sizeof(double)*(k*2));

        k *= 2;
    }

    for(int i=0; i < length; i++) {
        data[i] *= inv_norm;
    }

    delete [] tmp;
   }

   /* Inverse 2d Haar wavelet transform */
   void WaveletDenoising::InverseHaar2D(double* data, int width, int height)
   {
    for(int i=0; i < width; i++) {
        InverseHaar1D(&data[i*height], height);
    }

    Transpose(data, height, width);

    for(int i=0; i < height; i++) {
        InverseHaar1D(&data[i*width], width);
    }

    Transpose(data, width, height);
   }

   /* Image denoising by soft-thresholding */
   void WaveletDenoising::WaveletThresholdDenoising(int width, int height, double* src,    double* des, double threshold)
   {
    int i, j, x, y;

    this->ForwardHaar2D(src, width, height);

    double mi = src[0*width+0]; /* find min value */
    double ma = src[0*width+0]; /* find max value */

    for (y=0; y<height; y++)
    {
        for (x=0; x<width; x++)
        {
            if (mi > src[y*width+x])
                mi = src[y*width+x];
            if (ma < src[y*width+x])
                ma = src[y*width+x];
        }
    }

    /* soft-thresholding */
    for (y=0; y<height; y++)
    {
        for (x=0; x<width; x++)
        {
            if (src[y*width+x] < threshold) 
                src[y*width+x] = 0;
            else if (src[y*width+x] > threshold)
                src[y*width+x] = src[y*width+x] - threshold;
            else if (src[y*width+x] < -threshold)
                src[y*width+x] = src[y*width+x] + threshold;


        }
    }

    this->InverseHaar2D(src, width, height);

    for (y=0; y<height; y++)
    {
        for (x=0; x<width; x++)
        {
            // for normalized:
            src[y*width+x] = ((src[y*width+x] - mi) / (ma - mi))*255;

            double temp =  src[y*width+x];
            if (temp < 0) temp = 0;
            else if (temp >255) temp = 255;
            else temp = temp;
            des[y*width+x] = (BYTE) src[y*width+x];
        }
    }
   }

これどうやってするの ？