Changes for LZ framework from Quentin Riffard, LBNL. No functionality…

… difference for testNEST.
NESTCollaboration · Oct 1, 2019 · 18c3d09 · 18c3d09
1 parent b90777c
commit 18c3d09
Show file tree

Hide file tree

Showing 2 changed files with 36 additions and 9 deletions.
diff --git a/NEST.cpp b/NEST.cpp
@@ -814,7 +814,7 @@ vector<double> NESTcalc::GetS2(int Ne, double truthPos[3], double smearPos[3],
     double rho = fdetector->get_p_bar() * 1e5 /
                  (fdetector->get_T_Kelvin() * 8.314) * MOLAR_MASS * 1e-6;
     double driftVelocity_gas =
-        SetDriftVelocity_MagBoltz(rho, fdetector->get_E_gas() * 1000.);
+        GetDriftVelocity_MagBoltz(rho, fdetector->get_E_gas() * 1000.);
     double dt_gas = gasGap / driftVelocity_gas;
     double sigmaDLg = 10. * sqrt(2. * Diff_Long_Gas * dt_gas * 1e-6);
     double sigmaDTg = 10. * sqrt(2. * Diff_Tran_Gas * dt_gas * 1e-6);
@@ -1116,7 +1116,16 @@ vector<double> NESTcalc::GetSpike(int Nph, double dx, double dy, double dz,
 }
 
 double NESTcalc::SetDensity(double Kelvin,
-                            double bara) {  // currently only for fixed pressure
+                            double bara) {
+  bool inGas = false;
+  double density = GetDensity(Kelvin, bara, inGas);
+  fdetector->set_inGas(inGas);
+
+  return density;
+}
+
+double NESTcalc::GetDensity(double Kelvin,
+                            double bara, bool &inGas) {  // currently only for fixed pressure
                                             // (saturated vapor pressure); will
                                             // add pressure dependence later
 
@@ -1137,7 +1146,7 @@ double NESTcalc::SetDensity(double Kelvin,
         bara * 1e5 / (Kelvin * 8.314);  // ideal gas law approximation, mol/m^3
     density *= MOLAR_MASS * 1e-6;
     cerr << "\nWARNING: GAS PHASE. IS THAT WHAT YOU WANTED?\n";
-    fdetector->set_inGas(true);
+    inGas = true;
     return density;  // in g/cm^3
   }
 
@@ -1155,11 +1164,18 @@ double NESTcalc::SetDensity(double Kelvin,
                                   // zunzun fit to NIST data; will add gas later
 }
 
-double NESTcalc::SetDriftVelocity(double Kelvin, double Density,
-                                  double eField) {  // for liquid and solid only
+double NESTcalc::SetDriftVelocity(double Kelvin, double Density, double eField){
+  return GetDriftVelocity(Kelvin, Density, eField, fdetector->get_inGas());
+}
 
-  if (fdetector->get_inGas()) return SetDriftVelocity_MagBoltz(Density, eField);
+double NESTcalc::GetDriftVelocity(double Kelvin, double Density, double eField, bool inGas){
+  if (inGas) return GetDriftVelocity_MagBoltz(Density, eField);
+  else return GetDriftVelocity_Liquid(Kelvin, Density, eField);
+}
 
+double NESTcalc::GetDriftVelocity_Liquid(double Kelvin, double Density,
+                                  double eField) {  // for liquid and solid only
+
   double speed =
       0.0;  // returns drift speed in mm/usec. based on Fig. 14 arXiv:1712.08607
   int i, j;
@@ -1247,7 +1263,7 @@ double NESTcalc::SetDriftVelocity(double Kelvin, double Density,
   return speed;
 }
 
-double NESTcalc::SetDriftVelocity_MagBoltz(
+double NESTcalc::GetDriftVelocity_MagBoltz(
     double density, double efieldinput)  // Nichole Barry UCD 2011
 {
   density *= NEST_AVO / MOLAR_MASS;
@@ -1292,7 +1308,7 @@ vector<double> NESTcalc::SetDriftVelocity_NonUniform(double rho, double zStep,
                                                 fdetector->get_E_gas() /
                                                     (EPS_LIQ / EPS_GAS) * 1e3);
         else  // if gate == TopDrift properly set, shouldn't happen
-          driftTime += zStep / SetDriftVelocity_MagBoltz(
+          driftTime += zStep / GetDriftVelocity_MagBoltz(
                                    rho, fdetector->get_E_gas() * 1e3);
       } else
         driftTime +=

diff --git a/include/NEST.hh b/include/NEST.hh
@@ -225,7 +225,15 @@ class NESTcalc {
   // Gives one the drift velocity as a function of temperature and electric
   // field in liquid or solid. If density implies gas, kicks calculation down to
   // the next function below
-  double SetDriftVelocity_MagBoltz(double D, double F);
+  static double GetDriftVelocity(double T, double D, double F, bool inGas);
+  // Gives one the drift velocity as a function of temperature and electric
+  // field in liquid or solid. If density implies gas, kicks calculation down to
+  // the next function below
+  static double GetDriftVelocity_Liquid(double T, double D, double F);
+  // Gives one the drift velocity as a function of temperature and electric
+  // field in liquid or solid. If density implies gas, kicks calculation down to
+  // the next function below
+  static double GetDriftVelocity_MagBoltz(double D, double F);
   // Gas electron drift speed for S2 gas gap in 2-phase TPCs or the whole
   // detector for all gas. Based on simple fits to complicated MagBoltz software
   // output.
@@ -237,6 +245,9 @@ class NESTcalc {
   double SetDensity(double T, double P);
   // A simple, approximate but good, density is returned for solid, liquid, or
   // gaseous xenon, as a function of temperature and pressure
+  static double GetDensity(double T, double P, bool &inGas);
+  // A simple, approximate but good, density is returned for solid, liquid, or
+  // gaseous xenon, as a function of temperature and pressure
   std::vector<double> xyResolution(double xPos_mm, double yPos_mm,
                                    double A_top);
   // Utilizing a dependence on radius and the size of the S2 signal, takes MC