Dealing with lot of images and multiplications

Question

With some basic knowledge of Python and referring a lot of sources, I have written the code below. But it takes half an hour for execution. How can I reduce the time? I read about vectorization but not understanding how exactly I can use it here.

In this, I have to read 2D skeleton images( Size 1980x1080) and Depth images (size 512x424). As I have to do mapping of both the images of different sizes and form a 3D Skelton, depth_to_xyz_and_rgb(uu, vv, dep) function does that. Here the skeleton Image is generated from the OpenPose Library of Facebook. The skeleton image mainly has 25 key-points and those 25 key-points are joined together to form a 2D skeleton. The main task is to calculate Gait Parameters which could be found out if we know the 3D co-ordinate of each joint. Initially, I was generating 3D Point cloud to view if the code generates it properly or not and as the code was generating it properly and I didn't need the point cloud file for further processing, I have removed that part of code.

Now here I am mainly saving the 2D skeleton location of 25 key-points of a skeleton and 3D values of 25 key-points. I have added a sample of generated excel file as a sample in the link below.

import math
import sys
from PIL import Image
import numpy as np
import scipy.io as sio
import os



scalingFactor = 5000.0
fx_d=365.3768
fy_d=365.3768
cx_d=253.6238
cy_d=211.5918
fx_rgb=1054.8082
fy_rgb=1054.8082
cx_rgb=965.6725
cy_rgb=552.0879



RR = np.array([
    [0.99991, -0.013167,-0.0020807],
    [0.013164,0.99991,-0.0011972],
    [-0.0020963,0.0011697,1]
])
TT = np.array([ 0.052428,0.0006748,0.000098668 ])
extrinsics=np.array([[.99991,-0.013167,-0.0020807,0.052428],[0.013164,0.99991,-0.0011972,0.0006748],[-0.0020963,0.0011697,1,0.000098668],[0,0,0,1]])


path = 'G:\\SENDA\\Proband_172\\GAIT\\RG1_color\\mat\\'
file_lists = os.listdir(path)
path3 = 'G:\\SENDA\\Proband_123\\GAIT\\RG1_Depth\\selected\\'
included_extensions = ['bmp']
file_lists3 = [fn for fn in os.listdir(path3)
              if any(fn.endswith(ext) for ext in included_extensions)]

path4 = 'G:\\SENDA\\Proband_123\\GAIT\\RG1_results\\skeleton\\'
file_lists4 = os.listdir(path4) 


path6 = 'G:\\SENDA\\Proband_123\\GAIT\\RG1_results\\'
file_lists6 = os.listdir(path6)


def init_maxvalue():
    for center in centers3:
        if center==(0,0):
            continue
        min_xy.append(sys.float_info.max)
        min_vex.append((0,0))
        p_vex.append((0,0,0))

def depth_rgb_registration(rgb,depth):

#    f=open("RecordAll.xls",'w') # Taking second argument i.e the depth image name
    init_maxvalue()
    rgb = Image.open(rgb)
    depth = Image.open(depth).convert('L') # convert image to monochrome 
    if rgb.mode != "RGB":
            raise Exception("Color image is not in RGB format")


    for v in range(depth.size[0]):
        for u in range(depth.size[1]):
            try:
                (p,x,y)=depth_to_xyz_and_rgb(v,u,depth) # this gives p = [pcx, pcy,pcz] 
                #aligned(:,:,0) = p
            except:
                continue

            if (x > rgb.size[0]-1 or y > rgb.size[1]-1 or x < 1 or y < 1 or np.isnan(x) or np.isnan(y)):
                continue
            x = round(x)
            y = round(y)
            color=rgb.getpixel((x,y))
            #print(color)
            min_distance((x,y),p)

            if color==(0,0,0):
                p[0]=0
                p[1]=0
                p[2]=0
                continue
            #if.write("%f .%f %f \n"%(p[0],p[1],p[2]))
            points.append(" %f %f %f %d %d %d 0\n"%(p[0],p[1],p[2],255,0,0))
        i=0
        x=[]
        y=[]
        z=[]

    for val in min_vex:
        f.write(str(val)+' '+str(p_vex[i])+'');

        points.append(" %f %f %f %d %d %d 0\n"%(p_vex[i][0],p_vex[i][1],p_vex[i][2],0,255,0))
        x.append(p_vex[i][0])
        y.append(p_vex[i][1])
        z.append(p_vex[i][2])
        i=i+1
    else:
        f.write("\n")
#    f.close()

def min_distance(val,p):
    i=0
    for center in centers3:
        if center==(0,0):
            continue
        temp=math.sqrt(math.pow(center[0]-val[0],2)+math.pow(center[1]-val[1],2))
        if temp<min_xy[i]:
            min_xy[i]=temp
            min_vex[i]=val
            p_vex[i]=p
        i=i+1

def depth_to_xyz_and_rgb(uu , vv,dep):

    # get z value in meters
    pcz =dep.getpixel((uu, vv))
    if pcz==60:
        return

    pcx = (uu - cx_d) * pcz / fx_d
    pcy = (vv - cy_d) * pcz / fy_d

    # apply extrinsic calibration
    P3D = np.array( [pcx , pcy , pcz] )
    P3Dp = np.dot(RR , P3D) - TT

    # rgb indexes that P3D should match
    uup = P3Dp[0] * fx_rgb / P3Dp[2] + cx_rgb
    vvp = P3Dp[1] * fy_rgb / P3Dp[2] + cy_rgb

    # return a point in the point cloud and its corresponding color indices
    return P3D , uup , vvp

if __name__ == '__main__':


    f=open("Proband_123_RG1.xls",'w')    
    for idx, list1 in enumerate(file_lists4):
    # Iterate through items of list2
        for i in range(len(file_lists3)):

            if list1.split('.')[0] == file_lists3[i].split('.')[0]:


                rgb = os.path.join(path4, list1)
                depth = os.path.join(path3,sorted(file_lists3)[i])
                m = sorted(file_lists)[idx]    
                mat2 = sio.loadmat(os.path.join(path,m))
                abc= list1.split('.')[0]
                centers2 = mat2['b2']
                centers3 =np.array(centers2).tolist()
                min_xy=[]
                min_vex=[]
                p_vex=[]
                image_list = []
                image_list2 = []
                points = []
                depth_rgb_registration(rgb,depth)
    f.close

The below part of the code takes maximum time for execution because it does Rotation and Translation for almost every pixel considering first depth image. In my images, RGB images have a lot of zeros, so if this process can be reversed then it could reduce the time I think. I just got an idea and I am working on it. Link for algorithm

def depth_to_xyz_and_rgb(uu , vv,dep):

    # get z value in meters
    pcz =dep.getpixel((uu, vv))
    if pcz==60:
        return

    pcx = (uu - cx_d) * pcz / fx_d
    pcy = (vv - cy_d) * pcz / fy_d

    # apply extrinsic calibration
    P3D = np.array( [pcx , pcy , pcz] )
    P3Dp = np.dot(RR , P3D) - TT

    # rgb indexes that P3D should match
    uup = P3Dp[0] * fx_rgb / P3Dp[2] + cx_rgb
    vvp = P3Dp[1] * fy_rgb / P3Dp[2] + cy_rgb

    # return a point in the point cloud and its corresponding color indices
    return P3D , uup , vvp

• Link for images and mat files Images and Mat files

Answer 1

There are several things that make it is needlessly difficult to figure out what this code is supposed to do.

• Using lots of global variables (f, min_xy, min_vex, etc.) for passing data in to and out of functions makes it difficult to see when a variable's value might be set or changed.
• Many variables with similar, non-descriptive names (list1, list3, center2, center3, etc.) make it more difficult to figure out what the variables purpose is.
• No docstrings describing what a function does or how it is used.

Presuming your system isn't too old, it's likely that the processor has multiple cores. So one way to speed things up would be to use the multiprocessing library. First, create a list of (rgb_file, depth_file) pairs. Then use a process pool to process the list.

Something like:

from multiprocessing import Pool

from collections import default dict
from pathlib import Path


def process_pair(rgb, depth):
    # code to process one image goes here
    ...


if __name__ == '__main__':

    base = Path('G:/SENDA/Proband_123/GAIT')
    skeleton_path = base / 'RG1_results/skeleton'
    image_path = base / 'RG1_Depth/selected'

    rgb_depth_pairs = defaultdict(list)

    for skeleton in skeleton_path.iterdir():
        rgb_depth_pairs[skeleton.stem].append(skeleton)

    included_extensions = ['bmp']
    images = [fn for ext in included_extensions for fn in image_path.glob(f'*.{ext}')]

    for image in images:
        rgb_depth_pairs[image.stem].append(image)

    rgb_depth_pairs = [item for item in rgb_depth_pairs.items() if len(item)==2]


    with Pool() as p:
        p.starmap_async(process_pairs, rgb_depth_pairs)

(Note: This code has not been tested or otherwise debugged)