Reference Manual

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