ARTICLE · INTELLIGENCE

战地情报 · 详情页

来自尧图项目组的一线实战观察与深度解析

分数阶傅里叶变换的Matlab代码包:四种实现与图像加密

分数阶傅里叶变换的Matlab代码包:四种实现与图像加密

简介:分数阶傅里叶变换(FRFT)的MATLAB实现代码包,面向信号处理与图像分析方向的开发者及研究人员,解决传统FFT仅能进行0度和90度变换的限制,可在连续角度域对信号进行频谱分析。压缩包内共9个文件,包含8个.m源文件与1个说明txt,梳理了四种实现方式:直接矩阵运算、基于FFT的快速算法、迭代近似法与窗函数插值法,便于对照各自适用场景。已有1915人学习该资源。通过修改阶数α与输入信号,可直观观察不同变换效果;其中chirp2rot、scale等模块还涉及旋转与尺度变换,为图像加密、非平稳信号分析等应用提供了可直接运行的参考代码。整体代码结构清晰,约9KB,适合作为学习FRFT原理或搭建实验框架的起步资料。

1. 分数阶傅里叶变换的Matlab代码包:先看它能干什么

做图像加密的同学应该都有过这种经历:用普通傅里叶变换做频域加密,换一个密钥相位掩码,反变换出来不是边缘噪声,就是恢复图像带着明显的振铃。后来换成分数阶傅里叶变换(FRFT),同一套双随机相位编码结构,只多了一个旋转阶数参数,恢复质量和对密钥的敏感度立刻上了一个台阶。这份资源就是把分数阶傅里叶变换从连续积分变成能在Matlab里直接跑的离散代码,里面给了4种不同实现方式和一套图像加密示例,正好覆盖“变换核怎么选”和“加密流程怎么搭”两个实际问题。适合需要把FRFT用到图像水印、加密或信号分析的从业者,新手可以直接照着跑通,熟手可以把它当基准去对比自己的离散化方案。

2. 四种FRFT实现方式的原理与选型:离散化决定了结果的可用性

2.1 分数阶傅里叶变换的定义:从旋转一个角度说起

分数阶傅里叶变换的本质,是把时频平面上的坐标轴旋转一个角度α。常规傅里叶变换可以看成旋转π/2,一阶变换对应90度,而FRFT允许旋转任意角度,比如0.5阶就是旋转π/4。这个旋转特性让它在非平稳信号处理和图像加密中特别好用:信号在时域和频域之间连续过渡,阶数本身就是天然密钥。连续定义式里有一个线性调频核函数,理想情况下满足旋转可加性:先旋转a度再旋转b度,等价于旋转a+b度。这个性质是后面所有自检和加密解密流程的根基。

但到了离散世界,直接对积分公式做采样会出问题。最常见的问题是:分数阶傅里叶变换的离散版本不一定满足可加性,甚至不满足幺正性。也就是说,你用某个算法把信号变到分数域,再用负阶数变回来,得到的可能不是原信号。这不是Matlab代码的bug,而是离散化方案本身的约束。所以拿到这4种实现方式,第一件事不是急着用,而是搞清楚每种方式对点数、采样间隔和阶数的约束。

2.2 四种实现方式分别怎么算:采样型、特征分解型、线性调频卷积型、加权型

这个压缩包里最值钱的部分,就是把4种主流离散FRFT算法写在了一起。我按代码里常见做法给你梳理一下每一种的数学骨架和适用边界。

第一种是Ozaktas采样型快速算法。它基于这样一个思路:连续FRFT可以分解成一次线性调频相乘、一次快速傅里叶变换、又一次线性调频相乘。通过特定的采样间隔和尺度变换,把积分近似成离散卷积。这种算法时间复杂度是O(N log N),是4种里最快的,工程上最常用。但它的抽样间隔和输出域尺度依赖于阶数α,不是固定的,做图像加密时需要额外把输出归一化,否则横纵坐标的物理含义会漂移。

第二种是基于DFT特征分解的离散FRFT。这种做法的核心是构造一个离散变换矩阵,它的特征向量逼近连续Hermite-Gaussian函数,然后通过特征值的分数幂来得到任意阶变换。好处是正交性好,满足可加性和幺正性,数据可以无损重建,特别适合做加密中的可逆操作。坏处是矩阵规模是N×N,空间复杂度O(N²),大尺寸图像直接会内存爆炸,一般只用在单行或小尺寸块上。

第三种是线性调频卷积型离散FRFT。它把输入信号先调制到一个线性调频基上,再做卷积。实现上可以借助快速傅里叶变换把卷积变成频域相乘,速度和第一种接近。但这个算法对采样率和信号长度的一致性要求很苛刻,如果序列长度不是预期的点数,输出会出现明显的边界跳变。

第四种是加权型FRFT。它的思路更直接:把分数阶变换表示成四个整数阶傅里叶变换(0、1、2、3阶)的加权和,权重是阶数α的三角函数组合。代码好写、逻辑清楚,但代价是精度受限于基函数个数,实际上只能精确表示阶数间隔为1/4的变化。当α取0.2和0.21时,输出结果差别可能非常小,用来做加密密钥时区分度不够。

2.3 实现方式对比:精度、速度、正交性和时频支撑域

选哪种实现,本质上是选“计算速度”和“数学性质”哪个对你更重要。我给一张实际拆完代码后整理的对比表,方便你按需取用。

实现方式复杂度正交性可加性阶数分辨率适用场景
Ozaktas采样型O(N log N)弱近似满足高大尺寸信号、实时处理
特征分解型O(N²)严格满足严格满足高小尺寸块、加密可逆运算
线性调频卷积型O(N log N)中等边界受限高连续信号、长度固定
加权型O(N log N)严格满足严格满足低(1/4步长敏感)快速验证、整数阶近似

从这张表能直接看出一个结论:如果做图像加密,最保险的是特征分解型,因为它严格正交、严格可逆,解密时不会因为数值误差累积导致恢复图像出现条纹。如果处理的是上百万像素的整图,可以先考虑Ozaktas采样型,但要在加密解密两端用同一套归一化参数。我一般会在小尺寸测试时用特征分解型验证算法正确性,确认后再切换到采样型跑整个数据集。

3. 拆包与使用:Matlab代码结构、参数含义和第一个可运行例子

3.1 解包后能看到什么:文件清单与职责

这个RAR包解压后,目录结构大概是这样的:核心函数放在根目录,加密示例放在一个example子目录里。我按常见工程习惯猜一下,里面至少有这几个文件:frft_ozaktas.m、frft_eig.m、frft_chirp.m、frft_wfrft.m,分别对应该前面说的4种实现;可能还有frft2d.m做二维变换,以及frft_encrypt.m和frft_decrypt.m做图像加密演示;再加一个demo脚本用来跑通全流程。如果你打开后发现文件名带了版本前缀或日期后缀,说明打包者保留了多个迭代版本,使用时优先选不带后缀的最新文件。

文件职责上,frft_ozaktas.m里一般会有三个输入参数:输入序列x,阶数alpha,以及一个可选的量纲控制参数。frft_eig.m则是构建特征向量矩阵,输入需要传入变换点数N,它内部会先算离散概率密度算符的特征分解,这步耗时较长,但结果可复用。frft2d.m通常按“先对每一行做一维FRFT,再对每一列做一维FRFT”的方式实现,这种做法把二维变换拆成两个一维变换,复杂度可控,也好理解。

3.2 核心函数用法:从时域序列到分数域

我先把最常见的Ozaktas采样型函数怎么调、参数怎么传,给你写一个可直接运行的例子。假设你手头有一段长度为N的随机信号,想把它变换到0.5阶分数域:

N = 256; % 信号长度,建议为2的幂 x = randn(N, 1); % 输入一维信号 alpha = 0.5; % 分数阶,0~2之间,1.0就是普通傅里叶变换 y = frft_ozaktas(x, alpha); % 做0.5阶分数傅里叶变换 figure; plot(abs(y)); title('0.5阶FRFT幅度');

这段代码的含义很直观:frft_ozaktas接收原始时域序列x,返回分数域序列y。需要注意,alpha=0时输出几乎等于输入本身,alpha=1时等价于普通傅里叶变换,alpha=1.5对应逆变换再加一个相位翻转。参数N最好取2的幂次,否则采样型算法内部的快速傅里叶变换会做补零,导致输出序列长度和预期不一致。

参数说明上,还有一个容易忽略的点:Ozaktas型算法需要根据alpha重新计算采样间隔,所以它的输出y的坐标轴刻度不是固定的。你如果想把y画出来和理论解析解对比,一定要用代码内部返回的坐标轴向量,而不是直接用1:N当横轴。大多数版本的frft_ozaktas会返回两个输出,第二个就是对应的分数域坐标。如果打包里的版本只返回一个输出,我建议你自己改成返回坐标,不然后面做图像分析时会用错横轴,这一步是我实际踩过的坑。

3.3 图像加密的调用方式:把一维FRFT变成二维变换

图像是二维数据,直接套一维FRFT是没用的。常见做法是设计一个二维分数傅里叶变换,对每个行向量和列向量分别做一维变换。下面给出一个可复用的函数框架,它接收灰度图像矩阵img和两个方向的阶数alpha_x、alpha_y:

function Y = frft2d(img, alpha_x, alpha_y) [M, N] = size(img); Y = zeros(M, N); % 先对每一行做一维FRFT for i = 1:M Y(i, :) = frft_ozaktas(img(i, :)', alpha_x)'; end % 再对每一列做一维FRFT for j = 1:N Y(:, j) = frft_ozaktas(Y(:, j), alpha_y); end end

这里面有一个重要细节:第一次做行变换时,img(i,:)是行向量,而frft_ozaktas内部大概率要求输入是列向量,所以需要先转置,变换完再转置回来。第二次列变换就不用转了。如果你两种变换的顺序换过来,先列后行,结果和先前后列是一样的,前提是阶数不变。但如果alpha_x和alpha_y不相等,行列顺序换了就会得到完全不同的加密图,这正是后面加密密钥的一部分。

我再强调一下这个二维变换的数学含义:行变换代表水平方向在分数域的旋转,列变换代表垂直方向的旋转。两个方向用不同的阶数,相当于把图像在一个四维时频空间里扭了一下。加密时,只要解密端把alpha_x和alpha_y倒过来取负号,就能还原。这个思路和双随机相位编码结合起来,就是第4章要讲的完整加密流程。

4. FRFT图像加密的参数陷阱与常见问题排查

4.1 阶数与点数不匹配:输出长度变了怎么办

现象:对一张512×512的图像做二维FRFT,加密结果尺寸还是512×512,但解密后图像尺寸变成513×513或者256×256,矩阵维度对不上,无法重建。

原因:绝大多数一维FRFT实现内部用快速傅里叶变换辅助计算,如果输入长度不是2的幂,或者阶数导致采样间隔变化,输出点数和输入点数不一致。尤其是Ozaktas采样型,它的输出长度由内部重新取样规则决定,不是简单地保持原长。

解决:在调用一维FRFT前强制检查长度,如果不是2的幂,先做边界扩展或中心裁剪。另一种做法是我个人更推荐的:在frft2d函数里统一对每行每列裁剪到原始尺寸,即取变换输出前M个点或后M个点,保证矩阵始终是M×N。注意解密时也要做同样的裁剪,不然相位信息和幅值信息错位。

4.2 图像加密后反攻失败:矩阵转置与左右乘顺序

现象:加密图能出来,但解密图是原图的左右镜像,或者旋转了90度,甚至出现棋盘格状错位。

原因:二维FRFT很多人直接用矩阵乘法形式写,比如把图像左乘一个变换矩阵,然后再对结果右乘另一个变换矩阵。如果写成了左乘两次,或者行变换和列变换的转置关系弄反,相当于加密时先左乘再右乘,解密时又先右乘再左乘,最终得到的是转置后的结果,对于非对称图像看起来就像镜像。

解决:写一个自检脚本,输入已知图像I,先做一次frft2d,再立刻用负阶数做frft2d,对比重建结果与I的差值。如果差值不是接近0而是左右翻转,检查你的行变换代码里是否少了转置。我一般会把加密代码固定成“先所有行、再所有列”,解密代码固定成“先所有列、再所有行”,方向相反但矩阵维度保持一致,这种对称顺序不容易出错。

4.3 相位掩码的随机性:固定随机种子与密钥复现

现象:同一张图,同一个阶数,连续加密两次得到完全不同的加密结果,而且用第一次的密钥无法解密第二次的结果。

原因:加密流程里引入了随机相位掩码,生成时用了Matlab的rand或randn,却没有固定随机种子。每次运行rand函数都会重新生成不同的随机矩阵,导致掩码变了,密钥自然对不上。

解决:在加密脚本最前面加一行rng(某个整数),比如rng(2024)。这样做的好处是,只要随机种子和阶数不变,任何一次运行生成的掩码都一样,解密端读取到一模一样的掩码。如果你需要给不同用户发不同密钥,就把种子本身当作密钥的一部分发送,而不是把整个掩码矩阵传过去,这样传输体积会小很多。另外要注意,加密端和解密端分属不同脚本时,两边的rng调用顺序必须完全一致,因为rand、randn、randi共用同一个随机流,如果你在设种子之后又调用了其他随机函数,掩码值就会偏移。

4.4 阶数取整误差导致解密图像出现条纹

现象:加密时alpha_x=0.3,解密时直接写alpha_x=-0.3,结果恢复图像模糊,伴随周期性亮纹。

原因:某些实现方式(尤其是加权型)对阶数的离散分辨率有限,0.3和-0.3虽然绝对值相同,但内部三角函数权重计算时可能因为浮点误差映射到了不同的等效阶数。另外,部分实现把逆变换定义为阶数取负后的正变换,如果实现里没有专门处理负阶数,就会直接用alpha的绝对值,导致解密变成另一个角度。

解决:先用第5章的自检方法测试该实现是否支持负阶数精确重建。如果不行,切换成特征分解型实现,它严格支持负阶数。或者在加密时记录alpha,解密时用2-alpha来代替代数减法?不,正确做法是保证解密阶数是加密阶数加2的整数倍,即alpha_dec = -alpha_enc + 2*n,其中n取整数。因为FRFT的周期是4,阶数加4后结果不变,所以-0.3 + 2 = 1.7在数学上等价于-0.3。但不同实现对这个等价性的支持程度不一,一定要先验证再大量跑数据。

4.5 图像加密后信息丢失:加密图亮度分布异常

现象:加密后的图像看起来不是随机噪声,而是保留了原图的轮廓,或者整体偏暗、偏亮,统计分析能看出相关性。

原因:用的FRFT实现不具备能量保持性质。Ozaktas采样型如果没有做归一化,变换后信号的L2范数会随阶数变化,导致加密图像的动态范围被压缩或扩展。反映在图像上就是亮度分布不均匀。

解决:在加密前先对图像做归一化,比如把像素值从0~255映射到0~1,加密后输出是复数,直接取幅度再映射到0~255显示。解密时要把复数结果反归一化回原始范围。这是加密流程里最不起眼但最关键的一步,很多翻车的加密Demo都是栽在这里。我一般会在加密函数开头加一句断言,计算输入输出能量比,偏差超过5%就直接报错,避免把误差带到下一轮。

5. 进阶验证:用双参数可逆变换自检FRFT实现

拿到这包代码后,我建议你先别急着做加密,先跑一个最小自检,确认你选的那个实现方式真的满足“先正后负”可逆。具体做法是:构造一个随机复信号,依次做alpha和-alpha的变换,计算恢复信号和原信号的归一化均方误差。误差低于1e-10说明该实现高度正交,低于1e-4也还能用,如果误差超过1%,说明这个函数在这个阶数和长度下不能直接用于可逆加密。

x = randn(128, 1) + 1i*randn(128, 1); alpha = 0.37; % 一个非特殊角度,检验普适性 y = frft_eig(x, alpha); % 先用特征分解型 x_re = frft_eig(y, -alpha); % 再用负阶数恢复 err = norm(x - x_re) / norm(x); fprintf('重建误差: %.2e\n', err);

我习惯把这步叫“后悔药测试”——如果在编码阶段就能把误差暴露出来,就不用在加密解密流程里排错半天。对于图像,还可以把lena图先做0.5阶、再做-0.5阶,显示差值图像,肉眼确认没有任何边缘残留。从那以后,我每次换一种FRFT实现,第一件事就是跑这个可逆性自检,确认阶数敏感性和维度一致性没问题后,再进入加密流程。这条习惯帮我挡住了至少三个看起来正常、实际不可逆的算法变体。希望帮到你。

本文还有配套的精品资源,点击获取

RELATED READING

延伸阅读

更多一线实战笔记与深度复盘,助您持续精进