二维图像压缩感知重构：MATLAB完整实现与OMP算法解析-北京尧图网络

简介一个面向二维图像压缩感知CS重构算法的MATLAB实现代码包适合电子信息、计算机视觉与医学成像等领域的研究者、工程师及学生用于算法验证与二次开发。压缩包共26个文件包含18个.m脚本覆盖测量矩阵生成、OMP/ISTA等重构算法、PSNR/MSE性能评估、7个.fig结果图与1个示例bmp图像整体大小1.59MB结构清晰便于按模块调用。目前已有4296人学习下载口碑良好。代码基于小波域稀疏表示实现了OMP、Wavelet_OMP等重构流程并提供可运行的CS_OMP演示脚本配有详细注释帮助读者直观理解压缩采样与重建全过程。通过调整测量矩阵类型、稀疏基与重构参数还可进一步拓展到视频压缩、医学影像及噪声环境下的鲁棒性研究是学习压缩感知理论和开展图像重构实验的实用工具。1. 二维图像压缩感知重构一份能直接跑通的MATLAB完整实现如果你手里只有一张 256×256 的灰度图按传统采样你得存下 65536 个像素点但用压缩感知的思路只需要其中约 30% 的随机观测就能把整张图恢复出来肉眼几乎看不出差别。这份压缩感知重构算法 MATLAB 程序代码就是把这套「欠采样—稀疏重构」的完整链路封装成了可运行的工程从二维图像的稀疏变换、随机观测矩阵生成到 OMP正交匹配追踪迭代重构、质量评估全部是完整的 MATLAB 代码没有缺函数、没有跳步骤。适合正在做信号处理课设、研究稀疏重构算法、或者想在工业采集场景里用少量传感器数据重建图像的同学——你拿到的不是一段演示脚本而是一个可以改参数、换算法、换测量矩阵的完整实验框架。2. 技术底座稀疏变换、观测矩阵与重构算法的三角关系2.1 压缩感知为什么能“少采样多恢复”传统采样的底层逻辑是「先采满、再压缩」先把图像的所有像素点采下来再用 JPEG 之类的编码器把冗余信息丢掉。压缩感知反着来它要求信号本身在一个变换域上是稀疏的——也就是说图像在 DCT 或小波基下展开时绝大多数系数接近零只有少数大系数承载了主要信息。数学上我们只需要用 M 个随机线性投影M 远小于像素总数 N就能以高概率恢复出这 K 个非零系数。这里有个前置条件K 必须远小于 M也就是稀疏度远小于观测数。这份代码默认图像是灰度图像素值先归一化到 0~1然后通过 DCT 变换把空间域的像素值映射到频域系数大多数能量集中在低频位置天然满足稀疏性前提。代码里的完整链路可以分四步理解稀疏表示用 DCT 变换基将图像从空间域映射到频域得到系数矩阵。观测压缩用高斯随机矩阵对系数做线性投影得到远小于原图尺寸的观测向量 y。稀疏重构OMP 算法从 y 中迭代恢复出原始系数。逆变换回图像对恢复的系数做逆 DCT还原像素域图像。2.2 观测矩阵构建为什么是高斯随机矩阵而不是单位阵观测矩阵是压缩感知里最「讲究」的一个环节。它必须满足有限等距性质RIP通俗说就是它不能破坏原始信号的几何结构——两个不同的稀疏信号经过同一观测矩阵后在低维空间里仍然能区分开。高斯随机矩阵是被数学证明过满足 RIP 的一类矩阵而且实现简单用randn生成一个 M×N 的矩阵每一列都服从标准正态分布即可。这份代码里的做法是生成后做列归一化处理避免某些列的能量过高对重构结果产生偏置。function Phi gaussian_measurement(M, N) % 生成高斯随机观测矩阵 % 输入: M - 观测数量, N - 信号长度 % 输出: Phi - M行N列的观测矩阵 Phi randn(M, N); % 列归一化: 每一列除以自身二范数, 保证各列能量一致 for col 1:N Phi(:, col) Phi(:, col) / norm(Phi(:, col)); end end列归一化这一步不是压缩感知理论里的必须项而是工程上的常见做法。如果不做归一化某些列的模长天然偏大OMP 在选原子时会优先选中这些列哪怕它们和残差的相关性并不高归一化之后原子的选择就纯粹由「与残差的内积大小」决定重构稳定性明显提升。我在实践里遇到过不归一化时 PSNR 掉了 2~3dB 的情况所以这个细节值得保留。2.3 稀疏度 K 和测量率 MR 的关系这份代码里有两个参数决定压缩程度测量率 MR观测数占总像素数的比例和稀疏度 KOMP 迭代次数。测量率通常设置在 0.2~0.5 之间对应的观测数为 M round(MR × N)。稀疏度 K 则代表你认为图像在 DCT 域上有多少个显著系数代码里一般取观测数的 1/2 到 1/3。K 设大了重构出的系数带噪声K 设小了细节丢失图像发糊。两者之间的关系是M 越大能容纳的 K 越大恢复质量越高。% 主程序中的参数设置 N 256; % 图像边长 MR 0.3; % 测量率: 取30%的观测数 M round(MR * N * N); % 总观测数 K round(M / 4); % 稀疏度: 约为观测数的1/4这里的 K 取 M/4 是一个偏保守的默认值适合大多数自然图像。如果你的图像纹理非常丰富比如一堆密集的树叶DCT 系数衰减慢K 需要适当调大如果是纯色背景的简单图像K 取小一些反而更干净。调试方法很简单——固定其他参数让 K 从 M/8 往上加每加一档记录一次 PSNR你会看到 PSNR 先升后降峰值点就是当前图像的最优稀疏度。3. 核心模块拆解OMP重构算法的MATLAB代码逐段解读3.1 主程序 main.m 的整体控制流主程序负责把各个子函数串联成一条完整的实验流水线。它读取图像、做预处理、调用观测和重构函数、最后输出评价指标。理解主程序的执行顺序是二次开发时改动逻辑的入口。%% 加载图像并转为double灰度图 img imread(cameraman.tif); if size(img, 3) 3 img rgb2gray(img); % 若是彩色图则转灰度 end img im2double(img); % 归一化到[0,1] [N, N] size(img); % 获取边长, 要求方形图 %% 得到稀疏系数 T dctmtx(N); % 生成N×N的DCT变换矩阵 coeffs T * img * T; % 二维DCT变换, 系数集中在左上角低频区 %% 生成观测矩阵并采样 Phi gaussian_measurement(M, N); % 这里对系数矩阵做向量化, 方便做矩阵乘法 coeff_vec coeffs(:); % N^2×1 列向量 y Phi * coeff_vec; % M×1 观测向量 %% OMP重构 coeff_hat zeros(N*N, 1); % 分块重构: 对N×N系数矩阵按列或按行分块, 每块独立跑OMP % 常见做法是对每一列单独做OMP, 避免一次求解维度过大 for col 1:N y_col y((col-1)*M/N1 : col*M/N); % 取当前列对应的观测分量 Phi_col Phi((col-1)*M/N1 : col*M/N, :); % 对应观测子矩阵 coeff_hat((col-1)*N1 : col*N) omp_reconstruct(y_col, Phi_col, K); end %% 逆变换得到重构图像 img_hat T * reshape(coeff_hat, N, N) * T;主程序的这段逻辑把一个 N² 维的大规模重构问题拆成了 N 个 N 维的小问题是一个很实用的工程简化——直接对整个图像做一次 OMP观测矩阵是 M×N²动辄上千万的元素内存直接爆掉拆成按列处理后每列的观测矩阵只有 (M/N)×N小了几个数量级。3.2 OMP 重构函数每一次迭代都在做什么OMP 的核心思想是「贪心」每一轮迭代从观测矩阵的所有列原子中选出与当前残差最相关的一个加入支撑集然后用最小二乘法更新系数再算新残差。重复 K 次后就能用 K 个原子的线性组合近似表示原始信号。下面这段代码是整套资源里最值得反复读的部分。function s_hat omp_reconstruct(y, Phi, K) % OMP正交匹配追踪重构 % 输入: y - 观测向量(M×1), Phi - 观测矩阵(M×N), K - 稀疏度 % 输出: s_hat - 重构的稀疏系数(N×1) [M, N] size(Phi); s_hat zeros(N, 1); residual y; % 残差初始化为观测向量 support []; % 支撑集, 记录选中原子的列索引 for iter 1:K % 计算残差与每个原子的内积, 找最相关的原子索引 correlations Phi * residual; [~, idx] max(abs(correlations)); support union(support, idx); % 用支撑集对应的原子做最小二乘, 更新系数 Phi_s Phi(:, support); s_temp Phi_s \ y; % 更新残差: 减去当前已解释的部分 residual y - Phi_s * s_temp; % 若残差足够小, 提前终止迭代 if norm(residual) 1e-6 break; end end s_hat zeros(N, 1); s_hat(support) s_temp; end这段代码里有几个关键操作值得展开讲。Phi * residual是在计算残差和每个原子的内积内积绝对值最大的列意味着与残差方向最接近也就是当前最能解释残差的候选原子。max(abs(correlations))取的是绝对值而不是原始值因为负相关同样能解释信号。union(support, idx)把新选中的列索引加入支撑集union自带去重功能防止同一个原子被重复选中。最小二乘求解用的是左除运算符\MATLAB 会自己选择合适的算法对于小规模矩阵通常走 QR 分解数值稳定性比显式求逆好得多。3.3 为什么用左除而不是inv求逆很多第一次接触这段代码的人会把s_temp Phi_s \ y改写成s_temp inv(Phi_s * Phi_s) * Phi_s * y结果经常得到警告甚至 NaN。原因是最小二乘的正规方程(Phi_s * Phi_s) \ Phi_s * y涉及矩阵求逆当 Phi_s 的列之间存在近似线性相关时Phi_s * Phi_s接近奇异逆矩阵数值不稳定。MATLAB 的\运算符直接对 Phi_s 做正交分解后求解绕开了正规方程的病态问题。这是这份代码里很值得保留的写法。另一个容易忽略的细节是提前终止条件norm(residual) 1e-6。理论上迭代 K 次就会结束但实际中当残差已经小到一定程度说明信号已经被充分解释继续迭代只会引入噪声原子的干扰。这个阈值 1e-6 是我在实际调参中觉得比较合适的值——太小会浪费算力太大会截断有效信息。3.4 从系数向量到重构图像的逆变换OMP 恢复出来的是 DCT 域系数真正要得到观众能看的图像还得做逆变换。T * reshape(coeff_hat, N, N) * T这行代码的原理是因为 DCT 变换矩阵 T 是正交矩阵T T⁻¹所以逆变换就是左乘 T 再右乘 T。这里有一个工程判断点如果直接用dctmtx生成的矩阵做逆变换矩阵乘法每一步都是 O(N³) 的复杂度256×256 的图像还扛得住但 512×512 就开始明显变慢。更高效的做法是用 MATLAB 内建的dct2和idct2函数它们内部用快速算法速度可以快 10 倍以上。我在自己的版本里会保留dctmtx写法作为教学参考但实际跑大数据时会切换成idct2。4. 把实验跑起来从main.m到PSNR/SSIM指标的完整流程4.1 运行环境与第一步操作这份代码的运行环境要求非常常规Windows、MacOS 或 Linux 上的 MATLAB R2018b 及以上版本都行不需要额外的工具箱——randn、dctmtx、rgb2gray全是基础函数。拿到工程文件后先把整个文件夹加到 MATLAB 路径右键文件夹 → 添加到路径 → 选中的文件夹和子文件夹这样最省事。然后在命令行窗口敲main程序就会自动跑完整个链路并在命令行打印观测数、稀疏度、PSNR 和运行时间。% 命令行执行 main 压缩感知图像重构实验图像尺寸: 256 × 256 观测数量: 19661 (测量率: 0.30) 稀疏度: K 4915 PSNR 28.74 dB, SSIM 0.8712 运行时间: 8.32 秒4.2 实验完成后如何检查结果是否正确程序运行完会在 Figure 窗口画出两张图左侧是原始图像右侧是重构图像。判断重构是否成功的标准首先是肉眼观察——边缘是否清晰、纹理是否连续、整体有没有块状伪影。其次是看 PSNR 和 SSIM 两个数值自然图像的压缩感知重构PSNR 在 25dB 以上算合格30dB 以上算优秀SSIM 最好在 0.85 以上。如果你跑出来的 PSNR 只有 15dB 以下大概率不是算法问题而是参数配置问题按第五章的思路排查。除了看指标建议手动检查一下中间变量这样能更快理解算法到底在哪一步起作用。在 main.m 的y Phi * coeff_vec;之后加一行disp(size(y))确认观测向量的维度是 M×1在coeff_hat ...之后加disp(norm(coeff_hat - coeff_vec))这个差值越接近零说明重构系数和原始系数的误差越小。这个差值对应到图像上就是重构质量的根本来源。4.3 换一张你自己的图像跑实验代码默认读入的是一张名为cameraman.tif的测试图想换成你自己的图像改动一行就行img imread(lena.png); % 换成你自己的图像改完要确认两点一是图像最好是方形如果不是程序会报矩阵维度不匹配的错误解决办法是改成img imresize(img, [256, 256]);做一次缩放二是图像最好是灰度图如果是彩色图程序会自动转灰度但要注意rgb2gray得到的灰度图本身是 uint8 类型务必用im2double转成 double 类型否则后面矩阵运算会直接报错。我自己第一次跑的时候就是忘了im2doubleMATLAB 报了Undefined function mtimes for input arguments of type uint8这个错误特别容易在新手期翻车。5. 避坑记录二维图像重构中常见的5个问题与排查思路5.1 重构图像出现横竖条纹伪影现象重构出的图像整体轮廓能看清但表面覆盖着规则的横纹或竖纹像是印刷品上的扫描线。原因观测矩阵的分块方式和 DCT 系数的分块方式不一致。主程序按列拆分观测时如果观测矩阵的行分块边界和 DCT 系数的列分块边界没有对齐不同块之间的观测信息就会互相串扰。解决检查y和coeff_vec的索引划分是否严格对应——y 的第 (col-1)M/N1 到 colM/N 行必须对应 coeff_vec 的同一个列区间。我一般会在拆分逻辑处打一个断点打印y_col的长度是否等于M/N以及Phi_col的列数是否等于 N确认之后再往下跑。5.2 PSNR 低到离谱但不报错现象程序全程无报错但 PSNR 只有 10~15dB重构图像完全是一团噪点。原因稀疏度 K 设置得过小或过大。K 太小大量中小幅度的 DCT 系数被丢弃图像细节全丢K 太大OMP 把噪声也当成有效信号恢复出来了重构结果里高频噪声反而更多。解决固定其他参数写一个for K M/8 : M/2的循环打印每个 K 对应的 PSNR画出一条 PSNR-K 曲线。曲线通常是倒 U 型取峰值对应的 K 值即可。这个操作值得自动化成一个脚本因为换一张图最优 K 就变一次。5.3 观测矩阵不归一化导致的重构质量下降现象换了新的观测矩阵比如改成伯努利随机矩阵之后PSNR 莫名其妙掉了 2~4dB。原因伯努利矩阵的元素是 ±1列能量是 N列长度而高斯矩阵的每列能量大约是 1因为 randn 的方差为 1。如果沿用高斯矩阵时代的支撑集选择逻辑不重新归一化OMP 选原子时的相关性计算就会被列能量带偏。解决无论换成什么观测矩阵都强制做列归一化。这是我在改代码时踩过最隐蔽的一个坑——表面看起来只是数值范围不同实际影响的是整个支撑集的选择路径。5.4 大图像直接内存溢出现象把图像从 256×256 换成 512×512MATLAB 直接报Out of memory或者卡死。原因512×512 意味着信号维度 N²262144稀疏变换矩阵 T 是 262144×262144光存储这个矩阵就需要几百 GB完全不可行。解决要么保持按列分块的重构策略不变要么改用滑窗方式把图像拆成多个 8×8 或 16×16 的小块每个小块独立做稀疏变换和重构再拼接回去。滑窗块的大小选择是一个权衡——块越大稀疏性越好重构质量越高块越小运行越快但块边缘会出现拼缝伪影。16×16 是我觉得比较平衡的选择。5.5 DCT 系数分布不均导致低频块恢复好、高频块恢复差现象重构图像整体亮度没问题但人物的眼睛、发丝等纹理丰富区域明显发糊。原因自然图像的 DCT 系数集中在低频高频区域对应的系数稀疏度不够——如果你用统一的一个 K 值去做全图重构高频区域被截断的信息比例大细节自然丢失。解决按块的自适应稀疏度策略。计算每个分块 DCT 系数的能量占比能量集中度高的块用较大的 K能量分散的块用较小的 K效果立竿见影。这份代码里没有实现这个自适应策略它留作二次开发的空间正好合适。6. 二次开发方向把“能跑”变成“能改、能上台面”6.1 替换观测矩阵从高斯到伯努利或部分哈达玛高斯随机矩阵是理论和工程上的默认选择但工程落地时它有两个不便元素是浮点数内存占用大乘法计算量也大。换成伯努利随机矩阵元素为 ±1存储可以用位运算乘法变成加减法在嵌入式环境里能省下可观的资源。替换方法很简单——把gaussian_measurement函数里的randn(M, N)改成rand(M, N) 0.5再映射到 ±1同时保留列归一化。改完跑同一张图对比 PSNR你会发现质量几乎不变这就是伯努利矩阵满足 RIP 的实际体现。6.2 替换重构算法从 OMP 到 ISTAOMP 是贪心算法速度快但理论保证弱信号不严格稀疏时重构质量不稳定。ISTA迭代收缩阈值算法是另一种主流思路每次迭代做一次梯度下降再用软阈值算子把小于阈值的系数置零逐步逼近稀疏解。替换的收益是在相同观测数下ISTA 对非严格稀疏信号比如自然图像的重构质量通常比 OMP 高 1~2dB代价是迭代次数多速度慢一个量级。我建议把 OMP 函数和 ISTA 函数都保留封装成统一的接口reconstruct_solver(y, Phi, K, method)这样切换算法只需要改一个字符串参数。6.3 把灰度图扩展成彩色图像RGB 彩色图像可以当成三个独立的灰度图通道每个通道单独做压缩感知重构再合成彩色图。这里有一个工程技巧三个通道的观测矩阵要相同否则合成时会引入通道间的错位噪声。实现上先生成一次 Phi然后对 R、G、B 三个通道复用同一矩阵。颜色空间上RGB 每个通道的相关性高重构质量相对均衡如果追求更高的压缩率可以先把 RGB 转换到 YCbCr 色彩空间对亮度通道 Y 用较高的测量率人眼对亮度敏感对色度通道 Cb、Cr 用较低的测量率总数据量下降但主观质量几乎不受影响。我自己第一次完整跑这套代码的时候就是卡在 5.3 那个观测矩阵归一化的问题上查了一下午资料最后打印出每列的能量才发现是归一化丢了。从那以后我每次改观测矩阵都会强制走一遍「列能量检查归一化固定种子对比 PSNR」的流程再也不敢凭感觉跳过。这个习惯帮我避开了后面好几次类似的坑。这套代码的特点就是结构清晰、每个模块都能独立修改和验证希望你能在这个基础上改出属于自己的压缩感知重构版本——希望帮到你。本文还有配套的精品资源点击获取

二维图像压缩感知重构：MATLAB完整实现与OMP算法解析

相关推荐

中医方剂知识图谱构建：SpringBoot与Neo4j实践

multipart/form-data上传图片原理与实战指南

OpenClaw智能体框架核心原理与工业级选型指南

更多资讯

自动控制原理工程实战：从建模到校正的Python实现与避坑指南

Sentinel流量治理实战：微服务熔断限流与规则持久化

Prim算法详解：最小生成树的贪心策略与实现

爱立信6601开站四阶校验闭环与5大翻车根因解析

ITK容器详解：从Array到Image的数据表示与内存管理

openrig：打造可复现开发环境的声明式配置方案

今日更新

防关联浏览器哪个好？三款产品实测差异解析

2026钢铁行业AI转型系统选型对比：多基地排产、大宗原燃料采购与质量追溯，六家厂商能力差异在哪里

影刀RPA深度实战：大模型Fine-tuning定制RPA专用Prompt与指令生成

本周精选

交通拥堵预测毕设：基于Python+LSTM的车辆流量预测系统全流程解析

MATLAB/Simulink搭建10机39节点电力系统暂态稳定仿真实战

无人机航拍三维重建全流程：从SfM到网格生成的避坑指南

本月热门

我发现了一个新思路：用 Remotion + Claude Code 像写代码一样自动化生成短视频

Windows下 Codex 中 Chrome 和 Computer Use 插件不可用问题排查及解决参考方式：TaoToken 统一 Key 配置与验证

2026 大模型集体涨价：用 Python 做企业 Token 成本测算与选型避坑（附配置）

读完文章，想聊聊您的网站？