Element-wise matrix multiplication in NumPy

Question

I'm making my first real foray into Python and NumPy to do some image processing. I have an image loaded as a 3 dimensional NumPy Array, where axis 0 represents image bands, while axes 1 and 2 represent columns and rows of pixels. From this, I need to take the 3x1 matrix representing each pixel and perform a few operations which result in another 3x1 matrix, which will be used to build a results image.

My first approach (simplified and with random data) looks like this:

import numpy as np
import random

factor = np.random.rand(3,3)
input = np.random.rand(3,100,100)
results = np.zeros((3,100,100))

for x in range(100):
    for y in range(100):
        results[:,x,y] = np.dot(factor,input[:,x,y])

But this strikes me as inelegant and inefficient. Is there a way to do this in an element-wise fasion, e.g.:

results = np.dot(factor,input,ElementWiseOnAxis0)

In trying to find a solution to this problem I came across this question, which is obviously quite similar. However, the author was unable to solve the problem to their satisfaction. I am hoping that either something has changed since 2012, or my problem is sufficiently different from theirs to make it more easily solvable.

Solution

Numpy arrays use element-wise multiplication by default. Check out numpy.einsum, and numpy.tensordot. I think what you're looking for is something like this:

results = np.einsum('ij,jkl->ikl',factor,input)