diff --git a/src/loadRamses.cpp b/src/loadRamses.cpp index 0eed042..d15566a 100644 --- a/src/loadRamses.cpp +++ b/src/loadRamses.cpp @@ -194,7 +194,7 @@ typedef struct double unitLength; double aexp; double boxSize; - double hubble; + double unit_t; } InfoData; int readInfoFile(const char *basename, int outputId, InfoData& info) @@ -249,10 +249,7 @@ int readInfoFile(const char *basename, int outputId, InfoData& info) info.unitLength = infoMap["unit_l"]; info.aexp = infoMap["aexp"]; info.boxSize = infoMap["boxlen"]; - - static const double CM_IN_MPC = 3.08e24; - double unit_t = infoMap["unit_t"]; - info.hubble = info.aexp*info.aexp/unit_t / (1e5/CM_IN_MPC); + info.unit_t = infoMap["unit_t"]; return 1; } @@ -270,11 +267,12 @@ CosmoTool::PurePositionData *CosmoTool::loadRamsesPosition(const char *basename, if (!readInfoFile(basename, outputId, info)) return 0; + double hubble = info.aexp*info.aexp/info.unit_t / (1e5/CM_IN_MPC); double L0 = info.boxSize*info.unitLength/CM_IN_MPC/info.aexp; if (!quiet) { cout << "L0=" << L0 << " Mpc" << endl; - cout << "H=" << info.hubble << " km/s/Mpc" << endl; + cout << "H=" << hubble << " km/s/Mpc" << endl; } if (!quiet) @@ -328,7 +326,7 @@ CosmoTool::PurePositionData *CosmoTool::loadRamsesPosition(const char *basename, assert(gd->pos != 0); gd->NumPart = totPart; gd->BoxSize = L0*1000; - gd->hubble = info.hubble; + gd->hubble = hubble; if (!quiet) cout << " Total number part=" << totPart << endl @@ -438,9 +436,13 @@ CosmoTool::PhaseSpaceData *CosmoTool::loadRamsesPhase(const char *basename, int if (!readInfoFile(basename, outputId, info)) return 0; + double unit_vel = info.unitLength/info.unit_t; + double hubble = info.aexp*info.aexp/info.unit_t / (1e5/CM_IN_MPC); double L0 = info.boxSize*info.unitLength/CM_IN_MPC/info.aexp; - if (!quiet) + if (!quiet) { cout << "L0=" << L0 << " Mpc" << endl; + cout << "H=" << hubble << " km/s/Mpc" << endl; + } if (!quiet) cout << "Detecting number of files and particles..." << endl; @@ -495,7 +497,8 @@ CosmoTool::PhaseSpaceData *CosmoTool::loadRamsesPhase(const char *basename, int assert(gd->vel != 0); gd->NumPart = totPart; gd->BoxSize = L0*1000; - + gd->hubble = hubble; + if (!quiet) cout << " Total number part=" << totPart << endl << "Loading particles ..." << endl; @@ -551,7 +554,7 @@ CosmoTool::PhaseSpaceData *CosmoTool::loadRamsesPhase(const char *basename, int infile.beginCheckpoint(); for (uint32_t i = 0; i < nPar; i++) { - vel[i][k] = infile.readReal32()*gd->BoxSize/1000*100; + vel[i][k] = infile.readReal32()*unit_vel; } infile.endCheckpoint(); }