这次给的需求是图片的压缩感知并能够实现图片分类的需求
基本定义
借用之前看过的一个例子
如果我们要感知一个10x10像素的图像,最简单最直接的方式,就是直接在空间域里面感知,得到并储存图像的100个像素点,这样的话,你可以通过把这100个像素点排列为矩阵的方式,重建这副目标图像。这里获取这100个点的过程,就称之为sensing。
压缩感知指的是,当我们没法存储采集到100个像素点,但需要恢复100个像素点的目标图像。
在不借助信号先验的条件下,通过98个线性方程组,解100个未知数,理论上的解空间的自由度会有一个二维平面那么大,而这平面上的每个点都可以是我们恢复的图像,没有唯一解。但如果我们引入了一个信号固有的先验:比如真实的图像的解只能在某条线上。这就相当于我们在98个线性方程上,外加了一个linear constraint。由这个line和我们的解空间自由的平面的交点,就唯一确定了一个解:这个解既满足我们的98个线性方程组,又在我们的linear constraint上。这个例子其实就解了一个极其简化的compressed sensing问题。
工程实现
我们选择实现CoSAMP的算法来实现压缩感知
这样重建以后的图片含有的信息会更加 sensible
论文中的算法如下图所示
image
通过python实现CoSAMP的算法如下
#coding:utf-8
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
# DCT基作为稀疏基,重建算法为CoSaMP算法,图像按列进行处理
# 参考文献: D. Deedell andJ. Tropp, “COSAMP: Iterative Signal Recovery from
#Incomplete and Inaccurate Samples,” 2008.
#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
#导入集成库
import math
# 导入所需的第三方库文件
import numpy as np #对应numpy包
from PIL import Image #对应pillow包
#读取图像,并变成numpy类型的 array
im = np.array(Image.open('fish.jpg'))#图片大小256*256
#生成高斯随机测量矩阵
sampleRate=0.5 #采样率
Phi=np.random.randn(256*sampleRate,256)
# Phi=np.random.randn(256,256)
# u, s, vh = np.linalg.svd(Phi)
# Phi = u[:256*sampleRate,] #将测量矩阵正交化
#生成稀疏基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)
#CoSaMP算法函数
def cs_CoSaMP(y,D):
S=math.floor(y.shape[0]/4) #稀疏度
residual=y #初始化残差
pos_last=np.array([],dtype=np.int64)
result=np.zeros((256))
for j in range(S): #迭代次数
product=np.fabs(np.dot(D.T,residual))
pos_temp=np.argsort(product)
pos_temp=pos_temp[::-1]#反向,得到前面L个大的位置
pos_temp=pos_temp[0:2*S]#对应步骤3
pos=np.union1d(pos_temp,pos_last)
result_temp=np.zeros((256))
result_temp[pos]=np.dot(np.linalg.pinv(D[:,pos]),y)
pos_temp=np.argsort(np.fabs(result_temp))
pos_temp=pos_temp[::-1]#反向,得到前面L个大的位置
result[pos_temp[:S]]=result_temp[pos_temp[:S]]
pos_last=pos_temp
residual=y-np.dot(D,result)
return result
#重建
sparse_rec_1d=np.zeros((256,256)) # 初始化稀疏系数矩阵
Theta_1d=np.dot(Phi,mat_dct_1d) #测量矩阵乘上基矩阵
for i in range(256):
print('rebuilding',i,'th column...')
column_rec=cs_CoSaMP(img_cs_1d[:,i],Theta_1d) #利用CoSaMP算法计算稀疏系数
sparse_rec_1d[:,i]=column_rec;
img_rec=np.dot(mat_dct_1d,sparse_rec_1d) #稀疏系数乘上基矩阵
#显示重建后的图片
image2=Image.fromarray(img_rec)
image2.show()
网友评论