压缩感知重构算法之SP算法python实现
来源:互联网 发布:win10查看电脑mac地址 编辑:程序博客网 时间:2024/05/16 02:01
SP(subspace pursuit)算法是压缩感知中一种非常重要的贪婪算法,它有较快的计算速度和较好的重构概率,在实际中应用较多。本文给出了SP算法的python和matlab代码,以及完整的仿真过程。
参考文献:Dai W, Milenkovic O. Subspace pursuit for compressive sensing signal reconstruction[J]. Information Theory, IEEE Transactions on, 2009, 55(5): 2230-2249.
SP算法流程:
代码
要利用python实现,电脑必须安装以下程序
- python (本文用的python版本为3.5.1)
- numpy python包(本文用的版本为1.10.4)
- scipy python包(本文用的版本为0.17.0)
- pillow python包(本文用的版本为3.1.1)
另外需要下载lena图片放在和程序同一个目录下面
#coding:utf-8#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%# DCT基作为稀疏基,重建算法为SP算法 ,图像按列进行处理# 参考文献: W. Dai and O. Milenkovic, “Subspace Pursuit for Compressive# Sensing Signal Reconstruction,” 2009.#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%#导入集成库import math# 导入所需的第三方库文件import numpy as np #对应numpy包from PIL import Image #对应pillow包#读取图像,并变成numpy类型的 arrayim = np.array(Image.open('lena.bmp'))#图片大小256*256#生成高斯随机测量矩阵sampleRate=0.7 #采样率Phi=np.random.randn(256*sampleRate,256)#生成稀疏基DCT矩阵mat_dct_1d=np.zeros((256,256))v=range(256)for k in range(0,256): dct_1d=np.cos(np.dot(v,k*math.pi/256)) if k>0: dct_1d=dct_1d-np.mean(dct_1d) mat_dct_1d[:,k]=dct_1d/np.linalg.norm(dct_1d)#随机测量img_cs_1d=np.dot(Phi,im)#SP算法函数def cs_sp(y,D): K=math.floor(y.shape[0]/3) pos_last=np.array([],dtype=np.int64) result=np.zeros((256)) product=np.fabs(np.dot(D.T,y)) pos_temp=product.argsort() pos_temp=pos_temp[::-1]#反向,得到前面L个大的位置 pos_current=pos_temp[0:K]#初始化索引集 对应初始化步骤1 residual_current=y-np.dot(D[:,pos_current],np.dot(np.linalg.pinv(D[:,pos_current]),y))#初始化残差 对应初始化步骤2 while True: #迭代次数 product=np.fabs(np.dot(D.T,residual_current)) pos_temp=np.argsort(product) pos_temp=pos_temp[::-1]#反向,得到前面L个大的位置 pos=np.union1d(pos_current,pos_temp[0:K])#对应步骤1 pos_temp=np.argsort(np.fabs(np.dot(np.linalg.pinv(D[:,pos]),y)))#对应步骤2 pos_temp=pos_temp[::-1] pos_last=pos_temp[0:K]#对应步骤3 residual_last=y-np.dot(D[:,pos_last],np.dot(np.linalg.pinv(D[:,pos_last]),y))#更新残差 #对应步骤4 if np.linalg.norm(residual_last)>=np.linalg.norm(residual_current): #对应步骤5 pos_last=pos_current break residual_current=residual_last pos_current=pos_last result[pos_last[0:K]]=np.dot(np.linalg.pinv(D[:,pos_last[0:K]]),y) #对应输出步骤 return result#重建sparse_rec_1d=np.zeros((256,256)) # 初始化稀疏系数矩阵 Theta_1d=np.dot(Phi,mat_dct_1d) #测量矩阵乘上基矩阵for i in range(256): print('正在重建第',i,'列。。。') column_rec=cs_sp(img_cs_1d[:,i],Theta_1d) #利用SP算法计算稀疏系数 sparse_rec_1d[:,i]=column_rec; img_rec=np.dot(mat_dct_1d,sparse_rec_1d) #稀疏系数乘上基矩阵#显示重建后的图片image2=Image.fromarray(img_rec)image2.show()
2 0
- 压缩感知重构算法之SP算法python实现
- 压缩感知重构算法之OMP算法python实现
- 压缩感知重构算法之IHT算法python实现
- 压缩感知重构算法之CoSaMP算法python实现
- 压缩感知重构算法之OLS算法python实现
- 压缩感知重构算法之IRLS算法python实现
- 压缩感知重构算法——SP算法
- 压缩感知重构算法之子空间追踪(SP)
- 压缩感知重构算法之子空间追踪(SP)
- 压缩感知重构算法之子空间追踪(SP)
- 压缩感知——SP(subspace pursuit)重构算法前言翻译
- 压缩感知重构算法之正交匹配追踪(OMP)
- 压缩感知重构算法之正交匹配追踪(OMP)
- 压缩感知重构算法之迭代硬阈值(IHT)
- 压缩感知重构算法之正交匹配追踪(OMP)
- 压缩感知重构算法之迭代软阈值(IST)
- 压缩感知重构算法之正交匹配追踪(OMP)
- 压缩感知重构算法之压缩感知匹配追踪(CoSaMp)
- Android 省市县 三级联动(android-wheel的使用)
- do while
- iOS地图展示mapkit框架-基本使用
- 文章标题
- 22.5-5求有向图的分量图
- 压缩感知重构算法之SP算法python实现
- C++内存管理
- c语言 inline函数的总结
- Java实验:输出希腊字母
- linux 使用epoll实现网络通信
- HDU 3293 sort [Ad Hoc]
- CSS的未来:一些试验性CSS属性
- Xcode的"Could not insert new outlet connection"问题
- toString()方法