Fix LUT creation (#38)

- Enforcing seed set in pythia

Author: njacazio

examples/scripts/createO2tables.py

@@ -56,7 +56,8 @@ def run_cmd(cmd, comment=""):
         if content:
             verbose_msg("++", content.strip())
         if "Encountered error" in content:
-            msg("[WARNING] Error encountered runtime in", cmd, color=bcolors.BWARNING)
+            msg("[WARNING] Error encountered runtime in",
+                cmd, color=bcolors.BWARNING)
     except:
         fatal_msg("Error while running", cmd)
 
@@ -239,7 +240,8 @@ def write_to_runner(line, log_file=None, check_status=False):
                 if check_status:
                     f_run.write("\nReturnValue=$?\n")
                     f_run.write("if [[ $ReturnValue != 0 ]]; then\n")
-                    f_run.write(f"  echo \"Encountered error with command '{line.strip()}'\"\n")
+                    f_run.write(
+                        f"  echo \"Encountered error with command '{line.strip()}'\"\n")
                     f_run.write("  exit $ReturnValue\n")
                     f_run.write("fi\n")
 
@@ -270,7 +272,9 @@ def copy_and_link(file_name):
                 # Adjust configuration file
                 with open(generator_cfg, "a") as f_cfg:
                     # number of events and random seed
+                    f_cfg.write(f"\n\n\n#### Additional part ###\n\n\n\n")
                     f_cfg.write(f"Main:numberOfEvents {nevents}\n")
+                    f_cfg.write(f"Random:setSeed = on\n")
                     f_cfg.write(f"Random:seed = {run_number}\n")
                     # collision time spread [mm/c]
                     f_cfg.write("Beams:allowVertexSpread on \n")

src/lutWrite.its1.cc

@@ -25,5 +25,5 @@ void lutWrite_its1(const char *filename = "lutCovm.dat", int pdg = 211,
   // init FAT
   fatInit_its1(field, rmin);
   // write
-  lutWrite(filename, pdg, field, rmin);
+  lutWrite(filename, pdg, field);
 }

src/lutWrite.its2.cc

@@ -13,19 +13,21 @@ void fatInit_its2(float field = 0.5, float rmin = 100.) {
   fat.AddLayer((char *)"vertex", 0, 0); // dummy vertex for matrix calculation
 
   // new ideal Pixel properties?
-  Double_t x0 = 0.0050;
+  Double_t x0IB = 0.0035; // X/X0 of the inner barrel (first three layers)
+  Double_t x0OB = 0.008;  // X/X0 of the outer barrel
   Double_t resRPhi = 0.0006;
   Double_t resZ = 0.0006;
   Double_t eff = 0.98;
-  Double_t xrhoOB = 1.1646e-01; // 500 mum Si
+  Double_t xrhoIB = 1.1646e-02; // Surface density for 50 mum thick Si
+  Double_t xrhoOB = 1.1646e-01; // Surface density for 500 mum thick Si
 
-  fat.AddLayer((char *)"ddd1", 2.3, x0, xrhoOB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"ddd2", 3.1, x0, xrhoOB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"ddd3", 3.9, x0, xrhoOB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"ddd4", 19.4, x0, xrhoOB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"ddd5", 24.7, x0, xrhoOB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"ddd6", 35.3, x0, xrhoOB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"ddd7", 40.5, x0, xrhoOB, resRPhi, resZ, eff);
+  fat.AddLayer((char *)"ddd1", 2.3, x0IB, xrhoIB, resRPhi, resZ, eff);
+  fat.AddLayer((char *)"ddd2", 3.1, x0IB, xrhoIB, resRPhi, resZ, eff);
+  fat.AddLayer((char *)"ddd3", 3.9, x0IB, xrhoIB, resRPhi, resZ, eff);
+  fat.AddLayer((char *)"ddd4", 19.4, x0OB, xrhoOB, resRPhi, resZ, eff);
+  fat.AddLayer((char *)"ddd5", 24.7, x0OB, xrhoOB, resRPhi, resZ, eff);
+  fat.AddLayer((char *)"ddd6", 35.3, x0OB, xrhoOB, resRPhi, resZ, eff);
+  fat.AddLayer((char *)"ddd7", 40.5, x0OB, xrhoOB, resRPhi, resZ, eff);
 
   fat.AddTPC(0.1, 0.1); // TPC
   fat.SetAtLeastHits(4);
@@ -43,5 +45,5 @@ void lutWrite_its2(const char *filename = "lutCovm.dat", int pdg = 211,
   // init FAT
   fatInit_its2(field, rmin);
   // write
-  lutWrite(filename, pdg, field, rmin);
+  lutWrite(filename, pdg, field);
 }

src/lutWrite.its3.cc

@@ -5,28 +5,49 @@
 
 #include "lutWrite.cc"
 
+// Adds foam spacers between inner layers
+const bool add_foam = true;
+// Puts foam spacers mid-way betwen layers
+const bool foam_middle = false;
+
 void fatInit_its3(float field = 0.5, float rmin = 100.) {
   fat.SetBField(field);
   fat.SetdNdEtaCent(400.);
 
-  fat.AddLayer((char *)"bpipe", 1.6, 0.0022); // beam pipe
-  fat.AddLayer((char *)"vertex", 0, 0); // dummy vertex for matrix calculation
+  Double_t x0BP = 0.0014;     // 500 mum Be
+  Double_t xrhoBP = 9.24e-02; // 500 mum Be
+
+  // dummy vertex for matrix calculation
+  fat.AddLayer((char *)"vertex", 0.0, 0, 0);
+  fat.AddLayer((char *)"bpipe", 1.625, x0BP, xrhoBP); // 500 mum Be beam pipe
 
-  // new ideal Pixel properties?
   Double_t x0IB = 0.0005;       // X/X0 of the inner barrel (first three layers)
-  Double_t x0OB = 0.0035;       // X/X0 of the outer barrel
+  Double_t x0OB = 0.008;        // X/X0 of the outer barrel
   Double_t resRPhi = 0.0006;    // Resolution in Rphi
   Double_t resZ = 0.0006;       // Resolution in Z
   Double_t eff = 0.98;          // Efficiency
   Double_t xrhoIB = 1.1646e-02; // Surface density for 50 mum thick Si
   Double_t xrhoOB = 1.1646e-01; // Surface density for 500 mum thick Si
 
+  const Double_t x0Foam = 0.0008;   // X0 ~710cm for 0.6cm thick foam
+  float foam_radius = 9.370 * x0IB; // width of a chip
+
   fat.AddLayer((char *)"ddd1", 1.8, x0IB, xrhoIB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"foam1", 1.8 + 9.370 * x0IB, 0.0008); // Foam spacer
+  if (foam_middle)
+    foam_radius = (2.4 - 1.8) / 2;
+  if (add_foam)
+    fat.AddLayer((char *)"foam1", 1.8 + foam_radius, x0Foam); // Foam spacer
   fat.AddLayer((char *)"ddd2", 2.4, x0IB, xrhoIB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"foam2", 2.4 + 9.370 * x0IB, 0.0008); // Foam spacer
+  if (foam_middle)
+    foam_radius = (3.0 - 2.4) / 2;
+  if (add_foam)
+    fat.AddLayer((char *)"foam2", 2.4 + foam_radius, x0Foam); // Foam spacer
   fat.AddLayer((char *)"ddd3", 3.0, x0IB, xrhoIB, resRPhi, resZ, eff);
-  fat.AddLayer((char *)"foam3", 3.0 + 9.370 * x0IB, 0.0008); // Foam spacer
+  if (foam_middle)
+    foam_radius = (19.4 - 3.0) / 2;
+  if (add_foam)
+    fat.AddLayer((char *)"foam3", 3.0 + foam_radius, x0Foam); // Foam spacer
+  // Structural cylinder?
   fat.AddLayer((char *)"ddd4", 19.4, x0OB, xrhoOB, resRPhi, resZ, eff);
   fat.AddLayer((char *)"ddd5", 24.7, x0OB, xrhoOB, resRPhi, resZ, eff);
   fat.AddLayer((char *)"ddd6", 35.3, x0OB, xrhoOB, resRPhi, resZ, eff);
@@ -48,5 +69,5 @@ void lutWrite_its3(const char *filename = "lutCovm.dat", int pdg = 211,
   // init FAT
   fatInit_its3(field, rmin);
   // write
-  lutWrite(filename, pdg, field, rmin);
+  lutWrite(filename, pdg, field);
 }