]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - TRD/AliTRDQADataMakerRec.cxx
Fixes for bug #49914: Compilation breaks in trunk, and bug #48629: Trunk cannot read...
[u/mrichter/AliRoot.git] / TRD / AliTRDQADataMakerRec.cxx
index 2b95b0565bc46997e1dc780c0d13b46c12d2770d..d2caf5d5ee98f442877a1ff95c6af22eb808e09b 100644 (file)
 // --- AliRoot header files ---
 #include "AliESDEvent.h"
 #include "AliLog.h"
+#include "AliRawReader.h"
 #include "AliTRDcluster.h"
 #include "AliTRDQADataMakerRec.h"
 #include "AliTRDgeometry.h"
-#include "AliTRDdataArrayI.h"
-#include "AliTRDrawStreamTB.h"
+//#include "AliTRDdataArrayI.h"
+#include "AliTRDrawStream.h"
 
 #include "AliQAChecker.h"
 
@@ -50,7 +51,7 @@ ClassImp(AliTRDQADataMakerRec)
 
 //____________________________________________________________________________ 
   AliTRDQADataMakerRec::AliTRDQADataMakerRec() : 
-  AliQADataMakerRec(AliQA::GetDetName(AliQA::kTRD), "TRD Quality Assurance Data Maker")
+  AliQADataMakerRec(AliQAv1::GetDetName(AliQAv1::kTRD), "TRD Quality Assurance Data Maker")
 {
   //
   // Default constructor
@@ -83,7 +84,7 @@ AliTRDQADataMakerRec& AliTRDQADataMakerRec::operator=(const AliTRDQADataMakerRec
 }
 
 //____________________________________________________________________________ 
-void AliTRDQADataMakerRec::EndOfDetectorCycle(AliQA::TASKINDEX task, TObjArray * list)
+void AliTRDQADataMakerRec::EndOfDetectorCycle(AliQAv1::TASKINDEX_t task, TObjArray ** list)
 {
   //
   // Detector specific actions at end of cycle
@@ -93,155 +94,141 @@ void AliTRDQADataMakerRec::EndOfDetectorCycle(AliQA::TASKINDEX task, TObjArray *
   
   AliInfo("End of TRD cycle");
   
-  if (task == AliQA::kRECPOINTS) {
+  if (task == AliQAv1::kRECPOINTS) {
+    for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
+      TH1D * hist = new TH1D("fitHist", "", 200, -0.5, 199.5);
+      //list[specie]->Print();
     
-    TH1D *hist = new TH1D("fitHist", "", 200, -0.5, 199.5);
-    //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
-    for (Int_t i=0; i<540; i++) {
+      // fill detector map;
+      for(Int_t i = 0 ; i < 540 ; i++) {
+        Double_t v = ((TH1D*)list[specie]->At(0))->GetBinContent(i+1);
+        Int_t sm = i/30;
+        Int_t det = i%30;
+        TH2D *detMap = (TH2D*)list[specie]->At(87);
+        Int_t bin = detMap->FindBin(sm, det);
+        detMap->SetBinContent(bin, v);
+     }
+      // Rec points full chambers
+      for (Int_t i = 0 ; i < 540 ; i++) {
+        //AliInfo(Form("I = %d", i));
+        //TH1D *h = ((TH2D*)list[specie]->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[specie]->At(1))->FindBin(i,xvalue);
+          Double_t value =  ((TH2D*)list[specie]->At(1))->GetBinContent(bin);
+          hist->SetBinContent(b, value);
+        }
+        //AliInfo(Form("Sum = %d %f\n", i, hist->GetSum()));
+        if (hist->GetSum() < 100) 
+            continue; // chamber not present
       
-      //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);
+        hist->Fit("landau", "q0", "goff", 10, 180);
+        TF1 *fit = hist->GetFunction("landau");
+        ((TH1D*)list[specie]->At(12))->Fill(fit->GetParameter(1));
+        ((TH1D*)list[specie]->At(13))->Fill(fit->GetParameter(2));
       }
-      
-      //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));
-    }
-    
-    // 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
+      // 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[specie]->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.0;
+            for(Int_t det = i*30 ; det < (i+1)*30 ; det++) { // loop over detectors
+              Int_t bin = ((TH3D*)list[specie]->At(10))->FindBin(det,j,xvalue);
+              Double_t value =  ((TH3D*)list[specie]->At(10))->GetBinContent(bin);
+              svalue += value;
+            }
+            //AliInfo(Form("v = %f\n", value));
+            hist->SetBinContent(b, svalue);
+          }
        
-       //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);
-       }
+          if (hist->GetSum() < 100) 
+            continue;
+          //AliInfo(Form("fitting %d %d %f\n", i, j, hist->GetSum()));
        
-       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");
        
-       hist->Fit("landau", "q0", "goff", 10, 180);
-       TF1 *fit = hist->GetFunction("landau");
-       
-       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)));
+          TH1D *h1 = (TH1D*)list[specie]->At(14+18+i);
+          Int_t bin = h1->FindBin(j);
+          // AliInfo(Form("%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;
+      // time-bin by time-bin chamber by chamber
+      for(Int_t i = 0 ; i < 540 ; i++) {
+        //TH1D *test = ((TH3D*)list[specie]->At(10))->ProjectionZ(Form("ampTime_%d",i), i+1, i+1, 0, 35);     
+        //if (test->GetSum() < 100) continue;
       
-      //AliInfo(Form("fitting det = %d", i));
+        //AliInfo(Form("fitting det = %d", i));
       
-      for(Int_t j=0; j<35; j++) {
+        for(Int_t j = 0 ; j < 35 ; j++) {
+          //TH1D *h =  ((TH3D*)list[specie]->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[specie]->At(10))->FindBin(i,j,xvalue);
+            Double_t value =  ((TH3D*)list[specie]->At(10))->GetBinContent(bin);
+            //AliInfo(Form("v = %f\n", value));
+            hist->SetBinContent(b, value);
+          }
        
-       //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;
+          //AliInfo(Form("fitting %d %d %f\n", i, j, hist->GetSum()));
        
-       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");
        
-       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));
+          Int_t sm = i/30;
+          Int_t det = i%30;
+          TH2D *h2 = (TH2D*)list[specie]->At(14+sm);
+          Int_t bin = h2->FindBin(det,j);
+          // AliInfo(Form("%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;
     }
-    
-    if (hist) delete hist;
   }
-  
   //////////////////////////
   // 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) {
+  if (task == AliQAv1::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);
+    for (Int_t specie = 0 ; specie < AliRecoParam::kNSpecies ; specie++) {
+      for(Int_t type = 0 ; type < 2 ; type++) {
+        for(Int_t i = 0 ; i < knRatio ; i++) {
+          TH1D *ratio = (TH1D*)list[specie]->At(19 + 2*i + type);
+          TH1D *histN = (TH1D*)list[specie]->At(3 + 2*kN[i] + type);
+          TH1D *histD = (TH1D*)list[specie]->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[specie]->At(33);
+      TH1D *histD = (TH1D*)list[specie]->At(32);
+      TH1D *ratio = (TH1D*)list[specie]->At(34);
+      BuildRatio(ratio, histN, 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) ;    
+  AliQAChecker::Instance()->Run(AliQAv1::kTRD, task, list) ;    
 }
 
 //____________________________________________________________________________ 
@@ -251,7 +238,7 @@ void AliTRDQADataMakerRec::InitESDs()
   // Create ESDs histograms in ESDs subdir
   //
 
-  const Int_t kNhist = 36+5;
+  const Int_t kNhist = 36+5+4;
 
   TH1 *hist[kNhist];
   Int_t histoCounter = -1 ;
@@ -308,6 +295,29 @@ void AliTRDQADataMakerRec::InitESDs()
     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();
@@ -344,7 +354,7 @@ void AliTRDQADataMakerRec::InitRecPoints()
   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", 200, 0, 200);
+  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
@@ -431,9 +441,9 @@ void AliTRDQADataMakerRec::MakeESDs(AliESDEvent * esd)
   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();
 
@@ -494,16 +504,26 @@ void AliTRDQADataMakerRec::MakeESDs(AliESDEvent * esd)
     GetESDsData(28)->Fill(track->GetTRDBudget());
     GetESDsData(29)->Fill(track->GetTRDchi2());
     GetESDsData(30)->Fill(track->GetTRDTimBin(0));
-    GetESDsData(31)->Fill(track->GetTRDpidQuality());
+    GetESDsData(31)->Fill(track->GetTRDntrackletsPID());
+    
     
+    // 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 i=0; i<AliPID::kSPECIES; i++
-       GetESDsData(36+i)->Fill(track->GetTRDpid(i));
+      for(Int_t k=0; k<AliPID::kSPECIES; ++k
+       GetESDsData(36+k)->Fill(track->GetTRDpid(k));
     }
 
     // probabilities uniformity
-    if (track->GetTRDpidQuality() < 6) continue;
+    if (track->GetTRDntrackletsPID() < 6) continue;
     GetESDsData(35)->Fill(paramOut->GetZ()/paramOut->GetX());
     
     Int_t idx = 5 * sector + stack;
@@ -626,10 +646,10 @@ void AliTRDQADataMakerRec::MakeRaws(AliRawReader* rawReader)
   const Int_t kMCM = 16;
   //  const Int_t kADC = 22;
 
-  AliTRDrawStreamTB *raw = new AliTRDrawStreamTB(rawReader);
-
-  raw->SetRawVersion(3);
-  raw->Init();
+       rawReader->Reset() ; 
+  //AliTRDrawStreamBase::SetRawStreamVersion("TB");
+  AliTRDrawStreamBase *raw = AliTRDrawStreamBase::GetRawStream(rawReader);
+  AliDebug(2,Form("Stream version: %s", raw->IsA()->GetName()));
 
   while (raw->Next()) {