StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
franksPairCut.cxx
1 //#include "StHbtMaker/Infrastructure/StHbtPairResonanceInfo.hh"
2 #include "StHbtMaker/Infrastructure/StHbtTypes.hh"
3 #include "StHbtMaker/Cut/franksPairCut.h"
4 #include "SystemOfUnits.h"
5 #include <string>
6 #include <cstdio>
7 
8 #ifdef __ROOT__
9 ClassImp(franksPairCut)
10 #endif
11 
12 #define __2POWER16__ 65536
13 
14 //__________________
15  franksPairCut::franksPairCut() /* : mResonanceInfoOn(false) */ {
16  mQuality[0] = -1.; mQuality[1] = +1.;
17  mKt[0] = -1e9; mKt[1]= +1e9;
18  mPt[0] = -1e9; mPt[1]= +1e9;
19  mOpeningAngle[0] = -1e9; mOpeningAngle[1]= +1e9;
20  mRapidity[0] = -1e9; mRapidity[1]= +1e9;
21  mEta[0] = -1e9; mEta[1]= +1e9;
22  mQinv[0] = -1e9; mQinv[1]= +1e9;
23  mEntranceSeparation[0] = -1e9; mEntranceSeparation[1]= +1e9;
24  mDecayLength[0] = -1e9; mDecayLength[1]= 1e9;
25  mAngleToPrimaryVertex[0] = -1e9; mAngleToPrimaryVertex[1]= 1e9;
26  mDcaToPrimaryVertex[0] = -1e9; mDcaToPrimaryVertex[1]= 1e9;
27  mDcaOfDaughters[0] = -1e9; mDcaOfDaughters[1]= 1e9;
28  mNPairsPassed = mNPairsFailed = 0;
29  mIdenticalMother =0;
30 }
31 //__________________
32 franksPairCut::franksPairCut(const franksPairCut& c) : StHbtPairCut(c) /* : , mResonanceInfoOn(false) */ {
33 #ifdef STHBTDEBUG
34  cout << " franksPairCut::franksPairCut(const franksPairCut& c) " << endl;
35 #endif
36  mNPairsPassed = mNPairsFailed = 0;
37  mPrimaryVertex = c.mPrimaryVertex;
38  mQuality[0] = c.mQuality[0];
39  mQuality[1] = c.mQuality[1];
40  mKt[0] = c.mKt[0];
41  mKt[1] = c.mKt[1];
42  mPt[0] = c.mPt[0];
43  mPt[1] = c.mPt[1];
44  mOpeningAngle[0] = c.mOpeningAngle[0];
45  mOpeningAngle[1] = c.mOpeningAngle[1];
46  mQinv[0] = c.mQinv[0];
47  mQinv[1] = c.mQinv[1];
48  mRapidity[0] = c.mRapidity[0];
49  mRapidity[1] = c.mRapidity[1];
50  mEta[0] = c.mEta[0];
51  mEta[1] = c.mEta[1];
52  mEntranceSeparation[0] = c.mEntranceSeparation[0];
53  mEntranceSeparation[1] = c.mEntranceSeparation[1];
54  mDecayLength[0] = c.mDecayLength[0];
55  mDecayLength[1] = c.mDecayLength[1];
56  mAngleToPrimaryVertex[0] = c.mAngleToPrimaryVertex[0];
57  mAngleToPrimaryVertex[1] = c.mAngleToPrimaryVertex[1];
58  mDcaToPrimaryVertex[0] = c.mDcaToPrimaryVertex[0];
59  mDcaToPrimaryVertex[1] = c.mDcaToPrimaryVertex[1];
60  mDcaOfDaughters[0] = c.mDcaOfDaughters[0];
61  mDcaOfDaughters[1] = c.mDcaOfDaughters[1];
62  mIdenticalMother = c.mIdenticalMother;
63 }
64 //__________________
65 //franksPairCut::~franksPairCut(){
66 // /* no-op */
67 //}
68 //__________________
69 bool franksPairCut::Pass(const StHbtPair* pair){
70  if ( !(pair->quality() >= mQuality[0] && pair->quality() <= mQuality[1] ) )
71  return leave(false);
72  if ( !(pair->kT() >= mKt[0] && pair->kT() <= mKt[1] ) )
73  return leave(false);
74  if ( !(pair->fourMomentumSum().perp() >= mPt[0] && pair->fourMomentumSum().perp() <= mPt[1] ) )
75  return leave(false);
76  if ( !(pair->OpeningAngle() >= mOpeningAngle[0] && pair->OpeningAngle() <= mOpeningAngle[1] ) )
77  return leave(false);
78  if ( !(pair->fourMomentumSum().rapidity() >= mRapidity[0] && pair->fourMomentumSum().rapidity() <= mRapidity[1]) )
79  return leave(false);
80  if ( !(pair->fourMomentumSum().pseudoRapidity() >= mEta[0] && pair->fourMomentumSum().pseudoRapidity() <= mEta[1] ) )
81  return leave(false);
82  if ( !(fabs(pair->qInv()) >= mQinv[0] && fabs(pair->qInv()) <= mQinv[1]) )
83  return leave(false);
84  if ( !(pair->NominalTpcEntranceSeparation() >= mEntranceSeparation[0] && pair->NominalTpcEntranceSeparation() <= mEntranceSeparation[1] )
85  ) return leave(false);
86 
87  /*if (mResonanceInfoOn) {
88  //cout << " fix this " << endl;
89  pair->CalculateResonanceInfo(&mPrimaryVertex, 0.25*tesla);
90  //cout << " mDecayLength " << pair->ResonanceInfo()->mDecayLength << endl;
91  if ( !(pair->ResonanceInfo()->mDecayLength >= mDecayLength[0] &&
92  pair->ResonanceInfo()->mDecayLength <= mDecayLength[1]) )
93  return leave(false);
94  //cout << " mDcaOfDaughters " << pair->ResonanceInfo()->mDcaOfDaughters << endl;
95  if ( !(pair->ResonanceInfo()->mDcaOfDaughters >= mDcaOfDaughters[0] &&
96  pair->ResonanceInfo()->mDcaOfDaughters <= mDcaOfDaughters[1]) )
97  return leave(false);
98 
99  //cout << "mAngleToPrimaryVertex " << pair->ResonanceInfo()->mAngleToPrimaryVertex << endl;
100  if ( !(pair->ResonanceInfo()->mAngleToPrimaryVertex >= mAngleToPrimaryVertex[0] &&
101  pair->ResonanceInfo()->mAngleToPrimaryVertex <= mAngleToPrimaryVertex[1]) )
102  return leave(false);
103 
104  //cout << "mDcaToPrimaryVertex " << pair->ResonanceInfo()->mDcaToPrimaryVertex << endl;
105  if ( !(pair->ResonanceInfo()->mDcaToPrimaryVertex >= mDcaToPrimaryVertex[0] &&
106  pair->ResonanceInfo()->mDcaToPrimaryVertex <= mDcaToPrimaryVertex[1]) )
107  return leave(false);
108 
109  }*/
110 
111  if (mIdenticalMother)
112  if( !((int)(pair->track1()->TrackId()/__2POWER16__) != (int)(pair->track2()->TrackId()/__2POWER16__) ) )
113  return leave(false);
114 
115  return leave(true);
116 }
117 
118 //__________________
119 bool franksPairCut::leave(bool b){
120  b ? mNPairsPassed++ : mNPairsFailed++;
121  return b;
122 }
123 
124 //__________________
125 StHbtString franksPairCut::Report(){
126  string Stemp = "Franks Pair Cut have to finish cut and report; this is just a test - \n";
127  char Ctemp[100];
128  sprintf(Ctemp,"Number of pairs which passed:\t%ld Number which failed:\t%ld\n",mNPairsPassed,mNPairsFailed);
129  Stemp += Ctemp;
130  sprintf(Ctemp,"quality: %f -- %f\n",mQuality[0],mQuality[1]);
131  Stemp += Ctemp;
132  sprintf(Ctemp,"kT: %f -- %f\n",mKt[0],mKt[1]);
133  Stemp += Ctemp;
134  sprintf(Ctemp,"pT: %f -- %f\n",mPt[0],mPt[1]);
135  Stemp += Ctemp;
136  sprintf(Ctemp,"opening angle: %f -- %f\n",mOpeningAngle[0],mOpeningAngle[1]);
137  Stemp += Ctemp;
138  sprintf(Ctemp,"pair qinv: %f -- %f\n",mQinv[0],mQinv[1]);
139  Stemp += Ctemp;
140  sprintf(Ctemp,"pair rapidity: %f -- %f\n",mRapidity[0],mRapidity[1]);
141  Stemp += Ctemp;
142  sprintf(Ctemp,"pair eta: %f -- %f\n",mEta[0],mEta[1]);
143  Stemp += Ctemp;
144  sprintf(Ctemp,"EntranceSeparation: %f -- %f\n",mEntranceSeparation[0],mEntranceSeparation[1]);
145  Stemp += Ctemp;
146  sprintf(Ctemp,"DecayLength: %f -- %f\n",mDecayLength[0],mDecayLength[1]);
147  Stemp += Ctemp;
148  sprintf(Ctemp,"AngleToPrimaryVertex: %f -- %f\n",mAngleToPrimaryVertex[0],mAngleToPrimaryVertex[1]);
149  Stemp += Ctemp;
150  sprintf(Ctemp,"DcaToPrimaryVertex: %f -- %f\n",mDcaToPrimaryVertex[0],mDcaToPrimaryVertex[1]);
151  Stemp += Ctemp;
152  StHbtString returnThis = Stemp;
153  return returnThis;
154 }
155