Go to the documentation of this file.00001
00002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00003
00004 #include "FWCore/Framework/interface/Event.h"
00005 #include "FWCore/Utilities/interface/EDMException.h"
00006
00007
00008 #include <DataFormats/TrackReco/interface/Track.h>
00009 #include <DataFormats/METReco/interface/CaloMET.h>
00010 #include <DataFormats/METReco/interface/CaloMETFwd.h>
00011 #include <DataFormats/Math/interface/deltaPhi.h>
00012
00013
00014 #include <math.h>
00015
00016 #include "TLorentzVector.h"
00017
00018 #include "Alignment/CommonAlignmentProducer/interface/AlignmentTwoBodyDecayTrackSelector.h"
00019
00020 using namespace std;
00021 using namespace edm;
00022
00023
00024 AlignmentTwoBodyDecayTrackSelector::AlignmentTwoBodyDecayTrackSelector(const edm::ParameterSet & cfg) :
00025 theMissingETSource("met")
00026 {
00027 LogDebug("Alignment") << "> applying two body decay Trackfilter ...";
00028 theMassrangeSwitch = cfg.getParameter<bool>( "applyMassrangeFilter" );
00029 if (theMassrangeSwitch){
00030 theMinMass = cfg.getParameter<double>( "minXMass" );
00031 theMaxMass = cfg.getParameter<double>( "maxXMass" );
00032 theDaughterMass = cfg.getParameter<double>( "daughterMass" );
00033 LogDebug("Alignment") << "> Massrange min,max : " << theMinMass << "," << theMaxMass
00034 << "\n> Mass of daughter Particle : " << theDaughterMass;
00035
00036 }else{
00037 theMinMass = 0;
00038 theMaxMass = 0;
00039 theDaughterMass = 0;
00040 }
00041 theChargeSwitch = cfg.getParameter<bool>( "applyChargeFilter" );
00042 if(theChargeSwitch){
00043 theCharge = cfg.getParameter<int>( "charge" );
00044 theUnsignedSwitch = cfg.getParameter<bool>( "useUnsignedCharge" );
00045 if(theUnsignedSwitch)
00046 theCharge=std::abs(theCharge);
00047 LogDebug("Alignment") << "> Desired Charge, unsigned: "<<theCharge<<" , "<<theUnsignedSwitch;
00048 }else{
00049 theCharge =0;
00050 theUnsignedSwitch = true;
00051 }
00052 theMissingETSwitch = cfg.getParameter<bool>( "applyMissingETFilter" );
00053 if(theMissingETSwitch){
00054 theMissingETSource = cfg.getParameter<InputTag>( "missingETSource" );
00055 LogDebug("Alignment") << "> missing Et Source: "<< theMissingETSource;
00056 }
00057 theAcoplanarityFilterSwitch = cfg.getParameter<bool>( "applyAcoplanarityFilter" );
00058 if(theAcoplanarityFilterSwitch){
00059 theAcoplanarDistance = cfg.getParameter<double>( "acoplanarDistance" );
00060 LogDebug("Alignment") << "> Acoplanar Distance: "<<theAcoplanarDistance;
00061 }
00062
00063 }
00064
00065
00066
00067 AlignmentTwoBodyDecayTrackSelector::~AlignmentTwoBodyDecayTrackSelector()
00068 {}
00069
00070
00072 bool AlignmentTwoBodyDecayTrackSelector::useThisFilter()
00073 {
00074 return theMassrangeSwitch || theChargeSwitch || theAcoplanarityFilterSwitch;
00075 }
00076
00077
00078
00079 AlignmentTwoBodyDecayTrackSelector::Tracks
00080 AlignmentTwoBodyDecayTrackSelector::select(const Tracks& tracks, const edm::Event& iEvent)
00081 {
00082 Tracks result=tracks;
00083
00084 if(theMassrangeSwitch){
00085 if(theMissingETSwitch)
00086 result = checkMETMass(result,iEvent);
00087 else
00088 result = checkMass(result);
00089 }
00090 if(theChargeSwitch)
00091 result = checkCharge(result);
00092 if(theAcoplanarityFilterSwitch){
00093 if(theMissingETSwitch)
00094 result = checkMETAcoplanarity(result,iEvent);
00095 else
00096 result = checkAcoplanarity(result);
00097 }
00098 LogDebug("Alignment") << "> TwoBodyDecay tracks all,kept: " << tracks.size() << "," << result.size();
00099
00100
00101 return result;
00102
00103 }
00104
00106 AlignmentTwoBodyDecayTrackSelector::Tracks
00107 AlignmentTwoBodyDecayTrackSelector::checkMass(const Tracks& cands) const
00108 {
00109 Tracks result; result.clear();
00110
00111 if(cands.size() == 2){
00112
00113 TLorentzVector track0(cands.at(0)->px(),cands.at(0)->py(),cands.at(0)->pz(),
00114 sqrt((cands.at(0)->p()*cands.at(0)->p())+theDaughterMass*theDaughterMass));
00115 TLorentzVector track1(cands.at(1)->px(),cands.at(1)->py(),cands.at(1)->pz(),
00116 sqrt((cands.at(1)->p()*cands.at(1)->p())+theDaughterMass*theDaughterMass));
00117 TLorentzVector mother = track0+track1;
00118 if(mother.M() > theMinMass && mother.M() < theMaxMass)
00119 result = cands;
00120 LogDebug("Alignment") <<"> mass of mother: "<<mother.M()<<"GeV";
00121 }
00122 return result;
00123 }
00124
00126 AlignmentTwoBodyDecayTrackSelector::Tracks
00127 AlignmentTwoBodyDecayTrackSelector::checkMETMass(const Tracks& cands,const edm::Event& iEvent) const
00128 {
00129 Tracks result; result.clear();
00130 if(cands.size() == 1){
00131 Handle<reco::CaloMETCollection> missingET;
00132 iEvent.getByLabel(theMissingETSource ,missingET);
00133 if(missingET.isValid()){
00134
00135
00136
00137
00138 TLorentzVector track(cands.at(0)->px(),cands.at(0)->py(),cands.at(0)->pz(),
00139 sqrt((cands.at(0)->p()*cands.at(0)->p())+theDaughterMass*theDaughterMass));
00140 TLorentzVector met((*missingET).at(0).px(),(*missingET).at(0).py(),(*missingET).at(0).pz(),
00141 (*missingET).at(0).p());
00142 TLorentzVector motherSystem = track + met;
00143 if(motherSystem.M() > theMinMass && motherSystem.M() < theMaxMass)
00144 result = cands;
00145 LogDebug("Alignment") <<"> mass of motherSystem: "<<motherSystem.M()<<"GeV";
00146 }else
00147 LogError("Alignment")<<"@SUB=AlignmentTwoBodyDecayTrackSelector::checkMETMass"
00148 <<"> could not optain missingET Collection!";
00149 }
00150 return cands;
00151 }
00152
00154 AlignmentTwoBodyDecayTrackSelector::Tracks
00155 AlignmentTwoBodyDecayTrackSelector::checkCharge(const Tracks& cands) const
00156 {
00157 Tracks result; result.clear();
00158 int sumCharge = 0;
00159 for(Tracks::const_iterator it = cands.begin();it < cands.end();++it)
00160 sumCharge += (*it)->charge();
00161 if(theUnsignedSwitch)
00162 sumCharge = std::abs(sumCharge);
00163 if(sumCharge == theCharge)
00164 result = cands;
00165
00166 return result;
00167 }
00168
00170 AlignmentTwoBodyDecayTrackSelector::Tracks
00171 AlignmentTwoBodyDecayTrackSelector::checkAcoplanarity(const Tracks& cands) const
00172 {
00173 Tracks result; result.clear();
00174
00175 if(cands.size() == 2){
00176 LogDebug("Alignment") <<"> Acoplanarity: "<<fabs(fabs(deltaPhi(cands.at(0)->phi(),cands.at(1)->phi()))-M_PI)<<endl;
00177 if(fabs(fabs(deltaPhi(cands.at(0)->phi(),cands.at(1)->phi()))-M_PI)<theAcoplanarDistance)
00178 result = cands;
00179 }
00180 return result;
00181 }
00183 AlignmentTwoBodyDecayTrackSelector::Tracks
00184 AlignmentTwoBodyDecayTrackSelector::checkMETAcoplanarity(const Tracks& cands,const edm::Event& iEvent)const
00185 {
00186 Tracks result; result.clear();
00187 if(cands.size() == 1){
00188 Handle<reco::CaloMETCollection> missingET;
00189 iEvent.getByLabel(theMissingETSource ,missingET);
00190 if(missingET.isValid()){
00191
00192 LogDebug("Alignment") <<"> METAcoplanarity: "<<fabs(fabs(deltaPhi(cands.at(0)->phi(),(*missingET).at(0).phi()))-M_PI)<<endl;
00193 if(fabs(fabs(deltaPhi(cands.at(0)->phi(),(*missingET).at(0).phi()))-M_PI)<theAcoplanarDistance)
00194 result = cands;
00195
00196 }else
00197 LogError("Alignment")<<"@SUB=AlignmentTwoBodyDecayTrackSelector::checkMETAcoplanarity"
00198 <<"> could not optain missingET Collection!";
00199 }
00200 return result;
00201 }
00202
00203
00204
00206 void AlignmentTwoBodyDecayTrackSelector::printTracks(const Tracks& col) const
00207 {
00208 int count = 0;
00209 LogDebug("Alignment") << ">......................................";
00210 for(Tracks::const_iterator it = col.begin();it < col.end();++it,++count){
00211 LogDebug("Alignment")
00212 <<"> Track No. "<< count <<": p = ("<<(*it)->px()<<","<<(*it)->py()<<","<<(*it)->pz()<<")\n"
00213 <<"> pT = "<<(*it)->pt()<<" eta = "<<(*it)->eta()<<" charge = "<<(*it)->charge();
00214 }
00215 LogDebug("Alignment") << ">......................................";
00216 }
00217
00218