Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 11:59:02 

0001 //April 2015 : Removal of itrg1, itrg2, but addition of isect2, same is true in HOCalibVariables.h
0002 // -*- C++ -*-
0003 //
0004 // Package:    HOCalibAnalyzer
0005 // Class:      HOCalibAnalyzer
0006 //
0007 /**\class HOCalibAnalyzer HOCalibAnalyzer.cc Calibration/HOCalibAnalyzer/src/HOCalibAnalyzer.cc
0008 
0009  Description: <one line class summary>
0010 
0011  Implementation:
0012      <Notes on implementation>
0013 
0014 April 2015
0015 // Addition of these variables, ilumi (analyser), pileup (analyser), nprim
0016 
0017 
0018 */
0019 //
0020 // Original Author:  Gobinda Majumder
0021 //         Created:  Sat Jul  7 09:51:31 CEST 2007
0022 //
0023 //
0024 
0025 // system include files
0026 #include <memory>
0027 
0028 // user include files
0029 #include "FWCore/Framework/interface/Frameworkfwd.h"
0030 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0031 
0032 #include "FWCore/Framework/interface/Event.h"
0033 #include "FWCore/Framework/interface/MakerMacros.h"
0034 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0035 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0036 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0037 
0038 #include "FWCore/ServiceRegistry/interface/Service.h"
0039 #include "FWCore/Utilities/interface/InputTag.h"
0040 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0041 #include "DataFormats/HcalCalibObjects/interface/HOCalibVariables.h"
0042 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
0043 #include "DataFormats/Math/interface/angle_units.h"
0044 
0045 #include "TFile.h"
0046 #include "TH1F.h"
0047 #include "TH2F.h"
0048 #include "TTree.h"
0049 #include "TProfile.h"
0050 
0051 #include <iostream>
0052 #include <fstream>
0053 #include <iomanip>
0054 #include <string>
0055 #include <vector>
0056 
0057 //
0058 // class decleration
0059 //
0060 
0061 class HOCalibAnalyzer : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0062 public:
0063   explicit HOCalibAnalyzer(const edm::ParameterSet&);
0064   ~HOCalibAnalyzer() override;
0065 
0066   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0067 
0068 private:
0069   void beginJob() override {}
0070   void analyze(const edm::Event&, const edm::EventSetup&) override;
0071   void endJob() override {}
0072 
0073   static constexpr int netamx = 30;
0074   static constexpr int nphimx = 72;
0075   static constexpr int ringmx = 5;
0076   static constexpr int ncut = 14;
0077 
0078   const char* varcrit[3] = {"All", "steps", "n-1"};  // or opposite
0079 
0080   const double elosfact = (14.9 + 0.96 * std::fabs(std::log(8 * 2.8)) + 0.033 * 8 * (1.0 - std::pow(8, -0.33)));
0081 
0082   int getHOieta(int ij) { return (ij < netamx / 2) ? -netamx / 2 + ij : -netamx / 2 + ij + 1; }
0083   int invert_HOieta(int ieta) { return (ieta < 0) ? netamx / 2 + ieta : netamx / 2 + ieta - 1; }
0084 
0085   int mypow_2[31];
0086 
0087   const bool m_cosmic;
0088   const bool m_zeroField;
0089   const int m_bins;
0090   const double m_low;
0091   const double m_ahigh;
0092   const bool m_histFill;  //Stored signals of individual HO towers with default selection criteria
0093   const bool m_treeFill;  //Store rootuple without almost any selection criteria except a quality on muon
0094   const bool m_verbose;
0095 
0096   int ipass;
0097 
0098   TTree* T1;
0099 
0100   TH1F* ho_indenergy[netamx][nphimx];
0101 
0102   TH1F* muonnm;
0103   TH1F* muonmm;
0104   TH1F* muonth;
0105   TH1F* muonph;
0106   TH1F* muonch;
0107 
0108   TH1F* sel_muonnm;
0109   TH1F* sel_muonmm;
0110   TH1F* sel_muonth;
0111   TH1F* sel_muonph;
0112   TH1F* sel_muonch;
0113 
0114   TH2F* sig_eta_evt[3 * netamx][ncut];  //For individual eta
0115   TH2F* sigvsevt[3 * netamx][ncut];
0116   TH1F* variab[3 * netamx][ncut];
0117 
0118   TH2F* mu_projection[ncut + 1];
0119 
0120   unsigned ievt, hoflag;
0121   int irun, ilumi, nprim, isect, isect2, ndof, nmuon;
0122 
0123   float pileup, trkdr, trkdz, trkvx, trkvy, trkvz, trkmm, trkth, trkph, chisq, therr, pherr, hodx, hody, hoang, htime,
0124       hosig[9], hocorsig[18], hocro, hbhesig[9], caloen[3];
0125   float momatho, tkpt03, ecal03, hcal03;
0126   float tmphoang;
0127 
0128   int nevents[10];
0129 
0130   float ncount[ringmx][ncut + 10];
0131 
0132   edm::InputTag hoCalibVariableCollectionTag;
0133   const edm::EDGetTokenT<HOCalibVariableCollection> tok_ho_;
0134   const edm::EDGetTokenT<HORecHitCollection> tok_allho_;
0135   // ----------member data ---------------------------
0136 };
0137 
0138 //
0139 // constants, enums and typedefs
0140 //
0141 
0142 //
0143 // static data member definitions
0144 //
0145 
0146 //
0147 // constructors and destructor
0148 //
0149 
0150 HOCalibAnalyzer::HOCalibAnalyzer(const edm::ParameterSet& iConfig)
0151     : m_cosmic(iConfig.getUntrackedParameter<bool>("cosmic", true)),
0152       m_zeroField(iConfig.getUntrackedParameter<bool>("zeroField", false)),
0153       m_bins(iConfig.getUntrackedParameter<int>("HOSignalBins", 120)),
0154       m_low(iConfig.getUntrackedParameter<double>("lowerRange", -1.0)),
0155       m_ahigh(iConfig.getUntrackedParameter<double>("upperRange", 29.0)),
0156       m_histFill(iConfig.getUntrackedParameter<bool>("histFill", true)),
0157       m_treeFill(iConfig.getUntrackedParameter<bool>("treeFill", false)),
0158       m_verbose(iConfig.getUntrackedParameter<bool>("verbose", false)),
0159       tok_ho_(consumes<HOCalibVariableCollection>(iConfig.getParameter<edm::InputTag>("hoCalibVariableCollectionTag"))),
0160       tok_allho_(consumes<HORecHitCollection>(iConfig.getParameter<edm::InputTag>("hoInputTag"))) {
0161   // It is very likely you want the following in your configuration
0162   // hoCalibVariableCollectionTag = cms.InputTag('hoCalibProducer', 'HOCalibVariableCollection')
0163 
0164   usesResource(TFileService::kSharedResource);
0165 
0166   //now do what ever initialization is needed
0167   ipass = 0;
0168   for (int ij = 0; ij < 10; ij++) {
0169     nevents[ij] = 0;
0170   }
0171 
0172   edm::Service<TFileService> fs;
0173 
0174   T1 = fs->make<TTree>("T1", "HOSignal");
0175 
0176   T1->Branch("irun", &irun, "irun/I");
0177   T1->Branch("ievt", &ievt, "ievt/i");
0178 
0179   T1->Branch("isect", &isect, "isect/I");
0180   T1->Branch("isect2", &isect2, "isect2/I");
0181   T1->Branch("ndof", &ndof, "ndof/I");
0182   T1->Branch("nmuon", &nmuon, "nmuon/I");
0183 
0184   T1->Branch("ilumi", &ilumi, "ilumi/I");
0185   if (!m_cosmic) {
0186     T1->Branch("pileup", &pileup, "pileup/F");
0187     T1->Branch("nprim", &nprim, "nprim/I");
0188     T1->Branch("tkpt03", &tkpt03, " tkpt03/F");
0189     T1->Branch("ecal03", &ecal03, " ecal03/F");
0190     T1->Branch("hcal03", &hcal03, " hcal03/F");
0191   }
0192 
0193   T1->Branch("trkdr", &trkdr, "trkdr/F");
0194   T1->Branch("trkdz", &trkdz, "trkdz/F");
0195 
0196   T1->Branch("trkvx", &trkvx, "trkvx/F");
0197   T1->Branch("trkvy", &trkvy, "trkvy/F");
0198   T1->Branch("trkvz", &trkvz, "trkvz/F");
0199   T1->Branch("trkmm", &trkmm, "trkmm/F");
0200   T1->Branch("trkth", &trkth, "trkth/F");
0201   T1->Branch("trkph", &trkph, "trkph/F");
0202 
0203   T1->Branch("chisq", &chisq, "chisq/F");
0204   T1->Branch("therr", &therr, "therr/F");
0205   T1->Branch("pherr", &pherr, "pherr/F");
0206   T1->Branch("hodx", &hodx, "hodx/F");
0207   T1->Branch("hody", &hody, "hody/F");
0208   T1->Branch("hoang", &hoang, "hoang/F");
0209 
0210   T1->Branch("momatho", &momatho, "momatho/F");
0211   T1->Branch("hoflag", &hoflag, "hoflag/i");
0212   T1->Branch("htime", &htime, "htime/F");
0213   T1->Branch("hosig", hosig, "hosig[9]/F");
0214   T1->Branch("hocro", &hocro, "hocro/F");
0215   T1->Branch("hocorsig", hocorsig, "hocorsig[18]/F");
0216   T1->Branch("caloen", caloen, "caloen[3]/F");
0217 
0218   char name[200];
0219   char title[200];
0220 
0221   if (m_histFill) {
0222     for (int ij = 0; ij < netamx; ij++) {
0223       int ieta = getHOieta(ij);
0224       for (int jk = 0; jk < nphimx; jk++) {
0225         sprintf(name, "ho_indenergy_%i_%i", ij, jk);
0226         sprintf(title, "ho IndEnergy (GeV) i#eta=%i i#phi=%i", ieta, jk + 1);
0227         ho_indenergy[ij][jk] = fs->make<TH1F>(name, title, 1200, m_low, m_ahigh);
0228       }
0229     }
0230   }
0231 
0232   muonnm = fs->make<TH1F>("muonnm", "No of muon", 10, -0.5, 9.5);
0233   muonmm = fs->make<TH1F>("muonmm", "P_{mu}", 200, -100., 100.);
0234   muonth = fs->make<TH1F>("muonth", "{Theta}_{mu}", 180, 0., 180.);
0235   muonph = fs->make<TH1F>("muonph", "{Phi}_{mu}", 180, -180., 180.);
0236   muonch = fs->make<TH1F>("muonch", "{chi^2}/ndf", 100, 0., 1000.);
0237 
0238   sel_muonnm = fs->make<TH1F>("sel_muonnm", "No of muon(sel)", 10, -0.5, 9.5);
0239   sel_muonmm = fs->make<TH1F>("sel_muonmm", "P_{mu}(sel)", 200, -100., 100.);
0240   sel_muonth = fs->make<TH1F>("sel_muonth", "{Theta}_{mu}(sel)", 180, 0., 180.);
0241   sel_muonph = fs->make<TH1F>("sel_muonph", "{Phi}_{mu}(sel)", 180, -180., 180.);
0242   sel_muonch = fs->make<TH1F>("sel_muonch", "{chi^2}/ndf(sel)", 100, 0., 1000.);
0243 
0244   //if change order, change in iselect_wotime also and other efficiency numbers
0245   const char* varnam[ncut] = {"ndof",
0246                               "chisq",
0247                               "th",
0248                               "ph",
0249                               "therr",
0250                               "pherr",
0251                               "dircos",
0252                               "trkmm",
0253                               "nmuon",
0254                               "calo",
0255                               "trkiso",
0256                               "#phi-dir",
0257                               "#eta-dir",
0258                               "time"};
0259   int nbinxx[ncut] = {25, 60, 60, 60, 60, 60, 60, 120, 6, 60, 60, 120, 120, 60};
0260   double alowxx[ncut] = {5.5, 0., 0., -angle_units::piRadians, 0.0, 0.0, 0.0, 0., 0.5, 0.0, 0.0, -20., -32., -45.0};
0261   double ahghxx[ncut] = {
0262       30.5, 40., angle_units::piRadians, angle_units::piRadians, 0.8, 0.02, 0.5, 300., 6.5, 10.0, 24.0, 20.0, 32.0, 45.0};
0263 
0264   for (int kl = 0; kl < ncut; kl++) {
0265     for (int jk = 0; jk < 3; jk++) {
0266       for (int ij = 0; ij < netamx; ij++) {
0267         sprintf(name, "sigeta_%i_%i_%i", kl, jk, ij);
0268         sprintf(title, "sigeta %s %s i#eta=%i", varnam[kl], varcrit[jk], getHOieta(ij));
0269         sig_eta_evt[netamx * jk + ij][kl] =
0270             fs->make<TH2F>(name, title, nbinxx[kl], alowxx[kl], ahghxx[kl], m_bins, m_low, m_ahigh);
0271       }
0272     }
0273   }
0274 
0275   for (int kl = 0; kl < ncut; kl++) {
0276     for (int ij = 0; ij < ringmx * 3; ij++) {
0277       int iring = ij % ringmx - 2;
0278       int iset = ij / ringmx;
0279       sprintf(name, "sigring_%i_%i", kl, ij);
0280       sprintf(title, "Signal %s %s Ring%i", varnam[kl], varcrit[iset], iring);
0281       sigvsevt[ij][kl] = fs->make<TH2F>(name, title, nbinxx[kl], alowxx[kl], ahghxx[kl], m_bins, m_low, m_ahigh);
0282     }
0283   }
0284 
0285   for (int kl = 0; kl < ncut; kl++) {
0286     for (int ij = 0; ij < ringmx * 3; ij++) {
0287       int iring = ij % ringmx - 2;
0288       int iset = ij / ringmx;
0289       sprintf(name, "varring_%i_%i", kl, ij);
0290       sprintf(title, "%s %s Ring%i", varnam[kl], varcrit[iset], iring);
0291       variab[ij][kl] = fs->make<TH1F>(name, title, nbinxx[kl], alowxx[kl], ahghxx[kl]);
0292     }
0293   }
0294 
0295   for (int ij = 0; ij <= ncut; ij++) {
0296     sprintf(name, "mu_projection_%i", ij);
0297     if (ij == 0) {
0298       sprintf(title, "All projected muon");
0299     } else {
0300       sprintf(title, "Projected muon with selection %s", varnam[ij - 1]);
0301     }
0302     mu_projection[ij] =
0303         fs->make<TH2F>(name, title, netamx + 1, -netamx / 2 - 0.5, netamx / 2 + 0.5, nphimx, 0.5, nphimx + 0.5);
0304   }
0305 
0306   for (int ij = 0; ij < 31; ij++) {
0307     mypow_2[ij] = std::pow(2, ij);
0308   }
0309   for (int ij = 0; ij < ringmx; ij++) {
0310     for (int jk = 0; jk < ncut + 10; jk++) {
0311       ncount[ij][jk] = 0.0;
0312     }
0313   }
0314 }
0315 
0316 HOCalibAnalyzer::~HOCalibAnalyzer() {
0317   edm::LogVerbatim("HOCalibAnalyzer") << " Total events = " << std::setw(7) << nevents[0] << " " << std::setw(7)
0318                                       << nevents[1] << " " << std::setw(7) << nevents[2] << " " << std::setw(7)
0319                                       << nevents[3] << " " << std::setw(7) << nevents[4] << " " << std::setw(7)
0320                                       << nevents[5] << " Selected events # is " << ipass;
0321 }
0322 
0323 //
0324 // member functions
0325 //
0326 void HOCalibAnalyzer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0327   edm::ParameterSetDescription desc;
0328   desc.add<edm::InputTag>("hoCalibVariableCollectionTag",
0329                           edm::InputTag("hoCalibProducer", "HOCalibVariableCollection"));
0330   desc.add<edm::InputTag>("hoInputTag", edm::InputTag("horeco"));
0331   desc.addUntracked<bool>("cosmic", true);
0332   desc.addUntracked<bool>("zeroField", false);
0333   desc.addUntracked<int>("HOSignalBins", 120);
0334   desc.addUntracked<double>("lowerRange", -1.0);
0335   desc.addUntracked<double>("upperRange", 29.0);
0336   desc.addUntracked<bool>("histFill", true);
0337   desc.addUntracked<bool>("treeFill", false);
0338   desc.addUntracked<double>("sigma", 0.05);
0339   desc.addUntracked<bool>("verbose", false);
0340   descriptions.add("hoCalibAnalyzer", desc);
0341 }
0342 
0343 // ------------ method called to for each event  ------------
0344 void HOCalibAnalyzer::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0345   nevents[0]++;
0346 
0347   ievt = iEvent.id().event();
0348   ilumi = iEvent.luminosityBlock();
0349 
0350   const edm::Handle<HOCalibVariableCollection>& HOCalib = iEvent.getHandle(tok_ho_);
0351 
0352   if (nevents[0] % 20000 == 1) {
0353     edm::LogVerbatim("HOCalibAnalyzer") << "nmuon event # " << std::setw(7) << nevents[0] << " " << std::setw(7)
0354                                         << nevents[1] << " " << std::setw(7) << nevents[2] << " " << std::setw(7)
0355                                         << nevents[3] << " " << std::setw(7) << nevents[4] << " " << std::setw(7)
0356                                         << nevents[5];
0357     edm::LogVerbatim("HOCalibAnalyzer") << " Run # " << iEvent.id().run() << " Evt # " << iEvent.id().event() << " "
0358                                         << int(HOCalib.isValid()) << " " << ipass;
0359   }
0360 
0361   if (HOCalib.isValid()) {
0362     nevents[1]++;
0363     nmuon = (*HOCalib).size();
0364 
0365     for (HOCalibVariableCollection::const_iterator hoC = (*HOCalib).begin(); hoC != (*HOCalib).end(); hoC++) {
0366       trkdr = (*hoC).trkdr;
0367       trkdz = (*hoC).trkdz;
0368 
0369       trkvx = (*hoC).trkvx;
0370       trkvy = (*hoC).trkvy;
0371       trkvz = (*hoC).trkvz;
0372 
0373       trkmm = (*hoC).trkmm;
0374       trkth = (*hoC).trkth;
0375       trkph = (*hoC).trkph;
0376 
0377       ndof = static_cast<int>((*hoC).ndof);
0378       chisq = (*hoC).chisq;
0379       momatho = (*hoC).momatho;
0380 
0381       therr = (*hoC).therr;
0382       pherr = (*hoC).pherr;
0383       trkph = (*hoC).trkph;
0384 
0385       if (!m_cosmic) {
0386         nprim = (*hoC).nprim;
0387         pileup = (*hoC).pileup;
0388         tkpt03 = (*hoC).tkpt03;
0389         ecal03 = (*hoC).ecal03;
0390         hcal03 = (*hoC).hcal03;
0391       }
0392 
0393       isect = (*hoC).isect;
0394       isect2 = (*hoC).isect2;
0395       hodx = (*hoC).hodx;
0396       hody = (*hoC).hody;
0397       hoang = (*hoC).hoang;
0398 
0399       tmphoang = std::sin(trkth) - hoang;
0400 
0401       htime = (*hoC).htime;
0402       hoflag = (*hoC).hoflag;
0403       for (int ij = 0; ij < 9; ij++) {
0404         hosig[ij] = (*hoC).hosig[ij];
0405         if (m_verbose)
0406           edm::LogVerbatim("HOCalibAnalyzer") << "hosig " << ij << " " << hosig[ij];
0407       }
0408       for (int ij = 0; ij < 18; ij++) {
0409         hocorsig[ij] = (*hoC).hocorsig[ij];
0410         if (m_verbose)
0411           edm::LogVerbatim("HOCalibAnalyzer") << "hocorsig " << ij << " " << hocorsig[ij];
0412       }
0413       hocro = (*hoC).hocro;
0414       for (int ij = 0; ij < 3; ij++) {
0415         caloen[ij] = (*hoC).caloen[ij];
0416       }
0417 
0418       int ipsall = 0;
0419       int ips0 = 0;
0420       int ips1 = 0;
0421       int ips2 = 0;
0422       int ips3 = 0;
0423       int ips4 = 0;
0424       int ips5 = 0;
0425       int ips6 = 0;
0426       int ips7 = 0;
0427       int ips8 = 0;
0428       int ips9 = 0;
0429       int ips10 = 0;
0430       int ips11 = 0;
0431       int ips12 = 0;
0432       int ips13 = 0;
0433 
0434       nevents[2]++;
0435       bool isZSps = (hosig[4] < -99.0) ? false : true;
0436 
0437       if ((!m_cosmic) && std::fabs(trkmm) < momatho)
0438         continue;
0439 
0440       nevents[3]++;
0441       if (std::fabs(trkth - angle_units::piRadians / 2) < 0.000001)
0442         continue;  //22OCT07
0443       nevents[4]++;
0444 
0445       int ieta = int((std::abs(isect) % 10000) / 100.) - 50;  //an offset to acodate -ve eta values
0446       if (std::abs(ieta) >= 16)
0447         continue;
0448       nevents[5]++;
0449       int iphi = std::abs(isect) % 100;
0450 
0451       int iring = 0;
0452 
0453       int iring2 = iring + 2;
0454 
0455       double abshoang = (m_cosmic) ? std::fabs(hoang) : hoang;
0456 
0457       double elos = 1.0 / std::max(0.1, std::abs(static_cast<double>(hoang)));
0458 
0459       if (!m_zeroField)
0460         elos *= ((14.9 + 0.96 * std::fabs(log(momatho * 2.8)) + 0.033 * momatho * (1.0 - std::pow(momatho, -0.33))) /
0461                  elosfact);
0462 
0463       if (m_cosmic) {
0464         if (std::abs(ndof) >= 20 && std::abs(ndof) < 55) {
0465           ips0 = mypow_2[0];
0466           ipsall += ips0;
0467         }
0468         if (chisq > 0 && chisq < 12) {
0469           ips1 = mypow_2[1];
0470           ipsall += ips1;
0471         }  //18Jan2008
0472 
0473         if (trkth > 0.3 && trkth < angle_units::piRadians - 0.3) {
0474           ips2 = mypow_2[2];
0475           ipsall += ips2;
0476         }  //No nead for pp evt
0477         if (trkph > -angle_units::piRadians + 0.1 && trkph < -0.1) {
0478           ips3 = mypow_2[3];
0479           ipsall += ips3;
0480         }  //No nead for pp evt
0481 
0482         if (therr < 0.02) {
0483           ips4 = mypow_2[4];
0484           ipsall += ips4;
0485         }
0486         if (pherr < 0.0002) {
0487           ips5 = mypow_2[5];
0488           ipsall += ips5;
0489         }
0490         if (abshoang > 0.60 && abshoang < 1.0) {
0491           ips6 = mypow_2[6];
0492           ipsall += ips6;
0493         }
0494 
0495         if (m_zeroField || (std::fabs(momatho) > 5.0 && std::fabs(momatho) < 2000.0)) {
0496           ips7 = mypow_2[7];
0497           ipsall += ips7;
0498         }
0499 
0500         if (nmuon >= 1 && nmuon <= 3) {
0501           ips8 = mypow_2[8];
0502           ipsall += ips8;
0503         }
0504 
0505         // initially for: if (hodx>0 && hody>0) { }
0506         ips9 = mypow_2[9];
0507         ipsall += ips9;
0508 
0509         ips10 = mypow_2[10];
0510         ipsall += ips10;
0511 
0512         if (iring2 == 2) {
0513           if (std::fabs(hodx) < 100 && std::fabs(hodx) > 2 && std::fabs(hocorsig[8]) < 40 &&
0514               std::fabs(hocorsig[8]) > 2) {
0515             ips11 = mypow_2[11];
0516             ipsall += ips11;
0517           }
0518 
0519           if (std::fabs(hody) < 100 && std::fabs(hody) > 2 && std::fabs(hocorsig[9]) < 40 &&
0520               std::fabs(hocorsig[9]) > 2) {
0521             ips12 = mypow_2[12];
0522             ipsall += ips12;
0523           }
0524 
0525         } else {
0526           if (std::fabs(hodx) < 100 && std::fabs(hodx) > 2) {
0527             ips11 = mypow_2[11];
0528             ipsall += ips11;
0529           }
0530 
0531           if (std::fabs(hody) < 100 && std::fabs(hody) > 2) {
0532             ips12 = mypow_2[12];
0533             ipsall += ips12;
0534           }
0535         }
0536 
0537         if (m_zeroField) {
0538           if (iring2 == 0) {
0539             if (htime > -60 && htime < 60) {
0540               ips13 = mypow_2[13];
0541               ipsall += ips13;
0542             }
0543           }
0544           if (iring2 == 1) {
0545             if (htime > -60 && htime < 60) {
0546               ips13 = mypow_2[13];
0547               ipsall += ips13;
0548             }
0549           }
0550           if (iring2 == 2) {
0551             if (htime > -60 && htime < 60) {
0552               ips13 = mypow_2[13];
0553               ipsall += ips13;
0554             }
0555           }
0556           if (iring2 == 3) {
0557             if (htime > -60 && htime < 60) {
0558               ips13 = mypow_2[13];
0559               ipsall += ips13;
0560             }
0561           }
0562           if (iring2 == 4) {
0563             if (htime > -60 && htime < 60) {
0564               ips13 = mypow_2[13];
0565               ipsall += ips13;
0566             }
0567           }
0568         } else {
0569           if (htime > -100 && htime < 100) {
0570             ips13 = mypow_2[13];
0571             ipsall += ips13;
0572           }
0573         }
0574       } else {
0575         if (std::abs(ndof) >= 10 && std::abs(ndof) < 25) {
0576           ips0 = mypow_2[0];
0577           ipsall += ips0;
0578         }
0579         if (chisq > 0 && chisq < 10) {
0580           ips1 = mypow_2[1];
0581           ipsall += ips1;
0582         }  //18Jan2008
0583 
0584         if (std::fabs(trkth - angle_units::piRadians / 2) < 21.5) {
0585           ips2 = mypow_2[2];
0586           ipsall += ips2;
0587         }  //No nead for pp evt
0588         if (std::fabs(trkph + angle_units::piRadians / 2) < 21.5) {
0589           ips3 = mypow_2[3];
0590           ipsall += ips3;
0591         }  //No nead for pp evt
0592 
0593         if (therr < 0.00002) {
0594           ips4 = mypow_2[4];
0595           ipsall += ips4;
0596         }
0597         if (pherr < 0.000002) {
0598           ips5 = mypow_2[5];
0599           ipsall += ips5;
0600         }
0601         // earlier: if (abshoang >0.40 && abshoang <1.0) {ips6 = mypow_2[6];  ipsall +=ips6;}
0602         if (tmphoang < 0.065) {
0603           ips6 = mypow_2[6];
0604           ipsall += ips6;
0605         }
0606 
0607         if (std::fabs(momatho) < 250.0 && std::fabs(momatho) > 15.0) {
0608           if (iring2 == 2) {
0609             ips7 = mypow_2[7];
0610             ipsall += ips7;
0611           }
0612           if ((iring2 == 1 || iring2 == 3) && std::fabs(momatho) > 17.0) {
0613             ips7 = mypow_2[7];
0614             ipsall += ips7;
0615           }
0616           if ((iring2 == 0 || iring2 == 4) && std::fabs(momatho) > 20.0) {
0617             ips7 = mypow_2[7];
0618             ipsall += ips7;
0619           }
0620         }
0621 
0622         if (nmuon >= 1 && nmuon <= 3) {
0623           ips8 = mypow_2[8];
0624           ipsall += ips8;
0625         }
0626 
0627         if (ndof > 0 && caloen[0] < 15.0) {
0628           ips9 = mypow_2[9];
0629           ipsall += ips9;
0630         }  //5.0
0631         if (tkpt03 < 5.0) {
0632           ips10 = mypow_2[10];
0633           ipsall += ips10;
0634         }  //4.0
0635 
0636         if (iring2 == 2) {
0637           if (std::fabs(hodx) < 100 && std::fabs(hodx) > 2 && std::fabs(hocorsig[8]) < 40 &&
0638               std::fabs(hocorsig[8]) > 2) {
0639             ips11 = mypow_2[11];
0640             ipsall += ips11;
0641           }
0642 
0643           if (std::fabs(hody) < 100 && std::fabs(hody) > 2 && std::fabs(hocorsig[9]) < 40 &&
0644               std::fabs(hocorsig[9]) > 2) {
0645             ips12 = mypow_2[12];
0646             ipsall += ips12;
0647           }
0648 
0649         } else {
0650           if (std::fabs(hodx) < 100 && std::fabs(hodx) > 2) {
0651             ips11 = mypow_2[11];
0652             ipsall += ips11;
0653           }
0654 
0655           if (std::fabs(hody) < 100 && std::fabs(hody) > 2) {
0656             ips12 = mypow_2[12];
0657             ipsall += ips12;
0658           }
0659         }
0660 
0661         if (iring2 == 0) {
0662           if (htime > -20 && htime < 20) {
0663             ips13 = mypow_2[13];
0664             ipsall += ips13;
0665           }
0666         }
0667         if (iring2 == 1) {
0668           if (htime > -20 && htime < 20) {
0669             ips13 = mypow_2[13];
0670             ipsall += ips13;
0671           }
0672         }
0673         if (iring2 == 2) {
0674           if (htime > -30 && htime < 20) {
0675             ips13 = mypow_2[13];
0676             ipsall += ips13;
0677           }
0678         }
0679         if (iring2 == 3) {
0680           if (htime > -20 && htime < 20) {
0681             ips13 = mypow_2[13];
0682             ipsall += ips13;
0683           }
0684         }
0685         if (iring2 == 4) {
0686           if (htime > -20 && htime < 20) {
0687             ips13 = mypow_2[13];
0688             ipsall += ips13;
0689           }
0690         }
0691       }
0692       int tmpxet = invert_HOieta(ieta);
0693       double nomHOSig = hosig[4] / elos;
0694 
0695       if (ipsall - ips0 == mypow_2[ncut] - mypow_2[0] - 1) {
0696         if (isZSps) {
0697           sigvsevt[iring2][0]->Fill(std::abs(ndof), nomHOSig);
0698           sig_eta_evt[tmpxet][0]->Fill(std::abs(ndof), nomHOSig);
0699         }
0700         variab[iring2][0]->Fill(std::abs(ndof));
0701       }
0702       if (ipsall - ips1 == mypow_2[ncut] - mypow_2[1] - 1) {
0703         if (isZSps) {
0704           sigvsevt[iring2][1]->Fill(chisq, nomHOSig);
0705           sig_eta_evt[tmpxet][1]->Fill(chisq, nomHOSig);
0706         }
0707         variab[iring2][1]->Fill(chisq);
0708       }
0709       if (ipsall - ips2 == mypow_2[ncut] - mypow_2[2] - 1) {
0710         if (isZSps) {
0711           sigvsevt[iring2][2]->Fill(trkth, nomHOSig);
0712           sig_eta_evt[tmpxet][2]->Fill(trkth, nomHOSig);
0713         }
0714         variab[iring2][2]->Fill(trkth);
0715       }
0716       if (ipsall - ips3 == mypow_2[ncut] - mypow_2[3] - 1) {
0717         if (isZSps) {
0718           sigvsevt[iring2][3]->Fill(trkph, nomHOSig);
0719           sig_eta_evt[tmpxet][3]->Fill(trkph, nomHOSig);
0720         }
0721         variab[iring2][3]->Fill(trkph);
0722       }
0723       if (ipsall - ips4 == mypow_2[ncut] - mypow_2[4] - 1) {
0724         if (isZSps) {
0725           sigvsevt[iring2][4]->Fill(1000 * therr, nomHOSig);
0726           sig_eta_evt[tmpxet][4]->Fill(1000 * therr, nomHOSig);
0727         }
0728         variab[iring2][4]->Fill(1000 * therr);
0729       }
0730       if (ipsall - ips5 == mypow_2[ncut] - mypow_2[5] - 1) {
0731         if (isZSps) {
0732           sigvsevt[iring2][5]->Fill(1000 * pherr, nomHOSig);
0733           sig_eta_evt[tmpxet][5]->Fill(1000 * pherr, nomHOSig);
0734         }
0735         variab[iring2][5]->Fill(1000 * pherr);
0736       }
0737       if (ipsall - ips6 == mypow_2[ncut] - mypow_2[6] - 1) {
0738         if (isZSps) {
0739           sigvsevt[iring2][6]->Fill(tmphoang, (nomHOSig)*abshoang);
0740           sig_eta_evt[tmpxet][6]->Fill(tmphoang, (nomHOSig)*abshoang);
0741         }
0742         variab[iring2][6]->Fill(tmphoang);
0743       }
0744       if (ipsall - ips7 == mypow_2[ncut] - mypow_2[7] - 1) {
0745         if (isZSps) {
0746           sigvsevt[iring2][7]->Fill(std::fabs(trkmm), nomHOSig);
0747           sig_eta_evt[tmpxet][7]->Fill(std::fabs(trkmm), nomHOSig);
0748         }
0749         variab[iring2][7]->Fill(std::fabs(trkmm));
0750       }
0751       if (ipsall - ips8 == mypow_2[ncut] - mypow_2[8] - 1) {
0752         if (isZSps) {
0753           sigvsevt[iring2][8]->Fill(nmuon, nomHOSig);
0754           sig_eta_evt[tmpxet][8]->Fill(nmuon, nomHOSig);
0755         }
0756         variab[iring2][8]->Fill(nmuon);
0757       }
0758       if (!m_cosmic) {
0759         if (ipsall - ips9 == mypow_2[ncut] - mypow_2[9] - 1) {
0760           if (isZSps) {
0761             sigvsevt[iring2][9]->Fill(caloen[0], nomHOSig);
0762             sig_eta_evt[tmpxet][9]->Fill(caloen[0], nomHOSig);
0763           }
0764           variab[iring2][9]->Fill(caloen[0]);
0765         }
0766       }
0767       if (ipsall - ips10 == mypow_2[ncut] - mypow_2[10] - 1) {
0768         if (isZSps) {
0769           sigvsevt[iring2][10]->Fill(tkpt03, nomHOSig);
0770           sig_eta_evt[tmpxet][10]->Fill(tkpt03, nomHOSig);
0771         }
0772         variab[iring2][10]->Fill(tkpt03);
0773       }
0774       if (ipsall - ips11 == mypow_2[ncut] - mypow_2[11] - 1) {
0775         if (isZSps) {
0776           sigvsevt[iring2][11]->Fill(hodx, nomHOSig);
0777           sig_eta_evt[tmpxet][11]->Fill(hodx, nomHOSig);
0778         }
0779         variab[iring2][11]->Fill(hodx);
0780       }
0781       if (ipsall - ips12 == mypow_2[ncut] - mypow_2[12] - 1) {
0782         if (isZSps) {
0783           sigvsevt[iring2][12]->Fill(hody, nomHOSig);
0784           sig_eta_evt[tmpxet][12]->Fill(hody, nomHOSig);
0785         }
0786         variab[iring2][12]->Fill(hody);
0787       }
0788 
0789       if (ipsall - ips13 == mypow_2[ncut] - mypow_2[13] - 1) {
0790         if (isZSps) {
0791           sigvsevt[iring2][13]->Fill(htime, nomHOSig);
0792           sig_eta_evt[tmpxet][13]->Fill(htime, nomHOSig);
0793         }
0794         variab[iring2][13]->Fill(htime);
0795       }
0796 
0797       if (isZSps) {
0798         sigvsevt[iring2 + ringmx][0]->Fill(std::abs(ndof), nomHOSig);
0799         sig_eta_evt[netamx + tmpxet][0]->Fill(std::abs(ndof), nomHOSig);
0800       }
0801       variab[iring2 + 5][0]->Fill(std::abs(ndof));
0802 
0803       ncount[iring2][0]++;
0804       if (isZSps) {
0805         ncount[iring2][1]++;
0806       }
0807       if (ips0 > 0) {
0808         if (isZSps) {
0809           ncount[iring2][10]++;
0810           sigvsevt[iring2 + ringmx][1]->Fill(chisq, nomHOSig);
0811           sig_eta_evt[netamx + tmpxet][1]->Fill(chisq, nomHOSig);
0812         }
0813         variab[iring2 + ringmx][1]->Fill(chisq);
0814         mu_projection[1]->Fill(ieta, iphi);
0815         if (ips1 > 0) {
0816           if (isZSps) {
0817             ncount[iring2][11]++;
0818             sigvsevt[iring2 + ringmx][2]->Fill(trkth, nomHOSig);
0819             sig_eta_evt[netamx + tmpxet][2]->Fill(trkth, nomHOSig);
0820           }
0821           variab[iring2 + ringmx][2]->Fill(trkth);
0822           mu_projection[2]->Fill(ieta, iphi);
0823           if (ips2 > 0) {
0824             if (isZSps) {
0825               ncount[iring2][12]++;
0826               sigvsevt[iring2 + ringmx][3]->Fill(trkph, nomHOSig);
0827               sig_eta_evt[netamx + tmpxet][3]->Fill(trkph, nomHOSig);
0828             }
0829             variab[iring2 + ringmx][3]->Fill(trkph);
0830             mu_projection[3]->Fill(ieta, iphi);
0831             if (ips3 > 0) {
0832               if (isZSps) {
0833                 ncount[iring2][13]++;
0834                 sigvsevt[iring2 + ringmx][4]->Fill(1000 * therr, nomHOSig);
0835                 sig_eta_evt[netamx + tmpxet][4]->Fill(1000 * therr, nomHOSig);
0836               }
0837               variab[iring2 + ringmx][4]->Fill(1000 * therr);
0838               mu_projection[4]->Fill(ieta, iphi);
0839               if (ips4 > 0) {
0840                 if (isZSps) {
0841                   ncount[iring2][14]++;
0842                   sigvsevt[iring2 + ringmx][5]->Fill(1000 * pherr, nomHOSig);
0843                   sig_eta_evt[netamx + tmpxet][5]->Fill(1000 * pherr, nomHOSig);
0844                 }
0845                 variab[iring2 + ringmx][5]->Fill(1000 * pherr);
0846                 mu_projection[5]->Fill(ieta, iphi);
0847                 if (ips5 > 0) {
0848                   if (isZSps) {
0849                     ncount[iring2][15]++;
0850                     sigvsevt[iring2 + ringmx][6]->Fill(tmphoang, (nomHOSig)*abshoang);
0851                     sig_eta_evt[netamx + tmpxet][6]->Fill(tmphoang, (nomHOSig)*abshoang);
0852                   }
0853                   variab[iring2 + ringmx][6]->Fill(tmphoang);
0854                   mu_projection[6]->Fill(ieta, iphi);
0855                   if (ips6 > 0) {
0856                     if (isZSps) {
0857                       ncount[iring2][16]++;
0858                       sigvsevt[iring2 + ringmx][7]->Fill(std::fabs(trkmm), nomHOSig);
0859                       sig_eta_evt[netamx + tmpxet][7]->Fill(std::fabs(trkmm), nomHOSig);
0860                     }
0861                     variab[iring2 + ringmx][7]->Fill(std::fabs(trkmm));
0862                     mu_projection[7]->Fill(ieta, iphi);
0863                     if (ips7 > 0) {
0864                       ncount[iring2][4]++;  //Efficiency of Muon detection
0865                       if (isZSps) {
0866                         ncount[iring2][17]++;
0867                         sigvsevt[iring2 + ringmx][8]->Fill(nmuon, nomHOSig);
0868                         sig_eta_evt[netamx + tmpxet][8]->Fill(nmuon, nomHOSig);
0869                       }
0870                       variab[iring2 + ringmx][8]->Fill(nmuon);
0871                       mu_projection[8]->Fill(ieta, iphi);
0872                       if (ips8 > 0) {
0873                         if (!m_cosmic) {
0874                           if (isZSps) {
0875                             ncount[iring2][18]++;
0876                             sigvsevt[iring2 + ringmx][9]->Fill(caloen[0], nomHOSig);
0877                             sig_eta_evt[netamx + tmpxet][9]->Fill(caloen[0], nomHOSig);
0878                           }
0879                           variab[iring2 + ringmx][9]->Fill(caloen[0]);
0880                           mu_projection[9]->Fill(ieta, iphi);
0881                         }
0882                         if (ips9 > 0) {
0883                           if (isZSps) {
0884                             ncount[iring2][19]++;
0885                             sigvsevt[iring2 + ringmx][10]->Fill(tkpt03, nomHOSig);
0886                             sig_eta_evt[netamx + tmpxet][10]->Fill(tkpt03, nomHOSig);
0887                           }
0888                           variab[iring2 + ringmx][10]->Fill(tkpt03);
0889                           mu_projection[10]->Fill(ieta, iphi);
0890                           if (ips10 > 0) {
0891                             ncount[iring2][3]++;  //Efficiency of Muon detection
0892                             if (isZSps) {
0893                               ncount[iring2][20]++;
0894                               sigvsevt[iring2 + ringmx][11]->Fill(hodx, nomHOSig);
0895                               sig_eta_evt[netamx + tmpxet][11]->Fill(hodx, nomHOSig);
0896                             }
0897                             variab[iring2 + ringmx][11]->Fill(hodx);
0898                             mu_projection[11]->Fill(ieta, iphi);
0899 
0900                             if (ips11 > 0) {
0901                               if (isZSps) {
0902                                 ncount[iring2][21]++;
0903                                 sigvsevt[iring2 + ringmx][12]->Fill(hody, nomHOSig);
0904                                 sig_eta_evt[netamx + tmpxet][12]->Fill(hody, nomHOSig);
0905                               }
0906                               variab[iring2 + ringmx][12]->Fill(hody);
0907                               mu_projection[12]->Fill(ieta, iphi);
0908 
0909                               if (ips12 > 0) {
0910                                 ncount[iring2][2]++;  //Efficiency of Muon detection
0911                                 if (isZSps) {
0912                                   ncount[iring2][22]++;
0913                                   sigvsevt[iring2 + ringmx][13]->Fill(htime, nomHOSig);
0914                                   sig_eta_evt[tmpxet + ringmx][13]->Fill(htime, nomHOSig);
0915                                 }
0916                                 variab[iring2 + ringmx][13]->Fill(htime);
0917                                 mu_projection[13]->Fill(ieta, iphi);
0918 
0919                                 if (ips13 > 0) {
0920                                   if (isZSps) {
0921                                     ncount[iring2][23]++;
0922                                     mu_projection[14]->Fill(ieta, iphi);
0923                                   }
0924                                 }
0925                               }
0926                             }
0927                           }
0928                         }
0929                       }
0930                     }
0931                   }
0932                 }
0933               }
0934             }
0935           }
0936         }
0937       }
0938       if (isZSps) {
0939         sigvsevt[iring2 + 2 * ringmx][0]->Fill(std::abs(ndof), nomHOSig);
0940         sigvsevt[iring2 + 2 * ringmx][1]->Fill(chisq, nomHOSig);
0941         sigvsevt[iring2 + 2 * ringmx][2]->Fill(trkth, nomHOSig);
0942         sigvsevt[iring2 + 2 * ringmx][3]->Fill(trkph, nomHOSig);
0943         sigvsevt[iring2 + 2 * ringmx][4]->Fill(1000 * therr, nomHOSig);
0944         sigvsevt[iring2 + 2 * ringmx][5]->Fill(1000 * pherr, nomHOSig);
0945         if (abshoang > 0.01) {
0946           sigvsevt[iring2 + 2 * ringmx][6]->Fill(tmphoang, (nomHOSig)*abshoang);
0947         }
0948         sigvsevt[iring2 + 2 * ringmx][7]->Fill(std::fabs(trkmm), nomHOSig);
0949         sigvsevt[iring2 + 2 * ringmx][8]->Fill(nmuon, nomHOSig);
0950         if (!m_cosmic)
0951           sigvsevt[iring2 + 2 * ringmx][9]->Fill(caloen[0], nomHOSig);
0952         sigvsevt[iring2 + 2 * ringmx][10]->Fill(tkpt03, nomHOSig);
0953         sigvsevt[iring2 + 2 * ringmx][11]->Fill(hodx, nomHOSig);
0954         sigvsevt[iring2 + 2 * ringmx][12]->Fill(hody, nomHOSig);
0955         sigvsevt[iring2 + 2 * ringmx][13]->Fill(htime, nomHOSig);
0956 
0957         sig_eta_evt[2 * netamx + tmpxet][0]->Fill(std::abs(ndof), nomHOSig);
0958         sig_eta_evt[2 * netamx + tmpxet][1]->Fill(chisq, nomHOSig);
0959         sig_eta_evt[2 * netamx + tmpxet][2]->Fill(trkth, nomHOSig);
0960         sig_eta_evt[2 * netamx + tmpxet][3]->Fill(trkph, nomHOSig);
0961         sig_eta_evt[2 * netamx + tmpxet][4]->Fill(1000 * therr, nomHOSig);
0962         sig_eta_evt[2 * netamx + tmpxet][5]->Fill(1000 * pherr, nomHOSig);
0963         if (abshoang > 0.01) {
0964           sig_eta_evt[2 * netamx + tmpxet][6]->Fill(tmphoang, (nomHOSig)*abshoang);
0965         }
0966         sig_eta_evt[2 * netamx + tmpxet][7]->Fill(std::fabs(trkmm), nomHOSig);
0967         sig_eta_evt[2 * netamx + tmpxet][8]->Fill(nmuon, nomHOSig);
0968         if (!m_cosmic)
0969           sig_eta_evt[2 * netamx + tmpxet][9]->Fill(caloen[0], nomHOSig);
0970         sig_eta_evt[2 * netamx + tmpxet][10]->Fill(tkpt03, nomHOSig);
0971         sig_eta_evt[2 * netamx + tmpxet][11]->Fill(hodx, nomHOSig);
0972         sig_eta_evt[2 * netamx + tmpxet][12]->Fill(hody, nomHOSig);
0973         sig_eta_evt[2 * netamx + tmpxet][13]->Fill(htime, nomHOSig);
0974       }
0975 
0976       variab[iring2 + 2 * ringmx][0]->Fill(std::abs(ndof));
0977       variab[iring2 + 2 * ringmx][1]->Fill(chisq);
0978       variab[iring2 + 2 * ringmx][2]->Fill(trkth);
0979       variab[iring2 + 2 * ringmx][3]->Fill(trkph);
0980       variab[iring2 + 2 * ringmx][4]->Fill(1000 * therr);
0981       variab[iring2 + 2 * ringmx][5]->Fill(1000 * pherr);
0982       variab[iring2 + 2 * ringmx][6]->Fill(tmphoang);
0983       variab[iring2 + 2 * ringmx][7]->Fill(std::fabs(trkmm));
0984       variab[iring2 + 2 * ringmx][8]->Fill(nmuon);
0985       if (!m_cosmic)
0986         variab[iring2 + 2 * ringmx][9]->Fill(caloen[0]);
0987       variab[iring2 + 2 * ringmx][10]->Fill(tkpt03);
0988       variab[iring2 + 2 * ringmx][11]->Fill(hodx);
0989       variab[iring2 + 2 * ringmx][12]->Fill(hody);
0990       variab[iring2 + 2 * ringmx][13]->Fill(htime);
0991 
0992       muonnm->Fill(nmuon);
0993       muonmm->Fill(trkmm);
0994       muonth->Fill(trkth * 180 / angle_units::piRadians);
0995       muonph->Fill(trkph * 180 / angle_units::piRadians);
0996       muonch->Fill(chisq);
0997 
0998       int iselect = (ipsall == mypow_2[ncut] - 1) ? 1 : 0;
0999 
1000       if (iselect == 1) {
1001         ipass++;
1002         sel_muonnm->Fill(nmuon);
1003         sel_muonmm->Fill(trkmm);
1004         sel_muonth->Fill(trkth * 180 / angle_units::piRadians);
1005         sel_muonph->Fill(trkph * 180 / angle_units::piRadians);
1006         sel_muonch->Fill(chisq);
1007         if (m_histFill && tmpxet >= 0 && tmpxet < netamx && iphi >= 0 && iphi < nphimx) {
1008           ho_indenergy[tmpxet][iphi - 1]->Fill(nomHOSig);
1009         }
1010         if (m_treeFill) {
1011           T1->Fill();
1012         }
1013       }
1014     }  //close the for loop: (HOCalibVariableCollection::const_iterator hoC=(*HOCalib).begin(); hoC!=(*HOCalib).end(); hoC++){
1015   }  //end of the if loop (isCosMu)
1016 }
1017 
1018 //define this as a plug-in
1019 DEFINE_FWK_MODULE(HOCalibAnalyzer);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!