Maximum intensity projection from image stack

Question

I'm trying to recreate the function

max(array, [], 3)

From MatLab, which can take my 300x300px image stack of N images (I'm saying "Image" here because I'm processing images, really this is just a big double array), 300x300xN, and create a 300x300 array. What I think is happening in this function, if it were to operate inefficiently, is that it is parsing through each (x,y) point, then taking the maximum value from that point across the z-axis, then normalizing with maximum and minimum values of the entire array.

I've tried recreating this in python with

# Shape of dataset: (300, 300, 181)
# Type of dataset: <type 'numpy.ndarray'>
for x in range(numpy.size(self.dataset, 0)):
    for y in range(numpy.size(self.dataset, 1)):
        print "Point is", x, y
        # more would go here to find the maximum (x,y) value over Z axis in self.dataset

A very simple X,Y iterator. -- but not only does my IDE crash after a few milliseconds of running this code, but also it feels gross and inefficient.

Is there something I'm missing? I'm new to Python, and therefore the answer here isn't clear to me. Is there an existing function that does this operation?

Answer 1

import numpy as np
import matplotlib.pyplot as plt
from skimage import io

path = "test.tif"
IM = io.imread(path)
IM_MAX= np.max(IM, axis=0)
plt.imshow(IM_MAX)