另外参见俄罗斯写的http://www.codeproject.com/Articles/22243/Real-Time-Object-Tracker-in-C

Haar小波在图像处理和数字水印等方面应用较多，这里简单的介绍一下哈尔小波的基本原理以及其实现情况。

 一、Haar小波的基本原理

数学理论方面的东西我也不是很熟悉，这边主要用简单的例子来介绍下Haar小波的使用情况。

例如：有a=[8,7,6,9]四个数，并使用b[4]数组来保存结果.

则一级Haar小波变换的结果为:

b[0]=(a[0]+a[1])/2,                       b[2]=(a[0]-a[1])/2

b[1]=(a[2]+a[3])/2,                       b[3]=(a[2]-a[3])/2

即依次从数组中取两个数字，计算它们的和以及差，并将和一半和差的一半依次保存在数组的前半部分和后半部分。

例如：有a[8]，要进行一维Haar小波变换，结果保存在b[8]中

则一级Haar小波变换的结果为:

b[0]=(a[0]+a[1])/2,                        b[4]=(a[0]-a[1])/2

b[1]=(a[2]+a[3])/2,                        b[5]=(a[2]-a[3])/2

b[2]=(a[4]+a[5])/2,                        b[6]=(a[4-a[5]])/2

b[3]=(a[6]+a[7])/2,                        b[7]=(a[6]-a[7])/2

如果需要进行二级Haar小波变换的时候，只需要对b[0]-b[3]进行Haar小波变换.

对于二维的矩阵来讲，每一级Haar小波变换需要先后进行水平方向和竖直方向上的两次一维小波变换,行和列的先后次序对结果不影响。

二、Haar小波的实现

使用opencv来读取图片及像素，对图像的第一个8*8的矩阵做了一级小波变换

#include <cv.h>

#include <highgui.h>

#include <iostream>

using namespace std;

int main(){

IplImage* srcImg;

double  imgData[8][8];

int i,j;

srcImg=cvLoadImage(“lena.bmp”,0);

cout<<”原8*8数据”<<endl;

for( i=0;i<8;i++)

{

for( j=0;j<8;j++)

{

imgData[i][j]=cvGetReal2D(srcImg,i+256,j+16);

cout<<imgData[i][j]<<” “;

}

cout<<endl;

}

double tempData[8];

//行小波分解

for( i=0;i<8;i++)

{

for( j=0;j<4;j++)

{

double temp1=imgData[i][2*j];

double temp2=imgData[i][2*j+1];

tempData[j]=(temp1+temp2)/2;

tempData[j+4]=(temp1-temp2)/2;

}

for( j=0;j<8;j++)

imgData[i][j]=tempData[j];

}

//列小波分解

for( i=0;i<8;i++)

{

for( j=0;j<4;j++)

{

double temp1=imgData[2*j][i];

double temp2=imgData[2*j+1][i];

tempData[j]=(temp1+temp2)/2;

tempData[j+4]=(temp1-temp2)/2;

}

for( j=0;j<8;j++)

imgData[j][i]=tempData[j];

}

cout<<”1级小波分解数据”<<endl;

for( i=0;i<8;i++)

{

for( j=0;j<8;j++)

{

cout<<imgData[i][j]<<” “;

}

cout<<endl;

}

//列小波逆分解

for( i=0;i<8;i++)

{

for( j=0;j<4;j++)

{

double temp1=imgData[j][i];

double temp2=imgData[j+4][i];

tempData[2*j]=temp1+temp2;

tempData[2*j+1]=temp1-temp2;

}

for( j=0;j<8;j++)

{

imgData[j][i]=tempData[j];

}

}

//行小波逆分解

for( i=0;i<8;i++)

{

for( j=0;j<4;j++)

{

double temp1=imgData[i][j];

double temp2=imgData[i][j+4];

tempData[2*j]=temp1+temp2;

tempData[2*j+1]=temp1-temp2;

}

for( j=0;j<2*4;j++)

{

imgData[i][j]=tempData[j];

}

}

cout<<”1级小波逆分解数据”<<endl;

for( i=0;i<8;i++)

{

for( j=0;j<8;j++)

{

cout<<imgData[i][j]<<” “;

}

cout<<endl;

}

return 0;

}

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

///  小波变换
Mat WDT( const Mat &_src, const string _wname, const int _level )const
{
    int reValue = THID_ERR_NONE;
    Mat src = Mat_<float>(_src);
    Mat dst = Mat::zeros( src.rows, src.cols, src.type() ); 
    int N = src.rows;
    int D = src.cols;
 
    /// 高通低通滤波器
    Mat lowFilter; 
    Mat highFilter;
    wavelet( _wname, lowFilter, highFilter );
 
    /// 小波变换
    int t=1;
    int row = N;
    int col = D;
 
    while( t<=_level )
    {
        ///先进行行小波变换
        for( int i=0; i<row; i++ ) 
        {
            /// 取出src中要处理的数据的一行
            Mat oneRow = Mat::zeros( 1,col, src.type() );
            for ( int j=0; j<col; j++ )
            {
                oneRow.at<float>(0,j) = src.at<float>(i,j);
            }
            oneRow = waveletDecompose( oneRow, lowFilter, highFilter );
            /// 将src这一行置为oneRow中的数据
            for ( int j=0; j<col; j++ )
            {
                dst.at<float>(i,j) = oneRow.at<float>(0,j);
            }
        }
 
#if 0
        //normalize( dst, dst, 0, 255, NORM_MINMAX );
        IplImage dstImg1 = IplImage(dst); 
        cvSaveImage( "dst.jpg", &dstImg1 );
#endif
        /// 小波列变换
        for ( int j=0; j<col; j++ )
        {
            /// 取出src数据的一行输入
            Mat oneCol = Mat::zeros( row, 1, src.type() );
            for ( int i=0; i<row; i++ )
            {
                oneCol.at<float>(i,0) = dst.at<float>(i,j);
            }
            oneCol = ( waveletDecompose( oneCol.t(), lowFilter, highFilter ) ).t();
        
            for ( int i=0; i<row; i++ )
            {
                dst.at<float>(i,j) = oneCol.at<float>(i,0);
            }
        }
 
#if 0
        //normalize( dst, dst, 0, 255, NORM_MINMAX );
        IplImage dstImg2 = IplImage(dst); 
        cvSaveImage( "dst.jpg", &dstImg2 );
#endif
 
        /// 更新
        row /= 2;
        col /=2;
        t++;
        src = dst;
    }
 
    return dst;
}
 
///  小波逆变换
Mat IWDT( const Mat &_src, const string _wname, const int _level )const
{
    int reValue = THID_ERR_NONE;
    Mat src = Mat_<float>(_src);
    Mat dst = Mat::zeros( src.rows, src.cols, src.type() ); 
    int N = src.rows;
    int D = src.cols;
 
    /// 高通低通滤波器
    Mat lowFilter; 
    Mat highFilter;
    wavelet( _wname, lowFilter, highFilter );
 
    /// 小波变换
    int t=1;
    int row = N/std::pow( 2., _level-1);
    int col = D/std::pow(2., _level-1);
 
    while ( row<=N && col<=D )
    {
        /// 小波列逆变换
        for ( int j=0; j<col; j++ )
        {
            /// 取出src数据的一行输入
            Mat oneCol = Mat::zeros( row, 1, src.type() );
            for ( int i=0; i<row; i++ )
            {
                oneCol.at<float>(i,0) = src.at<float>(i,j);
            }
            oneCol = ( waveletReconstruct( oneCol.t(), lowFilter, highFilter ) ).t();
 
            for ( int i=0; i<row; i++ )
            {
                dst.at<float>(i,j) = oneCol.at<float>(i,0);
            }
        }
 
#if 0
        //normalize( dst, dst, 0, 255, NORM_MINMAX );
        IplImage dstImg2 = IplImage(dst); 
        cvSaveImage( "dst.jpg", &dstImg2 );
#endif
        ///行小波逆变换
        for( int i=0; i<row; i++ ) 
        {
            /// 取出src中要处理的数据的一行
            Mat oneRow = Mat::zeros( 1,col, src.type() );
            for ( int j=0; j<col; j++ )
            {
                oneRow.at<float>(0,j) = dst.at<float>(i,j);
            }
            oneRow = waveletReconstruct( oneRow, lowFilter, highFilter );
            /// 将src这一行置为oneRow中的数据
            for ( int j=0; j<col; j++ )
            {
                dst.at<float>(i,j) = oneRow.at<float>(0,j);
            }
        }
 
#if 0
        //normalize( dst, dst, 0, 255, NORM_MINMAX );
        IplImage dstImg1 = IplImage(dst); 
        cvSaveImage( "dst.jpg", &dstImg1 );
#endif
 
        row *= 2;
        col *= 2;
        src = dst;
    }
 
    return dst;
}
 
 
 
/// 调用函数
 
/// 生成不同类型的小波，现在只有haar，sym2
void wavelet( const string _wname, Mat &_lowFilter, Mat &_highFilter )const
{
    if ( _wname=="haar" || _wname=="db1" )
    {
        int N = 2;
        _lowFilter = Mat::zeros( 1, N, CV_32F );
        _highFilter = Mat::zeros( 1, N, CV_32F );
        
        _lowFilter.at<float>(0, 0) = 1/sqrtf(N); 
        _lowFilter.at<float>(0, 1) = 1/sqrtf(N); 
 
        _highFilter.at<float>(0, 0) = -1/sqrtf(N); 
        _highFilter.at<float>(0, 1) = 1/sqrtf(N); 
    }
    if ( _wname =="sym2" )
    {
        int N = 4;
        float h[] = {-0.483, 0.836, -0.224, -0.129 };
        float l[] = {-0.129, 0.224,    0.837, 0.483 };
 
        _lowFilter = Mat::zeros( 1, N, CV_32F );
        _highFilter = Mat::zeros( 1, N, CV_32F );
 
        for ( int i=0; i<N; i++ )
        {
            _lowFilter.at<float>(0, i) = l[i]; 
            _highFilter.at<float>(0, i) = h[i]; 
        }
 
    }
}
 
/// 小波分解
Mat waveletDecompose( const Mat &_src, const Mat &_lowFilter, const Mat &_highFilter )const
{
    assert( _src.rows==1 && _lowFilter.rows==1 && _highFilter.rows==1 );
    assert( _src.cols>=_lowFilter.cols && _src.cols>=_highFilter.cols );
    Mat &src = Mat_<float>(_src);
 
    int D = src.cols;
    
    Mat &lowFilter = Mat_<float>(_lowFilter);
    Mat &highFilter = Mat_<float>(_highFilter);
 
 
    /// 频域滤波，或时域卷积；ifft( fft(x) * fft(filter)) = cov(x,filter) 
    Mat dst1 = Mat::zeros( 1, D, src.type() );
    Mat dst2 = Mat::zeros( 1, D, src.type()  );
 
    filter2D( src, dst1, -1, lowFilter );
    filter2D( src, dst2, -1, highFilter );
 
 
    /// 下采样
    Mat downDst1 = Mat::zeros( 1, D/2, src.type() );
    Mat downDst2 = Mat::zeros( 1, D/2, src.type() );
 
    resize( dst1, downDst1, downDst1.size() );
    resize( dst2, downDst2, downDst2.size() );
 
 
    /// 数据拼接
    for ( int i=0; i<D/2; i++ )
    {
        src.at<float>(0, i) = downDst1.at<float>( 0, i );
        src.at<float>(0, i+D/2) = downDst2.at<float>( 0, i );
    }
 
    return src;
}
 
/// 小波重建
Mat waveletReconstruct( const Mat &_src, const Mat &_lowFilter, const Mat &_highFilter )const
{
    assert( _src.rows==1 && _lowFilter.rows==1 && _highFilter.rows==1 );
    assert( _src.cols>=_lowFilter.cols && _src.cols>=_highFilter.cols );
    Mat &src = Mat_<float>(_src);
 
    int D = src.cols;
 
    Mat &lowFilter = Mat_<float>(_lowFilter);
    Mat &highFilter = Mat_<float>(_highFilter);
 
    /// 插值;
    Mat Up1 = Mat::zeros( 1, D, src.type() );
    Mat Up2 = Mat::zeros( 1, D, src.type() );
 
    /// 插值为0
    //for ( int i=0, cnt=1; i<D/2; i++,cnt+=2 )
    //{
    //    Up1.at<float>( 0, cnt ) = src.at<float>( 0, i );     ///< 前一半
    //    Up2.at<float>( 0, cnt ) = src.at<float>( 0, i+D/2 ); ///< 后一半
    //}
 
    /// 线性插值
    Mat roi1( src, Rect(0, 0, D/2, 1) );
    Mat roi2( src, Rect(D/2, 0, D/2, 1) );
    resize( roi1, Up1, Up1.size(), 0, 0, INTER_CUBIC );
    resize( roi2, Up2, Up2.size(), 0, 0, INTER_CUBIC );
 
    /// 前一半低通，后一半高通
    Mat dst1 = Mat::zeros( 1, D, src.type() );
    Mat dst2= Mat::zeros( 1, D, src.type() );
    filter2D( Up1, dst1, -1, lowFilter );
    filter2D( Up2, dst2, -1, highFilter );
 
    /// 结果相加
    dst1 = dst1 + dst2;
 
    return dst1;
 
}

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

*****************************************************************************************************************

其他代码实现

*****************************************************************************************************************

// 哈尔小波.cpp : 定义控制台应用程序的入口点。

//

#include “stdafx.h”

#include <iostream>

#include “cv.h”

#include “highgui.h”

#include <ctype.h>

// 二维离散小波变换（单通道浮点图像）

void DWT(IplImage *pImage, int nLayer)

{

// 执行条件

if (pImage)

{

if (pImage->nChannels == 1 &&

pImage->depth == IPL_DEPTH_32F &&

((pImage->width >> nLayer) << nLayer) == pImage->width &&

((pImage->height >> nLayer) << nLayer) == pImage->height)

{

int     i, x, y, n;

float   fValue   = 0;

float   fRadius  = sqrt(2.0f);

int     nWidth   = pImage->width;

int     nHeight  = pImage->height;

int     nHalfW   = nWidth / 2;

int     nHalfH   = nHeight / 2;

float **pData    = new float*[pImage->height];

float  *pRow     = new float[pImage->width];

float  *pColumn  = new float[pImage->height];

for (i = 0; i < pImage->height; i++)

{

pData[i] = (float*) (pImage->imageData + pImage->widthStep * i);

}

// 多层小波变换

for (n = 0; n < nLayer; n++, nWidth /= 2, nHeight /= 2, nHalfW /= 2, nHalfH /= 2)

{

// 水平变换

for (y = 0; y < nHeight; y++)

{

// 奇偶分离

memcpy(pRow, pData[y], sizeof(float) * nWidth);

for (i = 0; i < nHalfW; i++)

{

x = i * 2;

pData[y][i] = pRow[x];

pData[y][nHalfW + i] = pRow[x + 1];

}

// 提升小波变换

for (i = 0; i < nHalfW - 1; i++)

{

fValue = (pData[y][i] + pData[y][i + 1]) / 2;

pData[y][nHalfW + i] -= fValue;

}

fValue = (pData[y][nHalfW - 1] + pData[y][nHalfW - 2]) / 2;

pData[y][nWidth - 1] -= fValue;

fValue = (pData[y][nHalfW] + pData[y][nHalfW + 1]) / 4;

pData[y][0] += fValue;

for (i = 1; i < nHalfW; i++)

{

fValue = (pData[y][nHalfW + i] + pData[y][nHalfW + i - 1]) / 4;

pData[y][i] += fValue;

}

// 频带系数

for (i = 0; i < nHalfW; i++)

{

pData[y][i] *= fRadius;

pData[y][nHalfW + i] /= fRadius;

}

}

// 垂直变换

for (x = 0; x < nWidth; x++)

{

// 奇偶分离

for (i = 0; i < nHalfH; i++)

{

y = i * 2;

pColumn[i] = pData[y][x];

pColumn[nHalfH + i] = pData[y + 1][x];

}

for (i = 0; i < nHeight; i++)

{

pData[i][x] = pColumn[i];

}

// 提升小波变换

for (i = 0; i < nHalfH - 1; i++)

{

fValue = (pData[i][x] + pData[i + 1][x]) / 2;

pData[nHalfH + i][x] -= fValue;

}

fValue = (pData[nHalfH - 1][x] + pData[nHalfH - 2][x]) / 2;

pData[nHeight - 1][x] -= fValue;

fValue = (pData[nHalfH][x] + pData[nHalfH + 1][x]) / 4;

pData[0][x] += fValue;

for (i = 1; i < nHalfH; i++)

{

fValue = (pData[nHalfH + i][x] + pData[nHalfH + i - 1][x]) / 4;

pData[i][x] += fValue;

}

// 频带系数

for (i = 0; i < nHalfH; i++)

{

pData[i][x] *= fRadius;

pData[nHalfH + i][x] /= fRadius;

}

}

}

delete[] pData;

delete[] pRow;

delete[] pColumn;

}

}

}

// 二维离散小波恢复（单通道浮点图像）

void IDWT(IplImage *pImage, int nLayer)

{

// 执行条件

if (pImage)

{

if (pImage->nChannels == 1 &&

pImage->depth == IPL_DEPTH_32F &&

((pImage->width >> nLayer) << nLayer) == pImage->width &&

((pImage->height >> nLayer) << nLayer) == pImage->height)

{

int     i, x, y, n;

float   fValue   = 0;

float   fRadius  = sqrt(2.0f);

int     nWidth   = pImage->width >> (nLayer - 1);

int     nHeight  = pImage->height >> (nLayer - 1);

int     nHalfW   = nWidth / 2;

int     nHalfH   = nHeight / 2;

float **pData    = new float*[pImage->height];

float  *pRow     = new float[pImage->width];

float  *pColumn  = new float[pImage->height];

for (i = 0; i < pImage->height; i++)

{

pData[i] = (float*) (pImage->imageData + pImage->widthStep * i);

}

// 多层小波恢复

for (n = 0; n < nLayer; n++, nWidth *= 2, nHeight *= 2, nHalfW *= 2, nHalfH *= 2)

{

// 垂直恢复

for (x = 0; x < nWidth; x++)

{

// 频带系数

for (i = 0; i < nHalfH; i++)

{

pData[i][x] /= fRadius;

pData[nHalfH + i][x] *= fRadius;

}

// 提升小波恢复

fValue = (pData[nHalfH][x] + pData[nHalfH + 1][x]) / 4;

pData[0][x] -= fValue;

for (i = 1; i < nHalfH; i++)

{

fValue = (pData[nHalfH + i][x] + pData[nHalfH + i - 1][x]) / 4;

pData[i][x] -= fValue;

}

for (i = 0; i < nHalfH - 1; i++)

{

fValue = (pData[i][x] + pData[i + 1][x]) / 2;

pData[nHalfH + i][x] += fValue;

}

fValue = (pData[nHalfH - 1][x] + pData[nHalfH - 2][x]) / 2;

pData[nHeight - 1][x] += fValue;

// 奇偶合并

for (i = 0; i < nHalfH; i++)

{

y = i * 2;

pColumn[y] = pData[i][x];

pColumn[y + 1] = pData[nHalfH + i][x];

}

for (i = 0; i < nHeight; i++)

{

pData[i][x] = pColumn[i];

}

}

// 水平恢复

for (y = 0; y < nHeight; y++)

{

// 频带系数

for (i = 0; i < nHalfW; i++)

{

pData[y][i] /= fRadius;

pData[y][nHalfW + i] *= fRadius;

}

// 提升小波恢复

fValue = (pData[y][nHalfW] + pData[y][nHalfW + 1]) / 4;

pData[y][0] -= fValue;

for (i = 1; i < nHalfW; i++)

{

fValue = (pData[y][nHalfW + i] + pData[y][nHalfW + i - 1]) / 4;

pData[y][i] -= fValue;

}

for (i = 0; i < nHalfW - 1; i++)

{

fValue = (pData[y][i] + pData[y][i + 1]) / 2;

pData[y][nHalfW + i] += fValue;

}

fValue = (pData[y][nHalfW - 1] + pData[y][nHalfW - 2]) / 2;

pData[y][nWidth - 1] += fValue;

// 奇偶合并

for (i = 0; i < nHalfW; i++)

{

x = i * 2;

pRow[x] = pData[y][i];

pRow[x + 1] = pData[y][nHalfW + i];

}

memcpy(pData[y], pRow, sizeof(float) * nWidth);

}

}

delete[] pData;

delete[] pRow;

delete[] pColumn;

}

}

}

int _tmain(int argc, _TCHAR* argv[])

{

// 小波变换层数

int nLayer = 1;

// 输入彩色图像

IplImage *pSrc = cvLoadImage(“lena.bmp”, CV_LOAD_IMAGE_COLOR);

// 计算小波图象大小

CvSize size = cvGetSize(pSrc);

if ((pSrc->width >> nLayer) << nLayer != pSrc->width)

{

size.width = ((pSrc->width >> nLayer) + 1) << nLayer;

}

if ((pSrc->height >> nLayer) << nLayer != pSrc->height)

{

size.height = ((pSrc->height >> nLayer) + 1) << nLayer;

}

// 创建小波图象

IplImage *pWavelet = cvCreateImage(size, IPL_DEPTH_32F, pSrc->nChannels);

if (pWavelet)

{

// 小波图象赋值

cvSetImageROI(pWavelet, cvRect(0, 0, pSrc->width, pSrc->height));

cvConvertScale(pSrc, pWavelet, 1, -128);

cvResetImageROI(pWavelet);

// 彩色图像小波变换

IplImage *pImage = cvCreateImage(cvGetSize(pWavelet), IPL_DEPTH_32F, 1);

if (pImage)

{

for (int i = 1; i <= pWavelet->nChannels; i++)

{

cvSetImageCOI(pWavelet, i);

cvCopy(pWavelet, pImage, NULL);

// 二维离散小波变换

DWT(pImage, nLayer);

// 二维离散小波恢复

// IDWT(pImage, nLayer);

cvCopy(pImage, pWavelet, NULL);

}

cvSetImageCOI(pWavelet, 0);

cvReleaseImage(&pImage);

}

// 小波变换图象

cvSetImageROI(pWavelet, cvRect(0, 0, pSrc->width, pSrc->height));

cvConvertScale(pWavelet, pSrc, 1, 128);

cvNamedWindow(“pWavelet”,CV_WINDOW_AUTOSIZE);

cvShowImage(“pWavelet”,pWavelet);

cvNamedWindow(“pSrc”,CV_WINDOW_AUTOSIZE);

cvShowImage(“pSrc”,pSrc);

cvWaitKey(0);

cvResetImageROI(pWavelet); // 本行代码有点多余，但有利用养成良好的编程习惯

cvReleaseImage(&pWavelet);

}

// 显示图像pSrc

// …

cvReleaseImage(&pSrc);

return 0;

}

************************

————————————-Codeproject Haar ————————————————————————————————-

void Haar::transrows(char** dest, char** sour, unsigned int w, unsigned int h) const

{

unsigned int w2 = w / 2;

__m64 m00FF;

m00FF.m64_u64 = 0x00FF00FF00FF00FF;

for (unsigned int y = 0; y < h; y++) {

__m64 *mlo = (__m64 *) & dest[y][0];

__m64 *mhi = (__m64 *) & dest[y][w2];

__m64 *msour = (__m64 *) & sour[y][0];

for (unsigned int k = 0; k < w2 / 8; k++) {   //k<w2/8   k=8*k

__m64 even = _mm_packs_pu16(_mm_and_si64(*msour, m00FF), _mm_and_si64(*(msour + 1), m00FF));       //even coeffs

__m64 odd = _mm_packs_pu16(_mm_srli_pi16(*msour, 8), _mm_srli_pi16(*(msour + 1), 8));              //odd coeffs

addsub(even, odd, mlo++, mhi++);

msour += 2;

}

if (w2 % 8) {

for (unsigned int k = w2 - (w2 % 8); k < w2; k++) {

dest[y][k] = char(((int)sour[y][2*k] + (int)sour[y][2*k+1]) / 2);

dest[y][k+w2] = char(((int)sour[y][2*k] - (int)sour[y][2*k+1]) / 2);

}

}

}

_mm_empty();

}

void Haar::transcols(char** dest, char** sour, unsigned int w, unsigned int h) const

{

unsigned int h2 = h / 2;

for (unsigned int k = 0; k < h2; k++) {

__m64 *mlo = (__m64 *) & dest[k][0];

__m64 *mhi = (__m64 *) & dest[k+h2][0];

__m64 *even = (__m64 *) & sour[2*k][0];

__m64 *odd = (__m64 *) & sour[2*k+1][0];

for (unsigned int x = 0; x < w / 8; x++) {

addsub(*even, *odd, mlo, mhi);

even++;

odd++;

mlo++;

mhi++;

}

}

_mm_empty();

//odd remainder

for (unsigned int x = w - (w % 8); x < w; x++) {

for (unsigned int k = 0; k < h2; k++) {

dest[k][x] = char(((int)sour[2*k][x] + (int)sour[2*k+1][x]) / 2);

dest[k+h2][x] = char(((int)sour[2*k][x] - (int)sour[2*k+1][x]) / 2);

}

}

}

****************************************************************************************************************

哈尔小波变换示例

1.哈尔基函数
  最简单的基函数是哈尔基函数(Haar basis function)。哈尔基函数在1909年提出，它是由一组分段常值函数组成的函数集。这个函数集定义在半开区间[0,1)上，每一个分段常值函数的数 值在一个小范围里是1，其他地方为0，现以图像为例并使用线性代数中的矢量空间来说明哈尔基函数。
  如果一幅图像仅由2^0=1个像素组成，这幅图像在整个[0,1) 区间中就是一个常值函数。

用q_00(x)表示这个常值函数，用V0表示由这个常值函数生成的矢量空间，

即V0：q_00(x)=1(0<=x<1)或0(其它),它是构成矢量空间V0的基。

如果一幅图像由2^1=2个像素组成，这幅图像在[0,1) 区间中有两个等间隔的子区间：[0,1/2) 和[1/2,1) ，每一个区间中各有1个常值函数，分别用q_10(x) q_11(x)表示。用V1表示由2个子区间中的常值函数生成的矢量空间，即V1: q_10(x)  =1(0<=x<1/2)或0(其它); q_11(x)=1(1/2<=x<1)或0(其它).这2个常值函数就是构成矢量空间V1的基。

如果一幅图像由2^2= 4个像素组成，这幅图像在[0,1)区间中被分成4个等间隔的子区间：[0,1/4),[1/4,1/2),[1/2,3/4)和[3/4,1)，它们的 常值函数分别用q_20(x) q_21(x) q_22(x) q_23(x)表示，用V2表示由4个子区间中的常值函数生成的矢量空间，即V2:q_20(x)=1(0<=x<1/4)或0(其它); q_21(x)=1(1/4<=x<1/2)或0(其它); q_22(x) =1(1/2<=x<3/4)或0(其它); q_23(x)=1(3/4<=x<1)或0(其它).这4个常值函数就是构成矢量空间V2的基。

我们可以按照这种方法继续定义基函数和由它生成的矢量空间。总而言之，为了表示矢量空间中的矢量，每一个矢量空间V_i都需要定义一个基(basis)。 为生成矢量空间V_i而定义的基函数也叫做尺度函数,这种函数通常用符号q_ij(x)表示。哈尔基函数定义为q(x)=1(0<=x<1) 或0(其它),哈尔尺度函数q_ij(x)定义为q_ij(x)=q(2^i·x-j),j=0,1,…, 2^i-1其中,i为尺度因子，改变i使函数图形缩小或者放大,j为平移参数，改变j使函数沿x轴方向平移。
  空间矢量V_i定义为V_i=span{ q_ij(x)},j=0,1,…, 2^i-1.由于定义了基和矢量空间，我们就可以把由2^i个像素组成的一维图像看成为矢量空间V_i中的矢量。由于这些矢量都是在单位区间[0,1)上 定义的函数，所以在V_i矢量空间中的每一个矢量也被包含在V_i+1矢量空间中，即V0∈V1∈…∈V_i∈V_i+1.

2.哈尔小波函数
  小波函数通常用Ψ_ij(x)表示。与哈尔基函数相对应的小波称为基本哈尔小波函数,并由下式定义:

Ψ(x)=1(0<=x<1/2)或 -1(1/2<=x<1)或0(其它); 
哈尔小波尺度函数Ψ_ij(x)定义为
Ψ_ij(x)= Ψ(2^i·x-j),j=0,1,…, 2^i-1.用小波函数构成的矢量空间Wi表示，Wi=span{Ψ_ij(x)}, j=0,1,…, 2^i-1
  根据哈尔小波函数的定义，可以写出生成W0，W1和W2等矢量空间的小波函数。
生成矢量空间W0的哈尔小波：Ψ_00(x)=1(0<=x<1/2)或-1(1/2<=x<1)或0(其它)
生成矢量空间W1的哈尔小波：Ψ_10(x)=1(0<=x<1/4)或-1(1/4<=x<1/2)或0(其它);
Ψ_11(x)=1(1/2<=x<3/4)或-1(3/4<=x<1/2)或0(其它).
生成矢量空间W2的哈尔小波：Ψ_00(x)=1(0<=x<1/2)或-1(1/2<=x<1)或0(其它)
Ψ_20(x)=1(0<=x<1/8)或-1(1/8<=x<2/8)或0(其它);
Ψ_21(x)=1(2/8<=x<3/8)或-1(3/8<=x<4/8)或0(其它);
Ψ_22(x)=1(4/8<=x<5/8)或-1(5/8<=x<6/8)或0(其它);
Ψ_00(x)=1(6/8<=x<7/8)或-1(7/8<=x<1)或0(其它).
哈尔基和哈尔小波分别使用下面两个式子进行规范化
q_ij(x)=2^(i/2)·q(2^i·x-j),
Ψ_ij(x)= 2^(i/2)·Ψ(2^i·x-j),
3.哈尔基的结构
  使用哈尔基函数 q_ij(x)和哈尔小波函数Ψ_ij(x)生成的矢量空间Vi和Wi具有下面的性质，V_i+1=V_i⊕W_i。这就是说，在矢量空间V_i+1中, 矢量空间W_i中的小波可用来表示一个函数在矢量空间V_i 中不能表示的部分。因此，小波可定义为生成矢量空间W_i的一组线性无关的函数Ψ_ij(x)的集合。
4.哈尔小波变换
  小波变换的基本思想是用一组小波函数或者基函数表示一个函数或者信号，例如图像信号。为了理解什么是小波变换，下面用一个具体的例子来说明小波变换的过程。假设有一幅分辨率只有4个像素的一维图像，对应的像素值分别为[9 7 3 5]
  计算它的哈尔小波变换系数:
(1).求均值(averaging)。计算相邻像素对的平均值，得到一幅分辨率比较低的新图像，它的像素数目变成了2个，即新的图像的分辨率是原来的1/2，相应的像素值为：[8 4]
(2). 求差值(diqqerencing)。很明显，用2个像素表示这幅图像时，图像的信息已经部分丢失。为了能够从由2个像素组成的图像重构出由4个像素组成 的原始图像，就需要存储一些图像的细节系数,以便在重构时找回丢失的信息。方法是把像素对的第一个像素值减去这个像素对的平均值，或者使用这个像素对的差 值除以2。在这个例子中，第一个细节系数是(9-8)=1，因为计算得到的平均值是8,它比9小1而比7大1，存储这个细节系数就可以恢复原始图像的前两 个像素值。使用同样的方法，第二个细节系数是(3-4)=-1，存储这个细节系数就可以恢复后2个像素值。因此，原始图像就可以用下面的两个平均值和两个 细节系数表示[8 4 1 -1]
(3).重复第1，2步，把由第一步分解得到的图像进一步分解成分辨率更低的图像[6 2 1 -1].
  由此可见，通过上述分解就把由4像素组成的一幅图像用一个平均像素值和三个细节系数表示，这个过程就叫做哈尔小波变换,也称哈尔小波分解(Haar wavelet decomposition)。这个概念可以推广到使用其他小波基的变换。在这个例子中我们可以看到，①变换过程中没有丢失信息，因为能够从所记录的数据 中重构出原始图像。②对这个给定的变换，我们可以从所记录的数据中重构出各种分辨率的图像。例如，在分辨率为1的图像基础上重构出分辨率为2的图像，在分 辨率为2的图像基础上重构出分辨率为4的图像。③通过变换之后产生的细节系数的幅度值比较小，这就为图像压缩提供了一种途径，例如去掉一些微不足道的细节 系数而不影响对重构图像的理解。
  求均值和差值的过程实际上就是一维小波变换的过程，现在用数学方法重新描述小波变换的过程。
(1) 用V2中的哈尔基表示
图像I(x)=[9 7 3 5]有2^j =2^2=4像素，因此可以用生成矢量空间V2中的基函数的线性组合表示，
I(x) =9 q_20(x)+ 7q_21(x)+3q_22(x)+5q_23(x)

(2)用V1和W1中的函数表示
根据V2=V1⊕W1,I(x)可表示成
I(x)=8q_10(x)+4q_11(x)+Ψ_10(x)-Ψ_11(x)

(3)用V0,W0和W1中的函数表示
V2=V0⊕W0⊕W1,I(x)可表示成
I(x) =6 q_00(x)+2Ψ_00(x)+Ψ_10(x)-Ψ_11(x)

5.规范化算法
  规范化的小波变换与非规范化的小波变换相比，唯一的差别是在规范化的小波变换中需要乘一个规范化的系数。下面用一个例子说明。对函数f(x)= [2,5,8,9,7,4, -1, 1]作哈尔小波变换。哈尔小波变换实际上是使用求均值和差值的方法进行分解。我们把f(x)看成是矢量空间V3中的一个向量，尺度因子i= 3，因此最多可分解为3层.分解过程如下。
步骤1：c=(2 5,8 9,7,4,-1,1)/√2=(7,17,11,0, -3, -1,3, -2)/√2
解释：
平均值（7,17,11,0）/2 
差值（-3，-1,3，-2）/2（与平均值的差值相对于直接差值的1/2  即 A- (A+B)/2 = (A-B)/2  ,我们将除以2提出来）
所以步骤一的结果是（7,17,11,0）/2  并上（-3，-1,3，-2）/2 = (7,17,11,0, -3, -1,3, -2)/2
步骤2：c=[(7+17)/√2,(11+0)/√2,(7-17)/√2,(11-0)/√2,-3,-1,3,-2]/√2
=(24/√2,11/√2,-10/√2,11/√2,-3,-1,3,-2)/√2
步骤3:c=[(24+11)/2,(24-11)/2 ,-10/√2,11/√2,-3,-1,3,-2]/√2
=(35/2,13/2,-10/√2,11/√2,-3,-1,3,-2)/√2
这个结果实际上是后面哈尔转换N=8转换矩阵相乘的结果。
N = 8，

Haar小波变换基本原理相关推荐

  1. SSE2实现HAAR小波变换（dwt2与idwt2）

    wiki链接:http://en.wikipedia.org/wiki/Haar_wavelet 可用SSE2实现HAAR小波变换,达到实时,关于HAAR小波的介绍可参考以上维基链接 参考MATLAB ...

  2. Haar小波变换的快速实现

    Haar小波变换的快速实现 2014年3月12日renjihe发表评论阅读评论 先举个例子,有a=[100,12,43,39]四个数,并使用b[4]数组来保存结果. 一级Haar小波变换的结果为: b ...

  3. 使用CUDA计算Haar小波变换

    在<Haar小波变换的快速实现>一文里我们提到了Haar小波变换的计算,在这里我们使用CUDA实现文中提到的计算方式. 01 __global__ void 02 _cuda_haar(f ...

  4. haar小波变换学习笔记

    本篇很大一部分内容借鉴了篇末所引用的优质博客~~ 小波可以认为是一个带通滤波器,只允许频率和小波基函数频率相近的信号通过.小波变换的基本思想是用一组小波函数和基函数表示一个函数或者信号. haar小波 ...

  5. 二维图像haar小波变换的分解与重构

    二维图像haar小波变换的分解与重构 二维离散小波的理论推导和一维小波类似,但是以其尺度函数生成的尺度函数集作为标准正交基的尺度空间Vi的正交补空间Wi不能直接得到,而是可以证明,正交补空间Wi是由三 ...

  6. 图像Haar小波变换

    说起小波变换就需要提起傅里叶变换.傅里叶变换就是把波进行分解,可以认为任意一个周期波都可以有足够多的正弦(余弦)波组成,这里足够多的正弦波对应的频率不同,把这些足够的正弦波放在频域中,就是傅里叶变换, ...

  7. 哈尔(Haar)小波变换的原理及opencv源代码

    1. 小波分析 小波分析是对傅里叶变换的继承,总结和重大突破.小波分析的优势在于可以同时进行时频域分析,比傅里叶变换更适合处理非平稳信号. 小波分析所用的波称为小波,小波的能量有限,有限长且会衰减,集 ...

  8. 一维的Haar小波变换

    本文转载自:http://blog.csdn.net/liulina603/article/details/8649339 小波变换的基本思想是用一组小波函数或者基函数表示一个函数或者信号,例如图像信 ...

  9. Haar小波变换代码实现

    代码1:以图像的形式显示. # include<opencv2/opencv.hpp> # include<iostream> using namespace std; usi ...

最新文章

  1. MATLAB_10-模式识别_
  2. 深度强化学习入门介绍
  3. tsp 选边 matlab,【转载】蚁群算法TSP(旅行商问题)通用matlab程序
  4. vue调用顺序(初学版) index.html → main.js → app.vue → index.js → components/组件 测试
  5. 计算机一级查询记录,技巧查看电脑中使用过的记录痕迹的详细教程
  6. 2021年，不平凡的一年~
  7. 设置域用户登录主目录
  8. java 获取本机的IP和hostname
  9. css中的.clearfix是什么意思？
  10. 异常处理 Exceptions
  11. 系列4—BabeLua常见问题
  12. 《Qt 实战一二三》
  13. 壳聚糖/纳米金水凝胶/纳米木质素/掺杂二硫化钼/微米级Ag2O2掺杂壳聚糖水凝胶的制备研究
  14. java 打印 边距_缩小边距 – Java打印
  15. CHD-5.3.6集群上sqoop安装
  16. Arun Jaitley：要健康最好让银行保持增长势头
  17. android 美团拆包,如何进行拆包合包操作？
  18. 数据传输完整性_基于IBIS模型的FPGA信号完整性仿真验证方法
  19. python股票网格交易_不知道哪位交易员可以解释下网格交易法？
  20. Hashmap底层源码

热门文章

  1. 哪些耳机音质不错？比苹果耳机音质好的耳机推荐
  2. bash 剪切文件_linux 剪切命令
  3. python Django session/cookie
  4. 20190328学习笔记 - JSP 中的 tag 文件
  5. ramos一键处理多合一_ramos一键处理多合一_Primo Ramdisk驱RAMOS一键制作生成工具多合一 V3.2.5版...
  6. Vue中html全屏背景色（height：100%不能生效时）
  7. python用turtle画月亮_使用Python turtle画表白分形树
  8. Guns 技术文档 v5.1
  9. linux操作系统的文件名最大长度为,linux和windows文件名长度限制
  10. MySQL数据库创建索引的方法和好处