//================================================================== // // TTopFindModule.cc // // Version 0.0: P. Koehn 5/5/01 // //================================================================== #include "TF1.h" #include #include #include #include "Stntuple/loop/TStnAna.hh" #include "Stntuple/obj/TStnHeaderBlock.hh" #include "TTopFindModule.hh" //_____________________________________________________________________________ TTopFindModule::TTopFindModule(TFile* file,const char* name, const char* title): TStnModule(name,title){ // // create a subdirectory "dir" in this file histFile = file; histDir = name; dir = histFile->mkdir(histDir); // // Initialize cuts init(); } //_____________________________________________________________________________ TTopFindModule::TTopFindModule(const char* name, const char* title): TStnModule(name,title){ // // create a subdirectory "dir" in this file histFile = new TFile("TEleCand.hbk","RECREATE","test file"); histDir = name; dir = histFile->mkdir(histDir); // // Initialize cuts init(); } void TTopFindModule::init() { // totEvents = 0; passEvents = 0; // // Default is to make histograms enableHist(true); // // Define the default cuts for this module enableNumberLeptonCut(true); setNumberLeptonCut(1); enableMaxNumberJetCut(false); setMaxNumberJetCut(1); enableMinNumberJetCut(true); setMinNumberJetCut(1); enableMetCut(true); setMetCut(20.0); printf("Cut Settings:\n"); printf(" Number Lepton Cut: enabled %l value %d\n", isEnableNumberLeptonCut(),getNumberLeptonCut()); } //_____________________________________________________________________________ TTopFindModule::~TTopFindModule() { } //_____________________________________________________________________________ void TTopFindModule::BookHistograms() { // // Delete the old histos // GetHistogramList()->Delete(); histFile->cd(histDir); // make the "dir" directory the current directory //----------------------------------------------------------------------------- // electron histograms //----------------------------------------------------------------------------- fHist.fEle.EtHist = new TH1F("fEleEtHist","EleCand Et",50,0,100); fHist.fEle.EoverPHist = new TH1F("fEleEoverPHist","EleCand E over P",50,0,5); fHist.fEle.DelxHist = new TH1F("fEleDelxHist","EleCand DeltaX",50,-20,20); fHist.fEle.DelzHist = new TH1F("fEleDelzHist","EleCand DeltaZ",50,-20,20); fHist.fEle.LshrHist = new TH1F("fEleLshrHist","EleCand Lshr",50,0.0,5.0); fHist.fEle.HademHist = new TH1F("fEleHademHist","EleCand HadEm",50,0,1.0); fHist.fEle.IsolHist = new TH1F("fEleIsolHist","EleCand Iso",50,0,2); fHist.fEle.TisolHist = new TH1F("fEleTisolHist","EleCand Track Iso",50,0,100); fHist.fEle.ChisHist = new TH1F("fEleChisHist","EleCand ChiStrip",50,0,100); fHist.fEle.DelZvrtHist = new TH1F("fEleDelZvrtHist","EleCand Delta(Zvtx-Zele)",50,-10,10); fHist.fEle.ZvrtHist = new TH1F("fEleZvrtHist","EleCand Event Zvertex",50,-100,100); fHist.fEle.DetHist = new TH1F("fEleDetHist","EleCand Detector",3,0,2); //----------------------------------------------------------------------------- // jet histograms //----------------------------------------------------------------------------- fHist.fJet.EtHist = new TH1F("fJetEtHist","JetCand Et",50,0,200); fHist.fJet.EvetaHist = new TH1F("fJetEvetaHist","JetCand Evt Eta",50,-5,5); fHist.fJet.DtetaHist = new TH1F("fJetDtetaHist","JetCand Det Eta",50,-5,5); fHist.fJet.EmfHist = new TH1F("fJetEmfHist","JetCand EMF",50,0.0,1.0); //----------------------------------------------------------------------------- // muon histograms //----------------------------------------------------------------------------- fHist.fMuon.fN = new TH1F("fMuonfN", "Number of muons",21,-0.5,20.5); fHist.fMuon.fDetCode = new TH1F("fMuonfDetCode", "DetCode",21,-0.5,20.5); fHist.fMuon.fCmuCmpPt = new TH1F("fMuonfCmuCmpPt", "CmuCmpPt",100,0,100); fHist.fMuon.fCmuCmpPhi = new TH1F("fMuonfCmuCmpPhi", "CmuCmpPhi",210,0,6.3); fHist.fMuon.fCmuCmpEta = new TH1F("fMuonfCmuCmpEta", "CmuCmpEta",150,-1.5,1.5); fHist.fMuon.fCmuCmpEmEnergy = new TH1F("fMuonfCmuCmpEmEnergy", "CmuCmpEmEnergy",100,0,5); fHist.fMuon.fCmuCmpHadEnergy = new TH1F("fMuonfCmuCmpHadEnergy", "CmuCmpHadEnergy",100,0,10); fHist.fMuon.fCmuCmpHadEm = new TH1F("fMuonfCmuCmpHadEm", "CmuCmpHadEm",50,0,10); fHist.fMuon.fCmxPt = new TH1F("fMuonfCmxPt", "CmxPt",100,0,100); fHist.fMuon.fCmxPhi = new TH1F("fMuonfCmxPhi", "CmxPhi",210,0,6.3); fHist.fMuon.fCmxEta = new TH1F("fMuonfCmxEta", "CmxEta",150,-1.5,1.5); fHist.fMuon.fCmxEmEnergy = new TH1F("fMuonfCmxEmEnergy", "CmxEmEnergy",100,0,5); fHist.fMuon.fCmxHadEnergy = new TH1F("fMuonfCmxHadEnergy", "",100,0,10); fHist.fMuon.fCmxHadEm = new TH1F("fMuonfCmxHadEm", "",50,0,10); fHist.fMuon.fISO = new TH1F("fMuonfISO", "ISO",40,0,0.2); fHist.fMuon.fTISO = new TH1F("fMuonfTISO", "Track ISO",40,0,0.2); fHist.fMuon.fd0 = new TH1F("fMuonfd0", "D0",50,-0.05,0.05); fHist.fMuon.fzv = new TH1F("fMuonfzv", "zv",100,-100,100); fHist.fMuon.fztrk = new TH1F("fMuonfztrk", "ztrk",100,-100,100); fHist.fMuon.fCmuCmpDx = new TH1F("fMuonfCmuCmpDx", "CmuCmpDx",50,-10,10); fHist.fMuon.fCmxDx = new TH1F("fMuonfCmxdx", "Cmxdx",50,-10,10); //----------------------------------------------------------------------------- // missing Et histograms //----------------------------------------------------------------------------- fHist.fMet.MetHist = new TH1F("fMetMetHist","Met",50 ,0,100); fHist.fMet.PhiHist = new TH1F("fMetPhiHist","MetPhi",128 ,0,6.3); fHist.fMet.SumEtHist = new TH1F("fMetSumEtHist", "SumEt",50,0,200); //----------------------------------------------------------------------------- // other histograms //----------------------------------------------------------------------------- fHist.fSkim = new TProfile("obj_skim", "Object Skim Rate",21,-0.5,20.5); fHist.fIsoVsMetHist = new TH2F("fIsoVsMetHist","Ele Iso vs Met", 50 ,0,100,40,0.0,0.2); fHist.fWEleEtHist = new TH1F("fWEleEtHist","EleCand Et, MET>20",50,0,100); fHist.fWEleMtHist = new TH1F("fWEleMtHist","EleCand TransMass, MET>20",25,0,100); fHist.fEtVsMetHist = new TH2F("fEtVsMetHist","Ele Et vs Met", 50 ,0,100,50,0,100); fHist.fZMass = new TH1F("fZMass","Z Mass",50,65.,115.); } //_____________________________________________________________________________ int TTopFindModule::BeginJob() { printf("----- begin job: ---- %s\n",GetName()); printf(" Nlepton cut enabled %d value %d\n", isEnableNumberLeptonCut(),getNumberLeptonCut()); printf(" Nmaxjet cut enabled %d value %d\n", isEnableMaxNumberJetCut(),getMaxNumberJetCut()); printf(" Nminjet cut enabled %d value %d\n", isEnableMinNumberJetCut(),getMinNumberJetCut()); printf(" NMissEt cut enabled %d value %f\n", isEnableMetCut(),getMetCut()); totEvents = 0; passEvents = 0; // register the data blocks pEleBlock = (TStnElectronBlock*) fAna->RegisterDataBlock("ElectronBlock","TStnElectronBlock"); pJetBlock = (TStnJetBlock*) fAna->RegisterDataBlock("JetBlock","TStnJetBlock"); pMetBlock = (TStnMetBlock*) fAna->RegisterDataBlock("MetBlock","TStnMetBlock"); // pMuonBlock = (TStnMuonBlock*) fAna->RegisterDataBlock("MuonBlock","TStnMuonBlock"); // book histograms if (isEnableHist()) BookHistograms(); if (! pEleBlock) { printf(" >>> branch *** %s *** doesn't exist \n","EleBlock"); fEnabled = 0; } if (! pJetBlock) { printf(" >>> branch *** %s *** doesn't exist \n","JetBlock"); fEnabled = 0; } // if (! pMuonBlock) { // printf(" >>> branch *** %s *** doesn't exist \n","MuonBlock"); // fEnabled = 0; // } return 0; } //_____________________________________________________________________________ int TTopFindModule::BeginRun() { return 0; } //================================== // Event(int ientry) //================================== int TTopFindModule::Event(int ientry) { bool pass = false; if (fDebugLevel > 0) printf(" ---- ientry = %7i\n",ientry); totEvents++; // // Get the header info TStnHeaderBlock* header = GetHeaderBlock(); int iev = header->EventNumber(); int irun = header->RunNumber(); // // Met Selection //=================== pMetBlock->GetEntry(ientry); // // Get access to the candidate modules (ele,muo,jet): //--------------------------------------------------- TStnAna * TestAna = TStnModule::GetAna(); TEleCand* eleCand = (TEleCand*) TestAna->TStnAna::GetModule("EleC"); TJetCand* jetCand = (TJetCand*) TestAna->TStnAna::GetModule("JetC"); // TZBosFnd* zBosCand = (TZBosFnd*) TestAna->TStnAna::GetModule("ZBosFnd"); if(zBosCand->isZboson()) { TLorentzVector* Zb = zBosCand->Zmomentum; double Zmass = Zb->M(); fHist.fZMass->Fill(Zmass); } // Get the basic event info int nelectrons = eleCand->getNgoodElectron(); int nmuons = 0; int nleptons = nelectrons + nmuons; int njets = jetCand->getNgoodJet(); njets = jetCand->getNgoodJet(); // printf("njets %d nele %d \n",njets,nelectrons); // // Lets make our cuts for this event if( ((!isEnableNumberLeptonCut()) || (nleptons >= getNumberLeptonCut())) && ((!isEnableMaxNumberJetCut()) || (njets <= getMaxNumberJetCut())) && ((!isEnableMinNumberJetCut()) || (njets >= getMinNumberJetCut())) && ((!isEnableMetCut()) || (pMetBlock->fMet[0] >= getMetCut())) ) { pass = true; passEvents++; } // // Some debug printout if (fDebugLevel > 0) { printf(" nele %d nlep %d njet %d pass %d\n", nelectrons,nleptons,njets,pass); } //================================================= // filter on events passing all cuts //================================================= if (pass) { SetPassed(1); } else { SetPassed(0); } //================================================ // Histogram the Electrons in events passing cuts= //================================================ if (pass) { if (fDebugLevel > 0) { printf(" Good event run %d event %d\n",irun,iev); } // // If hists are enabled, make these plots for all good electron candidates if (isEnableHist()) { fHist.fMet.MetHist->Fill(pMetBlock->fMet[0]); fHist.fMet.PhiHist->Fill(pMetBlock->fMetphi[0]); fHist.fMet.SumEtHist->Fill(pMetBlock->fSumet); } // // Now loop over the electrons for (int i=0; igetgoodElectron(i); // // If hists are enabled, make these plots for all good electron candidates if (isEnableHist()) { fHist.fEle.EtHist->Fill(electron->Et()); fHist.fEle.EoverPHist->Fill(electron->EoverP()); fHist.fEle.HademHist->Fill(electron->HadEm()); fHist.fEle.IsolHist->Fill(electron->Iso()/electron->Et()); fHist.fEle.TisolHist->Fill(electron->Tiso()); fHist.fEle.DetHist->Fill(electron->DetectorCode()); fHist.fEle.DelxHist->Fill(electron->DelX()); fHist.fEle.DelzHist->Fill(electron->DelZ()); fHist.fEle.LshrHist->Fill(electron->Lshr()); fHist.fEle.ChisHist->Fill(electron->Chi2Strip()); fHist.fEle.DelZvrtHist->Fill(electron->fZtrk - electron->fZv); fHist.fEle.ZvrtHist->Fill(electron->fZv); // // Caluate the transverse mass Float_t phiLepton = cdfphi(electron->fPhi); float delphi = fabs(phiLepton - pMetBlock->fMetphi[0]); if (delphi > 3.14159) { delphi = 2.0*3.14159 - delphi; } float TransverseMass = sqrt(2.0*pMetBlock->fMet[0]*electron->Et() *(1.0-cos(delphi))); // printf("tm %f et %f met %f elephi %f metphi %f delphi %f \n", // TransverseMass,electron->Et(),pMetBlock->fMet[0], // phiLepton,pMetBlock->fMetphi[0],delphi); // // Histos for W's fHist.fIsoVsMetHist->Fill(pMetBlock->fMet[0], electron->Iso()/electron->Et()); fHist.fEtVsMetHist->Fill(pMetBlock->fMet[0],electron->Et()); if (pMetBlock->fMet[0] > 20.0) { fHist.fWEleEtHist->Fill(electron->Et()); fHist.fWEleMtHist->Fill(TransverseMass); } } } // // Loop over the jets and make cuts for (int i=0; igetgoodJet(i); if (isEnableHist()) { fHist.fJet.EtHist->Fill(jet->Et()); fHist.fJet.EvetaHist->Fill(jet->fEveta); fHist.fJet.DtetaHist->Fill(jet->DetEta()); fHist.fJet.EmfHist->Fill(jet->EmFract()); } } // } if (ientry%1000==0) printf(" Events processed %d passing %d \n", totEvents,passEvents); return 0; } //================================== // EndJob() //================================== int TTopFindModule::EndJob() { printf("----- end job: ---- %s\n",GetName()); printf(" Events processed %d passed %d \n", totEvents,passEvents); return 0; } Float_t TTopFindModule::cdfphi(Float_t phi) { while (phi > 2.0*3.14159) { phi = phi - 2.0*3.14159; } while (phi < 0) { phi = phi + 2.0*3.14159; } return phi; }