Index: trunk/modules/MadGraphCDFConeJetFinder.cc =================================================================== --- trunk/modules/MadGraphCDFConeJetFinder.cc (revision 11) +++ trunk/modules/MadGraphCDFConeJetFinder.cc (revision 11) @@ -0,0 +1,111 @@ + +#include "modules/MadGraphCDFConeJetFinder.h" + +#include "ExRootAnalysis/ExRootClasses.h" +#include "ExRootAnalysis/ExRootFactory.h" +#include "ExRootAnalysis/ExRootCandidate.h" + +#include "CDFCones/JetCluAlgorithm.hh" +#include "CDFCones/MidPointAlgorithm.hh" + +#include "TString.h" +#include "TLorentzVector.h" + +#include +#include + +using namespace std; + +//------------------------------------------------------------------------------ + +MadGraphCDFConeJetFinder::MadGraphCDFConeJetFinder() : + fJetAlgo(0), fItInputArray(0) +{ + +} + +//------------------------------------------------------------------------------ + +MadGraphCDFConeJetFinder::~MadGraphCDFConeJetFinder() +{ + +} + +//------------------------------------------------------------------------------ + +void MadGraphCDFConeJetFinder::Init() +{ + + // define MidPoint algorithm + + double seedThreshold = GetDouble("SeedThreshold", 1.0); + double coneRadius = GetDouble("ConeRadius", 0.5); + double coneAreaFraction = GetDouble("ConeAreaFraction", 1.0); + int maxPairSize = GetInt("MaxPairSize", 2); + int maxIterations = GetInt("MaxIterations", 100); + double overlapThreshold = GetDouble("OverlapThreshold", 0.75); + + fJetAlgo = new MidPointAlgorithm(seedThreshold, coneRadius, coneAreaFraction, + maxPairSize, maxIterations, overlapThreshold); + + // import array with output from filter/classifier module + + fInputArray = ImportArray(GetString("InputArray", "selection/candidates")); + fItInputArray = fInputArray->MakeIterator(); + + // create output arrays + + fOutputArray = ExportArray("candidates"); + +} + +//------------------------------------------------------------------------------ + +void MadGraphCDFConeJetFinder::Finish() +{ + if(fJetAlgo) delete fJetAlgo; + if(fItInputArray) delete fItInputArray; +} + +//------------------------------------------------------------------------------ + +void MadGraphCDFConeJetFinder::Process() +{ + ExRootCandidate *candidate; + TLorentzVector momentum; + LorentzVector jetMomentum; + Int_t entry; + + ExRootFactory *factory = GetFactory(); + + fTowersList.clear(); + + // loop over all particles in event and select stable ones + fItInputArray->Reset(); + while((candidate = static_cast(fItInputArray->Next()))) + { + momentum = candidate->GetP4(); + fTowersList.push_back(PhysicsTower(LorentzVector(momentum.Px(), momentum.Py(), + momentum.Pz(), momentum.E()))); + } + + // construct jets from a list of stable particles + fJetsList.clear(); + fJetAlgo->run(fTowersList, fJetsList); + + // loop over all jets and export them + vector::iterator itJet; + for(itJet = fJetsList.begin(), entry = 1; itJet != fJetsList.end(); ++itJet, ++entry) + { + jetMomentum = itJet->fourVector; + + momentum.SetPxPyPzE(jetMomentum.px, jetMomentum.py, jetMomentum.pz, jetMomentum.E); + + candidate = factory->NewCandidate(); + + candidate->SetP4(momentum); + candidate->SetName(Form("jet_{%d}", entry )); + + fOutputArray->Add(candidate); + } +} Index: trunk/modules/MadGraphCDFConeJetFinder.h =================================================================== --- trunk/modules/MadGraphCDFConeJetFinder.h (revision 11) +++ trunk/modules/MadGraphCDFConeJetFinder.h (revision 11) @@ -0,0 +1,43 @@ +#ifndef MadGraphCDFConeJetFinder_h +#define MadGraphCDFConeJetFinder_h + +#include "ExRootAnalysis/ExRootModule.h" + +#include "CDFCones/PhysicsTower.hh" +#include "CDFCones/Cluster.hh" + +#include + +class TObjArray; +class TIterator; + +class MidPointAlgorithm; + +class MadGraphCDFConeJetFinder: public ExRootModule +{ +public: + + MadGraphCDFConeJetFinder(); + ~MadGraphCDFConeJetFinder(); + + void Init(); + void Process(); + void Finish(); + +private: + + std::vector fTowersList; + std::vector fJetsList; + + MidPointAlgorithm *fJetAlgo; //! + + TIterator *fItInputArray; //! + + const TObjArray *fInputArray; //! + + TObjArray *fOutputArray; //! + + ClassDef(MadGraphCDFConeJetFinder, 1) +}; + +#endif Index: trunk/modules/MadGraphConeJetFinder.cc =================================================================== --- trunk/modules/MadGraphConeJetFinder.cc (revision 2) +++ (revision ) @@ -1,111 +1,0 @@ - -#include "modules/MadGraphConeJetFinder.h" - -#include "ExRootAnalysis/ExRootClasses.h" -#include "ExRootAnalysis/ExRootFactory.h" -#include "ExRootAnalysis/ExRootCandidate.h" - -#include "CDFCones/JetCluAlgorithm.hh" -#include "CDFCones/MidPointAlgorithm.hh" - -#include "TString.h" -#include "TLorentzVector.h" - -#include -#include - -using namespace std; - -//------------------------------------------------------------------------------ - -MadGraphConeJetFinder::MadGraphConeJetFinder() : - fJetAlgo(0), fItInputArray(0) -{ - -} - -//------------------------------------------------------------------------------ - -MadGraphConeJetFinder::~MadGraphConeJetFinder() -{ - -} - -//------------------------------------------------------------------------------ - -void MadGraphConeJetFinder::Init() -{ - - // define MidPoint algorithm - - double seedThreshold = GetDouble("SeedThreshold", 1.0); - double coneRadius = GetDouble("ConeRadius", 0.5); - double coneAreaFraction = GetDouble("ConeAreaFraction", 1.0); - int maxPairSize = GetInt("MaxPairSize", 2); - int maxIterations = GetInt("MaxIterations", 100); - double overlapThreshold = GetDouble("OverlapThreshold", 0.75); - - fJetAlgo = new MidPointAlgorithm(seedThreshold, coneRadius, coneAreaFraction, - maxPairSize, maxIterations, overlapThreshold); - - // import array with output from filter/classifier module - - fInputArray = ImportArray(GetString("InputArray", "selection/candidates")); - fItInputArray = fInputArray->MakeIterator(); - - // create output arrays - - fOutputArray = ExportArray("candidates"); - -} - -//------------------------------------------------------------------------------ - -void MadGraphConeJetFinder::Finish() -{ - if(fJetAlgo) delete fJetAlgo; - if(fItInputArray) delete fItInputArray; -} - -//------------------------------------------------------------------------------ - -void MadGraphConeJetFinder::Process() -{ - ExRootCandidate *candidate; - TLorentzVector momentum; - LorentzVector jetMomentum; - Int_t entry; - - ExRootFactory *factory = GetFactory(); - - fTowersList.clear(); - - // loop over all particles in event and select stable ones - fItInputArray->Reset(); - while((candidate = static_cast(fItInputArray->Next()))) - { - momentum = candidate->GetP4(); - fTowersList.push_back(PhysicsTower(LorentzVector(momentum.Px(), momentum.Py(), - momentum.Pz(), momentum.E()))); - } - - // construct jets from a list of stable particles - fJetsList.clear(); - fJetAlgo->run(fTowersList, fJetsList); - - // loop over all jets and export them - vector::iterator itJet; - for(itJet = fJetsList.begin(), entry = 1; itJet != fJetsList.end(); ++itJet, ++entry) - { - jetMomentum = itJet->fourVector; - - momentum.SetPxPyPzE(jetMomentum.px, jetMomentum.py, jetMomentum.pz, jetMomentum.E); - - candidate = factory->NewCandidate(); - - candidate->SetP4(momentum); - candidate->SetName(Form("jet_{%d}", entry )); - - fOutputArray->Add(candidate); - } -} Index: trunk/modules/MadGraphConeJetFinder.h =================================================================== --- trunk/modules/MadGraphConeJetFinder.h (revision 2) +++ (revision ) @@ -1,43 +1,0 @@ -#ifndef MadGraphConeJetFinder_h -#define MadGraphConeJetFinder_h - -#include "ExRootAnalysis/ExRootModule.h" - -#include "CDFCones/PhysicsTower.hh" -#include "CDFCones/Cluster.hh" - -#include - -class TObjArray; -class TIterator; - -class MidPointAlgorithm; - -class MadGraphConeJetFinder: public ExRootModule -{ -public: - - MadGraphConeJetFinder(); - ~MadGraphConeJetFinder(); - - void Init(); - void Process(); - void Finish(); - -private: - - std::vector fTowersList; - std::vector fJetsList; - - MidPointAlgorithm *fJetAlgo; //! - - TIterator *fItInputArray; //! - - const TObjArray *fInputArray; //! - - TObjArray *fOutputArray; //! - - ClassDef(MadGraphConeJetFinder, 1) -}; - -#endif Index: trunk/modules/MadGraphPartonSelector.cc =================================================================== --- trunk/modules/MadGraphPartonSelector.cc (revision 2) +++ trunk/modules/MadGraphPartonSelector.cc (revision 11) @@ -33,14 +33,17 @@ void InsertParticleID(Int_t pid); - void InsertAncestorID(Int_t pid); + void InsertExclAncestorID(Int_t pid); + void InsertInclAncestorID(Int_t pid); private: Bool_t hasBadAncestor(ExRootGenParticle *object); + Bool_t hasGoodAncestor(ExRootGenParticle *object); TClonesArray *fBranchParticle; set< Int_t > fParticleIDSet; - set< Int_t > fAncestorIDSet; + set< Int_t > fExclAncestorIDSet; + set< Int_t > fInclAncestorIDSet; }; @@ -61,7 +64,14 @@ //------------------------------------------------------------------------------ -void MadGraphPartonClassifier::InsertAncestorID(Int_t pid) -{ - fAncestorIDSet.insert(pid); +void MadGraphPartonClassifier::InsertExclAncestorID(Int_t pid) +{ + fExclAncestorIDSet.insert(pid); +} + +//------------------------------------------------------------------------------ + +void MadGraphPartonClassifier::InsertInclAncestorID(Int_t pid) +{ + fInclAncestorIDSet.insert(pid); } @@ -76,4 +86,10 @@ set< Int_t >::const_iterator itParticleIDSet; + // keep all particles if there is no pid in the list + if(fExclAncestorIDSet.empty()) + { + return kFALSE; + } + for(i = 0; i < kMaxAncestors; ++i) { @@ -90,7 +106,45 @@ // skip particles with pid included in list - itAncestorIDSet = fAncestorIDSet.find(pidAbs); - - if(itAncestorIDSet != fAncestorIDSet.end()) return kTRUE; + itAncestorIDSet = fExclAncestorIDSet.find(pidAbs); + + if(itAncestorIDSet != fExclAncestorIDSet.end()) return kTRUE; + } + + return kFALSE; +} + +//------------------------------------------------------------------------------ + +Bool_t MadGraphPartonClassifier::hasGoodAncestor(ExRootGenParticle *object) +{ + const int kMaxAncestors = 10; + Int_t i, pidAbs; + ExRootGenParticle *particle = object; + set< Int_t >::const_iterator itAncestorIDSet; + set< Int_t >::const_iterator itParticleIDSet; + + // keep all particles if there is no pid in the list + if(fInclAncestorIDSet.empty()) + { + return kTRUE; + } + + for(i = 0; i < kMaxAncestors; ++i) + { + + if(particle->M1 < 0 || particle->M2 > -1) return kFALSE; + +// if(particle->PID == 21) return kFALSE; + + particle = static_cast(fBranchParticle->At(particle->M1)); + + if(particle->PID == 21) return kFALSE; + + pidAbs = TMath::Abs(particle->PID); + + // keep particles with pid included in list + itAncestorIDSet = fInclAncestorIDSet.find(pidAbs); + + if(itAncestorIDSet != fInclAncestorIDSet.end()) return kTRUE; } @@ -145,4 +199,6 @@ if(hasBadAncestor(particle)) return -1; + if(!hasGoodAncestor(particle)) return -1; + return 0; } @@ -197,5 +253,16 @@ { pid = param[i].GetInt(); - fClassifier->InsertAncestorID(pid); + fClassifier->InsertExclAncestorID(pid); + } + + // read ancestor IDs from configuration file and setup classifier + + param = GetParam("IncludedAncestorIDs"); + sizeParam = param.GetSize(); + + for(i = 0; i < sizeParam; ++i) + { + pid = param[i].GetInt(); + fClassifier->InsertInclAncestorID(pid); } @@ -209,5 +276,5 @@ void MadGraphPartonSelector::Finish() -{ +{ if(fFilter) delete fFilter; if(fClassifier) delete fClassifier; Index: trunk/modules/MadGraphSISConeJetFinder.cc =================================================================== --- trunk/modules/MadGraphSISConeJetFinder.cc (revision 11) +++ trunk/modules/MadGraphSISConeJetFinder.cc (revision 11) @@ -0,0 +1,106 @@ + +#include "modules/MadGraphSISConeJetFinder.h" + +#include "ExRootAnalysis/ExRootClasses.h" +#include "ExRootAnalysis/ExRootFactory.h" +#include "ExRootAnalysis/ExRootCandidate.h" + +#include "SISCone/momentum.h" +#include "SISCone/siscone.h" + +#include "TString.h" +#include "TLorentzVector.h" + +#include +#include + +using namespace std; +using namespace siscone; + +//------------------------------------------------------------------------------ + +MadGraphSISConeJetFinder::MadGraphSISConeJetFinder() : + fJetAlgo(0), fItInputArray(0) +{ + +} + +//------------------------------------------------------------------------------ + +MadGraphSISConeJetFinder::~MadGraphSISConeJetFinder() +{ + +} + +//------------------------------------------------------------------------------ + +void MadGraphSISConeJetFinder::Init() +{ + + // define MidPoint algorithm + + fConeRadius = GetDouble("ConeRadius", 0.5); + fConeAreaFraction = GetDouble("ConeAreaFraction", 1.0); + fMaxIterations = GetInt("MaxIterations", 100); + + fJetAlgo = new Csiscone; + + // import array with output from filter/classifier module + + fInputArray = ImportArray(GetString("InputArray", "selection/candidates")); + fItInputArray = fInputArray->MakeIterator(); + + // create output arrays + + fOutputArray = ExportArray("candidates"); + +} + +//------------------------------------------------------------------------------ + +void MadGraphSISConeJetFinder::Finish() +{ + if(fJetAlgo) delete fJetAlgo; + if(fItInputArray) delete fItInputArray; +} + +//------------------------------------------------------------------------------ + +void MadGraphSISConeJetFinder::Process() +{ + ExRootCandidate *candidate; + TLorentzVector momentum; + Cmomentum jetMomentum; + Int_t entry; + + ExRootFactory *factory = GetFactory(); + + fParticlesList.clear(); + + // loop over all particles in event and select stable ones + fItInputArray->Reset(); + while((candidate = static_cast(fItInputArray->Next()))) + { + momentum = candidate->GetP4(); + fParticlesList.push_back(Cmomentum(momentum.Px(), momentum.Py(), momentum.Pz(), momentum.E())); + } + + // construct jets from a list of stable particles + fJetAlgo->compute_jets(fParticlesList, fConeRadius, fConeAreaFraction, fMaxIterations); + + // loop over all jets and export them + vector::iterator itJet; + for(itJet = fJetAlgo->jets.begin(), entry = 1; itJet != fJetAlgo->jets.end(); ++itJet, ++entry) + { + jetMomentum = itJet->v; + + momentum.SetPxPyPzE(jetMomentum.px, jetMomentum.py, jetMomentum.pz, jetMomentum.E); + + candidate = factory->NewCandidate(); + + candidate->SetP4(momentum); + candidate->SetName(Form("jet_{%d}", entry )); + + fOutputArray->Add(candidate); + } +} Index: trunk/modules/MadGraphSISConeJetFinder.h =================================================================== --- trunk/modules/MadGraphSISConeJetFinder.h (revision 11) +++ trunk/modules/MadGraphSISConeJetFinder.h (revision 11) @@ -0,0 +1,49 @@ +#ifndef MadGraphSISConeJetFinder_h +#define MadGraphSISConeJetFinder_h + +#include "ExRootAnalysis/ExRootModule.h" + +#include "SISCone/momentum.h" + +#include + +namespace siscone { + class Csiscone; +} + +class TObjArray; +class TIterator; + +class MidPointAlgorithm; + +class MadGraphSISConeJetFinder: public ExRootModule +{ +public: + + MadGraphSISConeJetFinder(); + ~MadGraphSISConeJetFinder(); + + void Init(); + void Process(); + void Finish(); + +private: + + std::vector fParticlesList; + + siscone::Csiscone *fJetAlgo; //! + + double fConeRadius; + double fConeAreaFraction; + int fMaxIterations; + + TIterator *fItInputArray; //! + + const TObjArray *fInputArray; //! + + TObjArray *fOutputArray; //! + + ClassDef(MadGraphSISConeJetFinder, 1) +}; + +#endif Index: trunk/modules/MadGraphShowerPartonSelector.cc =================================================================== --- trunk/modules/MadGraphShowerPartonSelector.cc (revision 2) +++ trunk/modules/MadGraphShowerPartonSelector.cc (revision 11) @@ -35,5 +35,6 @@ void SetEtaMax(Double_t eta); void InsertParticleID(Int_t pid); - void InsertAncestorID(Int_t pid); + void InsertExclAncestorID(Int_t pid); + void InsertInclAncestorID(Int_t pid); void SetHadronizationInfo(Bool_t info); @@ -41,4 +42,5 @@ Bool_t hasBadAncestor(ExRootGenParticle *object); + Bool_t hasGoodAncestor(ExRootGenParticle *object); Double_t fEtaMax; @@ -47,5 +49,6 @@ set< Int_t > fParticleIDSet; - set< Int_t > fAncestorIDSet; + set< Int_t > fExclAncestorIDSet; + set< Int_t > fInclAncestorIDSet; }; @@ -73,7 +76,14 @@ //------------------------------------------------------------------------------ -void MadGraphShowerPartonClassifier::InsertAncestorID(Int_t pid) -{ - fAncestorIDSet.insert(pid); +void MadGraphShowerPartonClassifier::InsertExclAncestorID(Int_t pid) +{ + fExclAncestorIDSet.insert(pid); +} + +//------------------------------------------------------------------------------ + +void MadGraphShowerPartonClassifier::InsertInclAncestorID(Int_t pid) +{ + fInclAncestorIDSet.insert(pid); } @@ -96,10 +106,16 @@ if(particle->PID == 21) return kFALSE; + // keep all particles if there is no pid in the list + if(fExclAncestorIDSet.empty()) + { + return kFALSE; + } + pidAbs = TMath::Abs(particle->PID); // skip particles with pid included in list - itAncestorIDSet = fAncestorIDSet.find(pidAbs); - - if(itAncestorIDSet != fAncestorIDSet.end()) return kTRUE; + itAncestorIDSet = fExclAncestorIDSet.find(pidAbs); + + if(itAncestorIDSet != fExclAncestorIDSet.end()) return kTRUE; if(particle->M2 > -1) return kFALSE; @@ -115,7 +131,59 @@ // skip particles with pid included in list - itAncestorIDSet = fAncestorIDSet.find(pidAbs); - - if(itAncestorIDSet != fAncestorIDSet.end()) return kTRUE; + itAncestorIDSet = fExclAncestorIDSet.find(pidAbs); + + if(itAncestorIDSet != fExclAncestorIDSet.end()) return kTRUE; + if(particle->M2 > -1) return kFALSE; + } + + return kFALSE; +} + +//------------------------------------------------------------------------------ + +Bool_t MadGraphShowerPartonClassifier::hasGoodAncestor(ExRootGenParticle *object) +{ + const int kMaxAncestors = 10; + Int_t i, pidAbs; + ExRootGenParticle *particle = object; + set< Int_t >::const_iterator itAncestorIDSet; + + for(i = 0; i < kMaxAncestors && particle->Status != 3; ++i) + { + if(particle->M1 < 0) return kFALSE; + + particle = static_cast(fBranchParticle->At(particle->M1)); + } + + if(particle->PID == 21) return kFALSE; + + // keep all particles if there is no pid in the list + if(fInclAncestorIDSet.empty()) + { + return kTRUE; + } + + pidAbs = TMath::Abs(particle->PID); + + // keep particles with pid included in list + itAncestorIDSet = fInclAncestorIDSet.find(pidAbs); + + if(itAncestorIDSet != fInclAncestorIDSet.end()) return kTRUE; + if(particle->M2 > -1) return kFALSE; + + for(i = 0; i < kMaxAncestors; ++i) + { + if(particle->M1 < 0) return kFALSE; + + if(particle->PID == 21) return kFALSE; + + particle = static_cast(fBranchParticle->At(particle->M1)); + + pidAbs = TMath::Abs(particle->PID); + + // keep particles with pid included in list + itAncestorIDSet = fInclAncestorIDSet.find(pidAbs); + + if(itAncestorIDSet != fInclAncestorIDSet.end()) return kTRUE; if(particle->M2 > -1) return kFALSE; } @@ -159,4 +227,6 @@ if(hasBadAncestor(particle)) return -1; + if(!hasGoodAncestor(particle)) return -1; + return 0; } @@ -214,5 +284,16 @@ { pid = param[i].GetInt(); - fClassifier->InsertAncestorID(pid); + fClassifier->InsertExclAncestorID(pid); + } + + // read ancestor IDs from configuration file and setup classifier + + param = GetParam("IncludedAncestorIDs"); + sizeParam = param.GetSize(); + + for(i = 0; i < sizeParam; ++i) + { + pid = param[i].GetInt(); + fClassifier->InsertInclAncestorID(pid); } Index: trunk/modules/ModulesLinkDef.h =================================================================== --- trunk/modules/ModulesLinkDef.h (revision 2) +++ trunk/modules/ModulesLinkDef.h (revision 11) @@ -7,5 +7,6 @@ #include "modules/MadGraphMatchingTreeWriter.h" #include "modules/MadGraphKtJetFinder.h" -#include "modules/MadGraphConeJetFinder.h" +#include "modules/MadGraphCDFConeJetFinder.h" +#include "modules/MadGraphSISConeJetFinder.h" #include "modules/MadGraphIsolatedLeptonFinder.h" #include "modules/PythiaFix.h" @@ -29,5 +30,6 @@ #pragma link C++ class MadGraphMatchingTreeWriter+; #pragma link C++ class MadGraphKtJetFinder+; -#pragma link C++ class MadGraphConeJetFinder+; +#pragma link C++ class MadGraphCDFConeJetFinder+; +#pragma link C++ class MadGraphSISConeJetFinder+; #pragma link C++ class MadGraphIsolatedLeptonFinder+; #pragma link C++ class PythiaFix+;