StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
StHbtExampleII.C
1 // example how to access memberfunctions from a correlation function
2 // ((MinvCorrFctn*) ((StHbtAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(0))->CorrFctn(0))->Difference()->Draw()
3 
4 // ((MinvCorrFctnY_vs_Pt*) ((StHbtAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(0))->CorrFctn(2))->Difference()->Draw("colz")
5 
6 // ((MinvCorrFctnArmenteros*) ((StHbtAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(0))->CorrFctn(1))->Difference()->Draw("colz")
7 
8 // examples to access member function from cutMonitors (here the member functions return a pointer to a histogram)
9 
10 // for Like Sign Analysis
11 // ((MinvLikeSignCorrFctn*) ((StHbtLikeSignAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(1))->CorrFctn(0))->Numerator()->Draw()
12 // ((MinvLikeSignCorrFctn*) ((StHbtLikeSignAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(1))->CorrFctn(0))->MixedEventDenominator()->Draw()
13 // ((MinvLikeSignCorrFctn*) ((StHbtLikeSignAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(1))->CorrFctn(0))->PositiveDenominator()->Draw()
14 // ((MinvLikeSignCorrFctn*) ((StHbtLikeSignAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(1))->CorrFctn(0))->NegativeDenominator()->Draw()
15 // ((MinvLikeSignCorrFctn*) ((StHbtLikeSignAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(1))->CorrFctn(0))->MixedEventDifference()->Draw()
16 // ((MinvLikeSignCorrFctn*) ((StHbtLikeSignAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(1))->CorrFctn(0))->LikeSignDifference()->Draw()
17 // ((trackCutMonitor_P_vs_Dedx*) ((franksTrackCut*) ((StHbtLikeSignAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(1))->FirstParticleCut())->PassMonitor(0))->Histo()->Draw()
18 // ((trackCutMonitor_P_vs_Dedx*) ((franksTrackCut*) ((StHbtLikeSignAnalysis*)((StHbtMaker*) chain->GetMaker("HBT"))->HbtManager()->Analysis(1))->FirstParticleCut())->FailMonitor(0))->Histo()->Draw()
19 
20 #define STEVENT 0
21 #define MCEVENT 0
22 #define ASSOCIATION 0
23 #define STEVENTREADER 0
24 #define MCEVENTREADER 0
25 #define ASSOCIATIONREADER 0
26 #define XDFREADER 0
27 #define GSTARTXTREADER 0
28 
29 #define ASCIIWRITER 0
30 #define ASCIIREADER 0
31 
32 #define BINARYWRITER 0
33 #define BINARYREADER 1
34 
35 #if ASCIIREADER
36 #define STEVENT 0
37 #define MCEVENT 0
38 #define ASSOCIATION 0
39 #define STEVENTREADER 0
40 #define MCEVENTREADER 0
41 #define ASSOCIATIONREADER 0
42 #endif
43 
44 #if BINARYREADER
45 #define STEVENT 0
46 #define MCEVENT 0
47 #define ASSOCIATION 0
48 #define STEVENTREADER 0
49 #define MCEVENTREADER 0
50 #define ASSOCIATIONREADER 0
51 #endif
52 
53 #define ANALYSIS
54 
55 // NOTE - chain needs to be declared global so for StHbtEventReader
56 class StChain;
57 StChain *chain=0;
58 
59 // keep pointers to Correlation Functions global, so you can have access to them...
60 class QinvCorrFctn;
61 QinvCorrFctn* QinvCF;
62 
63 class QvecCorrFctn;
64 QvecCorrFctn* QvecCF;
65 
66 class MinvCorrFctn;
67 MinvCorrFctn* MinvCF;
68 MinvCorrFctn* MinvCF2;
69 MinvCorrFctn* MinvCFrho;
70 
72 MinvCorrFctnM_vs_Pt* MinvCFM_vs_Pt;
73 
75 MinvCorrFctnY_vs_Pt* MinvCFY_vs_Pt;
76 
77 class MinvCorrFctnM_vs_P;
78 MinvCorrFctnM_vs_P* MinvCFM_vs_P;
79 
80 class p1_vs_p2CorrFctn;
81 p1_vs_p2CorrFctn* p1_vs_p2_CF;
82 
83 class YCorrFctn;
84 YCorrFctn* YdiffCF;
85 
86 class pDiffCorrFctn;
87 pDiffCorrFctn* pDiffCF;
88 
89 class p_vs_angleCorrFctn;
90 p_vs_angleCorrFctn* p_vs_angleCF;
91 
92 class phiMeas_vs_phiCalcCorrFctn;
93 phiMeas_vs_phiCalcCorrFctn* phiMeas_vs_phiCalcCF;
94 
95 // keep pointers to Analysis global, so you can have access to themm ...
96 class StHbtAnalysis;
97 StHbtAnalysis* phiAnal;
98 StHbtAnalysis* phiAnal2;
99 StHbtAnalysis* phiAnal3;
100 StHbtAnalysis* phiAnalCopy;
101 StHbtAnalysis* deltaAnal;
102 StHbtAnalysis* rhoAnal;
103 StHbtAnalysis* lambdaAnal;
104 
105 // File-scope stuff needed by setFiles, nextFile. Someone ambitious
106 // can clean this up by putting it all into a nice clean class.
107 Int_t usePath = 0;
108 Int_t nFile = 0;
109 TString thePath;
110 TString theFileName;
111 TString originalPath;
112 //class StChain;
113 //StChain *chain=0;
114 
115 TBrowser *b=0;
116 const char *venusFile ="set*geant.root";
117 const char *venusPath ="/disk00001/star/auau200/venus412/default/b0_3/year_1b/hadronic_on/tfs_4/";
118 const char *dstFile ="/disk00001/star/auau200/two_photon/starlight/twogam/year_1b/hadronic_on/tfs/ric0022_01_14552evts.dst.root";
119 const char *xdfFile ="/afs/rhic.bnl.gov/star/data/samples/psc0054_07_40evts_dst.xdf";
120 const char *mdcFile ="/disk00001/star/auau200/venus412/default/b0_3/year_1b/hadronic_on/tss/psc0081_07_40evts.root";
121 const char *geantFile ="/disk00000/star/auau200/hijing135/jetq_off/b0_3/year_1b/hadronic_on/tfsr/set0041_01_53evts.geant.root";
122 const char *fileList[] = {dstFile,xdfFile,mdcFile,0};
123 
124 void wait(int n=1) {
125  for ( int i=0; i<n*1e6; i++) { /*no-op*/ }
126 }
127 void mess(const char* c="alive") {
128  for ( int i=0; i<10; i++) { cout << c << endl; }
129 }
130 
131 
132 void StHbtExampleQQ(const Int_t nevents, const Char_t **fileList, const char*);
133 
134 
135 //==========================================================================================
136 //==========================================================================================
137 void StHbtExampleII(const Int_t nevents=9999,
138  const Char_t *path=venusPath,
139  const Char_t *file=venusFile,
140  const char* histoFile="test.histo.root")
141 {
142  const char *fileListQQ[]={0,0};
143  if (path[0]=='-') {
144  fileListQQ[0]=file;
145  } else {
146  fileListQQ[0] = gSystem->ConcatFileName(path,file);
147  }
148  StHbtExampleQQ(nevents,fileListQQ,histoFile);
149 }
150 //==========================================================================================
151 //==========================================================================================
152 void StHbtExampleQQ(const Int_t nevents, const Char_t **fileList, const char* histoFile)
153 {
154 
155 // Dynamically link needed shared libs
156 gSystem->Load("St_base");
157 gSystem->Load("StChain");
158 gSystem->Load("St_Tables");
159 gSystem->Load("StMagF");
160 gSystem->Load("StUtilities"); // new addition 22jul99
161 gSystem->Load("StTreeMaker");
162 gSystem->Load("StIOMaker");
163 gSystem->Load("StarClassLibrary");
164 #if XDFREADER
165 gSystem->Load("xdf2root");
166 gSystem->Load("St_xdfin_Maker");
167 #endif
168 gSystem->Load("StEvent");
169 gSystem->Load("StEventMaker");
170 gSystem->Load("StMcEvent");
171 gSystem->Load("StMcEventMaker");
172 gSystem->Load("StAssociationMaker");
173 gSystem->Load("StMcAnalysisMaker");
174 
175 gSystem->Load("StHbtMaker");
176 gSystem->Load("global_Tables");
177 
178 cout << " loading done " << endl;
179 
180 chain = new StChain("StChain");
181 chain->SetDebug();
182 
183 
184 StFile *setFiles= new StFile();
185 for (int ifil=0; fileList[ifil]; ifil++)
186 setFiles->AddFile(fileList[ifil]);
187 
188 
189 // ********************************
190 // Now we add Makers to the chain...
191 // ********************************
192 
193 // *************
194 // file handling
195 // *************
196 #if STEVENT
197 StIOMaker* ioMaker = new StIOMaker("IO","r",setFiles,"bfcTree");
198 ioMaker->SetDebug();
199 
200 ioMaker->SetIOMode("r");
201 ioMaker->SetDebug();
202 ioMaker->SetBranch("*",0,"0"); //deactivate all branches
203 ioMaker->SetBranch("dstBranch",0,"r"); //activate EventBranch
204 #endif
205 #if MCEVENT
206 StIOMaker* ioMaker = new StIOMaker("IO","r",setFiles,"bfcTree");
207 ioMaker->SetDebug();
208 
209 ioMaker->SetIOMode("r");
210 ioMaker->SetDebug();
211 ioMaker->SetBranch("*",0,"0"); //deactivate all branches
212 ioMaker->SetBranch("geantBranch",0,"r"); //activate EventBranch
213 #endif
214 #if ASSOCIATION
215 StIOMaker* ioMaker = new StIOMaker("IO","r",setFiles,"bfcTree");
216 ioMaker->SetDebug();
217 
218 ioMaker->SetIOMode("r");
219 ioMaker->SetDebug();
220 ioMaker->SetBranch("*",0,"0"); //deactivate all branches
221 cout << " files open" << endl;
222 ioMaker->SetBranch("geantBranch",0,"r"); //activate EventBranch
223 cout << " files open" << endl;
224 ioMaker->SetBranch("dstBranch",0,"r"); //activate EventBranch
225 cout << " files open" << endl;
226 #endif
227 #if XDFREADER
228 StIOMaker* ioMaker = new StIOMaker("IO","r",setFiles,"bfcTree");
229 ioMaker->SetDebug();
230 ioMaker->SetBranch("*",0,"r"); //activate EventBranch
231 #endif
232 
233 // ***********
234 // Event Maker
235 // ***********
236 #if STEVENT
237 StEventMaker* eventMaker = new StEventMaker("events","title");
238  eventMaker->doLoadTpcHits=0;
239  eventMaker->doLoadFtpcHits=0;
240  eventMaker->doLoadSvtHits=0;
241 cout << "Just instantiated StEventMaker... lets go StHbtMaker!" << endl;
242 #endif
243 #if MCEVENT
244 StMcEventMaker* mcEventMaker = new StMcEventMaker; // Make an instance...
245 cout << "StMcEventMaker instantiated"<<endl;
246 #endif
247 #if ASSOCIATION
248 StEventMaker* eventMaker = new StEventMaker("events","title");
249 cout << "Just instantiated StEventMaker... lets go StHbtMaker!" << endl;
250 StMcEventMaker* mcEventMaker = new StMcEventMaker; // Make an instance...
251 cout << "StMcEventMaker instantiated"<<endl;
252 StAssociationMaker* associationMaker = new StAssociationMaker; // Make an instance...
253 cout << "StAssociationMaker instantiated"<<endl;
254 #endif
255 
256 // *********
257 // Hbt Maker
258 // *********
259 
260 StHbtMaker* hbtMaker = new StHbtMaker("HBT","title");
261 cout << "StHbtMaker instantiated"<<endl;
262 //StHbtTagMaker* hbtTagMaker = new StHbtTagMaker("HBTTAG");
263 //cout << "StHbtMaker instantiated"<<endl;
264 // -------------- set up of hbt stuff ----- //
265 cout << "StHbtMaker::Init - setting up Reader and Analyses..." << endl;
266 
267 StHbtManager* TheManager = hbtMaker->HbtManager();
268 
269 
270 // ***********************
271 // setup HBT event readers
272 // ***********************
273 #if STEVENTREADER
274 // *****************************************
275 // set up StHbtMcEventReader as Event Reader
276 // *****************************************
278 Reader->SetTheEventMaker(eventMaker); // gotta tell the reader where it should read from
279 TheManager->SetEventReader(Reader);
280 cout << "READER SET UP.... " << endl;
281 #endif
282 #if MCEVENTREADER
283 // *****************************************
284 // set up StHbtMcEventReader as Event Reader
285 // *****************************************
287 Reader->SetTheMcEventMaker(mcEventMaker); // gotta tell the reader where it should read from
288 TheManager->SetEventReader(Reader);
289 cout << "READER SET UP.... " << endl;
290 #endif
291 #if ASSOCIATIONREADER
292 // ********************************************
293 // set up StHbtAssociationReader as Event Reader
294 // ********************************************
295 // StHbtAssociationReader* Reader = new StHbtAssociationReader();
297 Reader->SetTheEventMaker(eventMaker); // gotta tell the reader where it should read from
298 Reader->SetTheMcEventMaker(mcEventMaker); // gotta tell the reader where it should read from
299 Reader->SetTheAssociationMaker(associationMaker); // gotta tell the reader where it should read from
300 Reader->SetTheMcEventMaker(mcEventMaker); // gotta tell the reader where it should read from
301 TheManager->SetEventReader(Reader);
302 // define cuts for the association maker
303 StMcParameterDB* parameterDB = StMcParameterDB::instance();
304 parameterDB->setXCut(0.1);
305 parameterDB->setZCut(0.2);
306 parameterDB->setReqCommonHits(10);
307 cout << "READER SET UP.... " << endl;
308 #endif
309 
310 #if ASCIIWRITER
311  StHbtAsciiReader* asciiWriter = new StHbtAsciiReader;
312  asciiWriter->SetFileName("test1.asc");
313  TheManager->AddEventWriter(asciiWriter);
314 // set up the front loaded cuts
315  mikesEventCut* asciiFrontLoadedEventCut = new mikesEventCut;
316  asciiFrontLoadedEventCut->SetEventMult(0,100000); // selected multiplicity range
317  asciiFrontLoadedEventCut->SetVertZPos(-35.0,35.0); // selected range of vertex z-position
318  asciiWriter->SetEventCut(asciiFrontLoadedEventCut);
319  franksTrackCut* asciiFrontLoadedTrackCut = new franksTrackCut;
320  asciiFrontLoadedTrackCut->SetNSigmaPion(-1000.0,1000.0); // number of Sigma in TPC dEdx away from nominal pion dEdx
321  asciiFrontLoadedTrackCut->SetNSigmaKaon(-3.,3.); // number of Sigma in TPC dEdx away from nominal kaon dEdx
322  asciiFrontLoadedTrackCut->SetNSigmaProton(-1000.0,1000.0); // number of Sigma in TPC dEdx away from nominal proton dEdx
323  asciiFrontLoadedTrackCut->SetNHits(5,1000); // range on number of TPC hits on the track
324  asciiFrontLoadedTrackCut->SetP(0.0,5.0); // range in P
325  asciiFrontLoadedTrackCut->SetPt(0.0,5.0); // range in Pt
326  asciiFrontLoadedTrackCut->SetRapidity(-1.5,1.5); // range in rapidity
327  asciiFrontLoadedTrackCut->SetDCA(0.0,2.); // range in Distance of Closest Approach to primary vertex
328  asciiFrontLoadedTrackCut->SetCharge(0); // no cut on charge
329  asciiFrontLoadedTrackCut->SetMass(0.494); // kaon mass
330  //asciiWriter->SetTrackCut(asciiFrontLoadedTrackCut);
331  // set up a microDstWriter
332  cout << "WRITER SET UP.... " << endl;
333 #endif
334 #if ASCIIREADER
335  // set up a microDstWriter
336  StHbtAsciiReader* Reader = new StHbtAsciiReader;
337  Reader->SetFileName(*fileList);
338  //Reader->SetFileName("test1.asc");
339  TheManager->SetEventReader(Reader);
340  cout << "READER SET UP.... " << endl;
341 #endif
342 #if BINARYWRITER
343  // set up a microDstWriter
344  StHbtBinaryReader* binaryWriter = new StHbtBinaryReader(ioMaker); // retrieve filename from ioMaker
345  //StHbtBinaryReader* binaryWriter = new StHbtBinaryReader(); // specify filename
346  binaryWriter->SetFileName("test1.bin");
347  TheManager->AddEventWriter(binaryWriter);
348  // set up the front loaded cuts
349  mikesEventCut* binaryFrontLoadedEventCut = new mikesEventCut;
350  binaryFrontLoadedEventCut->SetEventMult(0,100000); // selected multiplicity range
351  binaryFrontLoadedEventCut->SetVertZPos(-35.0,35.0); // selected range of vertex z-position
352  binaryWriter->SetEventCut(binaryFrontLoadedEventCut);
353  franksTrackCut* binaryFrontLoadedTrackCut = new franksTrackCut;
354  binaryFrontLoadedTrackCut->SetNSigmaPion(-1000.0,1000.0); // number of Sigma in TPC dEdx away from nominal pion dEdx
355  binaryFrontLoadedTrackCut->SetNSigmaKaon(-3.,3.); // number of Sigma in TPC dEdx away from nominal kaon dEdx
356  binaryFrontLoadedTrackCut->SetNSigmaProton(-1000.0,1000.0); // number of Sigma in TPC dEdx away from nominal proton dEdx
357  binaryFrontLoadedTrackCut->SetNHits(5,1000); // range on number of TPC hits on the track
358  binaryFrontLoadedTrackCut->SetP(0.0,5.0); // range in P
359  binaryFrontLoadedTrackCut->SetPt(0.0,5.0); // range in Pt
360  binaryFrontLoadedTrackCut->SetRapidity(-1.5,1.5); // range in rapidity
361  binaryFrontLoadedTrackCut->SetDCA(0.0,2.); // range in Distance of Closest Approach to primary vertex
362  binaryFrontLoadedTrackCut->SetCharge(0); // no cut on charge
363  binaryFrontLoadedTrackCut->SetMass(0.494); // kaon mass
364  binaryWriter->SetTrackCut(binaryFrontLoadedTrackCut);
365  cout << "WRITER SET UP.... " << endl;
366 #endif
367 #if BINARYREADER
368  // set up a microDstWriter
369  StHbtBinaryReader* Reader = new StHbtBinaryReader(0,*fileList,0);
370  //Reader->SetFileName("/star/u2d/laue/MDC3/MicroDst/rcf0110.kaon.microDst");
371  //Reader->AddFileList("test.lis");
372  TheManager->SetEventReader(Reader);
373  cout << "READER SET UP.... " << endl;
374 #endif
375 #if XDFREADER
376  // set up StHbtXDFReader as Event Reader
377  StHbtXDFReader* Reader = new StHbtXDFReader("evgen","evgen/particle");
378  Reader->AddAcceptedParticle(-321); // kaon+
379  Reader->AddAcceptedParticle(+321); // kaon-
380  Reader->AddAcceptedMother(333); // phi
381 
382  //Reader->AddAcceptedParticle(+211); // pion+
383  //Reader->AddAcceptedParticle(+2212); // proton
384  //Reader->AddAcceptedMother(2224); // delta++
385 
386  TheManager->SetEventReader(Reader);
387  cout << "READER SET UP.... " << endl;
388 #endif
389 #if GSTARTXTREADER
391  // Reader->SetFileName("/star/rcf/pwg/hbt/GstarTextFiles/blended_events.txt");
392  Reader->SetFileName("/star/rcf/data06/evgen/auau200/hbt/default/midcentral/year_2000/hadronic_on/gen/evgen.41.evt");
393  TheManager->SetEventReader(Reader);
394 #endif
395 
396 
397  // define example particle cut and cut monitors to use in the analyses
398  // example particle cut
399  franksTrackCut* aTrackCut = new franksTrackCut; // use "frank's" particle cut object
400  aTrackCut->SetNSigmaPion(3.0,1000.0); // number of Sigma in TPC dEdx away from nominal pion dEdx
401  aTrackCut->SetNSigmaKaon(-3.,3.); // number of Sigma in TPC dEdx away from nominal kaon dEdx
402  aTrackCut->SetNSigmaProton(-1000.,-1.0); // number of Sigma in TPC dEdx away from nominal proton dEdx
403  aTrackCut->SetNHits(5,50); // range on number of TPC hits on the track
404  aTrackCut->SetP(0.2,0.8); // range in P
405  aTrackCut->SetPt(0.0,5.0); // range in Pt
406  aTrackCut->SetRapidity(-1.5,1.5); // range in rapidity
407  aTrackCut->SetDCA(0.0,2.); // range in Distance of Closest Approach to primary vertex
408  aTrackCut->SetCharge(+1); // want positive kaons
409  aTrackCut->SetMass(0.494); // kaon mass
410  // define example track cut monitor
411  trackCutMonitor_P_vs_Dedx* aDedxMoniPos = new trackCutMonitor_P_vs_Dedx(+1,"P_vs_Dedx +","Momentum (GeV/c) vs Energy loss (a.u.)",100,0.,1.2,100,0.,1e-5);
412  trackCutMonitor_P_vs_Dedx* aDedxMoniNeg = new trackCutMonitor_P_vs_Dedx(-1,"P_vs_Dedx -","Momentum (GeV/c) vs Energy loss (a.u.)",100,0.,1.2,100,0.,1e-5);
413 
414  // now, we define another analysis that runs simultaneously with the previous one.
415  // this one looks at K+K- correlations (so NONidentical particles) in invariant mass
416 
417  // ****************************************** //
418  // * franks phi analysis - by Frank Laue, OSU //
419  // ****************************************** //
420  // 0) now define an analysis...
421  // StHbtAnalysis* phiAnal = new StHbtAnalysis;
422  phiAnal = new StHbtAnalysis;
423  // 1) set the Event cuts for the analysis
424  mikesEventCut* phiEvcut = new mikesEventCut; // use "mike's" event cut object
425  phiEvcut->SetEventMult(0,100000); // selected multiplicity range
426  phiEvcut->SetVertZPos(-35.0,35.0); // selected range of vertex z-position
427  eventCutMonitor_Mult* multMoniPass = new eventCutMonitor_Mult();
428  eventCutMonitor_Mult* multMoniFail = new eventCutMonitor_Mult();
429  phiEvcut->AddCutMonitor(multMoniPass, multMoniFail);
430  phiAnal->SetEventCut(phiEvcut); // this is the event cut object for this analsys
431  // 2) set the Track (particle) cuts for the analysis
432  franksTrackCut* kaonTrkcut = new franksTrackCut( *aTrackCut ); // copy from example
433  // new particle cut moni
434  trackCutMonitor_P_vs_Dedx* dedxMoniPosPass = new trackCutMonitor_P_vs_Dedx( *aDedxMoniPos);
435  trackCutMonitor_P_vs_Dedx* dedxMoniPosFail = new trackCutMonitor_P_vs_Dedx( *aDedxMoniPos);
436  kaonTrkcut->AddCutMonitor( dedxMoniPosPass, dedxMoniPosFail);
437  // new particle cut moni
438  trackCutMonitor_P* pMoni1 = new trackCutMonitor_P;
439  trackCutMonitor_P* pMoni2 = new trackCutMonitor_P("P","momentum (GeV/c)",20, 0., 4.);
440  kaonTrkcut->AddCutMonitor( pMoni1, pMoni2);
441  // new particle cut moni
442  trackCutMonitor_Pt* ptMoni1 = new trackCutMonitor_Pt;
443  trackCutMonitor_Pt* ptMoni2 = new trackCutMonitor_Pt;
444  kaonTrkcut->AddCutMonitor( ptMoni1, ptMoni2);
445  // new particle cut moni
448  kaonTrkcut->AddCutMonitor( yptMoni1, yptMoni2);
449  // new particle cut moni
450  trackCutMonitor_DCA* DCAMoni1 = new trackCutMonitor_DCA("DCA","DCA (cm)",100,0.,10.);
451  trackCutMonitor_DCA* DCAMoni2 = new trackCutMonitor_DCA("DCA","DCA (cm)",100,0.,10.);
452  kaonTrkcut->AddCutMonitor( DCAMoni1, DCAMoni2 );
453  phiAnal->SetFirstParticleCut(kaonTrkcut); // this is the track cut for the "first" particle
454 
455  // copy second particle cut from first particle cut
456  franksTrackCut* antikaonTrkcut = new franksTrackCut( *((franksTrackCut*)phiAnal->FirstParticleCut()) );
457  antikaonTrkcut->SetCharge(-1);
458  phiAnal->SetSecondParticleCut(antikaonTrkcut); // this is the track cut for the "first" particle
459  // new particle cut
460  trackCutMonitor_P_vs_Dedx* dedxMoniNegPass = new trackCutMonitor_P_vs_Dedx( *aDedxMoniNeg);
461  trackCutMonitor_P_vs_Dedx* dedxMoniNegFail = new trackCutMonitor_P_vs_Dedx( *aDedxMoniNeg);
462  antikaonTrkcut->AddCutMonitor( dedxMoniNegPass, dedxMoniNegFail);
463 
464  // 3) set the Pair cuts for the analysis
465  mikesPairCut* phiPairCut = new mikesPairCut; // use "mike's" pair cut object
466  //franksPairCut* phiPairCut = new franksPairCut; // use "frank's" pair cut object
467  phiAnal->SetPairCut(phiPairCut); // this is the pair cut for this analysis
468  // 4) set the number of events to mix (per event)
469  phiAnal->SetNumEventsToMix(3);
470  // ********************************************************************
471  // 5) now set up the correlation functions that this analysis will make
472  // ********************************************************************
473  // define example Minv correlation function
474  MinvCorrFctn* MinvCF = new MinvCorrFctn("Minv",100,0.95,1.25);
475  phiAnal->AddCorrFctn(MinvCF); // adds the just-defined correlation function to the analysis
476 
477  MinvCorrFctnArmenteros* MinvCFArm = new MinvCorrFctnArmenteros("ArmenterosPodolanski",200,-1.,1.,300,0.,.3);
478  phiAnal->AddCorrFctn(MinvCFArm); // adds the just-defined correlation function to the analysis
479 
480  MinvCFY_vs_Pt = new MinvCorrFctnY_vs_Pt("MinvCF dn/d(Y_vs_Pt)",31, -1.5, 1.5, 20, 0., 2. ); // defines a Minv function
481  phiAnal->AddCorrFctn(MinvCFY_vs_Pt); // adds the just-defined correlation function to the analysis
482 
483  MinvCFM_vs_Pt = new MinvCorrFctnM_vs_Pt("MinvCM dn/d(M_vs_Pt)",100,0.98, 1.18, 40, 0., 2. ); // defines a Minv function
484  phiAnal->AddCorrFctn(MinvCFM_vs_Pt); // adds the just-defined correlation function to the analysis
485 
486  /*
487  MinvCFM_vs_P = new MinvCorrFctnM_vs_P("Mass (GeV/c^2) vs Momentum (GeV/c)",50,0.98,1.08,20,0.,2.);
488  phiAnal->AddCorrFctn(MinvCFM_vs_P); // adds the just-defined correlation function to the analysis
489 
490  p1_vs_p2CF = new p1_vs_p2CorrFctn("K+ momentum (GeV/c^) vs K- momentum (GeV/c)",20,0.,2., 20,0.,2.);
491  phiAnal->AddCorrFctn(p1_vs_p2CF); // adds the just-defined correlation function to the analysis
492 
493  YCF = new YCorrFctn("Rapidity",20,-2.,2.);
494  phiAnal->AddCorrFctn(YCF); // adds the just-defined correlation function to the analysis
495 
496  pDiffCF = new pDiffCorrFctn(" p1-p2 (GeV/c)",100,0.,.1);
497  phiAnal->AddCorrFctn(pDiffCF); // adds the just-defined correlation function to the analysis
498 
499  p_vs_angleCF = new p_vs_angleCorrFctn("K+/- momentum (GeV/c^) vs opening angle",20,0.,2., 20,0.,2.*3.1415927);
500  phiAnal->AddCorrFctn(p_vs_angleCF); // adds the just-defined correlation function to the analysis
501 
502  phiMeas_vs_phiCalcCF = new phiMeas_vs_phiCalcCorrFctn(" phiMeas vs phiCalc",20,-1.,1.,20,-1.,1.);
503  phiAnal->AddCorrFctn(phiMeas_vs_phiCalcCF); // adds the just-defined correlation function to the analysis
504  */
505 
506  // now add as many more correlation functions to the Analysis as you like..
507  // 6) add the Analysis to the AnalysisCollection
508  TheManager->AddAnalysis(phiAnal);
509 
510 
511  // ***************************************************** //
512  // * franks phi ---> e+ e- analysis - by Frank Laue, OSU //
513  // ***************************************************** //
514  // 0) now define an analysis...
515  StHbtAnalysis* eeAnal = new StHbtLikeSignAnalysis;
516  // 1) set the Event cuts for the analysis
517  mikesEventCut* eeEvcut = new mikesEventCut; // use "mike's" event cut object
518  eeEvcut->SetEventMult(0,100000); // selected multiplicity range
519  eeEvcut->SetVertZPos(-35.0,35.0); // selected range of vertex z-position
520  eeAnal->SetEventCut(eeEvcut); // this is the event cut object for this analsys
521  // 2) set the Track (particle) cuts for the analysis
522  franksTrackCut* electronTrkCut = new franksTrackCut( *aTrackCut ); // copy from example
523  electronTrkCut->SetNSigmaElectron(-1e2,+1e2); // number of Sigma in TPC dEdx away from nominal pion dEdx
524  electronTrkCut->SetNSigmaPion(-1e5,+1e5); // number of Sigma in TPC dEdx away from nominal pion dEdx
525  electronTrkCut->SetNSigmaKaon(-1e5,+1e5); // number of Sigma in TPC dEdx away from nominal kaon dEdx
526  electronTrkCut->SetNSigmaProton(-1e5,+1e5); // number of Sigma in TPC dEdx away from nominal proton dEdx
527 #ifdef MCEVENT
528  electronTrkCut->SetNSigmaElectron(-3,+3); // number of Sigma in TPC dEdx away from nominal pion dEdx
529  electronTrkCut->SetNSigmaPion(-1e5,+1e5); // number of Sigma in TPC dEdx away from nominal pion dEdx
530  electronTrkCut->SetNSigmaKaon(-1e5,+1e5); // number of Sigma in TPC dEdx away from nominal kaon dEdx
531  electronTrkCut->SetNSigmaProton(-1e5,+1e5); // number of Sigma in TPC dEdx away from nominal proton dEdx
532 #endif
533  electronTrkCut->SetNHits(5,100); // range on number of TPC hits on the track
534  electronTrkCut->SetP(0.,2.0); // range in P
535  electronTrkCut->SetPt(0.0,2.0); // range in Pt
536  electronTrkCut->SetRapidity(-1.5,1.5); // range in rapidity
537  electronTrkCut->SetDCA(2.0,50); // range in Distance of Closest Approach to primary vertex
538  electronTrkCut->SetCharge(+1); // want positive kaons
539  electronTrkCut->SetMass(0.000511); // kaon mass
540  // new particle cut moni
541  trackCutMonitor_P_vs_Dedx* dedxMoniPosPass = new trackCutMonitor_P_vs_Dedx(+1,"P_vs_Dedx +","Momentum (GeV/c) vs Energy loss (a.u.)",100,0.,5.0,100,0.,1e-5);
542  trackCutMonitor_P_vs_Dedx* dedxMoniPosFail = new trackCutMonitor_P_vs_Dedx(-1,"P_vs_Dedx -","Momentum (GeV/c) vs Energy loss (a.u.)",100,0.,5.0,100,0.,1e-5);
543  electronTrkCut->AddCutMonitor( dedxMoniPosPass, dedxMoniPosFail);
544  // new particle cut moni
545  trackCutMonitor_P* pMoni1 = new trackCutMonitor_P;
546  trackCutMonitor_P* pMoni2 = new trackCutMonitor_P("P","momentum (GeV/c)",20, 0., 4.);
547  electronTrkCut->AddCutMonitor( pMoni1, pMoni2);
548  // new particle cut moni
549  trackCutMonitor_Pt* ptMoni1 = new trackCutMonitor_Pt;
550  trackCutMonitor_Pt* ptMoni2 = new trackCutMonitor_Pt;
551  electronTrkCut->AddCutMonitor( ptMoni1, ptMoni2);
552  // new particle cut moni
555  electronTrkCut->AddCutMonitor( yptMoni1, yptMoni2);
556  // new particle cut moni
557  trackCutMonitor_DCA* DCAMoni1 = new trackCutMonitor_DCA("DCA","DCA (cm)",100,0.,10.);
558  trackCutMonitor_DCA* DCAMoni2 = new trackCutMonitor_DCA("DCA","DCA (cm)",100,0.,10.);
559  electronTrkCut->AddCutMonitor( DCAMoni1, DCAMoni2 );
560  eeAnal->SetFirstParticleCut(electronTrkCut); // this is the track cut for the "first" particle
561 
562  // copy second particle cut from first particle cut
563  franksTrackCut* antiElectronTrkCut = new franksTrackCut( *((franksTrackCut*)eeAnal->FirstParticleCut()) );
564  antiElectronTrkCut->SetCharge(-1);
565  eeAnal->SetSecondParticleCut(antiElectronTrkCut); // this is the track cut for the "first" particle
566  // new particle cut
567  trackCutMonitor_P_vs_Dedx* dedxMoniNegPass = new trackCutMonitor_P_vs_Dedx( *aDedxMoniNeg);
568  trackCutMonitor_P_vs_Dedx* dedxMoniNegFail = new trackCutMonitor_P_vs_Dedx( *aDedxMoniNeg);
569  antiElectronTrkCut->AddCutMonitor( dedxMoniNegPass, dedxMoniNegFail);
570 
571  // 3) set the Pair cuts for the analysis
572  mikesPairCut* eePairCut = new mikesPairCut; // use "mike's" pair cut object
573  //franksPairCut* eePairCut = new franksPairCut; // use "frank's" pair cut object
574  eeAnal->SetPairCut(eePairCut); // this is the pair cut for this analysis
575  // 4) set the number of events to mix (per event)
576  eeAnal->SetNumEventsToMix(2);
577  // ********************************************************************
578  // 5) now set up the correlation functions that this analysis will make
579  // ********************************************************************
580  // define example MinvLikeSign correlation function
581  MinvLikeSignCorrFctn* MinvLikeSignCFconv = new MinvLikeSignCorrFctn("MinvLikeSign",500,0.,.25);
582  eeAnal->AddCorrFctn(MinvLikeSignCFconv); // adds the just-defined correlation function to the analysis
583  MinvLikeSignCorrFctn* MinvLikeSignCF = new MinvLikeSignCorrFctn("MinvLikeSign",100,0.95,1.15);
584  eeAnal->AddCorrFctn(MinvLikeSignCF); // adds the just-defined correlation function to the analysis
585 
586  // now add as many more correlation functions to the Analysis as you like..
587  // 6) add the Analysis to the AnalysisCollection
588  TheManager->AddAnalysis(eeAnal);
589 
590  // ********************************************* //
591  // * franks lambda analysis - by Frank Laue, OSU //
592  // ********************************************* //
593  //StHbtAnalysis* lambdaAnal = new StHbtAnalysis( *phiAnal );
594  lambdaAnal = new StHbtAnalysis( *phiAnal );
595  ((franksTrackCut*)lambdaAnal->FirstParticleCut())->SetNSigmaPion(-1000.0,1000.0); // proton
596  ((franksTrackCut*)lambdaAnal->FirstParticleCut())->SetNSigmaKaon(-1000.0,1000.0);
597  ((franksTrackCut*)lambdaAnal->FirstParticleCut())->SetNSigmaProton(-3.0,3.0);
598  ((franksTrackCut*)lambdaAnal->FirstParticleCut())->SetP(0.,2.);
599  ((franksTrackCut*)lambdaAnal->FirstParticleCut())->SetPt(0.,2.);
600  ((franksTrackCut*)lambdaAnal->FirstParticleCut())->SetCharge(-1);
601  ((franksTrackCut*)lambdaAnal->FirstParticleCut())->SetMass(0.938);
602  ((franksTrackCut*)lambdaAnal->FirstParticleCut())->SetDCA(8.0,16.);
603  ((franksTrackCut*)lambdaAnal->SecondParticleCut())->SetNSigmaPion(-3.0,3.0); //pion
604  ((franksTrackCut*)lambdaAnal->SecondParticleCut())->SetNSigmaKaon(-1000.0,1000.0);
605  ((franksTrackCut*)lambdaAnal->SecondParticleCut())->SetNSigmaProton(-1000.0,1000.0);
606  ((franksTrackCut*)lambdaAnal->SecondParticleCut())->SetP(0.,2.);
607  ((franksTrackCut*)lambdaAnal->SecondParticleCut())->SetPt(0.,2.);
608  ((franksTrackCut*)lambdaAnal->SecondParticleCut())->SetCharge(-1);
609  ((franksTrackCut*)lambdaAnal->SecondParticleCut())->SetMass(0.139);
610  ((franksTrackCut*)lambdaAnal->SecondParticleCut())->SetDCA(8.0,16.);
611  //TheManager->AddAnalysis(lambdaAnal);
612 
613  // ********************************************* //
614  // * franks lambda analysis - by Frank Laue, OSU //
615  // ********************************************* //
616  StHbtAnalysis* lambdaAnal2 = new StHbtAnalysis( *lambdaAnal );
617  delete ((franksPairCut*)lambdaAnal2->PairCut());
618  franksPairCut* lambdaAnal2PairCut = new franksPairCut; // use "frank's" pair cut object
619  lambdaAnal2->SetPairCut(lambdaAnal2PairCut); // this is the pair cut for this analysis
620  //TheManager->AddAnalysis(lambdaAnal2);
621 
622  // ****************************************** //
623  // * franks rho analysis - by Frank Laue, OSU //
624  // ****************************************** //
625  // 0) now define an analysis...
626  // StHbtAnalysis* rhoAnal = new StHbtAnalysis;
627  rhoAnal = new StHbtAnalysis;
628  // 1) set the Event cuts for the analysis
629  mikesEventCut* rhoEvcut = new mikesEventCut; // use "mike's" event cut object
630  rhoEvcut->SetEventMult(0,100000); // selected multiplicity range
631  rhoEvcut->SetVertZPos(-35.0,35.0); // selected range of vertex z-position
632  rhoAnal->SetEventCut(rhoEvcut); // this is the event cut object for this analsys
633  // 2) set the Track (particle) cuts for the analysis
634  franksTrackCut* piPosTrkcut = new franksTrackCut; // use "frank's" particle cut object
635  piPosTrkcut->SetNSigmaPion(-2.0,+2.0); // number of Sigma in TPC dEdx away from nominal pion dEdx
636  piPosTrkcut->SetNSigmaKaon(-1000.,+1000.); // number of Sigma in TPC dEdx away from nominal kaon dEdx
637  piPosTrkcut->SetNSigmaProton(-1000.,+1000.0); // number of Sigma in TPC dEdx away from nominal proton dEdx
638  piPosTrkcut->SetNHits(10,50); // range on number of TPC hits on the track
639  piPosTrkcut->SetP(0.1,0.5); // range in P
640  piPosTrkcut->SetPt(0.1,2.0); // range in Pt
641  piPosTrkcut->SetRapidity(-1.5,1.5); // range in rapidity
642  piPosTrkcut->SetDCA(0.0,2.); // range in Distance of Closest Approach to primary vertex
643  piPosTrkcut->SetCharge(+1); // want positive kaons
644  piPosTrkcut->SetMass(0.139); // kaon mass
645  rhoAnal->SetFirstParticleCut(piPosTrkcut); // this is the track cut for the "first" particle
646  // copy second particle cut from first particle cut
647  franksTrackCut* piNegTrkCut = new franksTrackCut( *((franksTrackCut*)rhoAnal->FirstParticleCut()) );
648  piNegTrkCut->SetCharge(-1);
649  rhoAnal->SetSecondParticleCut(piNegTrkCut); // this is the track cut for the "first" particle
650  // 3) set the Pair cuts for the analysis
651  mikesPairCut* rhoPairCut = new mikesPairCut; // use "frank's" pair cut object
652  //franksPairCut* rhoPairCut = new franksPairCut; // use "frank's" pair cut object
653  //rhoPairCut->SetPDiff(0.,1000.); // difference of pairs in momentum
654  // rhoPairCut->SetAngle(0.,180.); // opening angle
655  rhoAnal->SetPairCut(rhoPairCut); // this is the pair cut for this analysis
656  // 4) set the number of events to mix (per event)
657  rhoAnal->SetNumEventsToMix(10);
658  // ********************************************************************
659  // 5) now set up the correlation functions that this analysis will make
660  // ********************************************************************
661  // define example Minv correlation function
662  MinvCorrFctn* MinvCFrho = new MinvCorrFctn("Minv",100,0.4,.9);
663  // rhoAnal->AddCorrFctn(MinvCFrho); // adds the just-defined correlation function to the analysis
664 
665  // now add as many more correlation functions to the Analysis as you like..
666  // 6) add the Analysis to the AnalysisCollection
667  //TheManager->AddAnalysis(rhoAnal);
668 
669  // ******************************************** //
670  // * franks delta analysis - by Frank Laue, OSU //
671  // ******************************************** //
672  // 0) now define an analysis...
673  //StHbtAnalysis* deltaAnal = new StHbtAnalysis;
674  deltaAnal = new StHbtAnalysis;
675  // 1) set the Event cuts for the analysis
676  mikesEventCut* deltaEvcut = new mikesEventCut; // use "mike's" event cut object
677  deltaEvcut->SetEventMult(0,100000); // selected multiplicity range
678  deltaEvcut->SetVertZPos(-35.0,35.0); // selected range of vertex z-position
679  deltaAnal->SetEventCut(deltaEvcut); // this is the event cut object for this analsys
680  // 2) set the Track (particle) cuts for the analysis
681  franksTrackCut* pionTrkcut = new franksTrackCut; // use "frank's" particle cut object
682  pionTrkcut->SetNSigmaPion(-2.0,+2.0); // number of Sigma in TPC dEdx away from nominal pion dEdx
683  pionTrkcut->SetNSigmaKaon(-1000.,+1000.); // number of Sigma in TPC dEdx away from nominal kaon dEdx
684  pionTrkcut->SetNSigmaProton(-1000.,+1000.0); // number of Sigma in TPC dEdx away from nominal proton dEdx
685  pionTrkcut->SetNHits(10,50); // range on number of TPC hits on the track
686  pionTrkcut->SetP(0.1,0.5); // range in P
687  pionTrkcut->SetPt(0.1,2.0); // range in Pt
688  pionTrkcut->SetRapidity(-1.5,1.5); // range in rapidity
689  pionTrkcut->SetDCA(0.0,2.); // range in Distance of Closest Approach to primary vertex
690  pionTrkcut->SetCharge(+1); // want positive kaons
691  pionTrkcut->SetMass(0.139); // kaon mass
692  deltaAnal->SetFirstParticleCut(pionTrkcut); // this is the track cut for the "first" particle
693  // copy second particle cut from first particle cut
694  franksTrackCut* protonTrkcut = new franksTrackCut( *((franksTrackCut*)deltaAnal->FirstParticleCut()) );
695  protonTrkcut->SetNSigmaPion(-1000.0,+1000.0); // number of Sigma in TPC dEdx away from nominal pion dEdx
696  protonTrkcut->SetNSigmaProton(-2.,+2.0); // number of Sigma in TPC dEdx away from nominal proton dEdx
697  protonTrkcut->SetMass(0.938); // kaon mass
698  deltaAnal->SetSecondParticleCut(protonTrkcut); // this is the track cut for the "first" particle
699  // 3) set the Pair cuts for the analysis
700  mikesPairCut* deltaPairCut = new mikesPairCut; // use "frank's" pair cut object
701  //franksPairCut* deltaPairCut = new franksPairCut; // use "frank's" pair cut object
702  //deltaPairCut->SetPDiff(0.,1000.); // difference of pairs in momentum
703  // deltaPairCut->SetAngle(0.,180.); // opening angle
704  deltaAnal->SetPairCut(deltaPairCut); // this is the pair cut for this analysis
705  // 4) set the number of events to mix (per event)
706  deltaAnal->SetNumEventsToMix(10);
707  // ********************************************************************
708  // 5) now set up the correlation functions that this analysis will make
709  // ********************************************************************
710  // define example Minv correlation function
711  MinvCorrFctn* MinvCFdelta = new MinvCorrFctn("Minv",100,1.0,1.4);
712  deltaAnal->AddCorrFctn(MinvCFdelta); // adds the just-defined correlation function to the analysis
713 
714  // now add as many more correlation functions to the Analysis as you like..
715  // 6) add the Analysis to the AnalysisCollection
716  //TheManager->AddAnalysis(deltaAnal);
717 
718 
719 
720  // ------------------ end of setting up hbt stuff ------------------ //
721 
722  chain->Init(); // This should call the Init() method in ALL makers
723  chain->PrintInfo();
724 
725  // exit();
726  for (Int_t iev=0;iev<nevents; iev++) {
727  cout << "StHbtExample -- Working on eventNumber " << iev << endl;
728  chain->Clear();
729  int iret = chain->Make(iev); // This should call the Make() method in ALL makers
730  if (iret) {
731  cout << "Bad return code!" << endl;
732  break;
733  }
734  } // Event Loop
735  chain->Finish(); // This should call the Finish() method in ALL makers
736 
737  cout << " End of Analysis " << endl;
738  // write some histos into a file
739  TFile histoOutput(histoFile,"recreate");
740  MinvCF->Numerator()->Write();
741  MinvCF->Denominator()->Write();
742  MinvCF->Difference()->Write();
743  MinvCFY_vs_Pt->Numerator()->Write();
744  MinvCFY_vs_Pt->Denominator()->Write();
745  MinvCFY_vs_Pt->Difference()->Write();
746  MinvCFM_vs_Pt->Numerator()->Write();
747  MinvCFM_vs_Pt->Denominator()->Write();
748  MinvCFM_vs_Pt->Difference()->Write();
749  MinvCFArm->Numerator()->Write();
750  MinvCFArm->Denominator()->Write();
751  MinvCFArm->Difference()->Write();
752 ((eventCutMonitor_Mult*)phiAnal->EventCut()->PassMonitor(0))->Histo()->Write();
753  histoOutput->Write();
754  histoOutput->Close();
755 }
756 
757 
758 
Definition: StTree.h:125
virtual void SetIOMode(Option_t *iomode="w")
number of transactions
Definition: StIOInterFace.h:35
virtual void Clear(Option_t *option="")
User defined functions.
Definition: StChain.cxx:77
virtual Int_t Finish()
Definition: StChain.cxx:85
Filling of all StMcEvent classes from g2t tables Transform all the data in the g2t tables into the co...
virtual Int_t Make()
Definition: StChain.cxx:110