LCFIPlus  0.6.5
VertexSelector.h
Go to the documentation of this file.
1 // VertexSelector.h
2 
3 #ifndef VertexSelector_h
4 #define VertexSelector_h 1
5 
6 #include "lcfiplus.h"
7 #include <vector>
8 #include <cmath>
9 
10 using namespace std;
11 using namespace lcfiplus;
12 
13 namespace lcfiplus {
14 
16  public:
17  bool rejectdist;
18  double minpos;
19  double maxpos;
22 
23  bool rejectk0;
24  double k0width;
25  double k0dirdot;
26 
27  bool rejectl0;
28  double l0width;
29  double l0dirdot;
30 
31  bool rejectconv;
32  double convmass;
33  double convdist;
34  double convdirdot;
35 
37  rejectdist(false), minpos(0.), maxpos(1e+300),
38  rejectdistnegative(false), rejectdistor(false),
39  rejectk0(false), k0width(0.), k0dirdot(0.),
40  rejectl0(false), l0width(0.), l0dirdot(0.),
41  rejectconv(false), convmass(0.), convdist(0.), convdirdot(0.) {
42  }
43 
44  void setV0Tight() { // for rejecting high purity v0 only
45  rejectk0 = true;
46  k0width = 0.005;
47  k0dirdot = 0.999;
48  rejectl0 = true;
49  l0width = 0.005;
50  l0dirdot = 0.99995;
51  rejectconv = true;
52  convmass = 0.005;
53  convdist = 9;
54  convdirdot = 0.99995;
55  }
56  void setV0Loose() { // for rejecting high purity v0 only
57  rejectk0 = true;
58  k0width = 0.010;
59  k0dirdot = 0.999;
60  rejectl0 = true;
61  l0width = 0.010;
62  l0dirdot = 0.999;
63  rejectconv = true;
64  convmass = 0.010;
65  convdist = 9;
66  convdirdot = 0.999;
67  }
68  void setNoV0Cut() {
69  rejectk0 = false;
70  rejectl0 = false;
71  rejectconv = false;
72  }
73 
74 };
75 
77  public:
78  // const full version
79  vector<const Vertex*> operator () (VertexVec& vertices, VertexSelectorConfig& config, vector<const Track*>& residualTracks, bool addTracks) {
80  vector<const Vertex*> ret;
81 
82  for (unsigned int i=0; i<vertices.size(); i++) {
83  if (vertices[i]->isPrimary()) {
84  //cout << "Primary vertex found in Vertex Selector!!!!!!!" << endl;
85  continue;
86  }
87 
88  if (passesCut(vertices[i], config)) {
89  ret.push_back(vertices[i]);
90  if (!addTracks) {
91  // TODO: possible use of STL algorithm set_difference - but need to sort
92  for (TrackVecIte itt = vertices[i]->getTracks().begin(); itt != vertices[i]->getTracks().end(); itt++) {
93  vector<const Track*>::iterator itt2 = remove_if(residualTracks.begin(), residualTracks.end(), bind2nd(equal_to<const Track*>(), *itt));
94  residualTracks.erase(itt2, residualTracks.end());
95  }
96  }
97  } else if (addTracks) {
98  residualTracks.insert(residualTracks.end(),vertices[i]->getTracks().begin(), vertices[i]->getTracks().end());
99  }
100  }
101 
102  return ret;
103  }
104  // non-const full version
105  vector<Vertex*> operator () (const vector<Vertex*>& vertices, VertexSelectorConfig& config, vector<const Track*>& residualTracks, bool addTracks) {
106  vector<Vertex*> ret;
107 
108  for (unsigned int i=0; i<vertices.size(); i++) {
109  if (passesCut(vertices[i], config)) {
110  ret.push_back(vertices[i]);
111  if (!addTracks) {
112  // TODO: possible use of STL algorithm set_difference - but need to sort
113  for (TrackVecIte itt = vertices[i]->getTracks().begin(); itt != vertices[i]->getTracks().end(); itt++) {
114  vector<const Track*>::iterator itt2 = remove_if(residualTracks.begin(), residualTracks.end(), bind2nd(equal_to<const Track*>(), *itt));
115  residualTracks.erase(itt2, residualTracks.end());
116  }
117  }
118  } else if (addTracks) {
119  residualTracks.insert(residualTracks.end(),vertices[i]->getTracks().begin(), vertices[i]->getTracks().end());
120  }
121  }
122 
123  return ret;
124  }
125  // const simple version
126  vector<const Vertex*> operator () (VertexVec& vertices, VertexSelectorConfig& config) {
127  vector<const Vertex*> ret;
128 
129  for (unsigned int i=0; i<vertices.size(); i++) {
130  if (passesCut(vertices[i], config)) {
131  ret.push_back(vertices[i]);
132  }
133  }
134 
135  return ret;
136  }
137  // non-const very simple version
138  void operator () (vector<Vertex*>& vertices, VertexSelectorConfig& config) {
139 
140  for (vector<Vertex*>::iterator it = vertices.begin(); it!=vertices.end();) {
141  if (!passesCut(*it, config)) {
142  it = vertices.erase(it);
143  continue;
144  }
145  it ++;
146  }
147  }
148 
149  // true if vtx is V0
150  bool isV0(const Vertex* vtx, const VertexSelectorConfig& cfg, const Vertex* primary) {
151  // two-track selection
152  if ( vtx->getTracks().size() !=2 )
153  return false;
154 
155  // charge selection
156  if ( vtx->getTracks()[0]->getCharge() * vtx->getTracks()[1]->getCharge() != -1 )
157  return false;
158 
159  // obtain momentum vectors at the V0 vertex (IMPORTANT)
160  TVector3 mom1 = vtx->getTracks()[0]->momentumAtVertex(vtx);
161  TVector3 mom2 = vtx->getTracks()[1]->momentumAtVertex(vtx);
162 
163  // V0 direction selection
164  TVector3 posip;
165  if (primary) posip = primary->getPos();
166  double v0dirdot = (mom1+mom2).Unit().Dot( (vtx->getPos()-posip).Unit() );
167 
168  // K0 selection
169  if (cfg.rejectk0) {
170  TLorentzVector lvec1;
171  TLorentzVector lvec2;
172  lvec1.SetVectM( mom1, 0.1396 );
173  lvec2.SetVectM( mom2, 0.1396 );
174  double mass = (lvec1+lvec2).M();
175  double k0mass = .498;
176  if (v0dirdot > cfg.k0dirdot && fabs(mass-k0mass) < cfg.k0width )
177  return true;
178  }
179 
180  // Lambda selection
181  if (cfg.rejectl0) {
182  // assume higher momentum particle is the proton
183  TLorentzVector protonForLambda;
184  TLorentzVector pionForLambda;
185  if (mom1.Mag() > mom2.Mag()) {
186  protonForLambda.SetVectM( mom1, 0.9383 );
187  pionForLambda.SetVectM( mom2, 0.1396 );
188  } else {
189  protonForLambda.SetVectM( mom2, 0.9383 );
190  pionForLambda.SetVectM( mom1, 0.1396 );
191  }
192 
193  double lambdaMass = (protonForLambda+pionForLambda).M();
194  if (v0dirdot > cfg.l0dirdot && fabs(lambdaMass - 1.1157) < cfg.l0width)
195  return true;
196  }
197 
198  // check photon conversion; photon mass corrected geometrically
199  if (cfg.rejectconv) {
200  double ang1 = atan( vtx->getTracks()[0]->getTanLambda());
201  double ang2 = atan( vtx->getTracks()[1]->getTanLambda());
202  double convmassCorr = sqrt( vtx->getTracks()[0]->Vect().Mag()*vtx->getTracks()[1]->Vect().Mag()*(1-cos(ang1-ang2)) );
203  if (v0dirdot > cfg.convdirdot && vtx->getPos().Mag()> cfg.convdist && convmassCorr < cfg.convmass)
204  return true;
205  }
206 
207  return false;
208  }
209 
210  // true if passes cut
211  bool passesCut(const Vertex* vtx, const VertexSelectorConfig& cfg, const Vertex* primary=0) {
212 
213  // position selection
214  if (cfg.rejectdist) {
215  bool cond = (vtx->getPos().Mag() > cfg.minpos && vtx->getPos().Mag() < cfg.maxpos);
216 
217  if (cfg.rejectdistnegative)cond = !cond;
218  if (cfg.rejectdistor && cond)return true;
219  if (!cfg.rejectdistor && !cond)return false;
220  }
221 
222  // v0 selection
223  if (cfg.rejectk0 || cfg.rejectl0 || cfg.rejectconv) {
224  if (isV0(vtx,cfg,primary))return false;
225 // if (isV0(vtx,cfg,primary)) {cout << "rejected by v0 selection" << endl;return false;}
226  }
227 
228  return true;
229  }
230 
231  //c-tor / d-tor
234 };
235 }
236 
237 #endif //TrackSelector_h
vector< const Track * >::const_iterator TrackVecIte
Definition: lcfiplus.h:83
double minpos
Definition: VertexSelector.h:18
double convmass
Definition: VertexSelector.h:32
Definition: lcfiplus.h:771
void setV0Loose()
Definition: VertexSelector.h:56
bool rejectdistnegative
Definition: VertexSelector.h:20
bool isV0(const Vertex *vtx, const VertexSelectorConfig &cfg, const Vertex *primary)
Definition: VertexSelector.h:150
void setNoV0Cut()
Definition: VertexSelector.h:68
double l0width
Definition: VertexSelector.h:28
bool passesCut(const Vertex *vtx, const VertexSelectorConfig &cfg, const Vertex *primary=0)
Definition: VertexSelector.h:211
VertexSelectorConfig()
Definition: VertexSelector.h:36
double maxpos
Definition: VertexSelector.h:19
TVector3 getPos() const
Definition: lcfiplus.h:817
~VertexSelector()
Definition: VertexSelector.h:233
bool rejectk0
Definition: VertexSelector.h:23
Definition: VertexSelector.h:76
Definition: VertexSelector.h:15
double convdirdot
Definition: VertexSelector.h:34
VertexSelector()
Definition: VertexSelector.h:232
bool rejectl0
Definition: VertexSelector.h:27
bool rejectdist
Definition: VertexSelector.h:17
void setV0Tight()
Definition: VertexSelector.h:44
double l0dirdot
Definition: VertexSelector.h:29
bool rejectconv
Definition: VertexSelector.h:31
double k0width
Definition: VertexSelector.h:24
double k0dirdot
Definition: VertexSelector.h:25
bool rejectdistor
Definition: VertexSelector.h:21
const vector< const Track * > & getTracks() const
Definition: lcfiplus.h:823
double convdist
Definition: VertexSelector.h:33
const vector< const Vertex * > VertexVec
Definition: lcfiplus.h:79