SVD image compression

ith rank-1 image, and ith partial sum

1 :

component_0001.png using_0001.png

2 :

component_0002.png using_0002.png

3 :

component_0003.png using_0003.png

4 :

component_0004.png using_0004.png

5 :

component_0005.png using_0005.png

6 :

component_0006.png using_0006.png

7 :

component_0007.png using_0007.png

8 :

component_0008.png using_0008.png

9 :

component_0009.png using_0009.png

10 :

component_0010.png using_0010.png

11 :

component_0011.png using_0011.png

12 :

component_0012.png using_0012.png

16 :

component_0016.png using_0016.png

20 :

component_0020.png using_0020.png

32 :

component_0032.png using_0032.png

64 :

component_0064.png using_0064.png

125 :

component_0125.png using_0125.png

250 :

component_0250.png using_0250.png 
# compress a photo using SVD

from numpy import *
from scipy.misc import imread, imsave

a = imread('438_538_s16_group_33pc.jpg')
h,w = a.shape
a = array(a,dtype=float)
u,s,vt = linalg.svd(a)

S = zeros((h,w))
fill_diagonal(S,s)
# check USV^T really = A
print 'check\n',( dot(u,dot(S,vt)) - a ).round(4)

aa = zeros_like(a)  # to accumulate partial sum of rank-1 matrices

rst1 = open('svdimages/components.rst','w') # rst to display rank-1s
print >> rst1,'SVD image compression'
print >> rst1,'-------------------------------------------------------'
print >> rst1,'ith rank-1 image, and ith partial sum\n'

# choose some partial sums to write out as images
#ifigs = set(range(1,20)+[ 2**i for i in range(20) if 2**i<min(h,w) ])
ifigs = [1,2,3,4,5,6,7,8,9,10,11,12,16,20,32,64,125,250]

def normalize(a):
        ra = a - min(0,a.min()) # shift up if there are negative values
        ra *= 255/ra.max()
        return array(ra,dtype=uint8)

for i in range(min(h,w)):
        rank1 = outer(s[i]*u[:,i],vt[i,:])
        aa += rank1
        print i,linalg.norm(aa-a) # how far off from A is the partial sum

        if i+1 in ifigs:
                imagename = 'component_'+str(i+1).zfill(4)+'.png'
                imsave('svdimages/'+imagename,normalize(rank1))
                print >> rst1,i+1,':\n'
                print >> rst1,'.. image::',imagename
                print >> rst1,'\t:width: 45%'

                print 'doing figure',i+1
                imsave('svdimages/using_'+str(i+1).zfill(4)+'.png',normalize(aa))
                imagename = 'using_'+str(i+1).zfill(4)+'.png'
                print >> rst1,'.. image::',imagename
                print >> rst1,'\t:width: 45%'
                print >> rst1,'\n'
        print >> rst1,'\n'

#imagename = 'using_all.png'
#print >> rst2,'All:\n'
#print >> rst2,'.. image::',imagename
#imsave('svdimages/'+imagename,a)

# copy this code into the rst
py = open('compress.py').readlines()
print >> rst1,'\n.. code:: python\n'
for line in py:
        print >> rst1,'\t'+line[:-1]


rst1.close()
rst2.close()