How to perform operations on every pixel of an image faster than using two for loops?

Question

I am using the following method to perform operations on every pixel of an image but it is too slow. Roughly takes around 110-120s on my machine.

 for i, j in product(xrange(15, width - 15), xrange(15, height - 15)):
        # finding the avg of 15x15 window
        temp = image.crop((i - 7, j - 7, i + 8, j + 8))
        N = numpy.mean(list(temp.getdata()))

        # calling the functions
        avg0_3,avg0_5,avg0_7,avg0_9,avg0_11,avg0_13,avg0_15 = angle_0(7, 7, temp)
        avg15_3,avg15_5,avg15_7,avg15_9,avg15_11,avg15_13,avg15_15 = angle_15(7, 7, temp)
        avg30_3,avg30_5,avg30_7,avg30_9,avg30_11,avg30_13,avg30_15 = angle_30(7, 7, temp)
        avg45_3,avg45_5,avg45_7,avg45_9,avg45_11,avg45_13,avg45_15 = angle_45(7, 7, temp)
        avg60_3,avg60_5,avg60_7,avg60_9,avg60_11,avg60_13,avg60_15 = angle_60(7, 7, temp)
        avg75_3,avg75_5,avg75_7,avg75_9,avg75_11,avg75_13,avg75_15 = angle_75(7, 7, temp)
        avg90_3,avg90_5,avg90_7,avg90_9,avg90_11,avg90_13,avg90_15 = angle_90(7, 7, temp)

        avg105_3,avg105_5,avg105_7,avg105_9,avg105_11,avg105_13,avg105_15 = angle_105(7, 7, temp)
        avg120_3,avg120_5,avg120_7,avg120_9,avg120_11,avg120_13,avg120_15 = angle_120(7, 7, temp)
        avg135_3,avg135_5,avg135_7,avg135_9,avg135_11,avg135_13,avg135_15 = angle_135(7, 7, temp)
        avg150_3,avg150_5,avg150_7,avg150_9,avg150_11,avg150_13,avg150_15 = angle_150(7, 7, temp)
        avg165_3,avg165_5,avg165_7,avg165_9,avg165_11,avg165_13,avg165_15 = angle_165(7, 7, temp)

        # largest grey level lines (L3,L5,L7,L9,L11,L13,L15)
        L3 = max(avg0_3, avg15_3, avg30_3, avg45_3, avg60_3, avg75_3, avg90_3, avg105_3, avg120_3, avg135_3, avg150_3, avg165_3)
        L5 = max(avg0_5, avg15_5, avg30_5, avg45_5, avg60_5, avg75_5, avg90_5, avg105_5, avg120_5, avg135_5, avg150_5,avg165_5)
        L7 = max(avg0_7, avg15_7, avg30_7, avg45_7, avg60_7, avg75_7, avg90_7, avg105_7, avg120_7, avg135_7, avg150_7,avg165_7)
        L9 = max(avg0_9, avg15_9, avg30_9, avg45_9, avg60_9, avg75_9, avg90_9, avg105_9, avg120_9, avg135_9, avg150_9,avg165_9)
        L11 = max(avg0_11, avg15_11, avg30_11, avg45_11, avg60_11, avg75_11, avg90_11, avg105_11, avg120_11, avg135_11, avg150_11,avg165_11)
        L13 = max(avg0_13, avg15_13, avg30_13, avg45_13, avg60_13, avg75_13, avg90_13, avg105_13, avg120_13, avg135_13, avg150_13,avg165_13)
        L15 = max(avg0_15, avg15_15, avg30_15, avg45_15, avg60_15, avg75_15, avg90_15, avg105_15, avg120_15, avg135_15, avg150_15,avg165_15)

        '''
        # largest grey level orthognal line
        L2 = max(avgorth0, avgorth15, avgorth30, avgorth45, avgorth60, avgorth75, avgorth90, avgorth105, avgorth120,
             avgorth135, avgorth150, avgorth165)
        strength2 = L2 - N
        '''
        # line strengths of lines (L3,L5,L7,L9,L11,L13,L15)

        strength3 = L3 - N
        strength5 = L5 - N
        strength7 = L7 - N
        strength9 = L9 - N
        strength11 = L11 - N
        strength13 = L13 - N
        strength15 = L15 - N


        S3.append(strength3)
        S5.append(strength5)
        S7.append(strength7)
        S9.append(strength9)
        S11.append(strength11)
        S13.append(strength13)
        S15.append(strength15)
        #S2.append(strength2)
        p = image.getpixel((i,j))
        I.append(p)
        R = (strength3 + strength5 + strength7 + strength9 + strength11 + strength13 + strength15 + p) * 0.125
        R_comb.append(R)
        result.putpixel((i, j), R)

def angle_0(i, j, image):


sum3 = image.getpixel(((i - 1), j)) + image.getpixel((i, j)) + image.getpixel(((i + 1), j))
sum5 = image.getpixel(((i - 2), j)) + sum3 + image.getpixel(((i + 2), j))
sum7 = image.getpixel(((i - 3), j)) + sum5 + image.getpixel(((i + 3), j))
sum9 = image.getpixel(((i - 4), j)) + sum7 + image.getpixel(((i + 4), j))
sum11= image.getpixel(((i - 5), j)) + sum9 + image.getpixel(((i + 5), j))
sum13= image.getpixel(((i - 6), j)) + sum11+ image.getpixel(((i + 6), j))
sum15= image.getpixel(((i - 7), j)) + sum13+ image.getpixel(((i + 7), j))

avg_sum3 = sum3 / 3
avg_sum5 = sum5 / 5
avg_sum7 = sum7 / 7
avg_sum9 = sum9 / 9
avg_sum11 = sum11 / 11
avg_sum13 = sum13 / 13
avg_sum15 = sum15 / 15

return avg_sum3,avg_sum5,avg_sum7,avg_sum9,avg_sum11,avg_sum13,avg_sum15

what is more efficient way to do this? Keep in mind, my operations require that I need the coordinates of pixels because I also need pixel values some of the pixels that are at a position relative to i,j

Answer 1

This is a simple solution that you can adapt to your problem. The program loads, processes, and saves the output for a series of files using multithreading. Basically, you would split processing of your images into multiple threads that would, ideally, process batches images in parallel. Here I split processing of some .txt files named Matrix_#.txt.

Depending on the details of your problem this may or may not speed up your code - it may even slow it down, but you won't really know until you try it.

#!/usr/bin/python2.7

import time
import math
import numpy as np
import threading

def process_matrix(pathM, nameM, n0, nf):
    """ Operations applied by each thread
        for it's range of files """
    # pathM is the path to files ending in /
    # or \ depending on the OS. Files have a
    # common name nameM with .txt extension
    # and numbers that go from 0 to nf
    for i in range(n0,nf+1):
        # Load as numpy arrays
        in_name = pathM + nameM + str(i) + '.txt'
        temp = np.loadtxt(in_name)
        # Initialize array with results
        res = np.zeros(3)
        # Send to function that processes the data
        res = processed(temp)
        # Save the output
        out_name =  pathM + 'Out_' + nameM + str(i) + '.txt'
        np.savetxt(out_name, res)

def processed(M):
    """ Function with example operations
        on the data """
    # Add this to simulate processing time
    time.sleep(20)
    # Actual simple operations 
    res = np.zeros(3)
    res[0] = np.mean(M)
    res[1] = np.amax(M)
    res[2] = np.amin(M)
    return res

# Time it - at least to optimize thread number
start = time.time()

# Array with Thread objects
threads = []
# Number of threads
num_threads = 10
# Number of files to process
num_files = 10
# "Increment" for each thread
dnth = math.floor(num_files/num_threads)

# Distribute the processing across threads
for ith in range(0,num_threads):
    # Lower limit for file names for a given thread
    low_lim = int(ith*dnth+1)
    # Upper limit for file names (last thread gets all until 
    # the end)
    up_lim = int(num_files) if (ith == num_threads-1) else int((ith+1)*dnth)
    # Start a thread and append the resulting Thread object to the
    # threads list
    thread = threading.Thread(target=process_matrix, 
            args=('/Users/atru/Research/stack/multi_matrix/', 'Matrix_', 
            low_lim, up_lim))
    threads.append(thread)
    thread.start()

# This allows the program to wait until all
# threads finish the execution (i.e. reach this point)
for t in threads:
    t.join()

# Measure execution time and print it
end = time.time()
print(end - start)

This strategy worked great for me once, the processing was speed up 10 times with 14 threads (the system had 6-8 processors so the number was OK). I double checked the performance just today. This particular program has an almost 10 times speed up with 10 threads when ran on 10 relatively small files on a 2 core, 4 thread processor. Since the program is pausing instead of processing this may not be the best scenario, when running you should also check how it performs with more reasonable number of threads like 2 or 4.

It may also not work as intended, due to various overheads and RAM issues (if your one file is half of your RAM then if two threads open two files your RAM is getting full)- but since your processing takes 2 min it seems similar to my problem, where this was successful.

Let me know if you have questions.