Changeset 1315 in svn for trunk/modules/FastJetFinder.cc

Timestamp:

Nov 4, 2013, 11:59:27 AM (11 years ago)

Author:

Pavel Demin

Message:

add eta ranges for rho

File:

• trunk/modules/FastJetFinder.cc (modified) (9 diffs)

trunk/modules/FastJetFinder.cc

@@ -30 +30 @@
 #include "TLorentzVector.h"
 
-#include <algorithm> 
+#include <algorithm>
 #include <stdexcept>
 #include <iostream>
@@ -70 +70 @@
 {
   JetDefinition::Plugin *plugin = NULL;
+
+  // read eta ranges
+
+  ExRootConfParam param = GetParam("RhoEtaRange");
+  Long_t i, size;
+
+  fEtaRangeMap.clear();
+  size = param.GetSize();
+  for(i = 0; i < size/2; ++i)
+  {
+    fEtaRangeMap[param[i*2].GetDouble()] = param[i*2 + 1].GetDouble();
+  }
 
   // define algorithm
@@ -90 +102 @@
   fAreaAlgorithm = GetInt("AreaAlgorithm", 0);
   fComputeRho = GetBool("ComputeRho", false);
-  fRhoEtaMax = GetDouble("RhoEtaMax", 5.0);
   // - ghost based areas -
   fGhostEtaMax = GetDouble("GhostEtaMax", 5.0);
@@ -126 +137 @@
   switch(fJetAlgorithm)
   {
-    case 1: 
+    case 1:
       plugin = new fastjet::CDFJetCluPlugin(fSeedThreshold, fConeRadius, fAdjacencyCut, fMaxIterations, fIratch, fOverlapThreshold);
       fDefinition = new fastjet::JetDefinition(plugin);
@@ -149 +160 @@
       break;
   }
-  
+
   fPlugin = plugin;
 
   ClusterSequence::print_banner();
-  
+
   // import input array
 
@@ -187 +198 @@
   vector<PseudoJet> inputList, outputList;
   ClusterSequence *sequence;
+  map< Double_t, Double_t >::iterator itEtaRangeMap;
 
   DelphesFactory *factory = GetFactory();
@@ -196 +208 @@
   number = 0;
   while((candidate = static_cast<Candidate*>(fItInputArray->Next())))
-  {   
+  {
     momentum = candidate->Momentum;
     jet = PseudoJet(momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E());
@@ -203 +215 @@
     ++number;
   }
-   
-  // construct jets 
+
+  // construct jets
   if(fAreaDefinition)
   {
@@ -212 +224 @@
   {
     sequence = new ClusterSequence(inputList, *fDefinition);
-  } 
+  }
 
   // compute rho and store it
   if(fComputeRho && fAreaDefinition)
   {
-    Selector select_rapidity = SelectorAbsRapMax(fRhoEtaMax);
-    JetMedianBackgroundEstimator estimator(select_rapidity, *fDefinition, *fAreaDefinition);
-    estimator.set_particles(inputList);
-    rho = estimator.rho();
-
-    candidate = factory->NewCandidate();
-    candidate->Momentum.SetPtEtaPhiE(rho, 0.0, 0.0, rho); 
-    fRhoOutputArray->Add(candidate);
-  }
-  
+    for(itEtaRangeMap = fEtaRangeMap.begin(); itEtaRangeMap != fEtaRangeMap.end(); ++itEtaRangeMap)
+    {
+      Selector select_rapidity = SelectorAbsRapRange(itEtaRangeMap->first, itEtaRangeMap->second);
+      JetMedianBackgroundEstimator estimator(select_rapidity, *fDefinition, *fAreaDefinition);
+      estimator.set_particles(inputList);
+      rho = estimator.rho();
+
+      candidate = factory->NewCandidate();
+      candidate->Momentum.SetPtEtaPhiE(rho, 0.0, 0.0, rho);
+      candidate->Edges[0] = itEtaRangeMap->first;
+      candidate->Edges[1] = itEtaRangeMap->second;
+      fRhoOutputArray->Add(candidate);
+    }
+  }
+
   outputList.clear();
   outputList = sorted_by_pt(sequence->inclusive_jets(fJetPTMin));