1 |
// C/C++ headers |
2 |
// |
3 |
#include <fstream> |
4 |
#include <string.h> |
5 |
#include <iostream> |
6 |
#include <stdio.h> |
7 |
#include <ctype.h> |
8 |
// |
9 |
// ROOT headers |
10 |
// |
11 |
|
12 |
#include <TTree.h> |
13 |
|
14 |
#include <TSystem.h> |
15 |
#include <TStopwatch.h> |
16 |
#include <TSystemDirectory.h> |
17 |
#include <TGraphAsymmErrors.h> |
18 |
#include <TString.h> |
19 |
#include <TFile.h> |
20 |
#include <TClass.h> |
21 |
#include <TGraph.h> |
22 |
#include <TLine.h> |
23 |
#include <TH1F.h> |
24 |
#include <TMath.h> |
25 |
#include <TRandom3.h> |
26 |
#include <TF1.h> |
27 |
// |
28 |
// RunInfo header |
29 |
// |
30 |
#include <RunInfo.h> |
31 |
#include <PamLevel2.h> |
32 |
#include <TClonesArray.h> |
33 |
|
34 |
//PamVMC headers |
35 |
#include "PamVMCDigMgr.h" |
36 |
#include "PamVMCTrkHit.h" |
37 |
#include "PamVMCTrkF77GpsSpe.h" |
38 |
#include "PamVMCPrimaryGenerator.h" |
39 |
using namespace std; |
40 |
//This program is standalone Digitizer for PamVMC |
41 |
|
42 |
|
43 |
void Usage (); |
44 |
|
45 |
TClonesArray* FindColl(TTree*t, TString brname, TString type = "PamVMCDetectorHit"); |
46 |
void RegisterPDG(); |
47 |
//Here user should define all nuclei and particles that are not in PDG tabe of ROOT |
48 |
//otherwise application will crash on TOF digitization whrer charge information needed; |
49 |
|
50 |
void UnCompressCalo(TClonesArray* cal_comp, TClonesArray & cal_u); |
51 |
|
52 |
int main(int argc, char * argv[] ){ |
53 |
|
54 |
if(argc!=3){ |
55 |
Usage(); |
56 |
return 0; |
57 |
} |
58 |
|
59 |
TString DataPATH = gSystem->Getenv("PWD"); |
60 |
TString inputfilename(argv[1]); |
61 |
if(!inputfilename.Contains(".root")){ |
62 |
cerr<<"Error! Not root file on input!"<<endl; |
63 |
return 0; |
64 |
} |
65 |
|
66 |
TString outputfilename(argv[2]); |
67 |
|
68 |
TStopwatch timer; |
69 |
timer.Start(); |
70 |
TFile* rfin = new TFile(inputfilename, "UPDATE"); |
71 |
if(!rfin){ |
72 |
cerr<<"Error! Can't read root file!"<<endl; |
73 |
return 0; |
74 |
} |
75 |
|
76 |
TTree::SetMaxTreeSize(1000*Long64_t(2000000000)); |
77 |
|
78 |
TTree * hits_tree = (TTree*)rfin->Get("hits"); |
79 |
if(!hits_tree){ |
80 |
cerr<<"Tree hits are not found in file "<<rfin->GetName()<<endl; |
81 |
TString treename = rfin->GetName(); |
82 |
treename.Replace(treename.Index(".root",5),5,""); |
83 |
hits_tree = (TTree*)rfin->Get(treename); |
84 |
if(!hits_tree) return 0; |
85 |
} |
86 |
|
87 |
TObjArray * frndArr = new TObjArray(); |
88 |
hits_tree->SetBranchAddress("RNDM",&frndArr); |
89 |
|
90 |
PamVMCRndMgr* rnd_mgr = PamVMCRndMgr::Instance_ext(); |
91 |
cout<<"dmdr"<<endl; |
92 |
PamVMCDigMgr* d_mgr = PamVMCDigMgr::Instance_ext(); |
93 |
cout<<"dmdr"<<endl; |
94 |
|
95 |
//AC CARD: C1D1, C2D1, CAT: TOP1, CAS: SID1 |
96 |
PamVMCAcDig* ac = new PamVMCAcDig(); |
97 |
ac->RegisterCollections("C1D1", FindColl(hits_tree,"C1D1")); |
98 |
ac->RegisterCollections("C2D1", FindColl(hits_tree,"C2D1")); |
99 |
ac->RegisterCollections("TOP1", FindColl(hits_tree,"TOP1")); |
100 |
ac->RegisterCollections("SID1", FindColl(hits_tree,"SID1")); |
101 |
|
102 |
//TOF: S11Y, S12X, S21X, S22Y, S31Y, S32X |
103 |
PamVMCTofDig* tof = new PamVMCTofDig(); |
104 |
tof->RegisterCollections("S11Y", FindColl(hits_tree,"S11Y")); |
105 |
tof->RegisterCollections("S12X", FindColl(hits_tree,"S12X")); |
106 |
tof->RegisterCollections("S21X", FindColl(hits_tree,"S21X")); |
107 |
tof->RegisterCollections("S22Y", FindColl(hits_tree,"S22Y")); |
108 |
tof->RegisterCollections("S31Y", FindColl(hits_tree,"S31Y")); |
109 |
tof->RegisterCollections("S32X", FindColl(hits_tree,"S32X")); |
110 |
|
111 |
|
112 |
//Calorimeter CAST |
113 |
#define NSTRIPS 4224 |
114 |
PamVMCCaloDig* cal = new PamVMCCaloDig(); |
115 |
TClonesArray *cal_u = new TClonesArray("PamVMCDetectorHit",NSTRIPS); //uncompressed array |
116 |
TClonesArray* cal_comp = FindColl(hits_tree,"CAST"); // compressed array from file |
117 |
cal->RegisterCollections("CAST", cal_u); |
118 |
|
119 |
//Tracker GPSSPE |
120 |
PamVMCTrkDig* trk = new PamVMCTrkDig(); |
121 |
TBranch *gps = hits_tree->GetBranch("GPSSPE"); |
122 |
trk->RegisterCollections("TSPA", FindColl(hits_tree,"TSPA")); |
123 |
pGPSSPEHits * hits = new pGPSSPEHits(); |
124 |
if(gps){ |
125 |
hits_tree->SetBranchAddress("GPSSPE",&hits); |
126 |
hits_tree->SetBranchStatus("GPSSPE",1); |
127 |
cout<<"TRKCOLL="<<hits<<endl; |
128 |
trk->RegisterTrkData(hits); |
129 |
} |
130 |
|
131 |
|
132 |
//S4 |
133 |
PamVMCS4Dig* s4 = new PamVMCS4Dig(); |
134 |
s4->RegisterCollections("S4", FindColl(hits_tree,"S4")); |
135 |
|
136 |
//ND: NDTI |
137 |
PamVMCNDDig* nd = new PamVMCNDDig(); |
138 |
nd->RegisterCollections("NDTI", FindColl(hits_tree,"NDTI")); |
139 |
|
140 |
|
141 |
//Primaries |
142 |
TClonesArray *prims = FindColl(hits_tree,"PRIM","PamVMCPrimary"); |
143 |
|
144 |
//Initialization |
145 |
d_mgr->SetDIG("Tof",tof); |
146 |
d_mgr->SetDIG("AC",ac); |
147 |
d_mgr->SetDIG("CAST",cal); |
148 |
d_mgr->SetDIG("TSPA",trk); |
149 |
d_mgr->SetDIG("S4",s4); |
150 |
d_mgr->SetDIG("NDTI",nd); |
151 |
|
152 |
d_mgr->Initialize( new TRandom3(0)); //Load Calibrations and set dummy random |
153 |
|
154 |
PamVMCRawMgr* r_mgr = PamVMCRawMgr::Instance(); |
155 |
r_mgr->CreateOutFile(outputfilename.Data()); |
156 |
|
157 |
r_mgr->WriteToFile(); |
158 |
|
159 |
Int_t nev = hits_tree->GetEntries(); |
160 |
cout<<"TOTAL NUMBER OF EVENTS:"<<nev<<endl; |
161 |
RegisterPDG(); |
162 |
//Event loop |
163 |
d_mgr->PrintCollections(); |
164 |
for(Int_t i = 0; i<nev; i++){ |
165 |
hits_tree->GetEntry(i); |
166 |
cout<<"--> Setting Random"<<endl; |
167 |
for(Int_t j = 0; j<frndArr->GetEntries(); j++){ |
168 |
TRandom3 * rnd = (TRandom3*)frndArr->At(j); |
169 |
TString nm = rnd->GetName(); |
170 |
if ( nm=="AC" ) ac->SetRandom(rnd); |
171 |
if ( nm=="Calo" ) cal->SetRandom(rnd); |
172 |
if ( nm=="Tracker" ) trk->SetRandom(rnd); |
173 |
if ( nm=="S4" ) s4->SetRandom(rnd); |
174 |
if ( nm=="Tof" ) tof->SetRandom(rnd); |
175 |
if ( nm=="ND" ) nd->SetRandom(rnd); |
176 |
} |
177 |
|
178 |
cout<<"--> Digitizing evnt"<<endl; |
179 |
UnCompressCalo(cal_comp, *cal_u); |
180 |
d_mgr->Digitize(((PamVMCPrimary*)prims->At(0))->fID,((PamVMCPrimary*)prims->At(0))->fPDG); |
181 |
cout<<"--> Writing evnt"<<endl; |
182 |
r_mgr->WriteEvent(); |
183 |
rnd_mgr->ClearColl(); |
184 |
} |
185 |
|
186 |
r_mgr->FinishRun(); |
187 |
d_mgr->PrintCollections(); |
188 |
timer.Stop(); |
189 |
gDirectory->Delete("tmp.root"); |
190 |
rfin->Close(); |
191 |
return 0; |
192 |
} |
193 |
|
194 |
TClonesArray* FindColl(TTree*t, TString brname, TString type){ |
195 |
if(t->GetBranch(brname)){ |
196 |
TClonesArray* coll = 0;//new TClonesArray(type,4224); |
197 |
t->GetBranch(brname)->SetAddress(&coll); |
198 |
cout<<"coll="<<coll<<endl; |
199 |
return coll; |
200 |
} |
201 |
return 0; |
202 |
} |
203 |
|
204 |
void RegisterPDG(){ |
205 |
TDatabasePDG::Instance()->AddParticle("H-2","H-2",1.875613,kTRUE,0., 3. /* |e|/3 charge*/,"Ion",1000010020); |
206 |
TDatabasePDG::Instance()->AddParticle("H-3","H-3",2.808921,kTRUE,0., 3. /* |e|/3 charge*/,"Ion",1000010030); |
207 |
TDatabasePDG::Instance()->AddParticle("He-4","He-4",3.727379,kTRUE,0., 6. /* |e|/3 charge*/,"Ion",1000020040); |
208 |
TDatabasePDG::Instance()->AddParticle("He-3","He-3",2.808391,kTRUE,0., 6. /* |e|/3 charge*/,"Ion",1000020030); |
209 |
|
210 |
Double_t gev = 0.93146;// conversion a.m.u. to gev; |
211 |
|
212 |
//nuclei and other stuff |
213 |
TDatabasePDG::Instance()->AddParticle("Li-6","Li-6",6.015122*gev,kTRUE,0., 9. /* |e|/3 charge*/,"Ion",1000030060); |
214 |
TDatabasePDG::Instance()->AddParticle("Li-7","Li-7",7.016004*gev,kTRUE,0., 9. /* |e|/3 charge*/,"Ion",1000030070); |
215 |
TDatabasePDG::Instance()->AddParticle("Be-7","Be-7",7.01692983*gev,kTRUE,0., 12. /* |e|/3 charge*/,"Ion",1000040070); |
216 |
TDatabasePDG::Instance()->AddParticle("Be-9","Be-9",9.0121822*gev,kTRUE,0., 12. /* |e|/3 charge*/,"Ion",1000040090); |
217 |
TDatabasePDG::Instance()->AddParticle("Be-10","Be-10",10.013533818*gev,kTRUE,0., 12. /* |e|/3 charge*/,"Ion",1000040100); |
218 |
TDatabasePDG::Instance()->AddParticle("B-10","B-10",10.0129370*gev,kTRUE,0., 15. /* |e|/3 charge*/,"Ion",1000050100); |
219 |
TDatabasePDG::Instance()->AddParticle("B-11","B-11",11.0093054*gev,kTRUE,0., 15. /* |e|/3 charge*/,"Ion",1000050110); |
220 |
TDatabasePDG::Instance()->AddParticle("C-12","C-12",12.*gev,kTRUE,0., 18. /* |e|/3 charge*/,"Ion",1000060120); |
221 |
TDatabasePDG::Instance()->AddParticle("N-14","N-14",14.0030740048*gev,kTRUE,0., 21. /* |e|/3 charge*/,"Ion",1000070140); |
222 |
TDatabasePDG::Instance()->AddParticle("O-16","O-16",15.99491461956*gev,kTRUE,0., 23. /* |e|/3 charge*/,"Ion",1000080160); |
223 |
} |
224 |
|
225 |
void UnCompressCalo(TClonesArray* cal_comp, TClonesArray & cal_u){ |
226 |
if(!cal_comp) return; |
227 |
PamVMCDetectorHit* hit = NULL; |
228 |
PamVMCDetectorHit* newhit = NULL; |
229 |
cal_u.Clear("C"); |
230 |
for(Int_t j = 0; j<cal_comp->GetEntries(); j++){ |
231 |
hit=((PamVMCDetectorHit*)cal_comp->At(j)); |
232 |
newhit = (PamVMCDetectorHit*)cal_u.New(hit->GetPOS()); |
233 |
*newhit = *hit; |
234 |
} |
235 |
} |
236 |
|
237 |
void Usage(){ |
238 |
cout<<"NAME"<<endl; |
239 |
cout<<" vmc_dig - a standalone digitizer for PamVMC application"<<endl; |
240 |
cout<<"SYNOPSIS"<<endl; |
241 |
cout<<" vmc_dig INPROOTFILE OUTRAWFILE"<<endl; |
242 |
cout<<"DISCRIPTION"<<endl; |
243 |
cout<<" This tool convert hits collections of PamVMC in raw format "<<endl; |
244 |
cout<<" There are only 2 parameters: input root/output raw filenames."<<endl; |
245 |
cout<<" Please note that your input file should contain RNDM branch"<<endl; |
246 |
cout<<" ( produced with <all_detectors/> option in <save_mode> section )"<<endl; |
247 |
} |