图像频率域分析之傅里叶变换_研究傅里叶时域图像时频率是不变的吗-程序员宅基地

技术标签: 计算机视觉  图像处理  Computer Vision  

Overview

欢迎访问 持续更新:https://cgabc.xyz/posts/bcb5a7d3/

傅里叶变换基础

傅里叶级数

法国数学家傅里叶发现,任何周期函数都可以用正弦函数和余弦函数构成的无穷级数来表示(选择正弦函数与余弦函数作为基函数是因为它们是正交的),即 任何周期信号都可以表示成一系列正弦信号的叠加

  • 三角形式

f ( t ) = a 0 2 + ∑ k = 1 + ∞ [ a k c o s ( n ω t ) + b k s i n ( n ω t ) ] , a 0 2 = 1 T ∫ − T 2 T 2 f ( t ) d t f(t) = \frac{a_0}{2} + \sum_{k=1}^{+\infty} \big[ a_k cos (n \omega t) + b_k sin (n \omega t) \big], \quad \frac{a_0}{2} = \frac{1}{T} \int_{-\frac{T}{2}}^{\frac{T}{2}} f(t) dt f(t)=2a0+k=1+[akcos(t)+bksin(t)],2a0=T12T2Tf(t)dt

  • 复指数形式

f ( t ) = 1 T ∑ n = − ∞ + ∞ [ ∫ − T 2 T 2 f ( τ ) e − j ω n τ d τ ] e j ω n t f(t) = \frac{1}{T} \sum_{n=-\infty}^{+\infty} [ \int_{-\frac{T}{2}}^{\frac{T}{2}} f(\tau)e^{-j\omega_n\tau} d\tau ] e^{j\omega_nt} f(t)=T1n=+[2T2Tf(τ)ejωnτdτ]ejωnt

基波角频率 ω = 2 π T \omega = \frac{2\pi}{T} ω=T2π T T T f ( t ) f(t) f(t) 的周期, j j j 为虚数单位

傅里叶积分

复指数形式

f ( t ) = 1 2 π ∫ − ∞ + ∞ [ ∫ − ∞ + ∞ f ( τ ) e − j ω τ d τ ] e j ω t d ω f(t) = \frac{1}{2\pi} \int_{-\infty}^{+\infty} [ \int_{-\infty}^{+\infty} f(\tau)e^{-j\omega\tau} d\tau ] e^{j\omega t} d\omega f(t)=2π1+[+f(τ)eτdτ]etdω

傅里叶变换

一维连续傅里叶变换

正变换

F ( ω ) = ∫ − ∞ + ∞ f ( t ) e − j ω t d t F(\omega) = \int_{-\infty}^{+\infty} f(t) e^{-j\omega t} dt F(ω)=+f(t)etdt

逆变换

f ( t ) = 1 2 π ∫ − ∞ + ∞ F ( ω ) e j ω t d ω f(t) = \frac{1}{2\pi} \int_{-\infty}^{+\infty} F(\omega) e^{j\omega t} d\omega f(t)=2π1+F(ω)etdω

一维离散傅里叶变换

正变换

F ( u ) = ∑ x = 0 M − 1 f ( x ) e − j 2 π u x M F(u) = \sum_{x=0}^{M-1} f(x) e^{-j2\pi \frac{ux}{M}} F(u)=x=0M1f(x)ej2πMux

F ( 0 ) = ∑ x = 0 M − 1 f ( x ) F(0) = \sum_{x=0}^{M-1} f(x) F(0)=x=0M1f(x)

反变换

f ( x ) = 1 M ∑ u = 0 M − 1 F ( u ) e j 2 π u x M f(x) = \frac{1}{M} \sum_{u=0}^{M-1} F(u) e^{j2\pi \frac{ux}{M}} f(x)=M1u=0M1F(u)ej2πMux

对于反变换式前的系数 1 M \frac{1}{M} M1 ,也可放在正变换中,只要保证正变换与反变换之前的系数乘积为 1 M \frac{1}{M} M1 即可。

二维离散傅里叶变换

正变换

二维离散傅里叶变换:

F [ f ( x , y ) ] = F ( u , v ) = ∑ x = 0 M − 1 ∑ y = 0 N − 1 f ( x , y ) e − j 2 π ( u x M + v y N ) F[f(x,y)] = F(u,v) = \sum_{x=0}^{M-1} \sum_{y=0}^{N-1} f(x,y) e^{-j2\pi(\frac{ux}{M}+\frac{vy}{N})} F[f(x,y)]=F(u,v)=x=0M1y=0N1f(x,y)ej2π(Mux+Nvy)

( u , v ) (u,v) (u,v) 等于 ( 0 , 0 ) (0,0) (0,0) 时,直流分量 为:

F ( 0 , 0 ) = ∑ x = 0 M − 1 ∑ y = 0 N − 1 f ( x , y ) F(0,0) = \sum_{x=0}^{M-1} \sum_{y=0}^{N-1} f(x,y) F(0,0)=x=0M1y=0N1f(x,y)

幅度谱 为:

A ( u , v ) = ∣ F ( u , v ) ∣ = R e a l ( u , v ) 2 + I m a g e ( u , u ) 2 A(u,v) = |F(u,v)| = \sqrt{Real(u,v)^2 + Image(u,u)^2} A(u,v)=F(u,v)=Real(u,v)2+Image(u,u)2

功率谱 为:

P ( u , v ) = ∣ F ( u , v ) ∣ 2 = R e a l ( u , v ) 2 + I m a g e ( u , u ) 2 P(u,v) = |F(u,v)|^{2} = Real(u,v)^2 + Image(u,u)^2 P(u,v)=F(u,v)2=Real(u,v)2+Image(u,u)2

相位谱 为:

ϕ ( u , v ) = a r c t a n I m a g e ( u , v ) R e a l ( u , v ) \phi(u,v) = arctan \frac{Image(u,v)}{Real(u,v)} ϕ(u,v)=arctanReal(u,v)Image(u,v)

通过 幅度谱 和 相位谱,我们也能合成 其傅里叶变换(频谱):

F ( u , v ) = A ( u , v ) e j ϕ ( u , v ) = A ( c o s ϕ + j s i n ϕ ) (省略(u,v),应用 欧拉公式) = A c o s ϕ + j A s i n ϕ \begin{aligned} F(u,v) &= A(u,v)e^{j\phi(u,v)} \\ &= A( cos \phi + jsin \phi ) \quad \text{(省略(u,v),应用 欧拉公式)}\\ &= Acos\phi + jAsin\phi \end{aligned} F(u,v)=A(u,v)ejϕ(u,v)=A(cosϕ+jsinϕ)(省略(u,v),应用 欧拉公式)=Acosϕ+jAsinϕ

注意:

  • 上面式子中的 j j j虚数单位
  • R e a l ( u , v ) Real(u,v) Real(u,v) 为 复数的 实部
  • I m a g e ( u , v ) Image(u,v) Image(u,v) 为 复数的 虚部
反变换

f ( x , y ) = F − 1 ( u , v ) = 1 M N ∑ u = 0 M − 1 ∑ v = 0 N − 1 F ( u , v ) e j 2 π ( u x M + v y N ) f(x,y) = F^{-1}(u,v) = \frac{1}{MN} \sum_{u=0}^{M-1} \sum_{v=0}^{N-1} F(u,v) e^{j2\pi(\frac{ux}{M}+\frac{vy}{N})} f(x,y)=F1(u,v)=MN1u=0M1v=0N1F(u,v)ej2π(Mux+Nvy)

卷积

∫ − ∞ + ∞ f 1 ( τ ) f 2 ( t − τ ) d τ = f 1 ( t ) ∗ f 2 ( t ) \int_{-\infty}^{+\infty} f_1(\tau)f_2(t-\tau) d\tau = f_1(t) * f_2(t) +f1(τ)f2(tτ)dτ=f1(t)f2(t)

卷积定理

函数卷积的傅立叶变换是函数傅立叶变换的乘积

  • 时域卷积定理:时域内的卷积对应频域内的乘积

F [ f 1 ( t ) ∗ f 2 ( t ) ] = F 1 ( ω ) ⋅ F 2 ( ω ) F[f_1(t) * f_2(t)] = F_1(\omega) \cdot F_2(\omega) F[f1(t)f2(t)]=F1(ω)F2(ω)

  • 频域卷积定理:频域内的卷积对应时域内的乘积

F [ f 1 ( t ) ⋅ f 2 ( t ) ] = 1 2 π F 1 ( ω ) ∗ F 2 ( ω ) F[f_1(t) \cdot f_2(t)] = \frac{1}{2\pi} F_1(\omega) * F_2(\omega) F[f1(t)f2(t)]=2π1F1(ω)F2(ω)

数字图像DFT

借用知乎大神Heinrich的一张图,来个感性认识:

image_dft_plot.jpg

空间域和频域

  • 空间域:在图像处理中,时域可以理解为 空间域 或者 图像空间,处理对象为图像像元;
  • 频域:以 空间频率 为自变量描述图像的特征,可以将一幅图像像元值在空间上的变化分解为具有不同振幅、空间频率和相位的简振函数的线性叠加,图像中各种空间频率成分的组成和分布称为 图像频谱

空间域与频域可互相转换,对图像施行 二维离散傅立叶变换 或 小波变换 ,可以将图像由空间域转换到频域;通过 对应的反变换 又可转换回空间域图像,即人可以直接识别的图像。

图像频域滤波

二维数字图像的滤波主要分为 空间域滤波 和 频域滤波:

  • 空间域滤波: 用各种模板直接与图像进行 卷积运算,实现对图像的处理,这种方法直接对图像空间操作,操作简单

  • 频域滤波: 在实现某些图像处理的时候,频域的处理比空间域更简单;对于在空间域上的数字图像,根据 卷积定理 可以通过 傅立叶变换空域卷积滤波 变换为 频域滤波,然后再将频域滤波处理后的图像 反变换 回空间域

基本步骤

图像频域滤波步骤为(频谱图中心化):

  • 计算 原始图像 f ( x , y ) f(x,y) f(x,y) 的DFT,得到 频谱 F ( u , v ) F(u,v) F(u,v)
  • 中心化:将频谱 F ( u , v ) F(u,v) F(u,v) 的零频点移动到频谱图的中心位置
  • 计算 滤波器函数 H ( u , v ) H(u,v) H(u,v) F ( u , v ) F(u,v) F(u,v) 的乘积 G ( u , v ) = F ( u , v ) ⋅ H ( u , v ) G(u,v) = F(u,v) \cdot H(u,v) G(u,v)=F(u,v)H(u,v)
  • 反中心化:将频谱 G ( u , v ) G(u,v) G(u,v) 的零频点移回到频谱图的左上角位置
  • 计算上一步计算结果的 傅里叶反变换 g ( x , y ) g(x,y) g(x,y)
  • g ( x , y ) g(x,y) g(x,y) 的 实部 作为最终滤波后的结果图像

上面步骤是对 图像频谱 进行 中心变换;我们也可以先对 原始图像 进行 中心变换,再计算其 频谱图,滤波步骤如下(原始图中心化):

  • 原始图像 f ( x , y ) f(x,y) f(x,y) 中心变换: f ( x , y ) ⋅ ( − 1 ) ( x + y ) f(x,y) \cdot (-1)^{(x+y)} f(x,y)(1)(x+y)
  • 计算上一步计算结果的DFT,得到其 频谱 F ( u , v ) F(u,v) F(u,v)
  • 计算 滤波器函数 H ( u , v ) H(u,v) H(u,v) F ( u , v ) F(u,v) F(u,v) 的乘积 G ( u , v ) = F ( u , v ) ⋅ H ( u , v ) G(u,v) = F(u,v) \cdot H(u,v) G(u,v)=F(u,v)H(u,v)
  • 计算 G ( u , v ) G(u,v) G(u,v) 的 傅里叶反变换 g ( x , y ) g(x,y) g(x,y)
  • g ( x , y ) g(x,y) g(x,y) 的 实部
  • 上一步计算结果 乘以 ( − 1 ) ( x + y ) (-1)^{(x+y)} (1)(x+y) 作为最终滤波后的结果图像

滤波能否取得理想结果的关键取决于上面的 滤波器函数 H ( u , v ) H(u,v) H(u,v)

这时让我想到了《自动控制理论》中的 传递函数 G ( s ) G(s) G(s),定义为:初始条件为零的线性定常系统输出的拉普拉斯变换与输入的拉普拉斯变换之比。

下面以 控制论的思想 给出图像频域滤波的示意框图:

image_fft_flow.jpg

图像频率特性分析

频谱图上的每一个像素点都代表一个频率值,幅值由像素点亮度变码而得。对于一幅图像,图像信号的 频率特性 如下:

  • 直流分量 表示预想的平均灰度
  • 低频分量 代表了大面积背景区域和缓慢变化部分
  • 高频分量 代表了它的边缘、细节、跳跃部分以及颗粒噪声
  • 振幅 描述了图像灰度的亮度
  • 相位 决定了图像是什么样子

数字图像的二维离散傅立叶变换所得的结果的频域成分如下图所示,左上角是直流成分,变换结果四个角周围对应于低频成分,中央部分对应于高频部分

image_fft2_figure.png

为了便于观察,常常采取 换位 方法使直流成分出现在窗口的中央(中心化),变换后中心为低频,向外是高频。

在频域,可以很方便的实现 图像的锐化和模糊

  • 截取频率的低频分量,对其作傅立叶反变换,得到的就是模糊后的图像,即 低通滤波
  • 截取频率的高频分量,对其作傅立叶反变换,得到的就是锐化后的图像,即 高通滤波

图像滤波实践

下面,我们以 lena.bmp(点此下载) 图像进行滤波实践。

Python分析

(1)加载图像,并转换为 灰度图
image_gray.png

(2)对其 快速傅里叶变换,并经过 中心变换,得到 频率谱相位谱
image_frequency.png image_frequency_phase.png

(3)分别截取 频谱图低频部分(中间部分)高频分量(四周部分)
image_frequency_lf.png image_frequency_hf.png

(4)对以上处理过的频谱图分别进行 反中心化傅里叶反变换取实部,得到 低通滤波高通滤波 后的图像
image_back_lf.png image_back_hf.png

C++分析

使用 CImgFFTW库 对 lena图像进行傅里叶变换(源代码见文末),结果如下
lena_fftw.jpg

源代码

以上所有代码均存储在我的Github仓库:

参考资料

版权声明:本文为博主原创文章,遵循 CC 4.0 BY-SA 版权协议,转载请附上原文出处链接和本声明。
本文链接:https://blog.csdn.net/u011178262/article/details/84674412

智能推荐

【备忘录】 iPad air3 pro 10.5 原装键盘 失灵 维修_ipad10.5键盘套没反应-程序员宅基地

文章浏览阅读2.5k次。并不是苹果说的触点问题,而是软布导线被腐蚀断线问题。用铜箔胶带粘一下就好了。不过二手直接买一个咸鱼上也就三百多,也不贵。这个明显是苹果设计缺陷。直接丢了挺可惜的。拆的时候小心点。我的QQ: 13008312096,有空的话可以帮忙有偿代处理。参考链接如下:nullhttps://www.youtube.com/watch?v=buNYHzMZJdk修复好的样子,毫无违和感:..._ipad10.5键盘套没反应

gif图用photoshop快速抽帧,压缩,去底,修改时间_gif抽帧-程序员宅基地

文章浏览阅读1.3w次,点赞3次,收藏13次。网上下载的gif图很多都不透底,而且尺寸时间什么的都不合适,用ps简单修改一下就好了1.抽帧很多gif图帧数很多,抽帧可以大幅压缩gif大小1.导入gif到ps2.勾选动作和时间轴面板这时我们发现帧和图层的隐藏显示是相对应的这样的话我们操作图层就会打乱帧,我们要让所有的图层都显示,但又不干扰帧的显示3.这时我们就要把它转换为视频时间轴4.然后再转换帧>>转换为..._gif抽帧

c++ 优先队列_面试必知必会|堆和优先队列-程序员宅基地

文章浏览阅读235次。通过本文将了解到以下内容:优先队列的概念优先队列的实现优先队列的应用1.优先队列的概念优先队列是计算机科学中的一类抽象数据类型。优先队列中的每个元素都有各自的优先级,优先级最高的元素最先得到服务;优先级相同的元素按照其在优先队列中的顺序得到服务。优先队列至少需要支持下述操作:a.插入带优先级的元素b.取出具有最高优先级的元素c.查看最高优先级的元素。综合考虑插入和删除的性能 优先队列一般采用堆来实..._优先队列 c++面试

windows系统VMware安装Linux虚拟机、Linux配置JDK环境、Linux安装tomcat、Linux安装mysql、Linux通过SSH连接Navicat 数据库可视化工具_vm中linux下载安装mysql jdk-程序员宅基地

文章浏览阅读281次。windows系统VMware安装Linux虚拟机、配置JDK环境、安装tomcat、安装mysql、安装redis第一部分VMware安装Linux—————————————————————————————————xshell、VMware工具链接:解压之后安装VMware工具xshell不需要安装解压可直接使用..._vm中linux下载安装mysql jdk

Python 加密解密_咕咕加密v4-程序员宅基地

文章浏览阅读1k次。# -*- coding: utf-8 -*-import hashlib'''加密解密'''#MD5是最常见的摘要算法,速度很快,生成结果是固定的128 bit字节,# 通常用一个32位的16进制字符串表示。md5 = hashlib.md5()update = md5.update('hhhhhaaa')print(md5.hexdigest())#SHA1的结果是160..._咕咕加密v4

使用cesium primitive api绘制三维插值图,热力图_cesium 插值点-程序员宅基地

文章浏览阅读1.4k次,点赞3次,收藏13次。【代码】使用cesium primitive api绘制三维插值图。_cesium 插值点

随便推点

Ubuntu18.04安装教程——超详细的图文教程_ubuntu系统18.04-程序员宅基地

文章浏览阅读10w+次,点赞112次,收藏974次。Ubuntu18.04镜像_ubuntu系统18.04

二、Json对象、Json数组和Json字符串_json字符串数组-程序员宅基地

文章浏览阅读6.7k次,点赞3次,收藏17次。一、Json字符串和Json对象定义:1、Json字符串:所谓字符串:单引号或者双引号引起来,是一个String类型的字符串:如下:var person='{"name":"shily","sex":"女","age":"23"}';//json字符串console.log(person)console.log(person.name)console.log(typeof person) 2、Json对象:最显著的特征:对象的值可以用 “对象.属性” 进行访问,_json字符串数组

Linux系统100条命令:关于Ubuntu和 CentOS 7 相同功能的不同的终端操作命令_ubuntu 命令跟centos-程序员宅基地

文章浏览阅读718次。CentOS 7:ip link set interface_name up 或 ip link set interface_name down。Ubuntu:ifconfig interface_name up 或 ifconfig interface_name down。CentOS 7:编辑 /etc/sysconfig/network-scripts/ifcfg-eth0 文件。Ubuntu:编辑 /etc/network/interfaces 文件。_ubuntu 命令跟centos

windows10下VS2019编译jpegsrc.v9e.tar.gz为lib静态库(已验证)_jpeg library error vs2019-程序员宅基地

文章浏览阅读652次。jpegsr9e windows vs2019生成方法,以及库下载_jpeg library error vs2019

重磅?华为 Mate60/Pro 系列网速实测结果公布,最高 1205.57 Mbps_华为mate60pro+核实网络-程序员宅基地

文章浏览阅读647次。总的来说,华为Mate 60/Pro系列手机的高速网速表现引起了广泛的关注,这也是消费者对该系列手机购买热情高涨的一个重要因素。可以看出,华为Mate 60/Pro系列手机的网速表现非常出色,这也是消费者购买该系列手机的一个重要原因。此前,华为Mate 60 Pro的供应量已经增至1500万至1700万台,而最新消息称,华为Mate 60 Pro和Mate 60 Pro+的出货量甚至已上调至2000万台。目前,在中国市场上,手机竞争愈发激烈,不仅华为Mate 60系列,其他品牌的手机也都受到了高温的迎接。_华为mate60pro+核实网络

access查找出生日期年份_access怎样利用出生日期计算年龄呀!-程序员宅基地

文章浏览阅读7.1k次。公告: 为响应国家净网行动,部分内容已经删除,感谢读者理解。话题:access怎样利用出生日期计算年龄呀!回答:lt;%set rs = server.createobject("adodb.recordset") curid=request("id") sql = "UPDATE pany SET a_num=a_num+1,day_count=day_count+1 WHERE day_lda..._access出生年份表达式