Improving readibility of macro
[u/mrichter/AliRoot.git] / HLT / global / FlatESDConverter.C
1 /**
2  * >> Testing Macro to fill FlatESDEvent from ALiESDEvent <<
3  **
4  * Primary Authors : Sergey Gorbunov, Jochen Thaeder, Chiara Zampolli
5  *
6  * Usage:
7  *  aliroot -b -l -q LoadLibs.C FlatESDConverter.C++
8  *
9  **************************************************************************/
10
11 #if !defined(__CINT__) || defined(__MAKECINT__)
12 #include "AliESDEvent.h"
13 #include "AliESD.h"
14 #include "AliESDfriend.h"
15 #include "AliESDfriendTrack.h"
16 #include "../TPC/Rec/AliTPCseed.h"
17 #include <TFile.h>
18 #include <TTree.h>
19 #include <TSystem.h>
20 #include "./AliFlatESDEvent.h"
21 #include "./AliFlatESDTrack.h"
22 #include "./AliFlatTPCCluster.h"
23 #include "./AliFlatExternalTrackParam.h"
24 #include "Riostream.h"
25 #include "AliSysInfo.h"
26 #endif   
27
28 void FlatESDConverter(const char* filename="AliESDs.root", const char* filenameFriends="AliESDfriends.root",const char* filenameOut="out.dat", Bool_t useESDFriends = kTRUE, Bool_t useHLTtree = kFALSE,Int_t verbose = 0) {
29
30   // -- Convert AliESDEvent to AliFlatESDEvent
31
32   ofstream outFile(Form("%s",filenameOut), std::ifstream::binary | std::ifstream::out);
33   //ofstream outFile("outFlatESD.dat");
34
35   TFile *file    = new TFile(Form("%s", filename));
36   TTree *esdTree = useHLTtree? dynamic_cast<TTree*>(file->Get("HLTesdTree")) : dynamic_cast<TTree*>(file->Get("esdTree"));
37   
38   // -- Connect ESD
39   AliESDEvent *esd = new AliESDEvent;
40   esd->ReadFromTree(esdTree);
41   
42   // -- Connect ESD friend
43   AliESDfriend *esdFriend = NULL; 
44   if (useESDFriends && !esdTree->FindBranch("ESDfriend.")) {
45     esdTree->AddFriend("esdFriendTree", Form("%s", filenameFriends));
46     esdTree->SetBranchStatus("ESDfriend.", 1);
47     esdFriend = dynamic_cast<AliESDfriend*>((const_cast<AliESDEvent*>(esd))->FindListObject("AliESDfriend"));
48     if (esdFriend)
49       esdTree->SetBranchAddress("ESDfriend.", &esdFriend);
50   } // if (!esdTree->FindBranch("ESDfriend.")) {
51   
52   AliFlatESDEvent *flatEsd = NULL;
53
54   // -- Event Loop
55   for (Int_t idxEvent = 0; idxEvent < esdTree->GetEntries(); idxEvent++) {
56     Printf("Processing event nr %d", idxEvent);
57     esd->SaveAs("esdTemp.root");
58     TFile  fTmp = TFile("esdTemp.root");
59     Int_t sizeIn = fTmp.GetSize();
60
61     AliSysInfo::AddStamp("getEntry",0,0,idxEvent);
62
63     esdTree->GetEntry(idxEvent);
64     // -- Book memory for AliFlatESDEvent
65     // -- TEST >>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
66
67     Byte_t *mem = new Byte_t[AliFlatESDEvent::EstimateSize(esd, useESDFriends)];
68
69     flatEsd = reinterpret_cast<AliFlatESDEvent*>(mem);
70     new (flatEsd) AliFlatESDEvent;
71
72     AliSysInfo::AddStamp("DoEvent.Start",0,0,idxEvent);
73     // -- Fill AliFlatESDEvent
74
75     flatEsd->Fill(esd, useESDFriends);  
76         
77     AliSysInfo::AddStamp("DoEvent.Stop",sizeIn,flatEsd->GetSize(),idxEvent);
78
79     if (verbose) {
80       Printf("TEST: Event %d || Tracks %d | FRIEND Tracks %d || estimated size %llu || sizeof(AliFlatESDEvent) %llu", 
81              idxEvent, esd->GetNumberOfTracks(), esdFriend->GetNumberOfTracks(), 
82              AliFlatESDEvent::EstimateSize(esd, useESDFriends), flatEsd->GetSize());
83
84       AliFlatESDTrack *track = flatEsd->GetTracks();
85       for (Int_t idxTrack = 0; idxTrack < flatEsd->GetNumberOfTracks(); ++idxTrack) {      
86         AliESDtrack *esdTrack = esd->GetTrack(idxTrack);      
87         AliESDfriendTrack *friendTrack = esdFriend->GetTrack(idxTrack) ;
88         if (track && !esdTrack) {
89           Printf("ERROR THIS SHOULD NOT HAPPEN AT ALL !!! TRACK %d HAS NO ESD TRACK!!!", idxTrack);
90           return;
91         }
92         if (track) {
93           AliFlatExternalTrackParam* exp1 = track->GetTrackParamCp();
94           AliFlatExternalTrackParam* exp2 = track->GetTrackParamIp();
95           AliFlatExternalTrackParam* exp3 = track->GetTrackParamTPCInner();
96           AliFlatExternalTrackParam* exp4 = track->GetTrackParamOp();
97
98           Float_t alphaFLAT[4] = {0., 0., 0., 0.};
99           if (exp1) alphaFLAT[0] = exp1->GetAlpha();
100           if (exp2) alphaFLAT[1] = exp2->GetAlpha();
101           if (exp3) alphaFLAT[2] = exp3->GetAlpha();
102           if (exp4) alphaFLAT[3] = exp4->GetAlpha();
103
104           Float_t alphaOLD[4]  = {0., 0., 0., 0.};
105           if (esdTrack->GetConstrainedParam())  alphaOLD[0] = esdTrack->GetConstrainedParam()->GetAlpha();
106           if (esdTrack->GetInnerParam())        alphaOLD[1] = esdTrack->GetInnerParam()->GetAlpha();
107           if (esdTrack->GetTPCInnerParam())     alphaOLD[2] = esdTrack->GetTPCInnerParam()->GetAlpha();
108           if (esdTrack->GetOuterParam())        alphaOLD[3] = esdTrack->GetOuterParam()->GetAlpha();
109
110           Printf("TEST: FlatTrack %d > FlatExternalTrackParam > %p %p %p %p", idxTrack, exp1, exp2, exp3, exp4);
111           Printf("TEST: FlatTrack %d > Alpha %f %f %f %f", idxTrack, alphaFLAT[0], alphaFLAT[1], alphaFLAT[2], alphaFLAT[3]);
112           Printf("TEST: Old Track %d > Alpha %f %f %f %f", idxTrack, alphaOLD[0], alphaOLD[1], alphaOLD[2], alphaOLD[3]);
113           Printf("TEST: Diff      %d > Alpha %f %f %f %f", idxTrack, 
114                  alphaFLAT[0]-alphaOLD[0], alphaFLAT[1]-alphaOLD[1], alphaFLAT[2]-alphaOLD[2], alphaFLAT[3]-alphaOLD[3]);
115
116           Int_t nCl = track->GetNumberOfTPCClusters();
117           Printf("TEST: FlatTrack %d has %d FlatClusters", idxTrack, nCl);
118           if(nCl && useESDFriends && verbose > 1){
119             TObject* calibObject = NULL;
120             AliTPCseed* seed = NULL;
121             for (Int_t idx = 0; (calibObject = friendTrack->GetCalibObject(idx)); ++idx) {
122               if ((seed = dynamic_cast<AliTPCseed*>(calibObject))) break;
123             }
124             // -- Fill cluster
125             if (seed) {
126               Int_t idxRow2=0;
127               for (Int_t idxRow = 0; idxRow <  nCl; idxRow++){
128                 AliFlatTPCCluster * cl = track->GetTPCCluster(idxRow);
129                 cout << " idx fX fY fZ  fSigmaY2 fSigmaZ2 fCharge fQMax fPadRow" << endl;
130                 if(cl){
131                   cout << idxRow << " " << cl->GetX() << " " << cl->GetY() << " " << cl->GetZ() << " " << cl->GetSigmaY2() << " " << cl->GetSigmaZ2() << " " << cl->GetCharge() << " " << cl->GetQMax() << " " << cl->GetPadRow() << endl;
132                 }
133                 else{
134                   cout << idxRow << "---------------------------------" << endl << endl;                
135                 }
136                 AliTPCclusterMI* cl2 = NULL; 
137                 while(!cl2 && idxRow2<160){
138                   cl2 = seed->GetClusterPointer(idxRow2++);
139                 }
140                 if (cl2) {
141                   //cout<<" normalCl fX fY fZ fPadRow fSigmaY2 fSigmaZ2 fCharge fQMax" <<endl;
142                   cout << idxRow << " " << cl2->GetX() << " " << cl2->GetY() << " " << cl2->GetZ() << " " << cl2->GetSigmaY2() << " " << cl2->GetSigmaZ2() << " " << cl2->GetQ() << " " << cl2->GetMax() << " " << cl2->GetRow() << endl << endl;
143                 }
144                 else
145                   cout << idxRow << "---------------------------------" << endl << endl;        
146               }
147             }
148           }
149         
150         }
151         track = track->GetNextTrack();
152       }
153     }
154
155     outFile.write(reinterpret_cast<char*>(mem), flatEsd->GetSize());
156     delete[] mem; 
157
158   } // for (Int_t idxEvent = 1; idxEvent < 2; idxEvent++) {
159
160   outFile.close();
161
162   return;
163 }