StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
BPLCMSFrame3DCorrFctn_SIM.cxx
1 /***************************************************************************
2  *
3  * $Id: BPLCMSFrame3DCorrFctn_SIM.cxx,v 1.6 2003/09/02 17:58:20 perev Exp $
4  *
5  * Author: Mike Lisa, Ohio State, lisa@mps.ohio-state.edu
6  ***************************************************************************
7  *
8  * Description: part of STAR HBT Framework: StHbtMaker package
9  * 3D Bertsch-Pratt decomposition in the LCMS frame
10  * THIS IS A SIMULATION CORRELATION FUNCTION CLASS !!!
11  *
12  *
13  ***************************************************************************
14  *
15  * $Log: BPLCMSFrame3DCorrFctn_SIM.cxx,v $
16  * Revision 1.6 2003/09/02 17:58:20 perev
17  * gcc 3.2 updates + WarnOff
18  *
19  * Revision 1.5 2003/01/31 19:21:09 magestro
20  * Cleared up simple compiler warnings on i386_linux24
21  *
22  * Revision 1.4 2001/05/23 00:19:05 lisa
23  * Add in Smearing classes and methods needed for momentum resolution studies and correction
24  *
25  * Revision 1.3 2000/10/08 17:11:07 lisa
26  * now toggle between Num and Den at BEGINNING of BPLCMSFrame3DCorrFctn_SIM::AddMixedPair()
27  *
28  * Revision 1.2 2000/10/05 23:08:59 lisa
29  * Added kT-dependent radii to mixed-event simulator AND implemented AverageSeparation Cut and CorrFctn
30  *
31  * Revision 1.1 2000/09/14 18:36:53 lisa
32  * Added Qinv and ExitSep pair cuts and BPLCMSFrame3DCorrFctn_SIM CorrFctn
33  *
34  *
35  *
36  **************************************************************************/
37 
38 
39 #include "StHbtMaker/CorrFctn/BPLCMSFrame3DCorrFctn_SIM.h"
40 #include "StHbtMaker/Infrastructure/StHbtSmearPair.h"
41 
42 #include <cstdio>
43 
44 #ifdef __ROOT__
45 ClassImp(BPLCMSFrame3DCorrFctn_SIM)
46 #endif
47 
48 //____________________________
49 BPLCMSFrame3DCorrFctn_SIM::BPLCMSFrame3DCorrFctn_SIM(char* title, const int& nbins, const float& QLo, const float& QHi){
50 
51  // set some stuff...
52  mQinvNormLo = 0.15;
53  mQinvNormHi = 0.18;
54  mNumRealsNorm = 0;
55  mNumMixedNorm = 0;
56  mCorrection = 0; // pointer to Coulomb Correction object
57 
58  mPairCut = 0; // added Sept2000 - CorrFctn-specific PairCut
59 
60  mToggleNumDen = 0; // this toggles between a pair being given to numerator or denominator.
61 
62  mLambda=0.;
63  mRside2=25.0;
64  mRout2=25.0;
65  mRlong2=25.0;
66 
67  mRout_alpha = mRside_alpha = mRlong_alpha = 0; // set these zero by default
68 
69  // set up numerator
70  char TitNum[100] = "Num";
71  strcat(TitNum,title);
72  mNumerator = new StHbt3DHisto(TitNum,title,nbins,QLo,QHi,nbins,QLo,QHi,nbins,QLo,QHi);
73  // set up denominator
74  char TitDen[100] = "Den";
75  strcat(TitDen,title);
76  mDenominator = new StHbt3DHisto(TitDen,title,nbins,QLo,QHi,nbins,QLo,QHi,nbins,QLo,QHi);
77  // set up ratio
78  char TitRat[100] = "Rat";
79  strcat(TitRat,title);
80  mRatio = new StHbt3DHisto(TitRat,title,nbins,QLo,QHi,nbins,QLo,QHi,nbins,QLo,QHi);
81 
82  // to enable error bar calculation...
83  mNumerator->Sumw2();
84  mDenominator->Sumw2();
85  mRatio->Sumw2();
86 
87  // 20feb2001 - add some more histograms that are filled if there is any
88  // momentum smearing
89 
90  char titInv[100] = "QinvRes";
91  strcat(titInv,title);
92  mResolutionHistos[0] = new StHbt2DHisto(titInv,title,10,0.0,0.1,50,-0.025,0.025);
93  //
94  char titOut[100] = "QoutRes";
95  strcat(titOut,title);
96  mResolutionHistos[1] = new StHbt2DHisto(titOut,title,10,0.0,0.1,50,-0.025,0.025);
97  //
98  char titSid[100] = "QsidRes";
99  strcat(titSid,title);
100  mResolutionHistos[2] = new StHbt2DHisto(titSid,title,10,0.0,0.1,50,-0.025,0.025);
101  //
102  char titLon[100] = "QlonRes";
103  strcat(titLon,title);
104  mResolutionHistos[3] = new StHbt2DHisto(titLon,title,10,0.0,0.1,50,-0.025,0.025);
105 
106 }
107 
108 //____________________________
109 BPLCMSFrame3DCorrFctn_SIM::~BPLCMSFrame3DCorrFctn_SIM(){
110  delete mNumerator;
111  delete mDenominator;
112  delete mRatio;
113 }
114 //_________________________
115 void BPLCMSFrame3DCorrFctn_SIM::Finish(){
116 
117  mRatio->Divide(mNumerator,mDenominator);
118 }
119 
120 //____________________________
121 StHbtString BPLCMSFrame3DCorrFctn_SIM::Report(){
122  string stemp = "LCMS Frame Bertsch-Pratt 3D Correlation Function Report:\n";
123  char ctemp[100];
124  sprintf(ctemp," THIS IS A SIMULATION CORRELATION FUNCTION CLASS!!!\n");
125  stemp += ctemp;
126  sprintf(ctemp,"Number of entries in numerator:\t%E\n",mNumerator->GetEntries());
127  stemp += ctemp;
128  sprintf(ctemp,"Number of entries in denominator:\t%E\n",mDenominator->GetEntries());
129  stemp += ctemp;
130  sprintf(ctemp,"Number of entries in ratio:\t%E\n",mRatio->GetEntries());
131  stemp += ctemp;
132  sprintf(ctemp,"Normalization region in Qinv was:\t%E\t%E\n",mQinvNormLo,mQinvNormHi);
133  stemp += ctemp;
134  sprintf(ctemp,"Number of pairs in Normalization region was:\n");
135  stemp += ctemp;
136  sprintf(ctemp,"In numerator:\t%lu\t In denominator:\t%lu\n",mNumRealsNorm,mNumMixedNorm);
137  stemp += ctemp;
138  if (mCorrection)
139  {
140  float radius = mCorrection->GetRadius();
141  sprintf(ctemp,"Coulomb correction used radius of\t%E\n",radius);
142  }
143  else
144  {
145  sprintf(ctemp,"No Coulomb Correction applied to this CorrFctn\n");
146  }
147  stemp += ctemp;
148 
149  if (mPairCut){
150  sprintf(ctemp,"Here is the PairCut specific to this CorrFctn\n");
151  stemp += ctemp;
152  stemp += mPairCut->Report();
153  }
154  else{
155  sprintf(ctemp,"No PairCut specific to this CorrFctn\n");
156  stemp += ctemp;
157  }
158 
159  //
160  StHbtString returnThis = stemp;
161  return returnThis;
162 }
163 //____________________________
164 void BPLCMSFrame3DCorrFctn_SIM::AddRealPair(const StHbtPair* pair){
165 
166  // note that in this SIMULATION correlation function class, we do NOTHING with real
167  // pairs. The numerator is just a weighted denominator, filled in AddMixedPair (below)
168 
169  /* no-op */
170  return;
171 
172 }
173 
174 
175 
176 
177 //____________________________
178 void BPLCMSFrame3DCorrFctn_SIM::AddMixedPair(const StHbtPair* pair){
179 
180  mToggleNumDen = !mToggleNumDen; // toggle at the BEGINNING of call instead of end
181  // so that SAME pairs go into Num and Den as long as they pass cut
182 
183 
184  if (mPairCut){
185  if (!(mPairCut->Pass(pair))){
186 // // mike is checking stuff here
187 // double qOut = fabs(pair->qOutCMS());
188 // double qSide = fabs(pair->qSideCMS());
189 // double qLong = fabs(pair->qLongCMS());
190 // // if ((qSide > 0.035)||(qLong > 0.035)){
191 // if ((qSide > 0.0)||(qLong > 0.0)){
192 // cout << qOut << " " << qSide << " " << qLong << " " <<
193 // pair->NominalTpcEntranceSeparation() << " " <<
194 // pair->track1()->NominalTpcEntrancePoint() << " " <<
195 // pair->track2()->NominalTpcEntrancePoint() << endl;
196 // int ijunk;
197 // cout << "Enter a junk integer : ";
198 // cin >> ijunk;
199 // }
200 // // end of checking stuff
201  return;
202  }
203  }
204 
205  double weight=1.0;
206 
207  // these are the values we use to assign correlation WEIGHT
208  double qOut = fabs(pair->qOutCMS());
209  double qSide = fabs(pair->qSideCMS());
210  double qLong = fabs(pair->qLongCMS());
211 
212  // these are the values we BIN in - same as above if no smearing turned on...
213  double qOut_bin,qSide_bin,qLong_bin; // scope 'em!
214  if (mSmearPair){
215  mSmearPair->SetUnsmearedPair(pair);
216  qOut_bin = fabs(mSmearPair->SmearedPair().qOutCMS());
217  qSide_bin = fabs(mSmearPair->SmearedPair().qSideCMS());
218  qLong_bin = fabs(mSmearPair->SmearedPair().qLongCMS());
219  // fill resolution histos...
220  double Qinv = fabs(pair->qInv());
221  if (Qinv < 0.1){
222  double Qinv_smear = fabs(mSmearPair->SmearedPair().qInv());
223  mResolutionHistos[0]->Fill(Qinv,Qinv_smear-Qinv);
224  mResolutionHistos[1]->Fill(qOut,qOut_bin-qOut);
225  mResolutionHistos[2]->Fill(qSide,qSide_bin-qSide);
226  mResolutionHistos[3]->Fill(qLong,qLong_bin-qLong);
227  }
228  }
229  else{
230  qOut_bin = qOut;
231  qSide_bin = qSide;
232  qLong_bin = qLong;
233  }
234 
235 
236  // note the Coulomb bit below...
237  // for the *suppression* ("numerator" pairs), we use the TRUE momentum, as Nature would
238  // for the *correction* ("denominator" pairs), we use the SMEARED momentum, as an experimenter would
239 
240  if (mToggleNumDen){
241  if (mCorrection){
242  if (mSmearPair){
243  weight = mCorrection->CoulombCorrect(&(mSmearPair->SmearedPair()));
244  }
245  else{
246  weight = mCorrection->CoulombCorrect(pair);
247  }
248  }
249  mDenominator->Fill(qOut_bin,qSide_bin,qLong_bin,weight);
250  }
251  else{
252 
253  if (mCorrection) weight = mCorrection->CoulombCorrect(pair);
254 
255  // update 1oct2000 - now allow mT dependent evolution of Radii
256 
257  float mT2 = (pair->fourMomentumSum().mt2())/4.0; // divide by 4.0 because want mT of (p1+p2)/2
258  // and factor of two gets squared in mt2
259 
260  float Rout2,Rside2,Rlong2;
261 
262  if (mRout_alpha!=0){Rout2 = mRout2*::pow(mT2,mRout_alpha);}
263  else{Rout2 = mRout2;}
264 
265  if (mRside_alpha!=0){Rside2 = mRside2*::pow(mT2,mRside_alpha);}
266  else{Rside2 = mRside2;}
267 
268  if (mRlong_alpha!=0){Rlong2 = mRlong2*::pow(mT2,mRlong_alpha);}
269  else{Rlong2 = mRlong2;}
270 
271  double CorrWeight = 1.0 +
272  mLambda*exp((-qOut*qOut*Rout2 -qSide*qSide*Rside2 -qLong*qLong*Rlong2)/0.038936366329);
273  CorrWeight *= weight;
274  mNumerator->Fill(qOut_bin,qSide_bin,qLong_bin,CorrWeight);
275  }
276 
277 }
278 
279 
280 
BPLCMSFrame3DCorrFctn_SIM
Definition: BPLCMSFrame3DCorrFctn_SIM.h:37
Generated by   doxygen 1.8.5