Adding some further mother volumes to speed-up the overlap checking and particle...
[u/mrichter/AliRoot.git] / ITS / AliITSL2GConvertPointsV1.C
1 // 
2 // Macro to convert ITS local-coordinate points
3 // into globa lones
4 // 
5
6 Int_t AliITSL2GConvertPoints 
7 (const char* in_name = "ITS.RecPoints.root", 
8  const char* out_name = "ITS.Neural.PointsV1.root", Int_t nev = 0)
9 {
10         TStopwatch timer;
11
12         // Open output file
13         TFile *in = new TFile(in_name);
14         TFile *out = new TFile(out_name, "recreate");
15         
16         // Load event files
17         if (gAlice) {
18                 delete AliRunLoader::Instance();
19                 delete gAlice;
20                 gAlice=0;
21         } 
22         AliRunLoader* rl = AliRunLoader::Open("galice.root");
23     if (rl == 0x0) {
24                 cerr << "AliITSL2GConvertPoints.C : Can not open session." << endl;
25                 return 3;
26         }
27     Int_t retval = rl->LoadgAlice();
28         if (retval) {
29                 cerr << "AliITSL2GConvertPoints.C : LoadgAlice returned error" << endl;
30                 return 3;
31         }
32         gAlice=rl->GetAliRun();
33         AliITSLoader* gime = (AliITSLoader*)rl->GetLoader("ITSLoader");
34         if (gime == 0x0) {
35                 cerr << "AliITSL2GConvertPoints.C : can not get ITS loader" << endl;
36                 return 3;
37         }
38         AliITS *ITS = (AliITS*)gAlice->GetDetector("ITS");
39         if (!ITS) {
40                 cerr << "AliITSL2GConvertPoints.C : AliITS object not found on file" << endl;
41                 return 3;
42         }  // end if !ITS
43         AliITSgeom *geom = (AliITSgeom*)ITS->GetITSgeom();
44         if(!geom) {
45                 cerr << "AliITSL2GConvertPoints.C : AliITSgeom not found." << endl;
46                 return 4;
47         } // end if
48         
49         // Tree of recpoints
50         Int_t nModules = 0;
51         TTree *TR = (TTree*)in->Get(Form("Event%d/TreeR", nev));
52         if (!TR) {
53                 cout << "TreeR not found" << endl;
54                 return;
55         }
56         nModules = (Int_t)TR->GetEntries();
57         if (!nModules) {
58                 cout << "Empty TreeR!!!" << endl;
59                 return;
60         }
61
62         timer.Start();
63         
64         // Converts and stores the ITS points into global coordinate format
65         Int_t pos = 0;
66         AliITSRecPoint *local = 0;
67         AliITSNeuralPoint *global = 0;
68         TTree *TP = new TTree("TreeP", "Event points in global coords");
69         TP->Branch("pos", &pos, "pos/I");
70         TP->Branch("Points", "AliITSNeuralPoint", &global);
71         
72         TObjArray *localArray = 0;
73         TR->SetBranchAddress("ITSRecPoints", &localArray);
74         Int_t module, layer, i, j, count, index;
75         Double_t locPos[3], globPos[3], locErr[3][3], globErr[3][3];
76         
77         for(module = 0; module < nModules; module++) {
78                 TR->GetEvent(module);
79                 if (module > geom->GetLastSSD()) {
80                         cout << "Strange behavior: an entry greater than the last SSD index!" << endl;
81                         continue;
82                 }
83                 AliITSgeomMatrix *gm = geom->GetGeomMatrix(module);
84                 geom->GetModuleId(module, layer, i, j);
85                 count = (Int_t)localArray->GetEntriesFast();
86                 layer--;
87                 if (layer < 0 || layer > 5) {
88                         cout << "STRANGE layer value: " << layer << endl;
89                         continue;
90                 }
91                 for (index = 0; index < count; index++) {
92                         local = (AliITSRecPoint*)localArray->At(index);
93                         global = new AliITSNeuralPoint(local, gm);
94                         global->SetLayer(layer);
95                         global->SetModule(module);
96                         global->SetIndex(index);
97                         global->SetUser(-1);
98                         global->ConfMap(0.0, 0.0);
99                         TP->Fill();
100                         pos++;
101                 }
102         }
103
104         timer.Stop();
105         timer.Print();
106         cout << TP->GetEntries() << " points collected" << endl;
107         
108         out->cd();
109         out->mkdir(Form("Event%d", nev));
110         out->cd(Form("Event%d", nev));
111         TP->Write(Form("TreeP", nev));
112         out->Close();
113 }
114