egammaTools Namespace Reference

Typedefs
typedef math::XYZTLorentzVector	LorentzVector

Functions
ConversionInfo	arbitrateConversionPartnersbyR (const std::vector< ConversionInfo > &v_convCandidates)
double	ecalEta (const math::XYZVector &momentum, const math::XYZPoint &vertex)
double	ecalPhi (const MagneticField &magField, const math::XYZVector &momentum, const math::XYZPoint &vertex, const int charge)
ConversionInfo	findBestConversionMatch (const std::vector< ConversionInfo > &v_convCandidates)
ConversionInfo	getConversionInfo (const reco::Track *el_track, const reco::Track *candPartnerTk, const double bFieldAtOrigin)
ConversionInfo	getConversionInfo (const reco::GsfElectronCore &, const edm::Handle< reco::TrackCollection > &ctftracks_h, const edm::Handle< reco::GsfTrackCollection > &gsftracks_h, const double bFieldAtOrigin, const double minFracSharedHits=0.45)
ConversionInfo	getConversionInfo (const reco::GsfElectron &gsfElectron, const edm::Handle< reco::TrackCollection > &ctftracks_h, const edm::Handle< reco::GsfTrackCollection > &gsftracks_h, const double bFieldAtOrigin, const double minFracSharedHits=0.45)
ConversionInfo	getConversionInfo (const reco::GsfElectron &gsfElectron, const edm::Handle< reco::TrackCollection > &track_h, const double bFieldAtOrigin, const double minFracSharedHits=0.45)
std::pair< double, double >	getConversionInfo (LorentzVector trk1_p4, int trk1_q, float trk1_d0, LorentzVector trk2_p4, int trk2_q, float trk2_d0, float bFieldAtOrigin)
std::vector< ConversionInfo >	getConversionInfos (const reco::GsfElectronCore &, const edm::Handle< reco::TrackCollection > &ctftracks_h, const edm::Handle< reco::GsfTrackCollection > &gsftracks_h, const double bFieldAtOrigin, const double minFracSharedHits=0.45)
const reco::Track *	getElectronTrack (const reco::GsfElectron &, const float minFracSharedHits=0.45)
const reco::Track *	getElectronTrack (const reco::GsfElectronCore &, const float minFracSharedHits=0.45)
bool	isFromConversion (const ConversionInfo &, double maxAbsDist=0.02, double maxAbsDcot=0.02)
void	localEcalClusterCoordsEB (const reco::CaloCluster &bclus, const CaloGeometry &geom, float &etacry, float &phicry, int &ieta, int &iphi, float &thetatilt, float &phitilt)
void	localEcalClusterCoordsEE (const reco::CaloCluster &bclus, const CaloGeometry &geom, float &xcry, float &ycry, int &ix, int &iy, float &thetatilt, float &phitilt)

Typedef Documentation

typedef math::XYZTLorentzVector egammaTools::LorentzVector

Definition at line 10 of file ConversionFinder.cc.

Function Documentation

ConversionInfo egammaTools::arbitrateConversionPartnersbyR

(

const std::vector< ConversionInfo > &

v_convCandidates

)

Definition at line 304 of file ConversionFinder.cc.

References mps_fire::i, funct::pow(), dttmaxenums::R, mathSSE::sqrt(), and groupFilesInBlocks::temp.

Referenced by findBestConversionMatch().

                                                                                                 {

  if(v_convCandidates.size() == 1)
    return v_convCandidates.at(0);

  double R = sqrt(pow(v_convCandidates.at(0).dist,2) + pow(v_convCandidates.at(0).dcot,2));

  int iArbitrated = 0;
  int i = 0;

  for(auto const& temp : v_convCandidates) {
    double temp_R = sqrt(pow(temp.dist,2) + pow(temp.dcot,2));
    if(temp_R < R) {
      R = temp_R;
      iArbitrated = i;
    }
    ++i;
  }

  return v_convCandidates.at(iArbitrated);

 }

double egammaTools::ecalEta	(	const math::XYZVector &	momentum,
		const math::XYZPoint &	vertex
	)

Definition at line 71 of file ECALPositionCalculator.cc.

References ETA, etaBarrelEndcap, cmsBatch::log, Geom::pi(), R_ECAL, funct::tan(), theta(), and Z_Endcap.

Referenced by ContainmentCorrectionAnalyzer::analyze(), and PreshowerAndECALLinker::testLink().

{

   // Get kinematic variables

   float etaParticle = momentum.eta();
   float vZ = vertex.z();
   float vRho = vertex.Rho();

   if (etaParticle != 0.0)
   {
      float theta = 0.0;
      float zEcal = (R_ECAL-vRho)*sinh(etaParticle)+vZ;
      
      if(zEcal != 0.0) theta = atan(R_ECAL/zEcal);
      if(theta < 0.0) theta = theta + Geom::pi() ;

      float ETA = - log(tan(0.5*theta));
      
      if( fabs(ETA) > etaBarrelEndcap )
      {
         float Zend = Z_Endcap ;
     if(etaParticle < 0.0 )  Zend = - Zend;
     float Zlen = Zend - vZ ;
     float RR = Zlen/sinh(etaParticle);
     theta = atan((RR+vRho)/Zend);
     if(theta < 0.0) theta = theta+Geom::pi();
     ETA = -log(tan(0.5*theta));
      }
      return ETA;

    } else {
       edm::LogWarning("")  << "[EcalPositionFromTrack::etaTransformation] Warning: "
                << "Eta equals to zero, not correcting";
       return etaParticle;
    }

}

double egammaTools::ecalPhi	(	const MagneticField &	magField,
		const math::XYZVector &	momentum,
		const math::XYZPoint &	vertex,
		const int	charge
	)

Definition at line 16 of file ECALPositionCalculator.cc.

References ALCARECOTkAlJpsiMuMu_cff::charge, etaBarrelEndcap, MagneticField::inTesla(), Geom::pi(), Geom::twoPi(), and PV3DBase< T, PVType, FrameType >::z().

Referenced by EgammaHLTElectronDetaDphiProducer::calDEtaDPhiSCTrk(), MuonMETAlgo::correctMETforMuon(), and PreshowerAndECALLinker::testLink().

{

   // Get kinematic variables

   float ptParticle = momentum.Rho();
   float etaParticle = momentum.eta();
   float phiParticle = momentum.phi();
   float vRho = vertex.Rho();

   // Magnetic field

   const float RBARM = 1.357 ;      // was 1.31 , updated on 16122003
   const float ZENDM = 3.186 ;      // was 3.15 , updated on 16122003

   float rbend = RBARM-(vRho/100.0);    //Assumed vRho in cm
   float bend  = 0.3 * magField.inTesla(GlobalPoint(0.,0.,0.)).z() * rbend / 2.0; 
   float phi = 0.0;

   if( fabs(etaParticle) <=  etaBarrelEndcap)
   {
      if (fabs(bend/ptParticle) <= 1.) 
      {
         phi = phiParticle - asin(bend/ptParticle)*charge;
     if(phi >  Geom::pi()) phi = phi - Geom::twoPi();
     if(phi < -Geom::pi()) phi = phi + Geom::twoPi();
      } else {
         edm::LogWarning("") << "[EcalPositionFromTrack::phiTransformation] Warning: "
                     << "Too low Pt, giving up";
     return phiParticle;
      }

   } // end if in the barrel
  
   if(fabs(etaParticle) > etaBarrelEndcap)
   {
      float rHit = 0.0;
      rHit = ZENDM / sinh(fabs(etaParticle));
      if (fabs(((rHit-(vRho/100.0))/rbend)*bend/ptParticle) <= 1.0)
      {
         phi = phiParticle - asin(((rHit-(vRho/100.0)) / rbend)*bend/ptParticle)*charge;
     if(phi >  Geom::pi()) phi = phi - Geom::twoPi();
     if(phi < -Geom::pi()) phi = phi + Geom::twoPi();
      } else {
         edm::LogWarning("") << "[EcalPositionFromTrack::phiTransformation] Warning: "
                     << "Too low Pt, giving up";
     return phiParticle;
      }  

    } // end if in the endcap
  
   return phi;

}

ConversionInfo egammaTools::findBestConversionMatch

(

const std::vector< ConversionInfo > &

v_convCandidates

)

Definition at line 328 of file ConversionFinder.cc.

References arbitrateConversionPartnersbyR(), ConversionInfo::dcot, ConversionInfo::deltaMissingHits, ConversionInfo::dist, ConversionInfo::flag, mps_fire::i, funct::pow(), ConversionInfo::radiusOfConversion, mathSSE::sqrt(), and groupFilesInBlocks::temp.

Referenced by getConversionInfo().

 {
  using namespace std;

  if(v_convCandidates.empty())
    return   ConversionInfo{-9999.,-9999.,-9999.,
                math::XYZPoint(-9999.,-9999.,-9999),
                reco::TrackRef(), reco::GsfTrackRef(),
                -9999, -9999};


  if(v_convCandidates.size() == 1)
    return v_convCandidates.at(0);

  vector<ConversionInfo> v_0;
  vector<ConversionInfo> v_1;
  vector<ConversionInfo> v_2;
  vector<ConversionInfo> v_3;
  //loop over the candidates
  for(unsigned int i = 1; i < v_convCandidates.size(); i++) {
    ConversionInfo temp = v_convCandidates.at(i);

    if(temp.flag == 0) {
      bool isConv = false;
      if(fabs(temp.dist) < 0.02 &&
     fabs(temp.dcot) < 0.02 &&
     temp.deltaMissingHits < 3 &&
     temp.radiusOfConversion > -2)
    isConv = true;
      if(sqrt(pow(temp.dist,2) + pow(temp.dcot,2)) < 0.05  &&
     temp.deltaMissingHits < 2 &&
     temp.radiusOfConversion > -2)
    isConv = true;

      if(isConv)
    v_0.push_back(temp);
    }

    if(temp.flag == 1) {

      if(sqrt(pow(temp.dist,2) + pow(temp.dcot,2)) < 0.05  &&
     temp.deltaMissingHits < 2 &&
     temp.radiusOfConversion > -2)
    v_1.push_back(temp);
    }
    if(temp.flag == 2) {

      if(sqrt(pow(temp.dist,2) + pow(temp.dcot*temp.dcot,2)) < 0.05 &&
     temp.deltaMissingHits < 2 &&
     temp.radiusOfConversion > -2)
    v_2.push_back(temp);

    }
    if(temp.flag == 3) {

      if(sqrt(temp.dist*temp.dist + temp.dcot*temp.dcot) < 0.05
     && temp.deltaMissingHits < 2
     && temp.radiusOfConversion > -2)
    v_3.push_back(temp);

    }

  }//candidate conversion loop

  //now do some arbitration

  //give preference to conversion partners found in the CTF collection
  //using the electron's CTF track
  if(!v_0.empty())
    return arbitrateConversionPartnersbyR(v_0);

  if(!v_1.empty())
    return arbitrateConversionPartnersbyR(v_1);

  if(!v_2.empty())
    return arbitrateConversionPartnersbyR(v_2);

  if(!v_3.empty())
    return arbitrateConversionPartnersbyR(v_3);


  //if we get here, we didn't find a candidate conversion partner that
  //satisfied even the loose selections
  //return the the closest partner by R
  return arbitrateConversionPartnersbyR(v_convCandidates);

 }

ConversionInfo egammaTools::getConversionInfo	(	const reco::Track *	el_track,
		const reco::Track *	candPartnerTk,
		const double	bFieldAtOrigin
	)

Definition at line 242 of file ConversionFinder.cc.

References reco::TrackBase::charge(), funct::cos(), edmIntegrityCheck::d, reco::TrackBase::d0(), reco::TrackBase::dz(), reco::TrackBase::p(), funct::pow(), reco::TrackBase::pt(), reco::TrackBase::px(), reco::TrackBase::py(), reco::TrackBase::pz(), funct::sin(), mathSSE::sqrt(), and funct::tan().

                                                        {

  using namespace reco;

  //now calculate the conversion related information
  LorentzVector el_tk_p4(el_track->px(), el_track->py(), el_track->pz(), el_track->p());
  double elCurvature = -0.3*bFieldAtOrigin*(el_track->charge()/el_tk_p4.pt())/100.;
  double rEl = fabs(1./elCurvature);
  double xEl = -1*(1./elCurvature - el_track->d0())*sin(el_tk_p4.phi());
  double yEl = (1./elCurvature - el_track->d0())*cos(el_tk_p4.phi());


  LorentzVector cand_p4 = LorentzVector(candPartnerTk->px(), candPartnerTk->py(),candPartnerTk->pz(), candPartnerTk->p());
  double candCurvature = -0.3*bFieldAtOrigin*(candPartnerTk->charge()/cand_p4.pt())/100.;
  double rCand = fabs(1./candCurvature);
  double xCand = -1*(1./candCurvature - candPartnerTk->d0())*sin(cand_p4.phi());
  double yCand = (1./candCurvature - candPartnerTk->d0())*cos(cand_p4.phi());

  double d = sqrt(pow(xEl-xCand, 2) + pow(yEl-yCand , 2));
  double dist = d - (rEl + rCand);
  double dcot = 1./tan(el_tk_p4.theta()) - 1./tan(cand_p4.theta());

  //get the point of conversion
  double xa1 = xEl   + (xCand-xEl) * rEl/d;
  double xa2 = xCand + (xEl-xCand) * rCand/d;
  double ya1 = yEl   + (yCand-yEl) * rEl/d;
  double ya2 = yCand + (yEl-yCand) * rCand/d;

  double x=.5*(xa1+xa2);
  double y=.5*(ya1+ya2);
  double rconv = sqrt(pow(x,2) + pow(y,2));
  double z = el_track->dz() + rEl*el_track->pz()*TMath::ACos(1-pow(rconv,2)/(2.*pow(rEl,2)))/el_track->pt();

  math::XYZPoint convPoint(x, y, z);

  //now assign a sign to the radius of conversion
  float tempsign = el_track->px()*x + el_track->py()*y;
  tempsign = tempsign/fabs(tempsign);
  rconv = tempsign*rconv;

  //return an instance of ConversionInfo, but with a NULL track refs
  return ConversionInfo{dist, dcot, rconv, convPoint, TrackRef(), GsfTrackRef(), -9999, -9999};

}

ConversionInfo egammaTools::getConversionInfo	(	const reco::GsfElectronCore &	gsfElectron,
		const edm::Handle< reco::TrackCollection > &	ctftracks_h,
		const edm::Handle< reco::GsfTrackCollection > &	gsftracks_h,
		const double	bFieldAtOrigin,
		const double	minFracSharedHits = 0.45
	)

Definition at line 38 of file ConversionFinder.cc.

References findBestConversionMatch(), getConversionInfos(), and groupFilesInBlocks::temp.

Referenced by GsfElectronAlgo::createElectron(), AdHocNTupler::fill(), ZeeCandidateFilter::filter(), getConversionInfos(), and ElectronConversionRejectionVars::produce().

                                                            {

  std::vector<ConversionInfo> temp = getConversionInfos(gsfElectron,ctftracks_h,gsftracks_h,bFieldAtOrigin,minFracSharedHits) ;
  return findBestConversionMatch(temp);

}

ConversionInfo egammaTools::getConversionInfo	(	const reco::GsfElectron &	gsfElectron,
		const edm::Handle< reco::TrackCollection > &	ctftracks_h,
		const edm::Handle< reco::GsfTrackCollection > &	gsftracks_h,
		const double	bFieldAtOrigin,
		const double	minFracSharedHits = 0.45
	)

Definition at line 27 of file ConversionFinder.cc.

References reco::GsfElectron::core(), findBestConversionMatch(), getConversionInfos(), and groupFilesInBlocks::temp.

                                                                {

  std::vector<ConversionInfo> temp = getConversionInfos(*gsfElectron.core(),ctftracks_h,gsftracks_h,bFieldAtOrigin,minFracSharedHits) ;
  return findBestConversionMatch(temp);

}

ConversionInfo egammaTools::getConversionInfo	(	const reco::GsfElectron &	gsfElectron,
		const edm::Handle< reco::TrackCollection > &	track_h,
		const double	bFieldAtOrigin,
		const double	minFracSharedHits = 0.45
	)

Definition at line 453 of file ConversionFinder.cc.

References reco::TrackBase::charge(), reco::GsfElectron::closestCtfTrackRef(), funct::cos(), edmIntegrityCheck::d, reco::TrackBase::d0(), boostedElectronIsolation_cff::deltaR, PFRecoTauDiscriminationAgainstElectronDeadECAL_cfi::dR, reco::TrackBase::dz(), RemoveAddSevLevel::flag, getElectronTrack(), reco::TrackBase::hitPattern(), edm::Ref< C, T, F >::isNonnull(), edm::Ref< C, T, F >::key(), reco::HitPattern::MISSING_INNER_HITS, reco::HitPattern::numberOfLostHits(), reco::TrackBase::p(), funct::pow(), edm::Handle< T >::product(), reco::TrackBase::pt(), reco::TrackBase::px(), reco::TrackBase::py(), reco::TrackBase::pz(), reco::GsfElectron::shFracInnerHits(), funct::sin(), mathSSE::sqrt(), and funct::tan().

                                                           {


  using namespace reco;
  using namespace std;
  using namespace edm;


  const TrackCollection *ctftracks = track_h.product();
  const reco::TrackRef el_ctftrack = gsfElectron.closestCtfTrackRef();
  int ctfidx = -999.;
  int flag = -9999.;
  if(el_ctftrack.isNonnull() && gsfElectron.shFracInnerHits() > minFracSharedHits) {
    ctfidx = static_cast<int>(el_ctftrack.key());
    flag = 0;
  } else
    flag = 1;


  /*
    determine whether we're going to use the CTF track or the GSF track
    using the electron's CTF track to find the dist, dcot has been shown
    to reduce the inefficiency
  */
  const reco::Track* el_track = getElectronTrack(gsfElectron, minFracSharedHits);
  LorentzVector el_tk_p4(el_track->px(), el_track->py(), el_track->pz(), el_track->p());


  int tk_i = 0;
  double mindcot = 9999.;
  //make a null Track Ref
  TrackRef candCtfTrackRef = TrackRef() ;


  for(TrackCollection::const_iterator tk = ctftracks->begin();
      tk != ctftracks->end(); tk++, tk_i++) {
    //if the general Track is the same one as made by the electron, skip it
    if(tk_i == ctfidx)
      continue;

    LorentzVector tk_p4 = LorentzVector(tk->px(), tk->py(),tk->pz(), tk->p());

    //look only in a cone of 0.5
    double dR = deltaR(el_tk_p4, tk_p4);
    if(dR>0.5)
      continue;


    //require opp. sign -> Should we use the majority logic??
    if(tk->charge() + el_track->charge() != 0)
      continue;

    double dcot = fabs(1./tan(tk_p4.theta()) - 1./tan(el_tk_p4.theta()));
    if(dcot < mindcot) {
      mindcot = dcot;
      candCtfTrackRef = reco::TrackRef(track_h, tk_i);
    }
  }//track loop


  if(!candCtfTrackRef.isNonnull())
    return ConversionInfo{-9999.,-9999.,-9999.,
              math::XYZPoint(-9999.,-9999.,-9999),
              reco::TrackRef(), reco::GsfTrackRef(),
              -9999, -9999};



  //now calculate the conversion related information
  double elCurvature = -0.3*bFieldAtOrigin*(el_track->charge()/el_tk_p4.pt())/100.;
  double rEl = fabs(1./elCurvature);
  double xEl = -1*(1./elCurvature - el_track->d0())*sin(el_tk_p4.phi());
  double yEl = (1./elCurvature - el_track->d0())*cos(el_tk_p4.phi());


  LorentzVector cand_p4 = LorentzVector(candCtfTrackRef->px(), candCtfTrackRef->py(),candCtfTrackRef->pz(), candCtfTrackRef->p());
  double candCurvature = -0.3*bFieldAtOrigin*(candCtfTrackRef->charge()/cand_p4.pt())/100.;
  double rCand = fabs(1./candCurvature);
  double xCand = -1*(1./candCurvature - candCtfTrackRef->d0())*sin(cand_p4.phi());
  double yCand = (1./candCurvature - candCtfTrackRef->d0())*cos(cand_p4.phi());

  double d = sqrt(pow(xEl-xCand, 2) + pow(yEl-yCand , 2));
  double dist = d - (rEl + rCand);
  double dcot = 1./tan(el_tk_p4.theta()) - 1./tan(cand_p4.theta());

  //get the point of conversion
  double xa1 = xEl   + (xCand-xEl) * rEl/d;
  double xa2 = xCand + (xEl-xCand) * rCand/d;
  double ya1 = yEl   + (yCand-yEl) * rEl/d;
  double ya2 = yCand + (yEl-yCand) * rCand/d;

  double x=.5*(xa1+xa2);
  double y=.5*(ya1+ya2);
  double rconv = sqrt(pow(x,2) + pow(y,2));
  double z = el_track->dz() + rEl*el_track->pz()*TMath::ACos(1-pow(rconv,2)/(2.*pow(rEl,2)))/el_track->pt();

  math::XYZPoint convPoint(x, y, z);

  //now assign a sign to the radius of conversion
  float tempsign = el_track->px()*x + el_track->py()*y;
  tempsign = tempsign/fabs(tempsign);
  rconv = tempsign*rconv;

  int deltaMissingHits = -9999;

  deltaMissingHits = candCtfTrackRef->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS)
      - el_track->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);

  return ConversionInfo{dist, dcot, rconv, convPoint, candCtfTrackRef, GsfTrackRef(), deltaMissingHits, flag};

 }

std::pair<double, double> egammaTools::getConversionInfo	(	LorentzVector	trk1_p4,
		int	trk1_q,
		float	trk1_d0,
		LorentzVector	trk2_p4,
		int	trk2_q,
		float	trk2_d0,
		float	bFieldAtOrigin
	)

Definition at line 423 of file ConversionFinder.cc.

References funct::cos(), funct::pow(), funct::sin(), mathSSE::sqrt(), and funct::tan().

                                                        {


  double tk1Curvature = -0.3*bFieldAtOrigin*(trk1_q/trk1_p4.pt())/100.;
  double rTk1 = fabs(1./tk1Curvature);
  double xTk1 = -1.*(1./tk1Curvature - trk1_d0)*sin(trk1_p4.phi());
  double yTk1 = (1./tk1Curvature - trk1_d0)*cos(trk1_p4.phi());


  double tk2Curvature = -0.3*bFieldAtOrigin*(trk2_q/trk2_p4.pt())/100.;
  double rTk2 = fabs(1./tk2Curvature);
  double xTk2 = -1.*(1./tk2Curvature - trk2_d0)*sin(trk2_p4.phi());
  double yTk2 = (1./tk2Curvature - trk2_d0)*cos(trk2_p4.phi());


  double dist = sqrt(pow(xTk1-xTk2, 2) + pow(yTk1-yTk2 , 2));
  dist = dist - (rTk1 + rTk2);

  double dcot = 1./tan(trk1_p4.theta()) - 1./tan(trk2_p4.theta());

  return std::make_pair(dist, dcot);

 }

std::vector< ConversionInfo > egammaTools::getConversionInfos	(	const reco::GsfElectronCore &	gsfElectron,
		const edm::Handle< reco::TrackCollection > &	ctftracks_h,
		const edm::Handle< reco::GsfTrackCollection > &	gsftracks_h,
		const double	bFieldAtOrigin,
		const double	minFracSharedHits = 0.45
	)

Definition at line 51 of file ConversionFinder.cc.

References reco::GsfElectronCore::ctfGsfOverlap(), reco::GsfElectronCore::ctfTrack(), ConversionInfo::dcot, boostedElectronIsolation_cff::deltaR, ConversionInfo::dist, edm::Ref< C, T, F >::get(), getConversionInfo(), reco::GsfElectronCore::gsfTrack(), edm::HandleBase::id(), edm::Ref< C, T, F >::id(), edm::Ref< C, T, F >::isNonnull(), edm::Ref< C, T, F >::key(), reco::HitPattern::MISSING_INNER_HITS, ConversionInfo::pointOfConversion, edm::Handle< T >::product(), and ConversionInfo::radiusOfConversion.

Referenced by getConversionInfo().

                                                                {



  using namespace reco;
  using namespace std;
  using namespace edm;


  //get the track collections
  const TrackCollection    *ctftracks   = ctftracks_h.product();
  const GsfTrackCollection *gsftracks   = gsftracks_h.product();

  //get the references to the gsf and ctf tracks that are made
  //by the electron
  const reco::TrackRef    el_ctftrack   = gsfElectron.ctfTrack();
  const reco::GsfTrackRef& el_gsftrack  = gsfElectron.gsfTrack();

  //protect against the wrong collection being passed to the function
  if(el_ctftrack.isNonnull() && el_ctftrack.id() != ctftracks_h.id())
    throw cms::Exception("ConversionFinderError") << "ProductID of ctf track collection does not match ProductID of electron's CTF track! \n";
  if(el_gsftrack.isNonnull() && el_gsftrack.id() != gsftracks_h.id())
    throw cms::Exception("ConversionFinderError") << "ProductID of gsf track collection does not match ProductID of electron's GSF track! \n";

  //make p4s for the electron's tracks for use later
  LorentzVector el_ctftrack_p4;
  if(el_ctftrack.isNonnull() && gsfElectron.ctfGsfOverlap() > minFracSharedHits)
    el_ctftrack_p4 = LorentzVector(el_ctftrack->px(), el_ctftrack->py(), el_ctftrack->pz(), el_ctftrack->p());
  LorentzVector el_gsftrack_p4(el_gsftrack->px(), el_gsftrack->py(), el_gsftrack->pz(), el_gsftrack->p());

  //the electron's CTF track must share at least 45% of the inner hits
  //with the electron's GSF track
  int ctfidx = -999.;
  int gsfidx = -999.;
  if(el_ctftrack.isNonnull() && gsfElectron.ctfGsfOverlap() > minFracSharedHits)
    ctfidx = static_cast<int>(el_ctftrack.key());

  gsfidx = static_cast<int>(el_gsftrack.key());


  //these vectors are for those candidate partner tracks that pass our cuts
  vector<ConversionInfo> v_candidatePartners;
  //track indices required to make references
  int ctftk_i = 0;
  int gsftk_i = 0;


  //loop over the CTF tracks and try to find the partner track
  for(TrackCollection::const_iterator ctftk = ctftracks->begin();
      ctftk != ctftracks->end(); ctftk++, ctftk_i++) {

    if(ctftk_i == ctfidx)
      continue;

    //candidate track's p4
    LorentzVector ctftk_p4 = LorentzVector(ctftk->px(), ctftk->py(), ctftk->pz(), ctftk->p());

    //apply quality cuts to remove bad tracks
    if(ctftk->ptError()/ctftk->pt() > 0.05)
      continue;
    if(ctftk->numberOfValidHits() < 5)
      continue;

    if(el_ctftrack.isNonnull() &&
       gsfElectron.ctfGsfOverlap() > minFracSharedHits &&
       fabs(ctftk_p4.Pt() - el_ctftrack->pt())/el_ctftrack->pt() < 0.2)
      continue;

    //use the electron's CTF track, if not null, to search for the partner track
    //look only in a cone of 0.5 to save time, and require that the track is opp. sign
    if(el_ctftrack.isNonnull() && gsfElectron.ctfGsfOverlap() > minFracSharedHits &&
       deltaR(el_ctftrack_p4, ctftk_p4) < 0.5 &&
       (el_ctftrack->charge() + ctftk->charge() == 0) ) {

      ConversionInfo convInfo = getConversionInfo((const reco::Track*)(el_ctftrack.get()), &(*ctftk), bFieldAtOrigin);

      //need to add the track reference information for completeness
      //because the overloaded fnc above does not make a trackRef
      int deltaMissingHits = ctftk->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS)
          - el_ctftrack->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);


      v_candidatePartners.push_back({convInfo.dist,
                                     convInfo.dcot,
                                     convInfo.radiusOfConversion,
                                     convInfo.pointOfConversion,
                                     TrackRef(ctftracks_h, ctftk_i),
                                     GsfTrackRef() ,
                                     deltaMissingHits,
                                     0});

    }//using the electron's CTF track


    //now we check using the electron's gsf track
    if(deltaR(el_gsftrack_p4, ctftk_p4) < 0.5 &&
       (el_gsftrack->charge() + ctftk->charge() == 0) &&
       el_gsftrack->ptError()/el_gsftrack->pt() < 0.25) {

      int deltaMissingHits = ctftk->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS)
          - el_gsftrack->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);

      ConversionInfo convInfo = getConversionInfo((const reco::Track*)(el_gsftrack.get()), &(*ctftk), bFieldAtOrigin);

      v_candidatePartners.push_back({convInfo.dist,
                                     convInfo.dcot,
                                     convInfo.radiusOfConversion,
                                     convInfo.pointOfConversion,
                                     TrackRef(ctftracks_h, ctftk_i),
                                     GsfTrackRef(),
                                     deltaMissingHits,
                                     1});
    }//using the electron's GSF track

  }//loop over the CTF track collection


  //------------------------------------------------------ Loop over GSF collection ----------------------------------//
  for(GsfTrackCollection::const_iterator gsftk = gsftracks->begin();
      gsftk != gsftracks->end(); gsftk++, gsftk_i++) {

    //reject the electron's own gsfTrack
    if(gsfidx == gsftk_i)
      continue;

    LorentzVector gsftk_p4 = LorentzVector(gsftk->px(), gsftk->py(), gsftk->pz(), gsftk->p());

    //apply quality cuts to remove bad tracks
    if(gsftk->ptError()/gsftk->pt() > 0.5)
      continue;
    if(gsftk->numberOfValidHits() < 5)
      continue;

    if(fabs(gsftk->pt() - el_gsftrack->pt())/el_gsftrack->pt() < 0.25)
      continue;

    //try using the electron's CTF track first if it exists
    //look only in a cone of 0.5 around the electron's track
    //require opposite sign
    if(el_ctftrack.isNonnull() && gsfElectron.ctfGsfOverlap() > minFracSharedHits &&
       deltaR(el_ctftrack_p4, gsftk_p4) < 0.5 &&
       (el_ctftrack->charge() + gsftk->charge() == 0)) {

      int deltaMissingHits = gsftk->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS)
          - el_ctftrack->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);
      
      ConversionInfo convInfo = getConversionInfo((const reco::Track*)(el_ctftrack.get()), (const reco::Track*)(&(*gsftk)), bFieldAtOrigin);
      //fill the Ref info
      v_candidatePartners.push_back({convInfo.dist,
                                     convInfo.dcot,
                                     convInfo.radiusOfConversion,
                                     convInfo.pointOfConversion,
                                     TrackRef(),
                                     GsfTrackRef(gsftracks_h, gsftk_i),
                                     deltaMissingHits,
                                     2});
    }

    //use the electron's gsf track
    if(deltaR(el_gsftrack_p4, gsftk_p4) < 0.5 &&
       (el_gsftrack->charge() + gsftk->charge() == 0) &&
       (el_gsftrack->ptError()/el_gsftrack_p4.pt() < 0.5)) {
      ConversionInfo convInfo = getConversionInfo((const reco::Track*)(el_gsftrack.get()), (const reco::Track*)(&(*gsftk)), bFieldAtOrigin);
      //fill the Ref info

      int deltaMissingHits = gsftk->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS)
          - el_gsftrack->hitPattern().numberOfLostHits(reco::HitPattern::MISSING_INNER_HITS);


      v_candidatePartners.push_back({convInfo.dist,
                                     convInfo.dcot,
                                     convInfo.radiusOfConversion,
                                     convInfo.pointOfConversion,
                                     TrackRef(),
                                     GsfTrackRef(gsftracks_h, gsftk_i),
                                     deltaMissingHits,
                                     3});
    }
  }//loop over the gsf track collection


  return v_candidatePartners;

}

const reco::Track * egammaTools::getElectronTrack	(	const reco::GsfElectron &	electron,
		const float	minFracSharedHits = 0.45
	)

Definition at line 290 of file ConversionFinder.cc.

References reco::GsfElectron::closestCtfTrackRef(), edm::Ref< C, T, F >::get(), reco::GsfElectron::gsfTrack(), edm::Ref< C, T, F >::isNonnull(), and reco::GsfElectron::shFracInnerHits().

Referenced by getConversionInfo().

                                                                                                  {

  if(electron.closestCtfTrackRef().isNonnull() &&
     electron.shFracInnerHits() > minFracSharedHits)
    return (const reco::Track*)electron.closestCtfTrackRef().get();

  return (const reco::Track*)(electron.gsfTrack().get());
}

const reco::Track* egammaTools::getElectronTrack	(	const reco::GsfElectronCore &	,
		const float	minFracSharedHits = 0.45
	)

bool egammaTools::isFromConversion	(	const ConversionInfo &	convInfo,
		double	maxAbsDist = 0.02,
		double	maxAbsDcot = 0.02
	)

Definition at line 20 of file ConversionFinder.cc.

References funct::abs(), ConversionInfo::dcot, and ConversionInfo::dist.

{
  return (std::abs(convInfo.dist) < maxAbsDist) && (std::abs(convInfo.dcot) < maxAbsDcot);
}

void egammaTools::localEcalClusterCoordsEB	(	const reco::CaloCluster &	bclus,
		const CaloGeometry &	geom,
		float &	etacry,
		float &	phicry,
		int &	ieta,
		int &	iphi,
		float &	thetatilt,
		float &	phitilt
	)

Definition at line 14 of file EcalClusterLocal.cc.

References reco::deltaR(), egammaForCoreTracking_cff::depth, runTauDisplay::dr, DetId::Ecal, EcalBarrel, reco::CaloCluster::energy(), reco::tau::disc::Eta(), PV3DBase< T, PVType, FrameType >::eta(), plotBeamSpotDB::first, relativeConstraints::geom, TruncatedPyramid::getPhiAxis(), TruncatedPyramid::getPosition(), CaloGeometry::getSubdetectorGeometry(), TruncatedPyramid::getThetaAxis(), reco::CaloCluster::hitsAndFractions(), EBDetId::ieta(), EBDetId::iphi(), cmsBatch::log, colinearityKinematic::Phi, PV3DBase< T, PVType, FrameType >::phi(), Pi, reco::CaloCluster::position(), PV3DBase< T, PVType, FrameType >::theta(), and X0.

Referenced by EGEnergyCorrector::CorrectedEnergyWithError(), EGEnergyCorrector::CorrectedEnergyWithErrorV3(), EcalRegressionData::fill(), EGRegressionModifierV3::getSeedCrysCoord(), SuperClusterHelper::localCoordinates(), EGRegressionModifierV2::modifyObject(), and EGRegressionModifierV1::modifyObject().

{
  
  assert(bclus.hitsAndFractions().at(0).first.subdetId()==EcalBarrel);
  
  const CaloSubdetectorGeometry* geom=caloGeometry.getSubdetectorGeometry(DetId::Ecal,EcalBarrel);//EcalBarrel = 1
  
  const math::XYZPoint& position_ = bclus.position(); 
  double Theta = -position_.theta()+0.5*TMath::Pi();
  double Eta = position_.eta();
  double Phi = TVector2::Phi_mpi_pi(position_.phi());
  
  //Calculate expected depth of the maximum shower from energy (like in PositionCalc::Calculate_Location()):
  // The parameters X0 and T0 are hardcoded here because these values were used to calculate the corrections:
  const float X0 = 0.89; const float T0 = 7.4;
  double depth = X0 * (T0 + log(bclus.energy()));
  
  //find max energy crystal
  std::vector< std::pair<DetId, float> > crystals_vector = bclus.hitsAndFractions();
  float drmin = 999.;
  EBDetId crystalseed;
  //printf("starting loop over crystals, etot = %5f:\n",bclus.energy());
  for (unsigned int icry=0; icry!=crystals_vector.size(); ++icry) {    
    
    EBDetId crystal(crystals_vector[icry].first);
        
    auto cell=geom->getGeometry(crystal);
    const TruncatedPyramid* cpyr = dynamic_cast<const TruncatedPyramid*>(cell.get());
    GlobalPoint center_pos = cpyr->getPosition(depth);
    double EtaCentr = center_pos.eta();
    double PhiCentr = TVector2::Phi_mpi_pi(center_pos.phi());

    float dr = reco::deltaR(Eta,Phi,EtaCentr,PhiCentr);
    if (dr<drmin) {
      drmin = dr;
      crystalseed = crystal;
    }
    
  }
  
  ieta = crystalseed.ieta();
  iphi = crystalseed.iphi();
  
  // Get center cell position from shower depth
  auto cell=geom->getGeometry(crystalseed);
  const TruncatedPyramid* cpyr = dynamic_cast<const TruncatedPyramid*>(cell.get());

  thetatilt = cpyr->getThetaAxis();
  phitilt = cpyr->getPhiAxis();

  GlobalPoint center_pos = cpyr->getPosition(depth);
  
  double PhiCentr = TVector2::Phi_mpi_pi(center_pos.phi());
  double PhiWidth = (TMath::Pi()/180.);
  phicry = (TVector2::Phi_mpi_pi(Phi-PhiCentr))/PhiWidth;
  //Some flips to take into account ECAL barrel symmetries:
  if (ieta<0) phicry *= -1.;  
  
  double ThetaCentr = -center_pos.theta()+0.5*TMath::Pi();
  double ThetaWidth = (TMath::Pi()/180.)*TMath::Cos(ThetaCentr);
  etacry = (Theta-ThetaCentr)/ThetaWidth;    
  //flip to take into account ECAL barrel symmetries:
  if (ieta<0) etacry *= -1.;  
  
  return;

}

void egammaTools::localEcalClusterCoordsEE	(	const reco::CaloCluster &	bclus,
		const CaloGeometry &	geom,
		float &	xcry,
		float &	ycry,
		int &	ix,
		int &	iy,
		float &	thetatilt,
		float &	phitilt
	)

Definition at line 83 of file EcalClusterLocal.cc.

References Abs(), reco::deltaR(), egammaForCoreTracking_cff::depth, runTauDisplay::dr, DetId::Ecal, EcalEndcap, reco::CaloCluster::energy(), reco::tau::disc::Eta(), PV3DBase< T, PVType, FrameType >::eta(), reco::CaloCluster::eta(), plotBeamSpotDB::first, relativeConstraints::geom, TruncatedPyramid::getPhiAxis(), TruncatedPyramid::getPosition(), CaloGeometry::getSubdetectorGeometry(), TruncatedPyramid::getThetaAxis(), reco::CaloCluster::hitsAndFractions(), EEDetId::ix(), EEDetId::iy(), cmsBatch::log, colinearityKinematic::Phi, PV3DBase< T, PVType, FrameType >::phi(), reco::CaloCluster::position(), X, PV3DBase< T, PVType, FrameType >::x(), X0, DOFs::Y, and PV3DBase< T, PVType, FrameType >::y().

Referenced by EcalRegressionData::fill(), EGRegressionModifierV3::getSeedCrysCoord(), SuperClusterHelper::localCoordinates(), EGRegressionModifierV2::modifyObject(), and EGRegressionModifierV1::modifyObject().

{
    
  assert(bclus.hitsAndFractions().at(0).first.subdetId()==EcalEndcap);
  
  const CaloSubdetectorGeometry* geom=caloGeometry.getSubdetectorGeometry(DetId::Ecal,EcalEndcap);//EcalBarrel = 1
  
  const math::XYZPoint& position_ = bclus.position(); 
  //double Theta = -position_.theta()+0.5*TMath::Pi();
  double Eta = position_.eta();
  double Phi = TVector2::Phi_mpi_pi(position_.phi());
  double X = position_.x();
  double Y = position_.y();
  
  //Calculate expected depth of the maximum shower from energy (like in PositionCalc::Calculate_Location()):
  // The parameters X0 and T0 are hardcoded here because these values were used to calculate the corrections:
  const float X0 = 0.89; float T0 = 1.2;
  //different T0 value if outside of preshower coverage
  if (TMath::Abs(bclus.eta())<1.653) T0 = 3.1;
  
  double depth = X0 * (T0 + log(bclus.energy()));
  
  //find max energy crystal
  std::vector< std::pair<DetId, float> > crystals_vector = bclus.hitsAndFractions();
  float drmin = 999.;
  EEDetId crystalseed;
  //printf("starting loop over crystals, etot = %5f:\n",bclus.energy());
  for (unsigned int icry=0; icry!=crystals_vector.size(); ++icry) {    
    
    EEDetId crystal(crystals_vector[icry].first);
        
    auto cell=geom->getGeometry(crystal);
    const TruncatedPyramid* cpyr = dynamic_cast<const TruncatedPyramid*>(cell.get());
    GlobalPoint center_pos = cpyr->getPosition(depth);
    double EtaCentr = center_pos.eta();
    double PhiCentr = TVector2::Phi_mpi_pi(center_pos.phi());

    float dr = reco::deltaR(Eta,Phi,EtaCentr,PhiCentr);
    if (dr<drmin) {
      drmin = dr;
      crystalseed = crystal;
    }
    
  }
  
  ix = crystalseed.ix();
  iy = crystalseed.iy();
  
  // Get center cell position from shower depth
  auto cell=geom->getGeometry(crystalseed);
  const TruncatedPyramid* cpyr = dynamic_cast<const TruncatedPyramid*>(cell.get());

  thetatilt = cpyr->getThetaAxis();
  phitilt = cpyr->getPhiAxis();

  GlobalPoint center_pos = cpyr->getPosition(depth);
  
  double XCentr = center_pos.x();
  double XWidth = 2.59;
  xcry = (X-XCentr)/XWidth;
 
  
  double YCentr = center_pos.y();
  double YWidth = 2.59;
  ycry = (Y-YCentr)/YWidth;
  
  return;

}