I find myself wanting to run through a list of (x,y,z) coordinates of some data volume (here called “d”) to do some sort of processing on each voxel. What I have come up with is the following…

First, find the set of coordinates that match some criterion. For example, find all coordinates in “d” that are greater than the 70^th percentile:
[cc lang="python"]coords = array( nonzero( d > prctile( d, 70 ) ) ).transpose()[/cc]
Now that we have the list of coordinates, we can run through each coordinate and do some sort of processing on it:
[cc lang="python"]
for ii,coord in enumerate( coords ):
r = coord[0]
c = coord[1]

# more stuff here
[/cc]
Obviously if you are using a 3 dimensional volume “d” then you would use:
[cc lang="python"]
s = coord[0]
r = coord[2]
c = coord[2]
[/cc]