TRD/qaRec/AliTRDcheckDetector.cxx

   1 #include <TFile.h>
   2 #include <TH1F.h>
   3 #include <TGraph.h>
   4 #include <TMath.h>
   5 #include <TMap.h>
   6 #include <TObjArray.h>
   7 #include <TObject.h>
   8 #include <TObjString.h>
   9 #include <TProfile.h>
  10 #include <TROOT.h>
  11
  12 #include "AliTRDcluster.h"
  13 #include "AliESDHeader.h"
  14 #include "AliESDRun.h"
  15 #include "AliTRDseedV1.h"
  16 #include "AliTRDtrackV1.h"
  17 #include "TTreeStream.h"
  18
  19 #include "AliTRDtrackInfo/AliTRDtrackInfo.h"
  20 #include "AliTRDtrackInfo/AliTRDeventInfo.h"
  21 #include "AliTRDcheckDetector.h"
  22
  23 #include <cstdio>
  24
  25 ////////////////////////////////////////////////////////////////////////////
  26 //                                                                        //
  27 //  Reconstruction QA                                                     //
  28 //                                                                        //
  29 //  Task doing basic checks for tracking and detector performance         //
  30 //                                                                        //
  31 //  Authors:                                                              //
  32 //    Anton Andronic <A.Andronic@gsi.de>                                  //
  33 //    Markus Fasel <M.Fasel@gsi.de>                                       //
  34 //                                                                        //
  35 ////////////////////////////////////////////////////////////////////////////
  36
  37 //_______________________________________________________
  38 AliTRDcheckDetector::AliTRDcheckDetector():
  39   AliTRDrecoTask("DetChecker", "Basic Detector Checker")
  40   ,fEventInfo(0x0)
  41   ,fPHSdetector(0x0)
  42   ,fPHSsector(0x0)
  43   ,fQCLdetector(0x0)
  44   ,fQCLsector(0x0)
  45   ,fQTdetector(0x0)
  46   ,fQTsector(0x0)
  47   ,fTriggerNames(0x0)
  48 {
  49   //
  50   // Default constructor
  51   //
  52   DefineInput(1,AliTRDeventInfo::Class());
  53 }
  54
  55 //_______________________________________________________
  56 AliTRDcheckDetector::~AliTRDcheckDetector(){
  57   //
  58   // Destructor
  59   //
  60   if(fEventInfo) delete fEventInfo;
  61   if(fPHSdetector) delete fPHSdetector;
  62   if(fPHSsector) delete fPHSsector;
  63   if(fQCLdetector) delete fQCLdetector;
  64   if(fQCLsector) delete fQCLsector;
  65   if(fQTdetector) delete fQTdetector;
  66   if(fQTsector) delete fQTsector;
  67   if(fTriggerNames) delete fTriggerNames;
  68 }
  69
  70 //_______________________________________________________
  71 void AliTRDcheckDetector::ConnectInputData(Option_t *opt){
  72         //
  73         // Connect the Input data with the task
  74         //
  75         AliTRDrecoTask::ConnectInputData(opt);
  76         fEventInfo = dynamic_cast<AliTRDeventInfo *>(GetInputData(1));
  77 }
  78
  79 //_______________________________________________________
  80 void AliTRDcheckDetector::CreateOutputObjects(){
  81   //
  82   // Create Output Objects
  83   //
  84   fContainer = new TObjArray();
  85   // Register Histograms
  86   fContainer->Add(new TH1F("hNtrks", "Number of Tracks per event", 100, 0, 100));
  87   fContainer->Add(new TH1F("hEventsTriggerTracks", "Trigger Class (Tracks)", 20, 0, 20));
  88   fContainer->Add(new TH1F("hNcls", "Nr. of clusters per track", 30, 0, 180));
  89   fContainer->Add(new TH1F("hNtls", "Nr. tracklets per track", 10, 0, 10));
  90   fContainer->Add(new TH1F("hNclTls","Mean Number of clusters per tracklet", 30, 0, 30));
  91   fContainer->Add(new TH1F("hChi2", "Chi2", 200, 0, 20));
  92   fContainer->Add(new TH1F("hChi2n", "Norm. Chi2 (tracklets)", 50, 0, 5));
  93   fContainer->Add(new TH1F("hSM", "Track Counts in Supermodule", 18, 0., 18.));
  94   fContainer->Add(new TProfile("hPH", "Average PH", 31, -0.5, 30.5));
  95   fContainer->Add(new TH1F("hQcl1", "Cluster charge", 200, -1200, 1200));
  96   fContainer->Add(new TH1F("hQT", "Total Charge Deposit", 6000, 0, 6000));
  97   fContainer->Add(new TH1F("hEventsTrigger", "Trigger Class", 20, 0, 20));
  98   // Initialise the PHS, cluster charge and total charge deposit histograms for each detector
  99   fPHSdetector = new TObjArray();
 100   fPHSdetector->SetName("fPHSdetector");
 101   fQCLdetector = new TObjArray();
 102   fQCLdetector->SetName("fQCLdetector");
 103   fQTdetector = new TObjArray();
 104   fQTdetector->SetName("fQTdetector");
 105   for(Int_t idet = 0; idet < kNDetectors; idet++){
 106     fPHSdetector->Add(new TProfile(Form("hPHSdet%d", idet), Form("Average Pulse Height in Detector %d", idet), 31, -0.5, 30.5));
 107     fQCLdetector->Add(new TH1F(Form("hQd%d",idet),Form("Qd%d",idet), 100,0,5000));
 108     fQTdetector->Add(new TH1F(Form("hQTd%d", idet), Form("Qtotal%d", idet), 6000, 0, 6000));
 109   }
 110   fContainer->Add(fPHSdetector);
 111   fContainer->Add(fQCLdetector);
 112   fContainer->Add(fQTdetector);
 113
 114   // Initialise the PHS, cluster charge and total charge deposit histograms for each sector
 115   fPHSsector = new TObjArray();
 116   fPHSsector->SetName("fPHSsector");
 117   fQCLsector = new TObjArray();
 118   fQCLsector->SetName("fQCLsector");
 119   fQTsector = new TObjArray();
 120   fQTsector->SetName("fQTsector");
 121   for(Int_t isec = 0; isec < kNSectors; isec++){
 122     fPHSsector->Add(new TProfile(Form("hPHSsec%d", isec), Form("Average Pulse Height in Sector %d", isec), 31, -0.5, 30.5));
 123     fQCLsector->Add(new TH1F(Form("hQd%d",isec),Form("Qd%d",isec), 100,0,5000));
 124     fQTsector->Add(new TH1F(Form("hQTd%d", isec), Form("Qtotal%d", isec), 6000, 0, 6000));
 125   }
 126   fContainer->Add(fPHSsector);
 127   fContainer->Add(fQCLsector);
 128   fContainer->Add(fQTsector);
 129   fTriggerNames = new TMap();
 130   fTriggerNames->SetName("fTriggerNames");
 131   fContainer->Add(fTriggerNames);
 132 }
 133
 134 //_______________________________________________________
 135 void AliTRDcheckDetector::Exec(Option_t *){
 136   //
 137   // Execution function
 138   // Filling TRD quality histos
 139   //
 140   if(fEventInfo->GetEventHeader()->GetEventType() != 7) return;
 141   Int_t nTracks = 0;            // Count the number of tracks per event
 142   Int_t triggermask = fEventInfo->GetEventHeader()->GetTriggerMask();
 143   TString triggername =  fEventInfo->GetRunInfo()->GetTriggerClass(triggermask);
 144   dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTrigger))->Fill(triggermask);
 145   AliTRDtrackInfo *fTrackInfo = 0x0;
 146   AliTRDtrackV1 *fTRDtrack = 0x0;
 147   AliTRDseedV1 *fTracklet = 0x0;
 148   AliTRDcluster *fTRDcluster = 0x0;
 149   for(Int_t iti = 0; iti < fTracks->GetEntriesFast(); iti++){
 150     fTrackInfo = dynamic_cast<AliTRDtrackInfo *>(fTracks->UncheckedAt(iti));
 151     if(!fTrackInfo || !(fTRDtrack = fTrackInfo->GetTRDtrack())) continue;
 152     Int_t nclusters = fTRDtrack->GetNumberOfClusters();
 153     Int_t ntracklets = fTRDtrack->GetNumberOfTracklets();
 154     Float_t chi2 = fTRDtrack->GetChi2();
 155     // Fill Histograms
 156     dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNclustersHist))->Fill(nclusters);
 157     dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNtrackletsHist))->Fill(ntracklets);
 158     dynamic_cast<TH1F *>(fContainer->UncheckedAt(kChi2))->Fill(chi2);
 159     dynamic_cast<TH1F *>(fContainer->UncheckedAt(kChi2Normalized))->Fill(chi2/static_cast<Float_t>(ntracklets));
 160     // now loop over single tracklets
 161     if(ntracklets > 2){
 162       Int_t sector = -1;
 163       for(Int_t ilayer = 0; ilayer < kNLayers; ilayer++){
 164         if(!(fTracklet = fTRDtrack->GetTracklet(ilayer)) || !fTracklet->IsOK()) continue;
 165         Float_t nClustersTracklet = fTracklet->GetN();
 166         dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNclusterTrackletHist))->Fill(nClustersTracklet);
 167         if(nClustersTracklet){
 168           Float_t Qtot = 0;
 169           Int_t detector = -1;
 170           for(Int_t itb = 0; itb < kNTimebins; itb++){
 171             if(!(fTRDcluster = fTracklet->GetClusters(itb))) continue;
 172             Int_t localtime        = fTRDcluster->GetLocalTimeBin();
 173             Double_t clusterCharge = fTRDcluster->GetQ();
 174             detector               = fTRDcluster->GetDetector();
 175             sector                 = static_cast<Int_t>(detector/kNDetectorsSector);
 176             Double_t absolute_charge = TMath::Abs(clusterCharge);
 177             Qtot += absolute_charge;
 178             // Fill Histograms
 179             // 1st. PHS
 180             dynamic_cast<TProfile *>(fContainer->UncheckedAt(kPulseHeight))->Fill(localtime, absolute_charge);
 181             dynamic_cast<TProfile *>(fPHSdetector->UncheckedAt(detector))->Fill(localtime, absolute_charge);
 182             dynamic_cast<TProfile *>(fPHSsector->UncheckedAt(sector))->Fill(localtime, absolute_charge);
 183             // 2nd. cluster charge
 184             dynamic_cast<TH1F *>(fContainer->UncheckedAt(kClusterCharge))->Fill(absolute_charge);
 185             dynamic_cast<TH1F *>(fQCLdetector->UncheckedAt(detector))->Fill(absolute_charge);
 186             dynamic_cast<TH1F *>(fQCLsector->UncheckedAt(sector))->Fill(absolute_charge);
 187           }
 188           // Fill total cluster charge histograms
 189           dynamic_cast<TH1F *>(fContainer->UncheckedAt(kChargeDeposit))->Fill(Qtot);
 190           dynamic_cast<TH1F *>(fQTdetector->UncheckedAt(detector))->Fill(Qtot);
 191           dynamic_cast<TH1F *>(fQTsector->UncheckedAt(sector))->Fill(Qtot);
 192         }
 193       }
 194       dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNTracksSectorHist))->Fill(sector);
 195     }
 196     nTracks++;
 197   }
 198   if(nTracks) dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTriggerTracks))->Fill(triggermask);
 199   dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNTracksEventHist))->Fill(nTracks);
 200   if(!fTriggerNames->FindObject(Form("%d", triggermask)))
 201     fTriggerNames->Add(new TObjString(Form("%d", triggermask)), new TObjString(triggername));
 202   PostData(0, fContainer);
 203 }
 204
 205 //_______________________________________________________
 206 void AliTRDcheckDetector::Terminate(Option_t *){
 207   //
 208   // Terminate function
 209   //
 210 }
 211
 212 //_______________________________________________________
 213 Bool_t AliTRDcheckDetector::PostProcess(){
 214         //
 215         // Do Postprocessing (for the moment set the number of Reference histograms)
 216         //
 217
 218         TH1 * histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNTracksEventHist));
 219         histo->GetXaxis()->SetTitle("Number of Tracks");
 220         histo->GetYaxis()->SetTitle("Events");
 221         // Calculate percentage of events containing tracks
 222         Int_t entries = (Int_t)histo->GetEntries();
 223   Float_t percTracks = entries > 0 ? (entries - histo->GetBinContent(1))/Float_t(entries) : 0.;
 224         printf("Percentage of Events containing TRD Tracks: %6.3f %%\n", 1.E2*percTracks);
 225         histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNclustersHist));
 226         histo->GetXaxis()->SetTitle("Number of Clusters");
 227         histo->GetYaxis()->SetTitle("Entries");
 228         histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNtrackletsHist));
 229         histo->GetXaxis()->SetTitle("Number of Tracklets");
 230         histo->GetYaxis()->SetTitle("Entries");
 231         histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNclusterTrackletHist));
 232         histo->GetXaxis()->SetTitle("Number of Clusters");
 233         histo->GetYaxis()->SetTitle("Entries");
 234         histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kChi2));
 235         histo->GetXaxis()->SetTitle("#chi^2");
 236         histo->GetYaxis()->SetTitle("Entries");
 237         histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNTracksSectorHist));
 238         histo->GetXaxis()->SetTitle("Sector");
 239         histo->GetYaxis()->SetTitle("Number of Tracks");
 240         histo = dynamic_cast<TProfile *>(fContainer->UncheckedAt(kPulseHeight));
 241         histo->GetXaxis()->SetTitle("Time / 100ns");
 242         histo->GetYaxis()->SetTitle("Average Pulse Height (a. u.)");
 243         histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kClusterCharge));
 244         histo->GetXaxis()->SetTitle("Cluster Charge (a.u.)");
 245         histo->GetYaxis()->SetTitle("Entries");
 246         histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kChargeDeposit));
 247         histo->GetXaxis()->SetTitle("Charge Deposit (a.u.)");
 248         histo->GetYaxis()->SetTitle("Entries");
 249         // Calculate the Percentage of events containing tracks as function of the trigger cluster
 250         histo = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTrigger));
 251         TH1F *histoTracks = dynamic_cast<TH1F *>(fContainer->UncheckedAt(kNEventsTriggerTracks));
 252         histoTracks->Divide(histo);
 253         fTriggerNames = dynamic_cast<TMap *>(fContainer->UncheckedAt(17));
 254         TH1F *percentages = new TH1F("fTriggerInf", "Trigger Information", fTriggerNames->GetEntries(), 0, fTriggerNames->GetEntries());
 255         TObject *triggerclass = 0x0;
 256         Int_t ipt=0;
 257         for(Int_t ibin = 0; ibin < histo->GetNbinsX(); ibin++){
 258                 if((triggerclass = fTriggerNames->FindObject(Form("%d",ibin)))){
 259                         TObjString *label = dynamic_cast<TObjString *>(dynamic_cast<TPair *>(triggerclass)->Value());
 260                         //printf("Trigger Pattern for class %d: %s\n", ibin, label->String().Data());
 261                         percentages->SetBinContent(percentages->FindBin(ipt), histoTracks->GetBinContent(histoTracks->FindBin(ibin)));
 262                         percentages->GetXaxis()->SetBinLabel(percentages->FindBin(ipt),label->String().Data());
 263                         ipt++;
 264                 }
 265         }
 266         percentages->GetXaxis()->SetTitle("Trigger Cluster");
 267         percentages->GetYaxis()->SetTitle("%Events");
 268         percentages->GetYaxis()->SetRangeUser(0,1);
 269 //      percentages->SetMarkerColor(kBlue);
 270 //      percentages->SetMarkerStyle(22);
 271         fContainer->Add(percentages);
 272         fNRefFigures = 10;
 273         return kTRUE;
 274 }
 275
 276 //_______________________________________________________
 277 void AliTRDcheckDetector::GetRefFigure(Int_t ifig, Int_t &first, Int_t &last, Option_t *opt){
 278         //
 279         // Setting Reference Figures
 280         //
 281         opt = "pl";
 282         switch(ifig){
 283                 case 0: first = last = kNTracksEventHist;
 284                                                 break;
 285                 case 1: first = last = kNclustersHist;
 286                                                 break;
 287                 case 2: first = last = kNtrackletsHist;
 288                                                 break;
 289                 case 3: first = last = kNclusterTrackletHist;
 290                                                 break;
 291                 case 4: first = last = kChi2;
 292                                                 break;
 293                 case 5: first = last = kNTracksSectorHist;
 294                                                 break;
 295                 case 6: first = last = kPulseHeight;
 296                                                 break;
 297                 case 7: first = last = kClusterCharge;
 298                                                 break;
 299                 case 8: first = last = kChargeDeposit;
 300                                                 break;
 301                 case 9: first = last = kPostProcessing;
 302                                                 opt="bar";
 303                                                 break;
 304                 default: first = last = kNTracksEventHist;
 305                                                 break;
 306         };
 307 }
 308