特征变换（2）沃尔什-哈达玛变换

笔记->印象笔记（沃尔什变换篇）

MATLAB 源码：

clear all;
I=zeros(2.^8);
I(2.^7-2.^4+1:2.^7+2.^4,2.^7-2.^4+1:2.^7+2.^4)=ones(2*2.^4);  %读入数据
subplot(1,2,1);
colormap(gray(128));imagesc(I);
title('显示原始图像');
[m,n]=size(I);
for k=1:n
    wht(:,k)=hadamard(m)*I(:,k)/m;  %对每一列作walah-Hadamard变换
end
for j=1:m
%对进行列的walah-Hadamard变换后的系数进行walah-Hadamard变换
    wh(:,j)=hadamard(n)*wht(j,:)'/n;
end
wh=wh';
subplot(1,2,2);colormap(gray(128));imagesc(wh);
title('walah-Hadamard变换系数');

C++源码：

/*************************************************************************
 *
 * 函数名称：
 *   WALSH()
 *
 * 参数:
 *   double * dpf				- 指向时域值的指针
 *   double * dpF				- 指向频域值的指针
 *   r						－2的幂数
 *
 * 返回值:
 *   无。
 *
 * 说明:
 *   该函数用来实现一维快速沃尔什-哈达玛变换。
 *
 ***********************************************************************
*/

VOID WINAPI WALSH(double *dpf, double *dpF, int r)
{
	// 沃尔什-哈达玛变换点数
	LONG	lNum;
	
	// 快速沃尔什变换点数
	lNum = 1 << r;

	// 循环变量
	int		i,j,k;
	
	// 中间变量
	int		nTemp,m;
	
	double *X1,*X2,*X;
			
	// 分配运算所需的数组
	X1 = new double[lNum];
	X2 = new double[lNum];
	
	// 将时域点写入数组X1
	memcpy(X1, dpf, sizeof(double) * lNum);
		
	for(k = 0; k < r; k++)
	{
		for(j = 0; j < 1<<k; j++)
		{	
			// 按照蝶形运算图进行运算
			nTemp = 1 << (r-k);
			for(i = 0; i < nTemp / 2; i++)
			{
				m = j * nTemp;
				X2[i + m] = X1[i + m] + X1[i + m + nTemp / 2];
				X2[i + m + nTemp / 2] = X1[i + m] - X1[i + m + nTemp / 2];
			}
		}
		
		// 互换  
		X = X2;
		X2 = X1;
		X1 = X;
	}
	
	// 对系数做调整
	for(j = 0; j < lNum; j++)
	{
		m = 0;
		for(i = 0; i < r; i++)
		{
			if (j & (1<<i))
			{
				m += 1 << (r-i-1);
			}
		}

		dpF[j] = X1[m] / lNum;
	}
	
	// 释放内存
	delete X1;
	delete X2;
}

/*************************************************************************
 *
 * 函数名称：
 *   IWALSH()
 *
 * 参数:
 *   double * dpF				- 指向频域值的指针
 *   double * dpf				- 指向时域值的指针
 *   n						－2的幂数
 *
 * 返回值:
 *   无。
 *
 * 说明:
 *   该函数用来实现一维快速沃尔什-哈达玛反变换。
 *
 ***********************************************************************
 */

VOID WINAPI IWALSH(double *dpF, double *dpf, int n)
{
	// 变换点数
	LONG	lNum;
	
	// 循环变量
	int		i;
	
	// 计算变换点数
	lNum = 1 << n;
	
	// 用快速沃尔什－哈达玛变换进行反变换
	WALSH(dpF, dpf, n);
	
	// 对系数进行调整
	for(i = 0; i < lNum; i++)
	{
		dpf[i] *= lNum;
	}
}

/*************************************************************************
 *
 * 函数名称：
 *   DIBWalsh()
 *
 * 参数:
 *   CDib  *pDib       - 指向CDib类的指针
 *
 * 返回值:
 *   BOOL               - 成功返回TRUE，否则返回FALSE。
 *
 * 说明:
 *   该函数用来对图像进行二维快速沃尔什-哈达玛变换。
 *
 ***********************************************************************
 */

BOOL WINAPI DIBWalsh(CDib * pDib)
{
	// 指向源图像的指针
	unsigned char *lpSrc;		
	
	//图象的宽度和高度
	LONG    lWidth;
	LONG    lHeight;

	// 循环变量
	LONG	i;
	LONG	j;
	
	// 实际进行付立叶变换的宽度和高度
	LONG	lW = 1;
	LONG	lH = 1;
	
	int		wp = 0;
	int		hp = 0;
	
	// 中间变量
	double	dTemp;	

	//得到图象的宽度和高度
	CSize   SizeDim;
	SizeDim = pDib->GetDimensions();
	lWidth = SizeDim.cx;
	lHeight = SizeDim.cy;	
	
	//得到实际的图象存储大小
	CSize   SizeRealDim;
	SizeRealDim = pDib->GetDibSaveDim();

	// 图像每行的字节数
	LONG	lLineBytes;
	
	// 计算图像每行的字节数
	lLineBytes = SizeRealDim.cx;

	//图像数据的指针
	LPBYTE  lpDIBBits = pDib->m_lpImage;
		
	// 保证离散余弦变换的宽度和高度为2的整数次方
	while(lW * 2 <= lWidth)
	{
		lW = lW * 2;
		wp++;
	}
	
	while(lH * 2 <= lHeight)
	{
		lH = lH * 2;
		hp++;
	}
	
	// 分配内存
	double *dpf = new double[lW * lH];
	double *dpF = new double[lW * lH];
	
	// 时域赋值
	for(i = 0; i < lH; i++)
	{
		// 列
		for(j = 0; j < lW; j++)
		{
			// 指向DIBi行j列象素的指针
			lpSrc = lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;
			
			// 将象素值赋值给时域数组
			dpf[j + i * lW] = *(lpSrc);
		}
	}
	
	for(i = 0; i < lH; i++)
		// 对y方向进行沃尔什-哈达玛变换
		WALSH(dpf + lW * i, dpF + lW * i, wp);
		
	// 保存计算结果
	for(i = 0; i < lH; i++)
	{
		for(j = 0; j < lW; j++)
		{
			dpf[j * lH + i] = dpF[j + lW * i];
		}
	}
	
	for(j = 0; j < lW; j++)
		// 对x方向进行沃尔什-哈达玛变换
		WALSH(dpf + j * lH, dpF + j * lH, hp);
	
	// 行
	for(i = 0; i < lH; i++)
	{
		// 列
		for(j = 0; j < lW; j++)
		{
			// 计算频谱
			dTemp = fabs(dpF[j * lH + i] * 1000);
			
			if (dTemp > 255)
			{
				// 超过255直接设置为255
				dTemp = 255;
			}
			
			// 指向DIBi行j列象素的指针
			lpSrc = lpDIBBits + lLineBytes * (lHeight - 1 - i) + j;
			
			// 更新源图像
			* (lpSrc) = (BYTE)(dTemp);
		}
	}
	
	//释放内存
	delete dpf;
	delete dpF;

	// 返回
	return TRUE;
}