Update timestamps for new AMANDA simulation (17/02/2015)
[u/mrichter/AliRoot.git] / TPC / AliTPCCompareTracks.C
1 /// \file AliTPCCompareTracks.C
2
3 #ifndef __CINT__
4 #include "alles.h"
5 #include "AliComplexCluster.h"
6 //#include "AliTPCclusterM.h"
7
8 #include "AliTPCclusterMI.h"
9 #endif
10
11 Int_t AliTPCCompareTracks(Int_t eventn, Bool_t all = kFALSE) {
12    cerr<<"Comparing tracks...\n";
13    //CONNECT FILES
14    TFile *file=TFile::Open("galice.root");
15    if (!file->IsOpen()) {cerr<<"Can't open galice.root !\n"; return 1;}
16    //
17    TFile *ftracks=TFile::Open("AliTPCtracks.root","update");
18    if (!ftracks->IsOpen()){cerr<<"Can't open AliTPCtracks.root !\n"; return 3;}
19    //
20    AliTPCParam *param=(AliTPCParam *)file->Get("75x40_100x60_150x60");
21    if (!param) {cerr<<"TPC parameters have not been found !\n"; return 2;}   
22    //
23    TFile *fout=TFile::Open("AliTPCTracksDif.root","new"); 
24    if (!fout->IsOpen()) fout=TFile::Open("AliTPCClustersDif.root","recreate");  
25    //
26    //connect exact clusters 
27    file->cd();
28    AliTPCClustersArray *ce=new AliTPCClustersArray;
29    ce->Setup(param);
30    ce->SetClusterType("AliComplexCluster");
31    char  cname[100];
32    sprintf(cname,"TreeCExact_%d",eventn);
33    ce->ConnectTree(cname);
34    //
35    //connect reconstructed tracks
36    ftracks->cd();
37    sprintf(cname,"Seeds");
38    TTree * treetracks = (TTree*)ftracks->Get(cname);
39    TBranch * branchtracks = treetracks->GetBranch("seeds");   
40    //
41    //load seeds to the memory
42    Int_t trackmap[500000][4];  // map of tracks corresponding to given track number
43    memset(trackmap,0,sizeof(Int_t)*4*500000);
44    Int_t ntracks = treetracks->GetEntries();
45    TObjArray * trackarr= new TObjArray(ntracks);
46    Int_t nproces = TMath::Min(ntracks,4000);
47
48    //   for (Int_t itrack =0; itrack<ntracks; itrack++){
49    for (Int_t itrack =0; itrack<nproces; itrack++){
50      AliTPCseed * seed = new AliTPCseed;  
51      //
52      seed->fPoints = new TClonesArray("AliTPCTrackPoint",200);
53      branchtracks->SetAddress(&seed); 
54      branchtracks->GetEntry(itrack);
55
56      //if (seed->fRemoval>0 && (itrack%4) ) continue;
57      trackarr->AddLast(seed);    
58
59      //crete array with exact position information
60      seed->fEPoints = new TClonesArray("AliTPCTrackPointRef",1);
61      seed->fEPoints->ExpandCreateFast(200);
62
63      Int_t label = TMath::Abs(seed->GetLabel());
64      Int_t i;
65      if (label>500000) {
66        //too big track label
67      }else{
68        for (i=0;i<4;i++) {
69          if ( trackmap[label][i]==0) 
70            break;
71        }
72        if(i<4) trackmap[label][i]=itrack;
73      }
74    }
75
76    //add information about exact positions   
77    Int_t nrows=Int_t(ce->GetTree()->GetEntries());
78
79    for (Int_t n=0; n<nrows; n++) {
80        AliSegmentID *se=ce->LoadEntry(n);
81        Int_t sec,row;
82        param->AdjustSectorRow(se->GetID(),sec,row);
83        //
84        AliTPCClustersRow* clrowe = ce->GetRow(sec,row);
85        //
86        Int_t ncl=clrowe->GetArray()->GetEntriesFast();
87        const Int_t kNIS=param->GetNInnerSector(), kNOS=param->GetNOuterSector();
88        Int_t npads, sign;
89        if (sec < kNIS) {
90           npads = param->GetNPadsLow(row);
91           sign = (sec < kNIS/2) ? 1 : -1;
92        } else {
93           npads = param->GetNPadsUp(row);
94           sign = ((sec-kNIS) < kNOS/2) ? 1 : -1;
95        }
96        
97        Int_t trrow = row;
98        if (sec>=param->GetNInnerSector()) trrow += param->GetNRowLow(); 
99        
100        while (ncl--) {
101            AliComplexCluster *cl=(AliComplexCluster*)(clrowe->GetArray()->At(ncl));
102            Double_t x=param->GetPadRowRadii(sec,row);
103            Double_t y=cl->fY;
104            y = ( y + 0.5 -  0.5*npads) *param->GetPadPitchWidth(sec);
105            Double_t z=cl->fX*param->GetZWidth();
106            z = sign*(param->GetZLength() - z);
107            Float_t cs, sn, tmp;
108            param->AdjustCosSin(sec,cs,sn);
109            for (Int_t i=0;i<3;i++){
110              if (cl->fTracks[0]<500000) if (trackmap[cl->fTracks[0]][i]) {
111                AliTPCseed * seed = (AliTPCseed*)trackarr->At(trackmap[cl->fTracks[0]][i]);
112                TClonesArray * clarr = seed->fPoints;
113                if (!clarr){
114                  //printf("Seed %d without info",trackmap[cl->fTracks[0]][i]);
115                  continue;
116                }
117                AliTPCTrackPoint    * trcluster = (AliTPCTrackPoint*)(seed->fPoints->At(trrow));
118                AliTPCTrackPointRef * ecluster  = (AliTPCTrackPointRef*)(seed->fEPoints->At(trrow));
119                //
120                ecluster->GetCPoint() = trcluster->GetCPoint();
121                ecluster->GetTPoint() = trcluster->GetTPoint();
122                //
123                AliTPCExactPoint & epoint =  ecluster->GetExactPoint();
124                /*
125                  trcluster->fEZ = z;
126                  trcluster->fEY = y;
127                  trcluster->fAngleY = cl->fSigmaY2*param->GetPadPitchLength(sec);
128                  trcluster->fAngleZ = cl->fSigmaX2*param->GetPadPitchLength(sec);
129                  trcluster->fEAmp = cl->fQ;
130                  trcluster->fEPrim = cl->fMax;
131                */
132                epoint.fEZ = z;
133                epoint.fEY = y;
134                epoint.fEAngleY = cl->fSigmaY2*param->GetPadPitchLength(sec);
135                epoint.fEAngleZ = cl->fSigmaX2*param->GetPadPitchLength(sec);
136                epoint.fEAmp = cl->fQ;
137                epoint.fEPrim = cl->fMax;
138              }
139            } 
140        }
141        //       cr->ClearRow(sec,row);
142        ce->ClearRow(sec,row);
143    }
144
145    // make new tree - with tracks - exact position
146    fout->cd();
147    TTree * treenew = new TTree("Seedref","Seedref");
148    AliTPCseed *  ioseed  =   (AliTPCseed*)trackarr->At(0);   
149    TBranch * br = treenew->Branch("seeds","AliTPCseed",&ioseed,32000,99);
150
151    for (Int_t itrack =0; itrack<ntracks; itrack++){
152      ioseed  =  (AliTPCseed*) trackarr->At(itrack); 
153      br->SetAddress(&ioseed);
154      treenew->Fill();
155    }
156
157    treenew->Write();
158    printf("1\n");
159    delete ce;
160    printf("2\n");
161    //delete cr;
162    //   carray->GetTree()->Write();
163    printf("3\n");
164    //   delete carray;
165    printf("4\n");
166    //   cf->Close(); 
167    printf("5\n");
168    fout->Close(); 
169    printf("6\n");
170    ftracks->Close();
171    printf("7\n");
172    return 0;
173 }