[SciPy-user] 3D density calculation
Chris Lee
c.j.lee@tnw.utwente...
Mon Jun 25 09:12:09 CDT 2007
Thank you both for your reply. Sorry for my lack of response I was at a
conference and the webmail interface only gave me read only access :(
histogramdd looks like it might do it. I had thought about consecutive
1D histograms but realised that it would not work.
I had thought that I couldn't use a fixed grid because I require a large
total volume and also to be able to resolve small distances. However,
the David's idea looks pretty good so I will try that as well.
Thank you both for your help.
Cheers
Chris
David Huard wrote:
> 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
> <mailto: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
> <mailto: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 <mailto:SciPy-user@scipy.org>
> http://projects.scipy.org/mailman/listinfo/scipy-user
>
>
>
>
> _______________________________________________
> SciPy-user mailing list
> SciPy-user@scipy.org <mailto: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
>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: c.j.lee.vcf
Type: text/x-vcard
Size: 174 bytes
Desc: not available
Url : http://projects.scipy.org/pipermail/scipy-user/attachments/20070625/98a53efa/attachment.vcf
More information about the SciPy-user
mailing list