Added automatic construction of 2021 test beam geometry

MonteCarlo/geometry.mac

 
 /Detector/Optical     true
+/Detector/ReducedTree true
 /Detector/Overlaps    false
 /Detector/PI0         false
-/Detector/RunNumber 347
-/Detector/PrintDebugStatement TEST LINE
+/Detector/PrintDebugStatement TEST_LINE #Simply prints TEST_LINE to terminal. Used for debugging crashes
+/Detector/TestBeamYear 2021
+/Detector/RunNumber    843
+/Detector/ConfigFile Alignment_2018.xml
 
-/Detector/ForcePosition true
-#/Detector/ConfigFile Alignment_2018.xml
+# Force position means you set the position of detectors using the macro commands below
+# If true, the macro commands below will determine the geometry
+# If false, TestBeamYear and RunNumber set above will be used to simulate test beam geometry
+/Detector/ForcePosition false
+
+# World volue
 /Detector/SetWorldDimensions 2. 2. 260. m
 
-#EMCAL
+# EMCAL
 /Detector/ZDC/Add
+/Detector/ZDC/ZDCtype TB2021EM
 /Detector/ZDC/Position 0. 0. 100. mm
-/Detector/ZDC/FiberDiameters 1.5 0. 0. mm
-/Detector/ZDC/AbsorberDimensions 90. 180. 10.2 mm #run3
-/Detector/ZDC/nAbsorbers 11
-/Detector/ZDC/HousingThickness 5. mm
-/Detector/ZDC/GapThickness 2. mm
-/Detector/ZDC/SteelAbsorberHeight 290. mm
-/Detector/ZDC/Optical true
-/Detector/ZDC/CheckOverlaps false
-/Detector/ZDC/ReducedTree true
-/Detector/ZDC/HousingMaterial aluminum
-/Detector/ZDC/AbsorberMaterial pure
-
-/Detector/ZDC/Duplicate 1
+# /Detector/ZDC/FiberDiameters 1.5 0. 0. mm
+# /Detector/ZDC/AbsorberDimensions 90. 180. 10.2 mm #run3
+# /Detector/ZDC/nAbsorbers 11
+# /Detector/ZDC/HousingThickness 5. mm
+# /Detector/ZDC/GapThickness 2. mm
+# /Detector/ZDC/SteelAbsorberHeight 290. mm
+# /Detector/ZDC/Optical true
+# /Detector/ZDC/CheckOverlaps false
+# /Detector/ZDC/ReducedTree true
+# /Detector/ZDC/HousingMaterial aluminum
+# /Detector/ZDC/AbsorberMaterial pure
+
+# HAD1
+/Detector/ZDC/Add
+/Detector/ZDC/ZDCtype ATLASHAD
 /Detector/ZDC/Position 0. 0. 271. mm
 
+# HAD2
 /Detector/ZDC/Duplicate 2
 /Detector/ZDC/Position 0. 0. 423. mm
 
+# HAD3
 /Detector/ZDC/Duplicate 3
 /Detector/ZDC/Position 0. 0. 575. mm
 
+# Pan Flute RPD
 /Detector/RPD/Add
+/Detector/RPD/RPDtype TB2021PF
 /Detector/RPD/Position 0. 0. 185. mm
-/Detector/RPD/FiberDiameters .60 .66 .71 mm
-/Detector/RPD/HousingThickness 2. mm
-/Detector/RPD/FiberPitchX 1.2 mm
-/Detector/RPD/FiberPitchZ 2.0 mm
-/Detector/RPD/TileSize 9.6 mm
-/Detector/RPD/MinWallThickness 0.01 mm
-/Detector/RPD/FiberReadoutDistance 360 mm
-/Detector/RPD/RPDtype panflute
-/Detector/RPD/Optical true
-/Detector/RPD/CheckOverlaps false
-/Detector/RPD/ReducedTree FASTSCINTILLATIONRISETIME
-
-#Trigger
-#/Detector/TRI/Add
-#/Detector/TRI/Position 0. 0. -3183. mm
+# /Detector/RPD/FiberDiameters .60 .66 .71 mm
+# /Detector/RPD/HousingThickness 2. mm
+# /Detector/RPD/FiberPitchX 1.2 mm
+# /Detector/RPD/FiberPitchZ 2.0 mm
+# /Detector/RPD/TileSize 9.6 mm
+# /Detector/RPD/MinWallThickness 0.01 mm
+# /Detector/RPD/FiberReadoutDistance 360 mm
+# /Detector/RPD/Optical true
+# /Detector/RPD/CheckOverlaps false
+# /Detector/RPD/ReducedTree false
+
+# # Trigger
+/Detector/TRI/Add
+/Detector/TRI/Position 0. 0. -3183. mm

MonteCarlo/include/DetectorConstruction.hh

@@ -43,6 +43,8 @@
 #include "G4Cache.hh"
 #include "G4MagneticField.hh"
 
+#include "TQuaternion.h"
+
 #include <vector>
 
 
@@ -84,6 +86,21 @@ class Survey {
 
 };
 
+class Survey2021{
+
+public:
+  /** Name of the detector this survey entry is for */
+  G4String detector;
+  /** Position of the Detector during the survey*/
+  G4ThreeVector pos;
+  /** Table position during survey */
+  G4ThreeVector table;
+  /** Trigger table position during survey */
+  G4ThreeVector trig_table;
+  /** Rotation matrix for this detector*/
+  G4RotationMatrix rotation;
+};
+
 class Alignment {
 
  public:
@@ -107,6 +124,40 @@ class Alignment {
 
 };
 
+
+class Alignment2021 {
+
+ public:
+
+    /** Run number being analyzed **/
+    int runNumber;
+    /** Beam type - can be p or e in H2 2021 **/
+    std::string beam_type;
+    /** Beam energy **/
+    double beam_energy;
+    /** X position of the NIKHEF Table **/
+    double x_det_table;
+    /** Y position of the NIKHEF Table **/
+    double y_det_table;
+    /** X position of the Trigger Table **/
+    double x_trig_table;
+    /** Y position of the Trigger Table **/
+    double y_trig_table;
+    /** First detector met by the beam **/
+    std::string Det1;
+    /** Second detector met by the beam **/
+    std::string Det2;
+    /** Third detector met by the beam **/
+    std::string Det3;
+    /** Fourth detector met by the beam **/
+    std::string Det4;
+    /** Fifth detector met by the beam **/
+    std::string Det5;
+    /** Sixth detector met by the beam **/
+    std::string Det6;
+};
+
+
 /// Detector construction class to define materials and geometry.
 
 class DetectorConstruction : public G4VUserDetectorConstruction
@@ -117,22 +168,28 @@ public:
 
   virtual G4VPhysicalVolume* Construct();
   virtual G4VPhysicalVolume* ConstructWorldVolume(G4double x, G4double y, G4double z);
-  virtual G4VPhysicalVolume* ConstructSPSTestBeam();
+  virtual G4VPhysicalVolume* Construct2018SPSTestBeam();
+  virtual G4VPhysicalVolume* Construct2021SPSTestBeam();
   virtual G4VPhysicalVolume* ManualConstruction();
   virtual void ConstructSDandField();
 
   virtual void    LoadConfigurationFile ( G4String _inFile = "" );
   virtual void    LoadAlignmentFile     ( G4String _inFile = "" );
+  virtual void    LoadAlignmentFile2021 ( G4String _inFile = "" );
+  virtual void    LoadSurveyFile2021    ( G4String _inFile = "" );
   virtual void    AddZDC                ( G4ThreeVector* position = NULL );
   virtual void    AddRPD                ( G4ThreeVector* position = NULL );
   virtual void    AddTRI                ( G4ThreeVector* position = NULL );
   virtual void    SetConfigFileName     ( G4String _name ){ m_configFileName = _name; }
   virtual void    SetRunNumber          ( G4int _run ){ m_runNumber = _run; }
+  virtual void    SetTestYear           ( G4int _year){ m_year = _year; }
   inline  void    SetOverlapsFlag       ( G4bool arg ){CHECK_OVERLAPS = arg;}
   inline  void    SetOpticalFlag        ( G4bool arg ){OPTICAL = arg;}
+  inline  void    SetReducedTreeFlag    ( G4bool arg ){REDUCED_TREE = arg;}
   inline  void    ForcePosition         ( G4bool arg ){ForceDetectorPosition = arg;}
   inline  void    SetPI0Flag            ( G4bool arg ){PI0 = arg;}
 
+  virtual G4RotationMatrix QuaternionToRotationMatrix( TQuaternion Q );
   virtual Survey*    GetSurvey             ( G4String name );
   inline  Alignment  GetAlignment          ( ){ return *m_alignment;    }
   inline  G4bool     GetOverlapsFlag       ( ){ return CHECK_OVERLAPS;  }
@@ -153,6 +210,7 @@ public:
   inline  void SetZDCHousingThickness  ( G4double       arg ){ m_ZDCvec.at(currentZDC-1)->SetHousingThickness(arg);   }
   inline  void SetZDCGapThickness      ( G4double       arg ){ m_ZDCvec.at(currentZDC-1)->SetGapThickness(arg);       }
   inline  void SetZDCSteelAbsHeight    ( G4double       arg ){ m_ZDCvec.at(currentZDC-1)->SetSteelAbsHeight(arg);     }
+  inline  void SetZDCDetectorType      ( G4String       arg ){ m_ZDCvec.at(currentZDC-1)->SetDetectorType(arg);    }
   inline  void SetZDCOpticalFlag       ( G4bool         arg ){ m_ZDCvec.at(currentZDC-1)->SetOpticalFlag(arg);        }
   inline  void SetZDCOverlapsFlag      ( G4bool         arg ){ m_ZDCvec.at(currentZDC-1)->SetOverlapsFlag(arg);       }
   inline  void SetZDCReducedTreeFlag   ( G4bool         arg ){ m_ZDCvec.at(currentZDC-1)->SetReducedTreeFlag(arg);    }
@@ -198,9 +256,12 @@ protected:
   /* Configuration parser objects and storage */
   XMLSettingsReader*      m_XMLparser;
   Alignment* 	            m_alignment;
+  Alignment2021*          m_alignment2021;
   Survey*                 m_survey;
   G4int                   m_runNumber;
+  G4int                   m_year;
   std::vector < Survey* > m_surveyEntries;
+  std::vector < Survey2021* > m_surveyEntries2021;
 
   /* Number of each detector */
   std::vector< ModTypeZDC* > m_ZDCvec;
@@ -217,6 +278,7 @@ private:
   /* Run condition flags */
   G4bool m_GFlash;
   G4bool OPTICAL;
+  G4bool REDUCED_TREE;
   G4bool ForceDetectorPosition;
   G4bool CHECK_OVERLAPS;
   G4bool PI0;

MonteCarlo/include/DetectorMessenger.hh

@@ -62,11 +62,13 @@ class DetectorMessenger: public G4UImessenger{
     G4UIdirectory*             fZDCDir;
     G4UIdirectory*             fTRIDir;
 
+    G4UIcmdWithABool*          fReducedTreeCmd;
     G4UIcmdWithABool*          fOpticalCmd;
     G4UIcmdWithABool*          fOverlapsCmd;
     G4UIcmdWithABool*          fPI0Cmd;
     G4UIcmdWithABool*          fForcePositionCmd;
     G4UIcmdWithAnInteger*      fSetRunNumberCmd;
+    G4UIcmdWithAnInteger*      fSetTestYearCmd;
     G4UIcmdWithAString*        fConfigFileCmd;
     G4UIcmdWithAString*        fPrintDebugCmd;
     G4UIcmdWith3VectorAndUnit* fSetWorldDimensionsCmd;
@@ -81,6 +83,7 @@ class DetectorMessenger: public G4UImessenger{
     G4UIcmdWithADoubleAndUnit* fZDCHousingThicknessCmd;
     G4UIcmdWithADoubleAndUnit* fZDCGapThicknessCmd;
     G4UIcmdWithADoubleAndUnit* fZDCSteelAbsHeightCmd;
+    G4UIcmdWithAString*        fZDCTypeCmd;
     G4UIcmdWithABool*          fZDCOpticalFlagCmd;
     G4UIcmdWithABool*          fZDCOverlapsFlagCmd;
     G4UIcmdWithABool*          fZDCReducedTreeCmd;

MonteCarlo/include/ModTypeZDC.hh

@@ -62,6 +62,7 @@ public:
   inline  void  SetSteelAbsHeight    ( G4double       arg ){ STEEL_ABSORBER = true; m_SteelAbsHeight = arg; }
   inline  void  SetHousingThickness  ( G4double       arg ){ m_HousingThickness  = arg; }
   inline  void  SetGapThickness      ( G4double       arg ){ m_GapThickness      = arg; }
+  inline  void  SetDetectorType      ( G4String       arg ){ m_detType           = arg; }
   inline  void  SetOpticalFlag       ( G4bool         arg ){ OPTICAL             = arg; }
   inline  void  SetOverlapsFlag      ( G4bool         arg ){ CHECK_OVERLAPS      = arg; }
   inline  void  SetReducedTreeFlag   ( G4bool         arg ){ REDUCED_TREE        = arg; }
@@ -89,6 +90,7 @@ protected:
   G4double         m_GapThickness;
   G4double         m_SteelAbsHeight;
   G4double         m_topOfVolume;
+  G4String         m_detType;
   G4bool           OPTICAL;
   G4bool           CHECK_OVERLAPS;
   G4bool           STEEL_ABSORBER;

MonteCarlo/src/DetectorConstruction.cc

@@ -79,8 +79,11 @@
 
   # UI commands and their functions
   ## Commands for the world
+  - <b>/Detector/ReducedTree</b> true/false
+    - Sets global ReducedTree mode. Used for test beam construction
+
   - <b>/Detector/Optical</b> true/false
-    - Optical (Cherenkov) photons are killed in the world volume if false
+    - Sets global optical flag. Used for test beam construction. Detectors will only record optical photons if true
 
   - <b>/Detector/Overlaps</b> true/false
     - Check geometry for overlaps if true
@@ -148,6 +151,9 @@
   - <b>/Detector/ZDC/AbsorberMaterial</b> composite/pure
     - Chose between nickel tungsten composite or pure tungsten absorber material
 
+  - <b>/Detector/ZDC/ZDCtype</b> ATLASHAD/TB2021EM
+    - Chose between the the panflute prototype design and the CMS design
+
   ## RPD commands
   - <b>/Detector/RPD/FiberDiameters</b> core cladding buffer unit
     - Set the core, cladding, and buffer diameters of the optical fibers
@@ -167,8 +173,8 @@
   - <b>/Detector/RPD/FiberReadoutDistance</b> distance unit
     - Set the distance between the top of the top tile and the top of the fiber (PMT window)
 
-  - <b>/Detector/RPD/RPDtype</b> panflute/CMS
-    - Chose between the the panflute prototype design and the CMS design
+  - <b>/Detector/RPD/RPDtype</b> panflute/CMS/TB2021PF
+    - Chose between the the 2018 or 2021 panflute prototype design and the CMS tile design
 
 */
 
@@ -217,6 +223,8 @@
 #include "GFlashHitMaker.hh"
 #include "GFlashParticleBounds.hh"
 
+#include "Math/GenVector/Quaternion.h"
+
 #include <iostream>
 #include <stdio.h>
 
@@ -231,6 +239,7 @@ DetectorConstruction::DetectorConstruction(G4bool _GFlash)
   m_GFlash(_GFlash),
   OPTICAL(false),
   ForceDetectorPosition(false),
+  m_TRI(0),
   PI0(false)
 {
   new DetectorMessenger(this);
@@ -264,8 +273,18 @@ G4VPhysicalVolume* DetectorConstruction::Construct(){
   if( ForceDetectorPosition ){
     ManualConstruction();
   }else{
-    ConstructSPSTestBeam();
-    std::cout << "built test beam" << std::endl;
+    //Clear any existing geometry
+    for(ModTypeZDC *zdc : m_ZDCvec) delete zdc;
+    m_ZDCvec.clear();
+    currentZDC = -1;
+    for(ModTypeRPD *rpd : m_RPDvec) delete rpd;
+    m_RPDvec.clear();
+    currentRPD = -1;
+    if(m_TRI) delete m_TRI;
+    m_TRI = 0;
+
+    if(m_year == 2021) Construct2021SPSTestBeam();
+    if(m_year == 2018) Construct2018SPSTestBeam();
   }
   return m_physWorld;
 }
@@ -306,9 +325,9 @@ G4VPhysicalVolume* DetectorConstruction::ConstructWorldVolume(G4double x, G4doub
 
 G4VPhysicalVolume* DetectorConstruction::ManualConstruction(){
 
-  std::cout << "******************************************" << std::endl
-            << "        PLACING DETECTORS MANUALLY        " << std::endl
-            << "******************************************" << std::endl;
+  std::cout << "\n******************************************" << std::endl
+            << "*       PLACING DETECTORS MANUALLY       *" << std::endl
+            << "******************************************\n" << std::endl;
 
   G4ThreeVector* pos;
   char name[32];
@@ -329,7 +348,7 @@ G4VPhysicalVolume* DetectorConstruction::ManualConstruction(){
     printf( "RPD%d center = (%f,%f,%f)\n", rpd->GetModNum(), pos->x(), pos->y(), pos->z() );
     rpd->Construct();
   }
-if(isTrigger){
+  if(isTrigger){
     pos = m_TRI->GetPosition();
     printf( "Trigger%d center = (%f,%f,%f)\n", 1, pos->x(), pos->y(), pos->z() );
     m_TRI ->Construct();
@@ -340,7 +359,122 @@ if(isTrigger){
 
 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
 
-G4VPhysicalVolume* DetectorConstruction::ConstructSPSTestBeam(){
+G4VPhysicalVolume* DetectorConstruction::Construct2021SPSTestBeam(){
+
+  std::cout << "\n******************************************" << std::endl
+            << "*    BUILDING 2021 TESTBEAM, RUN " << std::left << setw(6) << m_runNumber << "  *" << std::endl
+            << "******************************************\n" << std::endl;
+
+  LoadAlignmentFile2021();
+  LoadSurveyFile2021();
+
+  G4double worldSizeX = 2000*mm;
+  G4double worldSizeY = 2000*mm;
+  G4double worldSizeZ = 32000*mm;
+
+  if(m_logicWorld){
+    m_solidWorld->SetXHalfLength(worldSizeX);
+    m_solidWorld->SetXHalfLength(worldSizeY);
+    m_solidWorld->SetXHalfLength(worldSizeZ);
+  }else{
+    ConstructWorldVolume( worldSizeX, worldSizeY, worldSizeZ );
+  }
+
+  vector< G4String > detectors = {m_alignment2021->Det1, m_alignment2021->Det2,
+                                  m_alignment2021->Det3, m_alignment2021->Det4,
+                                  m_alignment2021->Det5, m_alignment2021->Det6};
+
+  G4cout << "Run " << m_alignment2021->runNumber << " contains the following detectors:";
+  for( int det = 0; det < detectors.size(); det++){
+    if(detectors[det] != "none"){
+      G4cout << detectors[det] << ", ";
+      if(detectors[det] == "ZDC"){
+        detectors[det] = "HAD";
+      }
+    }
+  }
+  G4cout << G4endl;
+
+  std::vector< G4String > finishedHads;
+  G4ThreeVector* pos;
+  for( Survey2021* survey : m_surveyEntries2021 ){
+    for( G4String det : detectors){
+      if( survey->detector.contains(det) ){
+
+        /********** ATLAS run3 HAD module ************/
+        if( det.contains("HAD") ){
+          //Check if this had has already been constructed
+          G4bool finished = false;
+          for( G4String done : finishedHads ){
+            if(survey->detector == done) finished = true;
+          }
+          if(finished) break;
+
+          AddZDC( new G4ThreeVector(  m_alignment2021->x_det_table*mm - survey->table.x()*mm + survey->pos.x()*mm,
+                                      m_alignment2021->y_det_table*mm - survey->table.y()*mm + survey->pos.y()*mm,
+                                      survey->pos.z()*mm ) );
+          ModTypeZDC* zdc = m_ZDCvec.back();
+          zdc->SetDetectorType("ATLASHAD");
+          zdc->SetReducedTreeFlag(REDUCED_TREE);
+          zdc->SetOpticalFlag(OPTICAL);
+          pos = zdc->GetPosition();
+          printf( "ZDC%d center = (%f,%f,%f)\n", zdc->GetModNum(), pos->x(), pos->y(), pos->z() );
+          zdc->Construct();
+          finishedHads.push_back(survey->detector);
+
+        /********** 2021 Test Beam EM Prototype ************/
+        }else if( det.contains("UEM") ){
+          AddZDC( new G4ThreeVector(  m_alignment2021->x_det_table*mm - survey->table.x()*mm + survey->pos.x()*mm,
+                                      m_alignment2021->y_det_table*mm - survey->table.y()*mm + survey->pos.y()*mm,
+                                      survey->pos.z()*mm ) );
+          ModTypeZDC* zdc = m_ZDCvec.back();
+          zdc->SetDetectorType("TB2021EM");
+          zdc->SetReducedTreeFlag(REDUCED_TREE);
+          zdc->SetOpticalFlag(OPTICAL);
+          pos = zdc->GetPosition();
+          printf( "ZDC%d center = (%f,%f,%f)\n", zdc->GetModNum(), pos->x(), pos->y(), pos->z() );
+          zdc->Construct();
+
+        /********** 2021 Test Beam Pan Flute RPD ************/
+        }else if( det.contains("PFRPD") ){
+          AddRPD( new G4ThreeVector(  m_alignment2021->x_det_table*mm - survey->table.x()*mm + survey->pos.x()*mm,
+                                      m_alignment2021->y_det_table*mm - survey->table.y()*mm + survey->pos.y()*mm,
+                                      survey->pos.z()*mm ) );
+          ModTypeRPD* rpd = m_RPDvec.back();
+          rpd->SetDetectorType("TB2021PF");
+          rpd->SetReducedTreeFlag(REDUCED_TREE);
+          rpd->SetOpticalFlag(OPTICAL);
+          pos = rpd->GetPosition();
+          printf( "RPD%d center = (%f,%f,%f)\n", rpd->GetModNum(), pos->x(), pos->y(), pos->z() );
+          rpd->Construct();
+
+        }//end detector type checks
+      }//end detector loop
+    }//end if survey::detector.contains(det)
+
+    /********** 2021 Test Beam Trigger Paddles ************/
+    if( survey->detector.contains("TRI") ){
+      AddTRI( new G4ThreeVector(  m_alignment2021->x_trig_table*mm - survey->trig_table.x()*mm + survey->pos.x()*mm,
+                                  m_alignment2021->y_trig_table*mm - survey->trig_table.y()*mm + survey->pos.y()*mm,
+                                  survey->pos.z()*mm ) );
+      pos = m_TRI->GetPosition();
+      printf( "Trigger center = (%f,%f,%f)\n", pos->x(), pos->y(), pos->z() );
+      m_TRI->Construct();
+    }
+  }//end survey loop
+  G4cout << "2021 Test Beam Construction Completed!!!" << G4endl;
+  return m_physWorld;
+}
+
+//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
+
+G4VPhysicalVolume* DetectorConstruction::Construct2018SPSTestBeam(){
+
+
+  std::cout << "\n******************************************" << std::endl
+            << "*    BUILDING 2018 TESTBEAM, RUN " << std::left << setw(6) << m_runNumber << "  *" << std::endl
+            << "******************************************\n" << std::endl;
+
   // Create variables to be used in beamtest 2018 simulation
   G4ThreeVector zdc1Pos,zdc2Pos, rpdPos;
   bool ZDC1 = false, ZDC2 = false, RPD = false;
@@ -355,7 +489,13 @@ G4VPhysicalVolume* DetectorConstruction::ConstructSPSTestBeam(){
 
   //################################ World volume construction
 
+  if(m_logicWorld){
+    m_solidWorld->SetXHalfLength(worldSizeX);
+    m_solidWorld->SetXHalfLength(worldSizeY);
+    m_solidWorld->SetXHalfLength(worldSizeZ);
+  }else{
     ConstructWorldVolume( worldSizeX, worldSizeY, worldSizeZ );
+  }
 
   //################################ SURVEY/ALIGNMENT_SETUP
 
@@ -535,7 +675,8 @@ void DetectorConstruction::ConstructSDandField( ){
   for(ModTypeRPD* rpd : m_RPDvec){
     rpd->ConstructSDandField();
   }
-if(isTrigger){
+
+  if(isTrigger){
     m_TRI->ConstructSDandField();
   }
 }
@@ -596,6 +737,51 @@ void DetectorConstruction::LoadConfigurationFile( G4String _inFile  ){
     return;
 }
 
+/**
+ * @brief Reads the .xml configuration file and load characteristics for all the channels, immediately sorted into detectors objects
+ * @param _inFile
+ */
+void DetectorConstruction::LoadAlignmentFile2021( G4String _inFile ){
+    if( _inFile == "" ) _inFile = std::getenv("JZCaPA") + std::string("/Utils/Alignment_2021.xml");
+    m_XMLparser = new XMLSettingsReader();
+
+    if (!m_XMLparser->parseFile(_inFile)) {
+            std::cerr << " Data Reader could not parse file : " << _inFile << std::endl;
+            return;
+    }
+
+    m_alignment2021 = new Alignment2021();
+    m_alignment2021->runNumber = m_runNumber;
+
+    std::cout << "Loading .xml Alignment File..." << std::endl;
+    std::cout << "Found " << m_XMLparser->getBaseNodeCount("Alignment") << " alignment entries " << std::endl;
+    std::cout << "Retrieving the information for run " << m_runNumber << std::endl;
+
+    int run;
+    for (int i = 0; i < m_XMLparser->getBaseNodeCount("Alignment"); i++) {
+        m_XMLparser->getChildValue("Alignment",i,"run",run);
+        if(run != m_runNumber) continue;
+        std::cout << "Found Run Entry in Alignment file for run " << m_runNumber << std::endl;
+        m_XMLparser->getChildValue("Alignment",i,"beam_type",m_alignment2021->beam_type);
+        m_XMLparser->getChildValue("Alignment",i,"beam_energy",m_alignment2021->beam_energy);
+        m_XMLparser->getChildValue("Alignment",i,"x_table",m_alignment2021->x_det_table);
+        m_XMLparser->getChildValue("Alignment",i,"y_table",m_alignment2021->y_det_table);
+        m_XMLparser->getChildValue("Alignment",i,"x_trigger",m_alignment2021->x_trig_table);
+        m_XMLparser->getChildValue("Alignment",i,"y_trigger",m_alignment2021->y_trig_table);
+        m_XMLparser->getChildValue("Alignment",i,"Det1",m_alignment2021->Det1);
+        m_XMLparser->getChildValue("Alignment",i,"Det2",m_alignment2021->Det2);
+        m_XMLparser->getChildValue("Alignment",i,"Det3",m_alignment2021->Det3);
+        m_XMLparser->getChildValue("Alignment",i,"Det4",m_alignment2021->Det4);
+        m_XMLparser->getChildValue("Alignment",i,"Det5",m_alignment2021->Det5);
+        m_XMLparser->getChildValue("Alignment",i,"Det6",m_alignment2021->Det6);
+    }
+
+    if(m_alignment2021 == NULL) std::cout << "WARNING: ALIGNMENT NOT FOUND!!!" << std::endl;
+    delete m_XMLparser; m_XMLparser = NULL;
+    return;
+
+}
+
 /**
  * @brief Reads the .xml configuration file and load characteristics for all the channels, immediately sorted into detectors objects
  * @param _inFile
@@ -647,6 +833,144 @@ void DetectorConstruction::LoadAlignmentFile( G4String _inFile ){
     return;
 }
 
+
+/** @brief Loads calibrated timing information for DRS4 modules based on run number
+ *  @param _inFile Optional argument for loading a custom file
+ *
+ *  Uses run number to determine scan number and loads the appropriate DRS4 timing
+ *  information from the 2018 Testbeam. After loading, we hand each Channel a
+ *  pointer to its timing vector. Must be run after LoadConfigurationFile()
+ *
+ */
+void DetectorConstruction::LoadSurveyFile2021( G4String _inFile ){
+  if( _inFile == "" ) _inFile = std::getenv("JZCaPA") + std::string("/Utils/Survey_2021.xml");
+  m_XMLparser = new XMLSettingsReader();
+
+  if (!m_XMLparser->parseFile(_inFile)) {
+          std::cerr << " Data Reader could not parse file : " << _inFile << std::endl;
+          return;
+  }
+
+  int first_run, last_run;
+  double x,y,z,tx,ty,ttx,tty,qw,qi,qj,qk;
+  std::string names = "";
+
+  for (int i = 0; i < m_XMLparser->getBaseNodeCount("Survey"); i++) {
+      Survey2021 *buffer = new Survey2021();
+
+      m_XMLparser->getChildValue("Survey",i,"start_run",first_run);
+      m_XMLparser->getChildValue("Survey",i,"end_run",last_run);
+
+      //Discard entries for any channel that does not apply to our run
+      if(m_runNumber < first_run || m_runNumber > last_run) continue;
+
+      //If the entry applies, we store it in the vector
+      m_XMLparser->getChildValue("Survey",i,"detector",buffer->detector);
+      std::cout << "Parsing " << buffer->detector << " from run " << first_run << " to " << last_run << std::endl;
+      m_XMLparser->getChildValue("Survey",i,"x_pos",x);
+      m_XMLparser->getChildValue("Survey",i,"y_pos",y);
+      m_XMLparser->getChildValue("Survey",i,"z_pos",z);
+      buffer->pos.set(x,y,z);
+      m_XMLparser->getChildValue("Survey",i,"trigger_x",ttx);
+      m_XMLparser->getChildValue("Survey",i,"trigger_y",tty);
+      buffer->trig_table.set(ttx,tty,0.0);
+      m_XMLparser->getChildValue("Survey",i,"table_x",tx);
+      m_XMLparser->getChildValue("Survey",i,"table_y",ty);
+      buffer->table.set(tx,ty,0.0);
+      m_XMLparser->getChildValue("Survey",i,"w",qw);
+      m_XMLparser->getChildValue("Survey",i,"i",qi);
+      m_XMLparser->getChildValue("Survey",i,"j",qj);
+      m_XMLparser->getChildValue("Survey",i,"k",qk);
+      TQuaternion q(qw,qi,qj,qk);