nsct matlab 讲解,[转载]NSCT——Nonsubsampled contourlet 变换程序开发教程1

08年，被老板逼得走投无路，xx所得项目看来是实在躲不过去，只好硬着头皮上。开发一款图像处理软件，里面内容很多，不过核心就是NSCT变换。NSCT变换06年的新东西，现在跟风研究的热点。可惜的是作为一种新的变换滤波器长期以来只有matlab程序，没有VC程序。使很多算法开发不能应用到工程中，08暑假，大概花了我一个月吧，恩。。期间电脑还中了次病毒(可恶的机器狗第四代！！)，终于把这个软件搞得差不多了，现在回顾一下NSCT的程序包，下次有时间在写点别的。

Nonsubsampled Contourlet变换算法介绍：

对信号的稀疏表示是许多信号处理及应用的基础，2004年Minh N Do、Martin

Vetterli提出了一种能够较好表示二维信号的数学工具--Contourlet变换。Contourlet是用金字塔方向滤波器组(PDFB)来将图像分解成不同尺度下的方向子带的。根据PDFB的结构，PDFB是一个拉普拉斯金字塔滤波器Laplacian

Pyramid

(LP)和一个方向滤波器组的叠加。实验证明，Contourlet变换在图像降噪，纹理，形状的特征提取方面的性能比2-D离散小波变换有了明显的提高。

为了获得平移不变性，本章所用的Nonsubsampled

Contourlet变换(NSCT)是基于Nonsubsampled金字塔(NSP)和Nonsubsampled方向滤波器(NSDFB)的一种变换。首先由NSP对输入图像进行塔形分解，分解为高通和低通两个部分，然后由NSDFB将高频子带分解为多个方向子带，低频部分继续进行如上分解。NSCT是一种新型平移不变，多尺度，多方向性的快速变换。

1. Nonsubsampled Pyramid(NSP)：

Nonsubsampled Pyramid(NSP)和Contourlet的Laplacian

Pyramid(LP)多尺度分析特性不同。图像通过Nonsubsampled

Pyramid(NSP)进行多尺度分解，NSP去除了上采样和下采样，减少了采样在滤波器中的失真，获得了平移不变性。NSP为具有平移不变性滤波结构的NSCT多尺度分析，可以得到与LP分解一样的多尺度分析特性。图2.4(a)处分为3个尺度。

2. Nonsubsampled方向滤波器(NSDFB)

Nonsubsampled方向滤波器(NSDFB)是一个双通道的滤波器，将分布在同方向的奇异点合成NSCT的系数。方向滤波器(DFB)是Bamberger

and Smith提出的。其通过一个l层的树状结构的分解，有效的将信号分成了 个子带，其频带分割成为锲形。Nonsubsampled

DFB(NSDFB)为非采样，减少了采样在滤波器中的失真，获得了平移不变性。并且每个尺度下的方向子图的的大小都和原图同样大小，Contourlet变换为所有子带之和等于原图。NSCT有更多的细节得以保留，变换系数是冗余的。图2.4(b)为三个尺度下对图像频域的分割图，其中每个尺度的方向子带数目以2倍递增，以在1，2，3尺度下的方向子带数目分别为2，4，8个。

NSCT内部变量：

下面列出程序中主要的变量和结构体定义说明：

图像信息结构体:

typedef struct

{ int high;

//-----------------------------------------图像高度

int

width; //---------------------------------------图像宽度

double max;

//--------------------------------------最大像素值

double min;

//---------------------------------------最小像素值

double *data;

//--------------------------------------像素数据区域

}ImgData;

NSCT变换空间结构体：(用于存放NSCT变换滤波器，尺度、自带空间等基本信息)

typedef struct

{ unsigned int

*nlevels;//子带方向数

必须为偶数sublevels[0]为最高级方向数

sublevels[4]为最低级方向数

为负时为无效

unsigned

int

clevels; //------------------尺度层数

ImgData

***high;//------------------------高通子带区域

ImgData

low;//-----------------------------低通子带区域

ImgData

H1;//------------------------------高通滤波器组

ImgData H2;

//------------------------------高通滤波器组

ImgData G1;

//------------------------------低通滤波器组

ImgData G2;

//------------------------------低通滤波器组

ImgData

***filters_dec; //-----------------分解方向滤波器组

ImgData

***filters_rec; //-----------------重构方向滤波器组

BOOL sign;

}nsct_t;

尺度与分解层数滤波器：

const int

dlevels[3]={3,3,3};//--------------尺度为三，个尺度方向数为

const int

Qunx[4]={0,0,0,0};//--------------------------Q梅花采样矩阵

BMGImageStruct

SourceBmp;//---------------------------位图操作结构体

int pfilter_type,dfilter_type;

//滤波器类型选择,

两者都默认为,'maxflat'&'dmaxflat7'滤波器

int pralevels;

//------------图像金字塔层数(不能大于4，尺度在4以内足以满足一般图像处理需求)

nsct_t nsct;

//nsct-------数据计算区域

int shift[2];

//---------------------延迟分量

2.3.6

NSCT内部函数：

下面列出程序中主要的函数定义说明：

void ImgDataM_Var_3m(ImgData *x,float

*mean,float *var,float *mom);

//----------------------------------计算nsct子带系数的均值和三阶中心矩

void

ImgDataIFFT_2D(complex *

pCFData, complex * pCTData,

int nWidth, int nHeight);

//----------------------------------2D-FFT反变换

void

ImgDataFFT_2D(complex *

pCTData, int nWidth, int nHeight,

complex *

pCFData);

//----------------------------------2D-FFT反变换

void ImgDataFConv2(ImgData *out,ImgData

*img,ImgData *mask);

//----------------------------------频域卷积函数

void Nsct_Extend2(ImgData *in,int ru,int

rd,int cl,int cr,ImgData *out,const int flag=1);

//----------------------------------nsct镜像扩展函数

void Nsct_Zconv2(ImgData *in,ImgData

*out,ImgData *filter,const int *mup);

//----------------------------------时域滤波器卷积函数(快速算法)

void Nsct_Efilter2(ImgData *in,ImgData

*out,ImgData *filter);

//----------------------------------边缘扩展的2D滤波

void Nsct_Nssfbdec(ImgData * input,ImgData *

out1,ImgData * out2,ImgData * filter1,ImgData * filter2,const int *

Qx=Qunx);

//----------------------------------2通道非采样子带分解

void Nsct_Mirror2(ImgData

*tomirror);

//----------------------------------镜像变换

void ImgDataConv2(ImgData *out,ImgData

*img,ImgData *mask,ImgData *imgunext);

//----------------------------------相关卷积

void Nsct_Upsample2df(ImgData *h0,int

lev,ImgData *ho);

//----------------------------------频域2D上采样

void Nsct_Symext(ImgData *imglow,ImgData

*upedHi,ImgData *pp );

//----------------------------------2D平衡矩阵生成

void Nsct_Atrousc(ImgData * out,ImgData *

signal,ImgData *filter,int (*matrix)[2]);

//----------------------------------包含抽取的滤波器卷积

void Nsct_Nsdfbdec(ImgData *imghigh,int

currentlevel);

//----------------------------------非抽取方向子带分解

void Nsct_iTtransform();

//----------------------------------反变换主函数

void Nsct_New(HBITMAP hBitmap,const int

*Nlevels=dlevels,const unsigned int Level=3,const int

dfilter=1,const int pfilter=1);

//-----------------------------用来初始化原始目标图像，申请nsct计算空间和存储空间

void Nsct_Transform();

//-----------------------------变换主函数

void Nsct_NsfbDec(ImgData *imglow,ImgData

*imghigh,int i);

//-----------------------------nsct图像金字塔分解

void NsctDelete();

//-----------------------------nsct析构

void ImgDataFreeMemery(ImgData *image);

//-----------------------------内存空间释放

BOOL ImgDataAllocMemery(ImgData *

image);

//-----------------------------内存空间申请

一个示例程序，和显示的样图：