Merge branch 'ntupledumper' into 'master'

update ntuple dumper

See merge request faser/calypso!339

PhysicsAnalysis/NtupleDumper/scripts/faser_ntuple_maker.py

@@ -44,8 +44,17 @@ parser.add_argument("--isMC", action='store_true',
 parser.add_argument("--partial", action='store_true',
                     help="Allow partial input files")
 
+parser.add_argument("--trigFilt", action='store_true',
+                    help="apply trigger event filter")
+parser.add_argument("--scintFilt", action='store_true',
+                    help="apply scintillator event filter")
+parser.add_argument("--NoTrackFilt", action='store_true',
+                    help="Don't apply track event filter (default: do)")
+
 parser.add_argument("--unblind", action='store_true',
                     help="Don't apply signal blinding (default: do)")
+parser.add_argument("--onlyblind", action='store_true',
+                    help="Only store events that were blinded (will override unblind arg)")
  
 parser.add_argument("--fluka", action='store_true',
                     help="Add FLUKA weights to ntuple")
@@ -138,8 +147,12 @@ if filepath.is_dir():
     print(f"Run   = {runstr}")
     print(f"First = {firstseg}")
     print(f"Last  = {lastseg}")
-    print(f"Args  = {args.tag}")
+    print(f"Tag  = {args.tag}")
+    print(f"Trigger Filter  = {args.trigFilt}")
+    print(f"Scintillator Filter  = {args.scintFilt}")
+    print(f"Track Filter = {not args.NoTrackFilt}")
     print(f"Blind = {not args.unblind}")
+    print(f"OnlyBlinded = {args.onlyblind}")
 
     # Find any tags
     tagstr = firststem.replace(f"{firstfaser}-{firstshort}-{runstr}-{firstseg}", "")
@@ -243,7 +256,7 @@ if args.isMC:
         acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile))
 
 else:
-    acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile, DoBlinding=(not args.unblind)))
+    acc.merge(NtupleDumperAlgCfg(ConfigFlags, outfile, DoBlinding=(not args.unblind), OnlyBlinded=args.onlyblind, DoScintFilter = args.scintFilt, DoTrackFilter = (not args.NoTrackFilt), DoTrigFilter = args.trigFilt))
 
 if not args.verbose:
     from AthenaConfiguration.ComponentFactory import CompFactory

PhysicsAnalysis/NtupleDumper/src/NtupleDumperAlg.h

@@ -52,12 +52,13 @@ private:
 
   bool waveformHitOK(const xAOD::WaveformHit* hit) const;
   void clearTree() const;
-  void setNaN() const;
+  void clearTrackTruth() const;
   void addBranch(const std::string &name,float* var);
   void addBranch(const std::string &name,unsigned int* var);
   void addWaveBranches(const std::string &name, int nchannels, int first);
   void FillWaveBranches(const xAOD::WaveformHitContainer &wave) const;
   void addCalibratedBranches(const std::string &name, int nchannels, int first);
+  double radius(const Acts::Vector3 &position) const;
 
   ServiceHandle <ITHistSvc> m_histSvc;
 
@@ -101,17 +102,23 @@ private:
   const PreshowerID* m_preshowerHelper;
   const EcalID*      m_ecalHelper;
 
-  BooleanProperty m_useIFT            { this, "UseIFT", false, "Use IFT tracks" };
-  BooleanProperty m_doCalib           { this, "DoCalib", true, "Fill calibrated calorimeter quantities" };
-  BooleanProperty m_doBlinding        { this, "DoBlinding", true, "Blinding will not output events with Calo signal > 10 GeV e-" };
-  BooleanProperty m_useFlukaWeights   { this, "UseFlukaWeights", false, "Flag to weight events according to value stored in HepMC::GenEvent" };
-  BooleanProperty m_useGenieWeights   { this, "UseGenieWeights", false, "Flag to weight events according to Genie luminosity" };
-  IntegerProperty m_flukaCollisions   { this, "FlukaCollisions", 137130000, "Number of proton-proton collisions in FLUKA sample." };
-  DoubleProperty  m_flukaCrossSection { this, "FlukaCrossSection", 80.0, "Fluka p-p inelastic cross-section in millibarns." };
-  DoubleProperty  m_genieLuminosity   { this, "GenieLuminosity", 150.0, "Genie luminosity in inverse fb." };
-  DoubleProperty  m_minMomentum       { this, "MinMomentum", 50000.0, "Write out all truth particles with a momentum larger MinMomentum"};
+  BooleanProperty m_useIFT               { this, "UseIFT", false, "Use IFT tracks" };
+  BooleanProperty m_doCalib              { this, "DoCalib", true, "Fill calibrated calorimeter quantities" };
+  BooleanProperty m_doBlinding           { this, "DoBlinding", true, "Blinding will not output events with no veto signal adn Calo signal > 20 GeV e-" };
+  BooleanProperty m_onlyBlinded          { this, "OnlyBlinded", false, "Only events that would be blinded are saved" };
+  BooleanProperty m_doTrigFilter         { this, "DoTrigFilter", false, "Only events that pass trigger cuts are passed" };
+  BooleanProperty m_doScintFilter        { this, "DoScintFilter", false, "Only events that pass scintillator coincidence cuts are passed" };
+  BooleanProperty m_doTrackFilter        { this, "DoTrackFilter", true, "Only events that have >= 1 long track are passed, also non-colliding events with a track or calo signal are passed" };
+
+  BooleanProperty m_useFlukaWeights      { this, "UseFlukaWeights", false, "Flag to weight events according to value stored in HepMC::GenEvent" };
+  BooleanProperty m_useGenieWeights      { this, "UseGenieWeights", false, "Flag to weight events according to Genie luminosity" };
+  IntegerProperty m_flukaCollisions      { this, "FlukaCollisions", 137130000, "Number of proton-proton collisions in FLUKA sample." };
+  DoubleProperty  m_flukaCrossSection    { this, "FlukaCrossSection", 80.0, "Fluka p-p inelastic cross-section in millibarns." };
+  DoubleProperty  m_genieLuminosity      { this, "GenieLuminosity", 150.0, "Genie luminosity in inverse fb." };
+  DoubleProperty  m_minMomentum          { this, "MinMomentum", 50000.0, "Write out all truth particles with a momentum larger MinMomentum"};
 
   DoubleProperty m_blindingCaloThreshold {this, "BlindingCaloThreshold", 25000.0, "Blind events with Ecal energy above threshold (in MeV)"}; 
+  DoubleProperty m_caloMC_FudgeFactor    {this, "CaloFudgeFactorMC", 1.088, "Correct the MC energy calibration by this fudge factor"};
 
   double m_baseEventCrossSection {1.0};
   const double kfemtoBarnsPerMilliBarn {1.0e12};
@@ -150,7 +157,9 @@ private:
   mutable float m_wave_baseline_mean[15];
   mutable float m_wave_baseline_rms[15];
   mutable unsigned int m_wave_status[15];
-  
+ 
+  mutable unsigned int m_scintHit;
+ 
   mutable float m_calibrated_nMIP[15];
   mutable float m_calibrated_E_dep[15];
   mutable float m_calibrated_E_EM[15];
@@ -159,6 +168,8 @@ private:
   mutable float m_calo_total_E_dep;
   mutable float m_calo_total_E_EM;
 
+  mutable float m_calo_total_E_EM_fudged;
+
   mutable float m_preshower_total_nMIP;
   mutable float m_preshower_total_E_dep;
 
@@ -170,9 +181,6 @@ private:
   mutable float m_Preshower12_Edep;
   mutable float m_Preshower13_Edep;
 
-  mutable float m_MIP_sim_Edep_calo;
-  mutable float m_MIP_sim_Edep_preshower;
-
   mutable float m_clock_phase;
 
   mutable unsigned int m_station0Clusters;
@@ -247,6 +255,10 @@ private:
   mutable std::vector<double> m_yCalo;
   mutable std::vector<double> m_thetaxCalo;
   mutable std::vector<double> m_thetayCalo;
+  mutable std::vector<double> m_x_atMaxRadius;
+  mutable std::vector<double> m_y_atMaxRadius;
+  mutable std::vector<double> m_z_atMaxRadius;
+  mutable std::vector<double> m_r_atMaxRadius;
 
   mutable std::vector<int> m_t_pdg; // pdg code of the truth matched particle
   mutable std::vector<int> m_t_barcode; // barcode of the truth matched particle
@@ -324,9 +336,6 @@ private:
 
   mutable std::vector<int> m_truth_pdg; // pdg of first 10 truth particles 
 
-
-
-
   mutable double m_truthLeptonMomentum;
   mutable int    m_truthBarcode;
   mutable int    m_truthPdg;