diff --git a/sample/CMakeLists.txt b/sample/CMakeLists.txt index a333aff..606c8d1 100644 --- a/sample/CMakeLists.txt +++ b/sample/CMakeLists.txt @@ -77,3 +77,6 @@ if (Boost_FOUND) add_executable(simple3DFilter simple3DFilter.cpp) target_link_libraries(simple3DFilter ${tolink}) endif (Boost_FOUND) + +add_executable(gadgetToDensity gadgetToDensity.cpp) +target_link_libraries(gadgetToDensity ${tolink}) diff --git a/sample/gadgetToDensity.cpp b/sample/gadgetToDensity.cpp new file mode 100644 index 0000000..d17f263 --- /dev/null +++ b/sample/gadgetToDensity.cpp @@ -0,0 +1,62 @@ +#include +#include +#include +#include "cic.hpp" +#include "loadGadget.hpp" +#include "miniargs.hpp" +#include "yorick.hpp" + +using namespace std; +using namespace CosmoTool; + +int main(int argc, char **argv) +{ + uint32_t res; + char *fname; + int id; + + MiniArgDesc desc[] = { + { "SNAPSHOT", &fname, MINIARG_STRING }, + { "ID", &id, MINIARG_INT }, + { "RESOLUTION", &res, MINIARG_INT }, + { 0, 0, MINIARG_NULL } + }; + + if (!parseMiniArgs(argc, argv, desc)) + return 1; + + SimuData *p = loadGadgetMulti(fname, 0, 0); + double L0 = p->BoxSize/1000; + CICFilter filter(res, L0); + + delete p; + + try { + for (int cpuid=0;;cpuid++) { + p = loadGadgetMulti(fname, cpuid, NEED_POSITION); + for (uint32_t i = 0; i < p->NumPart; i++) + { + CICParticles a; + + a.mass = 1.0; + a.coords[0] = p->Pos[0][i]/1000; + a.coords[1] = p->Pos[1][i]/1000; + a.coords[2] = p->Pos[2][i]/1000; + filter.putParticles(&a, 1); + } + + delete p; + } + } catch (const NoSuchFileException& e) {} + + CICType *denField; + uint32_t Ntot; + filter.getDensityField(denField, Ntot); + + cout << "L0=" << L0 << endl; + cout << "Saving density field" << endl; + uint32_t dimList[] = { res, res, res}; + saveArray("densityField.nc", denField, dimList, 3); + + return 0; +} diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 8901027..611deaa 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -9,6 +9,7 @@ SET(CosmoTool_SRCS miniargs.cpp growthFactor.cpp cosmopower.cpp + cic.cpp ) IF(FOUND_NETCDF3)