StRoot  1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Groups Pages
bfcMixer_TpxAuAu39.C
1 //
3 // Macro for running chain with different inputs
4 //
5 // Owner: Yuri Fisyak
6 //
7 // $Id: bfcMixer_TpxAuAu39.C,v 1.3 2012/02/17 20:47:36 fisyak Exp $
8 //
10 
11 class StChain;
12 StChain *Chain=0;
13 class StBFChain;
14 StBFChain *chain1, *chain2, *chain3;
15 //_____________________________________________________________________
16 void bfcMixer_TpxAuAu39(const Int_t Nevents=100,
17  const Char_t *daqfile="/star/rcf/test/daq/2010/embed/auau39_prod/st_physics_adc_11100070_raw_1490001.daq",
18  const Char_t *tagfile="/star/rcf/test/daq/2010/embed/auau39_prod/st_physics_adc_11100070_raw_1490001.tags.root",
19  const Double_t pt_low=0.1,
20  const Double_t pt_high=5.0,
21  const Double_t eta_low=-1.5,
22  const Double_t eta_high=1.5,
23  const Double_t vzlow = -150.0,
24  const Double_t vzhigh = 150.0,
25  const Int_t pid=9,
26  const Double_t mult=100,
27  const std::vector<Int_t> triggers = 0,
28  const Char_t *prodName = "P10ihAuAu39",
29  const Char_t* type = "FlatPt"){
30  // production chains for P08ic - p+p, Au+Au 9 GeV and d+Au
31  TString prodP08iepp("DbV20081117 B2008a ITTF IAna ppOpt l3onl emcDY2 fpd ftpc trgd ZDCvtx NosvtIT NossdIT Corr4 OSpaceZ2 OGridLeak3D VFMCE -hitfilt");
32 // TString prodP08icpp("DbV20080712,pp2008,ITTF,OSpaceZ2,OGridLeak3D,beamLine,VFMCE,TpxClu -VFPPV -hitfilt");
33 // TString prodP08icAuAu9("DbV20080709 P2008 ITTF VFMCE -hitfilt");
34 // TString prodP08icAuAu200("DbV20070101 P2008 ITTF VFMCE -hitfilt");
35 // TString prodP08icdAu("DbV20080712 P2008 ITTF OSpaceZ2 OGridLeak3D beamLine, VFMCE TpxClu -VFMinuit -hitfilt");
36  TString prodP08iedAu("DbV20090213 P2008 ITTF OSpaceZ2 OGridLeak3D beamLine VFMCE TpxClu -VFMinuit -hitfilt");
37  TString prodP10iapp("DbV20091001 pp2009c TpcRS ITTF OSpaceZ2 OGridLeak3D beamLine, VFMCE TpcRS -VFMinuit -hitfilt");
38 
39  // BES Run10 chains
40  TString prodP10ihAuAu39("DbV20100909 P2010a,btof,BEmcChkStat,Corr4,OSpaceZ2,OGridLeak3D,VFMCE TpxClu -VFMinuit -hitfilt");
41  TString prodP10ihAuAu11("DbV20100821 P2010a,btof,BEmcChkStat,Corr4,OSpaceZ2,OGridLeak3D,VFMCE TpxClu -VFMinuit -hitfilt");
42  TString prodP10ihAuAu7("DbV20100821 P2010a,btof,BEmcChkStat,Corr4,OSpaceZ2,OGridLeak3D,VFMCE TpxClu -VFMinuit -hitfilt");
43 
44  TString geomP08ic("ry2008");
45  TString geomP10ih("ry2010");
46  TString chain1Opt("in,magF,tpcDb,NoDefault,TpxRaw,-ittf,NoOutput");
47  TString chain2Opt("NoInput,PrepEmbed,gen_T,geomT,sim_T,TpcRS,-ittf,-tpc_daq,nodefault");
48 // TString chain2Opt("NoInput,PrepEmbed,gen_T,geomT,sim_T,trs,-ittf,-tpc_daq,nodefault");
49  chain2Opt += " ";
50 
51  TString chain3Opt("");
52  if (prodName == "P08icpp") { chain3Opt = prodP08icpp; chain2Opt += geomP08ic; }
53  else if (prodName == "P08iepp") { chain3Opt = prodP08iepp; chain2Opt += geomP08ic; }
54  else if (prodName == "P08icAuAu9") { chain3Opt = prodP08icAuAu9; chain2Opt += geomP08ic; }
55  else if (prodName == "P08icdAu") { chain3Opt = prodP08icdAu; chain2Opt += geomP08ic; }
56  else if (prodName == "P08iedAu") { chain3Opt = prodP08iedAu; chain2Opt += geomP08ic; }
57  else if (prodName == "P08icAuAu200") { chain3Opt = prodP08icAuAu200; chain2Opt += geomP08ic; }
58  else if (prodName == "P10iapp") { chain3Opt = prodP10iapp; chain2Opt += geomP10ih; }
59  else if (prodName == "P10ihAuAu39") { chain3Opt = prodP10ihAuAu39; chain2Opt += geomP10ih; }
60  else if (prodName == "P10ihAuAu11") { chain3Opt = prodP10ihAuAu11; chain2Opt += geomP10ih; }
61  else if (prodName == "P10ihAuAu7") { chain3Opt = prodP10ihAuAu7; chain2Opt += geomP10ih; }
62  else {
63  cout << "Choice prodName " << prodName << " does not correspond to known chain. Processing impossible. " << endl;
64  return;
65  }
66  chain3Opt += ",Embedding,TpcMixer,GeantOut,MiniMcMk,McAna,-in,NoInput,useInTracker";
67  chain3Opt += ",";
68 
69  if (prodName == "P08icpp") { chain3Opt += geomP08ic; }
70  else if (prodName == "P08iepp") { chain3Opt += geomP08ic; }
71  else if (prodName == "P08icAuAu9") { chain3Opt += geomP08ic; }
72  else if (prodName == "P08icdAu") { chain3Opt += geomP08ic; }
73  else if (prodName == "P08iedAu") { chain3Opt += geomP08ic; }
74  else if (prodName == "P08icAuAu200") { chain3Opt += geomP08ic; }
75  else if (prodName == "P10iapp") { chain3Opt += geomP10ih; }
76  else if (prodName == "P10ihAuAu39") { chain3Opt += geomP10ih; }
77  else if (prodName == "P10ihAuAu11") { chain3Opt += geomP10ih; }
78  else if (prodName == "P10ihAuAu7") { chain3Opt += geomP10ih; }
79  else {
80  cout << "Choice prodName " << prodName << " does not correspond to known chain. Processing impossible. " << endl;
81  return;
82  }
83 
84  // Dynamically link some shared libs
85  gROOT->LoadMacro("bfc.C");
86  if (gClassTable->GetID("StBFChain") < 0) Load();
87  //______________Create the main chain object______________________________________
88  Chain = new StChain("Embedding");
89  //________________________________________________________________________________
90  bfc(-1,chain1Opt,daqfile);
91  chain1 = chain;
92  chain1->SetName("One");
93  Chain->cd();
94  //________________________________________________________________________________
95  bfc(-1,chain2Opt);
96  chain2 = chain;
97  chain2->SetName("Two");
98  Chain->cd();
99 #if 0
100  if (chain2->GetOption("TRS")){
101  StTrsMaker *trsMk = (StTrsMaker *) chain2->GetMaker("Trs");
102  if (! trsMk) {
103  cout << "Cannot find Trs in chain2" << endl;
104  return;
105  }
106  trsMk->setNormalFactor(1.32);
107  trsMk->SetMode(0);
108  }
109 #endif
110  //________________________________________________________________________________
111  // gSystem->Load("StFtpcMixerMaker");
112  // StFtpcMixerMaker *ftpcmixer = new StFtpcMixerMaker("FtpcMixer","daq","trs");
113  //________________________________________________________________________________
114  TString OutputFileName(gSystem->BaseName(daqfile));
115  OutputFileName.ReplaceAll("*","");
116  OutputFileName.ReplaceAll(".daq","");
117  // OutputFileName.Append("_emb.root");
118  OutputFileName.Append(".root");
119  bfc(-1,chain3Opt,0,OutputFileName);
120  chain3 = chain;
121  chain3->SetName("Three");
122  Chain->cd();
123  //________________________________________________________________________________
124  StTpcMixerMaker *mixer = (StTpcMixerMaker *) chain3->Maker("TpcMixer");
125  if( prodName == "P08icAuAu200")
126  {
127  mixer->SetInput("Input1","MixerEvent");
128  }
129  else
130  {
131  mixer->SetInput("Input1","TpxRaw/.data/Event");
132  }
133  // mixer->SetInput("Input2","Trs/.const/Event");
134 
135  if (chain2Opt.Contains("TpcRS",TString::kIgnoreCase)) {
136  mixer->SetInput("Input2","TpcRS/Event");
137  }else {
138  mixer->SetInput("Input2","Trs/.const/Event");
139  }
140 
141  Chain->cd();
142 
143  //............. begin of EMC embedding makers................
144 
145  //.............. Add BEmc stuff here ....................
146  gSystem->Load("StEmcSimulatorMaker");
147  gSystem->Load("StEmcMixerMaker");
148  gSystem->Load("StEEmcSimulatorMaker");
149 
150  StMcEventMaker* mcEventMaker = new StMcEventMaker();
151  StEmcSimulatorMaker *bemcSim = new StEmcSimulatorMaker();
152  StEmcMixerMaker *bemcMixer = new StEmcMixerMaker();
153  chain3->AddAfter("emcRaw",bemcMixer);
154  chain3->AddAfter("emcRaw",bemcSim);
155  chain3->AddAfter("emcRaw",mcEventMaker);
156  bemcMixer->SetDebug(0); // set it to 1 for more printouts
157  // note, Barrel slow sim is always ON, said Adam
158 
159  //........... Add EEmc Stuff ( Simu, and Mixer) here ..............
160  StEEmcFastMaker *eemcFastSim = new StEEmcFastMaker();
161  StEEmcMixerMaker *eemcMixer = new StEEmcMixerMaker();
162 
163  /* position B+E EMC makers in the chain
164  (order is reverse because 'After' is used - looks funny but is right)
165  */
166  chain3->AddAfter("emcRaw",eemcMixer);
167  chain3->AddAfter("emcRaw",eemcFastSim);
168 
169  eemcFastSim->SetEmbeddingMode();
170  // eemcFastSim->SetDebug();
171  // eemcMixer->SetDebug();
172 
173  bool useEndcapSlowSim = true;
174  if(useEndcapSlowSim) { // turn Endcap slow simu On/Off
175  StEEmcSlowMaker *slowSim=new StEEmcSlowMaker();
176  chain3->AddAfter("EEmcFastSim",slowSim);
177  slowSim->setEmbeddingMode();
178  }
179 
180 
181 
182  //________________________________________________________________________________
183  {
184  TDatime t;
185  gMessMgr->QAInfo() << Form("Run is started at Date/Time %i/%i",t.GetDate(),t.GetTime()) << endm;
186  }
187  gMessMgr->QAInfo() << Form("Run on %s in %s",gSystem->HostName(),gSystem->WorkingDirectory()) << endm;
188  gMessMgr->QAInfo() << Form("with %s", Chain->GetCVS()) << endm;
189  // embedded particle set
190  StPrepEmbedMaker *embMk = (StPrepEmbedMaker *) Chain->Maker("PrepEmbed");
191  if (! embMk) return;
192  cout << "bfcMixer: Setting PID: "<<pid<<endl;
193  embMk->SetTagFile(tagfile);
194  // pTlow,ptHigh,etaLow,etaHigh,phiLow,phiHigh
195  embMk->SetOpt( pt_low, pt_high, eta_low, eta_high, 0., 6.283185, type);
196  // pid, mult
197  embMk->SetPartOpt( pid,mult);
198 
199  // Default is no event selections
200  embMk->SetSkipMode(kFALSE);
201 
202  // Make trigger and z-vertex cuts (only if SkipMode is true)
203  // Trigger cut
204  // Can put multiple trigger id's
205  if ( !triggers.empty() ){
206  for(std::vector<Int_t>::iterator iter = triggers.begin(); iter != triggers.end(); iter++){
207  embMk->SetTrgOpt((*iter)) ;
208  }
209  }
210 
211  // z-vertex cuts
212  embMk->SetZVertexCut(vzlow, vzhigh) ;
213 
214  TAttr::SetDebug(0);
215  Chain->SetAttr(".Privilege",0,"*" ); //All makers are NOT priviliged
216  Chain->SetAttr(".Privilege",1,"StBFChain::*" ); //StBFChain is priviliged
217  Chain->SetAttr(".Privilege",1,"StIOInterFace::*" ); //All IO makers are priviliged
218  Chain->SetAttr(".Privilege",1,"St_geant_Maker::*"); //It is also IO maker
219  Chain->SetAttr(".Privilege",1,"StPrepEmbedMaker::*"); //It is also IO maker
220  // Chain->SetDEBUG(0);
221  if (Nevents < 0) return;
222  Int_t iInit = Chain->Init();
223  if (iInit >= kStEOF) {Chain->FatalErr(iInit,"on init"); return;}
224  StMaker *treeMk = Chain->GetMaker("outputStream");
225  Chain->EventLoop(Nevents,treeMk);
226  gMessMgr->QAInfo() << "Run completed " << endm;
227  gSystem->Exec("date");
228 }
229 
230 //$LOG$
Prepare GEANT Maker with input from embedding settings and DAQ event.
Filling of all StMcEvent classes from g2t tables Transform all the data in the g2t tables into the co...
Definition: Stypes.h:43
virtual Int_t Load()
Routine handling library loading depending on chain options.
Definition: StBFChain.cxx:117
void setEmbeddingMode(Bool_t e=true)
Sets all switches required to perform embedding.
Slow simulator for EEMC.