00001 #include "PhysicsTools/IsolationAlgos/interface/IsoDepositVetoFactory.h"
00002
00003 #include "DataFormats/RecoCandidate/interface/IsoDepositVetos.h"
00004 #include "PhysicsTools/IsolationAlgos/interface/EventDependentAbsVetos.h"
00005 #include <boost/regex.hpp>
00006
00007
00008 namespace reco { namespace isodeposit {
00009
00010 class SwitchingEcalVeto : public AbsVeto {
00011 public:
00012
00013 SwitchingEcalVeto(AbsVeto *veto, bool isBarrel) :
00014 veto_(veto), barrel_(isBarrel) {}
00015 virtual bool veto(double eta, double phi, float value) const {
00016 return (fabs(eta) < 1.479) == (barrel_) ? veto_->veto(eta,phi,value) : false;
00017 }
00018 virtual void centerOn(double eta, double phi) {
00019 veto_->centerOn(eta,phi);
00020 }
00021 private:
00022 std::auto_ptr<AbsVeto> veto_;
00023 bool barrel_;
00024 };
00025
00026 class NumCrystalVeto : public AbsVeto {
00027 public:
00028 NumCrystalVeto(Direction dir, double iR) : vetoDir_(dir), iR_(iR) {}
00029 virtual bool veto(double eta, double phi, float value) const {
00030 if( fabs(vetoDir_.eta()) < 1.479) {
00031 return ( vetoDir_.deltaR(Direction(eta,phi)) < 0.0174*iR_ );
00032 } else {
00033 return ( vetoDir_.deltaR(Direction(eta,phi)) < 0.00864*sinh(eta)*iR_ );
00034 }
00035 }
00036 virtual void centerOn(double eta, double phi) { vetoDir_ = Direction(eta,phi); }
00037 private:
00038 Direction vetoDir_; float iR_;
00039 };
00040
00041 class NumCrystalEtaPhiVeto : public AbsVeto {
00042 public:
00043 NumCrystalEtaPhiVeto(math::XYZVectorD dir, double iEta, double iPhi) :
00044 vetoDir_(dir.eta(),dir.phi()),
00045 iEta_(iEta),
00046 iPhi_(iPhi) {}
00047 NumCrystalEtaPhiVeto(Direction dir, double iEta, double iPhi) :
00048 vetoDir_(dir.eta(),dir.phi()),
00049 iEta_(iEta),
00050 iPhi_(iPhi) {}
00051 virtual bool veto(double eta, double phi, float value) const {
00052 double dPhi = phi - vetoDir_.phi();
00053 double dEta = eta - vetoDir_.eta();
00054 while( dPhi < -M_PI ) dPhi += 2*M_PI;
00055 while( dPhi >= M_PI ) dPhi -= 2*M_PI;
00056 if( fabs(vetoDir_.eta()) < 1.479) {
00057 return ( (fabs(dEta) < 0.0174*iEta_) && (fabs(dPhi) < 0.0174*iPhi_) );
00058 } else {
00059 return ( (fabs(dEta) < 0.00864*sinh(eta)*iEta_) &&
00060 (fabs(dPhi) < 0.00864*sinh(eta)*iPhi_) );
00061 }
00062 }
00063 virtual void centerOn(double eta, double phi) { vetoDir_ = Direction(eta,phi); }
00064 private:
00065 Direction vetoDir_;
00066 double iEta_, iPhi_;
00067 };
00068
00069 } }
00070
00071
00072 reco::isodeposit::AbsVeto *
00073 IsoDepositVetoFactory::make(const char *string) {
00074 reco::isodeposit::EventDependentAbsVeto * evdep = 0;
00075 std::auto_ptr<reco::isodeposit::AbsVeto> ret(make(string,evdep));
00076 if (evdep != 0) {
00077 throw cms::Exception("Configuration") << "The resulting AbsVeto depends on the edm::Event.\n"
00078 << "Please use the two-arguments IsoDepositVetoFactory::make.\n";
00079 }
00080 return ret.release();
00081 }
00082
00083 reco::isodeposit::AbsVeto *
00084 IsoDepositVetoFactory::make(const char *string, reco::isodeposit::EventDependentAbsVeto *&evdep) {
00085 using namespace reco::isodeposit;
00086 static boost::regex
00087 ecalSwitch("^Ecal(Barrel|Endcaps):(.*)"),
00088 threshold("Threshold\\((\\d+\\.\\d+)\\)"),
00089 thresholdtransverse("ThresholdFromTransverse\\((\\d+\\.\\d+)\\)"),
00090 absthreshold("AbsThreshold\\((\\d+\\.\\d+)\\)"),
00091 absthresholdtransverse("AbsThresholdFromTransverse\\((\\d+\\.\\d+)\\)"),
00092 cone("ConeVeto\\((\\d+\\.\\d+)\\)"),
00093 angleCone("AngleCone\\((\\d+\\.\\d+)\\)"),
00094 angleVeto("AngleVeto\\((\\d+\\.\\d+)\\)"),
00095 rectangularEtaPhiVeto("RectangularEtaPhiVeto\\(([+-]?\\d+\\.\\d+),([+-]?\\d+\\.\\d+),([+-]?\\d+\\.\\d+),([+-]?\\d+\\.\\d+)\\)"),
00096 numCrystal("NumCrystalVeto\\((\\d+\\.\\d+)\\)"),
00097 numCrystalEtaPhi("NumCrystalEtaPhiVeto\\((\\d+\\.\\d+),(\\d+\\.\\d+)\\)"),
00098 otherCandidates("OtherCandidatesByDR\\((\\w+:?\\w*:?\\w*),\\s*(\\d+\\.?|\\d*\\.\\d*)\\)"),
00099 number("^(\\d+\\.?|\\d*\\.\\d*)$");
00100 boost::cmatch match;
00101
00102 evdep = 0;
00103 if (regex_match(string, match, ecalSwitch)) {
00104 return new SwitchingEcalVeto(make(match[2].first), (match[1] == "Barrel") );
00105 } else if (regex_match(string, match, threshold)) {
00106 return new ThresholdVeto(atof(match[1].first));
00107 } else if (regex_match(string, match, thresholdtransverse)) {
00108 return new ThresholdVetoFromTransverse(atof(((std::string)match[1]).c_str()));
00109 } else if (regex_match(string, match, absthreshold)) {
00110 return new AbsThresholdVeto(atof(match[1].first));
00111 } else if (regex_match(string, match, absthresholdtransverse)) {
00112 return new AbsThresholdVetoFromTransverse(atof(((std::string)match[1]).c_str()));
00113 } else if (regex_match(string, match, cone)) {
00114 return new ConeVeto(Direction(), atof(match[1].first));
00115 } else if (regex_match(string, match, number)) {
00116 return new ConeVeto(Direction(), atof(match[1].first));
00117 } else if (regex_match(string, match, angleCone)) {
00118 return new AngleCone(Direction(), atof(match[1].first));
00119 } else if (regex_match(string, match, angleVeto)) {
00120 return new AngleConeVeto(Direction(), atof(match[1].first));
00121 } else if (regex_match(string, match, rectangularEtaPhiVeto)) {
00122 return new RectangularEtaPhiVeto(Direction(),
00123 atof(match[1].first), atof(match[2].first),
00124 atof(match[3].first), atof(match[4].first));
00125 } else if (regex_match(string, match, numCrystal)) {
00126 return new NumCrystalVeto(Direction(), atof(match[1].first));
00127 } else if (regex_match(string, match, numCrystalEtaPhi)) {
00128 return new NumCrystalEtaPhiVeto(Direction(),atof(match[1].first),atof(match[2].first));
00129 } else if (regex_match(string, match, otherCandidates)) {
00130 OtherCandidatesDeltaRVeto *ret = new OtherCandidatesDeltaRVeto(edm::InputTag(match[1]),
00131 atof(match[2].first));
00132 evdep = ret;
00133 return ret;
00134 } else {
00135 throw cms::Exception("Not Implemented") << "Veto " << string << " not implemented yet...";
00136 }
00137 }