To show how to use numba and blz we are going to downsample big images, the images we are using here are unconfortable to work with due to their big size on memory.
import numba
import blz
import numpy as np
from pylab import imshow
from scipy import misc
from time import time
from math import sqrt
from shutil import rmtree
For this tutorial I am going to use a chicken image, if you happen to have more than 4GiB of RAM I would recommend you to use The Garden Of Earthly Delights which you can find here.
%matplotlib inline
img = misc.imread('images/chicken.jpg')
imshow(img)
img.shape
(2000, 3008, 3)
The purpose of blz is to work with big data without the need of loading it into memory, so we are doing it wrong. Let's first see how can we convert the image to blz.
rmtree('images/chicken.blz', ignore_errors=True)
blz.barray(
img,
rootdir='images/chicken.blz',
bparams=blz.bparams(clevel=9, shuffle=True, cname='zlib')
)
barray((2000, 3008, 3), uint8) nbytes: 17.21 MB; cbytes: 15.90 MB; ratio: 1.08 bparams := bparams(clevel=9, shuffle=True, cname=zlib) rootdir := 'images/chicken.blz' [[[109 108 104] [112 111 107] [117 114 109] ..., [140 175 151] [136 171 147] [133 168 144]] [[104 103 98] [110 107 102] [114 111 106] ..., [138 172 147] [137 171 146] [138 172 147]] [[104 101 96] [111 106 102] [116 111 107] ..., [144 178 153] [143 175 151] [142 174 150]] ..., [[133 127 129] [138 132 134] [140 134 136] ..., [221 191 167] [224 194 170] [227 197 173]] [[129 125 126] [134 128 130] [135 129 131] ..., [224 192 167] [226 195 167] [228 197 169]] [[130 126 127] [134 130 131] [135 129 131] ..., [226 192 164] [225 191 163] [224 190 162]]]
There are some interesting parameters when we create it:
rootdir: this will tell blz where you want to store the blz, if this option is not given the blz will be created on memory so be careful with that.
bparams: this will tell blz how do you want to store the blz, let's go over the options:
Note that I have chosen zlib because it behaves better than blosclz when working with images.
blz is an amazing container that works exactly as a numpy array, let's see some examples:
#We open the image, this is a lazy operation so nothing is loaded into memory
img = blz.barray(rootdir='images/chicken.blz')
#Shape
img.shape
(2000, 3008, 3)
#First row with all RGB values
img[:, 0, :]
array([[109, 108, 104], [104, 103, 98], [104, 101, 96], ..., [133, 127, 129], [129, 125, 126], [130, 126, 127]], dtype=uint8)
#R values of the first row
img[:, 0, 0]
array([109, 104, 104, ..., 133, 129, 130], dtype=uint8)
#B values of the first two rows
img[:, :2, 2]
array([[104, 107], [ 98, 102], [ 96, 102], ..., [129, 134], [126, 130], [127, 131]], dtype=uint8)
In a barray, information is stored as rows, so it is better to process the image row by row than column by column. (Way better)
Now let's see how can we downsample the image without fully loading it into memory.
The actual downsampling function takes little pieces of the image and downsamples them one by one, therefore you don't need to load the full image, just the downsampled one. It takes an optional parameter to limit the cache used.
def downsample(orig, down_cell, cache_size=2**21):
c0, c1 = down_cell
#Let's calculate the matrix dimensions
pixel_size = orig[0, 0].nbytes
n = int(round(sqrt(cache_size/pixel_size), 0))
#How many complete matrices?
hor = int(orig.shape[1]) / n
ver = int(orig.shape[0]) / n
#Complete matrix dimensions
submatrix_n = round(n/float(c0))
submatrix_center_shape = (submatrix_n, submatrix_n, orig.shape[2])
submatrix_center = np.empty(submatrix_center_shape, dtype=orig.dtype)
#Bottom border matrix dimensions
ver_px = round((int(orig.shape[0]) % n) / c0, 0)
submatrix_bottom_shape = (ver_px, submatrix_n, orig.shape[2])
submatrix_bottom = np.empty(submatrix_bottom_shape, dtype=orig.dtype)
#Right border matrix dimensions
hor_px = round((int(orig.shape[1]) % n) / c1, 0)
submatrix_right_shape = (submatrix_n, hor_px, orig.shape[2])
submatrix_right = np.empty(submatrix_right_shape, dtype=orig.dtype)
#Corner matrix dimensions
submatrix_corner_shape = (ver_px, hor_px, orig.shape[2])
submatrix_corner = np.empty(submatrix_corner_shape, dtype=orig.dtype)
#We build the final container
final_shape = (submatrix_n * ver + ver_px,
submatrix_n * hor + hor_px,
orig.shape[2])
final = np.empty(final_shape, orig.dtype)
#Downsample the middle of the image
for i in xrange(ver):
for j in xrange(hor):
submatrix = orig[i*n:(i+1)*n, j*n:(j+1)*n, :]
mymean(submatrix, c0, c1, submatrix_center)
final[i*submatrix_n:(i+1)*submatrix_n,
j*submatrix_n:(j+1)*submatrix_n, :] = submatrix_center
#Downsample the right border
for i in range(ver):
submatrix = orig[i*n:(i+1)*n, hor*n:, :]
mymean(submatrix, c0, c1, submatrix_right)
final[i * submatrix_n:(i+1)*submatrix_n,
submatrix_n*hor:, :] = submatrix_right
#Downsample the bottom border
for j in range(hor):
submatrix = orig[ver*n:, j*n:(j+1)*n, :]
mymean(submatrix, c0, c1, submatrix_bottom)
final[submatrix_n*ver:,
j * submatrix_n:(j+1)*submatrix_n, :] = submatrix_bottom
#Downsample the corner
submatrix = orig[n*ver:, n*hor:, :]
mymean(submatrix, c0, c1, submatrix_corner)
final[submatrix_n*ver:, submatrix_n*hor:, :] = submatrix_corner
return final
Note that that code doesn't know anything about blz containers, as told before a blz container works pretty much as a numpy array. I call that magic.
But let's see how well does the code performs without numba. (HINT: Not good)
def mymean(src, p0, p1, y):
factor = 1/(1. * p0 * p1)
for z in range(y.shape[2]): # RGB
for i in range(y.shape[0]):
for j in range(y.shape[1]):
s = 0.
for k in range(p0):
for l in range(p1):
s += src[(p0*i)+k, (p1*j)+l, z] * factor
y[i, j, z] = s
t1 = time()
img_down = downsample(img, (16,16))
t2 = time()
elapsed1 = t2-t1
print elapsed1
imshow(img_down)
151.338282108
<matplotlib.image.AxesImage at 0x4d6bb90>
Now let's take the exact same function but with the numba decorator.
@numba.njit
def mymean(src, p0, p1, y):
factor = 1/(1. * p0 * p1)
for z in range(y.shape[2]): # RGB
for i in range(y.shape[0]):
for j in range(y.shape[1]):
s = 0.
for k in range(p0):
for l in range(p1):
s += src[(p0*i)+k, (p1*j)+l, z] * factor
y[i, j, z] = s
t1 = time()
img_down = downsample(img, (16,16))
t2 = time()
elapsed2 = t2-t1
elapsed2
print elapsed2
imshow(img_down)
1.05450105667
<matplotlib.image.AxesImage at 0x5b26a50>
Wow! That was fast.
print str(elapsed1/elapsed2) + ' times faster'
143.51648218 times faster
So we have worked with a "large" image, the larger the image the better performance blz has. I encourage you to try this same tutorial with this image although you should have at least 4GiB of RAM to do the first conversion from jpg to blz, after that you are good to go with any computer.