Fixing Effective C++ warnings
[u/mrichter/AliRoot.git] / PWG2 / FEMTOSCOPY / AliFemto / CorrFctn / AliFemtoQinvCorrFctn.cxx
1 /***************************************************************************
2  *
3  * $Id$
4  *
5  * Author: Mike Lisa, Ohio State, lisa@mps.ohio-state.edu
6  ***************************************************************************
7  *
8  * Description: part of STAR HBT Framework: AliFemtoMaker package
9  *   a simple Q-invariant correlation function           
10  *
11  ***************************************************************************
12  *
13  * $Log$
14  * Revision 1.1.1.1  2007/04/25 15:38:41  panos
15  * Importing the HBT code dir
16  *
17  * Revision 1.1.1.1  2007/03/07 10:14:49  mchojnacki
18  * First version on CVS
19  *
20  * Revision 1.4  2000/01/25 17:34:45  laue
21  * I. In order to run the stand alone version of the AliFemtoMaker the following
22  * changes have been done:
23  * a) all ClassDefs and ClassImps have been put into #ifdef __ROOT__ statements
24  * b) unnecessary includes of StMaker.h have been removed
25  * c) the subdirectory AliFemtoMaker/doc/Make has been created including everything
26  * needed for the stand alone version
27  *
28  * II. To reduce the amount of compiler warning
29  * a) some variables have been type casted
30  * b) some destructors have been declared as virtual
31  *
32  * Revision 1.3  1999/07/29 02:47:09  lisa
33  * 1) add OpeningAngle correlation function 2) add AliFemtoMcEventReader 3) make histos in CorrFctns do errors correctly
34  *
35  * Revision 1.2  1999/07/06 22:33:20  lisa
36  * Adjusted all to work in pro and new - dev itself is broken
37  *
38  * Revision 1.1.1.1  1999/06/29 16:02:57  lisa
39  * Installation of AliFemtoMaker
40  *
41  **************************************************************************/
42
43 #include "CorrFctn/AliFemtoQinvCorrFctn.h"
44 //#include "Infrastructure/AliFemtoHisto.h"
45 #include <cstdio>
46
47 #ifdef __ROOT__ 
48 ClassImp(AliFemtoQinvCorrFctn)
49 #endif
50
51 //____________________________
52 AliFemtoQinvCorrFctn::AliFemtoQinvCorrFctn(char* title, const int& nbins, const float& QinvLo, const float& QinvHi):
53   fNumerator(0),
54   fDenominator(0),
55   fRatio(0)
56 {
57   // set up numerator
58   //  title = "Num Qinv (MeV/c)";
59   char TitNum[100] = "Num";
60   strcat(TitNum,title);
61   fNumerator = new TH1D(TitNum,title,nbins,QinvLo,QinvHi);
62   // set up denominator
63   //title = "Den Qinv (MeV/c)";
64   char TitDen[100] = "Den";
65   strcat(TitDen,title);
66   fDenominator = new TH1D(TitDen,title,nbins,QinvLo,QinvHi);
67   // set up ratio
68   //title = "Ratio Qinv (MeV/c)";
69   char TitRat[100] = "Rat";
70   strcat(TitRat,title);
71   fRatio = new TH1D(TitRat,title,nbins,QinvLo,QinvHi);
72   // this next bit is unfortunately needed so that we can have many histos of same "title"
73   // it is neccessary if we typedef TH1D to TH1d (which we do)
74   //fNumerator->SetDirectory(0);
75   //fDenominator->SetDirectory(0);
76   //fRatio->SetDirectory(0);
77
78   // to enable error bar calculation...
79   fNumerator->Sumw2();
80   fDenominator->Sumw2();
81   fRatio->Sumw2();
82
83 }
84
85 //____________________________
86 AliFemtoQinvCorrFctn::AliFemtoQinvCorrFctn(const AliFemtoQinvCorrFctn& aCorrFctn) :
87   fNumerator(0),
88   fDenominator(0),
89   fRatio(0)
90 {
91   fNumerator = new TH1D(*aCorrFctn.fNumerator);
92   fDenominator = new TH1D(*aCorrFctn.fDenominator);
93   fRatio = new TH1D(*aCorrFctn.fRatio);
94 }
95 //____________________________
96 AliFemtoQinvCorrFctn::~AliFemtoQinvCorrFctn(){
97   delete fNumerator;
98   delete fDenominator;
99   delete fRatio;
100 }
101 //_________________________
102 AliFemtoQinvCorrFctn& AliFemtoQinvCorrFctn::operator=(const AliFemtoQinvCorrFctn& aCorrFctn)
103 {
104   if (this == &aCorrFctn)
105     return *this;
106
107   if (fNumerator) delete fNumerator;
108   fNumerator = new TH1D(*aCorrFctn.fNumerator);
109   if (fDenominator) delete fDenominator;
110   fDenominator = new TH1D(*aCorrFctn.fDenominator);
111   if (fRatio) delete fRatio;
112   fRatio = new TH1D(*aCorrFctn.fRatio);
113
114   return *this;
115 }
116
117 //_________________________
118 void AliFemtoQinvCorrFctn::Finish(){
119   // here is where we should normalize, fit, etc...
120   // we should NOT Draw() the histos (as I had done it below),
121   // since we want to insulate ourselves from root at this level
122   // of the code.  Do it instead at root command line with browser.
123   //  fNumerator->Draw();
124   //fDenominator->Draw();
125   //fRatio->Draw();
126   fRatio->Divide(fNumerator,fDenominator,1.0,1.0);
127
128 }
129
130 //____________________________
131 AliFemtoString AliFemtoQinvCorrFctn::Report(){
132   string stemp = "Qinv Correlation Function Report:\n";
133   char ctemp[100];
134   sprintf(ctemp,"Number of entries in numerator:\t%E\n",fNumerator->GetEntries());
135   stemp += ctemp;
136   sprintf(ctemp,"Number of entries in denominator:\t%E\n",fDenominator->GetEntries());
137   stemp += ctemp;
138   sprintf(ctemp,"Number of entries in ratio:\t%E\n",fRatio->GetEntries());
139   stemp += ctemp;
140   //  stemp += mCoulombWeight->Report();
141   AliFemtoString returnThis = stemp;
142   return returnThis;
143 }
144 //____________________________
145 void AliFemtoQinvCorrFctn::AddRealPair(const AliFemtoPair* pair){
146   double Qinv = fabs(pair->qInv());   // note - qInv() will be negative for identical pairs...
147   fNumerator->Fill(Qinv);
148   //  cout << "AliFemtoQinvCorrFctn::AddRealPair : " << pair->qInv() << " " << Qinv <<
149   //" " << pair->track1().FourMomentum() << " " << pair->track2().FourMomentum() << endl;
150 }
151 //____________________________
152 void AliFemtoQinvCorrFctn::AddMixedPair(const AliFemtoPair* pair){
153   double weight = 1.0;
154   double Qinv = fabs(pair->qInv());   // note - qInv() will be negative for identical pairs...
155   fDenominator->Fill(Qinv,weight);
156 }
157
158