5.14 图像频域处理

在图像处理和分析中，经常会将图像从图像空间转换到其他空间中，并利用这些空间的特点进行对转换后图像进行分析处理，然后再将处理后的图像转换到图像空间中，这称之为图像变换。在一些图像处理和分析中通过空间变换往往会取得更有效的结果。图像频域处理是指将图像从图像空间转换到频域空间进行处理的过程。最常用的频域转换是傅里叶变换。对于傅里叶变换的理论内容不是这里讨论的重点，读者可以从网上找到大量的资料。傅里叶变换的计算量较大，人们为了提高速度，提出了快速傅里叶变换，并得到了广泛的应用。下面我们首先来学习一下VTK中的快速傅里叶变换。

5.14.1 快速傅里叶变换

快速傅里叶变换(Fast Fourier Transform)，简称做FFT。它是根据离散傅氏变换的奇、偶、虚、实等特性，对离散傅立叶变换的算法进行改进获得的。傅里叶变换是可逆的，其逆变换为RFFT。FFT在数字图像处理中有着广泛的应用，例如数字图像频域滤波，去噪，增强等等。目前VTK中两变换都已经实现，对应的类分别为vtkImageFFT和vtkImageRFFT。下面代码演示了怎样对图像进行傅里叶变换，以及将傅里叶变换结果进行逆变换。

vtkImageFFT和vtkImageRFFT的输入为实数或者复数数据，输出都为复数数据。因此，vtkImageFFT与vtkImageRFFT的输出结果不能直接显示，因为VTK会将其当做彩色图像显示，需要通过vtkImageExtractComponents提取某一组分图像显示。

VTK频率域的图像处理步骤如下：

图5.35VTK频率域的图像处理步骤

下面代码说明了怎样对一副图像进行快速傅里叶变换和逆变换。

1: vtkSmartPointer<vtkJPEGReader>reader =

2: vtkSmartPointer<vtkJPEGReader>::New();

3: reader->SetFileName ( " lena2.jpg" );

4: reader->Update();

6: vtkSmartPointer<vtkImageFFT> fftFilter =

7: vtkSmartPointer<vtkImageFFT>::New();

8: fftFilter->SetInputConnection(reader->GetOutputPort());

9: fftFilter->Update();

10:

11: vtkSmartPointer<vtkImageExtractComponents> fftExtractReal =

12: vtkSmartPointer<vtkImageExtractComponents>::New();

13: fftExtractReal->SetInputConnection(fftFilter->GetOutputPort());

14: fftExtractReal->SetComponents(0);

15:

16: vtkSmartPointer<vtkImageCast> fftCastFilter =

17: vtkSmartPointer<vtkImageCast>::New();

18: fftCastFilter->SetInputConnection(fftExtractReal->GetOutputPort());

19: fftCastFilter->SetOutputScalarTypeToUnsignedChar();

20: fftCastFilter->Update();

21:

22: vtkSmartPointer<vtkImageRFFT> rfftFilter =

23: vtkSmartPointer<vtkImageRFFT>::New();

24: rfftFilter->SetInputConnection(fftFilter->GetOutputPort());

25: rfftFilter->Update();

26:

27: vtkSmartPointer<vtkImageExtractComponents> ifftExtractReal =

28: vtkSmartPointer<vtkImageExtractComponents>::New();

29: ifftExtractReal->SetInputConnection(rfftFilter->GetOutputPort());

30: ifftExtractReal->SetComponents(0);

31:

32: vtkSmartPointer<vtkImageCast> rfftCastFilter =

33: vtkSmartPointer<vtkImageCast>::New();

34: rfftCastFilter->SetInputConnection(ifftExtractReal->GetOutputPort());

35: rfftCastFilter->SetOutputScalarTypeToUnsignedChar();

36: rfftCastFilter->Update();

首先建立一个JPEG图像reader来读取一副灰度图像。然后定义一个vtkImageFFT指针，直接接收reader的输出即原图像数据作为输入进行二维快速傅里叶变换。其输出为一个像素类型为复数的vtkImageData数据，即每个像素值为两个组分（Component）：复数实部和虚部。因此如果直接显示这个vtkImageData，会发现是一个彩色图像。如果需要显示频域图像，需要通过vtkImageExtractComponents类来提取某一个组分图像来显示。上例中通过定义vtkImageExtractComponents类指针，利用函数SetComponents(0)指定提取实部图像显示；由于vtkImageActor类仅支持unsigned char数据类型的图像，利用vtkImageCast类的SetOutputScalarTypeToUnsignedChar()指定输出类型为unsigned char，将FFT结果图像转换为需要的类型。

对于逆变换的过程也是类似。首先定义vtkImageRFFT指针，并接收输入为vtkImageFFT指针的输出，调用Update执行函数完成快速傅里叶逆变换。vtkImageRFFT的输出同样为一副复数图像，通常不能直接显示或者进行其他操作。对于傅里叶逆变换的图像中虚数部分值为0，实数部分图像即为重建的原始图像。因此再次利用vtkImageExtractComponents提取实数部分图像，并通过vtkImageCast的SetOutputScalarTypeToUnsignedChar()将图像转换为unsigned char类型进行显示。下图为计算结果。

图5.36 快速傅立叶变换

图像通过FFT变换到频域后，可以通过改变不同的频率分量来对图像进行处理。图像频域处理的步骤是，先将空域图像通过FFT变换到频域图像，设计一个滤波器对不同频率的频域分量进行处理，将处理后的频域图像通过IFFT变换回空域中得到处理后的图像。滤波器的设计需要根据具体需求来变化，主要分为两种，低通滤波器和高通滤波器。下面我们来分别介绍。

5.14.2 低通滤波

5.14.2.1 理想低通滤波器

低通滤波是将频域图像中的高频部分滤除而通过低频部分。图像的边缘和噪声对应于频域图像中的高频部分，而低通滤波的作用即是减弱这部分的能量，从而达到图像平滑去噪的目的。最简单的低通滤波器是理想低通滤波器，基本思想是给定一个频率阈值，将高于该阈值的所有部分设置为0，而低于该频率的部分保持不变。理想是指该滤波器不能用电子元器件来实现，但是可以通过计算机来模拟。在VTK中定义了理想低通滤波器，下面我们来看下怎么使用该滤波器来对图像进行低通滤波。

1: vtkSmartPointer<vtkJPEGReader>reader =

2: vtkSmartPointer<vtkJPEGReader>::New();

3: reader->SetFileName("lena2.jpg");

4: reader->Update();

6: vtkSmartPointer<vtkImageFFT> fftFilter =

7: vtkSmartPointer<vtkImageFFT>::New();

8: fftFilter->SetInputConnection(reader->GetOutputPort());

9: fftFilter->Update();

10:

11: vtkSmartPointer<vtkImageIdealLowPass> lowPassFilter =

12: vtkSmartPointer<vtkImageIdealLowPass>::New();

13: lowPassFilter->SetInputConnection(fftFilter->GetOutputPort());

14: lowPassFilter->SetXCutOff(0.05);

15: lowPassFilter->SetYCutOff(0.05);

16: lowPassFilter->Update();

17:

18: vtkSmartPointer<vtkImageRFFT> rfftFilter =

19: vtkSmartPointer<vtkImageRFFT>::New();

20: rfftFilter->SetInputConnection(lowPassFilter->GetOutputPort());

21: rfftFilter->Update();

22:

23: vtkSmartPointer<vtkImageExtractComponents> ifftExtractReal =

24: vtkSmartPointer<vtkImageExtractComponents>::New();

25: ifftExtractReal->SetInputConnection(rfftFilter->GetOutputPort());

26: ifftExtractReal->SetComponents(0);

27:

28: vtkSmartPointer<vtkImageCast> castFilter =

29: vtkSmartPointer<vtkImageCast>::New();

30: castFilter->SetInputConnection(ifftExtractReal->GetOutputPort());

31: castFilter->SetOutputScalarTypeToUnsignedChar();

32: castFilter->Update();

首先读入一副jpg图像，通过vtkImageFFT将图像转换到频域空间。vtkImageIdealLowPass对频域图像做理想低通滤波，需要用户设置每个方向的截断频率，相应的设置函数SetXCutOff()和SetYCutOff()。执行完毕后，需要通过vtkImageRFFT将处理后的频域图像转换至空域图像。需要注意的是，转换后的图像每个像素都是一个复数，需要vtkImageExtractComponents将该图像的第一个分量提出出来显示，否则图像不能正确显示。由于傅里叶变换输入输出的数据类型都是double，为了方便显示，还需要将其转换为Unsigned char类型，这里vtkImageCast负责类型转换。下面是对图像做低通滤波的效果。从结果看，在过滤掉图像的高频部分后，图像变得模糊，丢失了许多细节，另外还可以看到图像会存在一定的振铃效应，这也是理想低通滤波的特点。

图5.37 理想低通滤波器

5.14.2.2 巴特沃斯低通滤波器

在实际中经常使用的低通滤波器是巴特沃斯滤波器。巴特沃斯滤波器对应的转移函数（可以看做是一个系数矩阵）是：

其中D(u,v)表示频域点(u,v)到频域图像原点的距离，称为截止频率，当D(u,v) = 时，H(u,v)=0.5，即对应的频域能量将为原来的一半。因为巴特沃斯低通滤波器在高低频间的过渡平滑，因此不会出现明显的振铃效应。VTK中实现巴特沃斯低通滤波器的类是vtkImageButterworthLowPass，其使用过程如下。

1: vtkSmartPointer<vtkImageButterworthLowPass> lowPassFilter =

2: vtkSmartPointer<vtkImageButterworthLowPass>::New();

3: lowPassFilter->SetInputConnection(fftFilter->GetOutputPort());

4: lowPassFilter->SetXCutOff(0.05);

5: lowPassFilter->SetYCutOff(0.05);

6: lowPassFilter->Update();

vtkImageButterworthLowPass与理想低通滤波器的使用一样，因此这里只复制了部分代码。为了便于比较，这里设置X和Y方向的截止频率时，与理想低通滤波器设置一致，下面是相应的执行结果。从结果来看，巴特沃斯低通滤波器产生的图像更为平滑，不会出现振铃现象。

图5.38 巴特沃斯滤波器

5.14.3 高通滤波

5.14.3.1 理想高通滤波器

高通滤波与低通滤波正好相反，是频域图像的高频部分通过而抑制低频部分。在图像中图像的边缘对应高频分量，因此高通滤波的效果是图像锐化。同样最简单的高通滤波器是理想高通滤波器。通过设置一个频率阈值，将高于该阈值的频率部分通过，而低于阈值的低频部分设置为0。下面来看一下VTK中理想高通滤波的实例。

1: vtkSmartPointer<vtkJPEGReader>reader =

2: vtkSmartPointer<vtkJPEGReader>::New();

3: reader->SetFileName("lena2.jpg");

4: reader->Update();

6: vtkSmartPointer<vtkImageFFT> fftFilter =

7: vtkSmartPointer<vtkImageFFT>::New();

8: fftFilter->SetInputConnection(reader->GetOutputPort());

9: fftFilter->Update();

10:

11: vtkSmartPointer<vtkImageIdealHighPass> highPassFilter =

12: vtkSmartPointer<vtkImageIdealHighPass>::New();

13: highPassFilter->SetInputConnection(fftFilter->GetOutputPort());

14: highPassFilter->SetXCutOff(0.1);

15: highPassFilter->SetYCutOff(0.1);

16: highPassFilter->Update();

17:

18: vtkSmartPointer<vtkImageRFFT> rfftFilter =

19: vtkSmartPointer<vtkImageRFFT>::New();

20: rfftFilter->SetInputConnection(highPassFilter->GetOutputPort());

21: rfftFilter->Update();

22:

23: vtkSmartPointer<vtkImageExtractComponents> ifftExtractReal =

24: vtkSmartPointer<vtkImageExtractComponents>::New();

25: ifftExtractReal->SetInputConnection(rfftFilter->GetOutputPort());

26: ifftExtractReal->SetComponents(0);

27:

28: vtkSmartPointer<vtkImageCast> castFilter =

29: vtkSmartPointer<vtkImageCast>::New();

30: castFilter->SetInputConnection(ifftExtractReal->GetOutputPort());

31: castFilter->SetOutputScalarTypeToUnsignedChar();

32: castFilter->Update();

同低通滤波一样，首先将读入图像通过vtkImageFFT转换到频域空间，定义vtkImageIdealHighPass对象，并通过SetXCutOff ()和SetYCutOff() 设置X和Y方向的截止频率。然后通过vtkImageRFFT将处理后的图像转换到空域中，得到高通滤波图像。为了显示的需要，还需要提取图像分量和数据类型的转换，这里与低通滤波一致，不再赘述。下面是理想高通滤波的执行结果。从结果看出高通滤波后图像得到锐化处理，图像中仅剩下边缘。

图5.39 理想高通滤波器

5.14.3.2 巴特沃斯高通滤波器

理想高通滤波器不能通过电子元器件来实现，而且存在振铃现象。在实际中最常使用的高通滤波器是巴特沃斯高通滤波器。该滤波器的转移函数是：

D(u,v)表示频域中点到频域平面的距离，是截止频率。当D(u,v)大于时，对应的H(u,v)逐渐接近1，从而使得高频部分得以通过；而当D(u,v)小于时，H(u,v)逐渐接近0，实现低频部分过滤。巴特沃斯高通滤波器在VTK中对应vtkImageButterworthHighPass类。下面代码说明了vtkImageButterworthHighPass对图像进行高通滤波。

1: vtkSmartPointer<vtkImageButterworthLowPass> lowPassFilter =

2: vtkSmartPointer<vtkImageButterworthLowPass>::New();

3: lowPassFilter->SetInputConnection(fftFilter->GetOutputPort());

4: lowPassFilter->SetXCutOff(0.05);

5: lowPassFilter->SetYCutOff(0.05);

6: lowPassFilter->Update();

vtkImageButterworthHighPass与理想高通滤波使用方法一致，这里没有给出傅里叶变换和逆变换代码。这里也需要设置X和Y轴的截止频率，为了便于比较，其截止频域与理想高通滤波设置一致。下图是执行结果。