CV HW1

main.py

# Jose Sotelo
# 2015
# Universite de Montreal
import matplotlib.pyplot as plt
import numpy as np

# x = np.linspace(0, 3*np.pi, 500)
# plt.plot(x, np.sin(x**2))
# plt.title('A simple chirp');
# plt.show()

DIR="/data/lisa/data/cifar10/cifar-10-batches-py/"

def unpickle(file):
    import cPickle
    fo = open(file, 'rb')
    dict = cPickle.load(fo)
    fo.close()
    return dict

def load_data(dir):
    #
    # dir: Directory with the python batches.
    #
    # train_X: numpy array of dimension (50000, 3072)
    # train_y: list of length 50000
    #
    from glob import glob
    files = glob(dir+"*data*")
    files = sorted(files)

    train_X = np.zeros((0, 3072))
    train_y = []

    test_X = np.zeros((0, 3072))
    test_y = []

    for file in files:
        data = unpickle(file)
        train_X = np.vstack((train_X, data['data']))
        train_y = train_y + data['labels']

    return (train_X, train_y)

raw_X,y = load_data(DIR)

def gray_image(X):
    X = X.reshape(50000,3,32,32).swapaxes(1,3).swapaxes(1,2)/float(255)
    X = (X*np.array([0.299, 0.587, 0.144])).sum(3)
    return X

gray_X = gray_image(raw_X)

categories = unpickle(DIR + "batches.meta")

def plot_image(index, X):
    plt.imshow(X[index,:,:], cmap=plt.get_cmap('gray'))
    plt.show()

plot_image(43, gray_X)

#######################

import matplotlib.pyplot as plt
import numpy as np
import scipy.misc as misc
from time import time

DIR = '/data/lisatmp3/sotelo/Results/'
image = misc.imread(DIR + 'nao_bw.jpg')
    
class Image:
    def get_value(self, x,y):
        if self.isFilter:
            x = x + self.center[0]
            y = y + self.center[1]
        if x < 0 or y < 0:
            return 0
        else:
            try:
                value = self.image[x,y]
            except:
                value = 0
            return value
    def __init__(self, image, isFilter = None):
        self.image = image
        self.isFilter = isFilter
        if isFilter:
            m,n = image.shape
            self.center = (m/2,n/2)
            self.size = np.max((m/2,n/2))
        
def convolve_patch(im1, im2, x,y ):
    result = 0
    for i in np.arange(-im2.size -1,im2.size + 1):
        for j in np.arange(-im2.size -1,im2.size + 1):
            result += im1.get_value(x - i, y - j) * im2.get_value(i,j)
    return result

def convolve(I1, I2):
    im1 = Image(I1)
    m,n = I1.shape
    im2 = Image(I2, isFilter = True)
    result = np.zeros(I1.shape)

    for i in xrange(m):
        for j in xrange(n):
            result[i,j] = convolve_patch(im1, im2, i,j)
    return result

def generate_filter(n, sigma = 1):
    filter_image = np.zeros((n,n))
    x = n/2
    y = n/2
    for i in xrange(n):
        for j in xrange(n):
            filter_image[i,j] = 0.5*np.exp(-0.5*((x-i)**2 + (y-j)**2)/sigma**2)/(np.pi*sigma)
    return filter_image

def compare(n):
    filter_image = generate_filter(n)
    plt.imshow(filter_image, cmap=plt.get_cmap('gray'))
    plt.axis('off')
    plt.savefig('filter_%s.png' % n, bbox_inches='tight')

    start = time()
    convolved_image = convolve(image, filter_image)
    t1 = time()-start

    plt.imshow(convolved_image, cmap=plt.get_cmap('gray'))
    plt.axis('off')
    plt.savefig('convolved_image_%s.png' % n, bbox_inches='tight')

    start = time()
    image_spectrum = np.fft.fft2(image)
    filter_spectrum = np.fft.fft2(filter_image, s = image.shape )
    convolved_spectrum = image_spectrum * filter_spectrum
    convolved_image_fft = np.fft.ifft2(convolved_spectrum).real
    t2 = time()-start

    plt.imshow(convolved_image_fft, cmap=plt.get_cmap('gray'))
    plt.axis('off')
    plt.savefig('convolved_image_fft_%s.png' % n, bbox_inches='tight')

    return t1, t2

times = [compare(n) for n in [3,5,11,15,21,25,31,35,41]]
t1, t2 = zip(*times)

t1 = np.array(t1)
t2 = np.array(t2)

plt.clf()
plt.plot(np.log(1+t1), '-sk')
plt.plot(np.log(1+t2), '-^k')
plt.title("Log runtime with respect to size of filter.")
plt.legend(['My implementation', 'Using FFT'], loc='upper left')
plt.savefig('runtime.png')