MATLAB在图像处理中的应用

图像处理概念、算法和MATLAB

FFT故事

我那天的最后一次会议刚刚结束不久,我开始认真考虑为本周写一篇博文。写什么?就在那时,我碰巧看到克利夫在就在昨天关于FFT。啊哈!我也会写这篇文章(为了纪念FFT周的国家博客)。

目录

是否重现?

我有三个小故事要讲,第一个故事的灵感来自克利夫的文章。克利夫描述了一个非常紧凑,优雅的分而治之思想的实现,这是快速傅立叶变换算法的核心。下面是Cleve教科书函数的关键代码片段ffttx公司:

%递归分治k=(0:n/2-1)';w=Ω。^k;u=ffttx(x(1:2:n-1));v=w.*ffttx(x(2:2:n));y=[u+v;u-v];

功能ffttx公司递归调用自身,两次。每次递归调用都是对原始输入值的一半进行调用。那么ffttx公司在最后做一点额外的计算来产生最终的结果。

像克里夫一样,我真的很喜欢这段代码。但这让我想起了我第一次学习数字信号处理时看到的所有FFT代码(这是1989年版的教科书样本离散时间信号处理作者:Oppenheim和Schafer):

这个代码在当时是典型的。重点是就地计算并尽可能少地使用乘法。有一件事你在这里看不到:显式递归!函数抖动不直接或间接地称呼自己。这在当时是很典型的。没有人愿意使用内存来复制两个半长数组,也没有人愿意承受许多递归函数调用的速度损失。你在这段代码中看到的另一个典型的东西是一个有趣的叫做位反转排序的东西,但这是另一天的话题(也许)。

无论如何,我的观点是,为当今计算机体系结构构建的现代FFT实现通常使用显式递归。随着分而治之的问题变得越来越小,将输入数组分割成两个或更多更小的副本所带来的额外成本被更易于缓存的执行所抵消。

蝴蝶与法律

我的第二个小故事的灵感来自克利夫发布的一行代码:

y=[u+v;u-v];

在这种情况下Ux(1)vx(2),那么上面的计算就是.

它的信号流程图如下所示(同样来自离散时间信号处理):

几十年来,数字信号处理专家一直把这个小图形称为“蝴蝶图”,因为它的出现。

嗯,为了在研究生院挣点外快,我写了几章的解决方案手册金宝搏官方网站离散时间信号处理. 其中一章涉及FFT。这种努力让我对画蝴蝶图和数数乘法完全厌倦了。一天,我把下面的图表贴在我的小隔间外面:

是的,数字信号处理的人甚至把上面的混乱称为“蝴蝶图”。我们是一群奇怪的人。

蝴蝶图间接地把我引向了我的第三个小故事。

用nan填充向量的最慢方法(有史以来?)

大约十年前,一个技术支持案例被转发给我。一位客户抱怨MA金宝appTLAB在计算一个特定向量的FFT时速度非常慢。其他长度相同的向量会走得很快,但不是这个。

所以我加载了包含向量的客户MAT文件,并开始对其运行FFTs。令我惊讶的是,它真的很慢!我挠头了。然后我看了看输出值,我更惊讶地发现它们都是nan!奇怪!有什么意义?这和速度问题有关吗?

我查看了客户提供的长输入向量,发现只有一个输入值是NaN。

灯泡!

再看看上面复杂的蝴蝶图。注意,每个输出值都依赖于每个输入值。对于任何长度的FFT计算,通常都是这样。对于NaN通过算术运算传播的方式,这意味着输入向量中只有一个NaN意味着快速傅里叶变换永远都会充满南。

事实证明,这与速度问题有关。NaN计算通常在CPU中优化程度较低的部分执行,甚至在软件中执行。大量的计算使这个客户的问题陷入困境。

FFT是一个永无止境的乐趣和神秘的来源!

我写了不少关于傅立叶变换理论和算法的博客文章。你可以在教室里看到他们傅里叶变换范畴.

如果你有更有趣的方法来创建一个矢量的南的想法,一定要把他们添加到这篇文章的评论。




与MATLAB®R2014a一起发布

|
  • 打印
  • 发送电子邮件

评论

如需留言,请点击在这里登录到您的MathWorks帐户或创建一个新帐户。