Hi Chris,

Have you tried numpy.histogramdd ? If its still too slow, I have a fortran implementation on the back burner. I could try to finish it quickly and send you a preliminary version.

Other thought: the kernel density estimator scipy.stats.gaussian_kde

David

2007/6/17, Bernhard Voigt <bernhard.voigt@gmail.com>:
Hi Chris!

you could try a grid of unit cells that cover your phase space (x,y,z,t). Count the number of photons per unit cell of your initial configuration and track photons leaving and entering a particular cell. A dictionary with a tuple of x,y,z,t coordinates obtained from integer division of the x,y,z,t coordinates could serve as keys.

Example for 2-D:

from numpy import *
# phase space in x,y
x = arange(-100,100.1,.1)
y = arange(-100,100.1,.1)
# cell dimension in both dimensions the same
GRID_WIDTH=7.5

# computes the grid key from x,y coordinates
def gridKey(x,y):
    '''return the a tuple of x,y integer divided by GRID_WIDHT'''
    return (int(x // GRID_WIDTH), int(y // GRID_WIDTH))

# setup your grid dictionary
gridLowX, gridHighX = gridKey(min(x), max(x))
gridLowY, gridHighY = gridKey(min(y), max(y))
keys = [(i,j) for i in xrange(gridLowX, gridHighX + 1) \
            for j in xrange(gridLowY, gridHighY + 1)]
grid = dict().fromkeys(keys, 0)

# random photons
photons = random.uniform(-100.,100., (100000,2))

# count photons in each grid cell
for p in photons:
    grid[gridKey(*p)] += 1

#########################################
# in your simulation you have to keep track of where your photons
# are going to...
# (the code below won't run, it's just an example)
#########################################
oldKey = gridKey(photon)
propagate(photon)  # changes x,y coordinates of photon
newKey = gridKey(photon)
if oldKey != newKey:
    grid[oldKey] -= 1
    grid[newKey] += 1

I hope this helps! Bernhard


On 6/15/07, Chris Lee < c.j.lee@tnw.utwente.nl> wrote:
Hi everyone,

I was hoping this list could point me in the direction of a more
efficient solution to a problem I have.

I have 4 vectors: x, y, z, and t that are about 1 million in length
that describe the positions of photons.  As my simulation progresses
it updates the positions so x, y, z, and t change by an unknown (and
unknowable) amount every update.

This worked very well for its original purpose but now I need to
calculate the photon density change over time.  Currently after each
update, I iterate over time slices, x slices, and y slices and then
make an histogram of z which I then stitch together to create a
density.  However, this becomes very slow as the photons spread out
in space and time.

Does anyone know how to take such a large vector set and return a
density efficiently?


_______________________________________________
SciPy-user mailing list
SciPy-user@scipy.org
http://projects.scipy.org/mailman/listinfo/scipy-user



_______________________________________________
SciPy-user mailing list
SciPy-user@scipy.org
http://projects.scipy.org/mailman/listinfo/scipy-user