Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:10:53

0001 #include "DQM/Physics/src/CentralityDQM.h"
0002 
0003 #include <memory>
0004 
0005 // DQM
0006 #include "DQMServices/Core/interface/DQMStore.h"
0007 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0008 
0009 // Framework
0010 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/Framework/interface/EventSetup.h"
0013 #include "FWCore/Framework/interface/LuminosityBlock.h"
0014 #include "FWCore/ServiceRegistry/interface/Service.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 
0017 // Centrality
0018 #include "FWCore/Framework/interface/EDAnalyzer.h"
0019 #include "DataFormats/VertexReco/interface/Vertex.h"
0020 
0021 using namespace edm;
0022 using namespace reco;
0023 using namespace hi;
0024 
0025 //
0026 // -- Constructor
0027 //
0028 CentralityDQM::CentralityDQM(const edm::ParameterSet& ps) {
0029   edm::LogInfo("CentralityDQM") << " Starting CentralityDQM "
0030                                 << "\n";
0031 
0032   centralityTag_ = ps.getParameter<InputTag>("centralitycollection");
0033   centralityToken = consumes<reco::Centrality>(centralityTag_);
0034 
0035   centralityBinTag_ = (ps.getParameter<edm::InputTag>("centralitybincollection"));
0036   centralityBinToken = consumes<int>(centralityBinTag_);
0037 
0038   vertexTag_ = ps.getParameter<InputTag>("vertexcollection");
0039   vertexToken = consumes<std::vector<reco::Vertex> >(vertexTag_);
0040 
0041   eventplaneTag_ = ps.getParameter<InputTag>("eventplanecollection");
0042   eventplaneToken = consumes<reco::EvtPlaneCollection>(eventplaneTag_);
0043 
0044   // just to initialize
0045 }
0046 
0047 //
0048 // -- Destructor
0049 //
0050 CentralityDQM::~CentralityDQM() {
0051   edm::LogInfo("CentralityDQM") << " Deleting CentralityDQM "
0052                                 << "\n";
0053 }
0054 
0055 //
0056 //  -- Book histograms
0057 //
0058 void CentralityDQM::bookHistograms(DQMStore::IBooker& bei, edm::Run const&, edm::EventSetup const&) {
0059   // void CentralityDQM::bookHistograms(DQMStore* bei){
0060 
0061   bei.setCurrentFolder("Physics/Centrality/");
0062 
0063   h_hiNpix = bei.book1D("h_hiNpix", "h_hiNpix", 750, 0, 75000);
0064   h_hiNpixelTracks = bei.book1D("h_hiNpixelTracks", "hiNpixelTracks", 500, 0, 5000);
0065   h_hiNtracks = bei.book1D("h_hiNtracks", "h_hiNtracks", 500, 0, 5000);
0066   h_hiNtracksPtCut = bei.book1D("h_hiNtracksPtCut", "h_hiNtracksPtCut", 500, 0, 5000);
0067   h_hiNtracksEtaCut = bei.book1D("h_hiNtracksEtaCut", "h_hiNtracksEtaCut", 500, 0, 5000);
0068   h_hiNtracksEtaPtCut = bei.book1D("h_hiNtracksEtaPtCut", "h_hiNtracksEtaPtCut", 500, 0, 5000);
0069 
0070   h_hiHF = bei.book1D("h_hiHF", "h_hiHF", 900, 0, 9000);
0071   h_hiHFplus = bei.book1D("h_hiHFplus", "h_hiHFplus", 900, 0, 9000);
0072   h_hiHFminus = bei.book1D("h_hiHFminus", "h_hiHFminus", 900, 0, 9000);
0073   h_hiHFplusEta4 = bei.book1D("h_hiHFplusEta4", "h_hiHFplusEta4", 900, 0, 9000);
0074   h_hiHFminusEta4 = bei.book1D("h_hiHFminusEta4", "h_hiHFminusEta4", 900, 0, 9000);
0075 
0076   h_hiHFhit = bei.book1D("h_hiHFhit", "h_hiHFhit", 3000, 0, 300000);
0077   h_hiHFhitPlus = bei.book1D("h_hiHFhitPlus", "h_hiHFhitPlus", 2000, 0, 200000);
0078   h_hiHFhitMinus = bei.book1D("h_hiHFhitMinus", "h_hiHFhitMinus", 2000, 0, 200000);
0079 
0080   h_hiEB = bei.book1D("h_hiEB", "h_hiEB", 600, 0, 6000);
0081   h_hiET = bei.book1D("h_hiET", "h_hiET", 600, 0, 6000);
0082   h_hiEE = bei.book1D("h_hiEE", "h_hiEE", 600, 0, 6000);
0083   h_hiEEplus = bei.book1D("h_hiEEplus", "h_hiEEplus", 600, 0, 6000);
0084   h_hiEEminus = bei.book1D("h_hiEEminus", "h_hiEEminus", 600, 0, 6000);
0085   h_hiZDC = bei.book1D("h_hiZDC", "h_hiZDC", 600, 0, 6000);
0086   h_hiZDCplus = bei.book1D("h_hiZDCplus", "h_hiZDCplus", 600, 0, 6000);
0087   h_hiZDCminus = bei.book1D("h_hiZDCminus", "h_hiZDCminus", 600, 0, 6000);
0088 
0089   h_vertex_x = bei.book1D("h_vertex_x", "h_vertex_x", 400, -4, 4);
0090   h_vertex_y = bei.book1D("h_vertex_y", "h_vertex_y", 400, -4, 4);
0091   h_vertex_z = bei.book1D("h_vertex_z", "h_vertex_z", 400, -40, 40);
0092 
0093   h_cent_bin = bei.book1D("h_cent_bin", "h_cent_bin", 200, 0, 200);
0094 
0095   Double_t psirange = 4;
0096   bei.setCurrentFolder("Physics/Centrality/EventPlane/");
0097   h_ep_HFm2 = bei.book1D("h_ep_HFm2", "h_ep_HFm2", 800, -psirange, psirange);
0098   h_ep_HFp2 = bei.book1D("h_ep_HFp2", "h_ep_HFp2", 800, -psirange, psirange);
0099   h_ep_trackmid2 = bei.book1D("h_ep_trackmid2", "h_ep_trackmid2", 800, -psirange, psirange);
0100   h_ep_trackm2 = bei.book1D("h_ep_trackm2", "h_ep_trackm2", 800, -psirange, psirange);
0101   h_ep_trackp2 = bei.book1D("h_ep_trackp2", "h_ep_trackp2", 800, -psirange, psirange);
0102   h_ep_HFm3 = bei.book1D("h_ep_HFm3", "h_ep_HFm3", 800, -psirange, psirange);
0103   h_ep_HFp3 = bei.book1D("h_ep_HFp3", "h_ep_HFp3", 800, -psirange, psirange);
0104   h_ep_trackmid3 = bei.book1D("h_ep_trackmid3", "h_ep_trackmid3", 800, -psirange, psirange);
0105 }
0106 
0107 //
0108 //  -- Analyze
0109 //
0110 void CentralityDQM::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0111   using namespace edm;
0112   edm::Handle<reco::Centrality> cent;
0113   iEvent.getByToken(centralityToken, cent);  //_centralitytag comes from the cfg
0114 
0115   edm::Handle<int> cbin;
0116   iEvent.getByToken(centralityBinToken, cbin);
0117 
0118   // as an inputTag and is
0119   //"hiCentrality"
0120   edm::Handle<reco::EvtPlaneCollection> ep;
0121   iEvent.getByToken(eventplaneToken, ep);
0122 
0123   if (cent.isValid()) {
0124     int hibin = -999;
0125     if (cbin.isValid())
0126       hibin = *cbin;
0127 
0128     //  std::cout<<  " ------------------------------------- "  << hibin << std::endl;
0129 
0130     h_cent_bin->Fill(hibin);
0131     h_hiNpix->Fill(cent->multiplicityPixel());
0132     h_hiNpixelTracks->Fill(cent->NpixelTracks());
0133     h_hiNtracks->Fill(cent->Ntracks());  //
0134 
0135     h_hiNtracksPtCut->Fill(cent->NtracksPtCut());
0136     h_hiNtracksEtaCut->Fill(cent->NtracksEtaCut());
0137     h_hiNtracksEtaPtCut->Fill(cent->NtracksEtaPtCut());
0138 
0139     h_hiHF->Fill(cent->EtHFtowerSum());
0140     h_hiHFplus->Fill(cent->EtHFtowerSumPlus());
0141     h_hiHFminus->Fill(cent->EtHFtowerSumMinus());
0142     h_hiHFplusEta4->Fill(cent->EtHFtruncatedPlus());
0143     h_hiHFminusEta4->Fill(cent->EtHFtruncatedMinus());
0144 
0145     h_hiHFhit->Fill(cent->EtHFhitSum());
0146     h_hiHFhitPlus->Fill(cent->EtHFhitSumPlus());
0147     h_hiHFhitMinus->Fill(cent->EtHFhitSumMinus());
0148 
0149     h_hiZDC->Fill(cent->zdcSum());
0150     h_hiZDCplus->Fill(cent->zdcSumPlus());
0151     h_hiZDCminus->Fill(cent->zdcSumMinus());
0152 
0153     h_hiEEplus->Fill(cent->EtEESumPlus());
0154     h_hiEEminus->Fill(cent->EtEESumMinus());
0155     h_hiEE->Fill(cent->EtEESum());
0156     h_hiEB->Fill(cent->EtEBSum());
0157     h_hiET->Fill(cent->EtMidRapiditySum());
0158 
0159     edm::Handle<std::vector<reco::Vertex> > vertex;
0160     iEvent.getByToken(vertexToken, vertex);
0161     h_vertex_x->Fill(vertex->begin()->x());
0162     h_vertex_y->Fill(vertex->begin()->y());
0163     h_vertex_z->Fill(vertex->begin()->z());
0164   }
0165 
0166   if (ep.isValid()) {
0167     EvtPlaneCollection::const_iterator rp = ep->begin();
0168 
0169     h_ep_HFm2->Fill(rp[HFm2].angle(0));
0170     h_ep_HFp2->Fill(rp[HFp2].angle(0));
0171     h_ep_trackmid2->Fill(rp[trackmid2].angle(0));
0172     h_ep_trackm2->Fill(rp[trackm2].angle(0));
0173     h_ep_trackp2->Fill(rp[trackp2].angle(0));
0174 
0175     h_ep_HFm3->Fill(rp[HFm3].angle(0));
0176     h_ep_HFp3->Fill(rp[HFp3].angle(0));
0177     h_ep_trackmid3->Fill(rp[trackmid3].angle(0));
0178   }
0179 }