StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
franksTrackCut.cxx
1 /***************************************************************************
2  *
3  *
4  *
5  * Author: Frank Laue, Ohio State, laue@mps.ohio-state.edu
6  ***************************************************************************
7  *
8  * Description: part of STAR HBT Framework: StHbtMaker package
9  * a simple particle cut that selects on phasespace, #hits, DCA, and PID
10  *
11  ***************************************************************************
12  *
13  *
14  **************************************************************************/
15 
16 #include "StHbtMaker/Cut/franksTrackCut.h"
17 #include <cstdio>
18 
19 #ifdef __ROOT__
20  ClassImp(franksTrackCut)
21 #endif
22 
23  //#define STHBTDEBUG
24 
25 franksTrackCut::franksTrackCut(){
26  mNTracksPassed = mNTracksFailed = 0;
27  mPidProbElectron[0] = -1e9; mPidProbElectron[1] = +1e9;
28  mPidProbPion[0] = -1e9; mPidProbPion[1] = +1e9;
29  mPidProbKaon[0] = -1e9; mPidProbKaon[1] = +1e9;
30  mPidProbProton[0] = -1e9; mPidProbProton[1] = +1e9;
31  mNSigmaElectron[0] = -1e9; mNSigmaElectron[1] = +1e9;
32  mNSigmaPion[0] = -1e9; mNSigmaPion[1] = +1e9;
33  mNSigmaKaon[0] = -1e9; mNSigmaKaon[1] = +1e9;
34  mNSigmaProton[0] = -1e9; mNSigmaProton[1] = +1e9;
35  mNSigmaAntiElectron[0] = 0.; mNSigmaAntiElectron[1] = 0.;
36  mNSigmaAntiPion[0] = 0.; mNSigmaAntiPion[1] = 0.;
37  mNSigmaAntiKaon[0] = 0.; mNSigmaAntiKaon[1] = 0.;
38  mNSigmaAntiProton[0] = 0.; mNSigmaAntiProton[1] = 0.;
39  mEta[0] = -1e9; mEta[1] = +1e9;
40  mRapidity[0] = -1e9; mRapidity[1] = +1e9;
41  mP[0] = -1e9; mP[1] = +1e9;
42  mPt[0] = -1e9; mPt[1] = +1e9;
43  mPx[0] = -1e9; mPx[1] = +1e9;
44  mPy[0] = -1e9; mPy[1] = +1e9;
45  mPz[0] = -1e9; mPz[1] = +1e9;
46  mDCA[0] = -1e9; mDCA[1] = +1e9;
47  mDCAGlobal[0] = -1e9; mDCAGlobal[1] = +1e9;
48  mNHits[0] = 0; mNHits[1] = 60;
49  mNdEdxHits[0] = 0; mNdEdxHits[1] = 60;
50 }
51 //------------------------------
52 franksTrackCut::franksTrackCut(franksTrackCut& c) : StHbtTrackCut(c) {
53  mNTracksPassed = mNTracksFailed = 0;
54  mCharge = c.mCharge;
55  mPidProbElectron[0] = c.mPidProbElectron[0];
56  mPidProbElectron[1] = c.mPidProbElectron[1];
57  mPidProbPion[0] = c.mPidProbPion[0];
58  mPidProbPion[1] = c.mPidProbPion[1];
59  mPidProbKaon[0] = c.mPidProbKaon[0];
60  mPidProbKaon[1] = c.mPidProbKaon[1];
61  mPidProbProton[0] = c.mPidProbProton[0];
62  mPidProbProton[1] = c.mPidProbProton[1];
63  mNSigmaElectron[0] = c.mNSigmaElectron[0];
64  mNSigmaElectron[1] = c.mNSigmaElectron[1];
65  mNSigmaPion[0] = c.mNSigmaPion[0];
66  mNSigmaPion[1] = c.mNSigmaPion[1];
67  mNSigmaKaon[0] = c.mNSigmaKaon[0];
68  mNSigmaKaon[1] = c.mNSigmaKaon[1];
69  mNSigmaProton[0] = c.mNSigmaProton[0];
70  mNSigmaProton[1] = c.mNSigmaProton[1];
71  mNSigmaAntiElectron[0] = c.mNSigmaAntiElectron[0];
72  mNSigmaAntiElectron[1] = c.mNSigmaAntiElectron[1];
73  mNSigmaAntiPion[0] = c.mNSigmaAntiPion[0];
74  mNSigmaAntiPion[1] = c.mNSigmaAntiPion[1];
75  mNSigmaAntiKaon[0] = c.mNSigmaAntiKaon[0];
76  mNSigmaAntiKaon[1] = c.mNSigmaAntiKaon[1];
77  mNSigmaAntiProton[0] = c.mNSigmaAntiProton[0];
78  mNSigmaAntiProton[1] = c.mNSigmaAntiProton[1];
79  mNHits[0] = c.mNHits[0];
80  mNHits[1] = c.mNHits[1];
81  mNdEdxHits[0] = c.mNdEdxHits[0];
82  mNdEdxHits[1] = c.mNdEdxHits[1];
83  mP[0] = c.mP[0]; mP[1] = c.mP[1];
84  mPt[0] = c.mPt[0]; mPt[1] = c.mPt[1];
85  mPx[0] = c.mPx[0]; mPx[1] = c.mPx[1];
86  mPy[0] = c.mPy[0]; mPy[1] = c.mPy[1];
87  mPz[0] = c.mPz[0]; mPz[1] = c.mPz[1];
88  mRapidity[0] = c.mRapidity[0];
89  mRapidity[1] = c.mRapidity[1];
90  mEta[0] = c.mEta[0];
91  mEta[1] = c.mEta[1];
92  mDCA[0] = c.mDCA[0];
93  mDCA[1] = c.mDCA[1];
94  mDCAGlobal[0] = c.mDCAGlobal[0];
95  mDCAGlobal[1] = c.mDCAGlobal[1];
96  mNTracksPassed=0;
97  mNTracksFailed=0;
98 #ifdef STHBTDEBUG
99  cout << " franksTrackCut::franksTrackCut(franksTrackCut& c) " << endl;
100 #endif
101 }
102 
103 //------------------------------
104 franksTrackCut::~franksTrackCut(){
105 }
106 
107 //------------------------------
108 bool franksTrackCut::Pass(const StHbtTrack* track){
109  // cout << " *** franksTrackCut::Pass(const StHbtTrack* track) " << endl;
110 
111  float TEnergy = ::sqrt(track->P().mag2()+mMass*mMass);
112  float TRapidity = 0.5*::log((TEnergy+track->P().z())/
113  (TEnergy-track->P().z()));
114 
115 #ifdef STHBTDEBUG
116  cout <<
117  track->NSigmaElectron() << " " <<
118  track->NSigmaPion() << " " <<
119  track->NSigmaKaon() << " " <<
120  track->NSigmaProton() << " " <<
121  track->DCAxy() << " " <<
122  track->DCAxyGlobal() << " " <<
123  track->NHits() << " " <<
124  track->NHitsDedx() << " " <<
125  track->P().mag() << " " <<
126  track->Pt() << " " <<
127  TRapidity << " " <<
128  track->Charge() << " ";
129 #endif
130 
131  bool goodPID = (
132  (track->PidProbElectron() >= mPidProbElectron[0]) &&
133  (track->PidProbElectron() <= mPidProbElectron[1]) &&
134  (track->PidProbPion() >= mPidProbPion[0]) &&
135  (track->PidProbPion() <= mPidProbPion[1]) &&
136  (track->PidProbKaon() >= mPidProbKaon[0]) &&
137  (track->PidProbKaon() <= mPidProbKaon[1]) &&
138  (track->PidProbProton() >= mPidProbProton[0]) &&
139  (track->PidProbProton() <= mPidProbProton[1]) &&
140  (track->NSigmaElectron() >= mNSigmaElectron[0]) &&
141  (track->NSigmaElectron() <= mNSigmaElectron[1]) &&
142  (track->NSigmaPion() >= mNSigmaPion[0]) &&
143  (track->NSigmaPion() <= mNSigmaPion[1]) &&
144  (track->NSigmaKaon() >= mNSigmaKaon[0]) &&
145  (track->NSigmaKaon() <= mNSigmaKaon[1]) &&
146  (track->NSigmaProton() >= mNSigmaProton[0]) &&
147  (track->NSigmaProton() <= mNSigmaProton[1]) &&
148  !( (track->NSigmaElectron() > mNSigmaAntiElectron[0]) &&
149  (track->NSigmaElectron() < mNSigmaAntiElectron[1]) ) &&
150  !( (track->NSigmaPion() > mNSigmaAntiPion[0]) &&
151  (track->NSigmaPion() < mNSigmaAntiPion[1]) ) &&
152  !( (track->NSigmaKaon() > mNSigmaAntiKaon[0]) &&
153  (track->NSigmaKaon() < mNSigmaAntiKaon[1]) ) &&
154  !( (track->NSigmaProton() > mNSigmaAntiProton[0]) &&
155  (track->NSigmaProton() < mNSigmaAntiProton[1]) ) &&
156  (track->Charge() == mCharge || mCharge==0 )
157  );
158 
159 #ifdef STHBTDEBUG
160  cout << mNSigmaElectron[0] << " << " << track->NSigmaElectron() << " << " << mNSigmaElectron[1] << " ";
161  cout << (track->NSigmaElectron() >= mNSigmaElectron[0]) << (track->NSigmaElectron() <= mNSigmaElectron[1]) << endl;
162 #endif
163 
164  bool goodTrack=( true &&
165  (track->DCAxy() >= mDCA[0]) &&
166  (track->DCAxy() <= mDCA[1]) &&
167  (track->DCAxyGlobal() >= mDCAGlobal[0]) &&
168  (track->DCAxyGlobal() <= mDCAGlobal[1]) &&
169  (track->NHits() >= mNHits[0]) &&
170  (track->NHits() <= mNHits[1]) &&
171  (track->NHitsDedx() >= mNdEdxHits[0]) &&
172  (track->NHitsDedx() <= mNdEdxHits[1]) &&
173  (track->P().mag() >= mP[0]) &&
174  (track->P().mag() <= mP[1]) &&
175  (track->Pt() >= mPt[0]) &&
176  (track->Pt() <= mPt[1]) &&
177  (track->P().x() >= mPx[0]) &&
178  (track->P().x() <= mPx[1]) &&
179  (track->P().y() >= mPy[0]) &&
180  (track->P().y() <= mPy[1]) &&
181  (track->P().z() >= mPz[0]) &&
182  (track->P().z() <= mPz[1]) &&
183  (track->P().pseudoRapidity() >= mEta[0]) &&
184  (track->P().pseudoRapidity() <= mEta[1]) &&
185  (TRapidity >= mRapidity[0]) &&
186  (TRapidity <= mRapidity[1])
187 );
188 
189 
190 #ifdef STHBTDEBUG
191  cout << " goodPID=" << goodPID << " ";
192  cout << " goodTrack=" << goodTrack << " ";
193  cout << endl;
194 #endif
195 
196  (goodTrack && goodPID) ? mNTracksPassed++ : mNTracksFailed++;
197 
198  return (goodPID && goodTrack);
199 }
200 //------------------------------
201 StHbtString franksTrackCut::Report(){
202  string Stemp;
203  char Ctemp[100];
204  sprintf(Ctemp,"\nParticle mass:\t%E",mMass);
205  Stemp+=Ctemp;
206  sprintf(Ctemp,"\nParticle charge:\t%d",mCharge);
207  Stemp+=Ctemp;
208  sprintf(Ctemp,"\nParticle Nsigma from pion:\t%E - %E",mNSigmaPion[0],mNSigmaPion[1]);
209  Stemp+=Ctemp;
210  sprintf(Ctemp,"\nParticle Nsigma from kaon:\t%E - %E",mNSigmaKaon[0],mNSigmaKaon[1]);
211  Stemp+=Ctemp;
212  sprintf(Ctemp,"\nParticle Nsigma from proton:\t%E - %E",mNSigmaProton[0],mNSigmaProton[1]);
213  Stemp+=Ctemp;
214  sprintf(Ctemp,"\nParticle PidProb from electron:\t%E - %E",mPidProbElectron[0],mPidProbElectron[1]);
215  Stemp+=Ctemp;
216  sprintf(Ctemp,"\nParticle PidProb from pion:\t%E - %E",mPidProbPion[0],mPidProbPion[1]);
217  Stemp+=Ctemp;
218  sprintf(Ctemp,"\nParticle PidProb from kaon:\t%E - %E",mPidProbKaon[0],mPidProbKaon[1]);
219  Stemp+=Ctemp;
220  sprintf(Ctemp,"\nParticle PidProb from proton:\t%E - %E",mPidProbProton[0],mPidProbProton[1]);
221  Stemp+=Ctemp;
222  sprintf(Ctemp,"\nParticle PidProb from electron:\t%E - %E",mPidProbElectron[0],mPidProbElectron[1]);
223  Stemp+=Ctemp;
224  sprintf(Ctemp,"\nParticle #hits:\t%d - %d",mNHits[0],mNHits[1]);
225  Stemp+=Ctemp;
226  sprintf(Ctemp,"\nParticle dEdx #hits:\t%d - %d",mNdEdxHits[0],mNdEdxHits[1]);
227  Stemp+=Ctemp;
228  sprintf(Ctemp,"\nParticle p:\t%E - %E",mP[0],mP[1]);
229  Stemp+=Ctemp;
230  sprintf(Ctemp,"\nParticle pT:\t%E - %E",mPt[0],mPt[1]);
231  Stemp+=Ctemp;
232  sprintf(Ctemp,"\nParticle rapidity:\t%E - %E",mRapidity[0],mRapidity[1]);
233  Stemp+=Ctemp;
234  sprintf(Ctemp,"\nParticle mEta:\t%E - %E",mEta[0],mEta[1]);
235  Stemp+=Ctemp;
236  sprintf(Ctemp,"\nParticle DCA:\t%E - %E",mDCA[0],mDCA[1]);
237  Stemp+=Ctemp;
238  sprintf(Ctemp,"\nParticle DCAGlobal:\t%E - %E",mDCAGlobal[0],mDCAGlobal[1]);
239  Stemp+=Ctemp;
240  sprintf(Ctemp,"\nNumber of tracks which passed:\t%ld Number which failed:\t%ld",mNTracksPassed,mNTracksFailed);
241  Stemp += Ctemp;
242  StHbtString returnThis = Stemp;
243  return returnThis;
244 }
245 
246 
247 std::ostringstream* franksTrackCut::finalReport() const{
248  std::ostringstream* tFinalReport = new std::ostringstream;
249  (*tFinalReport) << "_____ Track Cut _____ " << endl
250  << "Charge = " << mCharge << endl
251  << mNSigmaPion[0] << " < Sigma pion < "
252  << mNSigmaPion[1] << endl
253  << mNSigmaKaon[0] << " < Sigma Kaon < "
254  << mNSigmaKaon[1] << endl
255  << mNSigmaProton[0] << " < Sigma proton < "
256  << mNSigmaProton[1] << endl
257  << mNHits[0] << " < NHits < " << mNHits[1] << endl
258  << mNdEdxHits[0] << " < NHits < " << mNdEdxHits[1] << endl
259  << mPt[0] << " < pT < " << mPt[1] << endl
260  << mP[0] << " < p < " << mP[1] << endl
261  << mRapidity[0] << " < Y < " << mRapidity[1] << endl
262  << mDCA[0] << " < dca < " << mDCA[1] << endl
263  << " >>> N Tracks passed " << mNTracksPassed << endl
264  << " >>> N Tracks failed " << mNTracksFailed << endl
265  << ends;
266  return tFinalReport;
267 }
268