Added a lean CIC filter

This commit is contained in:
Guilhem Lavaux 2014-06-12 17:16:59 +02:00
parent c80e2833e0
commit cef0fb3786
5 changed files with 53 additions and 3 deletions

View File

@ -3,6 +3,7 @@ from libcpp cimport string as cppstring
import numpy as np import numpy as np
cimport numpy as np cimport numpy as np
from cpython cimport PyObject, Py_INCREF from cpython cimport PyObject, Py_INCREF
cimport cython
np.import_array() np.import_array()
@ -47,6 +48,22 @@ cdef extern from "cosmopower.hpp" namespace "CosmoTool":
void setNormalization(double) void setNormalization(double)
double power(double) double power(double)
cdef extern from "cic.hpp" namespace "CosmoTool":
ctypedef float CICType
ctypedef float Coordinates[3]
cdef cppclass CICParticles:
float mass
Coordinates coords
cdef cppclass CICFilter:
CICFilter(np.uint32_t resolution, double L) nogil
void resetMesh() nogil
void putParticles(CICParticles* particles, np.uint32_t N) nogil
void getDensityField(CICType *& field, np.uint32_t& res) nogil
cdef extern from "loadSimu.hpp" namespace "CosmoTool": cdef extern from "loadSimu.hpp" namespace "CosmoTool":
@ -400,3 +417,32 @@ cdef class CosmologyPower:
return out return out
else: else:
return self._compute(k) return self._compute(k)
@cython.boundscheck(False)
def leanCic(float[:,:] particles, float L, int Resolution):
cdef CICParticles p
cdef CICFilter *cic
cdef np.uint64_t i
cdef CICType *field
cdef np.uint32_t dummyRes
cic = new CICFilter(Resolution, L)
cic.resetMesh()
if particles.shape[1] != 3:
raise ValueError("Particles must be Nx3 array")
p.mass = 1
for i in xrange(particles.shape[0]):
p.coords[0] = particles[i,0]
p.coords[1] = particles[i,1]
p.coords[2] = particles[i,2]
cic.putParticles(&p, 1)
field = <CICType*>0
dummyRes = 0
cic.getDensityField(field, dummyRes)
del cic

View File

@ -748,7 +748,7 @@ cdef npx.uint64_t _mysum(int[:] jobs) nogil:
def spherical_projection(int Nside, def spherical_projection(int Nside,
npx.ndarray[DTYPE_t, ndim=3] density not None, npx.ndarray[DTYPE_t, ndim=3] density not None,
DTYPE_t min_distance, DTYPE_t min_distance,
DTYPE_t max_distance, int progress=1, int integrator_id=0, DTYPE_t[:] shifter = None, int booster=10): DTYPE_t max_distance, int progress=1, int integrator_id=0, DTYPE_t[:] shifter = None, int booster=1000):
import healpy as hp import healpy as hp
import progressbar as pb import progressbar as pb

View File

@ -19,7 +19,7 @@ snap_id=int(sys.argv[1])
astart=1/100. astart=1/100.
if doSimulation: if doSimulation:
s = ct.loadRamsesAll("/nethome/lavaux/remote2/borgsim2/", snap_id, doublePrecision=True) s = ct.loadRamsesAll("/nethome/lavaux/remote2/borgsim3/", snap_id, doublePrecision=True)
astart=s.getTime() astart=s.getTime()
L = s.getBoxsize() L = s.getBoxsize()

View File

@ -38,6 +38,8 @@ knowledge of the CeCILL license and that you accept its terms.
#include <inttypes.h> #include <inttypes.h>
#include "cic.hpp" #include "cic.hpp"
using namespace CosmoTool;
CICFilter::CICFilter(uint32_t N, double len) CICFilter::CICFilter(uint32_t N, double len)
{ {
spatialLen = len; spatialLen = len;

View File

@ -39,7 +39,7 @@ knowledge of the CeCILL license and that you accept its terms.
#include "config.hpp" #include "config.hpp"
#include <inttypes.h> #include <inttypes.h>
using namespace CosmoTool; namespace CosmoTool {
typedef float CICType; typedef float CICType;
@ -67,4 +67,6 @@ typedef float CICType;
uint32_t szGrid; uint32_t szGrid;
}; };
};
#endif #endif