CMSSW Reference Manual

00001 // system include files
00002 #include <memory>
00003 #include <string>
00004 
00005 // user include files
00006 #include "FWCore/Framework/interface/Frameworkfwd.h"
00007 #include "FWCore/Framework/interface/EDAnalyzer.h"
00008 
00009 #include "FWCore/Framework/interface/Event.h"
00010 #include "FWCore/Framework/interface/MakerMacros.h"
00011 
00012 #include "FWCore/ParameterSet/interface/ParameterSet.h"
00013 //needed for TFileService
00014 #include "PhysicsTools/UtilAlgos/interface/TFileService.h"
00015 #include "FWCore/ServiceRegistry/interface/Service.h"
00016 //needed for MessageLogger
00017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00018 
00019 #include "AnalysisDataFormats/TopObjects/interface/TtGenEvent.h"
00020 
00021 #include "DataFormats/PatCandidates/interface/Muon.h"
00022 #include "DataFormats/PatCandidates/interface/Jet.h"
00023 #include "DataFormats/PatCandidates/interface/MET.h"
00024 #include "DataFormats/PatCandidates/interface/Tau.h"
00025 #include "DataFormats/PatCandidates/interface/Electron.h"
00026 
00027 #include "TDirectory.h"
00028 #include "TH1F.h"
00029 #include "TF1.h"
00030 #include "TTree.h"
00031 
00032 #include <Math/VectorUtil.h>
00033 
00034 //
00035 // class declaration
00036 //
00037 
00038 class ResolutionCreator : public edm::EDAnalyzer {
00039    public:
00040       explicit ResolutionCreator(const edm::ParameterSet&);
00041       ~ResolutionCreator();
00042 
00043    private:
00044       virtual void beginJob(const edm::EventSetup&) ;
00045       virtual void analyze(const edm::Event&, const edm::EventSetup&);
00046       virtual void endJob() ;
00047 
00048       // ----------member data ---------------------------
00049                   std::string objectType_, labelName_;
00050                   std::vector<double> etabinVals_, pTbinVals_;
00051                 double minDR_;
00052       int matchingAlgo_;
00053       bool useMaxDist_;
00054       bool useDeltaR_;
00055       double maxDist_;
00056                 int ptnrbins, etanrbins;
00057                 int nrFilled;
00058 
00059                         //Histograms are booked in the beginJob() method
00060                   TF1 *fResPtEtaBin[10][20][20];
00061                 TF1 *fResEtaBin[10][20];
00062                 TH1F *hResPtEtaBin[10][20][20];
00063                 TH1F *hResEtaBin[10][20];
00064       TTree* tResVar;
00065 };
00066 
00067 
00068 //
00069 // constructors and destructor
00070 //
00071 ResolutionCreator::ResolutionCreator(const edm::ParameterSet& iConfig)
00072 {
00073   // input parameters
00074   objectType_   = iConfig.getParameter< std::string >           ("object");
00075   labelName_  = iConfig.getParameter< std::string >             ("label");
00076         if(objectType_ != "met"){
00077     etabinVals_ = iConfig.getParameter< std::vector<double> >   ("etabinValues");
00078   }
00079   pTbinVals_    = iConfig.getParameter< std::vector<double> >   ("pTbinValues");
00080   minDR_        = iConfig.getParameter< double > ("minMatchingDR");
00081 
00082         nrFilled = 0;
00083 
00084 }
00085 
00086 
00087 ResolutionCreator::~ResolutionCreator()
00088 {
00089 }
00090 
00091 
00092 //
00093 // member functions
00094 //
00095 
00096 // ------------ method called to for each event  ------------
00097 void
00098 ResolutionCreator::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup)
00099 {
00100    using namespace edm;
00101    using namespace std;
00102    using namespace pat;
00103    
00104          // Get the gen and cal object fourvector
00105    vector<reco::Particle *> p4gen, p4rec;
00106             
00107    Handle<TtGenEvent> genEvt;
00108    iEvent.getByLabel ("genEvt",genEvt);
00109 
00110    if(genEvt->particles().size()<10) return;
00111 
00112    if(objectType_ == "electron"){ 
00113      Handle<std::vector<pat::Electron> >  electrons; //to calculate the ResolutionCreator for the electrons, i used the TopElectron instead of the AOD information
00114      iEvent.getByLabel(labelName_,electrons);
00115      for(size_t e=0; e<electrons->size(); e++) { 
00116        for(size_t p=0; p<genEvt->particles().size(); p++){
00117          if( (abs(genEvt->particles()[p].pdgId()) == 11) && (ROOT::Math::VectorUtil::DeltaR(genEvt->particles()[p].p4(), (*electrons)[e].p4()) < minDR_) ) {
00118            p4gen.push_back(new reco::Particle(genEvt->particles()[p]));
00119                                  p4rec.push_back(new reco::Particle((pat::Electron)((*electrons)[e])));
00120                                  }
00121        }
00122      }
00123    }
00124    else if(objectType_ == "muon"){
00125      Handle<std::vector<pat::Muon> >  muons;
00126      iEvent.getByLabel(labelName_,muons);
00127      for(size_t m=0; m<muons->size(); m++) {      
00128        for(size_t p=0; p<genEvt->particles().size(); p++){
00129          if( (abs(genEvt->particles()[p].pdgId()) == 13) && (ROOT::Math::VectorUtil::DeltaR(genEvt->particles()[p].p4(), (*muons)[m].p4()) < minDR_) ) {
00130            p4gen.push_back(new reco::Particle(genEvt->particles()[p]));
00131            p4rec.push_back(new reco::Particle((pat::Muon)((*muons)[m])));
00132                                  }
00133        }
00134      }
00135    }
00136    else if(objectType_ == "lJets" ){
00137      Handle<std::vector<pat::Jet> > jets;
00138      iEvent.getByLabel(labelName_,jets);
00139                  if(jets->size()>=4) { 
00140        for(unsigned int j = 0; j<4; j++){      
00141          for(size_t p=0; p<genEvt->particles().size(); p++){
00142            if( (abs(genEvt->particles()[p].pdgId()) < 5) && (ROOT::Math::VectorUtil::DeltaR(genEvt->particles()[p].p4(), (*jets)[j].p4())< minDR_) ){
00143                   // std::cout <<"genEvt->particles()[p].status() "<<genEvt->particles()[p].status() << std::endl;
00144                                                 // const reco::Candidate *mom = genEvt->particles()[p].mother();
00145                                                 // std::cout <<"id mother particle "<<mom->pdgId()<< " status mother particle "<< mom->status()<< std::endl;
00146                                  p4gen.push_back(new reco::Particle(genEvt->particles()[p]));
00147                                  p4rec.push_back(new reco::Particle((pat::Jet)(*jets)[j]));
00148                                  }
00149                                  }
00150                          }
00151                  }
00152          }
00153          else if(objectType_ == "bJets" ){
00154      Handle<std::vector<pat::Jet> > jets;
00155      iEvent.getByLabel(labelName_,jets);
00156                  if(jets->size()>=4) { 
00157        for(unsigned int j = 0; j<4; j++){      
00158          for(size_t p=0; p<genEvt->particles().size(); p++){
00159                                          if( (abs(genEvt->particles()[p].pdgId()) == 5) && (ROOT::Math::VectorUtil::DeltaR(genEvt->particles()[p].p4(), (*jets)[j].p4())< minDR_) ) {
00160                   // std::cout <<"genEvt->particles()[p].status() "<<genEvt->particles()[p].status() << std::endl;
00161                                                 // const reco::Candidate *mom = genEvt->particles()[p].mother();
00162                                                 // std::cout <<"id mother particle "<<mom->pdgId()<< " status mother particle "<< mom->status()<< std::endl;
00163                                                  p4gen.push_back(new reco::Particle(genEvt->particles()[p]));
00164                                  p4rec.push_back(new reco::Particle((pat::Jet)(*jets)[j]));
00165                                  }
00166                                  }
00167        }
00168      }
00169    }
00170    else if(objectType_ == "met"){
00171      Handle<std::vector<pat::MET> >  mets;
00172      iEvent.getByLabel(labelName_,mets);
00173      if(mets->size()>=1) { 
00174        if( genEvt->isSemiLeptonic() && genEvt->singleNeutrino() != 0 && ROOT::Math::VectorUtil::DeltaR(genEvt->singleNeutrino()->p4(), (*mets)[0].p4()) < minDR_) {
00175          p4gen.push_back(new reco::Particle(0,genEvt->singleNeutrino()->p4(),math::XYZPoint()));
00176          p4rec.push_back(new reco::Particle((pat::MET)((*mets)[0])));
00177        }
00178      }
00179    } 
00180    else if(objectType_ == "tau"){
00181      Handle<std::vector<pat::Tau> > taus; 
00182      iEvent.getByLabel(labelName_,taus);
00183      for(std::vector<pat::Tau>::const_iterator tau = taus->begin(); tau != taus->end(); ++tau) {
00184        // find the tau (if any) that matches a MC tau from W
00185        reco::GenParticle genLepton = *(tau->genLepton());
00186        if( abs(genLepton.pdgId())==15 && genLepton.status()==2 &&
00187            genLepton.numberOfMothers()>0 &&
00188            abs(genLepton.mother(0)->pdgId())==15 &&
00189            genLepton.mother(0)->numberOfMothers()>0 &&
00190            abs(genLepton.mother(0)->mother(0)->pdgId())==24 &&
00191                                  ROOT::Math::VectorUtil::DeltaR(genLepton.p4(), tau->p4()) < minDR_  ) {
00192        }
00193        p4gen.push_back(new reco::Particle(genLepton));
00194        p4rec.push_back(new reco::Particle(*tau));
00195      }
00196    }
00197    // Fill the object's value
00198      for(unsigned m=0; m<p4gen.size(); m++){ 
00199        double Egen     = p4gen[m]->energy(); 
00200        double Thetagen = p4gen[m]->theta(); 
00201        double Phigen   = p4gen[m]->phi();
00202        double Etgen    = p4gen[m]->et();
00203        double Etagen   = p4gen[m]->eta();
00204        double Ecal     = p4rec[m]->energy(); 
00205        double Thetacal = p4rec[m]->theta();
00206        double Phical   = p4rec[m]->phi();
00207        double Etcal    = p4rec[m]->et();
00208        double Etacal   = p4rec[m]->eta();
00209        double phidiff  = Phical- Phigen;
00210        if(phidiff>3.14159)  phidiff = 2.*3.14159 - phidiff;
00211        if(phidiff<-3.14159) phidiff = -phidiff - 2.*3.14159;
00212    
00213        // find eta and et bin
00214        int etabin  =  0;
00215        if(etanrbins > 1){
00216          for(unsigned int b=0; b<etabinVals_.size()-1; b++) {
00217            if(fabs(Etacal) > etabinVals_[b]) etabin = b;
00218          }
00219        }
00220      
00221        int ptbin  =  0;
00222        for(unsigned int b=0; b<pTbinVals_.size()-1; b++) {
00223          if(p4rec[m]->pt() > pTbinVals_[b]) ptbin = b;
00224        }
00225      
00226        // calculate the resolution on "a", "b", "c" & "d" according to the definition (CMS-NOTE-2006-023):
00227        // p = a*|p_meas|*u_1 + b*u_2 + c*u_3
00228        // E(fit) = E_meas * d
00229        //
00230        // with u_1 = p/|p_meas|
00231        //      u_3 = (u_z x u_1)/|u_z x u_1|
00232        //      u_2 = (u_1 x u_3)/|u_1 x u_3|
00233        //
00234        // The initial parameters values are chosen like (a, b, c, d) = (1., 0., 0., 1.)
00235 
00236        // 1/ calculate the unitary vectors of the basis u_1, u_2, u_3
00237        ROOT::Math::SVector<double,3> pcalvec(p4rec[m]->px(),p4rec[m]->py(),p4rec[m]->pz());
00238        ROOT::Math::SVector<double,3> pgenvec(p4gen[m]->px(),p4gen[m]->py(),p4gen[m]->pz());
00239        
00240        ROOT::Math::SVector<double,3> u_z(0,0,1);
00241        ROOT::Math::SVector<double,3> u_1 = ROOT::Math::Unit(pcalvec);
00242        ROOT::Math::SVector<double,3> u_3 = ROOT::Math::Cross(u_z,u_1)/ROOT::Math::Mag(ROOT::Math::Cross(u_z,u_1));
00243        ROOT::Math::SVector<double,3> u_2 = ROOT::Math::Cross(u_1,u_3)/ROOT::Math::Mag(ROOT::Math::Cross(u_1,u_3));
00244        double acal = 1.;
00245        double bcal = 0.;
00246        double ccal = 0.;
00247        double dcal = 1.;
00248        double agen = ROOT::Math::Dot(pgenvec,u_1)/ROOT::Math::Mag(pcalvec);
00249        double bgen = ROOT::Math::Dot(pgenvec,u_2);
00250        double cgen = ROOT::Math::Dot(pgenvec,u_3);
00251        double dgen = Egen/Ecal;
00252                                 
00253        //fill histograms    
00254        ++nrFilled; 
00255        hResPtEtaBin[0][etabin][ptbin] -> Fill(acal-agen);
00256        hResPtEtaBin[1][etabin][ptbin] -> Fill(bcal-bgen);
00257        hResPtEtaBin[2][etabin][ptbin] -> Fill(ccal-cgen);
00258        hResPtEtaBin[3][etabin][ptbin] -> Fill(dcal-dgen);
00259        hResPtEtaBin[4][etabin][ptbin] -> Fill(Thetacal-Thetagen);
00260        hResPtEtaBin[5][etabin][ptbin] -> Fill(phidiff);
00261        hResPtEtaBin[6][etabin][ptbin] -> Fill(Etcal-Etgen);
00262        hResPtEtaBin[7][etabin][ptbin] -> Fill(Etacal-Etagen);
00263 
00264        delete p4gen[m];
00265        delete p4rec[m];
00266      }
00267                  
00268 }
00269 
00270 
00271 // ------------ method called once each job just before starting event loop  ------------
00272 void 
00273 ResolutionCreator::beginJob(const edm::EventSetup&)
00274 {
00275   edm::Service<TFileService> fs;
00276   if (!fs) throw edm::Exception(edm::errors::Configuration, "TFileService missing from configuration!");
00277 
00278         // input constants  
00279   TString         resObsName[8]         = {"_ares","_bres","_cres","_dres","_thres","_phres","_etres","_etares"};
00280   int             resObsNrBins          = 120;
00281   if( (objectType_ == "muon") || (objectType_ == "electron") ) resObsNrBins = 80;
00282   std::vector<double>  resObsMin, resObsMax;
00283   if(objectType_ == "electron"){ 
00284     resObsMin.push_back(-0.15);  resObsMin.push_back(-0.2);  resObsMin.push_back(-0.1);  resObsMin.push_back(-0.15);  resObsMin.push_back(-0.0012); resObsMin.push_back(-0.009);  resObsMin.push_back(-16);   resObsMin.push_back(-0.0012);   
00285     resObsMax.push_back( 0.15);  resObsMax.push_back( 0.2);  resObsMax.push_back( 0.1);  resObsMax.push_back( 0.15);  resObsMax.push_back( 0.0012); resObsMax.push_back( 0.009);  resObsMax.push_back( 16);   resObsMax.push_back( 0.0012);
00286   } else if(objectType_ == "muon"){
00287     resObsMin.push_back(-0.15);  resObsMin.push_back(-0.1);  resObsMin.push_back(-0.05);  resObsMin.push_back(-0.15);  resObsMin.push_back(-0.004);  resObsMin.push_back(-0.003);  resObsMin.push_back(-8);    resObsMin.push_back(-0.004);   
00288     resObsMax.push_back( 0.15);  resObsMax.push_back( 0.1);  resObsMax.push_back( 0.05);  resObsMax.push_back( 0.15);  resObsMax.push_back( 0.004);  resObsMax.push_back( 0.003);  resObsMax.push_back( 8);    resObsMax.push_back( 0.004);
00289   } else if(objectType_ == "tau"){ 
00290     resObsMin.push_back(-1.);    resObsMin.push_back(-10.);  resObsMin.push_back(-10);   resObsMin.push_back(-1.);   resObsMin.push_back(-0.1);    resObsMin.push_back(-0.1);    resObsMin.push_back(-80);   resObsMin.push_back(-0.1);   
00291     resObsMax.push_back( 1.);    resObsMax.push_back( 10.);  resObsMax.push_back( 10);   resObsMax.push_back( 1.);   resObsMax.push_back( 0.1);    resObsMax.push_back( 0.1);    resObsMax.push_back( 50);   resObsMax.push_back( 0.1);
00292   } else if(objectType_ == "lJets" || objectType_ == "bJets"){
00293     resObsMin.push_back(-1.);    resObsMin.push_back(-10.);  resObsMin.push_back(-10.);  resObsMin.push_back(-1.);   resObsMin.push_back(-0.4);    resObsMin.push_back(-0.6);    resObsMin.push_back( -80);  resObsMin.push_back(-0.6);   
00294     resObsMax.push_back( 1.);    resObsMax.push_back( 10.);  resObsMax.push_back( 10.);  resObsMax.push_back( 1.);   resObsMax.push_back( 0.4);    resObsMax.push_back( 0.6);    resObsMax.push_back( 80);   resObsMax.push_back( 0.6);
00295   } else{
00296     resObsMin.push_back(-2.);   resObsMin.push_back(-150.); resObsMin.push_back(-150.); resObsMin.push_back(-2.);   resObsMin.push_back(-6);      resObsMin.push_back(-6);      resObsMin.push_back( -180); resObsMin.push_back(-6);   
00297     resObsMax.push_back( 3.);   resObsMax.push_back( 150.); resObsMax.push_back( 150.); resObsMax.push_back( 3.);   resObsMax.push_back( 6);      resObsMax.push_back( 6);      resObsMax.push_back(  180); resObsMax.push_back( 6);
00298   }
00299   
00300   const char*   resObsVsPtFit[8]        = {     "[0]+[1]*exp(-[2]*x)",
00301                                           "[0]+[1]*exp(-[2]*x)",
00302                                           "[0]+[1]*exp(-[2]*x)",
00303                                           "[0]+[1]*exp(-[2]*x)",
00304                                                                                                                                                                 "[0]+[1]*exp(-[2]*x)",
00305                                                                                                                                                                 "[0]+[1]*exp(-[2]*x)",
00306                                                                                                                                                                 "pol1",
00307                                                                                                                                                                 "[0]+[1]*exp(-[2]*x)"
00308                                                                                                                                                         };
00309  
00310   ptnrbins        = pTbinVals_.size()-1;
00311   double *ptbins  = new double[pTbinVals_.size()];
00312   for(unsigned int b=0; b<pTbinVals_.size(); b++)  ptbins[b]  = pTbinVals_[b];
00313   double *etabins;
00314   if(objectType_ != "met"){
00315     etanrbins  = etabinVals_.size()-1;
00316     etabins    = new double[etabinVals_.size()];
00317     for(unsigned int b=0; b<etabinVals_.size(); b++) etabins[b] = etabinVals_[b];
00318   }else{
00319     etanrbins = 1;
00320     etabins    = new double[2];
00321     etabins[0] = 0; etabins[1] = 5.;
00322   }
00323         
00324 
00325   //define the histograms booked
00326   for(Int_t ro=0; ro<8; ro++) {
00327     for(Int_t etab=0; etab<etanrbins; etab++) { 
00328       for(Int_t ptb=0; ptb<ptnrbins; ptb++) {
00329         TString obsName = objectType_; obsName += resObsName[ro]; obsName += "_etabin"; obsName += etab; obsName += "_ptbin";
00330                                 obsName += ptb;
00331                                 hResPtEtaBin[ro][etab][ptb] = fs->make<TH1F>(obsName,obsName,resObsNrBins,resObsMin[ro],resObsMax[ro]);
00332         fResPtEtaBin[ro][etab][ptb] = fs->make<TF1>("F_"+obsName,"gaus");
00333       }
00334       TString obsName2 = objectType_; obsName2 += resObsName[ro]; obsName2 += "_etabin"; obsName2 += etab;
00335       hResEtaBin[ro][etab] = fs->make<TH1F>(obsName2,obsName2,ptnrbins,ptbins);
00336       fResEtaBin[ro][etab] = fs->make<TF1>("F_"+obsName2,resObsVsPtFit[ro],pTbinVals_[0],pTbinVals_[pTbinVals_.size()-1]);
00337     }
00338   }
00339         tResVar = fs->make< TTree >("tResVar","Resolution tree");
00340 
00341   delete [] etabins; 
00342   delete [] ptbins; 
00343 
00344 }
00345 
00346 // ------------ method called once each job just after ending the event loop  ------------
00347 void 
00348 ResolutionCreator::endJob() {
00349   TString         resObsName2[8]        = {"a","b","c","d","theta","phi","et","eta"};
00350   Int_t ro=0;
00351   Double_t pt=0.;
00352   Double_t eta=0.;
00353   Double_t value,error;
00354 
00355   tResVar->Branch("Pt",&pt,"Pt/D");
00356   tResVar->Branch("Eta",&eta,"Eta/D");
00357   tResVar->Branch("ro",&ro,"ro/I");
00358   tResVar->Branch("value",&value,"value/D");
00359   tResVar->Branch("error",&error,"error/D");
00360   
00361   for(ro=0; ro<8; ro++) {
00362     for(int etab=0; etab<etanrbins; etab++) {   
00363       //CD set eta at the center of the bin
00364       eta = etanrbins > 1 ? (etabinVals_[etab]+etabinVals_[etab+1])/2. : 2.5 ; 
00365       for(int ptb=0; ptb<ptnrbins; ptb++) {
00366                                 //CD set pt at the center of the bin
00367                                 pt = (pTbinVals_[ptb]+pTbinVals_[ptb+1])/2.; 
00368         double maxcontent = 0.;
00369                                 int maxbin = 0;
00370                                 for(int nb=1; nb<hResPtEtaBin[ro][etab][ptb]->GetNbinsX(); nb ++){
00371                                 if (hResPtEtaBin[ro][etab][ptb]->GetBinContent(nb)>maxcontent) {
00372                                 maxcontent = hResPtEtaBin[ro][etab][ptb]->GetBinContent(nb);
00373                                 maxbin = nb;
00374                                 }
00375                                 }
00376                                 int range = (int)(hResPtEtaBin[ro][etab][ptb]->GetNbinsX()/6); //in order that ~1/3 of X-axis range is fitted
00377                         fResPtEtaBin[ro][etab][ptb] -> SetRange(hResPtEtaBin[ro][etab][ptb]->GetBinCenter(maxbin-range),
00378                                 hResPtEtaBin[ro][etab][ptb]->GetBinCenter(maxbin+range));
00379                                 fResPtEtaBin[ro][etab][ptb] -> SetParameters(hResPtEtaBin[ro][etab][ptb] -> GetMaximum(),
00380                                                                                 hResPtEtaBin[ro][etab][ptb] -> GetMean(),
00381                                                                                                                                                                                                         hResPtEtaBin[ro][etab][ptb] -> GetRMS());
00382                                 hResPtEtaBin[ro][etab][ptb] -> Fit(fResPtEtaBin[ro][etab][ptb]->GetName(),"RQ");
00383                         hResEtaBin[ro][etab]        -> SetBinContent(ptb+1,fResPtEtaBin[ro][etab][ptb]->GetParameter(2));
00384                         hResEtaBin[ro][etab]        -> SetBinError(ptb+1,fResPtEtaBin[ro][etab][ptb]->GetParError(2));
00385                                 //CD: Fill the tree
00386                                 value = fResPtEtaBin[ro][etab][ptb]->GetParameter(2); //parameter value
00387                                 error = fResPtEtaBin[ro][etab][ptb]->GetParError(2);  //parameter error
00388                                 tResVar->Fill();
00389       }
00390       //CD: add a fake entry in pt=0 for the NN training
00391       // for that, use a linear extrapolation.
00392       pt = 0.;
00393       value =   ((pTbinVals_[0]+pTbinVals_[1])/2.)*(fResPtEtaBin[ro][etab][0]->GetParameter(2)-fResPtEtaBin[ro][etab][1]->GetParameter(2))/((pTbinVals_[2]-pTbinVals_[0])/2.)+fResPtEtaBin[ro][etab][0]->GetParameter(2);
00394       error = fResPtEtaBin[ro][etab][0]->GetParError(2)+fResPtEtaBin[ro][etab][1]->GetParError(2);
00395       tResVar->Fill();
00396       // standard fit
00397       hResEtaBin[ro][etab] -> Fit(fResEtaBin[ro][etab]->GetName(),"RQ");
00398     }
00399   } 
00400   if(objectType_ == "lJets" && nrFilled == 0) edm::LogProblem  ("SummaryError") << "No plots filled for light jets \n";    
00401   if(objectType_ == "bJets" && nrFilled == 0) edm::LogProblem  ("SummaryError") << "No plots filled for bjets \n";    
00402   if(objectType_ == "muon" && nrFilled == 0) edm::LogProblem  ("SummaryError") << "No plots filled for muons \n";    
00403   if(objectType_ == "electron" && nrFilled == 0) edm::LogProblem  ("SummaryError") << "No plots filled for electrons \n";    
00404   if(objectType_ == "tau" && nrFilled == 0) edm::LogProblem  ("SummaryError") << "No plots filled for taus \n";    
00405   if(objectType_ == "met" && nrFilled == 0) edm::LogProblem  ("SummaryError") << "No plots filled for met \n";    
00406         
00407   edm::LogVerbatim ("MainResults") << " \n\n";  
00408   edm::LogVerbatim ("MainResults") << " ----------------------------------------------";
00409   edm::LogVerbatim ("MainResults") << " ----------------------------------------------";
00410   edm::LogVerbatim ("MainResults") << " Resolutions on "<< objectType_ << " with nrfilled: "<<nrFilled;
00411   edm::LogVerbatim ("MainResults") << " ----------------------------------------------";
00412   edm::LogVerbatim ("MainResults") << " ----------------------------------------------";
00413   if(nrFilled != 0 && objectType_ != "met") {
00414         for(ro=0; ro<8; ro++) {
00415                         edm::LogVerbatim ("MainResults") << "-------------------- ";
00416         edm::LogVerbatim ("MainResults") << "\n Resolutions on " << resObsName2[ro] << "\n";
00417                 edm::LogVerbatim ("MainResults") << "-------------------- ";
00418                         for(int etab=0; etab<etanrbins; etab++) {       
00419                         if(nrFilled != 0 && ro != 6) {
00420                                         if(etab == 0){
00421                                                 edm::LogVerbatim   ("MainResults") << "if(fabs(eta)<"<<etabinVals_[etab+1] <<") res = " << 
00422                                                 fResEtaBin[ro][etab]->GetParameter(0) << "+" << fResEtaBin[ro][etab]->GetParameter(1) 
00423                                                 << "*exp(-(" <<fResEtaBin[ro][etab]->GetParameter(2) << "*pt));";  
00424                                         }else{ 
00425                                                 edm::LogVerbatim   ("MainResults") << "else if(fabs(eta)<"<<etabinVals_[etab+1] <<") res = " << 
00426                                                 fResEtaBin[ro][etab]->GetParameter(0) << "+" << fResEtaBin[ro][etab]->GetParameter(1) 
00427                                                 << "*exp(-(" <<fResEtaBin[ro][etab]->GetParameter(2) << "*pt));";                                       
00428                                         }
00429                                 }else if(nrFilled != 0 && ro == 6){
00430                                         if(etab == 0){
00431                                                 edm::LogVerbatim   ("MainResults") << "if(fabs(eta)<"<<etabinVals_[etab+1] <<") res = " << 
00432                                                 fResEtaBin[ro][etab]->GetParameter(0) << "+" << fResEtaBin[ro][etab]->GetParameter(1) 
00433                                                 << "*pt;";
00434                                         }else{                                          
00435                                                 edm::LogVerbatim   ("MainResults") << "else if(fabs(eta)<"<<etabinVals_[etab+1] <<") res = " << 
00436                                                 fResEtaBin[ro][etab]->GetParameter(0) << "+" << fResEtaBin[ro][etab]->GetParameter(1) 
00437                                                 << "*pt;";
00438 
00439                                         }
00440                                 }
00441                         }
00442                 }
00443         }else if(nrFilled != 0 && objectType_ == "met"){
00444         for(ro=0; ro<8; ro++) {
00445                         edm::LogVerbatim ("MainResults") << "-------------------- ";
00446         edm::LogVerbatim ("MainResults") << "\n Resolutions on " << resObsName2[ro] << "\n";
00447                 edm::LogVerbatim ("MainResults") << "-------------------- ";
00448                         if(nrFilled != 0 && ro != 6) {
00449                                         edm::LogVerbatim   ("MainResults") << "res = " <<
00450                                         fResEtaBin[ro][0]->GetParameter(0) << "+" << fResEtaBin[ro][0]->GetParameter(1) 
00451                                         << "*exp(-(" <<fResEtaBin[ro][0]->GetParameter(2) << "*pt));";                          
00452                         }else if(nrFilled != 0 && ro == 6){
00453                                         edm::LogVerbatim   ("MainResults") << "res = " << 
00454                                         fResEtaBin[ro][0]->GetParameter(0) << "+" << fResEtaBin[ro][0]->GetParameter(1) << "*pt;";
00455                         }
00456                 }
00457         }
00458 }
00459 
00460 //define this as a plug-in
00461 DEFINE_FWK_MODULE(ResolutionCreator);
ResolutionCreator.cc
