GIS数据处理-IDL自带函数的巴特沃斯（butterworth）高通/低通滤波

在之前的文章中，我曾经介绍过IDL实现巴特沃斯（butterworth）高通滤波和IDL实现巴特沃斯（butterworth）低通滤波，不过自己写的套上几重FOR循环，总感觉效率太低。今天介绍的是使用IDL自带的巴特沃斯（butterworth）高通/低通滤波。

下面是冈萨雷斯版本的巴特沃斯（butterworth）滤波器。

巴特沃斯（butterworth）低通滤波器

巴特沃斯（butterworth）高通滤波器

IDL自带巴特沃斯（butterworth）滤波器

\[\frac{1}{\sqrt{1+(\frac{\omega}{\omega_{c}})^{2N}}}\]

从公式以及IDL官方帮助文档里面可以看出这里实现的是一个低通的巴特沃斯。

示例代码

;AUTHOR:尹全超
PRO mybtworth
filter = BUTTERWORTH(500,500,CUTOFF=5,/ORIGIN)
img = image(filter)
END

效果图

这是一个500*500，截断频率为5的低通滤波器，具体函数原型为：

Result = BUTTERWORTH( X [, Y [, Z]] [, CUTOFF=value] [, ORDER=value] [, /ORIGIN] [, XDIM=value] [, YDIM=value] [, ZDIM=value] )

我看了下相关参数没有关于高通巴特沃斯滤波的实现，于是想了下面的办法。

更新：IDL中，可以通过BUTTERWORTH的ORDER参数控制高通和低通滤波，具体方法见：更新IDL实现巴特沃斯（butterworth）高通滤波，但是下面的方法仍有参考意义。

方法一：整体减去低通，这个用于在具体滤波的时候，乘以原始图像，然后再用原始图像减去低通的图像，剩下的就是高通。

方法二：可以对上面这个公式做一个推导，推导成高通滤波函数，最后结果如下：

高通 = sqrt(1 - 低通 * 低通)

示例代码

;Author: 尹全超
PRO higpassbutterworth
filter = BUTTERWORTH(500,500,CUTOFF=100,/ORIGIN)
filter = sqrt(1-filter*filter)
img = image(filter)
END

效果如图