CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_3/src/RecoEgamma/EgammaPhotonProducers/src/PhotonProducer.cc

Go to the documentation of this file.
00001 #include <iostream>
00002 #include <vector>
00003 #include <memory>
00004 
00005 // Framework
00006 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00007 #include "FWCore/Utilities/interface/Exception.h"
00008 
00009 #include "CommonTools/Utils/interface/StringToEnumValue.h"
00010 
00011 #include "Geometry/Records/interface/CaloGeometryRecord.h"
00012 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
00013 #include "Geometry/CaloTopology/interface/CaloTopology.h"
00014 #include "Geometry/CaloEventSetup/interface/CaloTopologyRecord.h"
00015 
00016 
00017 #include "DataFormats/VertexReco/interface/Vertex.h"
00018 #include "DataFormats/VertexReco/interface/VertexFwd.h"
00019 #include "DataFormats/EgammaReco/interface/ClusterShape.h"
00020 #include "DataFormats/EgammaCandidates/interface/PhotonCore.h"
00021 #include "DataFormats/EgammaCandidates/interface/Photon.h"
00022 #include "DataFormats/EgammaCandidates/interface/PhotonFwd.h"
00023 #include "DataFormats/EgammaCandidates/interface/Conversion.h"
00024 
00025 #include "DataFormats/EgammaReco/interface/ElectronSeed.h"
00026 #include "RecoCaloTools/Selectors/interface/CaloConeSelector.h"
00027 
00028 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgoRcd.h"
00029 #include "RecoLocalCalo/EcalRecAlgos/interface/EcalSeverityLevelAlgo.h"
00030 
00031 #include "RecoEgamma/EgammaPhotonProducers/interface/PhotonProducer.h"
00032 #include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaTowerIsolation.h"
00033 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterFunctionBaseClass.h" 
00034 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterFunctionFactory.h" 
00035 #include "RecoEcal/EgammaCoreTools/plugins/EcalClusterCrackCorrection.h"
00036 #include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaHadTower.h"
00037 
00038 PhotonProducer::PhotonProducer(const edm::ParameterSet& config) : 
00039 
00040   conf_(config)
00041 {
00042 
00043   // use onfiguration file to setup input/output collection names
00044 
00045   photonCoreProducer_   = conf_.getParameter<edm::InputTag>("photonCoreProducer");
00046   barrelEcalHits_   = conf_.getParameter<edm::InputTag>("barrelEcalHits");
00047   endcapEcalHits_   = conf_.getParameter<edm::InputTag>("endcapEcalHits");
00048   vertexProducer_   = conf_.getParameter<std::string>("primaryVertexProducer");
00049   hcalTowers_ = conf_.getParameter<edm::InputTag>("hcalTowers");
00050   hOverEConeSize_   = conf_.getParameter<double>("hOverEConeSize");
00051   highEt_        = conf_.getParameter<double>("highEt");
00052   // R9 value to decide converted/unconverted
00053   minR9Barrel_        = conf_.getParameter<double>("minR9Barrel");
00054   minR9Endcap_        = conf_.getParameter<double>("minR9Endcap");
00055   usePrimaryVertex_   = conf_.getParameter<bool>("usePrimaryVertex");
00056   runMIPTagger_       = conf_.getParameter<bool>("runMIPTagger");
00057 
00058   candidateP4type_ = config.getParameter<std::string>("candidateP4type") ;
00059  
00060   edm::ParameterSet posCalcParameters = 
00061     config.getParameter<edm::ParameterSet>("posCalcParameters");
00062   posCalculator_ = PositionCalc(posCalcParameters);
00063 
00064 
00065   //AA
00066   //Flags and Severities to be excluded from photon calculations
00067   const std::vector<std::string> flagnamesEB = 
00068     config.getParameter<std::vector<std::string> >("RecHitFlagToBeExcludedEB");
00069 
00070   const std::vector<std::string> flagnamesEE =
00071     config.getParameter<std::vector<std::string> >("RecHitFlagToBeExcludedEE");
00072 
00073   flagsexclEB_= 
00074     StringToEnumValue<EcalRecHit::Flags>(flagnamesEB);
00075 
00076   flagsexclEE_=
00077     StringToEnumValue<EcalRecHit::Flags>(flagnamesEE);
00078 
00079   const std::vector<std::string> severitynamesEB = 
00080     config.getParameter<std::vector<std::string> >("RecHitSeverityToBeExcludedEB");
00081 
00082   severitiesexclEB_= 
00083     StringToEnumValue<EcalSeverityLevel::SeverityLevel>(severitynamesEB);
00084 
00085   const std::vector<std::string> severitynamesEE = 
00086     config.getParameter<std::vector<std::string> >("RecHitSeverityToBeExcludedEE");
00087 
00088   severitiesexclEE_= 
00089     StringToEnumValue<EcalSeverityLevel::SeverityLevel>(severitynamesEE);
00090 
00091   //AA
00092 
00093   //
00094 
00095   // Parameters for the position calculation:
00096   //  std::map<std::string,double> providedParameters;
00097   // providedParameters.insert(std::make_pair("LogWeighted",conf_.getParameter<bool>("posCalc_logweight")));
00098   //providedParameters.insert(std::make_pair("T0_barl",conf_.getParameter<double>("posCalc_t0_barl")));
00099   //providedParameters.insert(std::make_pair("T0_endc",conf_.getParameter<double>("posCalc_t0_endc")));
00100   //providedParameters.insert(std::make_pair("T0_endcPresh",conf_.getParameter<double>("posCalc_t0_endcPresh")));
00101   //providedParameters.insert(std::make_pair("W0",conf_.getParameter<double>("posCalc_w0")));
00102   //providedParameters.insert(std::make_pair("X0",conf_.getParameter<double>("posCalc_x0")));
00103   //posCalculator_ = PositionCalc(providedParameters);
00104   // cut values for pre-selection
00105   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("minSCEtBarrel")); 
00106   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("maxHoverEBarrel")); 
00107   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("ecalRecHitSumEtOffsetBarrel")); 
00108   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("ecalRecHitSumEtSlopeBarrel")); 
00109   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("hcalTowerSumEtOffsetBarrel"));
00110   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("hcalTowerSumEtSlopeBarrel"));
00111   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("nTrackSolidConeBarrel"));
00112   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("nTrackHollowConeBarrel"));     
00113   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("trackPtSumSolidConeBarrel"));     
00114   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("trackPtSumHollowConeBarrel"));     
00115   preselCutValuesBarrel_.push_back(conf_.getParameter<double>("sigmaIetaIetaCutBarrel"));     
00116   //  
00117   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("minSCEtEndcap")); 
00118   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("maxHoverEEndcap")); 
00119   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("ecalRecHitSumEtOffsetEndcap")); 
00120   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("ecalRecHitSumEtSlopeEndcap")); 
00121   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("hcalTowerSumEtOffsetEndcap"));
00122   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("hcalTowerSumEtSlopeEndcap"));
00123   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("nTrackSolidConeEndcap"));
00124   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("nTrackHollowConeEndcap"));     
00125   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("trackPtSumSolidConeEndcap"));     
00126   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("trackPtSumHollowConeEndcap"));     
00127   preselCutValuesEndcap_.push_back(conf_.getParameter<double>("sigmaIetaIetaCutEndcap"));     
00128   //
00129 
00130   // Register the product
00131   produces< reco::PhotonCollection >(PhotonCollection_);
00132 
00133 }
00134 
00135 PhotonProducer::~PhotonProducer() 
00136 {
00137 
00138   //delete energyCorrectionF;
00139 }
00140 
00141 
00142 
00143 void  PhotonProducer::beginRun (edm::Run& r, edm::EventSetup const & theEventSetup) {
00144 
00145     thePhotonIsolationCalculator_ = new PhotonIsolationCalculator();
00146     edm::ParameterSet isolationSumsCalculatorSet = conf_.getParameter<edm::ParameterSet>("isolationSumsCalculatorSet"); 
00147     thePhotonIsolationCalculator_->setup(isolationSumsCalculatorSet, flagsexclEB_, flagsexclEE_, severitiesexclEB_, severitiesexclEE_);
00148 
00149     thePhotonMIPHaloTagger_ = new PhotonMIPHaloTagger();
00150     edm::ParameterSet mipVariableSet = conf_.getParameter<edm::ParameterSet>("mipVariableSet"); 
00151     thePhotonMIPHaloTagger_->setup(mipVariableSet);
00152     thePhotonEnergyCorrector_ = new PhotonEnergyCorrector(conf_);
00153     thePhotonEnergyCorrector_ -> init(theEventSetup); 
00154 }
00155 
00156 void  PhotonProducer::endRun (edm::Run& r, edm::EventSetup const & theEventSetup) {
00157 
00158   delete thePhotonIsolationCalculator_;
00159   delete thePhotonMIPHaloTagger_;
00160   delete thePhotonEnergyCorrector_;
00161 }
00162 
00163 
00164 void PhotonProducer::produce(edm::Event& theEvent, const edm::EventSetup& theEventSetup) {
00165 
00166   using namespace edm;
00167   //  nEvt_++;
00168  
00169   reco::PhotonCollection outputPhotonCollection;
00170   std::auto_ptr< reco::PhotonCollection > outputPhotonCollection_p(new reco::PhotonCollection);
00171 
00172 
00173   // Get the PhotonCore collection
00174   bool validPhotonCoreHandle=true;
00175   Handle<reco::PhotonCoreCollection> photonCoreHandle;
00176   theEvent.getByLabel(photonCoreProducer_,photonCoreHandle);
00177   if (!photonCoreHandle.isValid()) {
00178     edm::LogError("PhotonProducer") << "Error! Can't get the product "<<photonCoreProducer_.label();
00179     validPhotonCoreHandle=false;
00180   }
00181 
00182  // Get EcalRecHits
00183   bool validEcalRecHits=true;
00184   Handle<EcalRecHitCollection> barrelHitHandle;
00185   EcalRecHitCollection barrelRecHits;
00186   theEvent.getByLabel(barrelEcalHits_, barrelHitHandle);
00187   if (!barrelHitHandle.isValid()) {
00188     edm::LogError("PhotonProducer") << "Error! Can't get the product "<<barrelEcalHits_.label();
00189     validEcalRecHits=false; 
00190   }
00191   if (  validEcalRecHits)  barrelRecHits = *(barrelHitHandle.product());
00192 
00193   
00194   Handle<EcalRecHitCollection> endcapHitHandle;
00195   theEvent.getByLabel(endcapEcalHits_, endcapHitHandle);
00196   EcalRecHitCollection endcapRecHits;
00197   if (!endcapHitHandle.isValid()) {
00198     edm::LogError("PhotonProducer") << "Error! Can't get the product "<<endcapEcalHits_.label();
00199     validEcalRecHits=false; 
00200   }
00201   if( validEcalRecHits) endcapRecHits = *(endcapHitHandle.product());
00202 
00203   //AA
00204   //Get the severity level object
00205   edm::ESHandle<EcalSeverityLevelAlgo> sevLv;
00206   theEventSetup.get<EcalSeverityLevelAlgoRcd>().get(sevLv);
00207   //
00208 
00209 
00210 // get Hcal towers collection 
00211   Handle<CaloTowerCollection> hcalTowersHandle;
00212   theEvent.getByLabel(hcalTowers_, hcalTowersHandle);
00213 
00214 
00215   // get the geometry from the event setup:
00216   theEventSetup.get<CaloGeometryRecord>().get(theCaloGeom_);
00217 
00218   //
00219   // update energy correction function
00220   //  energyCorrectionF->init(theEventSetup);  
00221 
00222   edm::ESHandle<CaloTopology> pTopology;
00223   theEventSetup.get<CaloTopologyRecord>().get(theCaloTopo_);
00224   const CaloTopology *topology = theCaloTopo_.product();
00225 
00226   // Get the primary event vertex
00227   Handle<reco::VertexCollection> vertexHandle;
00228   reco::VertexCollection vertexCollection;
00229   bool validVertex=true;
00230   if ( usePrimaryVertex_ ) {
00231     theEvent.getByLabel(vertexProducer_, vertexHandle);
00232     if (!vertexHandle.isValid()) {
00233       edm::LogWarning("PhotonProducer") << "Error! Can't get the product primary Vertex Collection "<< "\n";
00234       validVertex=false;
00235     }
00236     if (validVertex) vertexCollection = *(vertexHandle.product());
00237   }
00238   //  math::XYZPoint vtx(0.,0.,0.);
00239   //if (vertexCollection.size()>0) vtx = vertexCollection.begin()->position();
00240 
00241 
00242   int iSC=0; // index in photon collection
00243   // Loop over barrel and endcap SC collections and fill the  photon collection
00244   if ( validPhotonCoreHandle) 
00245     fillPhotonCollection(theEvent,
00246                          theEventSetup,
00247                          photonCoreHandle,
00248                          topology,
00249                          &barrelRecHits,
00250                          &endcapRecHits,
00251                          hcalTowersHandle,
00252                          //vtx,
00253                          vertexCollection,
00254                          outputPhotonCollection,
00255                          iSC,
00256                          sevLv.product());
00257  
00258 
00259   // put the product in the event
00260   edm::LogInfo("PhotonProducer") << " Put in the event " << iSC << " Photon Candidates \n";
00261   outputPhotonCollection_p->assign(outputPhotonCollection.begin(),outputPhotonCollection.end());
00262   theEvent.put( outputPhotonCollection_p, PhotonCollection_);
00263 
00264 }
00265 
00266 void PhotonProducer::fillPhotonCollection(edm::Event& evt,
00267                                           edm::EventSetup const & es,
00268                                           const edm::Handle<reco::PhotonCoreCollection> & photonCoreHandle,
00269                                           const CaloTopology* topology,
00270                                           const EcalRecHitCollection* ecalBarrelHits,
00271                                           const EcalRecHitCollection* ecalEndcapHits,
00272                                           const edm::Handle<CaloTowerCollection> & hcalTowersHandle, 
00273                                           // math::XYZPoint & vtx,
00274                                           reco::VertexCollection & vertexCollection,
00275                                           reco::PhotonCollection & outputPhotonCollection, int& iSC,
00276                                           const EcalSeverityLevelAlgo * sevLv) {
00277   
00278   const CaloGeometry* geometry = theCaloGeom_.product();
00279   const CaloSubdetectorGeometry* subDetGeometry =0 ;
00280   const CaloSubdetectorGeometry* geometryES = theCaloGeom_->getSubdetectorGeometry(DetId::Ecal, EcalPreshower);
00281   const EcalRecHitCollection* hits = 0 ;
00282   std::vector<double> preselCutValues;
00283   float minR9=0;
00284 
00285 
00286   std::vector<int> flags_, severitiesexcl_;
00287 
00288   for(unsigned int lSC=0; lSC < photonCoreHandle->size(); lSC++) {
00289 
00290     reco::PhotonCoreRef coreRef(reco::PhotonCoreRef(photonCoreHandle, lSC));
00291     reco::SuperClusterRef scRef=coreRef->superCluster();
00292     //    const reco::SuperCluster* pClus=&(*scRef);
00293     iSC++;
00294 
00295     int subdet = scRef->seed()->hitsAndFractions()[0].first.subdetId();
00296     subDetGeometry =  theCaloGeom_->getSubdetectorGeometry(DetId::Ecal, subdet);
00297 
00298     if (subdet==EcalBarrel) { 
00299       preselCutValues = preselCutValuesBarrel_;
00300       minR9 = minR9Barrel_;
00301       hits = ecalBarrelHits;
00302       flags_ = flagsexclEB_;
00303       severitiesexcl_ = severitiesexclEB_;
00304     } else if  (subdet==EcalEndcap)  { 
00305       preselCutValues = preselCutValuesEndcap_;
00306       minR9 = minR9Endcap_;
00307       hits = ecalEndcapHits;
00308       flags_ = flagsexclEE_;
00309       severitiesexcl_ = severitiesexclEE_;
00310     } else {
00311       edm::LogWarning("")<<"PhotonProducer: do not know if it is a barrel or endcap SuperCluster"; 
00312     }
00313 
00314     
00315     // SC energy preselection
00316     if (scRef->energy()/cosh(scRef->eta()) <= preselCutValues[0] ) continue;
00317     // calculate HoE
00318 
00319     const CaloTowerCollection* hcalTowersColl = hcalTowersHandle.product();
00320     EgammaTowerIsolation towerIso1(hOverEConeSize_,0.,0.,1,hcalTowersColl) ;  
00321     EgammaTowerIsolation towerIso2(hOverEConeSize_,0.,0.,2,hcalTowersColl) ;  
00322     double HoE1=towerIso1.getTowerESum(&(*scRef))/scRef->energy();
00323     double HoE2=towerIso2.getTowerESum(&(*scRef))/scRef->energy(); 
00324     
00325     EgammaHadTower towerIsoBehindClus(es); 
00326     towerIsoBehindClus.setTowerCollection(hcalTowersHandle.product());
00327     std::vector<CaloTowerDetId> TowersBehindClus =  towerIsoBehindClus.towersOf(*scRef);
00328     float hcalDepth1OverEcalBc = towerIsoBehindClus.getDepth1HcalESum(TowersBehindClus)/scRef->energy();
00329     float hcalDepth2OverEcalBc = towerIsoBehindClus.getDepth2HcalESum(TowersBehindClus)/scRef->energy();
00330     //    std::cout << " PhotonProducer calculation of HoE with towers in a cone " << HoE1  << "  " << HoE2 << std::endl;
00331     //std::cout << " PhotonProducer calcualtion of HoE with towers behind the BCs " << hcalDepth1OverEcalBc  << "  " << hcalDepth2OverEcalBc << std::endl;
00332 
00333 
00334     // recalculate position of seed BasicCluster taking shower depth for unconverted photon
00335     math::XYZPoint unconvPos = posCalculator_.Calculate_Location(scRef->seed()->hitsAndFractions(),hits,subDetGeometry,geometryES);
00336 
00337 
00338     float maxXtal =   EcalClusterTools::eMax( *(scRef->seed()), &(*hits) );
00339     //AA
00340     //Change these to consider severity level of hits
00341     float e1x5    =   EcalClusterTools::e1x5(  *(scRef->seed()), &(*hits), &(*topology), flags_, severitiesexcl_, sevLv);
00342     float e2x5    =   EcalClusterTools::e2x5Max(  *(scRef->seed()), &(*hits), &(*topology),flags_, severitiesexcl_, sevLv );    
00343     float e3x3    =   EcalClusterTools::e3x3(  *(scRef->seed()), &(*hits), &(*topology), flags_, severitiesexcl_, sevLv);
00344     float e5x5    =   EcalClusterTools::e5x5( *(scRef->seed()), &(*hits), &(*topology),flags_, severitiesexcl_, sevLv);   
00345     std::vector<float> cov =  EcalClusterTools::covariances( *(scRef->seed()), &(*hits), &(*topology), geometry,flags_, severitiesexcl_, sevLv);
00346     std::vector<float> locCov =  EcalClusterTools::localCovariances( *(scRef->seed()), &(*hits), &(*topology),flags_, severitiesexcl_, sevLv);
00347       
00348     float sigmaEtaEta = sqrt(cov[0]);
00349     float sigmaIetaIeta = sqrt(locCov[0]);
00350     float r9 =e3x3/(scRef->rawEnergy());
00351 
00352     // compute position of ECAL shower
00353     math::XYZPoint caloPosition;
00354     if (r9>minR9) {
00355       caloPosition = unconvPos;
00356     } else {
00357       caloPosition = scRef->position();
00358     }
00359 
00361     double photonEnergy=1.;
00362     math::XYZPoint vtx(0.,0.,0.);
00363     if (vertexCollection.size()>0) vtx = vertexCollection.begin()->position();
00364     // compute momentum vector of photon from primary vertex and cluster position
00365     math::XYZVector direction = caloPosition - vtx;
00366     //math::XYZVector momentum = direction.unit() * photonEnergy ;
00367     math::XYZVector momentum = direction.unit() ;
00368 
00369     // Create dummy candidate with unit momentum and zero energy to allow setting of all variables. The energy is set for last.
00370     math::XYZTLorentzVectorD p4(momentum.x(), momentum.y(), momentum.z(), photonEnergy );
00371     reco::Photon newCandidate(p4, caloPosition, coreRef, vtx);
00372     //std::cout << " standard p4 before " << newCandidate.p4() << " energy " << newCandidate.energy() <<  std::endl;
00373     //std::cout << " type " <<newCandidate.getCandidateP4type() <<  " standard p4 after " << newCandidate.p4() << " energy " << newCandidate.energy() << std::endl;
00374 
00375     // Calculate fiducial flags and isolation variable. Blocked are filled from the isolationCalculator
00376     reco::Photon::FiducialFlags fiducialFlags;
00377     reco::Photon::IsolationVariables isolVarR03, isolVarR04;
00378     thePhotonIsolationCalculator_-> calculate ( &newCandidate,evt,es,fiducialFlags,isolVarR04, isolVarR03);
00379     newCandidate.setFiducialVolumeFlags( fiducialFlags );
00380     newCandidate.setIsolationVariables(isolVarR04, isolVarR03 );
00381 
00382     
00384     reco::Photon::ShowerShape  showerShape;
00385     showerShape.e1x5= e1x5;
00386     showerShape.e2x5= e2x5;
00387     showerShape.e3x3= e3x3;
00388     showerShape.e5x5= e5x5;
00389     showerShape.maxEnergyXtal =  maxXtal;
00390     showerShape.sigmaEtaEta =    sigmaEtaEta;
00391     showerShape.sigmaIetaIeta =  sigmaIetaIeta;
00392     showerShape.hcalDepth1OverEcal = HoE1;
00393     showerShape.hcalDepth2OverEcal = HoE2;
00394     showerShape.hcalDepth1OverEcalBc = hcalDepth1OverEcalBc;
00395     showerShape.hcalDepth2OverEcalBc = hcalDepth2OverEcalBc;
00396     showerShape.hcalTowersBehindClusters =  TowersBehindClus;
00397     newCandidate.setShowerShapeVariables ( showerShape ); 
00398 
00401     // Photon candidate takes by default (set in photons_cfi.py)  a 4-momentum derived from the ecal photon-specific corrections. 
00402     thePhotonEnergyCorrector_->calculate(evt, newCandidate, subdet, vertexCollection,es);
00403     if ( candidateP4type_ == "fromEcalEnergy") {
00404       newCandidate.setP4( newCandidate.p4(reco::Photon::ecal_photons) );
00405       newCandidate.setCandidateP4type(reco::Photon::ecal_photons);
00406     } else if ( candidateP4type_ == "fromRegression") {
00407       newCandidate.setP4( newCandidate.p4(reco::Photon::regression1) );
00408       newCandidate.setCandidateP4type(reco::Photon::regression1);
00409     }
00410 
00411     //       std::cout << " final p4 " << newCandidate.p4() << " energy " << newCandidate.energy() <<  std::endl;
00412 
00413 
00414     // std::cout << " PhotonProducer from candidate HoE with towers in a cone " << newCandidate.hadronicOverEm()  << "  " <<  newCandidate.hadronicDepth1OverEm()  << " " <<  newCandidate.hadronicDepth2OverEm()  << std::endl;
00415     //    std::cout << " PhotonProducer from candidate  of HoE with towers behind the BCs " <<  newCandidate.hadTowOverEm()  << "  " << newCandidate.hadTowDepth1OverEm() << " " << newCandidate.hadTowDepth2OverEm() << std::endl;
00416 
00417 
00418   // fill MIP Vairables for Halo: Block for MIP are filled from PhotonMIPHaloTagger
00419    reco::Photon::MIPVariables mipVar ;
00420    if(subdet==EcalBarrel && runMIPTagger_ )
00421     {
00422   
00423      thePhotonMIPHaloTagger_-> MIPcalculate( &newCandidate,evt,es,mipVar);
00424     newCandidate.setMIPVariables(mipVar);
00425     }
00426 
00427 
00428 
00430     bool isLooseEM=true;
00431     if ( newCandidate.pt() < highEt_) { 
00432       if ( newCandidate.hadronicOverEm()                   >= preselCutValues[1] )                                            isLooseEM=false;
00433       if ( newCandidate.ecalRecHitSumEtConeDR04()          > preselCutValues[2]+ preselCutValues[3]*newCandidate.pt() )       isLooseEM=false;
00434       if ( newCandidate.hcalTowerSumEtConeDR04()           > preselCutValues[4]+ preselCutValues[5]*newCandidate.pt() )       isLooseEM=false;
00435       if ( newCandidate.nTrkSolidConeDR04()                > int(preselCutValues[6]) )                                        isLooseEM=false;
00436       if ( newCandidate.nTrkHollowConeDR04()               > int(preselCutValues[7]) )                                        isLooseEM=false;
00437       if ( newCandidate.trkSumPtSolidConeDR04()            > preselCutValues[8] )                                             isLooseEM=false;
00438       if ( newCandidate.trkSumPtHollowConeDR04()           > preselCutValues[9] )                                             isLooseEM=false;
00439       if ( newCandidate.sigmaIetaIeta()                    > preselCutValues[10] )                                            isLooseEM=false;
00440     } 
00441     
00442 
00443         
00444     if ( isLooseEM)  
00445       outputPhotonCollection.push_back(newCandidate);
00446       
00447         
00448   }
00449 }
00450