Adding simple example to load default debug streamer
[u/mrichter/AliRoot.git] / TPC / AliTPCcalibCosmic.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /*
17     Comments to be written here:
18     1. What do we calibrate.
19     2. How to interpret results
20     3. Simple example
21     4. Analysis using debug streamers.
22
23
24
25     3.Simple example
26     // To make cosmic scan the user interaction neccessary
27     //
28     .x ~/UliStyle.C
29     gSystem->Load("libANALYSIS");
30     gSystem->Load("libTPCcalib");
31     TFile fcalib("CalibObjects.root");
32     TObjArray * array = (TObjArray*)fcalib.Get("TPCCalib");
33     AliTPCcalibCosmic * cosmic = ( AliTPCcalibCosmic *)array->FindObject("cosmicTPC");
34     
35
36
37 */
38
39
40
41 #include "Riostream.h"
42 #include "TChain.h"
43 #include "TTree.h"
44 #include "TH1F.h"
45 #include "TH2F.h"
46 #include "TList.h"
47 #include "TMath.h"
48 #include "TCanvas.h"
49 #include "TFile.h"
50 #include "TF1.h"
51
52 #include "AliTPCclusterMI.h"
53 #include "AliTPCseed.h"
54 #include "AliESDVertex.h"
55 #include "AliESDEvent.h"
56 #include "AliESDfriend.h"
57 #include "AliESDInputHandler.h"
58 #include "AliAnalysisManager.h"
59
60 #include "AliTracker.h"
61 #include "AliMagFMaps.h"
62 #include "AliTPCCalROC.h"
63
64 #include "AliLog.h"
65
66 #include "AliTPCcalibCosmic.h"
67
68 #include "TTreeStream.h"
69 #include "AliTPCTracklet.h"
70
71 ClassImp(AliTPCcalibCosmic)
72
73
74 AliTPCcalibCosmic::AliTPCcalibCosmic() 
75   :AliTPCcalibBase(),
76    fGainMap(0),
77    fHistNTracks(0),
78    fClusters(0),
79    fModules(0),
80    fHistPt(0),
81    fDeDx(0),
82    fDeDxMIP(0),
83    fMIPvalue(1), 
84    fCutMaxD(5),        // maximal distance in rfi ditection
85    fCutTheta(0.03),    // maximal distan theta
86    fCutMinDir(-0.99)   // direction vector products
87 {  
88   AliInfo("Default Constructor");    
89 }
90
91
92 AliTPCcalibCosmic::AliTPCcalibCosmic(const Text_t *name, const Text_t *title) 
93   :AliTPCcalibBase(),
94    fGainMap(0),
95    fHistNTracks(0),
96    fClusters(0),
97    fModules(0),
98    fHistPt(0),
99    fDeDx(0),
100    fDeDxMIP(0),
101    fMIPvalue(1),
102    fCutMaxD(5),        // maximal distance in rfi ditection
103    fCutTheta(0.03),    // maximal distan theta
104    fCutMinDir(-0.99)   // direction vector products
105 {  
106   SetName(name);
107   SetTitle(title);
108   AliMagFMaps * field = new AliMagFMaps("dummy1", "dummy2",0,5,0);
109   AliTracker::SetFieldMap(field, kTRUE);  
110
111   fHistNTracks = new TH1F("ntracks","Number of Tracks per Event; number of tracks per event; number of tracks",501,-0.5,500.5);
112   fClusters = new TH1F("signal","Number of Clusters per track; number of clusters per track n_{cl}; counts",160,0,160);
113   fModules = new TH2F("sector","Acorde hits; z (cm); x(cm)",1200,-650,650,600,-700,700);
114   fHistPt = new TH1F("Pt","Pt distribution; p_{T} (GeV); counts",2000,0,50);
115   fDeDx = new TH2F("DeDx","dEdx; momentum p (GeV); TPC signal (a.u.)",500,0.01,100.,500,2.,1000);
116   BinLogX(fDeDx);
117   fDeDxMIP =  new TH1F("DeDxMIP","MIP region; TPC signal (a.u.);counts ",500,2.,1000);
118
119   AliInfo("Non Default Constructor");  
120   //
121 }
122
123 AliTPCcalibCosmic::~AliTPCcalibCosmic(){
124   //
125   //
126   //
127 }
128
129
130
131
132
133 void AliTPCcalibCosmic::Process(AliESDEvent *event) {
134   //
135   //
136   //
137   if (!event) {
138     Printf("ERROR: ESD not available");
139     return;
140   }  
141   AliESDfriend *ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
142   if (!ESDfriend) {
143    Printf("ERROR: ESDfriend not available");
144    return;
145   }
146    
147   FindPairs(event); // nearly everything takes place in find pairs...
148
149   if (GetDebugLevel()>20) printf("Hallo world: Im here and processing an event\n");
150   Int_t ntracks=event->GetNumberOfTracks(); 
151   fHistNTracks->Fill(ntracks);
152   if (ntracks==0) return;
153
154 }
155
156
157
158 void AliTPCcalibCosmic::Analyze() {
159
160   fMIPvalue = CalculateMIPvalue(fDeDxMIP);
161
162   return;
163
164 }
165
166
167
168 void AliTPCcalibCosmic::FindPairs(AliESDEvent *event) {
169   //
170   // Find cosmic pairs
171   // 
172   // Track0 is choosen in upper TPC part
173   // Track1 is choosen in lower TPC part
174   //
175   if (GetDebugLevel()>20) printf("Hallo world: Im here\n");
176   AliESDfriend *ESDfriend=static_cast<AliESDfriend*>(event->FindListObject("AliESDfriend"));
177   Int_t ntracks=event->GetNumberOfTracks(); 
178   TObjArray  tpcSeeds(ntracks);
179   if (ntracks==0) return;
180   Double_t vtxx[3]={0,0,0};
181   Double_t svtxx[3]={0.000001,0.000001,100.};
182   AliESDVertex vtx(vtxx,svtxx);
183   //
184   //track loop
185   //
186   for (Int_t i=0;i<ntracks;++i) {
187    AliESDtrack *track = event->GetTrack(i);
188    fClusters->Fill(track->GetTPCNcls()); 
189   
190    const AliExternalTrackParam * trackIn = track->GetInnerParam();
191    const AliExternalTrackParam * trackOut = track->GetOuterParam();
192    if (!trackIn) continue;
193    if (!trackOut) continue;
194    
195    AliESDfriendTrack *friendTrack = ESDfriend->GetTrack(i);
196    TObject *calibObject;
197    AliTPCseed *seed = 0;
198    for (Int_t l=0;(calibObject=friendTrack->GetCalibObject(l));++l) {
199      if ((seed=dynamic_cast<AliTPCseed*>(calibObject))) break;
200    }
201    if (seed) tpcSeeds.AddAt(seed,i);
202
203    Double_t meanP = 0.5*(trackIn->GetP() + trackOut->GetP());
204    if (seed && track->GetTPCNcls() > 80 + 60/(1+TMath::Exp(-meanP+5))) {
205      fDeDx->Fill(meanP, seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap));
206      //
207      if (meanP > 0.4 && meanP < 0.45) fDeDxMIP->Fill(seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap));
208      //
209      if (GetDebugLevel()>0&&meanP>0.2&&seed->CookdEdxNorm(0.0,0.45,0,0,159,fGainMap)>300) {
210        TFile *curfile = AliAnalysisManager::GetAnalysisManager()->GetTree()->GetCurrentFile();
211        if (curfile) printf(">>> p+ in file: %s \t event: %i \t Number of ESD tracks: %i \n", curfile->GetName(), (int)event->GetEventNumberInFile(), (int)ntracks);
212        if (track->GetOuterParam()->GetAlpha()<0) cout << " Polartiy: " << track->GetSign() << endl;
213      }
214
215    }
216
217   }
218
219   if (ntracks<2) return;
220   //
221   // Find pairs
222   //
223   for (Int_t i=0;i<ntracks;++i) {
224     AliESDtrack *track0 = event->GetTrack(i);     
225     // track0 - choosen upper part
226     if (!track0) continue;
227     if (!track0->GetOuterParam()) continue;
228     if (track0->GetOuterParam()->GetAlpha()<0) continue;
229     Double_t d1[3];
230     track0->GetDirection(d1);    
231     for (Int_t j=0;j<ntracks;++j) {
232       if (i==j) continue;
233       AliESDtrack *track1 = event->GetTrack(j);   
234       //track 1 lower part
235       if (!track1) continue;
236       if (!track1->GetOuterParam()) continue;
237       if (track1->GetOuterParam()->GetAlpha()>0) continue;
238       //
239       Double_t d2[3];
240       track1->GetDirection(d2);
241       
242       AliTPCseed * seed0 = (AliTPCseed*) tpcSeeds.At(i);
243       AliTPCseed * seed1 = (AliTPCseed*) tpcSeeds.At(j);
244       if (! seed0) continue;
245       if (! seed1) continue;
246       Float_t dedx0 = seed0->CookdEdxNorm(0.05,0.55,0,0,159,fGainMap);
247       Float_t dedx1 = seed1->CookdEdxNorm(0.05,0.55,0,0,159,fGainMap);
248       //
249       Float_t dedx0I = seed0->CookdEdxNorm(0.05,0.55,0,0,63,fGainMap);
250       Float_t dedx1I = seed1->CookdEdxNorm(0.05,0.55,0,0,63,fGainMap);
251       //
252       Float_t dedx0O = seed0->CookdEdxNorm(0.05,0.55,0,64,159,fGainMap);
253       Float_t dedx1O = seed1->CookdEdxNorm(0.05,0.55,0,64,159,fGainMap);
254       //
255       Float_t dir = (d1[0]*d2[0] + d1[1]*d2[1] + d1[2]*d2[2]);
256       Float_t d0  = track0->GetLinearD(0,0);
257       Float_t d1  = track1->GetLinearD(0,0);
258       //
259       // conservative cuts - convergence to be guarantied
260       // applying before track propagation
261       if (TMath::Abs(d0+d1)>fCutMaxD) continue;   // distance to the 0,0
262       if (dir>fCutMinDir) continue;               // direction vector product
263       Float_t bz = AliTracker::GetBz();
264       Float_t dvertex0[2];   //distance to 0,0
265       Float_t dvertex1[2];   //distance to 0,0 
266       track0->GetDZ(0,0,0,bz,dvertex0);
267       track1->GetDZ(0,0,0,bz,dvertex1);
268       if (TMath::Abs(dvertex0[1])>250) continue;
269       if (TMath::Abs(dvertex1[1])>250) continue;
270       //
271       //
272       //
273       Float_t dmax = TMath::Max(TMath::Abs(d0),TMath::Abs(d1));
274       AliExternalTrackParam param0(*track0);
275       AliExternalTrackParam param1(*track1);
276       //
277       // Propagate using Magnetic field and correct fo material budget
278       //
279       AliTracker::PropagateTrackTo(&param0,dmax+1,0.0005,3,kTRUE);
280       AliTracker::PropagateTrackTo(&param1,dmax+1,0.0005,3,kTRUE);
281       //
282       // Propagate rest to the 0,0 DCA - z should be ignored
283       //
284       Bool_t b0 = param0.PropagateToDCA(&vtx,bz,1000);
285       Bool_t b1 = param1.PropagateToDCA(&vtx,bz,1000);
286       //      
287       param0.GetDZ(0,0,0,bz,dvertex0);
288       param1.GetDZ(0,0,0,bz,dvertex1);
289       //
290       Double_t xyz0[3];//,pxyz0[3];
291       Double_t xyz1[3];//,pxyz1[3];
292       param0.GetXYZ(xyz0);
293       param1.GetXYZ(xyz1);
294       Bool_t isPair = IsPair(&param0,&param1);
295       //
296       if (isPair) FillAcordeHist(track0);
297       //
298       if (fStreamLevel>0){
299         TTreeSRedirector * cstream =  GetDebugStreamer();
300         printf("My stream=%p\n",(void*)cstream);
301         if (cstream) {
302           (*cstream) << "Track0" <<
303             "dir="<<dir<<               //  direction
304             "OK="<<isPair<<             //  will be accepted
305             "b0="<<b0<<                 //  propagate status
306             "b1="<<b1<<                 //  propagate status
307             "Orig0.=" << track0 <<      //  original track  0
308             "Orig1.=" << track1 <<      //  original track  1
309             "Tr0.="<<&param0<<          //  track propagated to the DCA 0,0
310             "Tr1.="<<&param1<<          //  track propagated to the DCA 0,0        
311             "v00="<<dvertex0[0]<<       //  distance using kalman
312             "v01="<<dvertex0[1]<<       // 
313             "v10="<<dvertex1[0]<<       //
314             "v11="<<dvertex1[1]<<       // 
315             "d0="<<d0<<                 //  linear distance to 0,0
316             "d1="<<d1<<                 //  linear distance to 0,0
317             //
318             "x00="<<xyz0[0]<<
319             "x01="<<xyz0[1]<<
320             "x02="<<xyz0[2]<<
321             //
322             "x10="<<xyz1[0]<<
323             "x11="<<xyz1[1]<<
324             "x12="<<xyz1[2]<<
325             //
326             "Seed0.=" << seed0 <<       //  original seed 0
327             "Seed1.=" << seed1 <<       //  original seed 1
328             //
329             "dedx0="<<dedx0<<           //  dedx0 - all
330             "dedx1="<<dedx1<<           //  dedx1 - all
331             //
332             "dedx0I="<<dedx0I<<         //  dedx0 - inner ROC
333             "dedx1I="<<dedx1I<<         //  dedx1 - inner ROC
334             //
335             "dedx0O="<<dedx0O<<         //  dedx0 - outer ROC
336             "dedx1O="<<dedx1O<<         //  dedx1 - outer ROC
337             "\n";
338         }
339       }      
340     }
341   }  
342 }    
343
344
345
346
347 void  AliTPCcalibCosmic::FillAcordeHist(AliESDtrack *upperTrack) {
348
349   // Pt cut to select straight tracks which can be easily propagated to ACORDE which is outside the magnetic field
350   if (upperTrack->Pt() < 10 || upperTrack->GetTPCNcls() < 80) return;
351     
352   const Double_t AcordePlane = 850.; // distance of the central Acorde detectors to the beam line at y =0
353   const Double_t roof = 210.5;       // distance from x =0 to end of magnet roof
354
355   Double_t r[3];
356   upperTrack->GetXYZ(r);
357   Double_t d[3];
358   upperTrack->GetDirection(d);
359   Double_t x,z;
360   z = r[2] + (d[2]/d[1])*(AcordePlane - r[1]);
361   x = r[0] + (d[0]/d[1])*(AcordePlane - r[1]);
362   
363   if (x > roof) {
364     x = r[0] + (d[0]/(d[0]+d[1]))*(AcordePlane+roof-r[0]-r[1]);
365     z = r[2] + (d[2]/(d[0]+d[1]))*(AcordePlane+roof-r[0]-r[1]);
366   }
367   if (x < -roof) {
368     x = r[0] + (d[0]/(d[1]-d[0]))*(AcordePlane+roof+r[0]-r[1]);       
369     z = r[2] + (d[2]/(d[1]-d[0]))*(AcordePlane+roof+r[0]-r[1]);
370   } 
371
372   fModules->Fill(z, x);
373  
374 }
375
376
377
378 Long64_t AliTPCcalibCosmic::Merge(TCollection *li) {
379
380   TIterator* iter = li->MakeIterator();
381   AliTPCcalibCosmic* cal = 0;
382
383   while ((cal = (AliTPCcalibCosmic*)iter->Next())) {
384     if (!cal->InheritsFrom(AliTPCcalibCosmic::Class())) {
385       Error("Merge","Attempt to add object of class %s to a %s", cal->ClassName(), this->ClassName());
386       return -1;
387     }
388     
389     fHistNTracks->Add(cal->GetHistNTracks());
390     fClusters->Add(cal-> GetHistClusters());
391     fModules->Add(cal->GetHistAcorde());
392     fHistPt->Add(cal->GetHistPt());
393     fDeDx->Add(cal->GetHistDeDx());
394     fDeDxMIP->Add(cal->GetHistMIP());
395   
396   }
397   
398   return 0;
399   
400 }
401
402
403 Bool_t  AliTPCcalibCosmic::IsPair(AliExternalTrackParam *tr0, AliExternalTrackParam *tr1){
404   //
405   //
406   /*
407   // 0. Same direction - OPOSITE  - cutDir +cutT    
408   TCut cutDir("cutDir","dir<-0.99")
409   // 1. 
410   TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03")
411   //
412   // 2. The same rphi 
413   TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5")
414   //
415   //
416   //
417   TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<1&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");  
418   // 1/Pt diff cut
419   */
420   const Double_t *p0 = tr0->GetParameter();
421   const Double_t *p1 = tr1->GetParameter();
422   if (TMath::Abs(p0[3]+p1[3])>fCutTheta) return kFALSE;
423   if (TMath::Abs(p0[0]+p1[0])>fCutMaxD)  return kFALSE;
424   Double_t d0[3], d1[3];
425   tr0->GetDirection(d0);    
426   tr1->GetDirection(d1);       
427   if (d0[0]*d1[0] + d0[1]*d1[1] + d0[2]*d1[2] >fCutMinDir) return kFALSE;
428   //
429   return kTRUE;  
430 }
431  
432
433
434 Double_t AliTPCcalibCosmic::CalculateMIPvalue(TH1F * hist) {
435
436   TF1 * funcDoubleGaus = new TF1("funcDoubleGaus", "gaus(0)+gaus(3)",0,1000);
437   funcDoubleGaus->SetParameters(hist->GetEntries()*0.75,hist->GetMean()/1.3,hist->GetMean()*0.10,
438                                 hist->GetEntries()*0.25,hist->GetMean()*1.3,hist->GetMean()*0.10);
439   hist->Fit(funcDoubleGaus);
440   Double_t MIPvalue = TMath::Min(funcDoubleGaus->GetParameter(1),funcDoubleGaus->GetParameter(4));
441
442   delete funcDoubleGaus;
443
444   return MIPvalue;
445
446 }
447
448
449
450
451 void AliTPCcalibCosmic::CalculateBetheParams(TH2F *hist, Double_t * initialParam) {
452
453   return;
454
455 }
456
457
458 void AliTPCcalibCosmic::BinLogX(TH1 *h) {
459
460   // Method for the correct logarithmic binning of histograms
461
462   TAxis *axis = h->GetXaxis();
463   int bins = axis->GetNbins();
464
465   Double_t from = axis->GetXmin();
466   Double_t to = axis->GetXmax();
467   Double_t *new_bins = new Double_t[bins + 1];
468    
469   new_bins[0] = from;
470   Double_t factor = pow(to/from, 1./bins);
471   
472   for (int i = 1; i <= bins; i++) {
473    new_bins[i] = factor * new_bins[i-1];
474   }
475   axis->Set(bins, new_bins);
476   delete new_bins;
477   
478 }
479
480 AliExternalTrackParam *AliTPCcalibCosmic::Invert(AliExternalTrackParam *input)
481 {
482   //
483   // Invert paramerameter  - not covariance yet
484   //
485   AliExternalTrackParam *output = new AliExternalTrackParam(*input);
486   Double_t * param = (Double_t*)output->GetParameter();
487   param[0]*=-1;
488   param[3]*=-1;
489   param[4]*=-1;
490   //
491   return output;
492 }
493
494
495 /*
496
497 void AliTPCcalibCosmic::dEdxCorrection(){
498   TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03");  // OK
499   TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5");     // OK
500   TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.2&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
501   TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>70");
502   TCut cutA=cutT+cutD+cutPt+cutN;
503
504
505  .x ~/UliStyle.C
506   gSystem->Load("libANALYSIS");
507   gSystem->Load("libTPCcalib");
508   gSystem->AddIncludePath("-I$ALICE_ROOT/TPC/macros");
509   gROOT->LoadMacro("$ALICE_ROOT/TPC/macros/AliXRDPROOFtoolkit.cxx+")
510  AliXRDPROOFtoolkit tool; 
511   TChain * chain = tool.MakeChain("cosmic.txt","Track0",0,1000000);
512   chain->Lookup();
513
514   .x ~/rootlogon.C
515    gSystem->Load("libSTAT.so");
516
517   Double_t chi2=0;
518   Int_t    npoints=0;
519   TVectorD fitParam;
520   TMatrixD covMatrix;
521   
522   chain->Draw("Tr0.fP[4]+Tr1.fP[4]","OK"+cutA);
523   
524   TString strFit;
525   strFit+="(Tr0.fP[1]/250)++";
526   strFit+="(Tr0.fP[1]/250)^2++";
527   strFit+="(Tr0.fP[3])++";
528   strFit+="(Tr0.fP[3])^2++";
529
530   TString * ptParam = TStatToolkit::FitPlane(chain,"Tr0.fP[4]+Tr1.fP[4]", strFit.Data(),cutA, chi2,npoints,fitParam,covMatrix) 
531
532 strFit+="(Tr0.fP[1]/250)++";
533 strFit+="(Tr0.fP[1]/250)^2++";
534 strFit+="(Tr0.fP[3])++";
535 strFit+="(Tr0.fP[3])^2++";
536 strFit+="(Tr0.fP[1]/250)^2*Tr0.fP[3]++";
537 strFit+="(Tr0.fP[1]/250)^2*Tr0.fP[3]^2++";
538 //
539
540 strFit+="sign(Tr0.fP[1])++"
541 strFit+="sign(Tr0.fP[1])*(1-abs(Tr0.fP[1]/250))"
542                                             
543 TString * thetaParam = TStatToolkit::FitPlane(chain,"Tr0.fP[3]+Tr1.fP[3]", strFit.Data(),cutA, chi2,npoints,fitParam,covMatrix)
544
545
546
547 (-0.009263+(Tr0.fP[1]/250)*(-0.009693)+(Tr0.fP[1]/250)^2*(0.009062)+(Tr0.fP[3])*(0.002256)+(Tr0.fP[3])^2*(-0.052775)+(Tr0.fP[1]/250)^2*Tr0.fP[3]*(-0.020627)+(Tr0.fP[1]/250)^2*Tr0.fP[3]^2*(0.049030))
548
549 */
550
551
552 /*
553 gSystem->Load("libANALYSIS");
554 gSystem->Load("libSTAT");
555 gSystem->Load("libTPCcalib");
556
557 TStatToolkit toolkit;
558 Double_t chi2;
559 TVectorD fitParam;
560 TMatrixD covMatrix;
561 Int_t npoints;
562 //
563 TCut cutT("cutT","abs(Tr1.fP[3]+Tr0.fP[3])<0.03");  // OK
564 TCut cutD("cutD","abs(Tr0.fP[0]+Tr1.fP[0])<5");     // OK
565 TCut cutPt("cutPt","abs(Tr1.fP[4]+Tr0.fP[4])<0.2&&abs(Tr0.fP[4])+abs(Tr1.fP[4])<10");
566 TCut cutN("cutN","min(Orig0.fTPCncls,Orig1.fTPCncls)>110");
567 TCut cutA="abs(norm-1)<0.3"+cutT+cutD+cutPt+cutN;
568
569
570
571 TTree * chain = Track0;
572
573
574 chain->SetAlias("norm","signalTot0.fElements[3]/signalTot1.fElements[3]");
575 //
576 chain->SetAlias("dr1","(signalTot0.fElements[1]/signalTot0.fElements[3])");
577 chain->SetAlias("dr2","(signalTot0.fElements[2]/signalTot0.fElements[3])");
578 chain->SetAlias("dr4","(signalTot0.fElements[4]/signalTot0.fElements[3])");
579 chain->SetAlias("dr5","(signalTot0.fElements[5]/signalTot0.fElements[3])");
580
581 TString fstring="";  
582 fstring+="dr1++";
583 fstring+="dr2++";
584 fstring+="dr4++";
585 fstring+="dr5++";
586 //
587 fstring+="dr1*dr2++";
588 fstring+="dr1*dr4++";
589 fstring+="dr1*dr5++";
590 fstring+="dr2*dr4++";
591 fstring+="dr2*dr5++";
592 fstring+="dr4*dr5++";
593
594
595
596 TString *strqdedx = toolkit.FitPlane(chain,"norm",fstring->Data(), cutA, chi2,npoints,fitParam,covMatrix,-1,0,200000);
597   
598 chain->SetAlias("corQT",strqdedx->Data());
599
600 */
601
602
603
604
605
606