2782c1b5ab727c85606e1b413e18fed23be25fdd
[u/mrichter/AliRoot.git] / ITS / AliITSHits2SDR.C
1 #if !defined(__CINT__) || defined(__MAKECINT__)
2 #include <iostream.h>
3 #include <TClonesArray.h>
4 #include <TFile.h>
5 #include <TParticle.h>
6 #include <TStopwatch.h>
7 #include <TTree.h>
8 #include <TFile.h>
9
10 #include "AliRun.h"
11 #include "AliHeader.h"
12 #include "AliITS.h"
13 #include "AliITSDetType.h"
14 #include "AliITSsimulationFastPoints.h"
15 #include "AliITSresponse.h"
16 #include "AliITSresponseSDD.h"
17 #include "AliITSreconstruction.h"
18 #include "AliRunDigitizer.h"
19 #include "AliITSsegmentationSDD.h"
20 #include "AliITSDigitizer.h"
21
22 #endif
23
24 TFile* AccessInpFile(TString inFile);
25 void writeAR(TFile * fin, TFile *fou);
26 void AliITSSD2D(TString inFile = "galice.root", TString outFile = "");
27 void AliITSH2FR2files (TString inFile, TString outfile);
28 void AliITSH2SD2Files(TFile *file2);
29 void copyTK(TString inFile, TString outFile);
30
31 void AliITSHits2SDR(TString inFile = "galice.root", TString File1 = "galiceS.root",TString File2 = "galiceD.root", TString File3 = "galiceR.root"){
32   ////////////////////////////////////////////////////////////////////
33   // This macro takes hits from inFile file and 
34   // produce fast RecPoints, Sdigits, Digits and slow Rec Points 
35   // for the ITS using the standard detector simulations. 
36   // The first argument is the input file which contains the tree with hits
37   // The second argument is the output file name onto which TreeS will  
38   // be written. 
39   // The third and fourth arguments are the file names onto which TreeD  
40   // and TreeR+TreeC will be written (respectively).
41   // The AliRun object and the KINE tree will be saved onto 
42   // these files as well. 
43   // This macro will process all of the events on input the root file.
44   // This macro can be compiled and it should be launched from an aliroot 
45   // session.
46   // WARNING: Fast Rec Points are written on a branch named ITSRecPointsF
47   // this macro should be used with aliroot (it does not load the libraries) 
48   ////////////////////////////////////////////////////////////////////
49  
50  
51  
52   // Produce Fast Rec Points  (on File3)
53   cout<<"Producing Fast rec points...\n";
54   AliITSH2FR2files(inFile,File3);
55
56   // TreeK is copied to File3 (for comparison purposes)
57  
58   copyTK(inFile,File3);
59  
60   // Produce SDigits
61   TFile *file = AccessInpFile(inFile);
62   if(!file){
63     cerr<<"Problems in accessing the input file\n";
64     return;
65   }
66   TFile *file2 = gAlice->InitTreeFile("S",File1);
67   file2->cd();
68   cout<<"Producing Summable digits ....\n";
69   AliITSH2SD2Files(file2);
70   // write the AliRun object to the output file
71   writeAR(file,file2);
72   delete gAlice;
73   gAlice = 0;
74   file->Close();
75   delete file;
76   file2 = 0;
77   
78   // Produce Digits 
79   cout<<"Producing Digits ...\n";
80   AliITSSD2D(File1,File2);
81
82   // Produce Slow Rec. Points 
83   TFile *fin = new TFile(File2);
84   gAlice = (AliRun*) fin->Get("gAlice");
85   const char *nulptr=0;
86   AliITSreconstruction *itsr = new AliITSreconstruction(nulptr);
87   itsr->SetOutputFile(File3);
88   itsr->Init();
89   itsr->Exec();
90   TFile *fou = gAlice->GetTreeRFile();
91   writeAR(fin,fou);
92   delete itsr;
93   
94 }
95
96 void AliITSH2SD2Files(TFile *file2){
97
98   // Function used to produce s-digits
99  
100   AliITS *ITS  = (AliITS*) gAlice->GetModule("ITS");
101   if (!ITS) {
102         cerr<<"AliITS object not found on file" << endl;
103         return;
104   }
105   if(!(ITS->GetITSgeom())){
106         cerr << " AliITSgeom not found. Can't digitize without it." << endl;
107         return ;
108   }
109
110   // for old files
111   AliITSresponseSDD *resp1 = (AliITSresponseSDD*)ITS->DetType(1)->GetResponseModel();
112   TDatime *ct0= new TDatime(2002,04,26,00,00,00),ct=file2->GetCreationDate();
113   if(ct0->GetDate()<ct.GetDate()){
114     resp1 = new AliITSresponseSDD();
115     ITS->SetResponseModel(1,resp1);
116     cout << "Changed response class for SDD: \n";
117     resp1->Print();
118   } // end if
119   // end mods for old files (<26/4/2002)
120   for (Int_t nevent=0; nevent<gAlice->TreeE()->GetEntries(); nevent++) {
121     gAlice->GetEvent(nevent);
122     gAlice->MakeTree("S",file2);
123     ITS->MakeBranch("S");
124     ITS->SetTreeAddress();   
125     ITS->Hits2SDigits();
126   }
127 }
128
129
130 TFile * AccessInpFile(TString inFile){
131
132   // Function used to open the input file and fetch the AliRun object
133
134   TFile *retfil = 0;
135   TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(inFile);
136   if (file) {file->Close(); delete file;}
137   file = new TFile(inFile);
138   if (!file->IsOpen()) {
139         cerr<<"Can't open "<<inFile<<" !" << endl;
140         return retfil;
141   } 
142
143   // Get AliRun object from file or return if not on file
144   if (gAlice) delete gAlice; gAlice = 0;
145   gAlice = (AliRun*)file->Get("gAlice");
146   if (!gAlice) {
147         cerr << "AliRun object not found on file"<< endl;
148         return retfil;
149   } 
150   return file;
151 }
152
153 void writeAR(TFile * fin, TFile *fou) {
154   TDirectory *current = gDirectory;
155   TTree *Te;
156   TTree *TeNew;
157   AliHeader *alhe = new AliHeader();
158   Te = (TTree*)fin->Get("TE");
159   Te->SetBranchAddress("Header",&alhe);
160   Te->SetBranchStatus("*",1);
161   fou->cd();
162   TeNew = Te->CloneTree();
163   TeNew->Write(0,TObject::kOverwrite);
164   gAlice->Write(0,TObject::kOverwrite);
165   current->cd();
166   delete alhe;
167   cout<<"AliRun object written to file\n";
168 }
169
170
171 void AliITSSD2D(TString inFile, TString outFile){
172   AliRunDigitizer * manager = new AliRunDigitizer(1,1);
173   char ftmp[50];
174   sprintf(ftmp,"%s",inFile.Data());
175   manager->SetInputStream(0,ftmp);
176   if(outFile != "")manager->SetOutputFile(outFile);
177   AliITSDigitizer *dITS  = new AliITSDigitizer(manager);
178   manager->Exec("");
179   TFile *file = (TFile*)gROOT->GetListOfFiles()->FindObject(inFile);
180   TFile *file2= new TFile(outFile,"update");
181   writeAR(file,file2);
182   delete dITS;
183   delete manager;
184   file = (TFile*)gROOT->GetListOfFiles()->FindObject(inFile);
185   if(file){
186     file->Close();
187     delete file;
188   }
189   file2->Close();
190   delete file2;
191 }
192
193 void AliITSH2FR2files (TString inFile, TString outfile)
194 {
195
196   TFile * file = AccessInpFile(inFile);
197   if(!file){
198     cerr<<"*** AliITSH2FR2files: Problems in accessing the input file\n";
199     return;
200   }
201  
202   // open output file and create TreeR on it
203   TFile * file2 = gAlice->InitTreeFile("R",outfile);
204   file2->cd();
205   // get ITS  
206   AliITS *ITS  = (AliITS*) gAlice->GetModule("ITS");
207   if (!ITS) {
208         cerr<<"ITSH2FR2files.C : AliITS object not found on file" << endl;
209         return;
210   }  // end if !ITS
211
212   // Set the simulation model
213  
214   for (Int_t i=0;i<3;i++) {
215     ITS->SetSimulationModel(i,new AliITSsimulationFastPoints());
216   }
217  
218
219   //
220   // Event Loop
221   //
222   Int_t nsignal=25;
223   Int_t size=-1;
224   TStopwatch timer;
225   for (Int_t nevent=0; nevent<gAlice->TreeE()->GetEntries(); nevent++) {
226     gAlice->GetEvent(nevent);
227     if(!gAlice->TreeR())gAlice->MakeTree("R",file2);
228     ITS->MakeBranch("RF"); 
229     ITS->SetTreeAddress();  
230     Int_t bgr_ev=nevent/nsignal;
231     ITS->HitsToFastRecPoints(nevent,bgr_ev,size," ","All"," ");
232   }
233
234   timer.Stop(); timer.Print();
235
236   delete gAlice;
237   gAlice = 0;
238   file->Close();
239 }
240
241 void copyTK(TString inFile, TString outFile) {
242
243   TDirectory *current = gDirectory;
244   TParticle *part = new TParticle();
245   TTree *Tk;
246   TTree *TkNew;
247   TFile *fin = AccessInpFile(inFile);
248   TFile *fou = new TFile(outFile,"update");
249   char tname[30];
250   for(Int_t event= 0; event < gAlice->TreeE()->GetEntries(); event++){
251     current->cd();
252     sprintf(tname,"TreeK%d",event);
253     Tk = (TTree*)fin->Get(tname);
254     if(!Tk){
255       cerr<<"Trying to access a non-existent TTree : "<<tname<<endl;
256     }
257     else {
258       Tk->SetBranchAddress("Particles",&part);
259       Tk->SetBranchStatus("*",1);
260       fou->cd();
261       TkNew = Tk->CloneTree();
262       TkNew->Write();
263     }
264   }
265   delete gAlice;
266   gAlice = 0;
267   fin->Close();
268   fou->Close();
269   delete fin;
270   delete fou;
271   current->cd();
272
273 }