00001 // $Id: StHbtExample.C,v 1.11 2006/08/15 21:43:03 jeromel Exp $ 00002 // $Log: StHbtExample.C,v $ 00003 // Revision 1.11 2006/08/15 21:43:03 jeromel 00004 // Fix rhic -> rhic.bnl.gov 00005 // 00006 // Revision 1.10 2005/08/31 15:05:08 fisyak 00007 // Add dependence StMagF vs StarMagField 00008 // 00009 // Revision 1.9 2000/04/13 21:46:22 kathy 00010 // remove loading of libtpc_Tables since l3Track table is now dst_track type from global 00011 // 00012 // Revision 1.8 2000/04/12 17:39:02 kathy 00013 // change to only load table libraries needed: lib*_Tables instead of all tables: St_Tables 00014 // 00015 // Revision 1.7 2000/03/20 17:50:40 kathy 00016 // fix all macros so that they set all branches on that are needed - otherwise won't work with soft links 00017 // 00018 // Revision 1.6 2000/01/19 15:52:46 kathy 00019 // change default input file to be the one in /afs/rhic.bnl.gov/star/data/samples 00020 // 00021 // Revision 1.5 1999/09/24 00:00:34 lisa 00022 // HBT example macro updated to use TrackCuts not ParticleCuts 00023 // 00024 // Revision 1.4 1999/07/27 21:09:18 lisa 00025 // do event loop with goto because Cint pukes on for loop 00026 // 00027 // Revision 1.3 1999/07/13 02:31:39 lisa 00028 // Update StHbtExample macro for new StHbtMaker 00029 // 00030 // Revision 1.2 1999/07/13 01:13:02 kathy 00031 // moved rch.C to obsolete, put in id,log,owner into HbtExample, removed loading of StRootEvent and changed default input file in bfcread.C and Example_readdst_qa_tables.C 00032 // 00033 00034 //====================================================================== 00035 // owner: Mike Lisa 00036 // what it does: 00037 // runs the StHbt framework. Runs two simultaneous correlation 00038 // analyses. The first one is pi-pi HBT, building two simultaneous 00039 // correlation functions. The second in K+K-, building an invariant 00040 // mass spectrum with background subtraction 00041 //======================================================================= 00042 00043 class StChain; 00044 StChain *chain=0; 00045 00046 // keep pointers to Correlation Functions global, so you can have access to them... 00047 class QinvCorrFctn; 00048 QinvCorrFctn* QinvCF; 00049 class QvecCorrFctn; 00050 QvecCorrFctn* QvecCF; 00051 class MinvCorrFctn; 00052 MinvCorrFctn* MinvCF; 00053 00054 void StHbtExample(Int_t nevents=1, 00055 const char *MainFile= 00056 "/afs/rhic.bnl.gov/star/data/samples/gstar.dst.root") 00057 { 00058 00059 // Dynamically link needed shared libs 00060 gSystem->Load("St_base"); 00061 gSystem->Load("StChain"); 00062 gSystem->Load("libglobal_Tables"); 00063 gSystem->Load("libsim_Tables"); 00064 gSystem->Load("libgen_Tables"); 00065 gSystem->Load("StUtilities"); // new addition 22jul99 00066 gSystem->Load("StAnalysisUtilities"); // needed by V0dstMaker 00067 gSystem->Load("StarMagField"); 00068 gSystem->Load("StMagF"); 00069 gSystem->Load("StIOMaker"); 00070 gSystem->Load("StarClassLibrary"); 00071 gSystem->Load("StEvent"); 00072 gSystem->Load("StEventMaker"); 00073 gSystem->Load("StHbtMaker"); 00074 gSystem->Load("StV0MiniDstMaker"); 00075 00076 cout << "Dynamic loading done" << endl; 00077 00078 chain = new StChain("StChain"); 00079 chain->SetDebug(); 00080 00081 00082 // Now we add Makers to the chain... 00083 00084 StIOMaker* ioMaker = new StIOMaker("IO","r",MainFile,"bfcTree"); 00085 ioMaker->SetDebug(); 00086 00087 ioMaker->SetIOMode("r"); 00088 ioMaker->SetDebug(); 00089 ioMaker->SetBranch("*",0,"0"); //deactivate all branches 00090 ioMaker->SetBranch("dstBranch",0,"r"); //activate EventBranch 00091 ioMaker->SetBranch("runcoBranch",0,"r"); //activate runcoBranch 00092 00093 00094 StEventMaker* eventMaker = new StEventMaker("events","title"); 00095 cout << "Just instantiated StEventMaker... lets go StHbtMaker!" << endl; 00096 00097 // UNCOMMENT THIS NEXT PART OUT IF YOU WANT V0's 00098 //StV0MiniDstMaker* v0dst = new StV0MiniDstMaker("v0dst"); 00099 //cout << "Just instantiated StV0MiniDstMaker... lets go StHbt!" << endl; 00100 //v0dst.SetV0VertexType(); //Set v0MiniDstMaker to find v0s not Xis 00101 //v0dst.SetOutputFile("muv0dst.root"); // Set V0MiniDStMaker output file 00102 00103 00104 00105 StHbtMaker* hbtMaker = new StHbtMaker("HBT","title"); 00106 cout << "StHbtMaker instantiated"<<endl; 00107 00108 00109 00110 /* -------------- set up of hbt stuff ----- */ 00111 cout << "StHbtMaker::Init - setting up Reader and Analyses..." << endl; 00112 00113 StHbtManager* TheManager = hbtMaker->HbtManager(); 00114 00115 // here, we instantiate the appropriate StHbtEventReader 00116 // for STAR analyses in root4star, we instantiate StStandardHbtEventReader 00117 StStandardHbtEventReader* Reader = new StStandardHbtEventReader; 00118 Reader->SetTheEventMaker(eventMaker); // gotta tell the reader where it should read from 00119 00120 // UNCOMMENT THIS NEXT LINE OUT IF YOU WANT V0's 00121 // Reader->SetTheV0Maker(v0dst); //Gotta tell the reader where to read the v0 stuff from 00122 00123 // here would be the palce to plug in any "front-loaded" Event or Particle Cuts... 00124 TheManager->SetEventReader(Reader); 00125 00126 cout << "READER SET UP.... " << endl; 00127 00128 // Hey kids! Let's make a microDST! 00129 // in StHbt we do this by instantiating and plugging in a StHbtEventReader as a writer! 00130 // the particular StHbtEventReader that we will use will write (and read) ASCII files 00131 // 00132 // StHbtAsciiReader* Writer = new StHbtAsciiReader; 00133 // Writer->SetFileName("FirstMicroDst.asc"); 00134 // TheManager->SetEventWriter(Writer); 00135 // cout << "WRITER SET UP.... " << endl; 00136 00137 // 0) now define an analysis... 00138 StHbtAnalysis* anal = new StHbtAnalysis; 00139 // 1) set the Event cuts for the analysis 00140 mikesEventCut* evcut = new mikesEventCut; // use "mike's" event cut object 00141 evcut->SetEventMult(0,10000); // selected multiplicity range 00142 evcut->SetVertZPos(-35.0,35.0); // selected range of vertex z-position 00143 anal->SetEventCut(evcut); // this is the event cut object for this analsys 00144 // 2) set the Track (particle) cuts for the analysis 00145 mikesTrackCut* trkcut = new mikesTrackCut; // use "mike's" particle cut object 00146 trkcut->SetNSigmaPion(-1.5,1.5); // number of Sigma in TPC dEdx away from nominal pion dEdx 00147 trkcut->SetNSigmaKaon(-1000.0,1000.0); // number of Sigma in TPC dEdx away from nominal kaon dEdx 00148 trkcut->SetNSigmaProton(-1000.0,1000.0); // number of Sigma in TPC dEdx away from nominal proton dEdx 00149 trkcut->SetNHits(5,50); // range on number of TPC hits on the track 00150 trkcut->SetPt(0.1,1.0); // range in Pt 00151 trkcut->SetRapidity(-1.0,1.0); // range in rapidity 00152 trkcut->SetDCA(0.0,0.5); // range in Distance of Closest Approach to primary vertex 00153 trkcut->SetCharge(-1); // want negative pions 00154 trkcut->SetMass(0.139); // pion mass 00155 anal->SetFirstParticleCut(trkcut); // this is the track cut for the "first" particle 00156 anal->SetSecondParticleCut(trkcut); // NOTE - it is also for the "second" particle -- i.e. identical particle HBT 00157 // 3) set the Pair cuts for the analysis 00158 mikesPairCut* paircut = new mikesPairCut; // use "mike's" pair cut object 00159 anal->SetPairCut(paircut); // this is the pair cut for this analysis 00160 // 4) set the number of events to mix (per event) 00161 anal->SetNumEventsToMix(5); 00162 // 5) now set up the correlation functions that this analysis will make 00163 // this particular analysis will have two: the first is a Q-invariant correlation function 00164 QinvCF = new QinvCorrFctn("mikesQinvCF",50,0.0,0.2); // defines a Qinv correlation function 00165 anal->AddCorrFctn(QinvCF); // adds the just-defined correlation function to the analysis 00166 // for this analysis, we will also (simultaneously) build a Q-vector correlation function 00167 QvecCF = new QvecCorrFctn("randysQvecCF",50,0.0,0.2); 00168 anal->AddCorrFctn(QvecCF); // adds the just-defined correlation function to the analysis 00169 00170 // now add as many more correlation functions to the Analysis as you like.. 00171 00172 // 6) add the Analysis to the AnalysisCollection 00173 TheManager->AddAnalysis(anal); 00174 00175 // now, we define another analysis that runs simultaneously with the previous one. 00176 // this one looks at K+K- correlations (so NONidentical particles) in invariant mass 00177 00178 /* ****************************************** */ 00179 /* * franks phi analysis - by Frank Laue, OSU */ 00180 /* ****************************************** */ 00181 // 0) now define an analysis... 00182 StHbtAnalysis* phiAnal = new StHbtAnalysis; 00183 // 1) set the Event cuts for the analysis 00184 mikesEventCut* phiEvcut = new mikesEventCut; // use "mike's" event cut object 00185 phiEvcut->SetEventMult(0,10000); // selected multiplicity range 00186 phiEvcut->SetVertZPos(-35.0,35.0); // selected range of vertex z-position 00187 phiAnal->SetEventCut(phiEvcut); // this is the event cut object for this analsys 00188 // 2) set the Track (particle) cuts for the analysis 00189 mikesTrackCut* kaonTrkcut = new mikesTrackCut; // use "mike's" particle cut object 00190 kaonTrkcut->SetNSigmaPion(3,1000.0); // number of Sigma in TPC dEdx away from nominal pion dEdx 00191 kaonTrkcut->SetNSigmaKaon(-3.0,3.0); // number of Sigma in TPC dEdx away from nominal kaon dEdx 00192 kaonTrkcut->SetNSigmaProton(-1000.,-1.0); // number of Sigma in TPC dEdx away from nominal proton dEdx 00193 kaonTrkcut->SetNHits(0,50); // range on number of TPC hits on the track 00194 kaonTrkcut->SetPt(0.1,2.0); // range in Pt 00195 kaonTrkcut->SetRapidity(-2.0,2.0); // range in rapidity 00196 kaonTrkcut->SetDCA(0.0,0.5); // range in Distance of Closest Approach to primary vertex 00197 kaonTrkcut->SetCharge(+1); // want positive kaons 00198 kaonTrkcut->SetMass(0.494); // kaon mass 00199 phiAnal->SetFirstParticleCut(kaonTrkcut); // this is the track cut for the "first" particle 00200 mikesTrackCut* antikaonTrkcut = new mikesTrackCut; // use "mike's" particle cut object 00201 antikaonTrkcut->SetNSigmaPion(3.0,1000.0); // number of Sigma in TPC dEdx away from nominal pion dEdx 00202 antikaonTrkcut->SetNSigmaKaon(-3.0,3.0); // number of Sigma in TPC dEdx away from nominal kaon dEdx 00203 antikaonTrkcut->SetNSigmaProton(-1000.0,-1.0); // number of Sigma in TPC dEdx away from nominal proton dEdx 00204 antikaonTrkcut->SetNHits(0,50); // range on number of TPC hits on the track 00205 antikaonTrkcut->SetPt(0.1,2.0); // range in Pt 00206 antikaonTrkcut->SetRapidity(-2.0,2.0); // range in rapidity 00207 antikaonTrkcut->SetDCA(0.0,0.5); // range in Distance of Closest Approach to primary vertex 00208 antikaonTrkcut->SetCharge(-1); // want negative kaons 00209 antikaonTrkcut->SetMass(0.494); // kaon mass 00210 phiAnal->SetSecondParticleCut(antikaonTrkcut); // this is the track cut for the "second" particle 00211 // 3) set the Pair cuts for the analysis 00212 mikesPairCut* phiPaircut = new mikesPairCut; // use "mike's" pair cut object 00213 phiAnal->SetPairCut(phiPaircut); // this is the pair cut for this analysis 00214 // 4) set the number of events to mix (per event) 00215 phiAnal->SetNumEventsToMix(5); 00216 // 5) now set up the correlation functions that this analysis will make 00217 MinvCF = new MinvCorrFctn("franksMinvCF",100,0.98,1.18); // defines a Minv correlation function 00218 phiAnal->AddCorrFctn(MinvCF); // adds the just-defined correlation function to the analysis 00219 // now add as many more correlation functions to the Analysis as you like.. 00220 // 6) add the Analysis to the AnalysisCollection 00221 TheManager->AddAnalysis(phiAnal); 00222 00223 00224 00225 /* ------------------ end of setting up hbt stuff ------------------ */ 00226 00227 00228 // now execute the chain member functions 00229 00230 if (chain->Init()){ // This should call the Init() method in ALL makers 00231 cout << "Initialization failed \n"; 00232 goto TheEnd; 00233 } 00234 chain->PrintInfo(); 00235 00236 00237 // Event loop 00238 int istat=0,iev=1; 00239 EventLoop: if (iev <= nevents && !istat) { 00240 cout << "StHbtExample -- Working on eventNumber " << iev << " of " << nevents << endl; 00241 chain->Clear(); 00242 istat = chain->Make(iev); 00243 if (istat) {cout << "Last event processed. Status = " << istat << endl;} 00244 iev++; goto EventLoop; 00245 } 00246 00247 // good old Cint can't even handle a for-loop 00248 // for (Int_t iev=1;iev<=nevents; iev++) { 00249 // chain->Clear(); 00250 // int iret = chain->Make(iev); // This should call the Make() method in ALL makers 00251 // if (iret) { 00252 // cout << "StHbtExample.C -- chain returned nonzero value " << iret 00253 // << " on event " << iev << endl; 00254 // break; 00255 // } 00256 // } // Event Loop 00257 00258 00259 00260 cout << "StHbtExample -- Done with event loop" << endl; 00261 00262 chain->Finish(); // This should call the Finish() method in ALL makers 00263 TheEnd: 00264 }
1.5.9