Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Snowmass version is different from zhenbinwu/CMSupgrade. #1

Open
misho104 opened this issue May 21, 2014 · 0 comments
Open

Snowmass version is different from zhenbinwu/CMSupgrade. #1

misho104 opened this issue May 21, 2014 · 0 comments
Labels

Comments

@misho104
Copy link
Contributor

I expected that the snowmass_version @09d0552 is identical to Delphes-3.0.10 in zhenbinwu/CMSUpgrade@f6564e4, but there are slight differences as follows.
Why inconsistent? Which are the true version? (or irrelevant?)

diff -ur CMSUpgrade/Delphes-3.0.10/modules/Calorimeter.cc Delphes3/modules/Calorimeter.cc
--- CMSUpgrade/Delphes-3.0.10/modules/Calorimeter.cc    2014-05-21 18:37:15.941871251 +0900
+++ Delphes3/modules/Calorimeter.cc 2013-09-13 21:54:12.000000000 +0900
@@ -4,8 +4,8 @@
  *  Fills calorimeter towers, performs calorimeter resolution smearing,
  *  preselects towers hit by photons and creates energy flow objects.
  *
- *  $Date: 2013-06-25 16:24:34 +0200 (Tue, 25 Jun 2013) $
- *  $Revision: 1145 $
+ *  $Date: 2013-09-04 17:20:22 +0200 (Wed, 04 Sep 2013) $
+ *  $Revision: 1280 $
  *
  *
  *  \author P. Demin - UCL, Louvain-la-Neuve
@@ -42,15 +42,13 @@
 Calorimeter::Calorimeter() :
   fECalResolutionFormula(0), fHCalResolutionFormula(0),
   fItParticleInputArray(0), fItTrackInputArray(0),
-  fTowerTrackArray(0), fItTowerTrackArray(0),
-  fTowerPhotonArray(0), fItTowerPhotonArray(0)
+  fTowerTrackArray(0), fItTowerTrackArray(0)
 {
   fECalResolutionFormula = new DelphesFormula;
   fHCalResolutionFormula = new DelphesFormula;
+
   fTowerTrackArray = new TObjArray;
   fItTowerTrackArray = fTowerTrackArray->MakeIterator();
-  fTowerPhotonArray = new TObjArray;
-  fItTowerPhotonArray = fTowerPhotonArray->MakeIterator();
 }

 //------------------------------------------------------------------------------
@@ -59,10 +57,9 @@
 {
   if(fECalResolutionFormula) delete fECalResolutionFormula;
   if(fHCalResolutionFormula) delete fHCalResolutionFormula;
+
   if(fTowerTrackArray) delete fTowerTrackArray;
   if(fItTowerTrackArray) delete fItTowerTrackArray;
-  if(fTowerPhotonArray) delete fTowerPhotonArray;
-  if(fItTowerPhotonArray) delete fItTowerPhotonArray;
 }

 //------------------------------------------------------------------------------
@@ -160,7 +157,7 @@

 void Calorimeter::Finish()
 {
-  vector< vector< Double_t>* >::iterator itPhiBin;
+  vector< vector< Double_t >* >::iterator itPhiBin;
   if(fItParticleInputArray) delete fItParticleInputArray;
   if(fItTrackInputArray) delete fItTrackInputArray;
   for(itPhiBin = fPhiBins.begin(); itPhiBin != fPhiBins.end(); ++itPhiBin)
@@ -192,8 +189,10 @@

   DelphesFactory *factory = GetFactory();
   fTowerHits.clear();
-  fECalFractions.clear();
-  fHCalFractions.clear();
+  fTowerECalFractions.clear();
+  fTowerHCalFractions.clear();
+  fTrackECalFractions.clear();
+  fTrackHCalFractions.clear();

   // loop over all particles
   fItParticleInputArray->Reset();
@@ -214,8 +213,8 @@
     ecalFraction = itFractionMap->second.first;
     hcalFraction = itFractionMap->second.second;

-    fECalFractions.push_back(ecalFraction);
-    fHCalFractions.push_back(hcalFraction);
+    fTowerECalFractions.push_back(ecalFraction);
+    fTowerHCalFractions.push_back(hcalFraction);

     if(ecalFraction < 1.0E-9 && hcalFraction < 1.0E-9) continue;

@@ -232,9 +231,8 @@
     if(itPhiBin == phiBins->begin() || itPhiBin == phiBins->end()) continue;
     phiBin = distance(phiBins->begin(), itPhiBin);

-    flags = (particle->Charge == 0);
-    flags |= (pdgCode == 22) << 1;
-    flags |= (pdgCode == 11) << 2;
+    flags = 0;
+    flags |= (pdgCode == 11 || pdgCode == 22) << 1;

     // make tower hit {16-bits for eta bin number, 16-bits for phi bin number, 8-bits for flags, 24-bits for particle number}
     towerHit = (Long64_t(etaBin) << 48) | (Long64_t(phiBin) << 32) | (Long64_t(flags) << 24) | Long64_t(number);
@@ -250,6 +248,20 @@
     const TLorentzVector &trackPosition = track->Position;
     ++number;

+    pdgCode = TMath::Abs(track->PID);
+
+    itFractionMap = fFractionMap.find(pdgCode);
+    if(itFractionMap == fFractionMap.end())
+    {
+      itFractionMap = fFractionMap.find(0);
+    }
+
+    ecalFraction = itFractionMap->second.first;
+    hcalFraction = itFractionMap->second.second;
+
+    fTrackECalFractions.push_back(ecalFraction);
+    fTrackHCalFractions.push_back(hcalFraction);
+
     // find eta bin [1, fEtaBins.size - 1]
     itEtaBin = lower_bound(fEtaBins.begin(), fEtaBins.end(), trackPosition.Eta());
     if(itEtaBin == fEtaBins.begin() || itEtaBin == fEtaBins.end()) continue;
@@ -263,8 +275,10 @@
     if(itPhiBin == phiBins->begin() || itPhiBin == phiBins->end()) continue;
     phiBin = distance(phiBins->begin(), itPhiBin);

+    flags = 1;
+
     // make tower hit {16-bits for eta bin number, 16-bits for phi bin number, 8-bits for flags, 24-bits for track number}
-    towerHit = (Long64_t(etaBin) << 48) | (Long64_t(phiBin) << 32) | (Long64_t(1) << 27) | Long64_t(number);
+    towerHit = (Long64_t(etaBin) << 48) | (Long64_t(phiBin) << 32) | (Long64_t(flags) << 24) | Long64_t(number);

     fTowerHits.push_back(towerHit);
   }
@@ -312,53 +326,48 @@
       fTowerECalEnergy = 0.0;
       fTowerHCalEnergy = 0.0;

-      fTowerECalNeutralEnergy = 0.0;
-      fTowerHCalNeutralEnergy = 0.0;
+      fTrackECalEnergy = 0.0;
+      fTrackHCalEnergy = 0.0;

-      fTowerNeutralHits = 0;
-      fTowerPhotonHits = 0;
-      fTowerElectronHits = 0;
       fTowerTrackHits = 0;
-      fTowerAllHits = 0;
+      fTowerPhotonHits = 0;

       fTowerTrackArray->Clear();
-      fTowerPhotonArray->Clear();
     }

     // check for track hits
-    if(flags & 8)
+    if(flags & 1)
     {
       ++fTowerTrackHits;
+
       track = static_cast<Candidate*>(fTrackInputArray->At(number));
+      momentum = track->Momentum;
+
+      ecalEnergy = momentum.E() * fTrackECalFractions[number];
+      hcalEnergy = momentum.E() * fTrackHCalFractions[number];
+
+      fTrackECalEnergy += ecalEnergy;
+      fTrackHCalEnergy += hcalEnergy;
+
       fTowerTrackArray->Add(track);
+
       continue;
     }

+    // check for photon and electron hits in current tower
+    if(flags & 2) ++fTowerPhotonHits;
+
     particle = static_cast<Candidate*>(fParticleInputArray->At(number));
     momentum = particle->Momentum;

     // fill current tower
-    ecalEnergy = momentum.E() * fECalFractions[number];
-    hcalEnergy = momentum.E() * fHCalFractions[number];
+    ecalEnergy = momentum.E() * fTowerECalFractions[number];
+    hcalEnergy = momentum.E() * fTowerHCalFractions[number];

     fTowerECalEnergy += ecalEnergy;
     fTowerHCalEnergy += hcalEnergy;

-    ++fTowerAllHits;
     fTower->AddCandidate(particle);
-
-    // check for neutral hits in current tower
-    if(flags & 1) ++fTowerNeutralHits;
-
-    // check for photon hits in current tower
-    if(flags & 2)
-    {
-      ++fTowerPhotonHits;
-      fTowerPhotonArray->Add(particle);
-    }
-
-    // check for electron hits in current tower
-    if(flags & 4) ++fTowerElectronHits;
   }

   // finalize last tower
@@ -369,21 +378,26 @@

 void Calorimeter::FinalizeTower()
 {
-  Candidate *particle, *track, *tower;
+  Candidate *track, *tower;
   Double_t energy, pt, eta, phi;
   Double_t ecalEnergy, hcalEnergy;
+  Double_t ecalSigma, hcalSigma;

   if(!fTower) return;

-//  ecalEnergy = gRandom->Gaus(fTowerECalEnergy, fECalResolutionFormula->Eval(0.0, fTowerEta, 0.0, fTowerECalEnergy));
+  ecalSigma = fECalResolutionFormula->Eval(0.0, fTowerEta, 0.0, fTowerECalEnergy);
+
+//  ecalEnergy = gRandom->Gaus(fTowerECalEnergy, ecalSigma);
 //  if(ecalEnergy < 0.0) ecalEnergy = 0.0;

-  ecalEnergy = LogNormal(fTowerECalEnergy, fECalResolutionFormula->Eval(0.0, fTowerEta, 0.0, fTowerECalEnergy));
+  ecalEnergy = LogNormal(fTowerECalEnergy, ecalSigma);
+
+  hcalSigma = fHCalResolutionFormula->Eval(0.0, fTowerEta, 0.0, fTowerHCalEnergy);

-//  hcalEnergy = gRandom->Gaus(fTowerHCalEnergy, fHCalResolutionFormula->Eval(0.0, fTowerEta, 0.0, fTowerHCalEnergy));
+//  hcalEnergy = gRandom->Gaus(fTowerHCalEnergy, hcalSigma);
 //  if(hcalEnergy < 0.0) hcalEnergy = 0.0;

-  hcalEnergy = LogNormal(fTowerHCalEnergy, fHCalResolutionFormula->Eval(0.0, fTowerEta, 0.0, fTowerHCalEnergy));
+  hcalEnergy = LogNormal(fTowerHCalEnergy, hcalSigma);

   energy = ecalEnergy + hcalEnergy;

@@ -408,8 +422,7 @@
   // fill calorimeter towers and photon candidates
   if(energy > 0.0)
   {
-    if((fTowerPhotonHits > 0 || fTowerElectronHits > 0) &&
-        fTowerTrackHits == 0)
+    if(fTowerPhotonHits > 0 && fTowerTrackHits == 0)
     {
       fPhotonOutputArray->Add(fTower);
     }
@@ -418,55 +431,35 @@
   }

   // fill energy flow candidates
-  if(fTowerTrackHits == fTowerAllHits)
+
+  // save all the tracks as energy flow tracks
+  fItTowerTrackArray->Reset();
+  while((track = static_cast<Candidate*>(fItTowerTrackArray->Next())))
   {
-    fItTowerTrackArray->Reset();
-    while((track = static_cast<Candidate*>(fItTowerTrackArray->Next())))
-    {
-      fEFlowTrackOutputArray->Add(track);
-    }
+    fEFlowTrackOutputArray->Add(track);
   }
-  else if(fTowerTrackHits > 0 &&
-          fTowerElectronHits == 0 &&
-          fTowerPhotonHits + fTowerTrackHits == fTowerAllHits)
-  {
-    fItTowerTrackArray->Reset();
-    while((track = static_cast<Candidate*>(fItTowerTrackArray->Next())))
-    {
-      fEFlowTrackOutputArray->Add(track);
-    }

-    if(ecalEnergy > 0.0)
-    {
-      DelphesFactory *factory = GetFactory();
+  ecalEnergy -= fTrackECalEnergy;
+  if(ecalEnergy < 0.0) ecalEnergy = 0.0;

-      // create new tower
-      tower = factory->NewCandidate();
+  hcalEnergy -= fTrackHCalEnergy;
+  if(hcalEnergy < 0.0) hcalEnergy = 0.0;

-      fItTowerPhotonArray->Reset();
-      while((particle = static_cast<Candidate*>(fItTowerPhotonArray->Next())))
-      {
-        tower->AddCandidate(particle);
-      }
+  energy = ecalEnergy + hcalEnergy;

-      pt = ecalEnergy / TMath::CosH(eta);
+  // save ECAL and/or HCAL energy excess as an energy flow tower
+  if(energy > 0.0)
+  {
+    // create new tower
+    tower = static_cast<Candidate*>(fTower->Clone());

-      tower->Position.SetPtEtaPhiE(1.0, eta, phi, 0.0);
-      tower->Momentum.SetPtEtaPhiE(pt, eta, phi, ecalEnergy);
-      tower->Eem = ecalEnergy;
-      tower->Ehad = 0.0;
-
-      tower->Edges[0] = fTowerEdges[0];
-      tower->Edges[1] = fTowerEdges[1];
-      tower->Edges[2] = fTowerEdges[2];
-      tower->Edges[3] = fTowerEdges[3];
+    pt = energy / TMath::CosH(eta);

-      fEFlowTowerOutputArray->Add(tower);
-    }
-  }
-  else if(energy > 0.0)
-  {
-    fEFlowTowerOutputArray->Add(fTower);
+    tower->Momentum.SetPtEtaPhiE(pt, eta, phi, energy);
+    tower->Eem = ecalEnergy;
+    tower->Ehad = hcalEnergy;
+
+    fEFlowTowerOutputArray->Add(tower);
   }
 }

diff -ur CMSUpgrade/Delphes-3.0.10/modules/Calorimeter.h Delphes3/modules/Calorimeter.h
--- CMSUpgrade/Delphes-3.0.10/modules/Calorimeter.h 2014-05-21 18:37:15.941871251 +0900
+++ Delphes3/modules/Calorimeter.h  2013-09-13 21:54:12.000000000 +0900
@@ -6,8 +6,8 @@
  *  Fills calorimeter towers, performs calorimeter resolution smearing,
  *  preselects towers hit by photons and creates energy flow objects.
  *
- *  $Date: 2013-06-25 16:24:34 +0200 (Tue, 25 Jun 2013) $
- *  $Revision: 1145 $
+ *  $Date: 2013-09-03 17:54:56 +0200 (Tue, 03 Sep 2013) $
+ *  $Revision: 1273 $
  *
  *
  *  \author P. Demin - UCL, Louvain-la-Neuve
@@ -43,8 +43,8 @@
   Candidate *fTower;
   Double_t fTowerEta, fTowerPhi, fTowerEdges[4];
   Double_t fTowerECalEnergy, fTowerHCalEnergy;
-  Double_t fTowerECalNeutralEnergy, fTowerHCalNeutralEnergy;
-  Int_t fTowerNeutralHits, fTowerPhotonHits, fTowerElectronHits, fTowerTrackHits, fTowerAllHits;
+  Double_t fTrackECalEnergy, fTrackHCalEnergy;
+  Int_t fTowerTrackHits, fTowerPhotonHits;

   TFractionMap fFractionMap; //!
   TBinMap fBinMap; //!
@@ -54,8 +54,11 @@

   std::vector < Long64_t > fTowerHits;

-  std::vector < Double_t > fECalFractions;
-  std::vector < Double_t > fHCalFractions;
+  std::vector < Double_t > fTowerECalFractions;
+  std::vector < Double_t > fTowerHCalFractions;
+
+  std::vector < Double_t > fTrackECalFractions;
+  std::vector < Double_t > fTrackHCalFractions;

   DelphesFormula *fECalResolutionFormula; //!
   DelphesFormula *fHCalResolutionFormula; //!
@@ -75,9 +78,6 @@
   TObjArray *fTowerTrackArray; //!
   TIterator *fItTowerTrackArray; //!

-  TObjArray *fTowerPhotonArray; //!
-  TIterator *fItTowerPhotonArray; //!
-
   void FinalizeTower();
   Double_t LogNormal(Double_t mean, Double_t sigma);
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

No branches or pull requests

1 participant