Changeset 61dccd3 in git

Timestamp:

Aug 24, 2021, 11:57:44 AM (3 years ago)

Author:

michele <michele.selvaggi@…>

Branches:

master

Children:

83e77ee

Parents:

1ad8eca

Message:

added charged deposit to calo, fix time calc in DR, add path length to tower

Files:

• classes/DelphesClasses.cc (modified) (3 diffs)
• classes/DelphesClasses.h (modified) (5 diffs)
• modules/DualReadoutCalorimeter.cc (modified) (9 diffs)
• modules/DualReadoutCalorimeter.h (modified) (1 diff)
• modules/SimpleCalorimeter.cc (modified) (2 diffs)
• modules/TreeWriter.cc (modified) (2 diffs)

classes/DelphesClasses.cc

@@ -213 +213 @@
   Flavor(0), FlavorAlgo(0), FlavorPhys(0),
   BTag(0), BTagAlgo(0), BTagPhys(0),
-  TauTag(0), TauWeight(0.0), Eem(0.0), Ehad(0.0),
+  TauTag(0), TauWeight(0.0), Eem(0.0), Ehad(0.0), Etrk(0.0),
   DeltaEta(0.0), DeltaPhi(0.0),
   Momentum(0.0, 0.0, 0.0, 0.0),
@@ -378 +378 @@
   object.Eem = Eem;
   object.Ehad = Ehad;
+  object.Etrk = Etrk;
   object.Edges[0] = Edges[0];
   object.Edges[1] = Edges[1];
@@ -512 +513 @@
   Eem = 0.0;
   Ehad = 0.0;
+  Etrk = 0.0;
   Edges[0] = 0.0;
   Edges[1] = 0.0;

classes/DelphesClasses.h

@@ -515 +515 @@
   Float_t Eem; // calorimeter tower electromagnetic energy
   Float_t Ehad; // calorimeter tower hadronic energy
+  Float_t Etrk; // total charged energy hitting tower
 
   Float_t Edges[4]; // calorimeter tower edges
@@ -525 +526 @@
   TLorentzVector P4() const;
 
-  ClassDef(Tower, 2)
+  ClassDef(Tower, 3)
 };
 
@@ -604 +605 @@
   Float_t Eem; // calorimeter tower electromagnetic energy
   Float_t Ehad; // calorimeter tower hadronic energy
+  Float_t Etrk; // total charged energy hitting tower
 
   Float_t Edges[4]; // calorimeter tower edges
@@ -609 +611 @@
   TRefArray Particles; // references to generated particles
 
-  ClassDef(ParticleFlowCandidate, 2)
+  ClassDef(ParticleFlowCandidate, 3)
 
 };
@@ -674 +676 @@
   Float_t Eem;
   Float_t Ehad;
+  Float_t Etrk;
 
   Float_t Edges[4];

modules/DualReadoutCalorimeter.cc

@@ -248 +248 @@
   number = -1;
   fTowerRmax=0.;
+
+  //cout<<"--------- new event ---------- "<<endl;
+
   while((particle = static_cast<Candidate*>(fItParticleInputArray->Next())))
   {
@@ -391 +394 @@
       fTowerTrackHits = 0;
       fTowerPhotonHits = 0;
+
+      fTowerTime = 0.0;
+      fTowerTimeWeight = 0.0;
 
       fECalTowerTrackArray->Clear();
@@ -455 +461 @@
     position = particle->Position;
 
+
     // fill current tower
     ecalEnergy = momentum.E() * fECalTowerFractions[number];
@@ -462 +469 @@
     fHCalTowerEnergy += hcalEnergy;
 
-    if(ecalEnergy > fTimingEnergyMin && fTower)
-    {
-      if (abs(particle->PID) != 11 || !fElectronsFromTrack)
-      {
-        fTower->ECalEnergyTimePairs.push_back(make_pair<Float_t, Float_t>(ecalEnergy, particle->Position.T()));
-      }
-    }
+    // assume combined timing measurements in ECAL/HCAL sections
+    fTowerTime += (ecalEnergy + hcalEnergy) * position.T(); //sigma_t ~ 1/sqrt(E)
+    fTowerTimeWeight += ecalEnergy + hcalEnergy;
+    //fTowerTime += (hcalEnergy) * position.T(); //sigma_t ~ 1/sqrt(E)
+    //fTowerTimeWeight += hcalEnergy;
+    //fTowerTime +=  position.T(); //sigma_t ~ 1/sqrt(E)
+    //fTowerTimeWeight += 1;
+
+    //cout<<" tower particle PID, pt, eta, phi, l, tof:  "<<particle->PID<<", "<<momentum.E()<<", "<<momentum.Eta()<<", "<<momentum.Phi()<<", "<<position.Vect().Mag()<<", "<<position.T()/2.99792458E2<<endl;
+    //cout<<" tower particle time, weight:  "<<fTowerTime/2.99792458E2<<", "<<fTowerTimeWeight<<endl;
 
     fTower->AddCandidate(particle);
@@ -484 +494 @@
 
   Candidate *track, *tower, *mother;
-  Double_t energy, pt, eta, phi, r;
+  Double_t energy, pt, eta, phi, r, time;
   Double_t ecalEnergy, hcalEnergy;
   Double_t ecalNeutralEnergy, hcalNeutralEnergy, neutralEnergy;
@@ -529 +539 @@
   hcalEnergy = LogNormal(fHCalTowerEnergy, hcalSigma);
 
+  time = (fTowerTimeWeight < 1.0E-09) ? 0.0 : fTowerTime / fTowerTimeWeight;
+
   ecalSigma = fECalResolutionFormula->Eval(0.0, fTowerEta, 0.0, ecalEnergy);
   hcalSigma = fHCalResolutionFormula->Eval(0.0, fTowerEta, 0.0, hcalEnergy);
@@ -535 +547 @@
   if(hcalEnergy < fHCalEnergyMin || hcalEnergy < fHCalEnergySignificanceMin*hcalSigma) hcalEnergy = 0.0;
 
-  //cout<<"Measured energy: "<<energy<<endl;
-
   if(fSmearTowerCenter)
   {
@@ -550 +560 @@
   pt = energy / TMath::CosH(eta);
 
-  // Time calculation for tower
-  fTower->NTimeHits = 0;
-  sumWeightedTime = 0.0;
-  sumWeight = 0.0;
-
-  for(size_t i = 0; i < fTower->ECalEnergyTimePairs.size(); ++i)
-  {
-    weight = TMath::Power((fTower->ECalEnergyTimePairs[i].first),2);
-    sumWeightedTime += weight * fTower->ECalEnergyTimePairs[i].second;
-    sumWeight += weight;
-    fTower->NTimeHits++;
-  }
-
   // check whether barrel or endcap tower
-  if (fTower->Position.Perp() < fTowerRmax && TMath::Abs(eta) > 0.)
+  if ((fTowerRmax - fTower->Position.Perp()) < 1.e-06 && TMath::Abs(eta) > 0.)
     r = fTower->Position.Z()/TMath::SinH(eta);
   else
     r = fTower->Position.Pt();
 
-  if(sumWeight > 0.0)
-  {
-    fTower->Position.SetPtEtaPhiE(r, eta, phi, sumWeightedTime/sumWeight);
-  }
-  else
-  {
-    fTower->Position.SetPtEtaPhiE(r, eta, phi, 999999.9);
-  }
+  fTower->Position.SetPtEtaPhiE(r, eta, phi, time);
+
+  fTower->Momentum.SetPtEtaPhiE(pt, eta, phi, energy);
+  fTower->L = fTower->Position.Vect().Mag();
+  //cout<<"   tower pt, eta, phi, l, tof:  "<<fTower->Momentum.E()<<", "<<fTower->Momentum.Eta()<<", "<<fTower->Momentum.Phi()<<", "<<fTower->L<<", "<<fTower->Position.T()/2.99792458E2<<endl;
 
   fTower->Momentum.SetPtEtaPhiE(pt, eta, phi, energy);
   fTower->Eem = ecalEnergy;
   fTower->Ehad = hcalEnergy;
-
+  fTower->Etrk = fTrackEnergy;
   fTower->Edges[0] = fTowerEdges[0];
   fTower->Edges[1] = fTowerEdges[1];
@@ -602 +596 @@
   neutralSigma = neutralEnergy / TMath::Sqrt(fTrackSigma*fTrackSigma + sigma*sigma);
 
-  //cout<<"trackEnergy: "<<fTrackEnergy<<", trackSigma: "<<fTrackSigma<<", Ntracks: "<<fTowerTrackArray->GetEntries()<<endl;
-
-  //cout<<"neutralEnergy: "<<neutralEnergy<<", neutralSigma: "<<neutralSigma<<", :fEnergyMin "<<fEnergyMin<<", fEnergySignificanceMin: "<<fEnergySignificanceMin<<endl;
-
-   // For now, if neutral excess is significant, simply create neutral EflowPhoton tower and clone each track into eflowtrack !!! -> Creating only photons !! EFlowNeutralHadron collection will be empy!!! TO BE FIXED
   if(neutralEnergy > fEnergyMin && neutralSigma > fEnergySignificanceMin)
   {
-
-    //cout<<"significant neutral excess found:"<<endl;
     // create new photon tower
     tower = static_cast<Candidate*>(fTower->Clone());
     pt = neutralEnergy / TMath::CosH(eta);
-    //cout<<"Creating tower with Pt, Eta, Phi, Energy: "<<pt<<","<<eta<<","<<phi<<","<<neutralEnergy<<endl;
     tower->Momentum.SetPtEtaPhiE(pt, eta, phi, neutralEnergy);
 

modules/DualReadoutCalorimeter.h

@@ -79 +79 @@
   Double_t fTrackSigma;
 
+  Double_t fTowerTime;
+  Double_t fTowerTimeWeight;
+
   Bool_t fSmearTowerCenter;
 

modules/SimpleCalorimeter.cc

@@ -451 +451 @@
 
   // check whether barrel or endcap tower
-  if (fTower->Position.Perp() < fTowerRmax && TMath::Abs(eta) > 0.)
+  if ((fTowerRmax - fTower->Position.Perp()) < 1.e-06 && TMath::Abs(eta) > 0.)
     r = fTower->Position.Z()/TMath::SinH(eta);
   else
@@ -458 +458 @@
   fTower->Position.SetPtEtaPhiE(r, eta, phi, time);
   fTower->Momentum.SetPtEtaPhiE(pt, eta, phi, energy);
+  fTower->L = fTower->Position.Vect().Mag();
 
   fTower->Eem = (!fIsEcal) ? 0 : energy;
   fTower->Ehad = (fIsEcal) ? 0 : energy;
+  fTower->Etrk = fTrackEnergy;
 
   fTower->Edges[0] = fTowerEdges[0];

modules/TreeWriter.cc

@@ -455 +455 @@
     entry->Eem = candidate->Eem;
     entry->Ehad = candidate->Ehad;
+    entry->Etrk = candidate->Etrk;
     entry->Edges[0] = candidate->Edges[0];
     entry->Edges[1] = candidate->Edges[1];
@@ -575 +576 @@
     entry->Eem = candidate->Eem;
     entry->Ehad = candidate->Ehad;
+    entry->Etrk = candidate->Etrk;
     entry->Edges[0] = candidate->Edges[0];
     entry->Edges[1] = candidate->Edges[1];