StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
OpeningAngleCorrFctn.cxx
1 /***************************************************************************
2  *
3  * Author: Mike Lisa, Ohio State, lisa@mps.ohio-state.edu
4  ***************************************************************************
5  *
6  * Description: part of STAR HBT Framework: StHbtMaker package
7  * a simple opening-angle correlation function used for studying 2-track cuts
8  *
9  ***************************************************************************
10  *
11  **************************************************************************/
12 
13 #include "StHbtMaker/CorrFctn/OpeningAngleCorrFctn.h"
14 #include <cstdio>
15 
16 #ifdef __ROOT__
17 ClassImp(OpeningAngleCorrFctn)
18 #endif
19 //____________________________
20 OpeningAngleCorrFctn:: OpeningAngleCorrFctn(char* title, const int& nbinsQ, const float& QLo, const float& QHi,
21  const int& nbinsAng, const float& AngLo, const float& AngHi){
22  // set up numeratorS
23  char Tit[100];
24  sprintf(Tit,"2D Num");
25  strcat(Tit,title);
26  mNumerator2D = new StHbt2DHisto(Tit,title,nbinsQ,QLo,QHi,nbinsAng,AngLo,AngHi);
27 
28  // set up denominatorS
29  sprintf(Tit,"2D Den");
30  strcat(Tit,title);
31  mDenominator2D = new StHbt2DHisto(Tit,title,nbinsQ,QLo,QHi,nbinsAng,AngLo,AngHi);
32 
33  // set up ratioS
34  sprintf(Tit,"2D Rat");
35  strcat(Tit,title);
36  mRatio2D = new StHbt2DHisto(Tit,title,nbinsQ,QLo,QHi,nbinsAng,AngLo,AngHi);
37 
38  // these histograms should have errors associated with them...
39  mNumerator2D->Sumw2();
40  mDenominator2D->Sumw2();
41  mRatio2D->Sumw2();
42 
43 }
44 
45 //____________________________
46 OpeningAngleCorrFctn::~OpeningAngleCorrFctn(){
47  delete mNumerator2D;
48  delete mDenominator2D;
49  delete mRatio2D;
50 }
51 //_________________________
52 void OpeningAngleCorrFctn::Finish(){
53  // here is where we should normalize, fit, etc...
54  // we should NOT Draw() the histos (as I had done it below),
55  // since we want to insulate ourselves from root at this level
56  // of the code. Do it instead at root command line with browser.
57  mRatio2D->Divide(mNumerator2D,mDenominator2D,1.0,1.0);
58 }
59 
60 //____________________________
61 StHbtString OpeningAngleCorrFctn::Report(){
62  string stemp = "Opening Angle Correlation Function Report:\n";
63  char ctemp[100];
64  sprintf(ctemp,"Number of entries in numerator:\t%E\n",
65  mNumerator2D->GetEntries());
66  stemp += ctemp;
67  sprintf(ctemp,"Number of entries in denominator:\t%E\n",
68  mDenominator2D->GetEntries());
69  stemp += ctemp;
70  StHbtString returnThis = stemp;
71  return returnThis;
72 }
73 //____________________________
74 void OpeningAngleCorrFctn::AddRealPair(const StHbtPair* pair){
75  StHbtThreeVector p1 = pair->track1()->FourMomentum().vect();
76  StHbtThreeVector p2 = pair->track2()->FourMomentum().vect();
77 
78  // get relative angle in plane transverse to z:
79  // StHbtThreeVector p1Tran = p1; p1Tran.setZ(0.0);
80  //StHbtThreeVector p2Tran = p2; p2Tran.setZ(0.0);
81  //double dAngTran = 57.296*acos((p1Tran.dot(p2Tran))/(p1Tran.mag()*p2Tran.mag()));
82  //mNumeratorTran->Fill(dAngTran);
83 
84  // get "absolute" relative angle
85  double dAngInv = 57.296*acos((p1.dot(p2))/(p1.mag()*p2.mag()));
86 
87  double Qinv = fabs(pair->qInv()); // note - qInv() will be negative for identical pairs...
88 
89  mNumerator2D->Fill(Qinv,dAngInv,1.0);
90 }
91 //____________________________
92 void OpeningAngleCorrFctn::AddMixedPair(const StHbtPair* pair){
93  StHbtThreeVector p1 = pair->track1()->FourMomentum().vect();
94  StHbtThreeVector p2 = pair->track2()->FourMomentum().vect();
95 
96  // get relative angle in plane transverse to z:
97  // StHbtThreeVector p1Tran = p1; p1Tran.setZ(0.0);
98 // StHbtThreeVector p2Tran = p2; p2Tran.setZ(0.0);
99 // double dAngTran = 57.296*acos((p1Tran.dot(p2Tran))/(p1Tran.mag()*p2Tran.mag()));
100 // mDenominatorTran->Fill(dAngTran);
101 
102  // get "absolute" relative angle
103  double dAngInv = 57.296*acos((p1.dot(p2))/(p1.mag()*p2.mag()));
104  double Qinv = fabs(pair->qInv()); // note - qInv() will be negative for identical pairs...
105 
106  mDenominator2D->Fill(Qinv,dAngInv,1.0);
107 }
108 
109 
110