CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
CaloRecoTauAlgorithm.cc
Go to the documentation of this file.
3 
5 
6 using namespace reco;
7 
8 CaloRecoTauAlgorithm::CaloRecoTauAlgorithm() : TransientTrackBuilder_(0),MagneticField_(0),chargedpi_mass_(0.13957018){}
9 CaloRecoTauAlgorithm::CaloRecoTauAlgorithm(const edm::ParameterSet& iConfig) : TransientTrackBuilder_(0),MagneticField_(0),chargedpi_mass_(0.13957018){
10  LeadTrack_minPt_ = iConfig.getParameter<double>("LeadTrack_minPt");
11  Track_minPt_ = iConfig.getParameter<double>("Track_minPt");
12  IsolationTrack_minPt_ = iConfig.getParameter<double>("IsolationTrack_minPt");
13  IsolationTrack_minHits_ = iConfig.getParameter<unsigned int>("IsolationTrack_minHits");
14  UseTrackLeadTrackDZconstraint_ = iConfig.getParameter<bool>("UseTrackLeadTrackDZconstraint");
15  TrackLeadTrack_maxDZ_ = iConfig.getParameter<double>("TrackLeadTrack_maxDZ");
16  ECALRecHit_minEt_ = iConfig.getParameter<double>("ECALRecHit_minEt");
17 
18  MatchingConeMetric_ = iConfig.getParameter<std::string>("MatchingConeMetric");
19  MatchingConeSizeFormula_ = iConfig.getParameter<std::string>("MatchingConeSizeFormula");
20  MatchingConeSize_min_ = iConfig.getParameter<double>("MatchingConeSize_min");
21  MatchingConeSize_max_ = iConfig.getParameter<double>("MatchingConeSize_max");
22  TrackerSignalConeMetric_ = iConfig.getParameter<std::string>("TrackerSignalConeMetric");
23  TrackerSignalConeSizeFormula_ = iConfig.getParameter<std::string>("TrackerSignalConeSizeFormula");
24  TrackerSignalConeSize_min_ = iConfig.getParameter<double>("TrackerSignalConeSize_min");
25  TrackerSignalConeSize_max_ = iConfig.getParameter<double>("TrackerSignalConeSize_max");
26  TrackerIsolConeMetric_ = iConfig.getParameter<std::string>("TrackerIsolConeMetric");
27  TrackerIsolConeSizeFormula_ = iConfig.getParameter<std::string>("TrackerIsolConeSizeFormula");
28  TrackerIsolConeSize_min_ = iConfig.getParameter<double>("TrackerIsolConeSize_min");
29  TrackerIsolConeSize_max_ = iConfig.getParameter<double>("TrackerIsolConeSize_max");
30  ECALSignalConeMetric_ = iConfig.getParameter<std::string>("ECALSignalConeMetric");
31  ECALSignalConeSizeFormula_ = iConfig.getParameter<std::string>("ECALSignalConeSizeFormula");
32  ECALSignalConeSize_min_ = iConfig.getParameter<double>("ECALSignalConeSize_min");
33  ECALSignalConeSize_max_ = iConfig.getParameter<double>("ECALSignalConeSize_max");
34  ECALIsolConeMetric_ = iConfig.getParameter<std::string>("ECALIsolConeMetric");
35  ECALIsolConeSizeFormula_ = iConfig.getParameter<std::string>("ECALIsolConeSizeFormula");
36  ECALIsolConeSize_min_ = iConfig.getParameter<double>("ECALIsolConeSize_min");
37  ECALIsolConeSize_max_ = iConfig.getParameter<double>("ECALIsolConeSize_max");
38 
39  EBRecHitsLabel_ = iConfig.getParameter<edm::InputTag>("EBRecHitsSource");
40  EERecHitsLabel_ = iConfig.getParameter<edm::InputTag>("EERecHitsSource");
41  ESRecHitsLabel_ = iConfig.getParameter<edm::InputTag>("ESRecHitsSource");
42 
43 
44 
45  AreaMetric_recoElements_maxabsEta_ = iConfig.getParameter<double>("AreaMetric_recoElements_maxabsEta");
46 
47  //Computing the TFormula
48  myTrackerSignalConeSizeTFormula=TauTagTools::computeConeSizeTFormula(TrackerSignalConeSizeFormula_,"Tracker signal cone size");
49  myTrackerIsolConeSizeTFormula=TauTagTools::computeConeSizeTFormula(TrackerIsolConeSizeFormula_,"Tracker isolation cone size");
50  myECALSignalConeSizeTFormula=TauTagTools::computeConeSizeTFormula(ECALSignalConeSizeFormula_,"ECAL signal cone size");
51  myECALIsolConeSizeTFormula=TauTagTools::computeConeSizeTFormula(ECALIsolConeSizeFormula_,"ECAL isolation cone size");
52  myMatchingConeSizeTFormula=TauTagTools::computeConeSizeTFormula(MatchingConeSizeFormula_,"Matching cone size");
53 
54  mySelectedDetId_.clear();
55 }
58 
60  CaloJetRef myCaloJet=(*myCaloTauTagInfoRef).calojetRef(); // catch a ref to the initial CaloJet
61  const std::vector<CaloTowerPtr> myCaloTowers=(*myCaloJet).getCaloConstituents();
62  JetBaseRef jetRef = myCaloTauTagInfoRef->jetRef();
63  CaloTau myCaloTau(std::numeric_limits<int>::quiet_NaN(),jetRef->p4()); // create the CaloTau with the corrected Lorentz-vector
64 
65  myCaloTau.setcaloTauTagInfoRef(myCaloTauTagInfoRef);
66 
67  TrackRefVector myTks=(*myCaloTauTagInfoRef).Tracks();
68 
69  CaloTauElementsOperators myCaloTauElementsOperators(myCaloTau);
70  double myMatchingConeSize=myCaloTauElementsOperators.computeConeSize(myMatchingConeSizeTFormula,MatchingConeSize_min_,MatchingConeSize_max_);
71  TrackRef myleadTk=myCaloTauElementsOperators.leadTk(MatchingConeMetric_,myMatchingConeSize,LeadTrack_minPt_);
72  double myCaloTau_refInnerPosition_x=0.;
73  double myCaloTau_refInnerPosition_y=0.;
74  double myCaloTau_refInnerPosition_z=0.;
75  if(myleadTk.isNonnull()){
76  myCaloTau.setleadTrack(myleadTk);
77  double myleadTkDZ=(*myleadTk).dz(myPV.position());
79  {
80  const TransientTrack myleadTransientTk=TransientTrackBuilder_->build(&(*myleadTk));
81  GlobalVector myCaloJetdir((*myCaloJet).px(),(*myCaloJet).py(),(*myCaloJet).pz());
82  if(IPTools::signedTransverseImpactParameter(myleadTransientTk,myCaloJetdir,myPV).first)
83  myCaloTau.setleadTracksignedSipt(IPTools::signedTransverseImpactParameter(myleadTransientTk,myCaloJetdir,myPV).second.significance());
84  }
85  if((*myleadTk).innerOk()){
86  myCaloTau_refInnerPosition_x=(*myleadTk).innerPosition().x();
87  myCaloTau_refInnerPosition_y=(*myleadTk).innerPosition().y();
88  myCaloTau_refInnerPosition_z=(*myleadTk).innerPosition().z();
89  }
90 
91  if(MagneticField_!=0){
92  math::XYZPoint mypropagleadTrackECALSurfContactPoint=TauTagTools::propagTrackECALSurfContactPoint(MagneticField_,myleadTk);
93  if(mypropagleadTrackECALSurfContactPoint.R()!=0.){
94  double myleadTrackHCAL3x3hottesthitDEta=0.;
95  double myleadTrackHCAL3x3hottesthitEt=0.;
96  double myleadTrackHCAL3x3hitsEtSum=0.;
97  edm::ESHandle<CaloGeometry> myCaloGeometry;
98  iSetup.get<CaloGeometryRecord>().get(myCaloGeometry);
99  const CaloSubdetectorGeometry* myCaloSubdetectorGeometry=(*myCaloGeometry).getSubdetectorGeometry(DetId::Calo,CaloTowerDetId::SubdetId);
100  CaloTowerDetId mypropagleadTrack_closestCaloTowerId((*myCaloSubdetectorGeometry).getClosestCell(GlobalPoint(mypropagleadTrackECALSurfContactPoint.x(),
101  mypropagleadTrackECALSurfContactPoint.y(),
102  mypropagleadTrackECALSurfContactPoint.z())));
103  std::vector<CaloTowerDetId> mypropagleadTrack_closestCaloTowerNeighbourIds=getCaloTowerneighbourDetIds(myCaloSubdetectorGeometry,mypropagleadTrack_closestCaloTowerId);
104  for(std::vector<CaloTowerPtr>::const_iterator iCaloTower=myCaloTowers.begin();iCaloTower!=myCaloTowers.end();iCaloTower++){
105  CaloTowerDetId iCaloTowerId((**iCaloTower).id());
106  bool CaloTower_inside3x3matrix=false;
107  if (iCaloTowerId==mypropagleadTrack_closestCaloTowerId) CaloTower_inside3x3matrix=true;
108  if (!CaloTower_inside3x3matrix){
109  for(std::vector<CaloTowerDetId>::const_iterator iCaloTowerDetId=mypropagleadTrack_closestCaloTowerNeighbourIds.begin();iCaloTowerDetId!=mypropagleadTrack_closestCaloTowerNeighbourIds.end();iCaloTowerDetId++){
110  if (iCaloTowerId==(*iCaloTowerDetId)){
111  CaloTower_inside3x3matrix=true;
112  break;
113  }
114  }
115  }
116  if (!CaloTower_inside3x3matrix) continue;
117  myleadTrackHCAL3x3hitsEtSum+=(**iCaloTower).hadEt();
118  if((**iCaloTower).hadEt()>=myleadTrackHCAL3x3hottesthitEt ){
119  if ((**iCaloTower).hadEt()!=myleadTrackHCAL3x3hottesthitEt ||
120  ((**iCaloTower).hadEt()==myleadTrackHCAL3x3hottesthitEt && fabs((**iCaloTower).eta()-mypropagleadTrackECALSurfContactPoint.Eta())<myleadTrackHCAL3x3hottesthitDEta)) myleadTrackHCAL3x3hottesthitDEta = fabs((**iCaloTower).eta()-mypropagleadTrackECALSurfContactPoint.Eta());
121  myleadTrackHCAL3x3hottesthitEt=(**iCaloTower).hadEt();
122  }
123  }
124  myCaloTau.setleadTrackHCAL3x3hitsEtSum(myleadTrackHCAL3x3hitsEtSum);
125  if (myleadTrackHCAL3x3hottesthitEt!=0.) myCaloTau.setleadTrackHCAL3x3hottesthitDEta(myleadTrackHCAL3x3hottesthitDEta);
126  }
127  }
128 
130  TrackRefVector myTksbis;
131  for (TrackRefVector::const_iterator iTrack=myTks.begin();iTrack!=myTks.end();++iTrack) {
132  if (fabs((**iTrack).dz(myPV.position())-myleadTkDZ)<=TrackLeadTrack_maxDZ_) myTksbis.push_back(*iTrack);
133  }
134  myTks=myTksbis;
135  }
136 
137 
138  double myTrackerSignalConeSize=myCaloTauElementsOperators.computeConeSize(myTrackerSignalConeSizeTFormula,TrackerSignalConeSize_min_,TrackerSignalConeSize_max_);
139  double myTrackerIsolConeSize=myCaloTauElementsOperators.computeConeSize(myTrackerIsolConeSizeTFormula,TrackerIsolConeSize_min_,TrackerIsolConeSize_max_);
140  double myECALSignalConeSize=myCaloTauElementsOperators.computeConeSize(myECALSignalConeSizeTFormula,ECALSignalConeSize_min_,ECALSignalConeSize_max_);
141  double myECALIsolConeSize=myCaloTauElementsOperators.computeConeSize(myECALIsolConeSizeTFormula,ECALIsolConeSize_min_,ECALIsolConeSize_max_);
142 
143  TrackRefVector mySignalTks;
144  if (UseTrackLeadTrackDZconstraint_) mySignalTks=myCaloTauElementsOperators.tracksInCone((*myleadTk).momentum(),TrackerSignalConeMetric_,myTrackerSignalConeSize,Track_minPt_,TrackLeadTrack_maxDZ_,myleadTkDZ, myPV);
145  else mySignalTks=myCaloTauElementsOperators.tracksInCone((*myleadTk).momentum(),TrackerSignalConeMetric_,myTrackerSignalConeSize,Track_minPt_);
146  myCaloTau.setsignalTracks(mySignalTks);
147 
148  // setting invariant mass of the signal Tracks system
149  math::XYZTLorentzVector mySignalTksInvariantMass(0.,0.,0.,0.);
150  if((int)(mySignalTks.size())!=0){
151  int mySignalTks_qsum=0;
152  for(int i=0;i<(int)mySignalTks.size();i++){
153  mySignalTks_qsum+=mySignalTks[i]->charge();
154  math::XYZTLorentzVector mychargedpicand_fromTk_LorentzVect(mySignalTks[i]->momentum().x(),mySignalTks[i]->momentum().y(),mySignalTks[i]->momentum().z(),sqrt(std::pow((double)mySignalTks[i]->momentum().r(),2)+std::pow(chargedpi_mass_,2)));
155  mySignalTksInvariantMass+=mychargedpicand_fromTk_LorentzVect;
156  }
157  myCaloTau.setCharge(mySignalTks_qsum);
158  }
159  myCaloTau.setsignalTracksInvariantMass(mySignalTksInvariantMass.mass());
160 
161  TrackRefVector myIsolTks;
162  if (UseTrackLeadTrackDZconstraint_) myIsolTks=myCaloTauElementsOperators.tracksInAnnulus((*myleadTk).momentum(),TrackerSignalConeMetric_,myTrackerSignalConeSize,TrackerIsolConeMetric_,myTrackerIsolConeSize,IsolationTrack_minPt_,TrackLeadTrack_maxDZ_,myleadTkDZ, myPV);
163  else myIsolTks=myCaloTauElementsOperators.tracksInAnnulus((*myleadTk).momentum(),TrackerSignalConeMetric_,myTrackerSignalConeSize,TrackerIsolConeMetric_,myTrackerIsolConeSize,IsolationTrack_minPt_);
165  myCaloTau.setisolationTracks(myIsolTks);
166 
167  // setting sum of Pt of the isolation annulus Tracks
168  float myIsolTks_Ptsum=0.;
169  for(int i=0;i<(int)myIsolTks.size();i++) myIsolTks_Ptsum+=myIsolTks[i]->pt();
170  myCaloTau.setisolationTracksPtSum(myIsolTks_Ptsum);
171 
172 
173  //getting the EcalRecHits. Just take them all
174  std::vector<std::pair<math::XYZPoint,float> > thePositionAndEnergyEcalRecHits;
175  mySelectedDetId_.clear();
176  // std::vector<CaloTowerPtr> theCaloTowers=myCaloJet->getCaloConstituents();
177  edm::ESHandle<CaloGeometry> theCaloGeometry;
178  iSetup.get<CaloGeometryRecord>().get(theCaloGeometry);
179  const CaloSubdetectorGeometry* theCaloSubdetectorGeometry;
183  iEvent.getByLabel(EBRecHitsLabel_,EBRecHits);
184  iEvent.getByLabel(EERecHitsLabel_,EERecHits);
185  iEvent.getByLabel(ESRecHitsLabel_,ESRecHits);
186  double maxDeltaR = 0.8;
187  math::XYZPoint myCaloJetdir((*myCaloJet).px(),(*myCaloJet).py(),(*myCaloJet).pz());
188 
189  for(EBRecHitCollection::const_iterator theRecHit = EBRecHits->begin();theRecHit != EBRecHits->end(); theRecHit++){
190  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalBarrel);
191  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(theRecHit->id());
192  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
193  if(ROOT::Math::VectorUtil::DeltaR(myCaloJetdir,theRecHitCell_XYZPoint) < maxDeltaR){
194  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
195  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
196  mySelectedDetId_.push_back(theRecHit->id());
197  }
198  }
199 
200 for(EERecHitCollection::const_iterator theRecHit = EERecHits->begin();theRecHit != EERecHits->end(); theRecHit++){
201  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalEndcap);
202  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(theRecHit->id());
203  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
204  if(ROOT::Math::VectorUtil::DeltaR(myCaloJetdir,theRecHitCell_XYZPoint) < maxDeltaR){
205  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
206  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
207  mySelectedDetId_.push_back(theRecHit->id());
208  }
209 }
210  for(ESRecHitCollection::const_iterator theRecHit = ESRecHits->begin();theRecHit != ESRecHits->end(); theRecHit++){
211  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalPreshower);
212  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(theRecHit->id());
213  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
214  if(ROOT::Math::VectorUtil::DeltaR(myCaloJetdir,theRecHitCell_XYZPoint) < maxDeltaR){
215  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
216  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
217  mySelectedDetId_.push_back(theRecHit->id());
218  }
219  }
220 
221  /*
222  for(std::vector<CaloTowerPtr>::const_iterator i_Tower=theCaloTowers.begin();i_Tower!=theCaloTowers.end();i_Tower++){
223  size_t numRecHits = (**i_Tower).constituentsSize();
224  for(size_t j=0;j<numRecHits;j++) {
225  DetId RecHitDetID=(**i_Tower).constituent(j);
226 
227 
228  DetId::Detector DetNum=RecHitDetID.det();
229  if(DetNum==DetId::Ecal){
230  if((EcalSubdetector)RecHitDetID.subdetId()==EcalBarrel){
231  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalBarrel);
232  EBDetId EcalID=RecHitDetID;
233  EBRecHitCollection::const_iterator theRecHit=EBRecHits->find(EcalID);
234  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(RecHitDetID);
235  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
236  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
237  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
238  }else if((EcalSubdetector)RecHitDetID.subdetId()==EcalEndcap){
239  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalEndcap);
240  EEDetId EcalID = RecHitDetID;
241  EERecHitCollection::const_iterator theRecHit=EERecHits->find(EcalID);
242  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(RecHitDetID);
243  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
244  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
245  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
246  }else if((EcalSubdetector)RecHitDetID.subdetId()==EcalPreshower){
247  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalPreshower);
248  ESDetId EcalID = RecHitDetID;
249  ESRecHitCollection::const_iterator theRecHit=ESRecHits->find(EcalID);
250  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(RecHitDetID);
251 
252  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
253  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
254  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
255  }
256  }
257  }
258  }
259  */
260 
261  // setting sum of Et of the isolation annulus ECAL RecHits
262  float myIsolEcalRecHits_EtSum=0.;
263 
264  std::vector< std::pair<math::XYZPoint,float> > myIsolPositionAndEnergyEcalRecHits=myCaloTauElementsOperators.EcalRecHitsInAnnulus((*myleadTk).momentum(),ECALSignalConeMetric_,myECALSignalConeSize,ECALIsolConeMetric_,myECALIsolConeSize,ECALRecHit_minEt_,thePositionAndEnergyEcalRecHits);
265  for(std::vector< std::pair<math::XYZPoint,float> >::const_iterator iEcalRecHit=myIsolPositionAndEnergyEcalRecHits.begin();iEcalRecHit!=myIsolPositionAndEnergyEcalRecHits.end();iEcalRecHit++){
266  myIsolEcalRecHits_EtSum+=(*iEcalRecHit).second*fabs(sin((*iEcalRecHit).first.theta()));
267  }
268  myCaloTau.setisolationECALhitsEtSum(myIsolEcalRecHits_EtSum);
269  }
270 
271  math::XYZTLorentzVector myTks_XYZTLorentzVect(0.,0.,0.,0.);
272  math::XYZTLorentzVector alternatLorentzVect(0.,0.,0.,0.);
273  for(TrackRefVector::iterator iTrack=myTks.begin();iTrack!=myTks.end();iTrack++) {
274  // build a charged pion candidate Lorentz vector from a Track
275  math::XYZTLorentzVector iChargedPionCand_XYZTLorentzVect((**iTrack).momentum().x(),(**iTrack).momentum().y(),(**iTrack).momentum().z(),sqrt(std::pow((double)(**iTrack).momentum().r(),2)+std::pow(chargedpi_mass_,2)));
276  myTks_XYZTLorentzVect+=iChargedPionCand_XYZTLorentzVect;
277  alternatLorentzVect+=iChargedPionCand_XYZTLorentzVect;
278  }
279  myCaloTau.setTracksInvariantMass(myTks_XYZTLorentzVect.mass());
280 
281  std::vector<BasicClusterRef> myneutralECALBasicClusters=(*myCaloTauTagInfoRef).neutralECALBasicClusters();
282  for(std::vector<BasicClusterRef>::const_iterator iBasicCluster=myneutralECALBasicClusters.begin();iBasicCluster!=myneutralECALBasicClusters.end();iBasicCluster++) {
283  // build a gamma candidate Lorentz vector from a neutral ECAL BasicCluster
284  double iGammaCand_px=(**iBasicCluster).energy()*sin((**iBasicCluster).position().theta())*cos((**iBasicCluster).position().phi());
285  double iGammaCand_py=(**iBasicCluster).energy()*sin((**iBasicCluster).position().theta())*sin((**iBasicCluster).position().phi());
286  double iGammaCand_pz=(**iBasicCluster).energy()*cos((**iBasicCluster).position().theta());
287  math::XYZTLorentzVector iGammaCand_XYZTLorentzVect(iGammaCand_px,iGammaCand_py,iGammaCand_pz,(**iBasicCluster).energy());
288  alternatLorentzVect+=iGammaCand_XYZTLorentzVect;
289  }
290  myCaloTau.setalternatLorentzVect(alternatLorentzVect);
291 
292 
293  myCaloTau.setVertex(math::XYZPoint(myCaloTau_refInnerPosition_x,myCaloTau_refInnerPosition_y,myCaloTau_refInnerPosition_z));
294 
295  // setting Et of the highest Et HCAL CaloTower
296  double mymaxEtHCALtower_Et=0.;
297  for(unsigned int iTower=0;iTower<myCaloTowers.size();iTower++){
298  if((*myCaloTowers[iTower]).hadEt()>=mymaxEtHCALtower_Et) mymaxEtHCALtower_Et=(*myCaloTowers[iTower]).hadEt();
299  }
300  myCaloTau.setmaximumHCALhitEt(mymaxEtHCALtower_Et);
301 
302  return myCaloTau;
303 }
304 
305 std::vector<CaloTowerDetId> CaloRecoTauAlgorithm::getCaloTowerneighbourDetIds(const CaloSubdetectorGeometry* myCaloSubdetectorGeometry,CaloTowerDetId myCaloTowerDetId){
306  CaloTowerTopology myCaloTowerTopology;
307  std::vector<CaloTowerDetId> myCaloTowerneighbourDetIds;
308  std::vector<DetId> northDetIds=myCaloTowerTopology.north(myCaloTowerDetId);
309  std::vector<DetId> westDetIds=myCaloTowerTopology.west(myCaloTowerDetId);
310  std::vector<DetId> northwestDetIds,southwestDetIds;
311  if (westDetIds.size()>0){
312  northwestDetIds=myCaloTowerTopology.north(westDetIds[0]);
313  southwestDetIds=myCaloTowerTopology.south(westDetIds[(int)westDetIds.size()-1]);
314  }
315  std::vector<DetId> southDetIds=myCaloTowerTopology.south(myCaloTowerDetId);
316  std::vector<DetId> eastDetIds=myCaloTowerTopology.east(myCaloTowerDetId);
317  std::vector<DetId> northeastDetIds,southeastDetIds;
318  if (eastDetIds.size()>0){
319  northeastDetIds=myCaloTowerTopology.north(eastDetIds[0]);
320  southeastDetIds=myCaloTowerTopology.south(eastDetIds[(int)eastDetIds.size()-1]);
321  }
322  std::vector<DetId> myneighbourDetIds=northDetIds;
323  myneighbourDetIds.insert(myneighbourDetIds.end(),westDetIds.begin(),westDetIds.end());
324  myneighbourDetIds.insert(myneighbourDetIds.end(),northwestDetIds.begin(),northwestDetIds.end());
325  myneighbourDetIds.insert(myneighbourDetIds.end(),southwestDetIds.begin(),southwestDetIds.end());
326  myneighbourDetIds.insert(myneighbourDetIds.end(),southDetIds.begin(),southDetIds.end());
327  myneighbourDetIds.insert(myneighbourDetIds.end(),eastDetIds.begin(),eastDetIds.end());
328  myneighbourDetIds.insert(myneighbourDetIds.end(),northeastDetIds.begin(),northeastDetIds.end());
329  myneighbourDetIds.insert(myneighbourDetIds.end(),southeastDetIds.begin(),southeastDetIds.end());
330  for(std::vector<DetId>::const_iterator iDetId=myneighbourDetIds.begin();iDetId!=myneighbourDetIds.end();iDetId++){
331  CaloTowerDetId iCaloTowerId(*iDetId);
332  myCaloTowerneighbourDetIds.push_back(iCaloTowerId);
333  }
334  return myCaloTowerneighbourDetIds;
335 }
336 
double computeConeSize(const TFormula &ConeSizeTFormula, double ConeSizeMin, double ConeSizeMax)
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
reco::TrackRefVector filteredTracksByNumTrkHits(reco::TrackRefVector theInitialTracks, int tkminTrackerHitsn)
Definition: TauTagTools.cc:68
std::vector< CaloTowerDetId > getCaloTowerneighbourDetIds(const CaloSubdetectorGeometry *, CaloTowerDetId)
virtual std::vector< DetId > south(const DetId &id) const
std::vector< DetId > mySelectedDetId_
const reco::TrackRef leadTk(std::string matchingConeMetric, double matchingConeSize, double ptTrackMin) const
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
std::pair< bool, Measurement1D > signedTransverseImpactParameter(const reco::TransientTrack &track, const GlobalVector &direction, const reco::Vertex &vertex)
Definition: IPTools.cc:50
math::XYZPoint propagTrackECALSurfContactPoint(const MagneticField *, reco::TrackRef)
Definition: TauTagTools.cc:184
std::vector< EcalRecHit >::const_iterator const_iterator
reco::TransientTrack build(const reco::Track *p) const
T y() const
Definition: PV3DBase.h:62
std::string TrackerSignalConeMetric_
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
const_iterator end() const
Termination of iteration.
Definition: RefVector.h:249
const Point & position() const
position
Definition: Vertex.h:93
std::string ECALIsolConeSizeFormula_
double double double z
const_iterator begin() const
Initialize an iterator over the RefVector.
Definition: RefVector.h:244
virtual const CaloCellGeometry * getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
std::string TrackerSignalConeSizeFormula_
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:30
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:250
const reco::TrackRefVector tracksInCone(const math::XYZVector &coneAxis, const std::string coneMetric, const double coneSize, const double ptTrackMin) const
int iEvent
Definition: GenABIO.cc:243
reco::CaloTau buildCaloTau(edm::Event &, const edm::EventSetup &, const reco::CaloTauTagInfoRef &, const reco::Vertex &)
std::string TrackerIsolConeMetric_
static const int SubdetId
T sqrt(T t)
Definition: SSEVec.h:46
T z() const
Definition: PV3DBase.h:63
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
const MagneticField * MagneticField_
virtual std::vector< DetId > north(const DetId &id) const
TFormula myTrackerSignalConeSizeTFormula
std::string TrackerIsolConeSizeFormula_
std::string ECALSignalConeSizeFormula_
TFormula computeConeSizeTFormula(const std::string &ConeSizeFormula, const char *errorMessage)
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:356
void setMagneticField(const MagneticField *)
std::string MatchingConeSizeFormula_
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:13
virtual std::vector< DetId > west(const DetId &id) const
const T & get() const
Definition: EventSetup.h:55
std::vector< std::pair< math::XYZPoint, float > > EcalRecHitsInAnnulus(const math::XYZVector &coneAxis, const std::string innerconeMetric, const double innerconeSize, const std::string outerconeMetric, const double outerconeSize, const double EcalRecHit_minEt, const std::vector< std::pair< math::XYZPoint, float > > &myEcalRecHits) const
unsigned int IsolationTrack_minHits_
void setcaloTauTagInfoRef(const CaloTauTagInfoRef)
Definition: CaloTau.cc:30
virtual std::vector< DetId > east(const DetId &id) const
void push_back(value_type const &ref)
Add a Ref&lt;C, T&gt; to the RefVector.
Definition: RefVector.h:64
size_type size() const
Size of the RefVector.
Definition: RefVector.h:89
Definition: DDAxes.h:10
const reco::TrackRefVector tracksInAnnulus(const math::XYZVector &coneAxis, const std::string innerconeMetric, const double innerconeSize, const std::string outerconeMetric, const double outerconeSize, const double ptTrackMin) const
T x() const
Definition: PV3DBase.h:61
virtual const LorentzVector & p4() const
four-momentum Lorentz vector
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
const TransientTrackBuilder * TransientTrackBuilder_
void setTransientTrackBuilder(const TransientTrackBuilder *)