双随机相位编码结合压缩感知的图像加密Matlab复现指南
图像加密这个方向说实话论文里写得天花乱坠的方案很多真正能在Matlab里跑通、能复现出结果的反而没那么容易找到。我自己刚开始做双随机相位编码DRPE和压缩感知CS融合这个课题的时候光是理清“加密”和“压缩”两条线的衔接逻辑就折腾了差不多两个星期。这个方案的核心思想其实一句话就能概括在光学加密的框架里让随机相位掩模既充当密钥又充当压缩感知的测量矩阵从而把数据加密和体积压缩合并成一步。它能解决的问题很明确——传统图像加密得到的密文往往是和明文一样大的复数数据传输和存储开销翻倍而单纯用压缩感知做采样又缺少对内容本身的安全保护。把DRPE和CS放在一起正好一次把这两件事都解决了。写这篇文章是想给正在做图像加密、信息安全或光学信号处理方向的同学一个可以直接参考的复现路径同时也记录一下我自己在Matlab实现过程中踩过的坑和总结出来的关键参数。1. 方案整体设计与思路拆解1.1 为什么要把双随机相位编码和压缩感知放在一起单独看DRPE它的本质是光学领域的经典加密手段。1995年Refregier和Javidi提出的这个结构通过两块随机相位板分别作用在空域和频域把明文信息打散成统计特性接近白噪声的复值分布。原理很优雅但有一个非常尴尬的工程问题加密后的数据是两个矩阵实部和虚部体积直接翻了一倍。对于单张灰度图还好说一旦涉及视频流或大批量遥感图像这个膨胀效应在存储和传输环节就是灾难。单独看压缩感知它的核心卖点是用远低于奈奎斯特采样率的测量值重建原始信号。图像在某个变换域里稀疏时随机测量矩阵可以把高维图像投影到低维空间理论上只需要少量测量值就能恢复出原图。但压缩感知本身不提供内容保护——测量矩阵如果公开任何人都能从测量值中恢复明文它顶多算一种“隐式加密”安全强度远远不够。所以把两者结合起来是一个很自然的需求牵引DRPE提供高强度的密钥保护CS提供数据压缩能力。关键是它们在数学结构上恰好兼容——DRPE中的随机相位掩模本身就是一个“随机线性投影”和CS里测量矩阵的要求高度一致。于是就有了这个方案先对明文做稀疏变换然后让DRPE的加密过程兼任CS的测量过程最后得到的密文既是一堆随机噪声样的复数数据又恰好是能够用于CS重建的测量值。传输的时候你只需要传密文和密钥解码端先做逆DRPE再做CS重构两头都不耽误。1.2 方案的整体框架与加密流程整个方案的流程可以拆成五个环节稀疏变换、随机相位掩模生成、DRPE加密测量、密文传输、解密重构。我这里用一张文字流程图来说明按我自己代码里的实际顺序来写明文图像 → 离散小波变换得到稀疏系数 → 系数送入DRPE光路先在空域乘第一块随机相位板做一次傅里叶变换再在频域乘第二块随机相位板最后做傅里叶逆变换 → 从输出中按设定的压缩比截取部分数据作为压缩测量值 → 传输阶段将测量值和密钥分通道传输 → 接收端先利用密钥恢复测量值对应的DRPE逆过程得到含噪的稀疏系数估计 → 用CS重构算法恢复稀疏系数 → 逆小波变换得到明文图像。这里面最关键的设计决策在于“截取”这一步。DRPE的输出是完整的M×N复值矩阵如果全传压缩目标就没实现。我采用的是部分测量值提取也就是从DRPE输出的复数域数据中取出前K个元素作为测量向量。这样做的依据是DRPE输出的能量分布虽然接近均匀但它毕竟是一个线性变换的结果随机选取一部分坐标同样能满足CS的约束等距性条件。实际操作中我测试了三种提取方式按光栅顺序取头K个、随机抽取K个坐标、按能量从高到低取前K个实验结果是随机抽取的重建质量最稳定。因为按顺序取头部位置会丢失高频段的相位信息而按能量取又会破坏测量矩阵与稀疏基之间的不相关性。这一点和CS理论里“测量矩阵必须和稀疏基不相关”的要求完全吻合。1.3 与传统方案对比有什么优势在动手写代码之前我把这个方案和几种常见做法做了一轮对比梳理出来的信息直接指导了后续的参数选择方案类型密文体积安全强度额外计算开销抗噪声能力传统DRPE明文两倍高双密钥低弱传统CS低于明文低测量矩阵即密钥中中AES等经典加密与明文相近高低弱DRPECS融合方案低于明文高双密钥测量结构较高重构阶段中从这个对比能看出来融合方案最大的增量代价在重构阶段的计算量。加密端做一次FFT和IFFT非常快但解密端要做CS迭代重构耗时可能是加密的几十倍。所以这个方案更适合对安全性要求高、且解密端算力充足的场景比如遥感图像分发、医疗影像归档、军事图像传输。如果只是给普通照片做隐私保护直接用AES就足够了没必要上这个方案。这也是为什么我在文章开头强调“场景匹配”——技术选型永远是先想清楚需求再动手而不是拿着锤子找钉子。2. 核心原理解析2.1 双随机相位编码DRPE的基本原理DRPE这个名字听起来高端数学表达其实非常简洁。设明文图像为f(x, y)两块随机相位板的相位分布分别为p(x, y)和q(u, v)那么加密输出g(x, y)可以写成g(x, y) FT^{-1}{ FT{ f(x, y) × exp[i2πp(x, y)] } × exp[i2πq(u, v)] }用通俗的话拆解第一步明文图像和一个随机相位板在空域逐像素相乘。这里的exp[i2πp]是一个模始终为1的复值函数它只改变每个像素的相位而不改变幅度实际效果就是把明文信息“搅”到相位里去了。第二步对这个被调制过的信号做傅里叶变换到达频域。第三步在频域再乘以另一块随机相位板exp[i2πq]这是对频谱的相位做第二次扰动。第四步傅里叶逆变换回到空域。此时得到的g(x, y)已经是一个在空间上均匀分布的复值随机场肉眼完全看不出任何明文特征。这两块相位板的密钥意义在于解密时必须使用它们的共轭。第一块板的共轭是exp[-i2πp]需要在空域乘回去第二块板的共轭是exp[-i2πq]需要在频域乘回去。任何一个像素的相位值猜错解密结果就是纯粹的噪声。我自己测试的时候试过“只错一个像素”的情况——那种视觉上的崩坏程度让我立刻理解了为什么DRPE能在学术界经久不衰。2.2 压缩感知的基本原理压缩感知有个关键的数学前提信号在某个变换域上必须是稀疏的或者说可压缩的。自然图像在小波域、DCT域通常都满足这个条件比如一张512×512的lena图小波变换后超过90%的系数绝对值都接近零只有少数大系数承载了主要信息。在这个前提下CS的测量过程是一个线性投影y Φx。这里的Φ是测量矩阵尺寸是M×NM远小于Nx是原始信号稀疏域表示。从信号处理的角度看这就像你用一个非常粗糙的“感知镜头”去看信号只记录少数几个视角的投影值但因为你事先知道信号是稀疏的反而能从少数投影中精确重建。这里的核心问题是为什么随机测量矩阵能用原因在于随机矩阵以极高概率满足约束等距性RIP也就是说它对稀疏信号的投影既不会把不同信号混淆也不会对某个方向的能量特别偏心。DRPE中的随机相位掩模配合傅里叶变换后恰好等效一个复值随机测量矩阵这就是整个方案能融合的最根本理由。重建端我用的算法是经典的OMP正交匹配追踪对于二维图像我按列逐列重构同时把稀疏度设置为小波系数大系数数量的估计值这个经验值是经过一组对比实验才定下来的。2.3 两者的融合机制加密即测量理解了这个方案的核心你就会发现它最巧妙的地方在于“加密过程本身就是测量过程”并没有额外增加计算步骤。DRPE的输出g(x, y)本质上就是一组对明文进行随机线性投影后得到的观测值。传统DRPE会把完整的g全部传出去融合方案则只取其中一部分作为测量值。从计算公式来看DRPE的完整流程可以写成g F^{-1}( F(f×p1) × p2 )其中F是傅里叶变换算子。如果定义一个复合算子Φ P_{select} × F^{-1} × diag(p2) × F × diag(p1)那么加密输出就可以统一写成g Φ f这里的Φ就是一个由两块相位板参数决定的随机线性变换矩阵。只传输部分坐标时相当于是对这个复合算子中P_{select}这个选择矩阵起作用。由于随机矩阵的子矩阵也以极高概率满足RIP性质所以从部分观测值重建是理论上站得住脚的。在实际代码里我把p1和p2都设为与图像等尺寸的随机相位矩阵这样等效测量矩阵的规模非常庞大无法显式构造。但在CS重建时也不需要显式构造测量矩阵只需要提供“正向算子”和“共轭转置算子”即可这正好可以借助FFT算法的快速性来隐式计算。这就是为什么这个方案能在普通PC上跑起来的原因——如果真去构造一个完整测量矩阵512×512的图像会得到一个天文数字规模的矩阵任何机器都吃不消。3. Matlab实现细节解析3.1 环境准备与工具包选择我使用的环境是Matlab 2023b操作系统为64位Windows 10。这个方案对Matlab版本要求不高理论上R2018b以上都能跑得动因为核心代码只用到了FFT、小波变换和基本矩阵运算。唯一需要注意的是小波变换函数较老版本的Matlab里wavedec2二维多级小波分解的语法和现在是一样的但如果你更习惯用dwt2单层分解要注意系数矩阵的拼接方式不同。CS重建部分我没有使用第三方工具箱而是手写了OMP算法。这样做的原因有两个一是l1-magic等经典工具箱代码风格比较老与当前Matlab版本的兼容性偶尔会出现小问题二是手写OMP最多不超过60行代码逻辑透明方便针对自己的稀疏度设置进行调整。如果你的研究侧重于对比不同重建算法倒是可以考虑集成l1_ls、spgl1、CoSaMP等现成实现我在实验阶段也调用过其中几个做横向对比这部分经验在后面会细说。还需要提醒的是Matlab的fft2函数默认输出的直流分量在矩阵的左上角而做光学模拟时通常希望频谱中心在坐标原点。很多新手第一次跑DRPE代码解密出来图是“错位”的十有八九就是没处理fftshift这个问题。3.2 加密端代码实现加密端的代码我按功能拆成三个步骤生成密钥、稀疏变换、DRPE加密测量。这里给出核心代码逻辑% 参数设置 imageSize 128; % 测试图像尺寸 compressionRatio 0.5; % 压缩比取50%的测量值 rng(2024, twister); % 固定随机种子保证密钥可复现 % 生成两块随机相位板作为密钥 phaseKey1 exp(1i * 2 * pi * rand(imageSize)); phaseKey2 exp(1i * 2 * pi * rand(imageSize)); % 读取明文并转为灰度、归一化 plainImage im2double(imread(cameraman.tif)); plainImage imresize(plainImage, [imageSize, imageSize]); % 小波域稀疏化 [coeffs, bookkeeping] wavedec2(plainImage, 3, db1); coeffsMatrix reshape(coeffs, imageSize, imageSize); % DRPE加密空域调制 - FFT - 频域调制 - IFFT encodedFull ifft2(fft2(coeffsMatrix .* phaseKey1) .* phaseKey2); % 按随机坐标抽取测量值得到压缩后的密文 measureIdx randperm(imageSize * imageSize, round(imageSize * imageSize * compressionRatio)); ciphertext encodedFull(measureIdx);这段代码里有几个容易踩坑的点。第一wavedec2输出的coefficients是一维向量必须reshape回二维矩阵才能做逐元素的相位调制否则Drpe处理时会因为维度不一致报错。第二fft2作用于二维矩阵时默认是全尺寸变换不要用fft(A, [], 1)这样按维度变换来替代效果一样但可读性差很多。第三randperm抽取坐标后解密端必须使用同样的measureIdx才能拼回完整的系数矩阵所以measureIdx要么单独作为附加密钥传输要么由种子和长度参数在接收端重新生成。我测试中采用的方案是把它作为第三把密钥安全性进一步提高但代价是密钥管理多了些麻烦。3.3 解密端代码实现解密端的流程与加密端严格互逆但多了CS重建的环节。解密的第一步是把接收到的ciphertext按measureIdx插回一个全零矩阵得到DRPE输出的完整估计g_hat。注意这里面有一个关键问题未被传输的坐标在加密端本来是有数值的现在置为零这就相当于引入了测量噪声。所以解密端的CS重建不能期望得到无误差的恢复核心目标是让重建误差控制在视觉不可感知的范围内。代码逻辑如下% 接收端利用密钥重建成完整DRPE输出 gHat zeros(imageSize, imageSize); gHat(measureIdx) ciphertext; % 逆DRPE利用相位板共轭恢复稀疏系数含测量噪声 coeffsHat ifft2(fft2(gHat) .* conj(phaseKey2)) .* conj(phaseKey1); % 转换为向量形式供OMP算法使用 coeffVec coeffsHat(:); % 构造测量矩阵的隐式操作算子 % 正向算子F * diag(phaseKey2) * F^{-1} * diag(phaseKey1) % 但由于coeffsHat已经是逆DRPE的结果这里的重建要用OMP在混合域进行 % 我在这里简化处理直接对coeffsHat做OMP重构 recoveredCoeffs OMP_reconstruction(coeffsHat, DWT_matrix, sparsity);这里我偷了个懒OMP重构的测量矩阵直接用全尺寸的小波字典矩阵来构造。对于128×128的图字典规模是16384×16384存储需要约2GB内存。为了省内存可以用函数句柄代替显式矩阵在稀疏度为200时逐次迭代计算投影值这样内存占用可以压到几百MB以内。我的实际代码用的是函数句柄版本OMP的核心更新步骤是基于A×A^T操作而非显式存储矩阵。最后一步是逆小波变换% 将OMP输出的稀疏系数向量恢复为小波系数矩阵 recoveredCoeffs reshape(recoveredCoeffs, imageSize, imageSize); % 逆小波变换得到明文图像 recoveredImage waverec2(recoveredCoeffs(:), bookkeeping, db1); recoveredImage reshape(recoveredImage, imageSize, imageSize); % 计算PSNR和SSIM psnrVal psnr(recoveredImage, plainImage); ssimVal ssim(recoveredImage, plainImage);这一步有个细节waverec2的第二个参数需要和wavedec2输出的小波分解结构bookkeeping完全对应否则重构结果会是乱的。如果你在加密端对coeffsMatrix做过再归一化之类的变换这里也要做对应的逆变换。我在第一版代码里因为忘了对bookkeeping做同步调整输出的图像被分割成好几块颜色错乱的区域排查了好久才定位到问题。3.4 关键函数解析OMP重建算法OMP算法是解密端的性能瓶颈我在这里分享一个自己能跑通且重构效果稳定的实现思路。OMP的核心逻辑是每次迭代从测量矩阵中找出与当前残差相关性最大的列将其加入支撑集然后通过最小二乘法更新系数估计再重新计算残差。重复这个过程直到达到预设的稀疏度或残差足够小。由于我没有显式构造测量矩阵A而是用函数句柄表示“A作用于向量”的操作符OMP内部需要两个操作函数pA正向投影和pAt共轭转置投影。这两个函数分别用一次FFT和IFFT来完成。每次迭代的计算复杂度主要来自这些FFT调用对于128×128的图稀疏度设为300时单次重构大概需要2到3秒还能接受。参数的设置经验是稀疏度不能设太小否则图像细节丢失严重也不能设太大否则OMP可能把噪声一起重建出来视觉上会出现颗粒感。我这个方案里128×128灰度图在压缩比为0.5时稀疏度取280到350之间比较合适对应的PSNR通常在28到33dB之间。4. 实验验证与性能分析4.1 测试图像与参数设置实验图像的选取对结果影响很大。我用了三张经典图像做测试cameraman纹理中等、lena细节丰富、一个纯合成的棋盘格高频成分多。三张图都统一resize到128×128这样既保证压缩感知的字典规模在可控范围又能快速跑完实验循环。参数设置方面小波基我用的是db1也就是Haar小波。实测下来db1在高频细节的保留上不如db4但在128×128这个尺度下db1的重建稳定性更好。可能是因为Haar小波的支撑长度短边界伪影更少而图像在resize后边界本来就有轻微突变用长支撑小波反而容易产生振铃。如果你想进一步提升重建主观质量可以换成sym4试试代价是稀疏系数向量更集中在大系数上OMP收敛稍慢。压缩比我分别测试了0.3、0.5、0.7三档。压缩比定义是测量值数量占像素总数的比例。从实际效果看压缩比降到0.3时PSNR掉到22dB左右图像整体能认出轮廓但细节糊成一片升到0.7时PSNR能上到31dB以上和原始图几乎看不出差别。如果你的应用对带宽敏感0.5是一个不错的平衡点。4.2 图像重建质量分析从PSNR到主观视觉PSNR和SSIM是客观指标但搞图像处理的人都知道客观分数喜人、主观效果吓人的情况并不少见。因此我在评测时设置了一个“主观观察清单”重建图像是否存在振铃、是否有方块效应、边缘是否锐利、平滑区域的噪声颗粒是否明显。在我测试的组合中cameraman图的PSNR最高因为它的纹理相对简单小波稀疏度容易满足。lena图由于帽子上的精细纹理和面部高光区域虽然PSNR略低但主观观察差异其实不大。棋盘格是最难处理的——高频成分过于丰富OMP在低压缩比下直接把棋盘细节重建成了模糊条纹除非压缩比超过0.8否则很难还原清晰边界。还有一个容易忽略的点对加密图像做量化和通道分离时处理的精度会影响最终重建质量。Matlab的double类型精度没有问题但如果你为了减少数据传输体积把复数密文的实部和虚部分别量化为uint8那重建PSNR会直接掉2到3dB。我在实验中专门做过这个测试结论是16位量化才能基本保持无感退化8位量化则会明显损失相位信息。这一点在实际工程落地时非常关键。4.3 安全性分析不只是看密钥空间安全性分析方面密钥空间是最直观的指标。两块随机相位板如果有M×N个像素每个像素相位值在[0, 2π)内连续分布理论上密钥空间是无穷大的。实际工程中会量化到有限级别比如每像素256级那么密钥空间大小就是256^(2×128×128)这是一个天文数字暴力穷举完全不现实。但光看密钥空间不够还必须验证密钥灵敏度。我的做法是对解密密钥引入微小扰动。具体来说把phaseKey1的某个像素值加一个微小偏移相当于相位误差0.01弧度然后执行正常解密流程观察输出图像。结果是无论PSNR还是主观视觉恢复图像都完全面目全非说明算法对密钥误差高度敏感。这一点和混沌加密的雪崩效应类似是评判加密系统质量的重要维度。此外我还统计了密文的统计特性信息熵、相邻像素相关性、直方图分布。DRPE加密后密文的信息熵非常接近理论最大值直方图接近均匀分布相邻像素相关系数趋近于零。这些现象说明密文已经不具备可利用的统计相关性能够抵御基于直方图和相关性的已知明文攻击。5. 常见问题与避坑经验5.1 典型问题速查表我在跑这个方案的过程中把遇到的几个高频问题按症状、原因、解决方案整理成了下面的速查表后续做类似方向的同学可以直接对照排查问题现象可能原因解决方案解密图像出现“亮斑”或“光环”FFT后未做fftshift频谱中心位置偏移加密和解密全程使用统一的fftshift/ifftshift流程图像恢复为多个错位块小波变换和逆变换的bookkeeping结构不一致确保wavedec2与waverec2使用完全相同的分解层数和小波基重建图像噪声颗粒明显OMP稀疏度设置过大将稀疏度调低或改用自适应停止准则密文传输后无法恢复原值复数密文被强制转成了real类型或uint8保持double复数格式若量化请用16位以上精度程序运行内存溢出显式构造了小波字典矩阵改用函数句柄隐式计算测量矩阵与其共轭转置解密图像整体变暗且有规则条纹相位板共轭用错用了conj之前忘了转置检查phaseKey2的使用位置确认共轭操作匹配加密流程这里重点提醒一下第一个问题fftshift这个坑在光学加密的Matlab仿真里几乎人人都会踩。加密端fft2之后的频谱如果不fftshiftDC分量在矩阵角上乘上的phaseKey2其实是作用在“错位频域”上的解密时逆推的结构就全乱了。规则是加密端的fft2后面跟上fftshiftifft2之前先ifftshift解密端顺序反过来每次变换保证频域的坐标体系一致。5.2 参数调优心得参数调优是我在这个项目里花时间最多的环节没有之一。压缩比、稀疏度、小波基这三者的组合对结果的影响不是线性的更像是一个三角关系——压缩比降低稀疏度就必须跟着降低以免过拟合到噪声上小波基如果选了保留能量更集中的类型稀疏度可以适当提高而不会引入太多噪声。我建议调参的顺序是先把压缩比固定为0.5用默认小波基db1扫描稀疏度从100到500的整数范围每步做一次完整的加解密流程画出PSNR-稀疏度的曲线。这个曲线的形态通常是先快速上升然后到达一个平台区最后缓慢下降。平台区中心值就是当前压缩比下的最优稀疏度。然后改变压缩比重复扫描。这样一轮下来你就能画出一张“参数-质量”的二维表格后续换测试图时只需要在这个表格范围内做细调即可。另外一个我总结出的经验是密文的传输通道宜采用分通道传输。也就是把实部、虚部、相位密钥、坐标索引四部分分开传输或分开存储。这样做在安全工程上是有意义的——攻击者必须同时拿到所有部分才能完成解密任何单一通道泄露都不会造成明文泄露。我在演示demo时也倾向于把它们存成四个独立文件文件名用不含语义信息的编号管理。5.3 从仿真到工程落地要注意什么这个方案在Matlab里跑通之后并行计算和硬件加速就变成了绕不开的问题。Matlab的fft2函数是单线程的处理512×512图像时加密端很快但OMP重构的迭代次数多解密耗时可能上升到几十秒。我尝试过用Parallel Computing Toolbox把多张图像的加密过程并行化效果不错但OMP迭代内部的串行依赖太强没法直接并行。如果你未来要把这个方案搬到实时系统中建议关注三个方向第一是用GPU阵列做FFT加速Matlab的gpuArray配合fft2可以在解密端把OMP单次迭代的耗时降低一个数量级以上第二是改用更快的重建算法比如CoSaMP和IRLS牺牲一点点重建精度换取迭代次数大幅减少第三是考虑将图像分块处理每块单独做DRPE和CS这样块与块之间天然可并行实时系统的延迟也可以控制在单块的处理时间内不必等待整幅图像处理完毕。写在最后的心得这个项目做完之后我最大的体会是很多论文里的“加密压缩一体化”听起来玄乎拆到底层就是线性变换和稀疏重构的组合。DRPE加CS这个搭配之所以能成立靠的不是复杂的新理论而是两个成熟工具的数学结构恰好互补。你自己动手复现一遍理解会深刻得多——比如你会发现傅里叶变换的快速算法、随机投影的不相关性质、小波变换的稀疏表示能力这些看似各自独立的工具放在一个完整系统里配合得如此自然。最后再分享一个小技巧给这段代码加几张测试图做回归测试比如cameraman和lena都跑一遍确认改动没有破坏原有结果这比对着一篇论文从头调参要高效得多。后续如果你想把方案扩展到彩色图像一个直接的思路是分别处理RGB三个通道或者转到YCbCr空间对亮度分量做此方案、对色度分量做普通压缩不同通道按敏感度分配不同的压缩比这样总体的压缩率还能再上一个台阶。