HLT/TPCLib/AliHLTTPCdEdxMonitoringComponent.cxx

   1 // $Id$
   2
   3 //**************************************************************************
   4 //* This file is property of and copyright by the ALICE HLT Project        *
   5 //* ALICE Experiment at CERN, All rights reserved.                         *
   6 //*                                                                        *
   7 //* Primary Authors: Per-Ivar Lønne <ploenne@cern.ch>                      *
   8 //*                  for The ALICE HLT Project.                            *
   9 //*                                                                        *
  10 //* Permission to use, copy, modify and distribute this software and its   *
  11 //* documentation strictly for non-commercial purposes is hereby granted   *
  12 //* without fee, provided that the above copyright notice appears in all   *
  13 //* copies and that both the copyright notice and this permission notice   *
  14 //* appear in the supporting documentation. The authors make no claims     *
  15 //* about the suitability of this software for any purpose. It is          *
  16 //* provided "as is" without express or implied warranty.                  *
  17 //**************************************************************************/
  18
  19 /// @file   AliHLTTPCdEdxMonitoringComponent.cxx
  20 /// @author Per-Ivar Lønne, Jochen Thaeder, Matthias Richter, Alexander Kalweit
  21 /// @date   21.08.2011
  22 /// @brief  Component for reading ESD from chain and produce a dEdx monitoring plot
  23 ///
  24
  25 #include "TSystem.h"
  26 #include "AliESDtrackCuts.h"
  27 #include <AliHLTDAQ.h>
  28 #include "AliESDEvent.h"
  29 #include "AliESDtrack.h"
  30 #include "AliLog.h"
  31 #include "TMap.h"
  32 #include "TObjString.h"
  33 #include "AliHLTDataTypes.h"
  34 #include "AliHLTComponentBenchmark.h"
  35 #include "TROOT.h"
  36 #include "TStyle.h"
  37 #include "TColor.h"
  38 #include "TH2F.h"
  39
  40 #include "AliHLTTPCdEdxMonitoringComponent.h"
  41
  42 using namespace std;
  43
  44 /** ROOT macro for the implementation of ROOT specific class methods */
  45 ClassImp( AliHLTTPCdEdxMonitoringComponent )
  46
  47 AliHLTTPCdEdxMonitoringComponent::AliHLTTPCdEdxMonitoringComponent()
  48   :
  49   AliHLTProcessor(),
  50   fESDTrackCuts(NULL),
  51   fHist(),
  52   fxbins(),
  53   fxmin(),
  54   fxmax(),
  55   fybins(),
  56   fymin(),
  57   fymax()
  58 {
  59   //Constructor
  60 }
  61
  62 AliHLTTPCdEdxMonitoringComponent::~AliHLTTPCdEdxMonitoringComponent()
  63 {
  64   //Destructor
  65 }
  66
  67 // ######################################################################### //
  68 const char* AliHLTTPCdEdxMonitoringComponent::GetComponentID()
  69 {
  70   // get component id
  71   return "TPCdEdxMonitoring";
  72 }
  73
  74 // ######################################################################### //
  75 void AliHLTTPCdEdxMonitoringComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list)
  76 {
  77   // get list of input data types
  78   list.push_back(kAliHLTDataTypeESDObject|kAliHLTDataOriginAny);
  79 }
  80
  81 // ######################################################################### //
  82 AliHLTComponentDataType AliHLTTPCdEdxMonitoringComponent::GetOutputDataType()
  83 {
  84   // get the output data size
  85   return kAliHLTDataTypeHistogram|kAliHLTDataOriginHLT;
  86 }
  87
  88 // ######################################################################### //
  89 void AliHLTTPCdEdxMonitoringComponent::GetOutputDataSize( unsigned long& constBase, Double_t& inputMultiplier )
  90 {
  91   // get output size estimator
  92   constBase = 4096;
  93   inputMultiplier = 0;
  94 }
  95
  96 // ######################################################################### //
  97 void AliHLTTPCdEdxMonitoringComponent::GetOCDBObjectDescription( TMap* const targetMap)
  98 {
  99   // Get a list of OCDB object description.
 100   // The list of objects is provided in a TMap
 101   // - key: complete OCDB path, e.g. GRP/GRP/Data
 102   // - value: short description why the object is needed
 103   // Key and value objects created inside this class go into ownership of
 104   // target TMap.
 105
 106     if (!targetMap) return;
 107     targetMap->Add(new TObjString("HLT/ConfigTPC/TPCdEdxMonitoring"),
 108     new TObjString("configuration object"));
 109     targetMap->Add(new TObjString("GRP/GRP/Data"),
 110     new TObjString("GRP object"));
 111
 112 }
 113
 114 // ######################################################################### //
 115 AliHLTComponent* AliHLTTPCdEdxMonitoringComponent::Spawn()
 116 {
 117   // Spawn function, return new class instance
 118   return new AliHLTTPCdEdxMonitoringComponent;
 119 }
 120
 121 // ######################################################################### //
 122 Int_t AliHLTTPCdEdxMonitoringComponent::DoInit( Int_t argc, const char** argv )
 123 {
 124   // init the component
 125   Int_t iResult=0;
 126   fESDTrackCuts = new AliESDtrackCuts("AliESDtrackCuts","HLT");
 127
 128   fxbins=Int_t(300);
 129   fxmin=Double_t(-2);
 130   fxmax=Double_t(2);
 131   fybins=Int_t(500);
 132   fymin=Double_t(0);
 133   fymax=Double_t(500);
 134
 135 if (!fESDTrackCuts)
 136     {
 137       iResult=-ENOMEM;
 138     }
 139
 140   if (iResult<0)
 141     {
 142       if (fESDTrackCuts)
 143         delete fESDTrackCuts;
 144       fESDTrackCuts = NULL;
 145     }
 146
 147   if (iResult>=0)
 148     {
 149       SetDefaultConfiguration();
 150       TString cdbPath="HLT/ConfigTPC/";
 151       cdbPath+=GetComponentID();
 152       iResult=ConfigureFromCDBTObjString(cdbPath);
 153
 154
 155       if (iResult>=0)
 156         {
 157           iResult=ConfigureFromArgumentString(argc, argv);
 158         }
 159     }
 160
 161   if (iResult>=0) {
 162     HLTInfo("ESD track cuts : %s",fESDTrackCuts->GetTitle() );
 163   }
 164   fHist = new TH2F("hHLT", "HLT", fxbins, fxmin, fxmax, fybins, fymin, fymax);
 165   fHist->GetXaxis()->SetTitle("momentum/charge #frac{p}{z}  (GeV/c)");
 166   fHist->GetYaxis()->SetTitle("dE/dx in TPC (a.u.)");
 167   Plotstyle();
 168   return iResult;
 169 }
 170
 171
 172 // ######################################################################### //
 173 Int_t AliHLTTPCdEdxMonitoringComponent::DoDeinit()
 174 {
 175   // component cleanup, delete all instances of helper classes here
 176   if (fESDTrackCuts)
 177     delete fESDTrackCuts;
 178   fESDTrackCuts = NULL;
 179   if (fHist)
 180     delete fHist;
 181   fHist=NULL;
 182   return 0;
 183 }
 184
 185
 186 // ######################################################################### //
 187 Int_t AliHLTTPCdEdxMonitoringComponent::DoEvent(const AliHLTComponentEventData& /*evtData*/,
 188                                               AliHLTComponentTriggerData& /*trigData*/)
 189 {
 190   // event processing function
 191   float sig;
 192   Double_t ptot;
 193
 194   // check if this is a data event, there are a couple of special events
 195   // which should be ignored for normal processing
 196   if (!IsDataEvent()) return 0;
 197
 198   const TObject* obj = GetFirstInputObject(kAliHLTAllDataTypes, "AliESDEvent");
 199
 200   // input objects are not supposed to be changed by the component, so they
 201   // are defined const. However, the implementation of AliESDEvent does not
 202   // support this and we need the const_cast
 203   AliESDEvent* esd = dynamic_cast<AliESDEvent*>(const_cast<TObject*>(obj));
 204   if (esd != NULL) {
 205     esd->GetStdContent();
 206     for (Int_t i = 0; i < esd->GetNumberOfTracks(); i++) {
 207       AliESDtrack* track = esd->GetTrack(i);
 208       sig=track->GetTPCsignal();
 209       if(!fESDTrackCuts->AcceptTrack(track)) continue;
 210       if (!track->GetInnerParam()) continue;
 211       ptot = track->GetInnerParam()->GetP()*track->GetSign();
 212       fHist->Fill(ptot, sig);
 213     }
 214   }
 215   // publish the histogram
 216   PushBack(fHist, kAliHLTDataTypeHistogram | kAliHLTDataOriginHLT);
 217
 218   return 0;
 219 }
 220
 221 // ######################################################################### //
 222 Int_t AliHLTTPCdEdxMonitoringComponent::ScanConfigurationArgument(Int_t argc, const char** argv)
 223 {
 224   // Scan configuration arguments
 225   // Return the number of processed arguments
 226   //        -EPROTO if argument format error (e.g. number expected but not found)
 227   //
 228   // The AliHLTComponent base class implements a parsing loop for argument strings and
 229   // arrays of strings which is invoked by ConfigureFromArgumentString/ConfigureFromCDBTObjString
 230   // The component needs to implement ScanConfigurationArgument in order to decode the arguments.
 231
 232   if (argc<=0) return 0;
 233   Int_t ii =0;
 234   TString argument=argv[ii];
 235
 236   if (argument.IsNull()) return 0;
 237
 238   if( !fESDTrackCuts){
 239     HLTError("No ESD track cuts availible");
 240     return -ENOMEM;
 241   }
 242
 243
 244   //**********************************//
 245   //        Histogram Binning         //
 246   //**********************************//
 247
 248
 249   // -xbins
 250   if (argument.CompareTo("-xbins")==0)
 251     {
 252       if (++ii>=argc) return -EPROTO;
 253       argument=argv[ii];
 254
 255       fxbins = argument.Atoi();
 256
 257       return 2;
 258     }
 259
 260   // -xmin
 261   if (argument.CompareTo("-xmin")==0)
 262     {
 263       if (++ii>=argc) return -EPROTO;
 264       argument=argv[ii];
 265
 266       fxmin = argument.Atoi();
 267
 268       return 2;
 269     }
 270
 271   if (argument.CompareTo("-xmax")==0)
 272     {
 273       if (++ii>=argc) return -EPROTO;
 274       argument=argv[ii];
 275
 276       fxmax = argument.Atoi();
 277
 278       return 2;
 279     }
 280
 281   // -xbins
 282   if (argument.CompareTo("-ybins")==0)
 283     {
 284       if (++ii>=argc) return -EPROTO;
 285       argument=argv[ii];
 286
 287       fybins = argument.Atoi();
 288
 289       return 2;
 290     }
 291
 292   // -xmin
 293   if (argument.CompareTo("-ymin")==0)
 294     {
 295       if (++ii>=argc) return -EPROTO;
 296       argument=argv[ii];
 297
 298       fymin = argument.Atoi();
 299
 300       return 2;
 301     }
 302
 303   if (argument.CompareTo("-ymax")==0)
 304     {
 305       if (++ii>=argc) return -EPROTO;
 306       argument=argv[ii];
 307
 308       fymax = argument.Atoi();
 309
 310       return 2;
 311     }
 312
 313
 314   //**********************************//
 315   //           Track Cuts             //
 316   //**********************************//
 317
 318  // -maxpt
 319   if (argument.CompareTo("-maxpt")==0) {
 320     if (++ii>=argc) return -EPROTO;
 321     argument=argv[ii];
 322
 323     Float_t minPt, maxPt;
 324     fESDTrackCuts->GetPtRange(minPt,maxPt);
 325     maxPt = argument.Atof();
 326     fESDTrackCuts->SetPtRange(minPt,maxPt);
 327
 328     TString title = fESDTrackCuts->GetTitle();
 329     if (!title.CompareTo("No track cuts")) title = "";
 330     else title += " && ";
 331     title += Form("p_t < %f", maxPt);
 332     fESDTrackCuts->SetTitle(title);
 333     return 2;
 334   }
 335
 336   // -minpt
 337   if (argument.CompareTo("-minpt")==0) {
 338     if (++ii>=argc) return -EPROTO;
 339     argument=argv[ii];
 340
 341     Float_t minPt, maxPt;
 342     fESDTrackCuts->GetPtRange(minPt,maxPt);
 343     minPt = argument.Atof();
 344     fESDTrackCuts->SetPtRange(minPt,maxPt);
 345
 346     TString title = fESDTrackCuts->GetTitle();
 347     if (!title.CompareTo("No track cuts")) title = "";
 348     else title += " && ";
 349     title += Form("p_t > %f", minPt);
 350     fESDTrackCuts->SetTitle(title);
 351     return 2;
 352   }
 353
 354   // -min-ldca
 355   // minimum longitudinal dca to vertex
 356   if (argument.CompareTo("-min-ldca")==0) {
 357     if (++ii>=argc) return -EPROTO;
 358     argument=argv[ii];
 359
 360     fESDTrackCuts->SetMinDCAToVertexZ(argument.Atof());
 361     TString title = fESDTrackCuts->GetTitle();
 362     if (!title.CompareTo("No track cuts")) title = "";
 363     else title += " && ";
 364     title += Form("DCAz > %f", argument.Atof());
 365     fESDTrackCuts->SetTitle(title);
 366     return 2;
 367   }
 368
 369   // -max-ldca
 370   // maximum longitudinal dca to vertex
 371   if (argument.CompareTo("-max-ldca")==0) {
 372     if (++ii>=argc) return -EPROTO;
 373     argument=argv[ii];
 374
 375     fESDTrackCuts->SetMaxDCAToVertexZ(argument.Atof());
 376     TString title = fESDTrackCuts->GetTitle();
 377     if (!title.CompareTo("No track cuts")) title = "";
 378     else title += " && ";
 379     title += Form("DCAz < %f", argument.Atof());
 380     fESDTrackCuts->SetTitle(title);
 381     return 2;
 382   }
 383
 384   // -min-tdca
 385   // minimum transverse dca to vertex
 386   if (argument.CompareTo("-min-tdca")==0) {
 387     if (++ii>=argc) return -EPROTO;
 388     argument=argv[ii];
 389
 390     fESDTrackCuts->SetMinDCAToVertexXY(argument.Atof());
 391     TString title = fESDTrackCuts->GetTitle();
 392     if (!title.CompareTo("No track cuts")) title = "";
 393     else title += " && ";
 394     title += Form("DCAr > %f", argument.Atof());
 395     fESDTrackCuts->SetTitle(title);
 396     return 2;
 397   }
 398
 399   // -max-tdca
 400   // maximum transverse dca to vertex
 401   if (argument.CompareTo("-max-tdca")==0) {
 402     if (++ii>=argc) return -EPROTO;
 403     argument=argv[ii];
 404
 405     fESDTrackCuts->SetMaxDCAToVertexXY(argument.Atof());
 406     TString title = fESDTrackCuts->GetTitle();
 407     if (!title.CompareTo("No track cuts")) title = "";
 408     else title += " && ";
 409     title += Form("DCAr < %f", argument.Atof());
 410     fESDTrackCuts->SetTitle(title);
 411     return 2;
 412   }
 413
 414   // -etarange
 415   // +/- eta
 416   if (argument.CompareTo("-etarange")==0) {
 417     if (++ii>=argc) return -EPROTO;
 418     argument=argv[ii];
 419     Float_t eta = argument.Atof();
 420
 421     fESDTrackCuts->SetEtaRange(-eta,eta);
 422     TString title = fESDTrackCuts->GetTitle();
 423     if (!title.CompareTo("No track cuts")) title = "";
 424     else title += " && ";
 425     title += Form("Eta[%f,%f]", argument.Atof(),argument.Atof());
 426     fESDTrackCuts->SetTitle(title);
 427     return 2;
 428   }
 429
 430   // -minNClsTPC
 431   // minimum clusters in TPC
 432   if (argument.CompareTo("-minNClsTPC")==0) {
 433     if (++ii>=argc) return -EPROTO;
 434     argument=argv[ii];
 435
 436     Int_t ncls;
 437     ncls = Int_t(argument.Atof());
 438     fESDTrackCuts->SetMinNClustersTPC(ncls);
 439
 440     TString title = fESDTrackCuts->GetTitle();
 441     if (!title.CompareTo("No track cuts")) title = "";
 442     else title += " && ";
 443     title += Form("minNClsTPC < %i", ncls);
 444     fESDTrackCuts->SetTitle(title);
 445     return 2;
 446   }
 447
 448
 449
 450   // unknown argument
 451   return -EINVAL;
 452
 453   return 0;
 454
 455 }
 456
 457 // ######################################################################### //
 458 Int_t AliHLTTPCdEdxMonitoringComponent::Reconfigure(const char* cdbEntry, const char* chainId)
 459 {
 460   // reconfigure the component from the specified CDB entry, or default CDB entry
 461   HLTInfo("reconfigure '%s' from entry %s", chainId, cdbEntry);
 462
 463   Int_t iResult=0;
 464   TString cdbPath;
 465   if (cdbEntry)
 466     {
 467       cdbPath=cdbEntry;
 468     } else {
 469       cdbPath="HLT/ConfigTPC/";
 470       cdbPath+=GetComponentID();
 471     }
 472
 473   iResult=ConfigureFromCDBTObjString(cdbPath); //// Or use return 0, and skip this line?
 474
 475   return iResult;
 476 }
 477
 478 // ######################################################################### //
 479 Int_t AliHLTTPCdEdxMonitoringComponent::ReadPreprocessorValues(const char* modules)
 480 {
 481   // read the preprocessor values for the detectors in the modules list
 482   Int_t iResult=0;
 483   TString detectors(modules!=NULL?modules:"");
 484   HLTInfo("read preprocessor values for detector(s): %s", detectors.IsNull()?"none":detectors.Data());
 485   return iResult; //Done differently in AliHLTMultiplicityCorrelationsComponent...
 486 }
 487
 488
 489
 490 // ######################################################################### //
 491 void AliHLTTPCdEdxMonitoringComponent::SetDefaultConfiguration()
 492 {
 493   if (fESDTrackCuts)
 494     {
 495       //fESDTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2010(kTRUE);
 496       fESDTrackCuts->SetEtaRange(-0.8,+0.8);
 497       fESDTrackCuts->SetPtRange(0.15,1e10);
 498       // 2011-10-28 investigation by Per Ivar and Alexander
 499       // the following cuts are not needed
 500       /*
 501       fESDTrackCuts->SetMaxCovDiagonalElements(2, 2, 0.5, 0.5, 2);  // BEWARE STANDARD VALUES ARE: 2, 2, 0.5, 0.5, 2
 502       fESDTrackCuts->SetMaxNsigmaToVertex(3);
 503       fESDTrackCuts->SetRequireSigmaToVertex(kTRUE);
 504       */
 505       fESDTrackCuts->SetAcceptKinkDaughters(kFALSE);
 506       fESDTrackCuts->SetMinNClustersTPC(70);
 507       fESDTrackCuts->SetMaxChi2PerClusterTPC(4);
 508       fESDTrackCuts->SetMaxDCAToVertexXY(3);
 509       fESDTrackCuts->SetMaxDCAToVertexZ(3);
 510       fESDTrackCuts->SetRequireTPCRefit(kTRUE);
 511       //fESDTrackCuts->SetRequireITSRefit(kTRUE); //Kills HLT simulated reconstructions?
 512       // 2011-10-28 investigation by Per Ivar and Alexander
 513       // those cuts remove all data points because the filling is different in the
 514       // HLT, needs further investigation
 515       /*
 516       fESDTrackCuts->SetClusterRequirementITS(AliESDtrackCuts::kSPD, AliESDtrackCuts::kAny); //TEMPORARY <-> REMOVE
 517       fESDTrackCuts->SetMinNClustersITS(3);
 518       */
 519
 520     }
 521  fxbins=300;
 522  fxmin=-2;
 523  fxmax=2;
 524  fybins=500;
 525  fymin=0;
 526  fymax=500;
 527   return;
 528 }
 529
 530
 531
 532 // ######################################################################### //
 533 void AliHLTTPCdEdxMonitoringComponent::Plotstyle()
 534 {
 535   gROOT->SetStyle("Plain");
 536
 537   gStyle->SetCanvasBorderMode(0);
 538   gStyle->SetCanvasColor(10);
 539   gStyle->SetCanvasDefH(550);
 540   gStyle->SetCanvasDefW(575);
 541   gStyle->SetPadBorderMode(0);
 542   gStyle->SetPadColor(10);
 543   gStyle->SetPadTickX(1);
 544   gStyle->SetPadTickY(1);
 545   gStyle->SetStatColor(10);
 546   gStyle->SetFrameFillColor(10);
 547   gStyle->SetPalette(1,0);
 548
 549   //
 550
 551   //gStyle->SetStatX(0.7);
 552   //gStyle->SetStatW(0.2);
 553   //gStyle->SetLabelOffset(1.2);
 554   //gStyle->SetLabelFont(72);
 555   //gStyle->SetLabelSize(0.6);
 556   //gStyle->SetTitleOffset(1.2);
 557   gStyle->SetTitleFontSize(0.04);
 558
 559
 560   gStyle->SetOptStat(10);
 561   gStyle->SetLineWidth(2);
 562   gStyle->SetMarkerSize(1.0);
 563   gStyle->SetTextSize(0.04);
 564   gStyle->SetTitleSize(0.04,"xyz");
 565   gStyle->SetLabelSize(0.04,"xyz");
 566   gStyle->SetLabelOffset(0.02,"xyz");
 567   gStyle->SetLabelFont(42,"xyz");
 568   gStyle->SetTitleOffset(1.3,"x");
 569   gStyle->SetTitleOffset(1.6,"y");
 570   gStyle->SetTitleOffset(1.6,"z");
 571   gStyle->SetTitleFont(42,"xyz");
 572   gStyle->SetTitleColor(1,"xyz");
 573   //gStyle->SetPadTopMargin(0.1);
 574   //gStyle->SetPadRightMargin(0.1);
 575   gStyle->SetPadBottomMargin(0.2);
 576   gStyle->SetPadLeftMargin(0.2);
 577
 578
 579   const Int_t NCont=255;
 580   const Int_t NRGBs = 5;
 581   Double_t stops[NRGBs] = { 0.00, 0.34, 0.61, 0.84, 1.00 };
 582   Double_t red[NRGBs]   = { 0.00, 0.00, 0.87, 1.00, 0.51 };
 583   Double_t green[NRGBs] = { 0.00, 0.81, 1.00, 0.20, 0.00 };
 584   Double_t blue[NRGBs]  = { 0.51, 1.00, 0.12, 0.00, 0.00 };
 585   TColor::CreateGradientColorTable(NRGBs, stops, red, green, blue, NCont);
 586   gStyle->SetNumberContours(NCont);
 587 }