843f0c3d3c3278c200ff163b02146eb29599167e
[u/mrichter/AliRoot.git] / HLT / TPCLib / AliHLTTPCDisplay.cxx
1 // @(#) $Id$
2 // Original: AliL3Display.cxx,v 1.26 2005/06/14 10:55:21 cvetan 
3
4 /** \class AliHLTTPCDisplay
5 <pre>
6 //_____________________________________________________________
7 // AliHLTTPCDisplay
8 //
9 // Display class for the HLT TPC events.
10 </pre>
11 */
12 // Author: Jochen Thaeder <mailto:thaeder@kip.uni-heidelberg.de>
13 //         Anders Vestbo <mailto:vestbo@fi.uib.no>      
14 //*-- Copyright &copy ALICE HLT Group 
15
16 #define TRACKHELIX 0
17 #define TRACKPOLYMARKER 0
18 #define BACKWARD 0
19
20 #include <TView.h>
21 #include <TPolyMarker3D.h>
22 #include <TPolyLine3D.h>
23 #include <TH2.h>
24 #include <TTree.h>
25 #include <TNode.h>
26 #include <TGeometry.h>
27 #include <TShape.h>
28 #include <TParticle.h>
29 #include <TFile.h>
30 #include <THelix.h>
31 #include <TStyle.h>
32 #include <TGraph.h>
33 #include <TMultiGraph.h>
34 #include <TAttText.h>
35 #include <TAxis.h>
36
37 #if TRACKHELIX
38 #include <THelix.h>
39 #endif
40
41 #ifdef use_aliroot
42 #include <TClonesArray.h>
43 #include <AliRun.h>
44 #include <AliSimDigits.h>
45 #include <AliTPCParam.h>
46 #endif
47
48 #include "AliHLTTPCLogging.h"
49 #include "AliHLTTPCDisplay.h"
50 #include "AliHLTTPCTransform.h"
51 #include "AliHLTTPCTrack.h"
52 #include "AliHLTTPCTrackArray.h"
53 #include "AliHLTTPCSpacePointData.h"
54 #include "AliHLTTPCMemHandler.h"
55 #include "AliHLTTPCDigitReaderPacked.h"
56
57 #if __GNUC__ == 3
58 using namespace std;
59 #endif
60
61 ClassImp(AliHLTTPCDisplay)
62
63 // #############################################################################
64 void AliHLTTPCDisplay::InitDisplay(Char_t *gfile) {
65     //constructor
66     memset(fClusters,0,36*6*sizeof(AliHLTTPCSpacePointData*));
67     memset(fNcl, 0, 36*6*sizeof(UInt_t)); 
68
69     fTracks = NULL;
70     fHistrawcl = NULL;
71     fHistraw = NULL;
72     fHistpad1 = NULL;
73     fHistpad2 = NULL;
74     fHistpad3 = NULL;
75     fHistallresidualsY = NULL;   
76     fHistallresidualsZ = NULL;
77     fHistcharge = NULL;
78     fGraphresidualsY = NULL;
79     fGraphresidualsZ = NULL;
80     fGraphresidualsYLength = NULL;
81     fGraphresidualsZLength = NULL;
82
83
84     fGeom = NULL;
85 // ---------------------------------------------------
86 // In order to be backward compatible
87 // ---------------------------------------------------
88 #if BACKWARD
89     //fc1 = NULL;
90 #endif 
91 // ---------------------------------------------------
92     fNPads = 0;
93     fNTimes = 0;
94     fMinHits = 0;
95     fPtThreshold = 0.;
96     fPad = -1;
97     fPadRow = 0;
98     fSlicePadRow = 0; 
99     fSelectTrack = -1;
100     fSelectTrackSlice = 0;
101     fSelectTrackSwitch = kFALSE;
102     fSelectCluster = 0;
103
104     fMinSlice = 0;
105     fMaxSlice = 35;
106     fSlicePair = kFALSE;
107
108     SetSliceArray();
109
110     fBackColor = 1; 
111     fLineColor = 0;
112     fKeepView = kFALSE;
113
114     fSwitch3DCluster = kFALSE;
115     fSwitch3DTracks = kFALSE;
116     fSwitch3DPadRow = kFALSE;
117     fSwitch3DGeometry = kFALSE;
118
119     AliHLTTPCTransform::SetBField(0.4);
120     LoadGeometrie(gfile);
121 }
122
123
124 // #############################################################################
125 AliHLTTPCDisplay::~AliHLTTPCDisplay() {
126     //destructor
127     if(fTracks) delete fTracks;
128     fTracks = NULL;
129 }
130
131 // #############################################################################
132 Bool_t AliHLTTPCDisplay::LoadGeometrie(Char_t *gfile) {
133     if (gfile) {
134         TFile *file = TFile::Open(gfile);
135         if(!file) {
136             LOG(AliHLTTPCLog::kError,"AliHLTTPCDisplay::AliHLTTPCDisplay","File Open") <<"Geometry file " << gfile << " does not exist!"<<ENDLOG;
137             return kFALSE;
138         }
139         
140         fGeom = (TGeometry*)file->Get("AliceGeom");
141
142         file->Close();
143         delete file;
144     }
145     return kTRUE;
146 }
147
148 // #############################################################################
149 //                 SETTER
150 // #############################################################################
151 void AliHLTTPCDisplay::SetHistPadRowAxis() {
152     // Set Axis range of Histogramm, due to variable NPads per padrow
153
154     fNPads = AliHLTTPCTransform::GetNPads(fPadRow);
155     fHistrawcl->SetAxisRange(0,fNPads);
156     fHistraw->SetAxisRange(0,fNPads);
157     fHistrawcl->SetAxisRange(0,fNTimes,"Y");
158     fHistraw->SetAxisRange(0,fNTimes,"Y");
159 }
160
161 void AliHLTTPCDisplay::SetSliceArray() {
162     Int_t slice=0;
163     Int_t minSlice = fMinSlice; 
164     Int_t maxSlice = fMaxSlice; 
165     Int_t realslice = 0;
166
167     for (slice=0;slice<=35;slice++){
168         fSliceArray[slice] = kFALSE;
169     }
170
171     // Single Slice, or Range
172     if (minSlice > maxSlice) maxSlice += 17;
173         
174     for (slice=minSlice;slice<=maxSlice;slice++){
175         realslice = slice % 18;
176         fSliceArray[realslice] = kTRUE;
177         fSliceArray[realslice+18] = kTRUE;
178     }
179
180     // Pair of Slices
181     if (fSlicePair) {
182         minSlice = fMinSlice + 9;
183         maxSlice = fMaxSlice + 9;
184         
185         if (minSlice > maxSlice) maxSlice += 17;
186         
187         for (slice=minSlice;slice<=maxSlice;slice++){
188             realslice = slice % 18;
189             fSliceArray[realslice] = kTRUE;     
190             fSliceArray[realslice+18] = kTRUE;
191         }
192     }
193 }
194
195 // #############################################################################
196 //                 SETUP
197 // #############################################################################
198 void AliHLTTPCDisplay::SetupCluster(Int_t slice, Int_t patch, UInt_t nofClusters, AliHLTTPCSpacePointData* data)  {  
199
200     if (data && slice>=0 && slice<36 && patch>=0 && patch<AliHLTTPCTransform::GetNPatches()) {
201         if (fClusters[slice][patch]!=NULL) {
202             delete(fClusters[slice][patch]);
203             fClusters[slice][patch]=NULL;
204         }
205         Int_t arraysize=nofClusters*sizeof(AliHLTTPCSpacePointData);
206         fClusters[slice][patch] = (AliHLTTPCSpacePointData*)new Byte_t[arraysize];
207         if (fClusters[slice][patch]) {
208             memcpy(fClusters[slice][patch], data, arraysize);
209             fNcl[slice][patch]=nofClusters;
210         } else {
211             fNcl[slice][patch]=nofClusters;
212             LOG(AliHLTTPCLog::kError,"AliHLTTPCDisplay::SetupCluster","memory allocation") << "memory allocation failed "<<ENDLOG; 
213         }
214     } else LOG(AliHLTTPCLog::kError,"AliHLTTPCDisplay::SetupCluster","argument check") << "invalid argument "<<ENDLOG; 
215 }
216
217 // #############################################################################
218 void AliHLTTPCDisplay::SetupTracks(AliHLTTPCTrackArray *tracks) {
219     fTracks=tracks;
220
221     // Set USED cluster
222     Int_t ntracks = fTracks->GetNTracks();
223
224     for(Int_t j=0; j<ntracks; j++) {    
225         AliHLTTPCTrack *gtrack = fTracks->GetCheckedTrack(j); 
226         if(!gtrack) continue;
227
228         Int_t nHits = gtrack->GetNHits();
229         UInt_t *hitnum = gtrack->GetHitNumbers();
230
231         for(Int_t h=0; h<nHits; h++){
232           
233             UInt_t id=hitnum[h];
234             Int_t slice = (id>>25) & 0x7f;
235             Int_t patch = (id>>22) & 0x7;
236             UInt_t pos = id&0x3fffff;       
237                 
238             AliHLTTPCSpacePointData *points = fClusters[slice][patch];
239             
240             if(!points) {
241                 LOG(AliHLTTPCLog::kError,"AliHLTTPCDisplay::Draw3D","Clusterarray") <<"No points at slice "<<slice<<" patch "<<patch<<" pos "<<pos<<ENDLOG;
242                 continue;
243             }
244  
245             if(pos>=fNcl[slice][patch]) {
246                 LOG(AliHLTTPCLog::kError,"AliHLTTPCDisplay::Draw3D","Clusterarray") <<"Pos is too large: pos "<<pos <<" ncl "<<fNcl[slice][patch]<<ENDLOG;
247                 continue;
248             }
249             points[pos].fUsed = kTRUE;
250             points[pos].fTrackN = j;
251         }
252     }
253 }
254
255 // #############################################################################
256 void AliHLTTPCDisplay::SetupHist(){
257
258     Int_t maxpads = 150;
259     fNTimes = AliHLTTPCTransform::GetNTimeBins();
260
261     if ( fHistraw ){
262         delete fHistraw;
263         fHistraw = NULL;
264     }
265     if ( fHistrawcl ){
266         delete fHistrawcl;
267         fHistrawcl = NULL;
268     }
269
270     if ( fHistpad1 ){
271         delete fHistpad1;
272         fHistpad1 = NULL;
273     }
274
275     if ( fHistpad2 ){
276         delete fHistpad2;
277         fHistpad2 = NULL;
278     }
279
280     if ( fHistpad3 ){
281         delete fHistpad3;
282         fHistpad3 = NULL;
283     }
284
285     if ( fHistallresidualsY){
286         delete fHistallresidualsY;
287         fHistallresidualsY = NULL;
288     }
289
290     if ( fHistallresidualsZ){
291         delete fHistallresidualsZ;
292         fHistallresidualsZ = NULL;
293     }
294     if ( fHistcharge){
295         delete fHistcharge;
296         fHistcharge = NULL;
297     }
298
299     // Setup the histograms
300     Int_t padbinning = maxpads*10;
301     fHistraw = new TH2F("fHistraw","Selected PadRow with found Clusters;Pad #;Timebin #",maxpads,0,maxpads-1,fNTimes,0,fNTimes-1);
302     fHistrawcl = new TH1F("fHistrawcl","",padbinning,0,maxpads-1);
303     fHistpad1 = new TH1F ("fHistpad1","Selected Pad -1;Timebin #",fNTimes,0,fNTimes-1);
304     fHistpad2 = new TH1F ("fHistpad2","Selected Pad;Timebin #",fNTimes,0,fNTimes-1); 
305     fHistpad3 = new TH1F ("fHistpad3","Selected Pad +1;Timebin #",fNTimes,0,fNTimes-1);
306     fHistallresidualsY = new TH1F ("fHistallresiduals","Y Residuals of all Tracks in selected slices;residuals",5000,0,100);
307     fHistallresidualsZ = new TH1F ("fHistallresiduals","Z Residuals of all Tracks in selected slices;residuals",5000,0,100);
308     fHistcharge = new TH1F ("fHistcharge","Cluster distribution per charge;charge;#cluster",5000,0,30000);
309
310     fHistraw->SetOption("COLZ"); 
311
312     fHistallresidualsY->SetTitleSize(0.03);
313     fHistallresidualsY->GetXaxis()->SetLabelSize(0.03);
314     fHistallresidualsY->GetXaxis()->SetTitleSize(0.03);
315     fHistallresidualsY->GetYaxis()->SetLabelSize(0.03);
316     fHistallresidualsY->GetYaxis()->SetTitleSize(0.03);
317
318     fHistallresidualsZ->SetTitleSize(0.03);
319     fHistallresidualsZ->GetXaxis()->SetLabelSize(0.03);
320     fHistallresidualsZ->GetXaxis()->SetTitleSize(0.03);
321     fHistallresidualsZ->GetYaxis()->SetLabelSize(0.03);
322     fHistallresidualsZ->GetYaxis()->SetTitleSize(0.03);
323
324     fHistcharge->SetTitleSize(0.03);
325     fHistcharge->GetXaxis()->SetLabelSize(0.03);
326     fHistcharge->GetXaxis()->SetTitleSize(0.03);
327     fHistcharge->GetYaxis()->SetLabelSize(0.03);
328     fHistcharge->GetYaxis()->SetTitleSize(0.03);
329
330     fHistraw->SetTitleSize(0.03);
331     fHistraw->GetXaxis()->SetLabelSize(0.03);
332     fHistraw->GetXaxis()->SetTitleSize(0.03);
333     fHistraw->GetYaxis()->SetLabelSize(0.03);
334     fHistraw->GetYaxis()->SetTitleSize(0.03);
335
336     fHistpad1->SetTitleSize(0.03);
337     fHistpad1->GetXaxis()->SetLabelSize(0.03);
338     fHistpad1->GetXaxis()->SetTitleSize(0.03);
339     fHistpad1->GetYaxis()->SetLabelSize(0.03);
340     fHistpad1->GetYaxis()->SetTitleSize(0.03);
341
342     fHistpad2->SetTitleSize(0.03);
343     fHistpad2->GetXaxis()->SetLabelSize(0.03);
344     fHistpad2->GetXaxis()->SetTitleSize(0.03);
345     fHistpad2->GetYaxis()->SetLabelSize(0.03);
346     fHistpad2->GetYaxis()->SetTitleSize(0.03);
347
348     fHistpad3->SetTitleSize(0.03);
349     fHistpad3->GetXaxis()->SetLabelSize(0.03);
350     fHistpad3->GetXaxis()->SetTitleSize(0.03);
351     fHistpad3->GetYaxis()->SetLabelSize(0.03);
352     fHistpad3->GetYaxis()->SetTitleSize(0.03);
353
354     gStyle->SetPalette(1);
355     
356     SetHistPadRowAxis();
357 }
358
359 // ####################################################################################################
360 void AliHLTTPCDisplay::FillPadRow(Int_t patch, ULong_t dataBlock, ULong_t dataLen){
361 #if defined(HAVE_ALIRAWDATA) && defined(HAVE_ALITPCRAWSTREAM_H)
362     AliHLTTPCDigitReaderPacked* fDigitReader = new AliHLTTPCDigitReaderPacked();
363     bool readValue = true;
364     Int_t rowOffset = 0;
365
366     // Initialize RAW DATA
367     Int_t firstRow = AliHLTTPCTransform::GetFirstRow(patch);
368     Int_t lastRow = AliHLTTPCTransform::GetLastRow(patch);
369
370     // Outer sector, patches 2, 3, 4, 5 -  start counting in patch 2 with row 0
371     if ( patch >= 2 ) rowOffset = AliHLTTPCTransform::GetFirstRow( 2 );
372
373     // Initialize block for reading packed data
374     void* tmpdataBlock = (void*) dataBlock;
375     fDigitReader->InitBlock(tmpdataBlock,dataLen,firstRow,lastRow);
376
377     readValue = fDigitReader->Next();
378
379     if (!readValue){    
380         LOG(AliHLTTPCLog::kError,"AliHLTTPCDisplay::FillPadRow","Read first value") << "No value in data block" << ENDLOG;
381         return;
382     }
383     
384     // FILL PADROW 3D --- Initialize the colorbins
385     if (fSwitch3DPadRow){
386         for (UInt_t ii=0;ii < 20;ii++){
387             fbinct[ii] = 0;
388             fcolorbin[ii] = 0;
389         }
390
391         // read number of entries in colorbin
392         while ( readValue ){ 
393
394             Int_t row = fDigitReader->GetRow() + rowOffset;
395             
396             if (row == fPadRow){    
397                 UInt_t charge = fDigitReader->GetSignal();
398                 
399                 for (UInt_t ii=0;ii < 19;ii++){
400                     if ( charge > (ii*15) && charge <= ((ii*15) + 15) ) fcolorbin[ii]++;
401                 }
402                 // larger than 19 * 15  
403                 if (charge > 285 ) fcolorbin[19]++;
404             }
405
406             // read next value
407             readValue = fDigitReader->Next();
408       
409             if(!readValue) break; //No more value
410         } 
411         //Initialize fpmarr[color][3*colorbin[ii]]  
412         fpmarr[0] = new Float_t[fcolorbin[0]*3]; 
413         fpmarr[1] = new Float_t[fcolorbin[1]*3]; 
414         fpmarr[2] = new Float_t[fcolorbin[2]*3]; 
415         fpmarr[3] = new Float_t[fcolorbin[3]*3]; 
416         fpmarr[4] = new Float_t[fcolorbin[4]*3];  
417         fpmarr[5] = new Float_t[fcolorbin[5]*3]; 
418         fpmarr[6] = new Float_t[fcolorbin[6]*3]; 
419         fpmarr[7] = new Float_t[fcolorbin[7]*3]; 
420         fpmarr[8] = new Float_t[fcolorbin[8]*3]; 
421         fpmarr[9] = new Float_t[fcolorbin[9]*3]; 
422         fpmarr[10] = new Float_t[fcolorbin[10]*3]; 
423         fpmarr[11] = new Float_t[fcolorbin[11]*3]; 
424         fpmarr[12] = new Float_t[fcolorbin[12]*3]; 
425         fpmarr[13] = new Float_t[fcolorbin[13]*3]; 
426         fpmarr[14] = new Float_t[fcolorbin[14]*3]; 
427         fpmarr[15] = new Float_t[fcolorbin[15]*3]; 
428         fpmarr[16] = new Float_t[fcolorbin[16]*3]; 
429         fpmarr[17] = new Float_t[fcolorbin[17]*3]; 
430         fpmarr[18] = new Float_t[fcolorbin[18]*3]; 
431         fpmarr[19] = new Float_t[fcolorbin[19]*3]; 
432         
433         // Rewind the raw reader and fill the polymarker3D
434         fDigitReader->InitBlock(tmpdataBlock,dataLen,firstRow,lastRow);
435         
436         readValue = fDigitReader->Next();
437     } // END if (fSwitch3DPadRow)
438
439     // -- Fill Raw Data
440     while ( readValue ){ 
441
442         Int_t row = fDigitReader->GetRow() + rowOffset;
443
444         // select padrow to fill in histogramm
445         if (row == fPadRow){    
446             UChar_t pad = fDigitReader->GetPad();
447             UShort_t time = fDigitReader->GetTime();
448             UInt_t charge = fDigitReader->GetSignal();
449             Float_t xyz[3];
450             fHistraw->Fill(pad,time,charge);
451
452             if (pad == (fPad-1) ) fHistpad1->Fill(time,charge);
453             if (pad == fPad) fHistpad2->Fill(time,charge);
454             if (pad == (fPad+1) ) fHistpad3->Fill(time,charge);
455
456             if (fSwitch3DPadRow) {
457                 // Transform raw coordinates to local coordinates
458                 AliHLTTPCTransform::RawHLT2Global(xyz, fSlicePadRow, fPadRow, pad, time);
459
460                 for (UInt_t ii=0;ii < 19;ii++){
461                     if ( charge > (ii*15) && charge <= ((ii*15) + 15) ){
462                         fpmarr[ii][fbinct[ii]] = xyz[0];
463                         fpmarr[ii][fbinct[ii]+1] = xyz[1];
464                         fpmarr[ii][fbinct[ii]+2] = xyz[2];
465                         fbinct[ii] += 3;
466                     }
467                 }
468                 // larger than 19 * 15
469                 if (charge > 285 ) {
470                     fpmarr[19][fbinct[19]] = xyz[0];
471                     fpmarr[19][fbinct[19]+1] = xyz[1];
472                     fpmarr[19][fbinct[19]+2] = xyz[2];
473                     fbinct[19] += 3;
474                 }
475             } // END if (fSwitch3DPadRow)
476         
477         }
478         
479         // read next value
480         readValue = fDigitReader->Next();
481       
482         //Check where to stop:
483         if(!readValue) break; //No more value
484     } 
485     
486     if ( fDigitReader )
487         delete fDigitReader;
488     fDigitReader = NULL;
489
490     AliHLTTPCSpacePointData *points = fClusters[fSlicePadRow][patch];
491     if(!points) return;
492     Int_t npoints = fNcl[fSlicePadRow][patch];
493     
494     Float_t xyz[3];
495     for(Int_t i=0; i<npoints; i++){
496         xyz[0] = points[i].fX;
497         xyz[1] = points[i].fY;
498         xyz[2] = points[i].fZ;
499         
500         Int_t clrow = AliHLTTPCTransform::GetPadRow(xyz[0]);
501         // select padrow to fill in histogramm
502         if (clrow == fPadRow){
503             AliHLTTPCTransform::LocHLT2Raw(xyz, fSlicePadRow, fPadRow);
504             fHistrawcl->Fill(xyz[1],xyz[2]);
505         }
506     }
507 #else //! if defined(HAVE_ALIRAWDATA) && defined(HAVE_ALITPCRAWSTREAM_H)
508     HLTFatal("DigitReaderPacked not available - check your build");
509 #endif //defined(HAVE_ALIRAWDATA) && defined(HAVE_ALITPCRAWSTREAM_H)
510 }
511
512 // #############################################################################
513 void AliHLTTPCDisplay::ResetHistPadRow(){  
514     fHistraw->Reset();   
515     fHistrawcl->Reset(); 
516     fHistpad1->Reset(); 
517     fHistpad2->Reset();  
518     fHistpad3->Reset(); 
519 }
520
521 // #############################################################################
522 void AliHLTTPCDisplay::ResetHistResiduals(){  
523     fHistallresidualsY->Reset();
524     fHistallresidualsZ->Reset();
525 }
526
527 // #############################################################################
528 void AliHLTTPCDisplay::ResetHistCharge(){  
529     fHistcharge->Reset();
530 }
531
532
533 // #############################################################################
534 //                 DRAWER
535 // #############################################################################
536 void AliHLTTPCDisplay::DrawGeomSector(Int_t sector) {  
537   Char_t fname[256];
538   Int_t realsector = sector;// % 18;
539   
540   if (realsector < 10){
541     sprintf(fname,"LS0%d",realsector);
542     fGeom->GetNode(fname)->SetLineColor(fLineColor);
543     fGeom->GetNode(fname)->Draw("same");
544     sprintf(fname,"US0%d",realsector);
545     fGeom->GetNode(fname)->SetLineColor(fLineColor); 
546     fGeom->GetNode(fname)->Draw("same");
547   }
548   else {
549     sprintf(fname,"LS%d",realsector);
550     fGeom->GetNode(fname)->SetLineColor(fLineColor);
551     fGeom->GetNode(fname)->Draw("same");
552     sprintf(fname,"US%d",realsector);
553     fGeom->GetNode(fname)->SetLineColor(fLineColor); 
554     fGeom->GetNode(fname)->Draw("same");
555   }   
556 }
557 // #############################################################################
558 void AliHLTTPCDisplay::DrawHistPadRow(){  
559     Char_t title[256];
560     sprintf(title,"Selected PadRow %d with found Clusters",fPadRow);
561
562     fHistraw->SetTitle(title);
563     fHistraw->SetStats(kFALSE);
564     fHistraw->Draw("COLZ");
565
566     fHistrawcl->SetStats(kFALSE);
567     fHistrawcl->SetMarkerStyle(28);
568     fHistrawcl->SetMarkerSize(2);
569     fHistrawcl->SetMarkerColor(1);
570     fHistrawcl->Draw("psame");
571 }
572
573 // #############################################################################
574 void AliHLTTPCDisplay::DrawHistPad1(){  
575     Char_t title[256];
576     sprintf(title,"Selected Pad %d",fPad -1);
577     fHistpad1->SetStats(kFALSE);
578     fHistpad1->SetTitle(title);
579     fHistpad1->Draw();
580 }
581
582 // #############################################################################
583 void AliHLTTPCDisplay::DrawHistPad2(){  
584     Char_t title[256];
585     sprintf(title,"Selected Pad %d",fPad);
586
587     fHistpad2->SetStats(kFALSE);
588     fHistpad2->SetTitle(title);
589     fHistpad2->Draw();
590 }
591
592 // #############################################################################
593 void AliHLTTPCDisplay::DrawHistPad3(){  
594     Char_t title[256];
595     sprintf(title,"Selected Pad %d",fPad +1);
596
597     fHistpad3->SetStats(kFALSE);
598     fHistpad3->SetTitle(title);
599     fHistpad3->Draw();
600 }
601
602 // #############################################################################
603 void AliHLTTPCDisplay::DrawHistResiduals(Bool_t ySwitch){  
604     if (fSwitch3DTracks){
605         if (ySwitch){
606             // Y Residual histogram for 1 track
607
608             if (fSelectTrackSwitch){
609                 Char_t title[256];
610                 sprintf(title,"Y Residuals of Track %d in Slice %d",fSelectTrack, fSelectTrackSlice );
611
612                 TMultiGraph *mgY = new TMultiGraph();
613
614
615 //              fGraphresidualsY->SetTitle(title);      
616                 fGraphresidualsY->GetXaxis()->SetTitle("padrow");       
617                 fGraphresidualsY->GetYaxis()->SetTitle("residuals");
618 //              fGraphresidualsY->Draw("A*");
619                 fGraphresidualsY->GetXaxis()->SetLabelSize(0.02);
620                 fGraphresidualsY->GetXaxis()->SetTitleSize(0.02);
621                 fGraphresidualsY->GetYaxis()->SetLabelSize(0.02);
622                 fGraphresidualsY->GetYaxis()->SetTitleSize(0.02);
623                 fGraphresidualsYLength->SetMarkerColor(2);
624                 fGraphresidualsYLength->SetMarkerStyle(5);
625                 fGraphresidualsY->SetMarkerColor(1);
626                 fGraphresidualsY->SetMarkerStyle(3);
627
628 //              fGraphresidualsY->Draw("A*");
629 //              fGraphresidualsYLength->Draw("*");
630
631                 mgY->Add(fGraphresidualsY);
632                 mgY->Add(fGraphresidualsYLength);
633                 mgY->SetTitle(title);
634 //              mgY->GetXaxis()->SetTitle("padrow");    
635 //              mgY->GetYaxis()->SetTitle("residuals");
636                 mgY->Draw("AP");
637             }
638             // Global residuals histogram
639             else{
640                 fHistallresidualsY->SetStats(kFALSE);
641                 fHistallresidualsY->Draw();
642             }
643         }
644         else {
645             // Z Residual histogram for 1 track
646             if (fSelectTrackSwitch){
647                 Char_t title[256];
648                 sprintf(title,"Z Residuals of Track %d in Slice %d",fSelectTrack, fSelectTrackSlice );
649
650                 TMultiGraph *mgZ = new TMultiGraph();
651
652                 fGraphresidualsZ->SetTitle(title);      
653                 fGraphresidualsZ->GetXaxis()->SetTitle("padrow");       
654                 fGraphresidualsZ->GetYaxis()->SetTitle("residuals");
655                 fGraphresidualsZ->GetXaxis()->SetLabelSize(0.02);
656                 fGraphresidualsZ->GetXaxis()->SetTitleSize(0.02);
657                 fGraphresidualsZ->GetYaxis()->SetLabelSize(0.02);
658                 fGraphresidualsZ->GetYaxis()->SetTitleSize(0.02);
659 //              fGraphresidualsZLength->Draw("F*");
660 //              fGraphresidualsZ->Draw("A*");
661         
662                 mgZ->Add(fGraphresidualsZ);
663 //              mgZ->Add(fGraphresidualsZLength);
664                 mgZ->SetTitle(title);
665                 mgZ->Draw("A*");
666             }
667             // Global residuals histogram
668             else{
669                 fHistallresidualsZ->SetStats(kFALSE);
670                 fHistallresidualsZ->Draw();
671             }
672         }
673     }
674 }
675
676 // #############################################################################
677 void AliHLTTPCDisplay::DrawHistCharge(){  
678     if (fSwitch3DCluster){
679 //      fHistcharge->SetStats(kFALSE);
680         fHistcharge->Draw();
681     }
682 }
683
684 // #############################################################################
685 void AliHLTTPCDisplay::Draw3D(){        
686     
687     TView *v = new TView(1);
688     v->SetRange(-800,-800,-800,800,800,800);
689
690     Float_t* etaRange = NULL;   // ------  STILL TO FIX
691     
692
693     //--------------------------------------------------------------------------------------------
694     // DRAW 3D CLUSTER
695     //--------------------------------------------------------------------------------------------
696     if (fSwitch3DCluster){
697         Int_t maxCharge = 0;
698
699         for (Int_t slice=0; slice <= 35; slice++){
700
701             Int_t currenttrack = -1;
702
703             if (fSelectCluster == 1 && fSelectTrackSwitch && slice == fSelectTrackSlice ){
704
705                 Int_t trackcounter = 0;
706                 Int_t ntracks = fTracks->GetNTracks();
707         
708                 for(Int_t j=0; j<ntracks; j++) {        
709
710                     AliHLTTPCTrack *gtrack = fTracks->GetCheckedTrack(j); 
711                     if(!gtrack) continue;
712
713                     Int_t nHits = gtrack->GetNHits();  // Number of associated hits to track
714                     Int_t tmpslice = gtrack->GetSector();
715
716                     // --- CHECK if track is should be drawn
717                     // select Single Track
718                     if(tmpslice != fSelectTrackSlice) continue;
719                         
720                     if (trackcounter != fSelectTrack){
721                         trackcounter++;  
722                         continue;
723                     }
724                     trackcounter++;
725
726                     if((fPtThreshold > 0) && (gtrack->GetPt()< fPtThreshold)) continue;
727                     if(nHits < fMinHits) continue;
728
729                     currenttrack = j;
730                     break;
731                 }
732             }
733             
734             if (!fSliceArray[slice]) continue;
735             
736             for(Int_t p=0;p<6;p++){
737
738                 AliHLTTPCSpacePointData *points = fClusters[slice][p];
739                 if(!points) continue;
740                 Int_t npoints = fNcl[slice][p];
741                 TPolyMarker3D *pmUsed = new TPolyMarker3D(1,6);
742                 TPolyMarker3D *pmUnused = new TPolyMarker3D(1,6);
743                 Int_t nUsedCluster = 0;
744                 Int_t nUnusedCluster = 0;
745
746                 Float_t xyz[3];
747                 for(Int_t i=0; i<npoints; i++){
748                     // Used  cluster only
749                     if (fSelectCluster == 1  && points[i].fUsed == kFALSE) continue; 
750                     // Unused cluster only
751                     if (fSelectCluster == 2  && points[i].fUsed == kTRUE) continue; 
752
753                     // if single track is selcted draw only cluster for this track
754                     if (fSelectCluster == 1 && fSelectTrackSwitch && points[i].fTrackN != currenttrack) continue;
755                     
756                     xyz[0] = points[i].fX;
757                     xyz[1] = points[i].fY;
758                     xyz[2] = points[i].fZ;
759                     
760                     if ( etaRange ){              
761                         // Do this before the transform, because the tracker also uses
762                         // local coordinates when using this limit to determine 
763                         // which clusters to use for tracking
764                         Double_t pointEta = AliHLTTPCTransform::GetEta( xyz );
765                         if ( pointEta<etaRange[0] || pointEta>etaRange[1] )
766                             continue;
767                     }
768
769                     AliHLTTPCTransform::Local2Global(xyz,slice);
770                  
771                     if (points[i].fUsed == kTRUE){
772                         pmUsed->SetPoint(nUsedCluster,xyz[0],xyz[1],xyz[2]);
773                         nUsedCluster++;
774                     }
775                     else {
776                         pmUnused->SetPoint(nUnusedCluster,xyz[0],xyz[1],xyz[2]);
777                         nUnusedCluster++;
778                     }
779
780                     // Fill Charge Histogram
781                     fHistcharge->Fill(points[i].fCharge);
782                     if ((Int_t)points[i].fCharge > maxCharge ) maxCharge = (Int_t) points[i].fCharge; 
783                 }
784                 pmUsed->SetMarkerSize(1);
785                 pmUsed->SetMarkerColor(3); 
786                 pmUsed->Draw("");
787
788                 pmUnused->SetMarkerSize(1);
789                 pmUnused->SetMarkerColor(2); 
790                 pmUnused->Draw("");
791             } // END - PATCH LOOP           
792         }  // END - SLICE LOOP
793         fHistcharge->SetAxisRange(0,maxCharge);
794     }   // END - DRAW 3D CLUSTER 
795
796     //--------------------------------------------------------------------------------------------
797     // DRAW 3D TRACKS
798     //--------------------------------------------------------------------------------------------
799     if (fSwitch3DTracks){
800
801         Int_t trackcounter = 0;
802         Int_t ntracks = fTracks->GetNTracks();
803         Double_t drawStep = 0.2;
804
805         Double_t maxResidualY = 0.;
806         Double_t maxResidualZ = 0.;
807
808         TPolyLine3D *line = new TPolyLine3D[ntracks];
809 #if TRACKHELIX
810         THelix *helix = new THelix[ntracks];
811 #endif
812         for(Int_t j=0; j<ntracks; j++) {        
813
814             AliHLTTPCTrack *gtrack = fTracks->GetCheckedTrack(j); 
815             if(!gtrack) continue;
816
817             Int_t nHits = gtrack->GetNHits();  // Number of associated hits to track
818             Int_t slice = gtrack->GetSector();
819
820             // --- CHECK if track is should be drawn
821             // select if slice should be displayed or not
822             if (!fSliceArray[slice]) continue;  
823             
824             // select Single Track
825             if (fSelectTrackSwitch){
826                 if(slice != fSelectTrackSlice) continue;
827         
828                 if (trackcounter != fSelectTrack){
829                     trackcounter++;  
830                     continue;
831                 }
832                 trackcounter++;
833             }
834     
835             if((fPtThreshold > 0) && (gtrack->GetPt()< fPtThreshold)) continue;
836             if(nHits < fMinHits) continue;
837             
838             TPolyMarker3D *pmL = new TPolyMarker3D(1,2);
839             TPolyMarker3D *pmF = new TPolyMarker3D(1,2);
840
841             Double_t radius = gtrack->GetRadius();      // radius
842             Double_t kappa = gtrack->GetKappa();        // curvature = 1/R , signed
843             Double_t lambda = atan( gtrack->GetTgl() ); // dipAngle lambda
844             Double_t phi0 = gtrack->GetPsi() + (gtrack->GetCharge() * AliHLTTPCTransform::PiHalf() ); // azimuthal angle of startingpoint, with respect to helix axis
845
846             Double_t xyzL[3];      // lastpoint of track
847             Double_t xyzF[3];      // firstpoint of track
848
849             xyzF[0] = gtrack->GetFirstPointX();
850             xyzF[1] = gtrack->GetFirstPointY();
851             xyzF[2] = gtrack->GetFirstPointZ();
852             pmF->SetPoint(0,xyzF[0],xyzF[1],xyzF[2]);
853
854             xyzL[0] = gtrack->GetLastPointX();
855             xyzL[1] = gtrack->GetLastPointY();
856             xyzL[2] = gtrack->GetLastPointZ();
857             pmL->SetPoint(0,xyzL[0],xyzL[1],xyzL[2]);
858
859             Double_t s = 0.;       // length of the track
860
861             // Calculate the length of the track. If it is to flat in in s,z plane use sxy, otherwise use sz
862             if (fabs(lambda) > 0.05){
863                 // length of track calculated out of z
864                 s = fabs( (xyzL[2] - xyzF[2]) / sin(lambda) ); // length of track calculated out of z
865             }
866             else {
867                 Double_t d =  (xyzL[0] - xyzF[0])*(xyzL[0] - xyzF[0]) + (xyzL[1] - xyzF[1])*(xyzL[1] - xyzF[1]);
868                 // length of track calculated out of xy
869                 s = fabs ( acos( 0.5 * (2 - (d / (radius*radius)))) / ( kappa * cos(lambda) ) );                
870             }
871             
872             Int_t nTrackPoints = 2 + (Int_t) floor(s / drawStep);
873
874 #if TRACKPOLYMARKER
875             TPolyMarker3D *pmT = new TPolyMarker3D(nTrackPoints,6);
876 #endif
877
878             Double_t *xT = new Double_t[nTrackPoints];
879             Double_t *yT = new Double_t[nTrackPoints];
880             Double_t *zT = new Double_t[nTrackPoints];
881
882             //Write Track Parameters for single track
883             if (fSelectTrackSwitch){
884                 fTrackParam.id = trackcounter - 1;
885                 fTrackParam.nHits = nHits;
886                 fTrackParam.charge = gtrack->GetCharge();
887                 fTrackParam.lambda = lambda;
888                 fTrackParam.kappa = kappa;
889                 fTrackParam.radius = radius;
890                 fTrackParam.slice = slice;
891                 fTrackParam.phi0 = phi0;
892                 fTrackParam.pt = gtrack->GetPt();
893                 fTrackParam.bfield = AliHLTTPCTransform::GetBFieldValue();
894                 fTrackParam.xyzF[0] = gtrack->GetFirstPointX();
895                 fTrackParam.xyzF[1] = gtrack->GetFirstPointY();
896                 fTrackParam.xyzF[2] = gtrack->GetFirstPointZ();
897                 fTrackParam.xyzL[0] = gtrack->GetLastPointX();
898                 fTrackParam.xyzL[1] = gtrack->GetLastPointY();
899                 fTrackParam.xyzL[2] = gtrack->GetLastPointZ();
900                 fTrackParam.psi = gtrack->GetPsi();
901                 fTrackParam.s = s;
902             }
903
904             Int_t trackPointCounter = 0;
905
906             for (Double_t ds = 0.; ds < s; ds = ds + drawStep){
907                 // FILL ARRAYS IN ORDER TO DRAW THE TRACKPOINTS, OUT OF THE PARAMETER
908                 xT[trackPointCounter] = xyzF[0] + radius * ( cos( phi0 + (ds*kappa*cos(lambda)) ) - cos(phi0) );
909                 yT[trackPointCounter] = xyzF[1] + radius * ( sin( phi0 + (ds*kappa*cos(lambda)) ) - sin(phi0) );
910                 zT[trackPointCounter] = xyzF[2] + ds * sin(lambda);
911 #if TRACKPOLYMARKER
912                 pmT->SetPoint(trackPointCounter,xT[trackPointCounter],yT[trackPointCounter],zT[trackPointCounter]);
913 #endif
914                 trackPointCounter++;
915             }
916
917             xT[trackPointCounter] = xyzF[0] + radius * ( cos( phi0 + (s*kappa*cos(lambda)) ) - cos(phi0) );
918             yT[trackPointCounter] = xyzF[1] + radius * ( sin( phi0 + (s*kappa*cos(lambda)) ) - sin(phi0) );
919             zT[trackPointCounter] = xyzF[2] + s * sin(lambda);
920 #if TRACKPOLYMARKER       
921             pmT->SetPoint(trackPointCounter,xT[trackPointCounter],yT[trackPointCounter],zT[trackPointCounter]);
922 #endif
923             // --- RESIDUALS ---
924             gtrack->Rotate(slice,kTRUE);
925             Int_t nRes = 0;  // number of resiudals
926             
927             UInt_t *hitnum = gtrack->GetHitNumbers();
928
929             Double_t *resY= new Double_t[nHits];
930             Double_t *resZ= new Double_t[nHits];
931
932             Double_t *resYLength= new Double_t[2*nHits];
933             Double_t *resZLength= new Double_t[2*nHits];
934
935             Double_t *padrows = new Double_t[nHits];
936             Double_t *padrowsLength = new Double_t[2*nHits];
937
938             for(Int_t h=0; h<nHits; h++){
939                 UInt_t id=hitnum[h];
940                 Int_t patch = (id>>22) & 0x7;
941                 UInt_t pos = id&0x3fffff; 
942
943                 AliHLTTPCSpacePointData *points = fClusters[slice][patch];
944
945                 Float_t xyzCtmp[3];    // cluster tmp
946                 Float_t xyzTtmp[3];    // track tmp
947
948                 xyzCtmp[0] = points[pos].fX;
949                 xyzCtmp[1] = points[pos].fY;
950                 xyzCtmp[2] = points[pos].fZ;
951
952                 Int_t padrow = AliHLTTPCTransform::GetPadRow(points[pos].fX);
953                 xyzTtmp[0] = gtrack->GetFirstPointX();
954                 if(gtrack->GetCrossingPoint(padrow,xyzTtmp)) {
955
956                     Float_t deltaY = ( xyzCtmp[1] - xyzTtmp[1] );
957                     Float_t deltaZ = ( xyzCtmp[2] - xyzTtmp[2] );
958 //                  Float_t residual = sqrt( deltaY*deltaY + deltaZ*deltaZ );
959                     
960                     padrows[nRes] = (Double_t) padrow;
961                     resY[nRes] = (Double_t) deltaY;
962                     resZ[nRes] = (Double_t) deltaZ;
963
964                     resYLength[(2*nRes)] = 0.5 * AliHLTTPCTransform::GetPadLength(padrow);
965                     resYLength[(2*nRes)+1] = -0.5 * AliHLTTPCTransform::GetPadLength(padrow);
966                     resZLength[nRes] = AliHLTTPCTransform::GetZLength();
967                     padrowsLength[(2*nRes)] = (Double_t) padrow;
968                     padrowsLength[(2*nRes)+1] = (Double_t) padrow;
969
970                     // FILL RESIDUALS HISTOGRAM
971                     fHistallresidualsY->Fill(resY[nRes]);
972                     fHistallresidualsZ->Fill(resZ[nRes]);
973                     if (resY[nRes] > maxResidualY ) maxResidualY = resY[nRes];
974                     if (resZ[nRes] > maxResidualZ ) maxResidualZ = resZ[nRes];
975                     nRes++;
976                 }
977             }
978
979             gtrack->Rotate(slice,kFALSE);
980             // --- RESIDUALS ---
981
982             // Draw last point of Track
983             pmL->SetMarkerSize(3);
984             pmL->SetMarkerColor(4); 
985             pmL->Draw();
986
987             // Draw first point of Track
988             pmF->SetMarkerSize(3);
989             pmF->SetMarkerColor(5); 
990             pmF->Draw();
991
992 #if TRACKPOLYMARKER
993             // Draw Track -- as polymarker
994             pmT->SetMarkerSize(3);
995             pmT->SetMarkerColor(3); 
996             pmT->Draw();
997 #endif
998             // Draw Track -- as line
999             TPolyLine3D *currentline = &(line[j]);
1000             currentline = new TPolyLine3D(nTrackPoints,xT,yT,zT,"");
1001             currentline->SetLineColor(4);   
1002             currentline->SetLineWidth(2);
1003             currentline->Draw("same");
1004
1005 #if TRACKHELIX
1006             // Draw Track -- as helix
1007             // works ok, execpt for very small dipangles -> track almost horizontal
1008             Double_t hrange[2];
1009             Double_t v0[3];
1010             Double_t omega;
1011             hrange[0] = xyzF[2];
1012             hrange[1] = xyzL[2];
1013             v0[0] = gtrack->GetPx();
1014             v0[1] = gtrack->GetPy();
1015             v0[2] = gtrack->GetPz();
1016             omega = AliHLTTPCTransform::GetBFieldValue() * gtrack->GetCharge();
1017
1018             THelix *currenthelix = &(helix[j]);
1019             currenthelix = new THelix(xyzF,v0,omega,hrange,kHelixZ,0);
1020             currenthelix->SetLineColor(6);   
1021             currenthelix->SetLineWidth(1);
1022             currenthelix->Draw("same");             
1023 #endif
1024
1025             //Residuals
1026             if ( fGraphresidualsY){
1027                 delete fGraphresidualsY;
1028                 fGraphresidualsY = NULL;
1029             }
1030
1031             if ( fGraphresidualsZ){
1032                 delete fGraphresidualsZ;
1033                 fGraphresidualsZ = NULL;
1034             }
1035             //Residuals
1036             if ( fGraphresidualsYLength){
1037                 delete fGraphresidualsYLength;
1038                 fGraphresidualsYLength = NULL;
1039             }
1040
1041             if ( fGraphresidualsZLength){
1042                 delete fGraphresidualsZLength;
1043                 fGraphresidualsZLength = NULL;
1044             }
1045
1046
1047
1048             // FILL Y RESIDUALS GRAPH
1049             fGraphresidualsY = new TGraph(nRes-1,padrows,resY);
1050             fGraphresidualsYLength = new TGraph((2*nRes)-2,padrowsLength,resYLength);
1051             // FILL Z RESIDUALS GRAPH
1052             fGraphresidualsZ = new TGraph(nRes-1,padrows,resZ);
1053             fGraphresidualsZLength = new TGraph(nRes-1,padrows,resZLength);
1054
1055             if (xT) delete xT;
1056             if (yT) delete yT;
1057             if (zT) delete zT;
1058
1059         } // END for tracks
1060
1061         fHistallresidualsY->SetAxisRange(0,maxResidualY);
1062         fHistallresidualsZ->SetAxisRange(0,maxResidualZ);
1063
1064     }   // END - DRAW 3D Tracks
1065
1066     //--------------------------------------------------------------------------------------------
1067     // DRAW 3D GEOMETRY
1068     //--------------------------------------------------------------------------------------------
1069     if (fSwitch3DGeometry){
1070
1071         for (Int_t slice=0; slice <= 17; slice++){
1072             if (!fSliceArray[slice]) continue;
1073             DrawGeomSector(slice);
1074         }
1075     }   // END - DRAW 3D GEOMETRY
1076     
1077     //--------------------------------------------------------------------------------------------
1078     // DRAW 3D PadRow
1079     //--------------------------------------------------------------------------------------------
1080     if (fSwitch3DPadRow && fSliceArray[fSlicePadRow]){
1081         Int_t markercolor = 51;
1082
1083         for (UInt_t ii=0;ii < 20;ii++){
1084             if (fcolorbin[ii]> 0){
1085                 
1086                 TPolyMarker3D *pm = new TPolyMarker3D(fcolorbin[ii], fpmarr[ii], 7 );
1087
1088                 pm->SetMarkerColor(markercolor); 
1089                 pm->Draw(""); 
1090             }
1091
1092             // in order to have the SetPalette(1), so called "pretty"
1093             if (ii % 2 == 0 ) markercolor += 2;
1094             else  markercolor += 3;
1095         }
1096     }
1097
1098     //--------------------------------------------------------------------------------------------
1099     // DRAW 3D 
1100     //--------------------------------------------------------------------------------------------
1101     v->ZoomView(0,4);
1102     v->Draw();   
1103 }
1104
1105 // ---------------------------------------------------
1106 // In order to be backward compatible
1107 // ---------------------------------------------------
1108 #if BACKWARD
1109 void AliHLTTPCDisplay::DisplayClusters(Bool_t x3don,Float_t* etaRange) {
1110     if (!fc1){
1111         fc1 = new TCanvas("c1","",900,900);
1112         fc1->cd();
1113     }
1114
1115     fSwitch3DTracks = kFALSE; 
1116     fSwitch3DCluster = kTRUE; 
1117     fSwitch3DPadRow = kFALSE; 
1118     fSwitch3DGeometry = kFALSE;
1119
1120     Draw3D();
1121 }
1122 // ---------------------------------------------------
1123 void AliHLTTPCDisplay::DisplayTracks(Int_t minhits,Bool_t x3don,Float_t thr) {
1124     if (!fc1){
1125         fc1 = new TCanvas("c1","",900,900);
1126         fc1->cd();
1127     }
1128
1129     fMinHits = minhits; 
1130     fPtThreshold = thr;
1131     fSwitch3DTracks = kTRUE; 
1132     fSwitch3DCluster = kFALSE; 
1133     fSwitch3DPadRow = kFALSE; 
1134     fSwitch3DGeometry = kFALSE;
1135
1136     Draw3D();
1137 }
1138 // ---------------------------------------------------
1139 void AliHLTTPCDisplay::DisplayAll(Int_t minhits,Bool_t clusterswitch,Bool_t trackswitch,Bool_t x3don, Float_t thr, Float_t* etaRange){
1140     if (!fc1){
1141         fc1 = new TCanvas("c1","",900,900);
1142         fc1->cd();
1143     }
1144
1145     fMinHits = minhits; 
1146     fPtThreshold = thr;
1147     fSwitch3DTracks = trackswitch; 
1148     fSwitch3DCluster = clusterswitch; 
1149     fSwitch3DPadRow = kFALSE; 
1150     fSwitch3DGeometry = kFALSE;
1151
1152     Draw3D();
1153 }
1154 #endif
1155 // ---------------------------------------------------