Fix bug due to shadowed variables

[u/mrichter/AliRoot.git] / TRD / AliTRDQADataMakerRec.cxx

diff --git a/TRD/AliTRDQADataMakerRec.cxx b/TRD/AliTRDQADataMakerRec.cxx
index ba26e0b465fb8d8d449eeef9a4accf8a50c4dcc9..d49a257283757c0737f541b665484c8024d1752b 100644 (file)
--- a/TRD/AliTRDQADataMakerRec.cxx
+++ b/TRD/AliTRDQADataMakerRec.cxx
@@ -83,66 +83,165 @@ AliTRDQADataMakerRec& AliTRDQADataMakerRec::operator=(const AliTRDQADataMakerRec
 }
 
 //____________________________________________________________________________ 
 }
 
 //____________________________________________________________________________ 

-void AliTRDQADataMakerRec::EndOfDetectorCycle(AliQA::TASKINDEX task, TObjArray * list)
+void AliTRDQADataMakerRec::EndOfDetectorCycle(AliQA::TASKINDEX_t task, TObjArray * list)

 {
   //
   // Detector specific actions at end of cycle
   //
 {
   //
   // Detector specific actions at end of cycle
   //

-
-  //AliInfo(Form("EndOfCycle", "Fitting RecPoints %d", task));
-
- 
+  //TStopwatch watch;
+  //watch.Start();
+  
+  AliInfo("End of TRD cycle");
+  

   if (task == AliQA::kRECPOINTS) {
   if (task == AliQA::kRECPOINTS) {

-
+    
+    TH1D *hist = new TH1D("fitHist", "", 200, -0.5, 199.5);

     //list->Print();
     
     //list->Print();
     

+    // fill detector map;
+    for(int i=0; i<540; i++) {
+      Double_t v = ((TH1D*)list->At(0))->GetBinContent(i+1);
+      Int_t sm = i/30;
+      Int_t det = i%30;
+
+      TH2D *detMap = (TH2D*)list->At(87);
+      Int_t bin = detMap->FindBin(sm, det);
+      detMap->SetBinContent(bin, v);
+    }
+
+

     // Rec points full chambers
     // Rec points full chambers

-    if (((TH2D*)list->At(1))->GetEntries() > 1e4) {
-      for (Int_t i=0; i<540; i++) {
-       
-       TH1D *h = ((TH2D*)list->At(1))->ProjectionY(Form("qaTRD_recPoints_amp_%d",i), i+1, i+1);
-       if (h->GetSum() < 100) continue; // chamber not present
-       
-       h->Fit("landau", "q0", "goff", 10, 180);
-       TF1 *fit = h->GetFunction("landau");
-       ((TH1D*)list->At(12))->Fill(fit->GetParameter(1));
-       ((TH1D*)list->At(13))->Fill(fit->GetParameter(2));
-       delete h;      
+    for (Int_t i=0; i<540; i++) {
+      
+      //AliInfo(Form("I = %d", i));
+
+      //TH1D *h = ((TH2D*)list->At(1))->ProjectionY(Form("qaTRD_recPoints_amp_%d",i), i+1, i+1);
+      hist->Reset();
+      for(Int_t b=1; b<hist->GetXaxis()->GetNbins()-1; b++) {
+       Double_t xvalue = hist->GetBinCenter(b);
+       Int_t bin = ((TH2D*)list->At(1))->FindBin(i,xvalue);
+       Double_t value =  ((TH2D*)list->At(1))->GetBinContent(bin);
+       hist->SetBinContent(b, value);

       }
       }

+      
+      //printf("Sum = %d %f\n", i, hist->GetSum());
+      if (hist->GetSum() < 100) continue; // chamber not present
+      
+      hist->Fit("landau", "q0", "goff", 10, 180);
+      TF1 *fit = hist->GetFunction("landau");
+      ((TH1D*)list->At(12))->Fill(fit->GetParameter(1));
+      ((TH1D*)list->At(13))->Fill(fit->GetParameter(2));

     }
     
     }
     

-    
-    if (((TH2D*)list->At(10))->GetEntries() > 1e5) {
-      for (Int_t i=0; i<540; i++) {
+    // time-bin by time-bin sm by sm
+    for(Int_t i=0; i<18; i++) { // loop over super-modules
+      
+      for(Int_t j=0; j<35; j++) { // loop over time bins
+       
+       //TH1D *h =  ((TH3D*)list->At(10))->ProjectionZ(Form("ampTime_%d",i), i+1, i+1, j+1, j+1);     
+       hist->Reset();
+       for(Int_t b=1; b<hist->GetXaxis()->GetNbins()-1; b++) {
+         Double_t xvalue = hist->GetBinCenter(b);
+         Double_t svalue = 0;
+         
+         for(Int_t det=i*30; det<(i+1)*30; det++) { // loop over detectors
+           Int_t bin = ((TH3D*)list->At(10))->FindBin(det,j,xvalue);
+           Double_t value =  ((TH3D*)list->At(10))->GetBinContent(bin);
+           svalue += value;
+         }
+         //printf("v = %f\n", value);
+         hist->SetBinContent(b, svalue);
+       }

        
        

-       TH1D *test = ((TH3D*)list->At(10))->ProjectionZ(Form("ampTime_%d",i), i+1, i+1, 0, 35);     
-       if (test->GetSum() < 100) continue;
+       if (hist->GetSum() < 100) continue;
+       //printf("fitting %d %d %f\n", i, j, hist->GetSum());

        
        

-       //AliInfo(Form("fitting det = %d", i));
+       hist->Fit("landau", "q0", "goff", 10, 180);
+       TF1 *fit = hist->GetFunction("landau");

        
        

-       for(Int_t j=0; j<35; j++) {
-         
-         TH1D *h =  ((TH3D*)list->At(10))->ProjectionZ(Form("ampTime_%d",i), i+1, i+1, j+1, j+1);     
-         if (h->GetSum() < 50) continue;
-         
-         h->Fit("landau", "q0", "goff", 10, 180);
-         TF1 *fit = h->GetFunction("landau");
-         
-         Int_t sm = i/18;
-         Int_t det = i%18;
-         TH2D *h2 = (TH2D*)list->At(14+sm);
-         Int_t bin = h2->FindBin(det,j);
-         // printf("%d %d %d\n", det, j, bin);
-         h2->SetBinContent(bin, fit->GetParameter(1));
+       TH1D *h1 = (TH1D*)list->At(14+18+i);
+       Int_t bin = h1->FindBin(j);
+       // printf("%d %d %d\n", det, j, bin);
+       h1->SetBinContent(bin, TMath::Abs(fit->GetParameter(1)));
+      }
+    }
+      
+
+    // time-bin by time-bin chamber by chamber
+    for (Int_t i=0; i<540; i++) {
+      
+      //TH1D *test = ((TH3D*)list->At(10))->ProjectionZ(Form("ampTime_%d",i), i+1, i+1, 0, 35);     
+      //if (test->GetSum() < 100) continue;
+      
+      //AliInfo(Form("fitting det = %d", i));
+      
+      for(Int_t j=0; j<35; j++) {
+       
+       //TH1D *h =  ((TH3D*)list->At(10))->ProjectionZ(Form("ampTime_%d",i), i+1, i+1, j+1, j+1);     
+       hist->Reset();
+       for(Int_t b=1; b<hist->GetXaxis()->GetNbins()-1; b++) {
+         Double_t xvalue = hist->GetBinCenter(b);
+         Int_t bin = ((TH3D*)list->At(10))->FindBin(i,j,xvalue);
+         Double_t value =  ((TH3D*)list->At(10))->GetBinContent(bin);
+         //printf("v = %f\n", value);
+         hist->SetBinContent(b, value);

        }
        }

+       
+       if (hist->GetSum() < 100) continue;
+       //printf("fitting %d %d %f\n", i, j, hist->GetSum());
+       
+       hist->Fit("landau", "q0", "goff", 10, 180);
+       TF1 *fit = hist->GetFunction("landau");
+       
+       Int_t sm = i/30;
+       Int_t det = i%30;
+       TH2D *h2 = (TH2D*)list->At(14+sm);
+       Int_t bin = h2->FindBin(det,j);
+       // printf("%d %d %d\n", det, j, bin);
+       h2->SetBinContent(bin, TMath::Abs(fit->GetParameter(1)));
+       h2->SetBinError(bin,fit->GetParError(1));

       }
     }
       }
     }

+    
+    if (hist) delete hist;

   }
   
   }
   

-  // call the checker
-  AliQAChecker::Instance()->Run(AliQA::kTRD, task, list) ;    
+  //////////////////////////
+  // const Int_t knbits = 6;
+  // const char *suf[knbits] = {"TPCi", "TPCo", "TPCz", "TRDo", "TRDr", "TRDz"};
+  //const char *sufRatio[4] = {"TRDrTRDo", "TRDoTPCo", "TRDrTPCo", "TRDzTPCo"};

 
 

+  if (task == AliQA::kESDS) {
+    
+    const Int_t knRatio = 4;
+    const Int_t kN[knRatio] = {4,3,4,5};
+    const Int_t kD[knRatio] = {3,1,1,3}; 
+    
+    // create ratios
+    for(Int_t type=0; type<2; type++) {
+      for(Int_t i=0; i<knRatio; i++) {
+
+       TH1D *ratio = (TH1D*)list->At(19 + 2*i + type);
+       TH1D *histN = (TH1D*)list->At(3 + 2*kN[i] + type);
+       TH1D *histD = (TH1D*)list->At(3 + 2*kD[i] + type);
+
+       BuildRatio(ratio, histN, histD);
+       //ratio->Reset();
+       //ratio->Add(histN);
+       //ratio->Divide(histD);
+      }
+    }

 
 

+    // ratio for the fraction of electrons per stack
+    TH1D *histN = (TH1D*)list->At(33);
+    TH1D *histD = (TH1D*)list->At(32);
+    TH1D *ratio = (TH1D*)list->At(34);
+    BuildRatio(ratio, histN, histD);
+  }
+
+  
+  // call the checker
+  AliQAChecker::Instance()->Run(AliQA::kTRD, task, list) ;    

 }
 
 //____________________________________________________________________________ 
 }
 
 //____________________________________________________________________________ 


@@ -152,30 +251,88 @@ void AliTRDQADataMakerRec::InitESDs()
   // Create ESDs histograms in ESDs subdir
   //
 
   // Create ESDs histograms in ESDs subdir
   //
 

-  const Int_t kNhist = 19;
+  const Int_t kNhist = 36+5+4;
+

   TH1 *hist[kNhist];
   Int_t histoCounter = -1 ;
 
   TH1 *hist[kNhist];
   Int_t histoCounter = -1 ;
 

-  hist[++histoCounter] = new TH1D("qaTRD_esd_ntracks", ":Number of tracks", 300, -0.5, 299.5);
-  hist[++histoCounter] = new TH1D("qaTRD_esd_sector", ":Sector", 18, -0.5, 17.7);
+  hist[++histoCounter] = new TH1D("qaTRD_esd_ntracks", ";Number of tracks", 300, -0.5, 299.5);
+  hist[++histoCounter] = new TH1D("qaTRD_esd_sector", ";Sector", 18, -0.5, 17.7);

   hist[++histoCounter] = new TH1D("qaTRD_esd_bits", ";Bits", 64, -0.5, 63.5);
 
   const Int_t knbits = 6;
   const char *suf[knbits] = {"TPCi", "TPCo", "TPCz", "TRDo", "TRDr", "TRDz"};
 
   hist[++histoCounter] = new TH1D("qaTRD_esd_bits", ";Bits", 64, -0.5, 63.5);
 
   const Int_t knbits = 6;
   const char *suf[knbits] = {"TPCi", "TPCo", "TPCz", "TRDo", "TRDr", "TRDz"};
 

-  for(Int_t i=0; i<knbits; i++) {
-    hist[++histoCounter] = new TH1D(Form("qaTRD_esd_pt%s",suf[i]), ";p_{T} (GeV/c);", 50, 0, 10);
+  // histo = 3
+  for(Int_t i=0; i<knbits; i++) { 
+    hist[++histoCounter] = new TH1D(Form("qaTRD_esd_pt%s",suf[i]), ";p_{T} (GeV/c);", 100, 0, 10);

     hist[++histoCounter] = new TH1D(Form("qaTRD_esd_trdz%s", suf[i]), ";z (cm)", 200, -400, 400); 
   }
 
     hist[++histoCounter] = new TH1D(Form("qaTRD_esd_trdz%s", suf[i]), ";z (cm)", 200, -400, 400); 
   }
 

-  hist[++histoCounter] = new TH1D("qaTRD_esd_clsTRDo", "TRDo;number of clusters", 130, -0.5, 129.5);;
-  hist[++histoCounter] = new TH1D("qaTRD_esd_clsTRDr", "TRDr;number of clusters", 130, -0.5, 129.5);;
-  hist[++histoCounter] = new TH1D("qaTRD_esd_clsTRDz", "TRDz;number of clusters", 130, -0.5, 129.5);;
+  hist[++histoCounter] = new TH1D("qaTRD_esd_clsTRDo", "TRDo;number of clusters", 180, -0.5, 179.5);;
+  hist[++histoCounter] = new TH1D("qaTRD_esd_clsTRDr", "TRDr;number of clusters", 180, -0.5, 179.5);;
+  hist[++histoCounter] = new TH1D("qaTRD_esd_clsTRDz", "TRDz;number of clusters", 180, -0.5, 179.5);;

   //hist[++histoCounter] = new TH1D("qaTRD_esd_clsRatio", ";cluster ratio", 100, 0., 1.3);;
 
   hist[++histoCounter] = new TH2D("qaTRD_esd_sigMom", ";momentum (GeV/c);signal", 100, 0, 5, 200, 0, 1e3);
 
   //hist[++histoCounter] = new TH1D("qaTRD_esd_clsRatio", ";cluster ratio", 100, 0., 1.3);;
 
   hist[++histoCounter] = new TH2D("qaTRD_esd_sigMom", ";momentum (GeV/c);signal", 100, 0, 5, 200, 0, 1e3);
 

-  for(Int_t i=0; i<=histoCounter; i++) {
+  // end of cycle plots (hist 19)
+  const char *sufRatio[4] = {"TRDrTRDo", "TRDoTPCo", "TRDrTPCo", "TRDzTPCo"};
+
+  for(int i=0; i<4; i++) {
+    hist[++histoCounter] = new TH1D(Form("qaTRD_esd_pt%s",sufRatio[i]), 
+                                   Form("Efficiency in Pt %s;p_{T};eff", sufRatio[i]),
+                                   100, 0, 10);
+
+    hist[++histoCounter] = new TH1D(Form("qaTRD_esd_trdz%s",sufRatio[i]), 
+                                   Form("Efficiency in Z %s;z (cm);eff", sufRatio[i]),
+                                   200, -400, 400);
+  }
+
+  // 27 - 31
+  hist[27] = new TH1D("qaTRD_esd_quality", ";quality", 120, 0, 12);
+  hist[28] = new TH1D("qaTRD_esd_budget", ";NN", 110, -1000, 100);
+  hist[29] = new TH1D("qaTRD_esd_chi2", ";chi2", 200, 0, 100);
+  hist[30] = new TH1D("qaTRD_esd_timeBin", ";time bin", 7, -0.5, 6.5);
+  hist[31] = new TH1D("qaTRD_esd_pidQuality", "pid Quality", 7, -0.5, 6.5);
+
+  // stack by stack electron identyfication
+  hist[32] = new TH1D("qaTRD_esd_tracksStack", "number of all tracks;stack", 90, -0.5, 89.5);
+  hist[33] = new TH1D("qaTRD_esd_electronStack", "number of electron tracks;stack", 90, -0.5, 89.5);
+  hist[34] = new TH1D("qaTRD_esd_elRatioStack", "fraction of electron tracks;stack", 90, -0.5, 89.5);
+  hist[35] = new TH1D("qaTRD_esd_thetaOut", ";tan(theta);", 100, -1, 1);
+  
+  const char *partType[5] = {"Electron", "Muon", "Pion", "Kaon", "Proton"}; 
+
+  for(Int_t i=0; i<AliPID::kSPECIES; i++)
+    hist[36+i] = new TH1D(Form("qaTRD_esd_pid%d",i),
+                         Form("%s;probability",partType[i]), 100, 0, 1);
+ 
+  // dE/dX vs momentum in three regions
+  const char *zoneName[4] = {"total charge", "ampilification range", "plateau", "TR range"};
+ 
+  // prepare the scale from 0.1 to 10 GeV
+  const Int_t nscalex= 50;
+  Double_t scalex[nscalex+1];
+  Double_t dd = (TMath::Log(10) - TMath::Log(0.5)) / nscalex;
+  for(Int_t ix=0; ix<nscalex+1; ix++) {
+    scalex[ix] = 0.5 * TMath::Exp(dd * ix);
+  }
+
+  const Int_t nscaley = 50;
+  Double_t scaley[nscaley+1];
+  for(Int_t iy=0; iy<nscaley+1; iy++) {
+    scaley[iy] = iy * (3e3/nscaley);
+  }
+    
+
+  for(Int_t i=0; i<4; i++) {
+    hist[41+i] = new TH2D(Form("qaTRD_esd_signalPzone_%d",i), 
+                         Form("%s;momentum (GeV/c);singal (a.u.)", zoneName[i]),
+                         nscalex, scalex, nscaley, scaley);
+  }
+
+  for(Int_t i=0; i<kNhist; i++) {

     //hist[i]->Sumw2();
     Add2ESDsList(hist[i], i);
   }
     //hist[i]->Sumw2();
     Add2ESDsList(hist[i], i);
   }


@@ -189,7 +346,7 @@ void AliTRDQADataMakerRec::InitRecPoints()
   // Create Reconstructed Points histograms in RecPoints subdir
   //
 
   // Create Reconstructed Points histograms in RecPoints subdir
   //
 

-  const Int_t kNhist = 14 + 18;
+  const Int_t kNhist = 14 + 4 * 18 + 2;

   TH1 *hist[kNhist];
 
   hist[0] = new TH1D("qaTRD_recPoints_det", ";Detector ID of the cluster", 540, -0.5, 539.5);
   TH1 *hist[kNhist];
 
   hist[0] = new TH1D("qaTRD_recPoints_det", ";Detector ID of the cluster", 540, -0.5, 539.5);


@@ -206,25 +363,52 @@ void AliTRDQADataMakerRec::InitRecPoints()
   hist[9] = new TH1D("qaTRD_recPoints_nCls", ";number of clusters", 500, -0.5, 499.5);
 
   hist[10] = new TH3D("qaTRD_recPoints_sigTime", ";chamber;time bin;signal", 
   hist[9] = new TH1D("qaTRD_recPoints_nCls", ";number of clusters", 500, -0.5, 499.5);
 
   hist[10] = new TH3D("qaTRD_recPoints_sigTime", ";chamber;time bin;signal", 

-                     540, -0.5, 539.5, 35, -0.5, 34.5, 100, 0, 200);
+                     540, -0.5, 539.5, 35, -0.5, 34.5, 200, -0.5, 199.5);

   hist[11] = new TProfile("qaTRD_recPoints_prf", ";distance;center of gravity"
                          , 120, -0.6, 0.6, -1.2, 1.2, "");
 
   hist[11] = new TProfile("qaTRD_recPoints_prf", ";distance;center of gravity"
                          , 120, -0.6, 0.6, -1.2, 1.2, "");
 

-  hist[12] = new TH1D("qaTRD_recPoints_ampMPV", ";amplitude MPV", 100, 0, 100);
-  hist[13] = new TH1D("qaTRD_recPoints_ampSigma", ";amplitude Sigma", 100, 0, 100); 
-
+  hist[12] = new TH1D("qaTRD_recPoints_ampMPV", ";amplitude MPV", 150, 0, 150);
+  hist[13] = new TH1D("qaTRD_recPoints_ampSigma", ";amplitude Sigma", 200, 0, 200); 
+  
+  // chamber by chamber

   for(Int_t i=0; i<18; i++) {
     hist[14+i] = new TH2D(Form("qaTRD_recPoints_sigTime_sm%d",i), Form("sm%d;det;time bin"), 
                        30, -0.5, 29.5, 35, -0.5, 34.5);
     hist[14+i]->SetMinimum(20);
     hist[14+i]->SetMaximum(40);
   }
   for(Int_t i=0; i<18; i++) {
     hist[14+i] = new TH2D(Form("qaTRD_recPoints_sigTime_sm%d",i), Form("sm%d;det;time bin"), 
                        30, -0.5, 29.5, 35, -0.5, 34.5);
     hist[14+i]->SetMinimum(20);
     hist[14+i]->SetMaximum(40);
   }

+ 
+  // time bin by time bin sm-by-sm
+  for(Int_t i=0; i<18; i++) {
+    hist[14+18+i] = new TH1D(Form("qaTRD_recPoints_sigTimeShape_sm%d", i), 
+                            Form("sm%d;time bin;signal"),
+                            35, -0.5, 34.5);
+
+    hist[14+18+i]->SetMaximum(120);    
+  }
+
+  // str = 50
+  for(Int_t i=0; i<18; i++) {
+    hist[50+i] = new TH1D(Form("qaTRD_recPoints_nCls_sm%d",i),
+                         Form("sm%d;time bin;number of clusters",i),
+                         35, -0.5, 34.5);
+  }
+
+  // str = 68
+  for(Int_t i=0; i<18; i++) {
+    hist[68+i] = new TH1D(Form("qaTRD_recPoints_totalCharge_sm%d", i),
+                         Form("sm%d;time bin;total charge", i),
+                         35, -0.5, 34.5);
+  }
+
+  hist[86] = new TH1D("qaTRD_recPoints_signal", ";amplitude", 200, -0.5, 199.5);
+  hist[87] = new TH2D("qaTRD_recPoints_detMap", ";sm;chamber", 18, -0.5, 17.5, 30, -0.5, 29.5);
+

 
   for(Int_t i=0; i<kNhist; i++) {
     //hist[i]->Sumw2();
     Add2RecPointsList(hist[i], i);
   }
 
   for(Int_t i=0; i<kNhist; i++) {
     //hist[i]->Sumw2();
     Add2RecPointsList(hist[i], i);
   }

-

 }
 
 //____________________________________________________________________________ 
 }
 
 //____________________________________________________________________________ 


@@ -270,9 +454,9 @@ void AliTRDQADataMakerRec::MakeESDs(AliESDEvent * esd)
   GetESDsData(0)->Fill(nTracks);
 
   // track loop
   GetESDsData(0)->Fill(nTracks);
 
   // track loop

-  for (Int_t i=0; i<nTracks; i++) {
+  for (Int_t iTrack = 0; iTrack<nTracks; iTrack++) {

 
 

-    AliESDtrack *track = esd->GetTrack(i);
+    AliESDtrack *track = esd->GetTrack(iTrack);

     const AliExternalTrackParam *paramOut = track->GetOuterParam();
     const AliExternalTrackParam *paramIn = track->GetInnerParam();
 
     const AliExternalTrackParam *paramOut = track->GetOuterParam();
     const AliExternalTrackParam *paramIn = track->GetInnerParam();
 


@@ -288,7 +472,7 @@ void AliTRDQADataMakerRec::MakeESDs(AliESDEvent * esd)
 
     // .. in the acceptance
     Int_t sector = GetSector(paramOut->GetAlpha());
 
     // .. in the acceptance
     Int_t sector = GetSector(paramOut->GetAlpha());

-    GetESDsData(1)->Fill(sector);
+    Int_t stack = GetStack(paramOut);

 
     UInt_t u = 1;
     UInt_t status = track->GetStatus();
 
     UInt_t u = 1;
     UInt_t status = track->GetStatus();


@@ -317,7 +501,7 @@ void AliTRDQADataMakerRec::MakeESDs(AliESDEvent * esd)
 
     // clusters
     for(Int_t b=0; b<3; b++) 
 
     // clusters
     for(Int_t b=0; b<3; b++) 

-      if (bit[3+b]) GetESDsData(b+15)->Fill(track->GetTRDncls());
+      if (bit[3+b]) GetESDsData(b+15)->Fill(track->GetTRDncls0());

 
     // refitted only
     if (!bit[4]) continue;
 
     // refitted only
     if (!bit[4]) continue;


@@ -325,9 +509,51 @@ void AliTRDQADataMakerRec::MakeESDs(AliESDEvent * esd)
     //fQuality->Fill(track->GetTRDQuality());
     //fBudget->Fill(track->GetTRDBudget());
     //fSignal->Fill(track->GetTRDsignal());
     //fQuality->Fill(track->GetTRDQuality());
     //fBudget->Fill(track->GetTRDBudget());
     //fSignal->Fill(track->GetTRDsignal());

-       
+
+    GetESDsData(1)->Fill(sector);

     GetESDsData(18)->Fill(track->GetP(), track->GetTRDsignal());
 
     GetESDsData(18)->Fill(track->GetP(), track->GetTRDsignal());
 

+    GetESDsData(27)->Fill(track->GetTRDQuality());
+    GetESDsData(28)->Fill(track->GetTRDBudget());
+    GetESDsData(29)->Fill(track->GetTRDchi2());
+    GetESDsData(30)->Fill(track->GetTRDTimBin(0));
+    GetESDsData(31)->Fill(track->GetTRDpidQuality());
+    
+    
+    // dedx
+    for(Int_t k=0; k<4; ++k) {
+      Double_t dedx = 0;
+      for(Int_t j=0; j<6; j++) {
+       dedx += track->GetTRDslice(j, k-1);
+      }
+      GetESDsData(41+k)->Fill(paramOut->GetP(), dedx/6.);
+    }
+
+    // probabilities
+    if (status & AliESDtrack::kTRDpid) {
+      for(Int_t k=0; k<AliPID::kSPECIES; ++k) 
+       GetESDsData(36+k)->Fill(track->GetTRDpid(k));
+    }
+
+    // probabilities uniformity
+    if (track->GetTRDpidQuality() < 6) continue;
+    GetESDsData(35)->Fill(paramOut->GetZ()/paramOut->GetX());
+    
+    Int_t idx = 5 * sector + stack;
+    GetESDsData(32)->Fill(idx); // all tracks
+    if (track->GetTRDpid(AliPID::kElectron) > 0.9) 
+      GetESDsData(33)->Fill(idx); // electrons only
+
+    
+
+    /*
+    hist[27] = new TH1D("qaTRD_esd_quality", ";quality", 120, 0, 12);
+    hist[28] = new TH1D("qaTRD_esd_budget", ";NN", 110, -1000, 100);
+    hist[29] = new TH1D("qaTRD_esd_chi2", ";chi2", 300, 0, 100);
+    hist[30] = new TH1D("qaTRD_esd_timeBin", 7, -0.5, 6.5);
+    hist[31] = new TH1D("qaTRD_esd_pidQuality", 7, -0.5, 6.5);
+    */
+

     /*
     // PID only
     if (status & AliESDtrack::kTRDpid) {
     /*
     // PID only
     if (status & AliESDtrack::kTRDpid) {


@@ -367,17 +593,33 @@ void AliTRDQADataMakerRec::MakeESDs(AliESDEvent * esd)
 }
 
 //______________________________________________________________________________
 }
 
 //______________________________________________________________________________

-Int_t AliTRDQADataMakerRec::GetSector(const Double_t alpha) const 
+Int_t AliTRDQADataMakerRec::GetSector(Double_t alpha) const 

 {
   //
   // Gets the sector number 
   //
 
   Double_t size = TMath::DegToRad() * 20.; // shall use TRDgeo
 {
   //
   // Gets the sector number 
   //
 
   Double_t size = TMath::DegToRad() * 20.; // shall use TRDgeo

-  Int_t sector = (Int_t)((alpha + TMath::Pi())/size);
+  if (alpha < 0) alpha += 2*TMath::Pi();
+  Int_t sector = (Int_t)(alpha/size);

   return sector;
 
 }
   return sector;
 
 }

+//______________________________________________________________________________
+
+Int_t AliTRDQADataMakerRec::GetStack(const AliExternalTrackParam *paramOut) const
+{
+  //
+  // calculates the stack the track is in
+  //
+  
+  const Double_t L = -0.9;
+  const Double_t W = (2*L)/5;
+
+  Double_t tan = paramOut->GetZ() / paramOut->GetX();
+  Double_t pos = (tan - L) / W;
+  return (Int_t) pos;
+}

 
 //______________________________________________________________________________
 Double_t AliTRDQADataMakerRec::GetExtZ(const AliExternalTrackParam *in) const 
 
 //______________________________________________________________________________
 Double_t AliTRDQADataMakerRec::GetExtZ(const AliExternalTrackParam *in) const 


@@ -417,10 +659,9 @@ void AliTRDQADataMakerRec::MakeRaws(AliRawReader* rawReader)
   const Int_t kMCM = 16;
   //  const Int_t kADC = 22;
 
   const Int_t kMCM = 16;
   //  const Int_t kADC = 22;
 

-  AliTRDrawStreamTB *raw = new AliTRDrawStreamTB(rawReader);
-
-  raw->SetRawVersion(3);
-  raw->Init();
+  //AliTRDrawStreamBase::SetRawStreamVersion("TB");
+  AliTRDrawStreamBase *raw = AliTRDrawStreamBase::GetRawStream(rawReader);
+  AliDebug(2,Form("Stream version: %s", raw->IsA()->GetName()));

 
   while (raw->Next()) {
 
 
   while (raw->Next()) {
 


@@ -445,6 +686,8 @@ void AliTRDQADataMakerRec::MakeRaws(AliRawReader* rawReader)
     GetRawsData(3)->Fill(sm);
     GetRawsData(4+sm)->Fill(index);
   }
     GetRawsData(3)->Fill(sm);
     GetRawsData(4+sm)->Fill(index);
   }

+
+  delete raw;

 }
 
 //____________________________________________________________________________
 }
 
 //____________________________________________________________________________


@@ -453,6 +696,8 @@ void AliTRDQADataMakerRec::MakeRecPoints(TTree * clustersTree)
   //  
   // Makes data from RecPoints
   // 
   //  
   // Makes data from RecPoints
   // 

+  
+  //  Info("MakeRecPoints", "making");

 
   Int_t nsize = Int_t(clustersTree->GetTotBytes() / (sizeof(AliTRDcluster))); 
   TObjArray *clusterArray = new TObjArray(nsize+1000); 
 
   Int_t nsize = Int_t(clustersTree->GetTotBytes() / (sizeof(AliTRDcluster))); 
   TObjArray *clusterArray = new TObjArray(nsize+1000); 


@@ -486,6 +731,7 @@ void AliTRDQADataMakerRec::MakeRecPoints(TTree * clustersTree)
       Int_t iDet = c->GetDetector();
       nDet[iDet]++;
       GetRecPointsData(0)->Fill(iDet);
       Int_t iDet = c->GetDetector();
       nDet[iDet]++;
       GetRecPointsData(0)->Fill(iDet);

+      GetRecPointsData(86)->Fill(c->GetQ());

       GetRecPointsData(1)->Fill(iDet, c->GetQ());
       GetRecPointsData(2)->Fill(c->GetNPads());
       if (c->GetNPads() < 6)
       GetRecPointsData(1)->Fill(iDet, c->GetQ());
       GetRecPointsData(2)->Fill(c->GetNPads());
       if (c->GetNPads() < 6)


@@ -496,6 +742,10 @@ void AliTRDQADataMakerRec::MakeRecPoints(TTree * clustersTree)
       GetRecPointsData(8)->Fill(c->GetPadTime());
 
       ((TH3D*)GetRecPointsData(10))->Fill(iDet, c->GetPadTime(), c->GetQ());
       GetRecPointsData(8)->Fill(c->GetPadTime());
 
       ((TH3D*)GetRecPointsData(10))->Fill(iDet, c->GetPadTime(), c->GetQ());

+      
+      Int_t iSM = iDet / 30;
+      GetRecPointsData(50+iSM)->Fill(c->GetPadTime());
+      GetRecPointsData(68+iSM)->Fill(c->GetPadTime(), c->GetQ());

 
       // PRF for 2pad
       //if (c->GetNPads() == 2) {
 
       // PRF for 2pad
       //if (c->GetNPads() == 2) {


@@ -539,5 +789,37 @@ Int_t AliTRDQADataMakerRec::CheckPointer(TObject *obj, const char *name)
 
   if (!obj) AliWarning(Form("null pointer: %s", name));
   return !!obj;
 
   if (!obj) AliWarning(Form("null pointer: %s", name));
   return !!obj;

+}
+//__________________________________________________________________________
+void AliTRDQADataMakerRec::BuildRatio(TH1D *ratio, TH1D *histN, TH1D*histD) {
+  //
+  // Calculate the ratio of two histograms 
+  // error are calculated assuming the histos have the same counts
+  //
+
+  // calclate
+
+  Int_t nbins = histN->GetXaxis()->GetNbins();
+  for(Int_t i=1; i<nbins+2; i++) {
+    
+    Double_t valueN = histN->GetBinContent(i);
+    Double_t valueD = histD->GetBinContent(i);
+    
+    if (valueD < 1) {
+      ratio->SetBinContent(i, 0);
+      ratio->SetBinError(i, 0);
+      continue;
+    }

 
 

+    Double_t eps = (valueN < valueD-valueN)? valueN : valueD-valueN;
+    
+    ratio->SetBinContent(i, valueN/valueD);
+    ratio->SetBinError(i, TMath::Sqrt(eps)/valueD);
+  }
+
+  // style
+  ratio->SetMinimum(-0.1);
+  ratio->SetMaximum(1.1);
+  ratio->SetMarkerStyle(20);

 }
 }

+//__________________________________________________________________________