Go to the documentation of this file.00001
00002
00003
00004
00005
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022 #include "SUSYBSMAnalysis/HSCP/interface/HSCParticleProducer.h"
00023
00024 HSCParticleProducer::HSCParticleProducer(const edm::ParameterSet& iConfig) {
00025 using namespace edm;
00026 using namespace std;
00027
00028
00029 Filter_ = iConfig.getParameter<bool> ("filter");
00030
00031
00032 m_trackTag = iConfig.getParameter<edm::InputTag>("tracks");
00033 m_muonsTag = iConfig.getParameter<edm::InputTag>("muons");
00034 m_trackIsoTag = iConfig.getParameter<edm::InputTag>("tracksIsolation");
00035
00036 useBetaFromTk = iConfig.getParameter<bool> ("useBetaFromTk" );
00037 useBetaFromMuon = iConfig.getParameter<bool> ("useBetaFromMuon");
00038 useBetaFromRpc = iConfig.getParameter<bool> ("useBetaFromRpc" );
00039 useBetaFromEcal = iConfig.getParameter<bool> ("useBetaFromEcal");
00040
00041
00042 minTkP = iConfig.getParameter<double> ("minTkP");
00043 maxTkChi2 = iConfig.getParameter<double> ("maxTkChi2");
00044 minTkHits = iConfig.getParameter<uint32_t>("minTkHits");
00045 minMuP = iConfig.getParameter<double> ("minMuP");
00046 minDR = iConfig.getParameter<double> ("minDR");
00047 maxInvPtDiff = iConfig.getParameter<double> ("maxInvPtDiff");
00048
00049 if(useBetaFromTk )beta_calculator_TK = new BetaCalculatorTK (iConfig);
00050 if(useBetaFromMuon)beta_calculator_MUON = new BetaCalculatorMUON(iConfig);
00051 if(useBetaFromRpc )beta_calculator_RPC = new BetaCalculatorRPC (iConfig);
00052 if(useBetaFromEcal)beta_calculator_ECAL = new BetaCalculatorECAL(iConfig);
00053
00054
00055 std::vector<edm::ParameterSet> SelectionParameters = iConfig.getParameter<std::vector<edm::ParameterSet> >("SelectionParameters");
00056 for(unsigned int i=0;i<SelectionParameters.size();i++){
00057 Selectors.push_back(new CandidateSelector(SelectionParameters[i]) );
00058 }
00059
00060
00061 produces<susybsm::HSCParticleCollection >();
00062 if(useBetaFromEcal)produces<susybsm::HSCPCaloInfoCollection >();
00063
00064 }
00065
00066 HSCParticleProducer::~HSCParticleProducer() {
00067
00068
00069 }
00070
00071
00072
00073
00074
00075
00076 bool
00077 HSCParticleProducer::filter(edm::Event& iEvent, const edm::EventSetup& iSetup) {
00078
00079 using namespace edm;
00080 using namespace reco;
00081 using namespace std;
00082 using namespace susybsm;
00083
00084
00085 edm::Handle<reco::MuonCollection> muonCollectionHandle;
00086 iEvent.getByLabel(m_muonsTag,muonCollectionHandle);
00087
00088
00089 edm::Handle<reco::TrackCollection> trackCollectionHandle;
00090 iEvent.getByLabel(m_trackTag,trackCollectionHandle);
00091
00092
00093 edm::Handle<reco::TrackCollection> trackIsoCollectionHandle;
00094 iEvent.getByLabel(m_trackIsoTag,trackIsoCollectionHandle);
00095
00096
00097
00098 susybsm::HSCParticleCollection* hscp = new susybsm::HSCParticleCollection;
00099 std::auto_ptr<susybsm::HSCParticleCollection> result(hscp);
00100
00101 susybsm::HSCPCaloInfoCollection* caloInfoColl = new susybsm::HSCPCaloInfoCollection;
00102 std::auto_ptr<susybsm::HSCPCaloInfoCollection> caloInfoCollaptr(caloInfoColl);
00103
00104
00105
00106 *hscp = getHSCPSeedCollection(trackCollectionHandle, muonCollectionHandle);
00107
00108
00109 for(susybsm::HSCParticleCollection::iterator hscpcandidate = hscp->begin(); hscpcandidate != hscp->end(); ++hscpcandidate) {
00110
00111 reco::TrackRef track = hscpcandidate->trackRef();
00112 if(track.isNull())continue;
00113 float dRMin=1000; int found = -1;
00114 for(unsigned int t=0; t<trackIsoCollectionHandle->size();t++) {
00115 reco::TrackRef Isotrack = reco::TrackRef( trackIsoCollectionHandle, t );
00116 if( fabs( (1.0/track->pt())-(1.0/Isotrack->pt())) > maxInvPtDiff) continue;
00117 float dR = deltaR(track->momentum(), Isotrack->momentum());
00118 if(dR <= minDR && dR < dRMin){ dRMin=dR; found = t;}
00119 }
00120 if(found>=0)hscpcandidate->setTrackIso(reco::TrackRef( trackIsoCollectionHandle, found ));
00121 }
00122
00123
00124 if(useBetaFromTk){
00125 for(susybsm::HSCParticleCollection::iterator hscpcandidate = hscp->begin(); hscpcandidate != hscp->end(); ++hscpcandidate) {
00126 beta_calculator_TK->addInfoToCandidate(*hscpcandidate, iEvent,iSetup);
00127 }}
00128
00129
00130 if(useBetaFromMuon){
00131 for(susybsm::HSCParticleCollection::iterator hscpcandidate = hscp->begin(); hscpcandidate != hscp->end(); ++hscpcandidate) {
00132 beta_calculator_MUON->addInfoToCandidate(*hscpcandidate, iEvent,iSetup);
00133 }}
00134
00135
00136 if(useBetaFromRpc){
00137 for(susybsm::HSCParticleCollection::iterator hscpcandidate = hscp->begin(); hscpcandidate != hscp->end(); ++hscpcandidate) {
00138 beta_calculator_RPC->addInfoToCandidate(*hscpcandidate, iEvent, iSetup);
00139 }}
00140
00141
00142
00143
00144
00145 if(useBetaFromEcal){
00146 int Index=0;
00147 caloInfoColl->resize(hscp->size());
00148 for(susybsm::HSCParticleCollection::iterator hscpcandidate = hscp->begin(); hscpcandidate != hscp->end(); ++hscpcandidate, Index++) {
00149 beta_calculator_ECAL->addInfoToCandidate(*hscpcandidate,trackCollectionHandle,iEvent,iSetup, (*caloInfoColl)[Index]);
00150 }}
00151
00152
00153 for(int i=0;i<(int)hscp->size();i++) {
00154 susybsm::HSCParticleCollection::iterator hscpcandidate = hscp->begin() + i;
00155 bool decision = false;
00156 for(unsigned int s=0;s<Selectors.size();s++){decision |= Selectors[s]->isSelected(*hscpcandidate);}
00157 if(!decision){
00158 hscp->erase(hscpcandidate);
00159 if(useBetaFromEcal)caloInfoColl->erase(caloInfoColl->begin() + i);
00160 i--;
00161 }
00162 }
00163 bool filterResult = !Filter_ || (Filter_ && hscp->size()>=1);
00164
00165
00166
00167
00168
00169 if(useBetaFromEcal){
00170 edm::OrphanHandle<susybsm::HSCPCaloInfoCollection> caloInfoHandle= iEvent.put(caloInfoCollaptr);
00171
00172 for(int i=0;i<(int)hscp->size();i++) {
00173 susybsm::HSCParticleCollection::iterator hscpcandidate = hscp->begin() + i;
00174 hscpcandidate->setCaloInfo(HSCPCaloInfoRef(caloInfoHandle,i));
00175 }
00176 }
00177
00178
00179
00180
00181
00182 edm::OrphanHandle<susybsm::HSCParticleCollection> putHandle = iEvent.put(result);
00183
00184
00185
00186
00187
00188
00189
00190 return filterResult;
00191 }
00192
00193
00194 void
00195 HSCParticleProducer::beginJob() {
00196 }
00197
00198
00199 void
00200 HSCParticleProducer::endJob() {
00201 }
00202
00203 std::vector<HSCParticle> HSCParticleProducer::getHSCPSeedCollection(edm::Handle<reco::TrackCollection>& trackCollectionHandle, edm::Handle<reco::MuonCollection>& muonCollectionHandle)
00204 {
00205 std::vector<HSCParticle> HSCPCollection;
00206
00207
00208 std::vector<reco::TrackRef> tracks;
00209 for(unsigned int i=0; i<trackCollectionHandle->size(); i++){
00210 TrackRef track = reco::TrackRef( trackCollectionHandle, i );
00211 if(track->p()<minTkP || (track->chi2()/track->ndof())>maxTkChi2 || track->found()<minTkHits)continue;
00212 tracks.push_back( track );
00213 }
00214
00215
00216 for(unsigned int m=0; m<muonCollectionHandle->size(); m++){
00217 reco::MuonRef muon = reco::MuonRef( muonCollectionHandle, m );
00218 if(muon->p()<minMuP )continue;
00219 TrackRef innertrack = muon->innerTrack();
00220 if(innertrack.isNull())continue;
00221
00222
00223
00224 float dRMin=1000; int found = -1;
00225 for(unsigned int t=0; t<tracks.size();t++) {
00226 reco::TrackRef track = tracks[t];
00227 if( fabs( (1.0/innertrack->pt())-(1.0/track->pt())) > maxInvPtDiff) continue;
00228 float dR = deltaR(innertrack->momentum(), track->momentum());
00229 if(dR <= minDR && dR < dRMin){ dRMin=dR; found = t;}
00230 }
00231
00232 HSCParticle candidate;
00233 candidate.setMuon(muon);
00234 if(found>=0){
00235
00236 candidate.setTrack(tracks[found]);
00237 tracks.erase(tracks.begin()+found);
00238 }
00239 HSCPCollection.push_back(candidate);
00240 }
00241
00242
00243 for(unsigned int m=0; m<muonCollectionHandle->size(); m++){
00244 reco::MuonRef muon = reco::MuonRef( muonCollectionHandle, m );
00245 if(muon->p()<minMuP)continue;
00246 TrackRef innertrack = muon->innerTrack();
00247 if(innertrack.isNonnull())continue;
00248
00249
00250 float dRMin=1000; int found = -1;
00251 for(unsigned int t=0; t<tracks.size();t++) {
00252 reco::TrackRef track = tracks[t];
00253 if( fabs( (1.0/muon->pt())-(1.0/track->pt())) > maxInvPtDiff) continue;
00254 float dR = deltaR(muon->momentum(), track->momentum());
00255 if(dR <= minDR && dR < dRMin){ dRMin=dR; found = t;}
00256 }
00257
00258 HSCParticle candidate;
00259 candidate.setMuon(muon);
00260 if(found>=0){
00261
00262 candidate.setTrack(tracks[found]);
00263 tracks.erase(tracks.begin()+found);
00264 }
00265 HSCPCollection.push_back(candidate);
00266 }
00267
00268
00269 for(unsigned int i=0; i<tracks.size(); i++){
00270 HSCParticle candidate;
00271 candidate.setTrack(tracks[i]);
00272 HSCPCollection.push_back(candidate);
00273 }
00274
00275 return HSCPCollection;
00276 }
00277
00278
00279 DEFINE_FWK_MODULE(HSCParticleProducer);
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298