test
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  edm::ESHandle<CaloTowerTopology> caloTowerTopology;
101  iSetup.get<HcalRecNumberingRecord>().get(caloTowerTopology);
102  CaloTowerDetId mypropagleadTrack_closestCaloTowerId((*myCaloSubdetectorGeometry).getClosestCell(GlobalPoint(mypropagleadTrackECALSurfContactPoint.x(),
103  mypropagleadTrackECALSurfContactPoint.y(),
104  mypropagleadTrackECALSurfContactPoint.z())));
105  std::vector<CaloTowerDetId> mypropagleadTrack_closestCaloTowerNeighbourIds=getCaloTowerneighbourDetIds(myCaloSubdetectorGeometry, *caloTowerTopology, mypropagleadTrack_closestCaloTowerId);
106  for(std::vector<CaloTowerPtr>::const_iterator iCaloTower=myCaloTowers.begin();iCaloTower!=myCaloTowers.end();iCaloTower++){
107  CaloTowerDetId iCaloTowerId((**iCaloTower).id());
108  bool CaloTower_inside3x3matrix=false;
109  if (iCaloTowerId==mypropagleadTrack_closestCaloTowerId) CaloTower_inside3x3matrix=true;
110  if (!CaloTower_inside3x3matrix){
111  for(std::vector<CaloTowerDetId>::const_iterator iCaloTowerDetId=mypropagleadTrack_closestCaloTowerNeighbourIds.begin();iCaloTowerDetId!=mypropagleadTrack_closestCaloTowerNeighbourIds.end();iCaloTowerDetId++){
112  if (iCaloTowerId==(*iCaloTowerDetId)){
113  CaloTower_inside3x3matrix=true;
114  break;
115  }
116  }
117  }
118  if (!CaloTower_inside3x3matrix) continue;
119  myleadTrackHCAL3x3hitsEtSum+=(**iCaloTower).hadEt();
120  if((**iCaloTower).hadEt()>=myleadTrackHCAL3x3hottesthitEt ){
121  if ((**iCaloTower).hadEt()!=myleadTrackHCAL3x3hottesthitEt ||
122  ((**iCaloTower).hadEt()==myleadTrackHCAL3x3hottesthitEt && fabs((**iCaloTower).eta()-mypropagleadTrackECALSurfContactPoint.Eta())<myleadTrackHCAL3x3hottesthitDEta)) myleadTrackHCAL3x3hottesthitDEta = fabs((**iCaloTower).eta()-mypropagleadTrackECALSurfContactPoint.Eta());
123  myleadTrackHCAL3x3hottesthitEt=(**iCaloTower).hadEt();
124  }
125  }
126  myCaloTau.setleadTrackHCAL3x3hitsEtSum(myleadTrackHCAL3x3hitsEtSum);
127  if (myleadTrackHCAL3x3hottesthitEt!=0.) myCaloTau.setleadTrackHCAL3x3hottesthitDEta(myleadTrackHCAL3x3hottesthitDEta);
128  }
129  }
130 
132  TrackRefVector myTksbis;
133  for (TrackRefVector::const_iterator iTrack=myTks.begin();iTrack!=myTks.end();++iTrack) {
134  if (fabs((**iTrack).dz(myPV.position())-myleadTkDZ)<=TrackLeadTrack_maxDZ_) myTksbis.push_back(*iTrack);
135  }
136  myTks=myTksbis;
137  }
138 
139 
140  double myTrackerSignalConeSize=myCaloTauElementsOperators.computeConeSize(myTrackerSignalConeSizeTFormula,TrackerSignalConeSize_min_,TrackerSignalConeSize_max_);
141  double myTrackerIsolConeSize=myCaloTauElementsOperators.computeConeSize(myTrackerIsolConeSizeTFormula,TrackerIsolConeSize_min_,TrackerIsolConeSize_max_);
142  double myECALSignalConeSize=myCaloTauElementsOperators.computeConeSize(myECALSignalConeSizeTFormula,ECALSignalConeSize_min_,ECALSignalConeSize_max_);
143  double myECALIsolConeSize=myCaloTauElementsOperators.computeConeSize(myECALIsolConeSizeTFormula,ECALIsolConeSize_min_,ECALIsolConeSize_max_);
144 
145  TrackRefVector mySignalTks;
146  if (UseTrackLeadTrackDZconstraint_) mySignalTks=myCaloTauElementsOperators.tracksInCone((*myleadTk).momentum(),TrackerSignalConeMetric_,myTrackerSignalConeSize,Track_minPt_,TrackLeadTrack_maxDZ_,myleadTkDZ, myPV);
147  else mySignalTks=myCaloTauElementsOperators.tracksInCone((*myleadTk).momentum(),TrackerSignalConeMetric_,myTrackerSignalConeSize,Track_minPt_);
148  myCaloTau.setsignalTracks(mySignalTks);
149 
150  // setting invariant mass of the signal Tracks system
151  math::XYZTLorentzVector mySignalTksInvariantMass(0.,0.,0.,0.);
152  if((int)(mySignalTks.size())!=0){
153  int mySignalTks_qsum=0;
154  for(int i=0;i<(int)mySignalTks.size();i++){
155  mySignalTks_qsum+=mySignalTks[i]->charge();
156  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)));
157  mySignalTksInvariantMass+=mychargedpicand_fromTk_LorentzVect;
158  }
159  myCaloTau.setCharge(mySignalTks_qsum);
160  }
161  myCaloTau.setsignalTracksInvariantMass(mySignalTksInvariantMass.mass());
162 
163  TrackRefVector myIsolTks;
164  if (UseTrackLeadTrackDZconstraint_) myIsolTks=myCaloTauElementsOperators.tracksInAnnulus((*myleadTk).momentum(),TrackerSignalConeMetric_,myTrackerSignalConeSize,TrackerIsolConeMetric_,myTrackerIsolConeSize,IsolationTrack_minPt_,TrackLeadTrack_maxDZ_,myleadTkDZ, myPV);
165  else myIsolTks=myCaloTauElementsOperators.tracksInAnnulus((*myleadTk).momentum(),TrackerSignalConeMetric_,myTrackerSignalConeSize,TrackerIsolConeMetric_,myTrackerIsolConeSize,IsolationTrack_minPt_);
167  myCaloTau.setisolationTracks(myIsolTks);
168 
169  // setting sum of Pt of the isolation annulus Tracks
170  float myIsolTks_Ptsum=0.;
171  for(int i=0;i<(int)myIsolTks.size();i++) myIsolTks_Ptsum+=myIsolTks[i]->pt();
172  myCaloTau.setisolationTracksPtSum(myIsolTks_Ptsum);
173 
174 
175  //getting the EcalRecHits. Just take them all
176  std::vector<std::pair<math::XYZPoint,float> > thePositionAndEnergyEcalRecHits;
177  mySelectedDetId_.clear();
178  // std::vector<CaloTowerPtr> theCaloTowers=myCaloJet->getCaloConstituents();
179  edm::ESHandle<CaloGeometry> theCaloGeometry;
180  iSetup.get<CaloGeometryRecord>().get(theCaloGeometry);
181  const CaloSubdetectorGeometry* theCaloSubdetectorGeometry;
185  iEvent.getByLabel(EBRecHitsLabel_,EBRecHits);
186  iEvent.getByLabel(EERecHitsLabel_,EERecHits);
187  iEvent.getByLabel(ESRecHitsLabel_,ESRecHits);
188  double maxDeltaR = 0.8;
189  math::XYZPoint myCaloJetdir((*myCaloJet).px(),(*myCaloJet).py(),(*myCaloJet).pz());
190 
191  for(EBRecHitCollection::const_iterator theRecHit = EBRecHits->begin();theRecHit != EBRecHits->end(); theRecHit++){
192  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalBarrel);
193  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(theRecHit->id());
194  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
195  if(ROOT::Math::VectorUtil::DeltaR(myCaloJetdir,theRecHitCell_XYZPoint) < maxDeltaR){
196  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
197  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
198  mySelectedDetId_.push_back(theRecHit->id());
199  }
200  }
201 
202 for(EERecHitCollection::const_iterator theRecHit = EERecHits->begin();theRecHit != EERecHits->end(); theRecHit++){
203  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalEndcap);
204  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(theRecHit->id());
205  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
206  if(ROOT::Math::VectorUtil::DeltaR(myCaloJetdir,theRecHitCell_XYZPoint) < maxDeltaR){
207  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
208  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
209  mySelectedDetId_.push_back(theRecHit->id());
210  }
211 }
212  for(ESRecHitCollection::const_iterator theRecHit = ESRecHits->begin();theRecHit != ESRecHits->end(); theRecHit++){
213  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalPreshower);
214  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(theRecHit->id());
215  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
216  if(ROOT::Math::VectorUtil::DeltaR(myCaloJetdir,theRecHitCell_XYZPoint) < maxDeltaR){
217  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
218  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
219  mySelectedDetId_.push_back(theRecHit->id());
220  }
221  }
222 
223  /*
224  for(std::vector<CaloTowerPtr>::const_iterator i_Tower=theCaloTowers.begin();i_Tower!=theCaloTowers.end();i_Tower++){
225  size_t numRecHits = (**i_Tower).constituentsSize();
226  for(size_t j=0;j<numRecHits;j++) {
227  DetId RecHitDetID=(**i_Tower).constituent(j);
228 
229 
230  DetId::Detector DetNum=RecHitDetID.det();
231  if(DetNum==DetId::Ecal){
232  if((EcalSubdetector)RecHitDetID.subdetId()==EcalBarrel){
233  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalBarrel);
234  EBDetId EcalID=RecHitDetID;
235  EBRecHitCollection::const_iterator theRecHit=EBRecHits->find(EcalID);
236  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(RecHitDetID);
237  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
238  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
239  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
240  }else if((EcalSubdetector)RecHitDetID.subdetId()==EcalEndcap){
241  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalEndcap);
242  EEDetId EcalID = RecHitDetID;
243  EERecHitCollection::const_iterator theRecHit=EERecHits->find(EcalID);
244  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(RecHitDetID);
245  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
246  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
247  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
248  }else if((EcalSubdetector)RecHitDetID.subdetId()==EcalPreshower){
249  theCaloSubdetectorGeometry = theCaloGeometry->getSubdetectorGeometry(DetId::Ecal,EcalPreshower);
250  ESDetId EcalID = RecHitDetID;
251  ESRecHitCollection::const_iterator theRecHit=ESRecHits->find(EcalID);
252  const CaloCellGeometry* theRecHitCell=theCaloSubdetectorGeometry->getGeometry(RecHitDetID);
253 
254  math::XYZPoint theRecHitCell_XYZPoint(theRecHitCell->getPosition().x(),theRecHitCell->getPosition().y(),theRecHitCell->getPosition().z());
255  std::pair<math::XYZPoint,float> thePositionAndEnergyEcalRecHit(theRecHitCell_XYZPoint,theRecHit->energy());
256  thePositionAndEnergyEcalRecHits.push_back(thePositionAndEnergyEcalRecHit);
257  }
258  }
259  }
260  }
261  */
262 
263  // setting sum of Et of the isolation annulus ECAL RecHits
264  float myIsolEcalRecHits_EtSum=0.;
265 
266  std::vector< std::pair<math::XYZPoint,float> > myIsolPositionAndEnergyEcalRecHits=myCaloTauElementsOperators.EcalRecHitsInAnnulus((*myleadTk).momentum(),ECALSignalConeMetric_,myECALSignalConeSize,ECALIsolConeMetric_,myECALIsolConeSize,ECALRecHit_minEt_,thePositionAndEnergyEcalRecHits);
267  for(std::vector< std::pair<math::XYZPoint,float> >::const_iterator iEcalRecHit=myIsolPositionAndEnergyEcalRecHits.begin();iEcalRecHit!=myIsolPositionAndEnergyEcalRecHits.end();iEcalRecHit++){
268  myIsolEcalRecHits_EtSum+=(*iEcalRecHit).second*fabs(sin((*iEcalRecHit).first.theta()));
269  }
270  myCaloTau.setisolationECALhitsEtSum(myIsolEcalRecHits_EtSum);
271  }
272 
273  math::XYZTLorentzVector myTks_XYZTLorentzVect(0.,0.,0.,0.);
274  math::XYZTLorentzVector alternatLorentzVect(0.,0.,0.,0.);
275  for(TrackRefVector::iterator iTrack=myTks.begin();iTrack!=myTks.end();iTrack++) {
276  // build a charged pion candidate Lorentz vector from a Track
277  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)));
278  myTks_XYZTLorentzVect+=iChargedPionCand_XYZTLorentzVect;
279  alternatLorentzVect+=iChargedPionCand_XYZTLorentzVect;
280  }
281  myCaloTau.setTracksInvariantMass(myTks_XYZTLorentzVect.mass());
282 
283  std::vector<BasicClusterRef> myneutralECALBasicClusters=(*myCaloTauTagInfoRef).neutralECALBasicClusters();
284  for(std::vector<BasicClusterRef>::const_iterator iBasicCluster=myneutralECALBasicClusters.begin();iBasicCluster!=myneutralECALBasicClusters.end();iBasicCluster++) {
285  // build a gamma candidate Lorentz vector from a neutral ECAL BasicCluster
286  double iGammaCand_px=(**iBasicCluster).energy()*sin((**iBasicCluster).position().theta())*cos((**iBasicCluster).position().phi());
287  double iGammaCand_py=(**iBasicCluster).energy()*sin((**iBasicCluster).position().theta())*sin((**iBasicCluster).position().phi());
288  double iGammaCand_pz=(**iBasicCluster).energy()*cos((**iBasicCluster).position().theta());
289  math::XYZTLorentzVector iGammaCand_XYZTLorentzVect(iGammaCand_px,iGammaCand_py,iGammaCand_pz,(**iBasicCluster).energy());
290  alternatLorentzVect+=iGammaCand_XYZTLorentzVect;
291  }
292  myCaloTau.setalternatLorentzVect(alternatLorentzVect);
293 
294 
295  myCaloTau.setVertex(math::XYZPoint(myCaloTau_refInnerPosition_x,myCaloTau_refInnerPosition_y,myCaloTau_refInnerPosition_z));
296 
297  // setting Et of the highest Et HCAL CaloTower
298  double mymaxEtHCALtower_Et=0.;
299  for(unsigned int iTower=0;iTower<myCaloTowers.size();iTower++){
300  if((*myCaloTowers[iTower]).hadEt()>=mymaxEtHCALtower_Et) mymaxEtHCALtower_Et=(*myCaloTowers[iTower]).hadEt();
301  }
302  myCaloTau.setmaximumHCALhitEt(mymaxEtHCALtower_Et);
303 
304  return myCaloTau;
305 }
306 
307 std::vector<CaloTowerDetId> CaloRecoTauAlgorithm::getCaloTowerneighbourDetIds(const CaloSubdetectorGeometry* myCaloSubdetectorGeometry, const CaloTowerTopology & myCaloTowerTopology, CaloTowerDetId myCaloTowerDetId){
308  std::vector<CaloTowerDetId> myCaloTowerneighbourDetIds;
309  std::vector<DetId> northDetIds=myCaloTowerTopology.north(myCaloTowerDetId);
310  std::vector<DetId> westDetIds=myCaloTowerTopology.west(myCaloTowerDetId);
311  std::vector<DetId> northwestDetIds,southwestDetIds;
312  if (westDetIds.size()>0){
313  northwestDetIds=myCaloTowerTopology.north(westDetIds[0]);
314  southwestDetIds=myCaloTowerTopology.south(westDetIds[(int)westDetIds.size()-1]);
315  }
316  std::vector<DetId> southDetIds=myCaloTowerTopology.south(myCaloTowerDetId);
317  std::vector<DetId> eastDetIds=myCaloTowerTopology.east(myCaloTowerDetId);
318  std::vector<DetId> northeastDetIds,southeastDetIds;
319  if (eastDetIds.size()>0){
320  northeastDetIds=myCaloTowerTopology.north(eastDetIds[0]);
321  southeastDetIds=myCaloTowerTopology.south(eastDetIds[(int)eastDetIds.size()-1]);
322  }
323  std::vector<DetId> myneighbourDetIds=northDetIds;
324  myneighbourDetIds.insert(myneighbourDetIds.end(),westDetIds.begin(),westDetIds.end());
325  myneighbourDetIds.insert(myneighbourDetIds.end(),northwestDetIds.begin(),northwestDetIds.end());
326  myneighbourDetIds.insert(myneighbourDetIds.end(),southwestDetIds.begin(),southwestDetIds.end());
327  myneighbourDetIds.insert(myneighbourDetIds.end(),southDetIds.begin(),southDetIds.end());
328  myneighbourDetIds.insert(myneighbourDetIds.end(),eastDetIds.begin(),eastDetIds.end());
329  myneighbourDetIds.insert(myneighbourDetIds.end(),northeastDetIds.begin(),northeastDetIds.end());
330  myneighbourDetIds.insert(myneighbourDetIds.end(),southeastDetIds.begin(),southeastDetIds.end());
331  for(std::vector<DetId>::const_iterator iDetId=myneighbourDetIds.begin();iDetId!=myneighbourDetIds.end();iDetId++){
332  CaloTowerDetId iCaloTowerId(*iDetId);
333  myCaloTowerneighbourDetIds.push_back(iCaloTowerId);
334  }
335  return myCaloTowerneighbourDetIds;
336 }
337 
double computeConeSize(const TFormula &ConeSizeTFormula, double ConeSizeMin, double ConeSizeMax)
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:252
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:63
std::string TrackerSignalConeMetric_
const_iterator end() const
Termination of iteration.
Definition: RefVector.h:255
const Point & position() const
position
Definition: Vertex.h:99
std::string ECALIsolConeSizeFormula_
const_iterator begin() const
Initialize an iterator over the RefVector.
Definition: RefVector.h:250
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:29
T x() const
Cartesian x coordinate.
const reco::TrackRefVector tracksInCone(const math::XYZVector &coneAxis, const std::string coneMetric, const double coneSize, const double ptTrackMin) const
int iEvent
Definition: GenABIO.cc:230
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:18
T z() const
Definition: PV3DBase.h:64
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_
std::vector< CaloTowerDetId > getCaloTowerneighbourDetIds(const CaloSubdetectorGeometry *, const CaloTowerTopology &, CaloTowerDetId)
TFormula computeConeSizeTFormula(const std::string &ConeSizeFormula, const char *errorMessage)
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:418
void setMagneticField(const MagneticField *)
std::string MatchingConeSizeFormula_
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
virtual std::vector< DetId > west(const DetId &id) const
const T & get() const
Definition: EventSetup.h:56
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
reco::TrackRefVector filteredTracksByNumTrkHits(const reco::TrackRefVector &theInitialTracks, int tkminTrackerHitsn)
Definition: TauTagTools.cc:68
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:62
size_type size() const
Size of the RefVector.
Definition: RefVector.h:99
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
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
T x() const
Definition: PV3DBase.h:62
virtual const LorentzVector & p4() const final
four-momentum Lorentz vector
Definition: LeafCandidate.h:99
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:40
const TransientTrackBuilder * TransientTrackBuilder_
void setTransientTrackBuilder(const TransientTrackBuilder *)