StRoot: StHbtMaker/Cut/franksPairCut.cxx Source File

00001 //#include "StHbtMaker/Infrastructure/StHbtPairResonanceInfo.hh"
00002 #include "StHbtMaker/Infrastructure/StHbtTypes.hh"
00003 #include "StHbtMaker/Cut/franksPairCut.h"
00004 #include "SystemOfUnits.h"
00005 #include <string>
00006 #include <cstdio>
00007 
00008 #ifdef __ROOT__
00009 ClassImp(franksPairCut)
00010 #endif
00011 
00012 #define __2POWER16__ 65536
00013 
00014 //__________________
00015   franksPairCut::franksPairCut() /* : mResonanceInfoOn(false) */ {
00016   mQuality[0] = -1.; mQuality[1] = +1.;
00017   mKt[0] = -1e9; mKt[1]= +1e9;
00018   mPt[0] = -1e9; mPt[1]= +1e9;
00019   mOpeningAngle[0] = -1e9; mOpeningAngle[1]= +1e9;
00020   mRapidity[0] = -1e9; mRapidity[1]= +1e9;
00021   mEta[0] = -1e9; mEta[1]= +1e9;
00022   mQinv[0] = -1e9; mQinv[1]= +1e9;
00023   mEntranceSeparation[0] = -1e9; mEntranceSeparation[1]= +1e9;
00024   mDecayLength[0] = -1e9; mDecayLength[1]= 1e9;
00025   mAngleToPrimaryVertex[0] = -1e9; mAngleToPrimaryVertex[1]= 1e9;
00026   mDcaToPrimaryVertex[0] = -1e9; mDcaToPrimaryVertex[1]= 1e9;
00027   mDcaOfDaughters[0] = -1e9; mDcaOfDaughters[1]= 1e9;
00028   mNPairsPassed = mNPairsFailed = 0;
00029   mIdenticalMother =0;
00030 }
00031 //__________________
00032 franksPairCut::franksPairCut(const franksPairCut& c) : StHbtPairCut(c) /* :  ,  mResonanceInfoOn(false) */ {
00033 #ifdef STHBTDEBUG
00034   cout << " franksPairCut::franksPairCut(const franksPairCut& c) " << endl;
00035 #endif
00036   mNPairsPassed = mNPairsFailed = 0;
00037   mPrimaryVertex = c.mPrimaryVertex;
00038   mQuality[0] = c.mQuality[0];
00039   mQuality[1] = c.mQuality[1];
00040   mKt[0] = c.mKt[0];
00041   mKt[1] = c.mKt[1];
00042   mPt[0] = c.mPt[0];
00043   mPt[1] = c.mPt[1];
00044   mOpeningAngle[0] = c.mOpeningAngle[0];
00045   mOpeningAngle[1] = c.mOpeningAngle[1];
00046   mQinv[0] = c.mQinv[0];
00047   mQinv[1] = c.mQinv[1];
00048   mRapidity[0] = c.mRapidity[0];
00049   mRapidity[1] = c.mRapidity[1];
00050   mEta[0] = c.mEta[0];
00051   mEta[1] = c.mEta[1];
00052   mEntranceSeparation[0] = c.mEntranceSeparation[0];
00053   mEntranceSeparation[1] = c.mEntranceSeparation[1];
00054   mDecayLength[0] = c.mDecayLength[0];
00055   mDecayLength[1] = c.mDecayLength[1];
00056   mAngleToPrimaryVertex[0] = c.mAngleToPrimaryVertex[0];
00057   mAngleToPrimaryVertex[1] = c.mAngleToPrimaryVertex[1];
00058   mDcaToPrimaryVertex[0] = c.mDcaToPrimaryVertex[0];
00059   mDcaToPrimaryVertex[1] = c.mDcaToPrimaryVertex[1];
00060   mDcaOfDaughters[0] = c.mDcaOfDaughters[0];
00061   mDcaOfDaughters[1] = c.mDcaOfDaughters[1];
00062   mIdenticalMother = c.mIdenticalMother;
00063 }
00064 //__________________
00065 //franksPairCut::~franksPairCut(){
00066 //  /* no-op */
00067 //}  
00068 //__________________
00069 bool franksPairCut::Pass(const StHbtPair* pair){
00070   if ( !(pair->quality() >= mQuality[0]  && pair->quality() <= mQuality[1] ) ) 
00071     return leave(false);
00072   if ( !(pair->kT() >= mKt[0]  && pair->kT() <= mKt[1] ) ) 
00073     return leave(false);
00074   if ( !(pair->fourMomentumSum().perp() >= mPt[0]  && pair->fourMomentumSum().perp() <= mPt[1] ) ) 
00075     return leave(false);
00076   if ( !(pair->OpeningAngle() >= mOpeningAngle[0]  && pair->OpeningAngle() <= mOpeningAngle[1] )  ) 
00077     return leave(false);
00078   if ( !(pair->fourMomentumSum().rapidity() >= mRapidity[0]  && pair->fourMomentumSum().rapidity() <= mRapidity[1]) ) 
00079     return leave(false);
00080   if ( !(pair->fourMomentumSum().pseudoRapidity() >= mEta[0]  && pair->fourMomentumSum().pseudoRapidity() <= mEta[1] )  ) 
00081     return leave(false);
00082   if ( !(fabs(pair->qInv()) >= mQinv[0]  && fabs(pair->qInv()) <= mQinv[1])  ) 
00083     return leave(false);
00084   if ( !(pair->NominalTpcEntranceSeparation() >= mEntranceSeparation[0]  && pair->NominalTpcEntranceSeparation() <= mEntranceSeparation[1] )
00085        ) return leave(false);
00086   
00087   /*if (mResonanceInfoOn) {
00088     //cout << " fix this " << endl; 
00089     pair->CalculateResonanceInfo(&mPrimaryVertex, 0.25*tesla);
00090     //cout << " mDecayLength "  << pair->ResonanceInfo()->mDecayLength << endl;
00091     if ( !(pair->ResonanceInfo()->mDecayLength >= mDecayLength[0]  && 
00092            pair->ResonanceInfo()->mDecayLength <= mDecayLength[1]) ) 
00093       return leave(false);
00094     //cout << " mDcaOfDaughters "  << pair->ResonanceInfo()->mDcaOfDaughters << endl;
00095     if ( !(pair->ResonanceInfo()->mDcaOfDaughters >= mDcaOfDaughters[0]  && 
00096            pair->ResonanceInfo()->mDcaOfDaughters <= mDcaOfDaughters[1]) ) 
00097       return leave(false);
00098     
00099     //cout << "mAngleToPrimaryVertex " << pair->ResonanceInfo()->mAngleToPrimaryVertex << endl;
00100     if ( !(pair->ResonanceInfo()->mAngleToPrimaryVertex >= mAngleToPrimaryVertex[0]  && 
00101            pair->ResonanceInfo()->mAngleToPrimaryVertex <= mAngleToPrimaryVertex[1]) ) 
00102       return leave(false);
00103     
00104     //cout << "mDcaToPrimaryVertex " << pair->ResonanceInfo()->mDcaToPrimaryVertex << endl;
00105     if ( !(pair->ResonanceInfo()->mDcaToPrimaryVertex >= mDcaToPrimaryVertex[0]  && 
00106            pair->ResonanceInfo()->mDcaToPrimaryVertex <= mDcaToPrimaryVertex[1]) ) 
00107       return leave(false);
00108     
00109            }*/
00110 
00111   if (mIdenticalMother)
00112     if(  !((int)(pair->track1()->TrackId()/__2POWER16__) != (int)(pair->track2()->TrackId()/__2POWER16__) ) ) 
00113       return leave(false);
00114   
00115   return leave(true);
00116 }
00117 
00118 //__________________
00119 bool franksPairCut::leave(bool b){
00120    b ? mNPairsPassed++ : mNPairsFailed++;
00121    return b;
00122 }
00123 
00124 //__________________
00125 StHbtString franksPairCut::Report(){
00126   string Stemp = "Franks Pair Cut have to finish cut and report; this is just a test - \n";
00127   char Ctemp[100];
00128   sprintf(Ctemp,"Number of pairs which passed:\t%ld  Number which failed:\t%ld\n",mNPairsPassed,mNPairsFailed);
00129   Stemp += Ctemp;
00130   sprintf(Ctemp,"quality: %f  -- %f\n",mQuality[0],mQuality[1]);
00131   Stemp += Ctemp;
00132   sprintf(Ctemp,"kT: %f  -- %f\n",mKt[0],mKt[1]);
00133   Stemp += Ctemp;
00134   sprintf(Ctemp,"pT: %f  -- %f\n",mPt[0],mPt[1]);
00135   Stemp += Ctemp;
00136   sprintf(Ctemp,"opening angle: %f  -- %f\n",mOpeningAngle[0],mOpeningAngle[1]);
00137   Stemp += Ctemp;
00138   sprintf(Ctemp,"pair qinv: %f  -- %f\n",mQinv[0],mQinv[1]);
00139   Stemp += Ctemp;
00140   sprintf(Ctemp,"pair rapidity: %f  -- %f\n",mRapidity[0],mRapidity[1]);
00141   Stemp += Ctemp;
00142   sprintf(Ctemp,"pair eta: %f  -- %f\n",mEta[0],mEta[1]);
00143   Stemp += Ctemp;
00144   sprintf(Ctemp,"EntranceSeparation: %f  -- %f\n",mEntranceSeparation[0],mEntranceSeparation[1]);
00145   Stemp += Ctemp;
00146   sprintf(Ctemp,"DecayLength: %f  -- %f\n",mDecayLength[0],mDecayLength[1]);
00147   Stemp += Ctemp;
00148   sprintf(Ctemp,"AngleToPrimaryVertex: %f  -- %f\n",mAngleToPrimaryVertex[0],mAngleToPrimaryVertex[1]);
00149   Stemp += Ctemp;
00150   sprintf(Ctemp,"DcaToPrimaryVertex: %f  -- %f\n",mDcaToPrimaryVertex[0],mDcaToPrimaryVertex[1]);
00151   Stemp += Ctemp;
00152   StHbtString returnThis = Stemp;
00153   return returnThis;
00154 }
00155