小波变换图像处理

应用小波变换对图像处理

图像金字塔，每一层的像素数为(2^i,2^j):图像金字塔基本操作 对一张图像不断的模糊之后向下采样,得到不同分辨率的图像,同时每次得到的 新的图像宽与高是原来图像的1/2，http://baike.baidu.com/link?url=uxUdGHjdQPia38hy6sElpAwpjkiyX0tUx0DH-I6JC2LrEWIG9ov8WKibm29zz60fiAE3rAiP0Jzd-YDDuk_VzK,

Matlab版本

clear all;close all;clc;tic;img=double(imread('yyt.jpg'));img=img(:,:,1);img=img(1:512,1:512);[m n]=size(img);[LL LH HL HH]=haar_dwt2D(img);img=[LL LH;HL HH];imgn=zeros(m,n);for i=0:m/2:m/2    for j=0:n/2:n/2        [LL LH HL HH]=haar_dwt2D(img(i+1:i+m/2,j+1:j+n/2)); %对一层分解后的四个图像分别再分解        imgn(i+1:i+m/2,j+1:j+n/2)=[LL LH;HL HH];      endendtoc;imshow(imgn)function [LL LH HL HH]=haar_dwt2D(img)    [m n]=size(img);    for i=1:m       %每一行进行分解        [L H]=haar_dwt(img(i,:));        img(i,:)=[L H];    end    for j=1:n       %每一列进行分解       [L H]=haar_dwt(img(:,j));       img(:,j)=[L H];    end    %本来分解不应该加mat2gray的，不过为了有好的显示效果就加上了    LL=mat2gray(img(1:m/2,1:n/2));          %行列都是低频      LH=mat2gray(img(1:m/2,n/2+1:n));        %行低频列高频    HL=mat2gray(img(m/2+1:m,1:n/2));        %行高频列低频  实现图像矩阵的归一化操作    HH=mat2gray(img(m/2+1:m,n/2+1:n));      %行列都是高频    endfunction [L H]=haar_dwt(f) %显然，我没有做边界处理，图片最好是2^n*2^n型的    n=length(f);    n=n/2;    L=zeros(1,n);   %低频分量    H=zeros(1,n);   %高频分量    for i=1:n        L(i)=(f(2*i-1)+f(2*i))/sqrt(2);        H(i)=(f(2*i-1)-f(2*i))/sqrt(2);    end    end

OPENCV版本

#include<stdio.h>#include "cv.h"#include "cxcore.h"#include "highgui.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 main( int argc, char** argv ) {// 小波变换层数int nLayer = 1;// 输入彩色图像IplImage *pSrc = cvLoadImage("G:\\me.jpg", 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);cvResetImageROI(pWavelet); // 本行代码有点多余，但有利用养成良好的编程习惯cvReleaseImage(&pWavelet);}// 显示图像pSrccvNamedWindow("pSrc", 1);   cvShowImage( "pSrc", pSrc);  cvWaitKey(0); //等待按键   cvDestroyWindow( "pSrc" );  cvReleaseImage(&pSrc);return 0;}