Changeset 1361 in svn

Timestamp:

Apr 8, 2014, 6:43:52 PM (11 years ago)

Author:

Pavel Demin

Message:

add HectorHit class and corresponding ROOT tree branch

Location:

trunk

Files:

• classes/ClassesLinkDef.h (modified) (2 diffs)
• classes/DelphesClasses.cc (modified) (1 diff)
• classes/DelphesClasses.h (modified) (1 diff)
• examples/LHCB1IR5_5TeV.tfs (added)
• examples/LHCB1IR5_7TeV.tfs (added)
• examples/LHCB2IR5_5TeV.tfs (added)
• examples/LHCB2IR5_7TeV.tfs (added)
• examples/delphes_card_CMS_PileUp_Hector.tcl (added)
• modules/Hector.cc (modified) (5 diffs)
• modules/Hector.h (modified) (1 diff)
• modules/TreeWriter.cc (modified) (22 diffs)
• modules/TreeWriter.h (modified) (1 diff)

trunk/classes/ClassesLinkDef.h

@@ -7 +7 @@
  *  $Revision: 1.2 $
  *
- *  
+ *
  *  \author P. Demin - UCL, Louvain-la-Neuve
  *
@@ -45 +45 @@
 #pragma link C++ class Track+;
 #pragma link C++ class Tower+;
+#pragma link C++ class HectorHit+;
 
 #pragma link C++ class Candidate+;

trunk/classes/DelphesClasses.cc

@@ -26 +26 @@
 CompBase *Track::fgCompare = CompPT<Track>::Instance();
 CompBase *Tower::fgCompare = CompE<Tower>::Instance();
+CompBase *HectorHit::fgCompare = CompE<HectorHit>::Instance();
 CompBase *Candidate::fgCompare = CompMomentumPt<Candidate>::Instance();
 

trunk/classes/DelphesClasses.h

@@ -390 +390 @@
 //---------------------------------------------------------------------------
 
+class HectorHit: public SortableObject
+{
+public:
+  Float_t E; // reconstructed energy [GeV]
+
+  Float_t Tx; // angle of the momentum in the horizontal (x,z) plane [urad]
+  Float_t Ty; // angle of the momentum in the verical (y,z) plane [urad]
+
+  Float_t T; // time of flight to the detector [s]
+
+  Float_t X; // horizontal distance to the beam [um]
+  Float_t Y; // vertical distance to the beam [um]
+  Float_t S; // distance to the interaction point [m]
+
+  static CompBase *fgCompare; //!
+  const CompBase *GetCompare() const { return fgCompare; }
+
+  ClassDef(HectorHit, 1)
+};
+
+//---------------------------------------------------------------------------
+
 class Candidate: public SortableObject
 {

trunk/modules/Hector.cc

@@ -63 +63 @@
   fBeamLineLength = GetDouble("BeamLineLength", 430.0);
   fDistance = GetDouble("Distance", 420.0);
+  fOffsetX = GetDouble("OffsetX", 0.0);
+  fOffsetS = GetDouble("OffsetS", 120.0);
   fSigmaE = GetDouble("SigmaE", 0.0);
   fSigmaX = GetDouble("SigmaX", 0.0);
   fSigmaY = GetDouble("SigmaY", 0.0);
+  fSigmaT = GetDouble("SigmaT", 0.0);
   fEtaMin = GetDouble("EtaMin", 5.0);
 
   fBeamLine = new H_BeamLine(fDirection, fBeamLineLength + 0.1);
-  fBeamLine->fill(GetString("BeamLineFile", "examples/LHCB1IR5_5TeV.tfs"), fDirection, "IP5" );
+  fBeamLine->fill(GetString("BeamLineFile", "examples/LHCB1IR5_5TeV.tfs"), fDirection, GetString("IPName", "IP5"));
   fBeamLine->offsetElements(120, -0.097*fDirection);
+  fBeamLine->offsetElements(fOffsetS, fOffsetX);
   fBeamLine->calcMatrix();
 
@@ -80 +84 @@
   // create output array
 
-  fOutputArray = ExportArray(GetString("OutputArray", "stableParticles"));
+  fOutputArray = ExportArray(GetString("OutputArray", "hits"));
 }
 
@@ -97 +101 @@
   Candidate *candidate, *mother;
   Double_t pz;
-  Double_t x, y, z, tx, ty;
+  Double_t x, y, z, tx, ty, theta;
+  Double_t distance, time;
+
+  const Double_t c_light = 2.99792458E8;
 
   fItInputArray->Reset();
@@ -118 +125 @@
     ty = 0.0;
 
+    theta = TMath::Hypot(TMath::ATan(candidateMomentum.Px()/pz), TMath::ATan(candidateMomentum.Py()/pz));
+    distance = (fDistance - 1.0E-3 * candidatePosition.Z())/TMath::Cos(theta);
+    time = gRandom->Gaus((distance + 1.0E-3 * candidatePosition.T())/c_light, fSigmaT);
+
     H_BeamParticle particle(candidate->Mass, candidate->Charge);
     particle.set4Momentum(candidateMomentum);
@@ -133 +144 @@
     mother = candidate;
     candidate = static_cast<Candidate*>(candidate->Clone());
-    candidate->Position.SetXYZT(particle.getX(), particle.getY(), fDistance, 0.0);
+    candidate->Position.SetXYZT(particle.getX(), particle.getY(), particle.getS(), time);
     candidate->Momentum.SetPxPyPzE(particle.getTX(), particle.getTY(), 0.0, particle.getE());
     candidate->AddCandidate(mother);

trunk/modules/Hector.h

@@ -34 +34 @@
 
   Int_t fDirection;
-  
+
   Double_t fBeamLineLength, fDistance;
-  Double_t fSigmaE, fSigmaX, fSigmaY;
+  Double_t fOffsetX, fOffsetS;
+  Double_t fSigmaE, fSigmaX, fSigmaY, fSigmaT;
   Double_t fEtaMin;
-  
+
   H_BeamLine *fBeamLine;
 

trunk/modules/TreeWriter.cc

@@ -67 +67 @@
   fClassMap[Rho::Class()] = &TreeWriter::ProcessRho;
   fClassMap[Weight::Class()] = &TreeWriter::ProcessWeight;
+  fClassMap[HectorHit::Class()] = &TreeWriter::ProcessHectorHit;
 
   TBranchMap::iterator itBranchMap;
@@ -95 +96 @@
       continue;
     }
-    
+
     itClassMap = fClassMap.find(branchClass);
     if(itClassMap == fClassMap.end())
@@ -160 +161 @@
   GenParticle *entry = 0;
   Double_t pt, signPz, cosTheta, eta, rapidity;
-  
-  const Double_t c_light = 2.99792458E8;
-  
+
+  const Double_t c_light = 2.99792458E8;
+
   // loop over all particles
   iterator.Reset();
@@ -220 +221 @@
   Candidate *candidate = 0;
   Vertex *entry = 0;
-  
+
   const Double_t c_light = 2.99792458E8;
 
@@ -248 +249 @@
   Double_t pt, signz, cosTheta, eta, rapidity;
   const Double_t c_light = 2.99792458E8;
-  
+
   // loop over all tracks
   iterator.Reset();
@@ -310 +311 @@
   Double_t pt, signPz, cosTheta, eta, rapidity;
   const Double_t c_light = 2.99792458E8;
-  
+
   // loop over all towers
   iterator.Reset();
@@ -317 +318 @@
     const TLorentzVector &momentum = candidate->Momentum;
     const TLorentzVector &position = candidate->Position;
-    
+
     pt = momentum.Pt();
     cosTheta = TMath::Abs(momentum.CosTheta());
@@ -339 +340 @@
     entry->Edges[2] = candidate->Edges[2];
     entry->Edges[3] = candidate->Edges[3];
-    
+
     entry->T = position.T()*1.0E-3/c_light;
-    
+
     FillParticles(candidate, &entry->Particles);
   }
@@ -355 +356 @@
   Double_t pt, signPz, cosTheta, eta, rapidity;
   const Double_t c_light = 2.99792458E8;
-  
+
   array->Sort();
 
@@ -365 +366 @@
     const TLorentzVector &momentum = candidate->Momentum;
     const TLorentzVector &position = candidate->Position;
-    
+
 
     pt = momentum.Pt();
@@ -379 +380 @@
     entry->PT = pt;
     entry->E = momentum.E();
-    
+
     entry->T = position.T()*1.0E-3/c_light;
-   
+
     entry->EhadOverEem = candidate->Eem > 0.0 ? candidate->Ehad/candidate->Eem : 999.9;
 
@@ -397 +398 @@
   Double_t pt, signPz, cosTheta, eta, rapidity;
   const Double_t c_light = 2.99792458E8;
-  
+
   array->Sort();
 
@@ -406 +407 @@
     const TLorentzVector &momentum = candidate->Momentum;
     const TLorentzVector &position = candidate->Position;
-    
+
     pt = momentum.Pt();
     cosTheta = TMath::Abs(momentum.CosTheta());
@@ -418 +419 @@
     entry->Phi = momentum.Phi();
     entry->PT = pt;
-    
+
     entry->T = position.T()*1.0E-3/c_light;
-    
+
     entry->Charge = candidate->Charge;
 
@@ -437 +438 @@
   Muon *entry = 0;
   Double_t pt, signPz, cosTheta, eta, rapidity;
-  
-  const Double_t c_light = 2.99792458E8;
-  
+
+  const Double_t c_light = 2.99792458E8;
+
   array->Sort();
 
@@ -448 +449 @@
     const TLorentzVector &momentum = candidate->Momentum;
     const TLorentzVector &position = candidate->Position;
-    
+
 
     pt = momentum.Pt();
@@ -466 +467 @@
 
     entry->T = position.T()*1.0E-3/c_light;
-    
+
    // cout<<entry->PT<<","<<entry->T<<endl;
 
@@ -485 +486 @@
   Double_t ecalEnergy, hcalEnergy;
   const Double_t c_light = 2.99792458E8;
-  
+
   array->Sort();
 
@@ -493 +494 @@
   {
     TIter itConstituents(candidate->GetCandidates());
-    
-    const TLorentzVector &momentum = candidate->Momentum;
-    const TLorentzVector &position = candidate->Position;
-   
+
+    const TLorentzVector &momentum = candidate->Momentum;
+    const TLorentzVector &position = candidate->Position;
+
     pt = momentum.Pt();
     cosTheta = TMath::Abs(momentum.CosTheta());
@@ -510 +511 @@
 
     entry->T = position.T()*1.0E-3/c_light;
-    
+
     entry->Mass = momentum.M();
 
@@ -533 +534 @@
 
     entry->EhadOverEem = ecalEnergy > 0.0 ? hcalEnergy/ecalEnergy : 999.9;
-  
+
     //---   Pile-Up Jet ID variables ----
 
@@ -632 +633 @@
 //------------------------------------------------------------------------------
 
+void TreeWriter::ProcessHectorHit(ExRootTreeBranch *branch, TObjArray *array)
+{
+  TIter iterator(array);
+  Candidate *candidate = 0;
+  HectorHit *entry = 0;
+
+  // loop over all roman pot hits
+  iterator.Reset();
+  while((candidate = static_cast<Candidate*>(iterator.Next())))
+  {
+    const TLorentzVector &position = candidate->Position;
+    const TLorentzVector &momentum = candidate->Momentum;
+
+    entry = static_cast<HectorHit*>(branch->NewEntry());
+
+    entry->E = momentum.E();
+
+    entry->Tx = momentum.Px();
+    entry->Ty = momentum.Py();
+
+    entry->T = position.T();
+
+    entry->X = position.X();
+    entry->Y = position.Y();
+    entry->S = position.Z();
+  }
+}
+
+//------------------------------------------------------------------------------
+
 void TreeWriter::Process()
 {

trunk/modules/TreeWriter.h

@@ -53 +53 @@
   void ProcessRho(ExRootTreeBranch *branch, TObjArray *array);
   void ProcessWeight(ExRootTreeBranch *branch, TObjArray *array);
+  void ProcessHectorHit(ExRootTreeBranch *branch, TObjArray *array);
 
 #ifndef __CINT__