小波变换笔记

实现haar小波

参考文章中讲解了haar小波的原理，比较易懂，实现了多级haar小波分解，本文参照它进行了练习，同时添加了小波重建。

# include<iostream>
# include<opencv2/opencv.hpp>

using namespace std;
using namespace cv;

int main()
{
	Mat src = imread("F:\\testdata\\pic\\16.jpg");
	Mat HSV;
	cvtColor(src, HSV, CV_BGR2HSV);
	vector<Mat> channels;
	split(HSV, channels);

	Mat img = channels.at(1);
        //上述是获取HSV中的饱和度图像，按需修改
	namedWindow("img", 0);
	imshow("img", img);
	int height = img.rows;
	int width = img.cols;
	int depth = 2;
	int depthcount = 1;

	Mat tmp = Mat::ones(Size(width, height), CV_32FC1);
	Mat wavelet = Mat::ones(Size(width, height), CV_32FC1);
	Mat imgtmp = img.clone();

	imgtmp.convertTo(imgtmp, CV_32FC1);

	//---------------------------------------小波分解-----------------------------//
	while (depthcount <= depth)
	{
		height = img.rows / pow(2,depthcount-1);
		width = img.cols / pow(2,depthcount-1);

		//水平方向变换
		for (int i = 0; i < height; ++i)
		{
			for (int j = 0; j < width / 2; ++j)
			{
				tmp.at<float>(i, j) = (imgtmp.at<float>(i, 2 * j) + imgtmp.at<float>(i, 2 * j + 1)) / 2;
				tmp.at<float>(i, j+width/2) = (imgtmp.at<float>(i, 2 * j) - imgtmp.at<float>(i, 2 * j + 1)) / 2;
			}
		}
		
		//垂直方向变换
		for (int i = 0; i < height / 2; ++i)
		{
			for (int j = 0; j < width; ++j)
			{
				wavelet.at<float>(i, j) = (tmp.at<float>(2 * i, j) + tmp.at<float>(2 * i + 1, j)) / 2;
				wavelet.at<float>(i+height/2, j) = (tmp.at<float>(2 * i, j) - tmp.at<float>(2 * i + 1, j)) / 2;
			}
		}
		/*
		//低通滤波，选用高斯低通滤波器
		Mat ROI = wavelet(Rect(0, 0, width, height));
		GaussianBlur(ROI, ROI, Size(7, 7), 0.5);
		Mat dst(wavelet, Rect(0, 0, width, height));
		ROI.copyTo(dst);
		*/
		imgtmp = wavelet;
		depthcount++;


	}

	
	//------------------------------------------小波重建--------------------------------------//
	while (depth > 0)
	{
		height = img.rows / pow(2, depth - 1);
		width = img.cols / pow(2, depth - 1);
		//列逆变换
		for (int i = 0; i < height/2; ++i)
		{
			for (int j = 0; j < width; ++j)
			{
				float value1 = imgtmp.at<float>(i, j);
				float value2 = imgtmp.at<float>(i+height/2, j);
				tmp.at<float>(2 * i, j) = value1 + value2;
				tmp.at<float>(2 * i+1, j) = value1 - value2;
			}
		}
		
		//行逆变换
		for (int i = 0; i < height; ++i)
		{
			for (int j = 0; j < width / 2; ++j)
			{
				float value1 = tmp.at<float>(i, j);
				float value2 = tmp.at<float>(i , j+width/2);
				wavelet.at<float>(i, 2*j) = value1 + value2;
				wavelet.at<float>(i , 2*j+1) = value1 - value2;
			}
		}

		Mat ROI=wavelet(Rect(0, 0, width, height));
		Mat dst(imgtmp, Rect(0, 0, width, height));
		ROI.copyTo(dst);
		depth--;
	}



	namedWindow("res", 0);
	imgtmp.convertTo(imgtmp, CV_8UC1);
	imshow("res", imgtmp);
	waitKey(0);
	return 0;
}

---

开源代码实现多种小波函数的小波变换

• C++ wavelet library
• 可以选择按照其文档中的配置方法，调用动态链接库，但是我的尝试不成功，出现了内存异常，所以本文直接将其源码放入自己的项目中，和自己的项目一起编译，可以实现同样的功能。需要注意的是，需要配置好FFTW这个依赖，方法见VS2015配置FFTW
• 函数说明
• 利用这个库中的函数可以得到小波分解的近似系数和细节系数

使用示例

重要的只有如下两行函数调用

dwt_2d(data, depth, name, output, flag, length);

idwt_2d(output, flag, name, res, length);

#include <iostream>
#include <fstream>
#include "wavelet2d.h"
#include <vector>
#include <string>
#include <cmath>
#include<opencv2\opencv.hpp>
using namespace std;
using namespace cv;

int main() {
  //读入图片，获取饱和度
  Mat src = imread("F:\\testdata\\poreImgV3\\10.jpg");
  resize(src, src, Size(450, 450));
  Mat HSV;
  cvtColor(src, HSV, CV_BGR2HSV);
  vector<Mat> channels;
  channels.resize(3);
  split(HSV, channels);
  Mat saturation = channels.at(1);
  namedWindow("s", 0);
  imshow("s", saturation);
  waitKey(0);
  //将饱和度图像数据存入vector<vector<double>>中,一行一行放置
  int row = saturation.rows;
  int col = saturation.cols;
  //分配大小
  vector<vector<double>> data;
  data.resize(row);
  for (int i = 0; i < data.size(); ++i)
  {
  	data[i].resize(col);
  }
  //放入数据
  for (int i = 0; i < row; ++i)
  {
  	for (int j = 0; j < col; ++j)
  	{
  		data[i][j] = saturation.at<uchar>(i, j);
  	}
  }

  int depth = 2;
  string name = "haar";
  vector<double> output;
  //output.resize(row*col);
  vector<double> flag;
  vector<int> length; 
  //小波分解
  dwt_2d(data, depth, name, output, flag, length);

  vector<vector<double>> res;
  res.resize(row);
  for (int i = 0; i < res.size(); ++i)
  {
  	res[i].resize(col);
  }
  //----------------------------------------------------------------------------------

  //小波重建
  idwt_2d(output, flag, name, res, length);
  //写回原图像中
  Mat recovery(src.size(), CV_8UC1);
  for (int i = 0; i < row; ++i)
  {
  	for (int j = 0; j < col; ++j)
  	{
  		recovery.at<uchar>(i, j) = (uchar)res[i][j];
  	}
  }
  imshow("res", recovery);
  waitKey(0);
  system("pause");
  return 0;
}