LCFIVertex  0.7.2
twotrackpid.cpp
1 #include <algo/inc/twotrackpid.h>
2 #include <string>
3 #include <vector>
4 #include <inc/jet.h>
5 #include <inc/track.h>
6 #include <inc/trackstate.h>
7 #include <zvtop/include/interactionpoint.h>
8 #include <zvtop/include/vertexfitterlsm.h>
9 #include <util/inc/string.h>
10 #include <util/inc/memorymanager.h>
11 
12 using std::string;
13 using namespace vertex_lcfi::util;
14 
15 namespace vertex_lcfi
16 {
17  using namespace ZVTOP;
18 
19  TwoTrackPid::TwoTrackPid()
20  {
21  _MaxGammaMass = 0.02;
22  _MinKsMass = 0.475;
23  _MaxKsMass = 0.525;
24  _Chi2Cut = 6.63;
25  _RPhiCut = 20;
26  _SignificanceCut =3;
27  _ParameterNames.push_back("MaxGammaMass");
28  _ParameterNames.push_back("MinKsMass");
29  _ParameterNames.push_back("MaxKsMass");
30  _ParameterNames.push_back("Chi2Cut");
31  _ParameterNames.push_back("RPhiCut");
32  _ParameterNames.push_back("SignificanceCut");
33 
34  }
35 
36  string TwoTrackPid::name() const
37  {
38  return _Name;
39  }
40 
41  std::vector<string> TwoTrackPid::parameterNames() const
42  {
43  return _ParameterNames;
44  }
45 
46  std::vector<string> TwoTrackPid::parameterValues() const
47  {
48  _ParameterValues.clear();
49  _ParameterValues.push_back(makeString(_MaxGammaMass));
50  _ParameterValues.push_back(makeString(_MinKsMass));
51  _ParameterValues.push_back(makeString(_MaxKsMass));
52  _ParameterValues.push_back(makeString(_Chi2Cut));
53  _ParameterValues.push_back(makeString(_RPhiCut));
54  _ParameterValues.push_back(makeString(_SignificanceCut));
55  return _ParameterValues;
56  }
57 
58  void TwoTrackPid::setStringParameter(const string & Parameter, const string & Value)
59  {
60  this->badParameter(Parameter);
61  }
62 
63  void TwoTrackPid::setDoubleParameter(const string & Parameter, const double Value)
64  {
65  if (Parameter == "MaxGammaMass")
66  {
67  _MaxGammaMass = Value;
68  return;
69  }
70  if (Parameter == "MinKsMass")
71  {
72  _MinKsMass = Value;
73  return;
74  }
75  if (Parameter == "MaxKsMass")
76  {
77  _MaxKsMass = Value;
78  return;
79  }
80  if (Parameter == "Chi2Cut")
81  {
82  _Chi2Cut = Value;
83  return;
84  }
85  if (Parameter == "RPhiCut")
86  {
87  _RPhiCut = Value;
88  return;
89  }
90  if (Parameter == "SignificanceCut")
91  {
92  _SignificanceCut = Value;
93  return;
94  }
95  this->badParameter(Parameter);
96  }
97 
98  void TwoTrackPid::setPointerParameter(const string & Parameter, void * Value)
99  {
100  this->badParameter(Parameter);
101  }
102 
103  std::map< PidCutType,std::vector<Track*> > TwoTrackPid::calculateFor(Jet* MyJet) const
104  {
105 
106  Vector3 Position;
107 
108  double chi2;
109 
110  double RPhiProjection;
111  InteractionPoint* IP = 0;
112  double momentummagnitude;
113  std::vector< TrackState* > AllTrackStates;
114  VertexFitterLSM Fitter;
115  double eeCalculatedM = 0;
116  double pipiCalculatedM = 0;
117  double electronmass2 = 0.00051*0.00051;
118  double pionmass2 = 0.1396*0.1396;
119 
120  double MaxGammaMass = _MaxGammaMass;
121  double MinKsMass = _MinKsMass;
122  double MaxKsMass = _MaxKsMass;
123  double Chi2Cut = _Chi2Cut;
124  double RPhiCut = _RPhiCut;
125  double significancecut = _SignificanceCut;
126 
127  //just a dummy for initialization
128  std::vector<Track*> Dummy;
129  std::map<PidCutType,std::vector<Track*> > ResultMap;
130 
131  ResultMap[Gamma] = Dummy;
132  ResultMap[KShort] = Dummy;
133 
134  //nested loop running over pairs of tracks.
135 
136 
137  if( MyJet->tracks().size()>0)
138  {
139  for(std::vector<Track*>::const_iterator iTrack= (MyJet->tracks().begin()); iTrack != --(MyJet->tracks().end()) ;++iTrack)
140  {
141 
142  for(std::vector<Track*>::const_iterator iTrack2 = iTrack+1 ; iTrack2 != (MyJet->tracks().end()) ;++iTrack2)
143  {
144  //check first track d0 significance
145 
146  if( (*iTrack)->significance(RPhi) > significancecut )
147  {
148 
149  AllTrackStates.push_back( (**iTrack).makeState() );
150 
151  //check second track d0 significance
152  if( (*iTrack2)->significance(RPhi) > significancecut )
153  {
154 
155  //chack that the charge of the tracks is opposite)
156  if (((*iTrack)->charge() + (*iTrack2)->charge()) == 0)
157  {
158  //vertex the tracks and see what happens.
159  AllTrackStates.push_back( (**iTrack2).makeState() );
160 
161  Fitter.fitVertex( AllTrackStates, IP, Position, chi2 );
162 
163  //the rphi projection
164  RPhiProjection = Position.mag(RPhi);
165 
166  //cuts on position of the vertecs and its significance.
167  if( RPhiProjection < RPhiCut && chi2< Chi2Cut )
168  {
169  Vector3 Totalmomentum;
170  Totalmomentum = Totalmomentum.add((**iTrack2).momentum());
171  Totalmomentum = Totalmomentum.add((**iTrack).momentum());
172  momentummagnitude = Totalmomentum.mag2();
173 
174  //here we are just using the standard e^2 =m^2 +p^2s
175  //first with the electron mass
176 
177  eeCalculatedM = pow(sqrt((**iTrack2).momentum().mag2()+electronmass2)+
178  sqrt((**iTrack).momentum().mag2()+electronmass2),2)- momentummagnitude;
179  if (eeCalculatedM >0)
180  {
181  eeCalculatedM = sqrt(eeCalculatedM);
182 
183  }
184  else
185  {
186  eeCalculatedM = 0;
187  }
188 
189  //then assume pion mass
190 
191  pipiCalculatedM = pow(sqrt((**iTrack2).momentum().mag2()+pionmass2)+
192  sqrt((**iTrack).momentum().mag2()+pionmass2),2)- momentummagnitude;
193  if (pipiCalculatedM >0)
194  {
195  pipiCalculatedM = sqrt(pipiCalculatedM);
196  }
197  else
198  {
199  pipiCalculatedM = 0;
200  }
201 
202 
203  //here we assume that a track will never end up in both
204  //given the different mass I think it is a fair assumption
205  if(eeCalculatedM < MaxGammaMass)
206  {
207 
208  //these algorithms basically put the tracks into the map in the case that they are not there already!
209  std::vector<Track*>::const_iterator iTrack3 = find(ResultMap[Gamma].begin(),ResultMap[Gamma].end(), (*iTrack));
210  if(iTrack3 == ResultMap[Gamma].end() )
211  {
212  ResultMap[Gamma].push_back(*iTrack);
213  }
214 
215 
216  std::vector<Track*>::const_iterator iTrack4 = find(ResultMap[Gamma].begin(),ResultMap[Gamma].end(), (*iTrack2));
217  if(iTrack4 == ResultMap[Gamma].end() )
218  {
219  ResultMap[Gamma].push_back(*iTrack2);
220  }
221  }
222  else if( MinKsMass < pipiCalculatedM && pipiCalculatedM < MaxKsMass )
223  {
224  std::vector<Track*>::const_iterator iTrack3 = find(ResultMap[KShort].begin(),ResultMap[KShort].end(), (*iTrack));
225  if(iTrack3 == ResultMap[KShort].end() )
226  {
227  ResultMap[KShort].push_back(*iTrack);
228 
229  }
230  std::vector<Track*>::const_iterator iTrack4 = find(ResultMap[KShort].begin(),ResultMap[KShort].end(), (*iTrack2));
231  if(iTrack4 == ResultMap[KShort].end() )
232  {
233  ResultMap[KShort].push_back(*iTrack2);
234 
235  }
236  }
237  }
238 
239  }
240  }
241  }
242 
243  AllTrackStates.clear();
244  }
245  }
246  }
247  return ResultMap;
248 
249  }
250 
251 }
252 
253 
254 
void setDoubleParameter(const string &Parameter, const double Value)
Set Double Parameter.
Definition: twotrackpid.cpp:63
const std::vector< Track * > & tracks() const
Tracks.
Definition: jet.cpp:23
void setStringParameter(const string &Parameter, const string &Value)
Set String Parameter.
Definition: twotrackpid.cpp:58
std::map< PidCutType, std::vector< Track * > > calculateFor(Jet *MyJet) const
Run the algorithm on a Jet.
void setPointerParameter(const string &Parameter, void *Value)
Set Pointer Parameter.
Definition: twotrackpid.cpp:98
string name() const
Name.
Definition: twotrackpid.cpp:36
std::vector< string > parameterValues() const
Parameter Values.
Definition: twotrackpid.cpp:46
std::vector< string > parameterNames() const
Parameter Names.
Definition: twotrackpid.cpp:41