CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_2_7_hltpatch2/src/PhysicsTools/PatAlgos/plugins/PATElectronProducer.cc

Go to the documentation of this file.
00001 //
00002 // $Id: PATElectronProducer.cc,v 1.57 2012/04/17 09:27:33 tjkim Exp $
00003 //
00004 #include "PhysicsTools/PatAlgos/plugins/PATElectronProducer.h"
00005 
00006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00007 #include "FWCore/ParameterSet/interface/FileInPath.h"
00008 
00009 #include "DataFormats/Common/interface/Association.h"
00010 #include "DataFormats/Common/interface/ValueMap.h"
00011 #include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
00012 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
00013 
00014 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
00015 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
00016 
00017 #include "PhysicsTools/PatUtils/interface/TrackerIsolationPt.h"
00018 #include "PhysicsTools/PatUtils/interface/CaloIsolationEnergy.h"
00019 
00020 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00021 #include "DataFormats/VertexReco/interface/Vertex.h"
00022 
00023 
00024 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
00025 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
00026 
00027 #include "TrackingTools/TransientTrack/interface/TransientTrackBuilder.h"
00028 #include "TrackingTools/Records/interface/TransientTrackRecord.h"
00029 #include "TrackingTools/TransientTrack/interface/TransientTrack.h"
00030 #include "TrackingTools/IPTools/interface/IPTools.h"
00031 
00032 #include "DataFormats/GsfTrackReco/interface/GsfTrackFwd.h"
00033 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
00034 
00035 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterLazyTools.h"
00036 #include "RecoEgamma/EgammaTools/interface/ConversionTools.h"
00037 
00038 #include <vector>
00039 #include <memory>
00040 
00041 
00042 using namespace pat;
00043 using namespace std;
00044 
00045 
00046 PATElectronProducer::PATElectronProducer(const edm::ParameterSet & iConfig) :
00047   isolator_(iConfig.exists("userIsolation") ? iConfig.getParameter<edm::ParameterSet>("userIsolation") : edm::ParameterSet(), false) ,
00048   useUserData_(iConfig.exists("userData"))
00049 {
00050 
00051   // general configurables
00052   electronSrc_      = iConfig.getParameter<edm::InputTag>( "electronSource" );
00053   embedGsfElectronCore_    = iConfig.getParameter<bool>         ( "embedGsfElectronCore" );
00054   embedGsfTrack_    = iConfig.getParameter<bool>         ( "embedGsfTrack" );
00055   embedSuperCluster_= iConfig.getParameter<bool>         ( "embedSuperCluster" );
00056   embedTrack_       = iConfig.getParameter<bool>         ( "embedTrack" );
00057 
00058   // pflow specific
00059   pfElecSrc_           = iConfig.getParameter<edm::InputTag>( "pfElectronSource" );
00060   pfCandidateMap_   = iConfig.getParameter<edm::InputTag>( "pfCandidateMap" );
00061   useParticleFlow_        = iConfig.getParameter<bool>( "useParticleFlow" );
00062   embedPFCandidate_   = iConfig.getParameter<bool>( "embedPFCandidate" );
00063 
00064   // MC matching configurables
00065   addGenMatch_      = iConfig.getParameter<bool>          ( "addGenMatch" );
00066   if (addGenMatch_) {
00067     embedGenMatch_ = iConfig.getParameter<bool>         ( "embedGenMatch" );
00068     if (iConfig.existsAs<edm::InputTag>("genParticleMatch")) {
00069       genMatchSrc_.push_back(iConfig.getParameter<edm::InputTag>( "genParticleMatch" ));
00070     } else {
00071       genMatchSrc_ = iConfig.getParameter<std::vector<edm::InputTag> >( "genParticleMatch" );
00072     }
00073   }
00074 
00075   // resolution configurables
00076   addResolutions_   = iConfig.getParameter<bool>         ( "addResolutions" );
00077   if (addResolutions_) {
00078     resolutionLoader_ = pat::helper::KinResolutionsLoader(iConfig.getParameter<edm::ParameterSet>("resolutions"));
00079   }
00080 
00081 
00082   // electron ID configurables
00083   addElecID_        = iConfig.getParameter<bool>         ( "addElectronID" );
00084   if (addElecID_) {
00085     // it might be a single electron ID
00086     if (iConfig.existsAs<edm::InputTag>("electronIDSource")) {
00087       elecIDSrcs_.push_back(NameTag("", iConfig.getParameter<edm::InputTag>("electronIDSource")));
00088     }
00089     // or there might be many of them
00090     if (iConfig.existsAs<edm::ParameterSet>("electronIDSources")) {
00091       // please don't configure me twice
00092       if (!elecIDSrcs_.empty()) throw cms::Exception("Configuration") <<
00093                                   "PATElectronProducer: you can't specify both 'electronIDSource' and 'electronIDSources'\n";
00094       // read the different electron ID names
00095       edm::ParameterSet idps = iConfig.getParameter<edm::ParameterSet>("electronIDSources");
00096       std::vector<std::string> names = idps.getParameterNamesForType<edm::InputTag>();
00097       for (std::vector<std::string>::const_iterator it = names.begin(), ed = names.end(); it != ed; ++it) {
00098         elecIDSrcs_.push_back(NameTag(*it, idps.getParameter<edm::InputTag>(*it)));
00099       }
00100     }
00101     // but in any case at least once
00102     if (elecIDSrcs_.empty()) throw cms::Exception("Configuration") <<
00103                                "PATElectronProducer: id addElectronID is true, you must specify either:\n" <<
00104                                "\tInputTag electronIDSource = <someTag>\n" << "or\n" <<
00105                                "\tPSet electronIDSources = { \n" <<
00106                                "\t\tInputTag <someName> = <someTag>   // as many as you want \n " <<
00107                                "\t}\n";
00108   }
00109 
00110   // construct resolution calculator
00111 
00112   //   // IsoDeposit configurables
00113   //   if (iConfig.exists("isoDeposits")) {
00114   //      edm::ParameterSet depconf = iConfig.getParameter<edm::ParameterSet>("isoDeposits");
00115   //      if (depconf.exists("tracker")) isoDepositLabels_.push_back(std::make_pair(TrackerIso, depconf.getParameter<edm::InputTag>("tracker")));
00116   //      if (depconf.exists("ecal"))    isoDepositLabels_.push_back(std::make_pair(ECalIso, depconf.getParameter<edm::InputTag>("ecal")));
00117   //      if (depconf.exists("hcal"))    isoDepositLabels_.push_back(std::make_pair(HCalIso, depconf.getParameter<edm::InputTag>("hcal")));
00118 
00119 
00120   //      if (depconf.exists("user")) {
00121   //         std::vector<edm::InputTag> userdeps = depconf.getParameter<std::vector<edm::InputTag> >("user");
00122   //         std::vector<edm::InputTag>::const_iterator it = userdeps.begin(), ed = userdeps.end();
00123   //         int key = UserBaseIso;
00124   //         for ( ; it != ed; ++it, ++key) {
00125   //             isoDepositLabels_.push_back(std::make_pair(IsolationKeys(key), *it));
00126   //         }
00127   //      }
00128   //   }
00129 
00130   // read isoDeposit labels, for direct embedding
00131   readIsolationLabels(iConfig, "isoDeposits", isoDepositLabels_);
00132 
00133   // read isolation value labels, for direct embedding
00134   readIsolationLabels(iConfig, "isolationValues", isolationValueLabels_);
00135 
00136   // read isolation value labels for non PF identified electron, for direct embedding
00137   readIsolationLabels(iConfig, "isolationValuesNoPFId", isolationValueLabelsNoPFId_);
00138 
00139   // Efficiency configurables
00140   addEfficiencies_ = iConfig.getParameter<bool>("addEfficiencies");
00141   if (addEfficiencies_) {
00142     efficiencyLoader_ = pat::helper::EfficiencyLoader(iConfig.getParameter<edm::ParameterSet>("efficiencies"));
00143   }
00144 
00145   // Check to see if the user wants to add user data
00146   if ( useUserData_ ) {
00147     userDataHelper_ = PATUserDataHelper<Electron>(iConfig.getParameter<edm::ParameterSet>("userData"));
00148   }
00149 
00150   // embed high level selection variables?
00151   embedHighLevelSelection_ = iConfig.getParameter<bool>("embedHighLevelSelection");
00152   if ( embedHighLevelSelection_ ) {
00153     beamLineSrc_ = iConfig.getParameter<edm::InputTag>("beamLineSrc");
00154     usePV_ = iConfig.getParameter<bool>("usePV");
00155     pvSrc_ = iConfig.getParameter<edm::InputTag>("pvSrc");
00156   }
00157 
00158 
00159   // produces vector of muons
00160   produces<std::vector<Electron> >();
00161 
00162 }
00163 
00164 
00165 PATElectronProducer::~PATElectronProducer() {
00166 }
00167 
00168 
00169 void PATElectronProducer::produce(edm::Event & iEvent, const edm::EventSetup & iSetup) {
00170 
00171   // Get the collection of electrons from the event
00172   edm::Handle<edm::View<reco::GsfElectron> > electrons;
00173   iEvent.getByLabel(electronSrc_, electrons);
00174 
00175   if (iEvent.isRealData()){
00176        addGenMatch_ = false;
00177        embedGenMatch_ = false;
00178    }
00179 
00180   // for additional mva variables
00181   edm::InputTag  reducedEBRecHitCollection(string("reducedEcalRecHitsEB"));
00182   edm::InputTag  reducedEERecHitCollection(string("reducedEcalRecHitsEE"));
00183   EcalClusterLazyTools lazyTools(iEvent, iSetup, reducedEBRecHitCollection, reducedEERecHitCollection);
00184 
00185   // for conversion veto selection  
00186   edm::Handle<reco::ConversionCollection> hConversions;
00187   iEvent.getByLabel("allConversions", hConversions);
00188 
00189   // Get the ESHandle for the transient track builder, if needed for
00190   // high level selection embedding
00191   edm::ESHandle<TransientTrackBuilder> trackBuilder;
00192 
00193   if (isolator_.enabled()) isolator_.beginEvent(iEvent,iSetup);
00194 
00195   if (efficiencyLoader_.enabled()) efficiencyLoader_.newEvent(iEvent);
00196   if (resolutionLoader_.enabled()) resolutionLoader_.newEvent(iEvent, iSetup);
00197 
00198   IsoDepositMaps deposits(isoDepositLabels_.size());
00199   for (size_t j = 0, nd = deposits.size(); j < nd; ++j) {
00200     iEvent.getByLabel(isoDepositLabels_[j].second, deposits[j]);
00201   }
00202 
00203   IsolationValueMaps isolationValues(isolationValueLabels_.size());
00204   for (size_t j = 0; j<isolationValueLabels_.size(); ++j) {
00205     iEvent.getByLabel(isolationValueLabels_[j].second, isolationValues[j]);
00206   }
00207 
00208   IsolationValueMaps isolationValuesNoPFId(isolationValueLabelsNoPFId_.size());
00209   for (size_t j = 0; j<isolationValueLabelsNoPFId_.size(); ++j) {
00210     iEvent.getByLabel(isolationValueLabelsNoPFId_[j].second, isolationValuesNoPFId[j]);
00211   }
00212 
00213   // prepare the MC matching
00214   GenAssociations  genMatches(genMatchSrc_.size());
00215   if (addGenMatch_) {
00216     for (size_t j = 0, nd = genMatchSrc_.size(); j < nd; ++j) {
00217       iEvent.getByLabel(genMatchSrc_[j], genMatches[j]);
00218     }
00219   }
00220 
00221   // prepare ID extraction
00222   std::vector<edm::Handle<edm::ValueMap<float> > > idhandles;
00223   std::vector<pat::Electron::IdPair>               ids;
00224   if (addElecID_) {
00225     idhandles.resize(elecIDSrcs_.size());
00226     ids.resize(elecIDSrcs_.size());
00227     for (size_t i = 0; i < elecIDSrcs_.size(); ++i) {
00228       iEvent.getByLabel(elecIDSrcs_[i].second, idhandles[i]);
00229       ids[i].first = elecIDSrcs_[i].first;
00230     }
00231   }
00232 
00233 
00234   // prepare the high level selection:
00235   // needs beamline
00236   reco::TrackBase::Point beamPoint(0,0,0);
00237   reco::Vertex primaryVertex;
00238   reco::BeamSpot beamSpot;
00239   bool beamSpotIsValid = false;
00240   bool primaryVertexIsValid = false;
00241 
00242   // Get the beamspot
00243   edm::Handle<reco::BeamSpot> beamSpotHandle;
00244   iEvent.getByLabel(beamLineSrc_, beamSpotHandle);
00245 
00246   // Get the primary vertex
00247   edm::Handle< std::vector<reco::Vertex> > pvHandle;
00248   iEvent.getByLabel( pvSrc_, pvHandle );
00249 
00250 
00251   if ( embedHighLevelSelection_ ) {
00252 
00253     // This is needed by the IPTools methods from the tracking group
00254     iSetup.get<TransientTrackRecord>().get("TransientTrackBuilder", trackBuilder);
00255 
00256     if ( ! usePV_ ) {
00257 
00258       if ( beamSpotHandle.isValid() ){
00259         beamSpot = *beamSpotHandle;
00260         beamSpotIsValid = true;
00261       } else{
00262         edm::LogError("DataNotAvailable")
00263           << "No beam spot available from EventSetup, not adding high level selection \n";
00264       }
00265 
00266       double x0 = beamSpot.x0();
00267       double y0 = beamSpot.y0();
00268       double z0 = beamSpot.z0();
00269 
00270       beamPoint = reco::TrackBase::Point ( x0, y0, z0 );
00271     } else {
00272       if ( pvHandle.isValid() && !pvHandle->empty() ) {
00273         primaryVertex = pvHandle->at(0);
00274         primaryVertexIsValid = true;
00275       } else {
00276         edm::LogError("DataNotAvailable")
00277           << "No primary vertex available from EventSetup, not adding high level selection \n";
00278       }
00279     }
00280   }
00281 
00282   std::vector<Electron> * patElectrons = new std::vector<Electron>();
00283 
00284   if( useParticleFlow_ ) {
00285     edm::Handle< reco::PFCandidateCollection >  pfElectrons;
00286     iEvent.getByLabel(pfElecSrc_, pfElectrons);
00287     unsigned index=0;
00288 
00289     for( reco::PFCandidateConstIterator i = pfElectrons->begin();
00290          i != pfElectrons->end(); ++i, ++index) {
00291 
00292       reco::PFCandidateRef pfRef(pfElectrons, index);
00293       reco::PFCandidatePtr ptrToPFElectron(pfElectrons,index);
00294 //       reco::CandidateBaseRef pfBaseRef( pfRef );
00295 
00296       reco::GsfTrackRef PfTk= i->gsfTrackRef();
00297 
00298       bool Matched=false;
00299       bool MatchedToAmbiguousGsfTrack=false;
00300       for (edm::View<reco::GsfElectron>::const_iterator itElectron = electrons->begin(); itElectron != electrons->end(); ++itElectron) {
00301         unsigned int idx = itElectron - electrons->begin();
00302         if (Matched || MatchedToAmbiguousGsfTrack) continue;
00303 
00304         reco::GsfTrackRef EgTk= itElectron->gsfTrack();
00305 
00306         if (itElectron->gsfTrack()==i->gsfTrackRef()){
00307           Matched=true;
00308         }
00309         else {
00310           for( reco::GsfTrackRefVector::const_iterator it = itElectron->ambiguousGsfTracksBegin() ;
00311                it!=itElectron->ambiguousGsfTracksEnd(); it++ ){
00312             MatchedToAmbiguousGsfTrack |= (bool)(i->gsfTrackRef()==(*it));
00313           }
00314         }
00315 
00316         if (Matched || MatchedToAmbiguousGsfTrack){
00317 
00318           // ptr needed for finding the matched gen particle
00319           reco::CandidatePtr ptrToGsfElectron(electrons,idx);
00320 
00321           // ref to base needed for the construction of the pat object
00322           const edm::RefToBase<reco::GsfElectron>& elecsRef = electrons->refAt(idx);
00323           Electron anElectron(elecsRef);
00324           anElectron.setPFCandidateRef( pfRef  );
00325 
00326           if( embedPFCandidate_ ) anElectron.embedPFCandidate();
00327 
00328           if ( useUserData_ ) {
00329             userDataHelper_.add( anElectron, iEvent, iSetup );
00330           }
00331 
00332           double ip3d = -999; // for mva variable
00333  
00334           // embed high level selection
00335           if ( embedHighLevelSelection_ ) {
00336             // get the global track
00337             reco::GsfTrackRef track = PfTk;
00338 
00339             // Make sure the collection it points to is there
00340             if ( track.isNonnull() && track.isAvailable() ) {
00341 
00342               reco::TransientTrack tt = trackBuilder->build(track);
00343               embedHighLevel( anElectron, 
00344                               track,
00345                               tt,
00346                               primaryVertex,
00347                               primaryVertexIsValid,
00348                               beamSpot,
00349                               beamSpotIsValid );
00350 
00351               std::pair<bool,Measurement1D> ip3dpv = IPTools::absoluteImpactParameter3D(tt, primaryVertex);
00352               ip3d = ip3dpv.second.value(); // for mva variable
00353 
00354               if ( !usePV_ ) {
00355                 double corr_d0 = track->dxy( beamPoint );
00356                 anElectron.setDB( corr_d0, -1.0 );
00357               } else {
00358                  std::pair<bool,Measurement1D> result = IPTools::absoluteTransverseImpactParameter(tt, primaryVertex);
00359                 double d0_corr = result.second.value();
00360                 double d0_err = result.second.error();
00361                 anElectron.setDB( d0_corr, d0_err );
00362               }
00363             }
00364           }
00365 
00366           //Electron Id
00367 
00368           if (addElecID_) {
00369             //STANDARD EL ID
00370             for (size_t i = 0; i < elecIDSrcs_.size(); ++i) {
00371               ids[i].second = (*idhandles[i])[elecsRef];
00372             }
00373             //SPECIFIC PF ID
00374             ids.push_back(std::make_pair("pf_evspi",pfRef->mva_e_pi()));
00375             ids.push_back(std::make_pair("pf_evsmu",pfRef->mva_e_mu()));
00376             anElectron.setElectronIDs(ids);
00377           }
00378 
00379           // add missing mva variables
00380           double r9 = lazyTools.e3x3( *( itElectron->superCluster()->seed())) / itElectron->superCluster()->rawEnergy() ;
00381           double sigmaIphiIphi;
00382           double sigmaIetaIphi;
00383           std::vector<float> vCov = lazyTools.localCovariances(*( itElectron->superCluster()->seed()));
00384           if( !isnan(vCov[2])) sigmaIphiIphi = sqrt(vCov[2]);
00385           else sigmaIphiIphi = 0;
00386           sigmaIetaIphi = vCov[1];
00387           anElectron.setMvaVariables( r9, sigmaIphiIphi, sigmaIetaIphi, ip3d);
00388 
00389           // set conversion veto selection
00390           bool passconversionveto = false;
00391           if( hConversions.isValid()){          
00392             // this is recommended method
00393             passconversionveto = !ConversionTools::hasMatchedConversion( *itElectron, hConversions, beamSpotHandle->position());
00394           }else{
00395             // use missing hits without vertex fit method
00396             passconversionveto = itElectron->gsfTrack()->trackerExpectedHitsInner().numberOfLostHits() < 1;
00397           }
00398 
00399           anElectron.setPassConversionVeto( passconversionveto );
00400 
00401 
00402 //        fillElectron(anElectron,elecsRef,pfBaseRef,
00403 //                     genMatches, deposits, isolationValues);
00404 
00405           //COLIN small warning !
00406           // we are currently choosing to take the 4-momentum of the PFCandidate;
00407           // the momentum of the GsfElectron is saved though
00408           // we must therefore match the GsfElectron.
00409           // because of this, we should not change the source of the electron matcher
00410           // to the collection of PFElectrons in the python configuration
00411           // I don't know what to do with the efficiencyLoader, since I don't know
00412           // what this class is for.
00413           fillElectron2( anElectron,
00414                          ptrToPFElectron,
00415                          ptrToGsfElectron,
00416                          ptrToGsfElectron,
00417                          genMatches, deposits, isolationValues );
00418 
00419           //COLIN need to use fillElectron2 in the non-pflow case as well, and to test it.
00420 
00421           patElectrons->push_back(anElectron);
00422         }
00423       }
00424       //if( !Matched && !MatchedToAmbiguousGsfTrack) std::cout << "!!!!A pf electron could not be matched to a gsf!!!!"  << std::endl;
00425     }
00426   }
00427 
00428   else{
00429     for (edm::View<reco::GsfElectron>::const_iterator itElectron = electrons->begin(); itElectron != electrons->end(); ++itElectron) {
00430       // construct the Electron from the ref -> save ref to original object
00431       //FIXME: looks like a lot of instances could be turned into const refs
00432       unsigned int idx = itElectron - electrons->begin();
00433       edm::RefToBase<reco::GsfElectron> elecsRef = electrons->refAt(idx);
00434       reco::CandidateBaseRef elecBaseRef(elecsRef);
00435       Electron anElectron(elecsRef);
00436 
00437       // Is this GsfElectron also identified as an e- in the particle flow?
00438       bool pfId = false;
00439 
00440       // first try to access PF electron collection
00441       edm::Handle<edm::ValueMap<reco::PFCandidatePtr> >ValMapH;
00442       bool valMapPresent = iEvent.getByLabel(pfCandidateMap_,ValMapH);
00443       if( valMapPresent ) {
00444         const edm::ValueMap<reco::PFCandidatePtr> & myValMap(*ValMapH); 
00445         
00446         // Get the PFCandidate
00447         const reco::PFCandidatePtr& pfElePtr(myValMap[elecsRef]);
00448         pfId= pfElePtr.isNonnull();
00449       }
00450       else {
00451         // PF electron collection not available. 
00452         // trying to access a PFCandidate collection, as supplied by the user
00453         edm::Handle< reco::PFCandidateCollection >  pfElectrons;
00454         iEvent.getByLabel(pfElecSrc_, pfElectrons);
00455         const reco::GsfTrackRef& trkRef = itElectron->gsfTrack();
00456         for( reco::PFCandidateConstIterator ie = pfElectrons->begin();
00457              ie != pfElectrons->end(); ++ie) {
00458           if(ie->particleId()!=reco::PFCandidate::e) continue;
00459           const reco::GsfTrackRef& pfTrkRef= ie->gsfTrackRef();
00460           if( trkRef == pfTrkRef ) {
00461             pfId = true;
00462             break;
00463           }
00464         }
00465       }
00466             
00467       // add resolution info
00468 
00469       // Isolation
00470       if (isolator_.enabled()) {
00471         isolator_.fill(*electrons, idx, isolatorTmpStorage_);
00472         typedef pat::helper::MultiIsolator::IsolationValuePairs IsolationValuePairs;
00473         // better to loop backwards, so the vector is resized less times
00474         for (IsolationValuePairs::const_reverse_iterator it = isolatorTmpStorage_.rbegin(), ed = isolatorTmpStorage_.rend(); it != ed; ++it) {
00475           anElectron.setIsolation(it->first, it->second);
00476         }
00477       }
00478 
00479       for (size_t j = 0, nd = deposits.size(); j < nd; ++j) {
00480         anElectron.setIsoDeposit(isoDepositLabels_[j].first, (*deposits[j])[elecsRef]);
00481       }
00482 
00483       // add electron ID info
00484       if (addElecID_) {
00485         for (size_t i = 0; i < elecIDSrcs_.size(); ++i) {
00486           ids[i].second = (*idhandles[i])[elecsRef];
00487         }
00488         anElectron.setElectronIDs(ids);
00489       }
00490 
00491 
00492       if ( useUserData_ ) {
00493         userDataHelper_.add( anElectron, iEvent, iSetup );
00494       }
00495 
00496 
00497       double ip3d = -999; //for mva variable
00498 
00499       // embed high level selection
00500       if ( embedHighLevelSelection_ ) {
00501         // get the global track
00502         reco::GsfTrackRef track = itElectron->gsfTrack();
00503 
00504         // Make sure the collection it points to is there
00505         if ( track.isNonnull() && track.isAvailable() ) {
00506 
00507           reco::TransientTrack tt = trackBuilder->build(track);
00508           embedHighLevel( anElectron, 
00509                           track,
00510                           tt,
00511                           primaryVertex,
00512                           primaryVertexIsValid,
00513                           beamSpot,
00514                           beamSpotIsValid );
00515 
00516           std::pair<bool,Measurement1D> ip3dpv = IPTools::absoluteImpactParameter3D(tt, primaryVertex);
00517           ip3d = ip3dpv.second.value(); // for mva variable
00518 
00519           if ( !usePV_ ) {
00520             double corr_d0 = track->dxy( beamPoint );
00521             anElectron.setDB( corr_d0, -1.0 );
00522           } else {
00523             std::pair<bool,Measurement1D> result = IPTools::absoluteTransverseImpactParameter(tt, primaryVertex);
00524             double d0_corr = result.second.value();
00525             double d0_err = result.second.error();
00526             anElectron.setDB( d0_corr, d0_err );
00527           }
00528         }
00529       }
00530 
00531       // add mva variables
00532       double r9 = lazyTools.e3x3( *( itElectron->superCluster()->seed())) / itElectron->superCluster()->rawEnergy() ;
00533       double sigmaIphiIphi;
00534       double sigmaIetaIphi;
00535       std::vector<float> vCov = lazyTools.localCovariances(*( itElectron->superCluster()->seed()));
00536       if( !isnan(vCov[2])) sigmaIphiIphi = sqrt(vCov[2]);
00537       else sigmaIphiIphi = 0;
00538       sigmaIetaIphi = vCov[1];
00539       anElectron.setMvaVariables( r9, sigmaIphiIphi, sigmaIetaIphi, ip3d);
00540 
00541       // set conversion veto selection
00542       bool passconversionveto = false;
00543       if( hConversions.isValid()){
00544         // this is recommended method 
00545         passconversionveto = !ConversionTools::hasMatchedConversion( *itElectron, hConversions, beamSpotHandle->position());
00546       }else{
00547         // use missing hits without vertex fit method
00548         passconversionveto = itElectron->gsfTrack()->trackerExpectedHitsInner().numberOfLostHits() < 1;
00549       }
00550       anElectron.setPassConversionVeto( passconversionveto );
00551 
00552       // add sel to selected
00553       fillElectron( anElectron, elecsRef,elecBaseRef,
00554                     genMatches, deposits, pfId, isolationValues, isolationValuesNoPFId);
00555       patElectrons->push_back(anElectron);
00556     }
00557   }
00558 
00559   // sort electrons in pt
00560   std::sort(patElectrons->begin(), patElectrons->end(), pTComparator_);
00561 
00562   // add the electrons to the event output
00563   std::auto_ptr<std::vector<Electron> > ptr(patElectrons);
00564   iEvent.put(ptr);
00565 
00566   // clean up
00567   if (isolator_.enabled()) isolator_.endEvent();
00568 
00569 }
00570 
00571 void PATElectronProducer::fillElectron(Electron& anElectron,
00572                                        const edm::RefToBase<reco::GsfElectron>& elecRef,
00573                                        const reco::CandidateBaseRef& baseRef,
00574                                        const GenAssociations& genMatches,
00575                                        const IsoDepositMaps& deposits,
00576                                        const bool pfId,
00577                                        const IsolationValueMaps& isolationValues,
00578                                        const IsolationValueMaps& isolationValuesNoPFId
00579                                        ) const {
00580 
00581   //COLIN: might want to use the PFCandidate 4-mom. Which one is in use now?
00582   //   if (useParticleFlow_)
00583   //     aMuon.setP4( aMuon.pfCandidateRef()->p4() );
00584 
00585   //COLIN:
00586   //In the embedding case, the reference cannot be used to look into a value map.
00587   //therefore, one has to had the PFCandidateRef to this function, which becomes a bit
00588   //too much specific.
00589 
00590   // in fact, this function needs a baseref or ptr for genmatch
00591   // and a baseref or ptr for isodeposits and isolationvalues.
00592   // baseref is not needed
00593   // the ptrForIsolation and ptrForMatching should be defined upstream.
00594 
00595   // is the concrete elecRef needed for the efficiency loader? what is this loader?
00596   // how can we make it compatible with the particle flow electrons?
00597 
00598   if (embedGsfElectronCore_) anElectron.embedGsfElectronCore();
00599   if (embedGsfTrack_) anElectron.embedGsfTrack();
00600   if (embedSuperCluster_) anElectron.embedSuperCluster();
00601   if (embedTrack_) anElectron.embedTrack();
00602 
00603   // store the match to the generated final state muons
00604   if (addGenMatch_) {
00605     for(size_t i = 0, n = genMatches.size(); i < n; ++i) {
00606       if(useParticleFlow_) {
00607         reco::GenParticleRef genElectron = (*genMatches[i])[anElectron.pfCandidateRef()];
00608         anElectron.addGenParticleRef(genElectron);
00609       }
00610       else {
00611         reco::GenParticleRef genElectron = (*genMatches[i])[elecRef];
00612         anElectron.addGenParticleRef(genElectron);
00613       }
00614     }
00615     if (embedGenMatch_) anElectron.embedGenParticle();
00616   }
00617 
00618   if (efficiencyLoader_.enabled()) {
00619     efficiencyLoader_.setEfficiencies( anElectron, elecRef );
00620   }
00621 
00622   if (resolutionLoader_.enabled()) {
00623     resolutionLoader_.setResolutions(anElectron);
00624   }
00625 
00626   for (size_t j = 0, nd = deposits.size(); j < nd; ++j) {
00627     if(useParticleFlow_) {
00628 
00629       reco::PFCandidateRef pfcandref =  anElectron.pfCandidateRef();
00630       assert(!pfcandref.isNull());
00631       reco::CandidatePtr source = pfcandref->sourceCandidatePtr(0);
00632       anElectron.setIsoDeposit(isoDepositLabels_[j].first,
00633                           (*deposits[j])[source]);
00634     }
00635     else
00636       anElectron.setIsoDeposit(isoDepositLabels_[j].first,
00637                           (*deposits[j])[elecRef]);
00638   }
00639 
00640   for (size_t j = 0; j<isolationValues.size(); ++j) {
00641     if(useParticleFlow_) {
00642       reco::CandidatePtr source = anElectron.pfCandidateRef()->sourceCandidatePtr(0);
00643       anElectron.setIsolation(isolationValueLabels_[j].first,
00644                          (*isolationValues[j])[source]);
00645     }
00646     else
00647       if(pfId){
00648         anElectron.setIsolation(isolationValueLabels_[j].first,(*isolationValues[j])[elecRef]);
00649       }else{
00650         anElectron.setIsolation(isolationValueLabelsNoPFId_[j].first,(*isolationValuesNoPFId[j])[elecRef]);
00651       }
00652   }
00653 
00654 }
00655 
00656 void PATElectronProducer::fillElectron2( Electron& anElectron,
00657                                          const reco::CandidatePtr& candPtrForIsolation,
00658                                          const reco::CandidatePtr& candPtrForGenMatch,
00659                                          const reco::CandidatePtr& candPtrForLoader,
00660                                          const GenAssociations& genMatches,
00661                                          const IsoDepositMaps& deposits,
00662                                          const IsolationValueMaps& isolationValues) const {
00663 
00664   //COLIN/Florian: use the PFCandidate 4-mom.
00665   anElectron.setEcalDrivenMomentum(anElectron.p4()) ;
00666   anElectron.setP4( anElectron.pfCandidateRef()->p4() );
00667 
00668 
00669   // is the concrete elecRef needed for the efficiency loader? what is this loader?
00670   // how can we make it compatible with the particle flow electrons?
00671 
00672   if (embedGsfElectronCore_) anElectron.embedGsfElectronCore();
00673   if (embedGsfTrack_) anElectron.embedGsfTrack();
00674   if (embedSuperCluster_) anElectron.embedSuperCluster();
00675   if (embedTrack_) anElectron.embedTrack();
00676 
00677   // store the match to the generated final state muons
00678 
00679   if (addGenMatch_) {
00680     for(size_t i = 0, n = genMatches.size(); i < n; ++i) {
00681       reco::GenParticleRef genElectron = (*genMatches[i])[candPtrForGenMatch];
00682       anElectron.addGenParticleRef(genElectron);
00683     }
00684     if (embedGenMatch_) anElectron.embedGenParticle();
00685   }
00686 
00687   //COLIN what's this? does it have to be GsfElectron specific?
00688   if (efficiencyLoader_.enabled()) {
00689     efficiencyLoader_.setEfficiencies( anElectron, candPtrForLoader );
00690   }
00691 
00692   if (resolutionLoader_.enabled()) {
00693     resolutionLoader_.setResolutions(anElectron);
00694   }
00695 
00696   for (size_t j = 0, nd = deposits.size(); j < nd; ++j) {
00697     if( isoDepositLabels_[j].first==pat::TrackIso ||
00698         isoDepositLabels_[j].first==pat::EcalIso ||
00699         isoDepositLabels_[j].first==pat::HcalIso || 
00700         deposits[j]->contains(candPtrForGenMatch.id())) {
00701       anElectron.setIsoDeposit(isoDepositLabels_[j].first,
00702                                (*deposits[j])[candPtrForGenMatch]);
00703     }
00704     else if (deposits[j]->contains(candPtrForIsolation.id())) {
00705       anElectron.setIsoDeposit(isoDepositLabels_[j].first,
00706                                (*deposits[j])[candPtrForIsolation]);
00707     }
00708     else {
00709       anElectron.setIsoDeposit(isoDepositLabels_[j].first,
00710                                (*deposits[j])[candPtrForIsolation->sourceCandidatePtr(0)]);
00711     }
00712   }
00713 
00714   for (size_t j = 0; j<isolationValues.size(); ++j) {
00715     if( isolationValueLabels_[j].first==pat::TrackIso ||
00716         isolationValueLabels_[j].first==pat::EcalIso ||
00717         isolationValueLabels_[j].first==pat::HcalIso || 
00718         isolationValues[j]->contains(candPtrForGenMatch.id())) {
00719       anElectron.setIsolation(isolationValueLabels_[j].first,
00720                               (*isolationValues[j])[candPtrForGenMatch]);
00721     }
00722     else if (isolationValues[j]->contains(candPtrForIsolation.id())) {
00723       anElectron.setIsolation(isolationValueLabels_[j].first,
00724                               (*isolationValues[j])[candPtrForIsolation]);
00725     }
00726     else {
00727       anElectron.setIsolation(isolationValueLabels_[j].first,
00728                               (*isolationValues[j])[candPtrForIsolation->sourceCandidatePtr(0)]);
00729     }    
00730   }
00731 }
00732 
00733 
00734 // ParameterSet description for module
00735 void PATElectronProducer::fillDescriptions(edm::ConfigurationDescriptions & descriptions)
00736 {
00737   edm::ParameterSetDescription iDesc;
00738   iDesc.setComment("PAT electron producer module");
00739 
00740   // input source
00741   iDesc.add<edm::InputTag>("pfCandidateMap", edm::InputTag("no default"))->setComment("input collection");
00742   iDesc.add<edm::InputTag>("electronSource", edm::InputTag("no default"))->setComment("input collection");
00743 
00744   // embedding
00745   iDesc.add<bool>("embedGsfElectronCore", true)->setComment("embed external gsf electron core");
00746   iDesc.add<bool>("embedGsfTrack", true)->setComment("embed external gsf track");
00747   iDesc.add<bool>("embedSuperCluster", true)->setComment("embed external super cluster");
00748   iDesc.add<bool>("embedTrack", false)->setComment("embed external track");
00749 
00750   // pf specific parameters
00751   iDesc.add<edm::InputTag>("pfElectronSource", edm::InputTag("pfElectrons"))->setComment("particle flow input collection");
00752   iDesc.add<bool>("useParticleFlow", false)->setComment("whether to use particle flow or not");
00753   iDesc.add<bool>("embedPFCandidate", false)->setComment("embed external particle flow object");
00754 
00755   // MC matching configurables
00756   iDesc.add<bool>("addGenMatch", true)->setComment("add MC matching");
00757   iDesc.add<bool>("embedGenMatch", false)->setComment("embed MC matched MC information");
00758   std::vector<edm::InputTag> emptySourceVector;
00759   iDesc.addNode( edm::ParameterDescription<edm::InputTag>("genParticleMatch", edm::InputTag(), true) xor
00760                  edm::ParameterDescription<std::vector<edm::InputTag> >("genParticleMatch", emptySourceVector, true)
00761                  )->setComment("input with MC match information");
00762 
00763   // electron ID configurables
00764   iDesc.add<bool>("addElectronID",true)->setComment("add electron ID variables");
00765   edm::ParameterSetDescription electronIDSourcesPSet;
00766   electronIDSourcesPSet.setAllowAnything();
00767   iDesc.addNode( edm::ParameterDescription<edm::InputTag>("electronIDSource", edm::InputTag(), true) xor
00768                  edm::ParameterDescription<edm::ParameterSetDescription>("electronIDSources", electronIDSourcesPSet, true)
00769                  )->setComment("input with electron ID variables");
00770 
00771 
00772   // IsoDeposit configurables
00773   edm::ParameterSetDescription isoDepositsPSet;
00774   isoDepositsPSet.addOptional<edm::InputTag>("tracker");
00775   isoDepositsPSet.addOptional<edm::InputTag>("ecal");
00776   isoDepositsPSet.addOptional<edm::InputTag>("hcal");
00777   isoDepositsPSet.addOptional<edm::InputTag>("pfAllParticles");
00778   isoDepositsPSet.addOptional<edm::InputTag>("pfChargedHadrons");
00779   isoDepositsPSet.addOptional<edm::InputTag>("pfChargedAll");
00780   isoDepositsPSet.addOptional<edm::InputTag>("pfPUChargedHadrons");
00781   isoDepositsPSet.addOptional<edm::InputTag>("pfNeutralHadrons");
00782   isoDepositsPSet.addOptional<edm::InputTag>("pfPhotons");
00783   isoDepositsPSet.addOptional<std::vector<edm::InputTag> >("user");
00784   iDesc.addOptional("isoDeposits", isoDepositsPSet);
00785 
00786   // isolation values configurables
00787   edm::ParameterSetDescription isolationValuesPSet;
00788   isolationValuesPSet.addOptional<edm::InputTag>("tracker");
00789   isolationValuesPSet.addOptional<edm::InputTag>("ecal");
00790   isolationValuesPSet.addOptional<edm::InputTag>("hcal");
00791   isolationValuesPSet.addOptional<edm::InputTag>("pfAllParticles");
00792   isolationValuesPSet.addOptional<edm::InputTag>("pfChargedHadrons");
00793   isolationValuesPSet.addOptional<edm::InputTag>("pfChargedAll");
00794   isolationValuesPSet.addOptional<edm::InputTag>("pfPUChargedHadrons");
00795   isolationValuesPSet.addOptional<edm::InputTag>("pfNeutralHadrons");
00796   isolationValuesPSet.addOptional<edm::InputTag>("pfPhotons");
00797   isolationValuesPSet.addOptional<std::vector<edm::InputTag> >("user");
00798   iDesc.addOptional("isolationValues", isolationValuesPSet);
00799 
00800   // isolation values configurables
00801   edm::ParameterSetDescription isolationValuesNoPFIdPSet;
00802   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("tracker");
00803   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("ecal");
00804   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("hcal");
00805   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("pfAllParticles");
00806   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("pfChargedHadrons");
00807   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("pfChargedAll");
00808   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("pfPUChargedHadrons");
00809   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("pfNeutralHadrons");
00810   isolationValuesNoPFIdPSet.addOptional<edm::InputTag>("pfPhotons");
00811   isolationValuesNoPFIdPSet.addOptional<std::vector<edm::InputTag> >("user");
00812   iDesc.addOptional("isolationValuesNoPFId", isolationValuesNoPFIdPSet);
00813 
00814   // Efficiency configurables
00815   edm::ParameterSetDescription efficienciesPSet;
00816   efficienciesPSet.setAllowAnything(); // TODO: the pat helper needs to implement a description.
00817   iDesc.add("efficiencies", efficienciesPSet);
00818   iDesc.add<bool>("addEfficiencies", false);
00819 
00820   // Check to see if the user wants to add user data
00821   edm::ParameterSetDescription userDataPSet;
00822   PATUserDataHelper<Electron>::fillDescription(userDataPSet);
00823   iDesc.addOptional("userData", userDataPSet);
00824 
00825   // electron shapes
00826   iDesc.add<bool>("addElectronShapes", true);
00827   iDesc.add<edm::InputTag>("reducedBarrelRecHitCollection", edm::InputTag("reducedEcalRecHitsEB"));
00828   iDesc.add<edm::InputTag>("reducedEndcapRecHitCollection", edm::InputTag("reducedEcalRecHitsEE"));
00829 
00830   edm::ParameterSetDescription isolationPSet;
00831   isolationPSet.setAllowAnything(); // TODO: the pat helper needs to implement a description.
00832   iDesc.add("userIsolation", isolationPSet);
00833 
00834   // Resolution configurables
00835   pat::helper::KinResolutionsLoader::fillDescription(iDesc);
00836 
00837   iDesc.add<bool>("embedHighLevelSelection", true)->setComment("embed high level selection");
00838   edm::ParameterSetDescription highLevelPSet;
00839   highLevelPSet.setAllowAnything();
00840   iDesc.addNode( edm::ParameterDescription<edm::InputTag>("beamLineSrc", edm::InputTag(), true)
00841                  )->setComment("input with high level selection");
00842   iDesc.addNode( edm::ParameterDescription<edm::InputTag>("pvSrc", edm::InputTag(), true)
00843                  )->setComment("input with high level selection");
00844   iDesc.addNode( edm::ParameterDescription<bool>("usePV", bool(), true)
00845                  )->setComment("input with high level selection, use primary vertex (true) or beam line (false)");
00846 
00847   descriptions.add("PATElectronProducer", iDesc);
00848 
00849 }
00850 
00851 
00852 
00853 void PATElectronProducer::readIsolationLabels( const edm::ParameterSet & iConfig,
00854                                                const char* psetName,
00855                                                IsolationLabels& labels) {
00856 
00857   labels.clear();
00858 
00859   if (iConfig.exists( psetName )) {
00860     edm::ParameterSet depconf
00861       = iConfig.getParameter<edm::ParameterSet>(psetName);
00862 
00863     if (depconf.exists("tracker")) labels.push_back(std::make_pair(pat::TrackIso, depconf.getParameter<edm::InputTag>("tracker")));
00864     if (depconf.exists("ecal"))    labels.push_back(std::make_pair(pat::EcalIso, depconf.getParameter<edm::InputTag>("ecal")));
00865     if (depconf.exists("hcal"))    labels.push_back(std::make_pair(pat::HcalIso, depconf.getParameter<edm::InputTag>("hcal")));
00866     if (depconf.exists("pfAllParticles"))  {
00867       labels.push_back(std::make_pair(pat::PfAllParticleIso, depconf.getParameter<edm::InputTag>("pfAllParticles")));
00868     }
00869     if (depconf.exists("pfChargedHadrons"))  {
00870       labels.push_back(std::make_pair(pat::PfChargedHadronIso, depconf.getParameter<edm::InputTag>("pfChargedHadrons")));
00871     }
00872     if (depconf.exists("pfChargedAll"))  {
00873       labels.push_back(std::make_pair(pat::PfChargedAllIso, depconf.getParameter<edm::InputTag>("pfChargedAll")));
00874     }
00875     if (depconf.exists("pfPUChargedHadrons"))  {
00876       labels.push_back(std::make_pair(pat::PfPUChargedHadronIso, depconf.getParameter<edm::InputTag>("pfPUChargedHadrons")));
00877     }
00878     if (depconf.exists("pfNeutralHadrons"))  {
00879       labels.push_back(std::make_pair(pat::PfNeutralHadronIso, depconf.getParameter<edm::InputTag>("pfNeutralHadrons")));
00880     }
00881     if (depconf.exists("pfPhotons")) {
00882       labels.push_back(std::make_pair(pat::PfGammaIso, depconf.getParameter<edm::InputTag>("pfPhotons")));
00883     }
00884     if (depconf.exists("user")) {
00885       std::vector<edm::InputTag> userdeps = depconf.getParameter<std::vector<edm::InputTag> >("user");
00886       std::vector<edm::InputTag>::const_iterator it = userdeps.begin(), ed = userdeps.end();
00887       int key = UserBaseIso;
00888       for ( ; it != ed; ++it, ++key) {
00889         labels.push_back(std::make_pair(IsolationKeys(key), *it));
00890       }
00891     }
00892   }
00893 
00894 
00895 }
00896 
00897 
00898 // embed various impact parameters with errors
00899 // embed high level selection
00900 void PATElectronProducer::embedHighLevel( pat::Electron & anElectron, 
00901                                           reco::GsfTrackRef track,
00902                                           reco::TransientTrack & tt,
00903                                           reco::Vertex & primaryVertex,
00904                                           bool primaryVertexIsValid,
00905                                           reco::BeamSpot & beamspot,
00906                                           bool beamspotIsValid
00907                                           )
00908 {
00909   // Correct to PV
00910 
00911   // PV2D
00912   std::pair<bool,Measurement1D> result =
00913     IPTools::signedTransverseImpactParameter(tt,
00914                                              GlobalVector(track->px(),
00915                                                           track->py(),
00916                                                           track->pz()),
00917                                              primaryVertex); 
00918   double d0_corr = result.second.value();
00919   double d0_err = primaryVertexIsValid ? result.second.error() : -1.0;
00920   anElectron.setDB( d0_corr, d0_err, pat::Electron::PV2D);
00921 
00922 
00923   // PV3D
00924   result =
00925     IPTools::signedImpactParameter3D(tt,
00926                                      GlobalVector(track->px(),
00927                                                   track->py(),
00928                                                   track->pz()),
00929                                      primaryVertex);
00930   d0_corr = result.second.value();
00931   d0_err = primaryVertexIsValid ? result.second.error() : -1.0;
00932   anElectron.setDB( d0_corr, d0_err, pat::Electron::PV3D);
00933   
00934 
00935   // Correct to beam spot
00936   // make a fake vertex out of beam spot
00937   reco::Vertex vBeamspot(beamspot.position(), beamspot.covariance3D());
00938   
00939   // BS2D
00940   result =
00941     IPTools::signedTransverseImpactParameter(tt,
00942                                              GlobalVector(track->px(),
00943                                                           track->py(),
00944                                                           track->pz()),
00945                                              vBeamspot);
00946   d0_corr = result.second.value();
00947   d0_err = beamspotIsValid ? result.second.error() : -1.0;
00948   anElectron.setDB( d0_corr, d0_err, pat::Electron::BS2D);
00949   
00950   // BS3D
00951   result =
00952     IPTools::signedImpactParameter3D(tt,
00953                                      GlobalVector(track->px(),
00954                                                   track->py(),
00955                                                   track->pz()),
00956                                      vBeamspot);
00957   d0_corr = result.second.value();
00958   d0_err = beamspotIsValid ? result.second.error() : -1.0;
00959   anElectron.setDB( d0_corr, d0_err, pat::Electron::BS3D);
00960 }
00961 
00962 #include "FWCore/Framework/interface/MakerMacros.h"
00963 
00964 DEFINE_FWK_MODULE(PATElectronProducer);