- Timestamp:
- Nov 4, 2013, 11:59:27 AM (11 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/modules/FastJetFinder.cc
r1122 r1315 30 30 #include "TLorentzVector.h" 31 31 32 #include <algorithm> 32 #include <algorithm> 33 33 #include <stdexcept> 34 34 #include <iostream> … … 70 70 { 71 71 JetDefinition::Plugin *plugin = NULL; 72 73 // read eta ranges 74 75 ExRootConfParam param = GetParam("RhoEtaRange"); 76 Long_t i, size; 77 78 fEtaRangeMap.clear(); 79 size = param.GetSize(); 80 for(i = 0; i < size/2; ++i) 81 { 82 fEtaRangeMap[param[i*2].GetDouble()] = param[i*2 + 1].GetDouble(); 83 } 72 84 73 85 // define algorithm … … 90 102 fAreaAlgorithm = GetInt("AreaAlgorithm", 0); 91 103 fComputeRho = GetBool("ComputeRho", false); 92 fRhoEtaMax = GetDouble("RhoEtaMax", 5.0);93 104 // - ghost based areas - 94 105 fGhostEtaMax = GetDouble("GhostEtaMax", 5.0); … … 126 137 switch(fJetAlgorithm) 127 138 { 128 case 1: 139 case 1: 129 140 plugin = new fastjet::CDFJetCluPlugin(fSeedThreshold, fConeRadius, fAdjacencyCut, fMaxIterations, fIratch, fOverlapThreshold); 130 141 fDefinition = new fastjet::JetDefinition(plugin); … … 149 160 break; 150 161 } 151 162 152 163 fPlugin = plugin; 153 164 154 165 ClusterSequence::print_banner(); 155 166 156 167 // import input array 157 168 … … 187 198 vector<PseudoJet> inputList, outputList; 188 199 ClusterSequence *sequence; 200 map< Double_t, Double_t >::iterator itEtaRangeMap; 189 201 190 202 DelphesFactory *factory = GetFactory(); … … 196 208 number = 0; 197 209 while((candidate = static_cast<Candidate*>(fItInputArray->Next()))) 198 { 210 { 199 211 momentum = candidate->Momentum; 200 212 jet = PseudoJet(momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E()); … … 203 215 ++number; 204 216 } 205 206 // construct jets 217 218 // construct jets 207 219 if(fAreaDefinition) 208 220 { … … 212 224 { 213 225 sequence = new ClusterSequence(inputList, *fDefinition); 214 } 226 } 215 227 216 228 // compute rho and store it 217 229 if(fComputeRho && fAreaDefinition) 218 230 { 219 Selector select_rapidity = SelectorAbsRapMax(fRhoEtaMax); 220 JetMedianBackgroundEstimator estimator(select_rapidity, *fDefinition, *fAreaDefinition); 221 estimator.set_particles(inputList); 222 rho = estimator.rho(); 223 224 candidate = factory->NewCandidate(); 225 candidate->Momentum.SetPtEtaPhiE(rho, 0.0, 0.0, rho); 226 fRhoOutputArray->Add(candidate); 227 } 228 231 for(itEtaRangeMap = fEtaRangeMap.begin(); itEtaRangeMap != fEtaRangeMap.end(); ++itEtaRangeMap) 232 { 233 Selector select_rapidity = SelectorAbsRapRange(itEtaRangeMap->first, itEtaRangeMap->second); 234 JetMedianBackgroundEstimator estimator(select_rapidity, *fDefinition, *fAreaDefinition); 235 estimator.set_particles(inputList); 236 rho = estimator.rho(); 237 238 candidate = factory->NewCandidate(); 239 candidate->Momentum.SetPtEtaPhiE(rho, 0.0, 0.0, rho); 240 candidate->Edges[0] = itEtaRangeMap->first; 241 candidate->Edges[1] = itEtaRangeMap->second; 242 fRhoOutputArray->Add(candidate); 243 } 244 } 245 229 246 outputList.clear(); 230 247 outputList = sorted_by_pt(sequence->inclusive_jets(fJetPTMin));
Note:
See TracChangeset
for help on using the changeset viewer.