hyperion-rt / hyperion

Hyperion Radiative Transfer Code
http://www.hyperion-rt.org
BSD 2-Clause "Simplified" License
52 stars 26 forks source link

Sampling #134

Closed bluescarni closed 9 years ago

bluescarni commented 9 years ago

This is the initial implementation of the sampling code. Example usage:

g = vh.voronoi_grid(sites_arr, np.array([[-1, 1.], [-1, 1], [-1, 1]]),verbose=True,n_samples=1000000)

This will sample at least 1000000 points in the domain. The number of samples for each Voronoi cell is at least 10 and otherwise proportional to the cell volume.

When constructed with a nonzero n_samples parameter the voronoi_grid object will expose two properties: samples and samples_idx. The first one is the list of coordinates of the sampling points, the second one the indices (in sparse format) of the sampling points for each cell.

Example:

>>> g.samples
array([[ 0.95339842,  0.79318593,  0.24375318],
       [ 0.95099843,  0.77655176,  0.23685462],
       [ 0.95147365,  0.78657782,  0.24113526],
       ..., 
       [ 0.50510705,  0.162278  ,  0.45528691],
       [ 0.5166702 ,  0.16003916,  0.46033259],
       [ 0.52101828,  0.14597512,  0.45877929]])

>>> g.samples_idx
array([      0,      10,      20, ..., 1852917, 1852927, 1852937], dtype=int32)

This means that, for the first cell, the sampling points are from index 0 to 10 in the first array, for the second cell from index 10 to 20, etc.

astrofrog commented 9 years ago

@bluescarni - could the minimum number of samples made to be a parameter in voronoi_grid? (defaulting to 10)

astrofrog commented 9 years ago

@bluescarni - could you create a new page in the docs in the advanced section that just describes how to use this? In particular it would be great if you can show a minimal example of how one would use this to find a density in each cell for example.

Using the output that you have at the moment, I think that the most efficient way would be something like (untested code):

def density(x, y, z):
    # stuff here
    return density_value

xyz = g.samples
idx = g.samples_idx

dens_all = density(xyz[0], xyz[1], xyz[2])

dens_average = np.add.reduceat(dens_all, idx[:-1]) / np.diff(idx)

Can you confirm that something like this works? If so, can you add an example like this to the docs?

astrofrog commented 9 years ago

@bluescarni - there was an issue with the Travis config which I just fixed in https://github.com/hyperion-rt/hyperion/pull/135 - can you rebase on the latest upstream master? (should be easy since no conflicts). Let me know if you need any instructions.

bluescarni commented 9 years ago

@astrofrog Thanks for the review, I am adding the configurable option for min samples and writing the docs. I will rebase before pushing the new iteration.

astrofrog commented 9 years ago

@bluescarni - apart from my comment above, this is good to go! Could you implement the small change I suggested today?

bluescarni commented 9 years ago

@astrofrog thanks for the review! I just committed the changes you suggested.

astrofrog commented 9 years ago

@bluescarni - thanks! I will go ahead and merge this. I think it would be nice to be able to access this functionality from the top-level VoronoiGrid class (the one users are normally exposed to) but that should be easy to add later. Let me know later if you have any ideas on how we could do this.

astrofrog commented 9 years ago

Actually I might open a quick PR now to add an interface in the top-level class to make it easier for Christine to use.

bluescarni commented 9 years ago

@astrofrog Thanks! Sounds good for the top level class interface, even though I have not much experience with it. I would imagine we can just expose the same interface as in the helper class and trickle down the construction parameters to voronoi_grid?

astrofrog commented 9 years ago

@bluescarni - yes, I'll open a PR shortly and will let you know so you can review it.

astrofrog commented 9 years ago

@bluescarni - see https://github.com/hyperion-rt/hyperion/pull/138