1 ////////////////////////////////////////////////////////////////////////////////
3 /// AliFemtoShareQualityCorrFctn - A correlation function that saves the ///
4 /// amount of sharing and splitting hits per pair as a function of qinv ///
5 /// Authors: Adam Kisiel kisiel@mps.ohio-state.edu ///
7 ////////////////////////////////////////////////////////////////////////////////
9 #include "AliFemtoShareQualityCorrFctn.h"
10 //#include "Infrastructure/AliFemtoHisto.hh"
14 ClassImp(AliFemtoShareQualityCorrFctn)
17 //____________________________
18 AliFemtoShareQualityCorrFctn::AliFemtoShareQualityCorrFctn(char* title, const int& nbins, const float& QinvLo, const float& QinvHi):
22 fQualityDenominator(0)
25 // title = "Num Qinv (MeV/c)";
26 char TitNum[100] = "NumShare";
28 fShareNumerator = new TH2D(TitNum,title,nbins,QinvLo,QinvHi,50,0.0,1.00001);
30 //title = "Den Qinv (MeV/c)";
31 char TitDen[100] = "DenShare";
33 fShareDenominator = new TH2D(TitDen,title,nbins,QinvLo,QinvHi,50,0.0,1.00001);
35 char Tit2Num[100] = "NumQuality";
36 strcat(Tit2Num,title);
37 fQualityNumerator = new TH2D(Tit2Num,title,nbins,QinvLo,QinvHi,75,-0.500001,1.000001);
39 //title = "Den Qinv (MeV/c)";
40 char Tit2Den[100] = "DenQuality";
41 strcat(Tit2Den,title);
42 fQualityDenominator = new TH2D(Tit2Den,title,nbins,QinvLo,QinvHi,75,-0.500001,1.000001);
44 //title = "Ratio Qinv (MeV/c)";
45 // this next bit is unfortunately needed so that we can have many histos of same "title"
46 // it is neccessary if we typedef TH2D to TH1d (which we do)
47 //mShareNumerator->SetDirectory(0);
48 //mShareDenominator->SetDirectory(0);
49 //mRatio->SetDirectory(0);
51 // to enable error bar calculation...
52 fShareNumerator->Sumw2();
53 fShareDenominator->Sumw2();
55 fQualityNumerator->Sumw2();
56 fQualityDenominator->Sumw2();
59 //____________________________
60 AliFemtoShareQualityCorrFctn::AliFemtoShareQualityCorrFctn(const AliFemtoShareQualityCorrFctn& aCorrFctn) :
64 fQualityDenominator(0)
66 if (aCorrFctn.fShareNumerator)
67 fShareNumerator = new TH2D(*aCorrFctn.fShareNumerator);
68 if (aCorrFctn.fShareDenominator)
69 fShareDenominator = new TH2D(*aCorrFctn.fShareDenominator);
70 if (aCorrFctn.fQualityNumerator)
71 fQualityNumerator = new TH2D(*aCorrFctn.fQualityNumerator);
72 if (aCorrFctn.fQualityDenominator)
73 fQualityDenominator = new TH2D(*aCorrFctn.fQualityDenominator);
75 //____________________________
76 AliFemtoShareQualityCorrFctn::~AliFemtoShareQualityCorrFctn(){
77 delete fShareNumerator;
78 delete fShareDenominator;
79 delete fQualityNumerator;
80 delete fQualityDenominator;
82 //_________________________
83 AliFemtoShareQualityCorrFctn& AliFemtoShareQualityCorrFctn::operator=(const AliFemtoShareQualityCorrFctn& aCorrFctn)
85 if (this == &aCorrFctn)
88 if (aCorrFctn.fShareNumerator)
89 fShareNumerator = new TH2D(*aCorrFctn.fShareNumerator);
92 if (aCorrFctn.fShareDenominator)
93 fShareDenominator = new TH2D(*aCorrFctn.fShareDenominator);
95 fShareDenominator = 0;
96 if (aCorrFctn.fQualityNumerator)
97 fQualityNumerator = new TH2D(*aCorrFctn.fQualityNumerator);
99 fQualityNumerator = 0;
100 if (aCorrFctn.fQualityDenominator)
101 fQualityDenominator = new TH2D(*aCorrFctn.fQualityDenominator);
103 fQualityDenominator = 0;
107 //_________________________
108 void AliFemtoShareQualityCorrFctn::Finish(){
109 // here is where we should normalize, fit, etc...
110 // we should NOT Draw() the histos (as I had done it below),
111 // since we want to insulate ourselves from root at this level
112 // of the code. Do it instead at root command line with browser.
113 // mShareNumerator->Draw();
114 //mShareDenominator->Draw();
119 //____________________________
120 AliFemtoString AliFemtoShareQualityCorrFctn::Report(){
121 string stemp = "Qinv Correlation Function Report:\n";
123 sprintf(ctemp,"Number of entries in numerator:\t%E\n",fShareNumerator->GetEntries());
125 sprintf(ctemp,"Number of entries in denominator:\t%E\n",fShareDenominator->GetEntries());
127 // stemp += mCoulombWeight->Report();
128 AliFemtoString returnThis = stemp;
131 //____________________________
132 void AliFemtoShareQualityCorrFctn::AddRealPair(const AliFemtoPair* pair){
133 double Qinv = fabs(pair->qInv()); // note - qInv() will be negative for identical pairs...
138 for (unsigned int imap=0; imap<pair->track1()->Track()->TPCclusters().GetNbits(); imap++) {
139 // If both have clusters in the same row
140 if (pair->track1()->Track()->TPCclusters().TestBitNumber(imap) &&
141 pair->track2()->Track()->TPCclusters().TestBitNumber(imap)) {
142 // Do they share it ?
143 if (pair->track1()->Track()->TPCsharing().TestBitNumber(imap) ||
144 pair->track2()->Track()->TPCsharing().TestBitNumber(imap))
147 cout << "Shared cluster in row " << imap << endl;
154 // Different hits on the same padrow
160 else if (pair->track1()->Track()->TPCclusters().TestBitNumber(imap) ||
161 pair->track2()->Track()->TPCclusters().TestBitNumber(imap)) {
162 // One track has a hit, the other does not
168 cout << "Qinv of the pair is " << Qinv << endl;
169 cout << "Clusters: " << endl;
170 for (unsigned int imap=0; imap<pair->track1()->Track()->TPCclusters().GetNbits(); imap++) {
172 if (pair->track1()->Track()->TPCclusters().TestBitNumber(imap)) cout << " 1 ";
174 if (pair->track2()->Track()->TPCclusters().TestBitNumber(imap)) cout << " 1 ";
177 if (pair->track1()->Track()->TPCsharing().TestBitNumber(imap)) cout << " S ";
179 if (pair->track2()->Track()->TPCsharing().TestBitNumber(imap)) cout << " S ";
185 Float_t hsmval = 0.0;
186 Float_t hsfval = 0.0;
194 cout << "Quality Sharity " << hsmval << " " << hsfval << " " << pair->track1()->Track() << " " << pair->track2()->Track() << endl;
197 fShareNumerator->Fill(Qinv, hsfval);
198 fQualityNumerator->Fill(Qinv, hsmval);
199 // cout << "AliFemtoShareQualityCorrFctn::AddRealPair : " << pair->qInv() << " " << Qinv <<
200 //" " << pair->track1().FourMomentum() << " " << pair->track2().FourMomentum() << endl;
202 //____________________________
203 void AliFemtoShareQualityCorrFctn::AddMixedPair(const AliFemtoPair* pair){
205 double Qinv = fabs(pair->qInv()); // note - qInv() will be negative for identical pairs...
210 for (unsigned int imap=0; imap<pair->track1()->Track()->TPCclusters().GetNbits(); imap++) {
211 // If both have clusters in the same row
212 if (pair->track1()->Track()->TPCclusters().TestBitNumber(imap) &&
213 pair->track2()->Track()->TPCclusters().TestBitNumber(imap)) {
214 // Do they share it ?
215 if (pair->track1()->Track()->TPCsharing().TestBitNumber(imap) ||
216 pair->track2()->Track()->TPCsharing().TestBitNumber(imap))
218 // cout << "A shared cluster !!!" << endl;
219 // cout << "imap idx1 idx2 "
221 // << tP1idx[imap] << " " << tP2idx[imap] << endl;
227 // Different hits on the same padrow
233 else if (pair->track1()->Track()->TPCclusters().TestBitNumber(imap) ||
234 pair->track2()->Track()->TPCclusters().TestBitNumber(imap)) {
235 // One track has a hit, the other does not
241 Float_t hsmval = 0.0;
242 Float_t hsfval = 0.0;
249 fShareDenominator->Fill(Qinv,hsfval,weight);
250 fQualityDenominator->Fill(Qinv,hsmval,weight);
254 void AliFemtoShareQualityCorrFctn::WriteHistos()
256 fShareNumerator->Write();
257 fShareDenominator->Write();
258 fQualityNumerator->Write();
259 fQualityDenominator->Write();