HLT/trigger/AliHLTTRDTriggerComponent.cxx

   1 // $Id$
   2 //**************************************************************************
   3 //* This file is property of and copyright by the ALICE HLT Project        *
   4 //* ALICE Experiment at CERN, All rights reserved.                         *
   5 //*                                                                        *
   6 //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *
   7 //*                  Jochen Thaeder <jochen@thaeder.de>                    *
   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   AliHLTTRDTriggerComponent.cxx
  20 /// @author Felix Rettig, Stefan Kirsch
  21 /// @date   2012-08-16
  22 /// @brief
  23
  24 #include <cstdlib>
  25 #include "TSystem.h"
  26 #include "TObjArray.h"
  27 #include "TObjString.h"
  28 #include "TFile.h"
  29 #include "TH1I.h"
  30 #include "TH2I.h"
  31 #include "AliESDtrack.h"
  32 #include "AliESDEvent.h"
  33 #include "AliHLTTRDDefinitions.h"
  34 #include "AliHLTTriggerDecision.h"
  35 #include "AliHLTDomainEntry.h"
  36 #include "AliHLTLogging.h"
  37 #include "AliHLTErrorGuard.h"
  38 #include "AliHLTGlobalBarrelTrack.h"
  39 #include "AliESDtrackCuts.h"
  40 //#include "AliHLTESDTrackCuts.h"
  41 #include "AliGeomManager.h"
  42 #include "AliHLTComponentBenchmark.h"
  43 #include "AliHLTTRDTriggerComponent.h"
  44 #include "AliTRDonlineTrackingDataContainer.h"
  45 #include "AliTRDpadPlane.h"
  46
  47 ClassImp(AliHLTTRDTriggerComponent)
  48
  49 #define LogError( ... ) { HLTError(__VA_ARGS__); if (fDebugLevel >= 1) { DbgLog("ERROR", __VA_ARGS__); } }
  50 #define LogInfo( ... ) { HLTInfo(__VA_ARGS__); if (fDebugLevel >= 1) { DbgLog("INFO", __VA_ARGS__); } }
  51 #define LogInspect( ... ) { HLTDebug(__VA_ARGS__); if (fDebugLevel >= 1) { DbgLog("INSPECT", __VA_ARGS__); } }
  52 #define LogDebug( ... ) { if (fDebugLevel >= 1) { HLTInfo(__VA_ARGS__); DbgLog("DEBUG", __VA_ARGS__); } }
  53
  54 AliHLTTRDTriggerComponent::AliHLTTRDTriggerComponent() :
  55   AliHLTTrigger(),
  56   fName(),
  57   fRefTrackSelectionEtaLimit(0.9),
  58   fRefTrackSelectionVertexXYLimit(20.),
  59   fRefTrackSelectionVertexZLimit(30.),
  60   fRefTrackSelectionPtThreshold(0.7),
  61   fMatchRatingThreshold(0.25),
  62   fElectronTriggerPtThresholdHSE(3.),
  63   fElectronTriggerPIDThresholdHSE(144),
  64   fElectronTriggerPtThresholdHQU(2.),
  65   fElectronTriggerPIDThresholdHQU(164),
  66   fApplyRefTrackCuts(kFALSE),
  67   fElectronTriggerOnL1TrgOnly(kFALSE),
  68   fHistoMode(1),
  69   fDebugLevel(0),
  70   fExtendedHistos(kFALSE),
  71   fEventRendering(kFALSE),
  72   fPushHistos(kFALSE),
  73   fWriteHistos(kFALSE),
  74   fEventId(fgkInvalidEventId),
  75   fRunNumber(-1),
  76   fChunkId(NULL),
  77   fSectorsWithData(0),
  78   fIsMinBiasEvent(kFALSE),
  79   fIsTRDElectronEvent(kFALSE),
  80   fESDtracksPresent(kFALSE),
  81   fHLTtracksPresent(kFALSE),
  82   fTRDGeometry(NULL),
  83   fEsdEvent(NULL),
  84   fTrackingData(NULL),
  85   fHLTTracks(NULL),
  86   fRefTrackCuts(NULL),
  87 #ifdef __TRDHLTDEBUG
  88   fEventDisplay(NULL),
  89   fBenchmark(NULL),
  90 #endif
  91   fHistArray(NULL),
  92   fHistMatchRating(NULL),
  93   fHistMatchRatingByPt(NULL),
  94   fHistMatchRatingByPid(NULL),
  95   fHistTrackPt(NULL),
  96   fHistTrackPtMatched(NULL),
  97   fHistTrackPtCorr(NULL),
  98   fHistTrackPid(NULL),
  99   fHistTrackPidMatched(NULL),
 100   fHistElectronCandidatePt(NULL),
 101   fHistElectronCandidateMatchedPt(NULL),
 102   fHistElectronCandidatePid(NULL),
 103   fHistElectronCandidateMatchedPid(NULL),
 104   fHistRefTrackPid(NULL),
 105   fHistMatchedRefTrackPid(NULL),
 106   fHistPIDvsTruncPID(NULL),
 107   fHistElectronFalsePIDvsTruncPID(NULL),
 108   fHistElectronConfirmedPIDvsTruncPID(NULL),
 109   fHistTrackMatchingCombinations(NULL),
 110   fHistElectronTriggerBaseMinBias(NULL),
 111   fHistElectronTriggerBaseTrdL1(NULL)
 112 {
 113 }
 114
 115 const char* AliHLTTRDTriggerComponent::fgkDefaultOCDBEntry = "HLT/ConfigHLT/TRDTrigger";
 116 const char* AliHLTTRDTriggerComponent::fgkTriggerDecisionElectronHSE = "TRD-ELECTRON-HSE";
 117 const char* AliHLTTRDTriggerComponent::fgkTriggerDecisionElectronHQU = "TRD-ELECTRON-HQU";
 118
 119 AliHLTTRDTriggerComponent::~AliHLTTRDTriggerComponent()
 120 {
 121 }
 122
 123 const char* AliHLTTRDTriggerComponent::GetTriggerName() const
 124 {
 125   if (!fName.IsNull())
 126     return fName.Data();
 127   else
 128     return "TRDTriggerComponent";
 129 }
 130
 131 AliHLTComponent* AliHLTTRDTriggerComponent::Spawn()
 132 {
 133   return new AliHLTTRDTriggerComponent;
 134 }
 135
 136 int AliHLTTRDTriggerComponent::DoInit(int argc, const char** argv)
 137 {
 138   int iResult = 0;
 139
 140   do {
 141
 142     fChunkId = new TString("XXXXX");
 143 #ifdef __TRDHLTDEBUG
 144     if (fChunkId){
 145       // chunk identification for debug output
 146       *fChunkId = gSystem->WorkingDirectory();
 147       fChunkId->Remove(0, fChunkId->Last('/') + 1);
 148       if (fChunkId->Contains("hlt_trd_comp"))
 149         *fChunkId = "L";
 150     } else {
 151       iResult = -ENOMEM;
 152       break;
 153     }
 154 #endif
 155
 156     AliGeomManager::LoadGeometry();
 157
 158     fTRDGeometry = new AliTRDgeometry();
 159     if (!fTRDGeometry) {
 160       iResult = -ENOMEM;
 161       break;
 162     }
 163
 164     fTrackingData = new AliTRDonlineTrackingDataContainer();
 165     if (!fTrackingData) {
 166       iResult = -ENOMEM;
 167       break;
 168     }
 169
 170     fHLTTracks = new vector<AliHLTGlobalBarrelTrack>;
 171     if (!fHLTTracks) {
 172       iResult = -ENOMEM;
 173       break;
 174     }
 175
 176     fRefTrackCuts = new AliESDtrackCuts("AliESDtrackCuts", "No track cuts");
 177     if (fRefTrackCuts){
 178       // fRefTrackCuts = AliESDtrackCuts::GetStandardITSTPCTrackCuts2011(kTRUE, 0);
 179       fRefTrackCuts = AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
 180       fRefTrackCuts->SetEtaRange(-0.8, 0.8);
 181     } else {
 182       iResult = -ENOMEM;
 183       break;
 184     }
 185
 186     fHistArray = new TObjArray(25);
 187     if(!fHistArray)
 188       return -ENOMEM;
 189     fHistArray->SetOwner(kTRUE);
 190
 191     fHistMatchRating = new TH1I("trdtrg_match_rating", "Match rating distribution HLT tracks vs. GTU track;match rating;abundance;",
 192                                 1./0.02, 0., 1.);
 193     fHistArray->AddLast(fHistMatchRating);
 194
 195     fHistMatchRatingByPt = new TH2I("trdtrg_match_rating_pt", "Match rating distribution HLT tracks vs. GTU track by p_{T};match rating m;p_{T}^{TRD,online} (GeV/c)",
 196                                     1./0.02, 0., 1., 20./0.02, 0., 20.);
 197     fHistArray->AddLast(fHistMatchRatingByPt);
 198
 199     fHistMatchRatingByPid = new TH2I("trdtrg_match_rating_pid", "Match rating distribution HLT tracks vs. GTU track by PID;match rating m;PID (a.u.)",
 200                                      1./0.02, 0., 1., 256, 0, 256);
 201     fHistArray->AddLast(fHistMatchRatingByPid);
 202
 203     fHistTrackPt = new TH1I("trdtrg_track_pt", "GTU track p_{T};GTU track p_{T} (GeV/c);abundance", 200, -20, 20.);
 204     fHistArray->AddLast(fHistTrackPt);
 205
 206     fHistTrackPtMatched = new TH1I("trdtrg_matched_track_pt", "matched GTU track p_{T};GTU track p_{T} (GeV/c);abundance", 200, -20, 20.);
 207     fHistArray->AddLast(fHistTrackPtMatched);
 208
 209     fHistTrackPtCorr = new TH2I("trdtrg_track_pt_corr", "HLT vs. GTU track p_{T};p_{T}^{HLT} (GeV/c);p_{T}^{GTU} (GeV/c)",
 210                400, -40., 40., 400, -40., 40.);
 211     fHistArray->AddLast(fHistTrackPtCorr);
 212
 213     fHistTrackPid = new TH1I("trdtrg_track_pid", "GTU track PID;GTU track PID (a.u.);abundance", 256, 0, 256);
 214     fHistArray->AddLast(fHistTrackPid);
 215
 216     fHistTrackPidMatched = new TH1I("trdtrg_matched_track_pid", "matched GTU track PID;GTU track PID (a.u.);abundance", 256, 0, 256);
 217     fHistArray->AddLast(fHistTrackPidMatched);
 218
 219     fHistTrackMatchingCombinations = new TH2I("trdtrg_matching_combinations", "HLT-GTU track matching combinations;set;combinations to check;",
 220                                               10, 0, 10, 2000, 0, 10000);
 221     fHistArray->AddLast(fHistTrackMatchingCombinations);
 222
 223     fHistElectronTriggerBaseMinBias = new TH1I("trdtrg_electron_trigger_base_mb", "min. bias base numbers for electron trigger analysis;set;abundance;",
 224                                               10, 0, 10);
 225     fHistElectronTriggerBaseMinBias->GetXaxis()->SetBinLabel(1, "min. bias events");
 226     fHistElectronTriggerBaseMinBias->GetXaxis()->SetBinLabel(2, "TRD L1 electron triggered");
 227     fHistElectronTriggerBaseMinBias->GetXaxis()->SetBinLabel(3, "TRD HLT electron triggered");
 228     fHistArray->AddLast(fHistElectronTriggerBaseMinBias);
 229
 230     fHistElectronTriggerBaseTrdL1 = new TH1I("trdtrg_electron_trigger_base_l1", "TRD L1 base numbers for electron trigger analysis;set;abundance;",
 231                                               10, 0, 10);
 232     fHistElectronTriggerBaseTrdL1->GetXaxis()->SetBinLabel(1, "TRD L1 electron triggered");
 233     fHistElectronTriggerBaseTrdL1->GetXaxis()->SetBinLabel(2, "TRD HLT electron triggered");
 234     fHistArray->AddLast(fHistElectronTriggerBaseTrdL1);
 235
 236     fHistElectronCandidatePt = new TH1I("trdtrg_electron_candidate_pt", "GTU electron candidate p_{T};GTU track p_{T} (GeV/c);abundance", 200, -20, 20.);
 237     fHistArray->AddLast(fHistElectronCandidatePt);
 238
 239     fHistElectronCandidateMatchedPt = new TH1I("trdtrg_electron_candidate_matched_pt", "matching GTU electron candidate p_{T};GTU track p_{T} (GeV/c);abundance", 200, -20, 20.);
 240     fHistArray->AddLast(fHistElectronCandidateMatchedPt);
 241
 242     fHistElectronCandidatePid = new TH1I("trdtrg_electron_candidate_pid", "GTU electron candidate PID;GTU track PID (a.u.);abundance", 256, 0, 256);
 243     fHistArray->AddLast(fHistElectronCandidatePid);
 244
 245     fHistElectronCandidateMatchedPid = new TH1I("trdtrg_electron_candidate_matched_pid", "matching GTU electron candidate PID;GTU track PID (a.u.);abundance", 256, 0, 256);
 246     fHistArray->AddLast(fHistElectronCandidateMatchedPid);
 247
 248   } while (0);
 249
 250   if (iResult < 0){
 251
 252     if (fHistArray) delete fHistArray;
 253     fHistArray = NULL;
 254
 255     if (fTRDGeometry) delete fTRDGeometry;
 256     fTRDGeometry = NULL;
 257
 258     if (fRefTrackCuts) delete fRefTrackCuts;
 259     fRefTrackCuts = NULL;
 260
 261     if (fHLTTracks) delete fHLTTracks;
 262     fHLTTracks = NULL;
 263
 264     if (fTrackingData) delete fTrackingData;
 265     fTrackingData = NULL;
 266
 267     if (fChunkId) delete fChunkId;
 268     fChunkId = NULL;
 269
 270   }
 271
 272   // check if the -triggername argument is used, the trigger name determines the following initialization
 273   vector<const char*> remainingArgs;
 274   for (int i = 0; i < argc; ++i) {
 275     if (strcmp(argv[i], "-triggername") == 0) {
 276       if (++i < argc){
 277         fName = argv[i];
 278       } else {
 279         LogError("invalid parameter for argument '-triggername', string expected");
 280         return -EINVAL;
 281       }
 282       continue;
 283     }
 284     remainingArgs.push_back(argv[i]);
 285   }
 286
 287   TString cdbPath;
 288   if (!fName.IsNull()) {
 289     cdbPath = "HLT/ConfigHLT/";
 290     cdbPath += fName;
 291   } else {
 292     cdbPath = fgkDefaultOCDBEntry;
 293   }
 294
 295   LogInfo("cdbPath = <%s>", cdbPath.Data());
 296   iResult = ConfigureFromCDBObject(cdbPath);
 297
 298   if (iResult >= 0 && argc > 0)
 299     iResult = ConfigureFromArgumentString(remainingArgs.size(), &(remainingArgs[0]));
 300
 301   return iResult;
 302 }
 303
 304 int AliHLTTRDTriggerComponent::DoDeinit()
 305 {
 306
 307 #ifdef __TRDHLTDEBUG
 308   if (fEventDisplay) delete fEventDisplay;
 309   fEventDisplay = NULL;
 310 #endif
 311
 312   if ((fHistoMode == 1) && (fWriteHistos)){
 313     TFile out("trg_out/trg_hists.root", "RECREATE");
 314     if (!out.IsZombie()) {
 315       out.cd();
 316       UInt_t numHists = fHistArray->GetEntries();
 317       for (UInt_t iHist = 0; iHist < numHists; ++iHist)
 318         if (fHistArray->At(iHist))
 319           fHistArray->At(iHist)->Write();
 320       out.Close();
 321     }
 322   }
 323
 324   if (fHistArray) delete fHistArray;
 325   fHistArray = NULL;
 326
 327   if (fRefTrackCuts) delete fRefTrackCuts;
 328   fRefTrackCuts = NULL;
 329
 330   if (!fHLTTracks) delete fHLTTracks;
 331   fHLTTracks = NULL;
 332
 333   if (fTrackingData) delete fTrackingData;
 334   fTrackingData = NULL;
 335
 336   if (fTRDGeometry) delete fTRDGeometry;
 337   fTRDGeometry = NULL;
 338
 339   if (fChunkId) delete fChunkId;
 340   fChunkId = NULL;
 341
 342   return 0;
 343 }
 344
 345 int AliHLTTRDTriggerComponent::Reconfigure(const char* cdbEntry, const char* /*chainId*/)
 346 {
 347   TString cdbPath;
 348   if (!cdbEntry || cdbEntry[0] == 0) {
 349     if (!fName.IsNull()) {
 350       cdbPath = "HLT/ConfigHLT/";
 351       cdbPath += fName;
 352     } else {
 353       cdbPath = fgkDefaultOCDBEntry;
 354     }
 355   } else
 356     cdbPath = cdbEntry;
 357
 358   return ConfigureFromCDBObject(cdbPath);
 359 }
 360
 361 int AliHLTTRDTriggerComponent::ReadPreprocessorValues(const char* /*modules*/)
 362 {
 363   return 0;
 364 }
 365
 366 Int_t AliHLTTRDTriggerComponent::ConfigureFromCDBObject(TString /*cdbPath*/)
 367 {
 368   return 0;
 369 }
 370
 371 int AliHLTTRDTriggerComponent::ScanConfigurationArgument(int argc, const char** argv)
 372 {
 373
 374   if (argc <= 0)
 375     return 0;
 376
 377   UShort_t iArg = 0;
 378   TString argument(argv[iArg]);
 379
 380   if (!argument.CompareTo("-match-sel-eta")){
 381     if (++iArg >= argc) return -EPROTO;
 382     argument = argv[iArg];
 383     fRefTrackSelectionEtaLimit = argument.Atof();
 384     LogInfo("ref track selection eta limit set to %.1f.", fRefTrackSelectionEtaLimit);
 385     return 2;
 386   }
 387
 388   if (!argument.CompareTo("-match-sel-vxy")){
 389     if (++iArg >= argc) return -EPROTO;
 390     argument = argv[iArg];
 391     fRefTrackSelectionVertexXYLimit = argument.Atof();
 392     LogInfo("ref track selection vertex xy limit set to %.1f.", fRefTrackSelectionVertexXYLimit);
 393     return 2;
 394   }
 395
 396   if (!argument.CompareTo("-match-sel-vz")){
 397     if (++iArg >= argc) return -EPROTO;
 398     argument = argv[iArg];
 399     fRefTrackSelectionVertexZLimit = argument.Atof();
 400     LogInfo("ref track selection vertex z limit set to %.1f.", fRefTrackSelectionVertexZLimit);
 401     return 2;
 402   }
 403
 404   if (!argument.CompareTo("-match-sel-pt")){
 405     if (++iArg >= argc) return -EPROTO;
 406     argument = argv[iArg];
 407     fRefTrackSelectionPtThreshold = argument.Atof();
 408     LogInfo("ref track selection pt threshold set to %.1f GeV/c.", fRefTrackSelectionPtThreshold);
 409     return 2;
 410   }
 411
 412   if (!argument.CompareTo("-match-rating")){
 413     if (++iArg >= argc) return -EPROTO;
 414     argument = argv[iArg];
 415     fMatchRatingThreshold = argument.Atof();
 416     LogInfo("match rating threshold set to %.2f GeV/c.", fMatchRatingThreshold);
 417     return 2;
 418   }
 419
 420   if (!argument.CompareTo("-histo-mode")){
 421     if (++iArg >= argc) return -EPROTO;
 422     argument = argv[iArg];
 423     fHistoMode = argument.Atoi();
 424     LogInfo("histo mode set to %d", fHistoMode);
 425     return 2;
 426   }
 427
 428   if (!argument.CompareTo("-trghse")){
 429     LogInfo("TRD HLT electron trigger HSE enabled.");
 430     return 1;
 431   }
 432
 433   if (!argument.CompareTo("-hse-pt")){
 434     if (++iArg >= argc) return -EPROTO;
 435     argument = argv[iArg];
 436     fElectronTriggerPtThresholdHSE = argument.Atof();
 437     LogInfo("pt threshold for HSE trigger set to %.1f GeV/c.", fElectronTriggerPtThresholdHSE);
 438     return 2;
 439   }
 440
 441   if (!argument.CompareTo("-hse-pid")){
 442     if (++iArg >= argc) return -EPROTO;
 443     argument = argv[iArg];
 444     fElectronTriggerPIDThresholdHSE = argument.Atof();
 445     LogInfo("PID threshold for HSE trigger set to %d.", fElectronTriggerPIDThresholdHSE);
 446     return 2;
 447   }
 448
 449   if (!argument.CompareTo("-trghqu")){
 450     LogInfo("TRD HLT electron trigger HQU enabled.");
 451     return 1;
 452   }
 453
 454   if (!argument.CompareTo("-hqu-pt")){
 455     if (++iArg >= argc) return -EPROTO;
 456     argument = argv[iArg];
 457     fElectronTriggerPtThresholdHQU = argument.Atof();
 458     LogInfo("pt threshold for HQU trigger set to %.1f GeV/c.", fElectronTriggerPtThresholdHQU);
 459     return 2;
 460   }
 461
 462   if (!argument.CompareTo("-hqu-pid")){
 463     if (++iArg >= argc) return -EPROTO;
 464     argument = argv[iArg];
 465     fElectronTriggerPIDThresholdHQU = argument.Atof();
 466     LogInfo("PID threshold for HQU trigger set to %.1f GeV/c.", fElectronTriggerPIDThresholdHQU);
 467     return 2;
 468   }
 469
 470   if (!argument.CompareTo("-l1-only")){
 471     LogInfo("evaluation of electron trigger only for events with TRD L1 electron trigger enabled.");
 472     fElectronTriggerOnL1TrgOnly = kTRUE;
 473     return 1;
 474   }
 475
 476   if (!argument.CompareTo("-histo-ext")){
 477     LogInfo("extended histogramming enabled.");
 478     fExtendedHistos = kTRUE; // enable extended histogramming, for debugging/tuning only!
 479
 480     if (!fHistRefTrackPid){
 481       fHistRefTrackPid = new TH2I("trdtrg_reftrack_pid", "TPC dE/dx by p_{T} for matching GTU tracks;p_{T} (GeV/c);dE/dx (a. u.)",
 482                                   500, 0., 10., 500, 0., 500);
 483       fHistArray->AddLast(fHistRefTrackPid);
 484     }
 485
 486     if (!fHistMatchedRefTrackPid){
 487       fHistMatchedRefTrackPid = new TH2I("trdtrg_matched_reftrack_pid", "TPC dE/dx by p_{T} for matching GTU tracks;p_{T} (GeV/c);dE/dx (a. u.)",
 488                                          500, 0., 10., 500, 0., 500);
 489       fHistArray->AddLast(fHistMatchedRefTrackPid);
 490     }
 491
 492     if (!fHistPIDvsTruncPID){
 493       fHistPIDvsTruncPID = new TH2I("trdtrg_pid_trunc_pid", "GTU track PID vs. truncated PID;PID (a.u.);truncated PID (a. u.)",
 494                                     256, 0., 256., 256, 0., 256.);
 495       fHistArray->AddLast(fHistPIDvsTruncPID);
 496     }
 497
 498     if (!fHistElectronFalsePIDvsTruncPID){
 499       fHistElectronFalsePIDvsTruncPID = new TH2I("trdtrg_electron_false_pid_trunc_pid", "false electron PID vs. truncated PID;PID (a.u.);truncated PID (a. u.)",
 500                                                  256, 0., 256., 256, 0., 256.);
 501       fHistArray->AddLast(fHistElectronFalsePIDvsTruncPID);
 502     }
 503
 504     if (!fHistElectronConfirmedPIDvsTruncPID){
 505       fHistElectronConfirmedPIDvsTruncPID = new TH2I("trdtrg_electron_confirmed_pid_trunc_pid", "confirmed electron PID vs. truncated PID;PID (a.u.);truncated PID (a. u.)",
 506                                                      256, 0., 256., 256, 0., 256.);
 507       fHistArray->AddLast(fHistElectronConfirmedPIDvsTruncPID);
 508     }
 509
 510     return 1;
 511   }
 512
 513   if (!argument.CompareTo("-ref-cuts")){
 514     LogInfo("ref track cuts for matching enabled.");
 515     fApplyRefTrackCuts = kTRUE;
 516     return 1;
 517   }
 518
 519   if (!argument.CompareTo("-push-histograms")){
 520     LogInfo("pushing of histograms enabled.");
 521     fPushHistos = kTRUE; // enable histogram pushing, for debugging/tuning only!
 522     return 1;
 523   }
 524
 525   if (!argument.CompareTo("-write-histograms")){
 526     LogInfo("writing of histograms enabled.");
 527     fWriteHistos = kTRUE; // enable histogram writing, for debugging/tuning only!
 528     return 1;
 529   }
 530
 531   if (!argument.CompareTo("-render")){
 532 #ifdef __TRDHLTDEBUG
 533     LogInfo("rendering of interesting events enabled.");
 534     if (!fEventDisplay)
 535       fEventDisplay = new AliTRDtrackingEventDisplay();
 536       LogInfo("event rendering activated. this is for debugging only!");
 537     fEventRendering = kTRUE; // enable event rendering, for debugging/tuning only!
 538 #endif
 539     return 1;
 540   }
 541
 542   if (!argument.CompareTo("-debug")){
 543     if (++iArg >= argc) return -EPROTO;
 544     argument = argv[iArg];
 545     fDebugLevel = argument.Atoi();
 546     LogInfo("debug level set to %d.", fDebugLevel);
 547     return 2;
 548   }
 549
 550   if (!argument.CompareTo("-ref-pt")){
 551     if (++iArg >= argc) return -EPROTO;
 552     argument = argv[iArg];
 553     Float_t minPt, maxPt;
 554     fRefTrackCuts->GetPtRange(minPt, maxPt);
 555     maxPt = argument.Atof();
 556     fRefTrackCuts->SetPtRange(minPt, maxPt);
 557     LogInfo("ref track pt range set to %.1f .. %.1f GeV/c.", minPt, maxPt);
 558     return 2;
 559   }
 560
 561   if (!argument.CompareTo("-ref-tdca")){
 562     if (++iArg >= argc) return -EPROTO;
 563     argument = argv[iArg];
 564     fRefTrackCuts->SetMaxDCAToVertexXY(argument.Atof());
 565     LogInfo("ref track DCA transverse threshold set to %.3f", argument.Atof());
 566     return 2;
 567   }
 568
 569   if (!argument.CompareTo("-ref-ldca")){
 570     if (++iArg >= argc) return -EPROTO;
 571     argument = argv[iArg];
 572     fRefTrackCuts->SetMaxDCAToVertexZ(argument.Atof());
 573     LogInfo("ref track longitudinal DCA threshold set to %.3f", argument.Atof());
 574     return 2;
 575   }
 576
 577   return 1;
 578 }
 579
 580 void AliHLTTRDTriggerComponent::ScanTriggerClasses(const char* firedTriggerClasses) {
 581
 582   fIsMinBiasEvent = kFALSE;
 583
 584   TString trg(firedTriggerClasses);
 585   if (trg.Index("CINT7WU-S-NOPF-ALL") >= 0)
 586     fIsMinBiasEvent = kTRUE;
 587
 588   if (trg.Index("CINT8WU-S-NOPF-ALL") >= 0)
 589     fIsMinBiasEvent = kTRUE;
 590
 591 }
 592
 593 Bool_t AliHLTTRDTriggerComponent::CheckRefTrackCuts(AliESDtrack* track){
 594
 595   // cuts by cut class
 596   if (fApplyRefTrackCuts)
 597     return fRefTrackCuts->AcceptTrack(track) ? kTRUE : kFALSE;
 598
 599   // simple custom cuts
 600   Float_t dcaToVertexXY, dcaToVertexZ;
 601   track->GetImpactParametersTPC(dcaToVertexXY, dcaToVertexZ);
 602   // LogDebug("IMPACT TPC %.4f  %.4f", dcaToVertexXY, dcaToVertexZ);
 603
 604   if ((dcaToVertexXY > fRefTrackSelectionVertexXYLimit) || (dcaToVertexZ >= fRefTrackSelectionVertexZLimit))
 605     return kFALSE;
 606
 607   if (TMath::Abs(track->Eta()) > fRefTrackSelectionEtaLimit)
 608     return kFALSE;
 609
 610   return kTRUE;
 611 }
 612
 613
 614 void AliHLTTRDTriggerComponent::DbgLog(const char* prefix, ...){
 615 #ifdef __TRDHLTDEBUG
 616   AliHLTEventID_t eventNumber = fEventId;
 617   printf("TRDHLTGM %s-%s-%s: [TRG] %s",
 618          (fRunNumber >= 0) ? Form("%06d", fRunNumber) : "XXXXXX",
 619          fChunkId->Data(),
 620          (eventNumber != fgkInvalidEventId) ? Form("%05llu", eventNumber) : "XXXXX",
 621          (strlen(prefix) > 0) ? Form("<%s> ", prefix) : "");
 622 #endif
 623   va_list args;
 624   va_start(args, prefix);
 625   char* fmt = va_arg(args, char*);
 626   vprintf(fmt, args);
 627   printf("\n");
 628   va_end(args);
 629 }
 630
 631 int AliHLTTRDTriggerComponent::PrepareESDData(){
 632
 633   int result = 0;
 634   fESDtracksPresent = kFALSE;
 635
 636   // check access to ESD event data
 637   const TObject* obj = GetFirstInputObject(kAliHLTAllDataTypes, "AliESDEvent");
 638   fEsdEvent = dynamic_cast<AliESDEvent*>(const_cast<TObject*>(obj));
 639
 640   if (fEsdEvent) {
 641     fEsdEvent->GetStdContent();
 642     fRunNumber = fEsdEvent->GetRunNumber();
 643     fESDtracksPresent = kTRUE;
 644     unsigned int numTracks = fEsdEvent->GetNumberOfTracks();
 645
 646     // process trigger classes
 647     ScanTriggerClasses(fEsdEvent->GetFiredTriggerClasses().Data());
 648
 649     LogDebug("ESD event data: %d ESD tracks [event num: %d] [minbias: %d (fired:%s)]",
 650                     numTracks, fEsdEvent->GetEventNumberInFile(), fIsMinBiasEvent, fEsdEvent->GetFiredTriggerClasses().Data());
 651
 652     // process ESD tracks
 653     if ((fDebugLevel >= 2) || fExtendedHistos){
 654       AliESDtrack* esdTrack;
 655       TString paramStr("");
 656       for (unsigned int iTrack = 0; iTrack < numTracks; ++iTrack) {
 657         esdTrack = fEsdEvent->GetTrack(iTrack);
 658
 659         if (fExtendedHistos){
 660           fHistRefTrackPid->Fill(esdTrack->Pt(), esdTrack->GetTPCsignal());
 661         }
 662
 663         if (fDebugLevel >= 2){
 664           paramStr = "";
 665           if (esdTrack){
 666
 667             if (esdTrack->GetInnerParam())
 668               paramStr += "I";
 669
 670             if (esdTrack->GetOuterParam())
 671               paramStr += "O";
 672
 673             LogDebug("ESD track %4d - pt: %+.2fGeV/c  [params: S%s]", iTrack, esdTrack->GetSignedPt(), paramStr.Data());
 674           } else
 675             LogError("ESD data for track %d invalid", iTrack);
 676         }
 677
 678       } // loop over ESD tracks
 679     }
 680
 681     fTrackingData->SetGtuPtMultiplierFromMagField(fEsdEvent->GetMagneticField()); // used for sign correction
 682
 683     result = 1;
 684
 685   }
 686
 687   return result;
 688
 689 }
 690
 691 int AliHLTTRDTriggerComponent::PrepareHLTData(){
 692
 693   int iResult = 0;
 694   UInt_t numHLTTracks = 0;
 695   fHLTtracksPresent = kFALSE;
 696
 697   if (!fHLTTracks){
 698     LogError("HLT track vector instance not available.");
 699     return 0;
 700   }
 701
 702   fHLTTracks->clear();
 703
 704   for (const AliHLTComponentBlockData *pBlock = GetFirstInputBlock(kAliHLTDataTypeTrack | kAliHLTDataOriginTPC);
 705        pBlock != NULL; pBlock = GetNextInputBlock()) {
 706     LogDebug("#hlttrk - data block with HLT raw tracks received");
 707     fHLTtracksPresent = kTRUE;
 708
 709     // vector<AliHLTGlobalBarrelTrack> hltTracks;
 710     if ((iResult = AliHLTGlobalBarrelTrack::ConvertTrackDataArray(reinterpret_cast<const AliHLTTracksData*>(pBlock->fPtr), pBlock->fSize, *fHLTTracks)) >= 0) {
 711       for (vector<AliHLTGlobalBarrelTrack>::iterator element = fHLTTracks->begin();
 712            element != fHLTTracks->end(); element++) {
 713
 714         numHLTTracks++;
 715         //## AliHLTGlobalBarrelTrack -> AliKalmanTrack -> AliExternalTrackParam
 716
 717       } // loop over HLT tracks
 718     }
 719   } // loop over data blocks
 720
 721   LogDebug("#hlttrk - %d HLT raw tracks found", numHLTTracks);
 722
 723   return iResult;
 724
 725 }
 726
 727 int AliHLTTRDTriggerComponent::PrepareTRDData() {
 728
 729   int result = 1;
 730
 731   fSectorsWithData = 0;
 732   fTrackingData->Clear();
 733   for (const AliHLTComponentBlockData* datablock = GetFirstInputBlock(AliHLTTRDDefinitions::fgkOnlineDataType);
 734        datablock != NULL;
 735        datablock = GetNextInputBlock())
 736     {
 737       fSectorsWithData |= datablock->fSpecification;
 738       fTrackingData->Decompress(datablock->fPtr, datablock->fSize, kTRUE);
 739     }
 740
 741   fTrackingData->SetGtuPtMultiplierFromMagField(GetBz()); // used for sign correction
 742
 743   fTrackingData->PrintSummary("trigger component");
 744
 745   return result;
 746
 747 }
 748
 749 int AliHLTTRDTriggerComponent::MatchTRDTracksESD(){
 750
 751   if (!fEsdEvent) {
 752     LogError("ESD event data not available in MatchTRDTracks().");
 753     return 0;
 754   }
 755
 756   if (fHistoMode == 0){
 757     UInt_t numHists = fHistArray->GetEntries();
 758     for (UInt_t iHist = 0; iHist < numHists; ++iHist)
 759       if (fHistArray->At(iHist))
 760         ((TH1*) (fHistArray->At(iHist)))->Reset();
 761   }
 762
 763   int result = 1;
 764   unsigned int numRefTracks = fEsdEvent->GetNumberOfTracks();
 765   unsigned int numGtuTracks;
 766   Int_t refTrackIndices[fkMaxRefTracksPerStack];
 767   UInt_t refTrackCount = 0;
 768   Bool_t isRelevant;
 769   Double_t distY, distZ;
 770   Double_t matchRating;
 771   Double_t bestMatchRating;
 772   Int_t bestMatchRefIndex;
 773   AliESDtrack* refTrack = NULL;
 774   UInt_t numComparisonsDone = 0;
 775   UInt_t numUnmatchedTracks = 0;
 776   UInt_t numMatchedTracks = 0;
 777   Double_t magField = fEsdEvent->GetMagneticField();
 778
 779   for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack) {
 780     numGtuTracks = fTrackingData->GetNumTracks(iStack);
 781     refTrackCount = 0;
 782
 783     // preselect ESD relevant ESD tracks
 784     for (UInt_t iRefTrack = 0; iRefTrack < numRefTracks; ++iRefTrack) {
 785       refTrack = fEsdEvent->GetTrack(iRefTrack);
 786       isRelevant = (refTrack->Pt() >= fRefTrackSelectionPtThreshold); //## implement properly
 787       if (isRelevant){
 788         if (refTrackCount < fkMaxRefTracksPerStack)
 789           refTrackIndices[refTrackCount++] = iRefTrack;
 790         else {
 791           LogError("number of HLT tracks exceeding limit of %d. Skipping some tracks.", fkMaxRefTracksPerStack);
 792           break;
 793         }
 794       }
 795     }
 796
 797     // try to match GTU track with ref tracks
 798     for (UInt_t iGtuTrack = 0; iGtuTrack < numGtuTracks; ++iGtuTrack) {
 799       bestMatchRating = 0.;
 800       bestMatchRefIndex = -1;
 801       Double_t gpt = fTrackingData->GetTrackPt(iStack, iGtuTrack);
 802       UShort_t gpid = fTrackingData->GetTrackPID(iStack, iGtuTrack);
 803       UShort_t layerMask = fTrackingData->GetTrackLayerMask(iStack, iGtuTrack);
 804
 805       Float_t trklLocalY[fkTRDLayers];
 806       Int_t trklBinZ[fkTRDLayers];
 807       for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){
 808         if ((layerMask >> iLayer) & 1){
 809           trklLocalY[iLayer] = fTrackingData->GetTrackTrackletLocalY(iStack, iGtuTrack, iLayer);
 810           trklBinZ[iLayer] = fTrackingData->GetTrackTrackletBinZ(iStack, iGtuTrack, iLayer);
 811         }
 812       }
 813
 814       for (UInt_t iRefTrack = 0; iRefTrack < refTrackCount; ++iRefTrack) {
 815         refTrack = fEsdEvent->GetTrack(refTrackIndices[iRefTrack]);
 816
 817 //      if (!CheckRefTrackCuts(refTrack))
 818 //        continue;
 819
 820         numComparisonsDone++;
 821
 822 //      if ((!refTrack->GetOuterParam()) && (!refTrack->GetInnerParam()))  // use tracks with TPC outer param only (all HLT tracks have it anyways)
 823 //        continue;
 824
 825         Int_t distRes = EstimateTrackDistance(refTrack, iStack, layerMask, trklLocalY, trklBinZ,
 826                                               magField, &distY, &distZ);
 827
 828         if (fDebugLevel >= 3){
 829           printf("CHECKMATCH = %i   distY %.2f   distZ %.2f   pt: %.2f   %.2f\n",
 830                 distRes, distY, distZ, gpt, refTrack->GetSignedPt());
 831         }
 832
 833         if (distRes == 0){
 834           matchRating = RateTrackMatch(distY, distZ, gpt, refTrack->GetSignedPt());
 835         } else {
 836           matchRating = 0.;
 837         }
 838
 839         if ((matchRating >= fMatchRatingThreshold) && (matchRating > bestMatchRating)) {
 840           bestMatchRefIndex = refTrackIndices[iRefTrack];
 841           bestMatchRating = matchRating;
 842         } else if (matchRating > bestMatchRating)
 843           bestMatchRating = matchRating;
 844
 845 //      DbgLog("", Form("#match - comparing GTU track %d in S%02d-%d with ref track %d: [gpt: %+5.2f rpt: %+5.2f] dy: %.1f dz: %.1f rating: %.1f",
 846 //                      iGtuTrack, iStack/5, iStack%5, refTrackIndices[iRefTrack],
 847 //                      fTrackingData->GetTrackPt(iStack, iGtuTrack), refTrack->GetSignedPt(),
 848 //                      distY, distZ, matchRating));
 849
 850       } // loop over ref tracks in stack
 851
 852       fHistMatchRating->Fill(bestMatchRating);
 853       fHistMatchRatingByPt->Fill(bestMatchRating, TMath::Abs(gpt));
 854       fHistMatchRatingByPid->Fill(bestMatchRating, gpid);
 855       fHistTrackPt->Fill(gpt);
 856       fHistTrackPid->Fill(gpid);
 857
 858       if (fExtendedHistos){
 859         fHistMatchedRefTrackPid->Fill(refTrack->Pt(), refTrack->GetTPCsignal());
 860       }
 861
 862       if (bestMatchRefIndex >= 0){
 863         // GTU track has matching reference track
 864         refTrack = fEsdEvent->GetTrack(bestMatchRefIndex);
 865         Double_t rpt = refTrack->GetSignedPt();
 866         LogDebug("#match - match found: rating %.2f, gpt: %+5.2f, rpt: %+5.2f (%i) -> diff: %.2f%%",
 867                             bestMatchRating,gpt, rpt, bestMatchRefIndex, (gpt - rpt)/rpt*100.);
 868         fTrackingData->SetTrackAddInfo(iStack, iGtuTrack, bestMatchRefIndex);
 869         if (fDebugLevel >= 3){
 870           LogDebug("#match-info rating: %.2f gpt: %+.2f  matchref: %d  rpt: %+.2f       gpid: %d  S%02d-%d %d   rpid: %.3f",
 871                           bestMatchRating, gpt, bestMatchRefIndex, rpt, gpid,
 872                           iStack/5, iStack%5, iGtuTrack, refTrack->GetTPCsignal());
 873         }
 874         numMatchedTracks++;
 875
 876         fHistTrackPtMatched->Fill(gpt);
 877         fHistTrackPtCorr->Fill(rpt, gpt);
 878         fHistTrackPidMatched->Fill(gpid);
 879
 880       } else {
 881         if (fDebugLevel >= 3){
 882           LogDebug("#match-info rating: %.2f gpt: %+.2f  no ref matching  gpid: %d  S%02d-%d %d",
 883                           bestMatchRating, gpt, gpid,
 884                           iStack/5, iStack%5, iGtuTrack);
 885         }
 886         numUnmatchedTracks++;
 887       }
 888
 889     } // loop over gtu tracks in stack
 890
 891   } // loop over stacks
 892
 893   LogInfo("#match - %d matched GTU tracks, %d unmatched",
 894                       numMatchedTracks, numUnmatchedTracks);
 895
 896   fHistTrackMatchingCombinations->Fill(0., (Double_t)(numRefTracks * fTrackingData->GetNumTracks()));  // index 0: full combinatorics
 897   fHistTrackMatchingCombinations->Fill(1., numComparisonsDone);  // index 1: track matching comparisons actually done
 898
 899   return result;
 900
 901 }
 902
 903 int AliHLTTRDTriggerComponent::MatchTRDTracksHLT(){
 904
 905   if (!fHLTTracks){
 906     LogError("HLT track data available.");
 907     return 0;
 908   }
 909
 910   Double_t magField = GetBz();
 911   unsigned int numGtuTracks;
 912   unsigned int numHLTTracks = 0;
 913   Bool_t isRelevant;
 914   Double_t distY, distZ;
 915   Double_t matchRating;
 916   Double_t bestMatchRating;
 917   Int_t bestMatchRefIndex;
 918   UInt_t numComparisonsDone = 0;
 919   UInt_t numUnmatchedTracks = 0;
 920   UInt_t numMatchedTracks = 0;
 921
 922   Bool_t hltTrackPreSel[fkMaxRefTracks];
 923   UInt_t iHltTrack;
 924
 925   for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack) {
 926     numGtuTracks = fTrackingData->GetNumTracks(iStack);
 927
 928     // preselect relevant HLT tracks
 929     iHltTrack = 0;
 930     for (vector<AliHLTGlobalBarrelTrack>::iterator element = fHLTTracks->begin();
 931          element != fHLTTracks->end(); element++) {
 932       numHLTTracks++;
 933       isRelevant = (element->Pt() >= fRefTrackSelectionPtThreshold); //## implement properly
 934
 935       //## use cuts here
 936       hltTrackPreSel[iHltTrack++] = isRelevant;
 937       if (iHltTrack >= fkMaxRefTracks)
 938         LogError("maximum number of HLT tracks exceeded.");
 939     } // loop over HLT tracks;
 940
 941     // search for matching HLT track for each GTU track
 942     for (UInt_t iGtuTrack = 0; iGtuTrack < numGtuTracks; ++iGtuTrack) {
 943       bestMatchRating = 0.;
 944       bestMatchRefIndex = -1;
 945       Double_t gpt = fTrackingData->GetTrackPt(iStack, iGtuTrack);
 946       UShort_t gpid = fTrackingData->GetTrackPID(iStack, iGtuTrack);
 947       UShort_t layerMask = fTrackingData->GetTrackLayerMask(iStack, iGtuTrack);
 948
 949       Float_t trklLocalY[fkTRDLayers];
 950       Int_t trklBinZ[fkTRDLayers];
 951       for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){
 952         if ((layerMask >> iLayer) & 1){
 953           trklLocalY[iLayer] = fTrackingData->GetTrackTrackletLocalY(iStack, iGtuTrack, iLayer);
 954           trklBinZ[iLayer] = fTrackingData->GetTrackTrackletBinZ(iStack, iGtuTrack, iLayer);
 955         }
 956       }
 957
 958       iHltTrack = 0;
 959       for (vector<AliHLTGlobalBarrelTrack>::iterator element = fHLTTracks->begin();
 960            element != fHLTTracks->end(); element++) {
 961         if (!hltTrackPreSel[iHltTrack]){
 962           iHltTrack++;
 963           continue;
 964         }
 965
 966         // compare GTU track and relevant HLT track
 967         numComparisonsDone++;
 968
 969         AliHLTGlobalBarrelTrack hltPar(*element);
 970         Int_t distRes = EstimateTrackDistance(&hltPar, iStack, layerMask, trklLocalY, trklBinZ,
 971                                                magField, &distY, &distZ);
 972
 973         if (fDebugLevel >= 3){
 974           printf("CHECKMATCH = %i   distY %.2f   distZ %.2f   pt: %.2f   %.2f\n",
 975                  distRes, distY, distZ, gpt, element->GetSignedPt());
 976         }
 977
 978         if (distRes == 0){
 979           matchRating = RateTrackMatch(distY, distZ, gpt, element->GetSignedPt());
 980         } else {
 981           matchRating = 0.;
 982         }
 983
 984         if ((matchRating >= fMatchRatingThreshold) && (matchRating > bestMatchRating)) {
 985           bestMatchRefIndex = iHltTrack;
 986           bestMatchRating = matchRating;
 987         } else if (matchRating > bestMatchRating)
 988           bestMatchRating = matchRating;
 989
 990
 991         iHltTrack++;
 992       } // loop over HLT tracks;
 993
 994       fHistMatchRating->Fill(bestMatchRating);
 995       fHistMatchRatingByPt->Fill(bestMatchRating, TMath::Abs(gpt));
 996       fHistMatchRatingByPid->Fill(bestMatchRating, gpid);
 997       fHistTrackPt->Fill(gpt);
 998       fHistTrackPid->Fill(gpid);
 999
1000       if (bestMatchRefIndex >= 0){
1001         // GTU track has matching reference track
1002         Double_t rpt = fHLTTracks->at(bestMatchRefIndex).GetSignedPt();
1003         LogDebug("#match - match found: rating %.2f, gpt: %+5.2f, rpt: %+5.2f (%i) -> diff: %.2f%%",
1004                             bestMatchRating,gpt, rpt, bestMatchRefIndex, (gpt - rpt)/rpt*100.);
1005         fTrackingData->SetTrackAddInfo(iStack, iGtuTrack, bestMatchRefIndex);
1006
1007 //      if (fExtendedHistos){
1008 //        fHistMatchedRefTrackPid->Fill(element->Pt(), element->GetTPCsignal());
1009 //      }
1010 //
1011 //      if (fDebugLevel >= 3){
1012 //        LogDebug("#match-info rating: %.2f gpt: %+.2f  matchref: %d  rpt: %+.2f       gpid: %d  S%02d-%d %d   rpid: %.3f",
1013 //                        bestMatchRating, gpt, bestMatchRefIndex, rpt, gpid,
1014 //                        iStack/5, iStack%5, iGtuTrack, refTrack->GetTPCsignal());
1015 //      }
1016         numMatchedTracks++;
1017
1018         fHistTrackPtMatched->Fill(gpt);
1019         fHistTrackPtCorr->Fill(rpt, gpt);
1020         fHistTrackPidMatched->Fill(gpid);
1021
1022       } else {
1023         if (fDebugLevel >= 3){
1024           LogDebug("#match-info rating: %.2f gpt: %+.2f  no ref matching  gpid: %d  S%02d-%d %d",
1025                           bestMatchRating, gpt, gpid,
1026                           iStack/5, iStack%5, iGtuTrack);
1027         }
1028         numUnmatchedTracks++;
1029       }
1030
1031     } // loop over gtu tracks
1032   } // loop over stacks
1033
1034   LogInfo("#match - %d matched GTU tracks, %d unmatched (%d comparisons)",
1035           numMatchedTracks, numUnmatchedTracks, numComparisonsDone);
1036
1037   fHistTrackMatchingCombinations->Fill(0., (Double_t)(numHLTTracks * fTrackingData->GetNumTracks()));  // index 0: full combinatorics
1038   fHistTrackMatchingCombinations->Fill(1., numComparisonsDone);  // index 1: track matching comparisons actually done
1039
1040   return kTRUE;
1041 }
1042
1043 int AliHLTTRDTriggerComponent::MatchTRDTracks(){
1044
1045   if (!fEsdEvent) {
1046     LogError("ESD event data not available in MatchTRDTracks().");
1047     return 0;
1048   }
1049
1050   if (fHistoMode == 0){
1051     UInt_t numHists = fHistArray->GetEntries();
1052     for (UInt_t iHist = 0; iHist < numHists; ++iHist)
1053       if (fHistArray->At(iHist))
1054         ((TH1*) (fHistArray->At(iHist)))->Reset();
1055   }
1056
1057   int result = 1;
1058   unsigned int numRefTracks = fEsdEvent->GetNumberOfTracks();
1059   unsigned int numGtuTracks;
1060   Int_t refTrackIndices[fkMaxRefTracksPerStack];
1061   UInt_t refTrackCount = 0;
1062   Bool_t isRelevant;
1063   Double_t distY, distZ;
1064   Double_t matchRating;
1065   Double_t bestMatchRating;
1066   Int_t bestMatchRefIndex;
1067   AliESDtrack* refTrack;
1068   UInt_t numComparisonsDone = 0;
1069   UInt_t numUnmatchedTracks = 0;
1070   UInt_t numMatchedTracks = 0;
1071
1072   for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack) {
1073     numGtuTracks = fTrackingData->GetNumTracks(iStack);
1074     refTrackCount = 0;
1075
1076     // preselect ESD relevant ESD tracks
1077     for (UInt_t iRefTrack = 0; iRefTrack < numRefTracks; ++iRefTrack) {
1078       refTrack = fEsdEvent->GetTrack(iRefTrack);
1079       isRelevant = (refTrack->Pt() >= fRefTrackSelectionPtThreshold);
1080       if (isRelevant){
1081         if (refTrackCount < fkMaxRefTracksPerStack)
1082           refTrackIndices[refTrackCount++] = iRefTrack;
1083         else {
1084           LogError("number of HLT tracks exceeding limit of %d. Skipping some tracks.", fkMaxRefTracksPerStack);
1085           break;
1086         }
1087       }
1088     }
1089
1090     // try to match GTU track with ref tracks
1091     for (UInt_t iGtuTrack = 0; iGtuTrack < numGtuTracks; ++iGtuTrack) {
1092       bestMatchRating = 0.;
1093       bestMatchRefIndex = -1;
1094       Double_t gpt = fTrackingData->GetTrackPt(iStack, iGtuTrack);
1095       UShort_t gpid = fTrackingData->GetTrackPID(iStack, iGtuTrack);
1096       UShort_t layerMask = fTrackingData->GetTrackLayerMask(iStack, iGtuTrack);
1097
1098       for (UInt_t iRefTrack = 0; iRefTrack < refTrackCount; ++iRefTrack) {
1099         refTrack = fEsdEvent->GetTrack(refTrackIndices[iRefTrack]);
1100
1101         if (!CheckRefTrackCuts(refTrack))
1102           continue;
1103
1104         numComparisonsDone++;
1105
1106 //      if ((!refTrack->GetOuterParam()) && (!refTrack->GetInnerParam()))  // use tracks with TPC outer param only (all HLT tracks have it anyways)
1107 //        continue;
1108
1109         Float_t trklLocalY[fkTRDLayers];
1110         Int_t trklBinZ[fkTRDLayers];
1111         for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){
1112           if ((layerMask >> iLayer) & 1){
1113             trklLocalY[iLayer] = fTrackingData->GetTrackTrackletLocalY(iStack, iGtuTrack, iLayer);
1114             trklBinZ[iLayer] = fTrackingData->GetTrackTrackletBinZ(iStack, iGtuTrack, iLayer);
1115           }
1116         }
1117
1118         Int_t distRes = EstimateTrackDistance(refTrack, iStack, layerMask, trklLocalY, trklBinZ,
1119                                               fEsdEvent->GetMagneticField(), &distY, &distZ);
1120
1121         if (distRes == 0){
1122           matchRating = RateTrackMatch(distY, distZ, gpt, refTrack->GetSignedPt());
1123         } else {
1124           matchRating = 0.;
1125         }
1126
1127         if ((matchRating >= fMatchRatingThreshold) && (matchRating > bestMatchRating)) {
1128           bestMatchRefIndex = refTrackIndices[iRefTrack];
1129           bestMatchRating = matchRating;
1130         } else if (matchRating > bestMatchRating)
1131           bestMatchRating = matchRating;
1132
1133 //      DbgLog("", Form("#match - comparing GTU track %d in S%02d-%d with ref track %d: [gpt: %+5.2f rpt: %+5.2f] dy: %.1f dz: %.1f rating: %.1f",
1134 //                      iGtuTrack, iStack/5, iStack%5, refTrackIndices[iRefTrack],
1135 //                      fTrackingData->GetTrackPt(iStack, iGtuTrack), refTrack->GetSignedPt(),
1136 //                      distY, distZ, matchRating));
1137
1138       } // loop over ref tracks in stack
1139
1140       fHistMatchRating->Fill(bestMatchRating);
1141       fHistMatchRatingByPt->Fill(bestMatchRating, TMath::Abs(gpt));
1142       fHistMatchRatingByPid->Fill(bestMatchRating, gpid);
1143       fHistTrackPt->Fill(gpt);
1144       fHistTrackPid->Fill(gpid);
1145
1146       if (bestMatchRefIndex >= 0){
1147         // GTU track has matching reference track
1148         refTrack = fEsdEvent->GetTrack(bestMatchRefIndex);
1149         Double_t rpt = refTrack->GetSignedPt();
1150         LogDebug("#match - match found: rating %.2f, gpt: %+5.2f, rpt: %+5.2f (%i) -> diff: %.2f%%",
1151                             bestMatchRating,gpt, rpt, bestMatchRefIndex, (gpt - rpt)/rpt*100.);
1152         fTrackingData->SetTrackAddInfo(iStack, iGtuTrack, bestMatchRefIndex);
1153
1154         if (fExtendedHistos){
1155           fHistMatchedRefTrackPid->Fill(refTrack->Pt(), refTrack->GetTPCsignal());
1156         }
1157
1158         if (fDebugLevel >= 3){
1159           LogDebug("#match-info rating: %.2f gpt: %+.2f  matchref: %d  rpt: %+.2f       gpid: %d  S%02d-%d %d   rpid: %.3f",
1160                           bestMatchRating, gpt, bestMatchRefIndex, rpt, gpid,
1161                           iStack/5, iStack%5, iGtuTrack, refTrack->GetTPCsignal());
1162         }
1163         numMatchedTracks++;
1164
1165         fHistTrackPtMatched->Fill(gpt);
1166         fHistTrackPtCorr->Fill(rpt, gpt);
1167         fHistTrackPidMatched->Fill(gpid);
1168
1169       } else {
1170         if (fDebugLevel >= 3){
1171           LogDebug("#match-info rating: %.2f gpt: %+.2f  no ref matching  gpid: %d  S%02d-%d %d",
1172                           bestMatchRating, gpt, gpid,
1173                           iStack/5, iStack%5, iGtuTrack);
1174         }
1175         numUnmatchedTracks++;
1176       }
1177
1178     } // loop over gtu tracks in stack
1179
1180   } // loop over stacks
1181
1182   LogInfo("#match - %d matched GTU tracks, %d unmatched",
1183                       numMatchedTracks, numUnmatchedTracks);
1184
1185   fHistTrackMatchingCombinations->Fill(0., (Double_t)(numRefTracks * fTrackingData->GetNumTracks()));  // index 0: full combinatorics
1186   fHistTrackMatchingCombinations->Fill(1., numComparisonsDone);  // index 1: track matching comparisons actually done
1187
1188   return result;
1189 }
1190
1191 void AliHLTTRDTriggerComponent::DumpTrackingData(){
1192
1193   if (fTrackingData->GetNumTracklets() + fTrackingData->GetNumTracks() == 0)
1194     return;
1195
1196   TString trklStr("");
1197   TString matchStr("");
1198   UShort_t layerMask;
1199
1200 //  for (UShort_t iSector = 0; iSector < 18; ++iSector){
1201 //    if (fTrackingData->GetSectorTrgWord(iSector) != ((UInt_t)0x2345352 ^ ((UInt_t)iSector + 34)))
1202 //      LogError("invalid sector trigger word in sector %02d: trg word is 0x%08x, should be 0x%08x",
1203 //             iSector, fTrackingData->GetSectorTrgWord(iSector), ((UInt_t)0x2345352 ^ ((UInt_t)iSector + 34)));
1204 //    for (UShort_t iStack = 0; iStack < 5; ++iStack){
1205 //      ULong64_t dwl = ((ULong64_t)0xaffe << 16) | (ULong64_t)iSector << 8 | iStack;
1206 //      ULong64_t dw = (((ULong64_t)0xbead << 16) | (ULong64_t)iSector << 8 | iStack) | ((ULong64_t)dwl << 32);
1207 //      if (fTrackingData->GetStackTrgWord(iSector, iStack) != dw)
1208 //      LogError("stack %02d-%d trg word is 0x%016llx, should be 0x%016llx",
1209 //               iSector, iStack, fTrackingData->GetStackTrgWord(iSector, iStack), dw);
1210 //    }
1211 //  }
1212
1213   for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack){
1214     for (Int_t iTrk = 0; iTrk < fTrackingData->GetNumTracks(iStack); ++iTrk){
1215
1216       layerMask = fTrackingData->GetTrackLayerMask(iStack, iTrk);
1217
1218       trklStr = Form("trkl: ");
1219       for (Short_t iLayer = 5; iLayer >= 0; --iLayer){
1220         if ((layerMask >> iLayer) & 1)
1221           trklStr += Form("0x%08x (%+8.3f)  ",
1222                           fTrackingData->GetTrackTrackletWord(iStack, iTrk, iLayer),
1223                           fTrackingData->GetTrackTrackletLocalY(iStack, iTrk, iLayer));
1224         else
1225           trklStr += "---------------------  ";
1226       } // loop over layers
1227       trklStr.Remove(trklStr.Length() - 2, 2);
1228
1229       if (fTrackingData->GetTrackAddInfo(iStack, iTrk) >= 0){
1230         Double_t rpt = (fESDtracksPresent) ?  fEsdEvent->GetTrack(fTrackingData->GetTrackAddInfo(iStack, iTrk))->GetSignedPt() :
1231           fHLTTracks->at(fTrackingData->GetTrackAddInfo(iStack, iTrk)).GetSignedPt();
1232         matchStr = Form("mpt: %+7.2f", rpt);
1233       } else
1234         matchStr = "unmatched";
1235
1236       if (fDebugLevel >= 3){
1237
1238         printf("###DOTDA EV%04llu  GTU TRACK - S%02d-%d  pt: %+7.2f  pid: %3d  lm: 0x%02x %s  %s\n",
1239                fEventId,
1240                iStack/5, iStack%5, fTrackingData->GetTrackPt(iStack, iTrk), fTrackingData->GetTrackPID(iStack, iTrk),
1241                layerMask, trklStr.Data(),
1242                matchStr.Data());
1243
1244         printf("###DOTDB EV%04llu  GTU TRACK - S%02d-%d  pt: %+7.2f  pid: %3d  lm: 0x%02x %s\n",
1245                fEventId,
1246                iStack/5, iStack%5, fTrackingData->GetTrackPt(iStack, iTrk), fTrackingData->GetTrackPID(iStack, iTrk),
1247                layerMask, trklStr.Data());
1248       }
1249
1250
1251       // paranoia checks
1252       for (Short_t iLayer = 5; iLayer >= 0; --iLayer){
1253         if (((layerMask >> iLayer) & 1) && (fTrackingData->GetTrackTrackletWord(iStack, iTrk, iLayer) == 0x10001000))
1254           LogError("invalid layer mask / tracklet value combination A");
1255
1256         if ((((layerMask >> iLayer) & 1) == 0) && (fTrackingData->GetTrackTrackletWord(iStack, iTrk, iLayer) != 0x10001000))
1257           LogError("invalid layer mask / tracklet value combination B");
1258       }
1259
1260     } // loop over tracks in stack
1261   } // loop over stacks
1262
1263 }
1264
1265 void AliHLTTRDTriggerComponent::AssignTrackInfo(TString* infoStr, const UInt_t stack, const UInt_t trackIndex, const char* flagStr) {
1266
1267   TString flags(flagStr);
1268   Bool_t appendRefData = kFALSE;
1269
1270   if (fTrackingData->GetTrackAddInfo(stack, trackIndex) >= 0){
1271     appendRefData = kTRUE;
1272     if (flags.First('M') < 0)
1273       flags += "M";
1274   }
1275
1276   *infoStr = Form("EXCH-TRK-INFO DS<%s> CH<%s> EV%05llu SEC%02d-%d-%d  gpt: %+6.1f  gpid: %3d  flags: %s",
1277                   "TRDHLTTRG",
1278                   fChunkId->Data(), fEventId,
1279                   stack/5, stack%5, trackIndex,
1280                   fTrackingData->GetTrackPt(stack, trackIndex),
1281                   fTrackingData->GetTrackPID(stack, trackIndex), flags.Data());
1282
1283   if (appendRefData){
1284     Double_t rpt = (fESDtracksPresent) ? fEsdEvent->GetTrack(fTrackingData->GetTrackAddInfo(stack, trackIndex))->GetSignedPt() :
1285       fHLTTracks->at(fTrackingData->GetTrackAddInfo(stack, trackIndex)).GetSignedPt();
1286     *infoStr += Form("  rpt: %+6.1f", rpt);
1287   }
1288
1289 }
1290
1291 #ifdef __TRDHLTDEBUG
1292 void AliHLTTRDTriggerComponent::RenderEvent(const Bool_t showGtuTracks, const Bool_t showTracklets, const Bool_t showRefTracks) {
1293
1294   if ((!fEventDisplay) || (!fEsdEvent))
1295     return;
1296
1297   LogDebug("rendering event");
1298
1299   const Float_t refTrackPtDisplayThreshold = 0.7;
1300   const Float_t trackPtEmphasizeThreshold = 1.8;
1301
1302   fEventDisplay->Reset();
1303   fEventDisplay->SetMagField(fEsdEvent->GetMagneticField());
1304   for (UInt_t iDet = 0; iDet < fkTRDStacks*fkTRDLayers; ++iDet)
1305     fEventDisplay->SetChamberState(iDet/fkTRDLayers, iDet%6, ((fSectorsWithData >> (iDet/30)) & 1) );
1306
1307   if (showTracklets){
1308     for (UShort_t iDet = 0; iDet < fkTRDStacks*fkTRDLayers; ++iDet){
1309       UInt_t numTrkl = fTrackingData->GetNumTracklets(iDet);
1310       for (UInt_t iTrkl = 0; iTrkl < numTrkl; ++iTrkl){
1311         fEventDisplay->AddTracklet(fTrackingData->GetTracklet(iDet, iTrkl));
1312       }
1313     }
1314   }
1315
1316   if (showGtuTracks){
1317     for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack){
1318       for (Int_t iTrk = 0; iTrk < fTrackingData->GetNumTracks(iStack); ++iTrk){
1319         AliESDTrdTrack* track = fTrackingData->GetTrack(iStack, iTrk);
1320         Int_t trkIndex = fEventDisplay->AddTrack(track, (TMath::Abs(track->Pt()) >= trackPtEmphasizeThreshold) ? (kMagenta + 2) : kMagenta, 1, 1);
1321         if (fTrackingData->GetTrackAddInfo(iStack, iTrk) >= 0)
1322           fEventDisplay->SetTrackTrackletStyle(trkIndex, kRed, 12);
1323         else
1324           fEventDisplay->SetTrackTrackletStyle(trkIndex, kViolet - 5, 12);
1325       }
1326     }
1327   }
1328
1329   unsigned int numRefTracks = fEsdEvent->GetNumberOfTracks();
1330   unsigned int numRefTracksRendered = 0;
1331   if (showRefTracks){
1332     // check for some marks for rendering
1333     UShort_t marks[40000];
1334     memset(marks, 0, sizeof(UShort_t)*40000);
1335     for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack){
1336       for (Int_t iTrk = 0; iTrk < fTrackingData->GetNumTracks(iStack); ++iTrk){
1337         if (fTrackingData->GetTrackAddInfo(iStack, iTrk) >= 0){
1338           marks[fTrackingData->GetTrackAddInfo(iStack, iTrk)] |= 1;
1339           if (
1340               (fTrackingData->GetTrackPID(iStack, iTrk) >= fElectronTriggerPIDThresholdHQU) &&
1341               (TMath::Abs(fTrackingData->GetTrackPt(iStack, iTrk)) >= fElectronTriggerPtThresholdHQU))
1342             marks[fTrackingData->GetTrackAddInfo(iStack, iTrk)] |= 2;
1343         }
1344       }
1345     }
1346
1347     // add to rendering
1348     for (unsigned int iTrack = 0; iTrack < numRefTracks; ++iTrack){
1349       AliESDtrack* track = fEsdEvent->GetTrack(iTrack);
1350       if ((track->Pt() >= refTrackPtDisplayThreshold) || (marks[iTrack] != 0)){
1351         Color_t color = (track->Pt() >= trackPtEmphasizeThreshold) ? (kGray + 2) : kGray;
1352         UShort_t width = 1;
1353         UShort_t style = 1;
1354         if (marks[iTrack] & 0x2){
1355           color = kRed + 1;
1356           width = 6;
1357         }
1358         if (!CheckRefTrackCuts(track))
1359           style = 2;
1360
1361         fEventDisplay->AddTrack(track, color, width, style);
1362         numRefTracksRendered++;
1363       }
1364     }
1365   }
1366
1367   if (!fEventDisplay->IsEmpty()){
1368     fEventDisplay->SetTitle("");
1369     fEventDisplay->SetSetupText("HLT", Form("%.1f#scale[0.5]{ }/#scale[0.5]{ }%.1f",
1370                                           refTrackPtDisplayThreshold, trackPtEmphasizeThreshold));
1371     fEventDisplay->SetBottomText(Form("%05i-%s-%05llu",
1372                                       fRunNumber, fChunkId->Data(), fEventId));
1373     fEventDisplay->SetBottomTextRight(Form("%d (%d) HLT tracks, %d tracklets, %d GTU tracks",
1374                                            numRefTracks,
1375                                            numRefTracksRendered,
1376                                            fTrackingData->GetNumTracklets(),
1377                                            fTrackingData->GetNumTracks()));
1378     fEventDisplay->SetLook(AliTRDtrackingEventDisplay::dmMediumLight);
1379     fEventDisplay->SetDisplayMode(AliTRDtrackingEventDisplay::dmFullXY);
1380     fEventDisplay->SaveToFile(Form("display/event-%s-%05llu.eps",
1381                                    fChunkId->Data(), fEventId));
1382   }
1383
1384 }
1385 #endif
1386
1387 Bool_t AliHLTTRDTriggerComponent::TRDElectronTrigger(const char *ident, const Double_t minPt, const UShort_t minPID){
1388
1389   LogDebug("Electron trigger processing (%s: pt>=%.1f, pid>=%d)...", ident, minPt, minPID);
1390
1391   UInt_t numTracks;
1392   Bool_t highPtElectronSeenGTU = kFALSE;
1393   Bool_t highPtElectronSeen = kFALSE;
1394   UInt_t truncMeanPID = 0;
1395   TString trackExchangeInfo("");
1396   TString flags("");
1397
1398   UInt_t trdSectorTrgContribs = 0;
1399   for (UShort_t iSector = 0; iSector < fkTRDSectors; ++iSector)
1400     trdSectorTrgContribs |= fTrackingData->GetSectorTrgContribs(iSector);
1401   if ((trdSectorTrgContribs >> 5) & 0x3)
1402     fIsTRDElectronEvent = kTRUE;
1403   else
1404     fIsTRDElectronEvent = kFALSE;
1405
1406   if (fIsMinBiasEvent)
1407     fHistElectronTriggerBaseMinBias->Fill(0.);
1408
1409   if (fIsTRDElectronEvent)
1410     fHistElectronTriggerBaseTrdL1->Fill(0.);
1411
1412   if ((fElectronTriggerOnL1TrgOnly) && (!fIsTRDElectronEvent)){
1413     // evaluate trigger for events with TRD L1 electron trigger fired
1414     DbgLog("skipping %s electron trigger evaluation for event, because no TRD trigger flag set (ctbs: 0x%02x)",
1415            ident, trdSectorTrgContribs);
1416     return 0;
1417   }
1418
1419   for (UShort_t iStack = 0; iStack < fkTRDStacks; ++iStack){
1420     numTracks = fTrackingData->GetNumTracks(iStack);
1421     for (UInt_t iTrk = 0; iTrk < numTracks; ++iTrk){
1422       Double_t gpt = fTrackingData->GetTrackPt(iStack, iTrk);
1423       Bool_t trdElectronCandidate = kFALSE;
1424       Bool_t hltElectronCandidate = kFALSE;
1425
1426       // re-evaluate GTU only decision for comparison
1427       if (
1428           (TMath::Abs(gpt) >= minPt) &&
1429           (fTrackingData->GetTrackPID(iStack, iTrk) >= minPID)
1430           ) {
1431         trdElectronCandidate = kTRUE;
1432         highPtElectronSeenGTU = kTRUE;
1433         fHistElectronCandidatePt->Fill(gpt);
1434         fHistElectronCandidatePid->Fill(fTrackingData->GetTrackPID(iStack, iTrk));
1435
1436         if (fExtendedHistos){
1437
1438           // idea to be checked
1439           truncMeanPID = 0;
1440           Short_t minLyr = -1;
1441           UShort_t minValue = 255;
1442           Short_t maxLyr = -1;
1443           UShort_t maxValue = 0;
1444           UShort_t lyrCount = 0;
1445           UInt_t layerMask = fTrackingData->GetTrackLayerMask(iStack, iTrk);
1446
1447           // scan for min & max values
1448           for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){
1449             if ((layerMask >> iLayer) & 1){
1450               if (fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer) < minValue){
1451                 minValue = fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer);
1452                 minLyr = iLayer;
1453               }
1454               if (fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer) > maxValue){
1455                 maxValue = fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer);
1456                 maxLyr = iLayer;
1457               }
1458             }
1459           } // loop over layers
1460
1461           // calculate trunc mean
1462           for (UShort_t iLayer = 0; iLayer < fkTRDLayers; ++iLayer){
1463             if (((layerMask >> iLayer) & 1) && (iLayer != minLyr) && (iLayer != maxLyr)){
1464               truncMeanPID += fTrackingData->GetTrackTrackletPID(iStack, iTrk, iLayer);
1465               lyrCount++;
1466             }
1467           } // loop over layers
1468           truncMeanPID = TMath::Nint((Double_t)(truncMeanPID)/(Double_t)(lyrCount));
1469
1470           fHistPIDvsTruncPID->Fill(fTrackingData->GetTrackPID(iStack, iTrk), truncMeanPID);
1471           if (fTrackingData->GetTrackAddInfo(iStack, iTrk) < 0)
1472             fHistElectronFalsePIDvsTruncPID->Fill(fTrackingData->GetTrackPID(iStack, iTrk), truncMeanPID);
1473
1474         }
1475
1476         LogInspect("#hlt-trd-trg - GTU flagged %s %s high-pt electron seen in S%02d-%d: pt: %+6.1f  pid: %d  [id: S%02d-%d-%d] [trunc-pid: %d]",
1477                    (fTrackingData->GetTrackAddInfo(iStack, iTrk) < 0) ? "unmatched" : "matched", ident,
1478                    iStack/5, iStack%5,
1479                    gpt, fTrackingData->GetTrackPID(iStack, iTrk),
1480                    iStack/5, iStack%5, iTrk,
1481                    truncMeanPID);
1482       }
1483
1484       // evaluate HLT-level trigger decision using match information
1485       if (
1486           (fTrackingData->GetTrackAddInfo(iStack, iTrk) >= 0) &&
1487           (TMath::Abs(gpt) >= minPt) &&
1488           (fTrackingData->GetTrackPID(iStack, iTrk) >= minPID)
1489           ) {
1490         hltElectronCandidate = kTRUE;
1491         highPtElectronSeen = kTRUE;
1492
1493         fHistElectronCandidateMatchedPt->Fill(gpt);
1494         fHistElectronCandidateMatchedPid->Fill(fTrackingData->GetTrackPID(iStack, iTrk));
1495
1496         if (fExtendedHistos){
1497           fHistElectronConfirmedPIDvsTruncPID->Fill(fTrackingData->GetTrackPID(iStack, iTrk), truncMeanPID);
1498         }
1499
1500         Double_t rpt = (fESDtracksPresent) ? fEsdEvent->GetTrack(fTrackingData->GetTrackAddInfo(iStack, iTrk))->GetSignedPt() :
1501           fHLTTracks->at(fTrackingData->GetTrackAddInfo(iStack, iTrk)).GetSignedPt();
1502         LogInspect("#hlt-trd-trg - HLT matched %s high-pt electron seen in S%02d-%d: gpt: %+6.1f  gpid: %d   rpt: %+6.1f  [id: S%02d-%d-%d]",
1503                    ident, iStack/5, iStack%5, fTrackingData->GetTrackPt(iStack, iTrk), fTrackingData->GetTrackPID(iStack, iTrk),
1504                    rpt, iStack/5, iStack%5, iTrk);
1505       }
1506
1507       // log output for subsequent offline analysis
1508       if ((fDebugLevel >= 3) && (trdElectronCandidate || hltElectronCandidate)){
1509         trackExchangeInfo = "";
1510         flags = "";
1511         if (trdElectronCandidate)
1512           flags += "G";
1513         if (hltElectronCandidate)
1514           flags += "H";
1515         AssignTrackInfo(&trackExchangeInfo, iStack, iTrk, flags.Data());
1516         LogDebug("%s\n", trackExchangeInfo.Data());
1517       }
1518
1519     } // loop over tracks in stack
1520   } // loop over stacks
1521
1522   if (highPtElectronSeenGTU || highPtElectronSeen){
1523     LogInspect("#hlt-trd-trg - event triggered by %s electron trigger (TRD L1: %d, TRD HLT: %d)",
1524                ident, highPtElectronSeenGTU, highPtElectronSeen);
1525   }
1526
1527   if (highPtElectronSeenGTU){
1528     if (fIsMinBiasEvent)
1529       fHistElectronTriggerBaseMinBias->Fill(1.);
1530   }
1531
1532   if (highPtElectronSeen){
1533     if (fIsMinBiasEvent)
1534       fHistElectronTriggerBaseMinBias->Fill(2.);
1535     if (fIsTRDElectronEvent)
1536       fHistElectronTriggerBaseTrdL1->Fill(1.);
1537   }
1538
1539   return (highPtElectronSeen) ? kTRUE : kFALSE;
1540 }
1541
1542
1543 int AliHLTTRDTriggerComponent::DoTrigger()
1544 {
1545
1546   fEsdEvent = NULL;
1547   fRunNumber = -1;
1548   int iResult = 0;
1549   UShort_t firedTriggers = 0;
1550
1551   const AliHLTComponentEventData* hltEventData = GetEventData();
1552   fEventId = hltEventData->fEventID;
1553   LogDebug("### START DoTrigger [event id: %llu, %d blocks, size: %d]",
1554                       fEventId, hltEventData->fBlockCnt, hltEventData->fStructSize);
1555
1556 //  LogDebug("### START DoTrigger [event id: %llu, %d blocks, size: %d]",
1557 //                    fEventId, hltEventData->fBlockCnt, hltEventData->fStructSize);
1558
1559   if (!IsDataEvent()) {  // process data events only
1560     IgnoreEvent();
1561     LogDebug("### END   DoTrigger [event id: %llu, %d blocks, size: %d] (skipped: no data event)",
1562              fEventId, hltEventData->fBlockCnt, hltEventData->fStructSize);
1563     return iResult;
1564   }
1565
1566   fTrackingData->SetLogPrefix(Form("TRDHLTGM XXXXXX-%05llu: [TRG] {TrkDat} ", fEventId));
1567
1568   do {
1569
1570     // access to TRD specific data from AliHLTTRDPreprocessorComponent
1571     if (!PrepareTRDData()){
1572       LogError("access to TRD data failed. Skipping event...");
1573       break;
1574     }
1575
1576     if (fTrackingData->GetNumTracks() + fTrackingData->GetNumTracklets() == 0) {
1577       LogDebug("no trigger-relevant TRD information, skipping further event processing");
1578       break;
1579     }
1580
1581     // access to ESD data
1582     if (!PrepareESDData()){
1583       LogInfo("access to ESD event data failed.");
1584     }
1585
1586     // access to alternative HLT data
1587     if (!fESDtracksPresent){
1588       if (!PrepareHLTData()){
1589         LogError("access to HLT event data failed.");
1590       }
1591     }
1592
1593     // match TRD and HLT tracks
1594     if (fESDtracksPresent){
1595       if (!MatchTRDTracksESD()){
1596         LogError("matching TRD tracks to ESD tracks failed. Skipping event...");
1597         break;
1598       }
1599     } else if (fHLTtracksPresent){
1600       if (!MatchTRDTracksHLT()){
1601         LogError("matching TRD tracks to HLT tracks failed. Skipping event...");
1602         break;
1603       }
1604     } else {
1605       LogError("No HLT track information available. Skipping event...");
1606       break;
1607     }
1608
1609 //    if (!MatchTRDTracks()){
1610 //      LogError("matching TRD tracks to TPC tracks failed. Skipping event...");
1611 //      break;
1612 //    }
1613
1614     if (fDebugLevel >= 1)
1615       DumpTrackingData();
1616
1617     // evaluate electron trigger conditions
1618     if (TRDElectronTrigger("HSE", fElectronTriggerPtThresholdHSE, fElectronTriggerPIDThresholdHSE))
1619       firedTriggers |= fkElectronTriggerHSE;
1620
1621     if (TRDElectronTrigger("HQU", fElectronTriggerPtThresholdHQU, fElectronTriggerPIDThresholdHQU))
1622       firedTriggers |= fkElectronTriggerHQU;
1623
1624     break;
1625
1626   } while (1);
1627
1628
1629   // trigger decision
1630   TString description("");
1631   if (firedTriggers & fkElectronTriggerHSE){
1632     if (description.Length() > 0)
1633       description += " ";
1634     description += fgkTriggerDecisionElectronHSE;
1635   }
1636
1637   if (firedTriggers & fkElectronTriggerHQU){
1638     if (description.Length() > 0)
1639       description += " ";
1640     description += fgkTriggerDecisionElectronHQU;
1641   }
1642
1643   SetDescription(description.Data());
1644   AliHLTTriggerDecision decision((firedTriggers) ? kTRUE : kFALSE,
1645                                  GetTriggerName(),
1646                                  GetReadoutList(),
1647                                  GetDescription()
1648                                  );
1649   TriggerEvent(&decision, kAliHLTDataTypeTObject | kAliHLTDataOriginOut);
1650
1651   if (firedTriggers){
1652     LogInspect("TRD HLT trigger fired for event: description: >%s<, flags: 0x%04x",
1653                description.Data(), firedTriggers);
1654 #ifdef __TRDHLTDEBUG
1655     if (fEventRendering)
1656       RenderEvent();
1657 #endif
1658   } else {
1659     LogDebug("TRD HLT trigger did not fire for event");
1660   }
1661
1662   if (fPushHistos){
1663     PushBack(fHistArray, (kAliHLTDataTypeTObjArray | kAliHLTDataOriginTRD), 0x3fffff);
1664   }
1665
1666   LogDebug("### END   DoTrigger [event id: %llu, %d blocks, size: %d]",
1667                       fEventId, hltEventData->fBlockCnt, hltEventData->fStructSize);
1668
1669   return iResult;
1670 }
1671
1672 Bool_t AliHLTTRDTriggerComponent::TrackPlaneIntersect(AliExternalTrackParam *trk, Double_t pnt[3], Double_t norm[3], Double_t mag){
1673
1674   UInt_t its = 0;
1675   Double_t r = 290.;
1676   Double_t dist = 99999, dist_prev = 99999;
1677   Double_t x[3] = {0., 0., 0.};
1678   Bool_t ret = kTRUE;
1679
1680   dist = (x[0] - pnt[0]) * norm[0] + (x[1] - pnt[1]) *norm[1] + (x[2] - pnt[2]) * norm[2];
1681
1682   while(TMath::Abs(dist) > 0.1) {
1683
1684     trk->GetXYZAt(r, mag, x);
1685
1686     if ((x[0] * x[0] + x[1] * x[1]) < 100.) {  // extrapolation to radius failed
1687       ret = kFALSE;
1688       break;
1689     }
1690
1691     //distance between current track position and plane
1692     dist_prev = TMath::Abs(dist);
1693     dist = (x[0] - pnt[0]) * norm[0] + (x[1] - pnt[1]) * norm[1];
1694     r -= dist;
1695     its++;
1696     if(TMath::Abs(dist) >= dist_prev ||
1697        (r > 380.) || (r < 100.)){
1698       ret = kFALSE;
1699       break;
1700     }
1701   }
1702
1703   for (Int_t i=0; i<3; i++){
1704     if(ret)
1705       pnt[i] = x[i];
1706     else
1707       pnt[i] = 0.;
1708   }
1709
1710   return kTRUE;
1711 }
1712
1713
1714 Double_t AliHLTTRDTriggerComponent::RateTrackMatch(Double_t distY, Double_t distZ, Double_t rpt, Double_t gpt){
1715
1716   // maximum limits for spatial distance
1717   if ((distY > 5.) || (distZ > 20.))
1718     return 0.;
1719
1720   // same pt sign required
1721   if ((rpt * gpt) < 0.)
1722     return 0.;
1723
1724   Double_t rating_distY = -0.1 * distY + 1.;
1725   Double_t rating_distZ = -0.025 * distZ + 1.;
1726   Double_t rating_ptDiff = 1. - TMath::Abs((TMath::Abs(rpt) > 0.000001) ? ((gpt-rpt)/rpt) : 0.);
1727
1728   if (rating_ptDiff <  0.)
1729     rating_ptDiff = 0.2;
1730
1731   Double_t total = rating_distY * rating_distZ * rating_ptDiff;
1732
1733 //  DbgLog("", Form("#matching: rating:   dy: %.3f   dz: %.3f   dpt: %.3f     -> total: %.3f",
1734 //                  rating_distY, rating_distZ, rating_ptDiff, total));
1735
1736   if (total > 1.)
1737     LogError("track match rating exceeds limit of 1.0: %.3f", total);
1738
1739   return total;
1740 }
1741
1742 Int_t AliHLTTRDTriggerComponent::EstimateTrackDistance(AliESDtrack *esd_track,
1743                                                         const UShort_t stack,
1744                                                         const UShort_t layerMask,
1745                                                         const Float_t trklLocalY[6], const Int_t trklBinZ[6],
1746                                                         Double_t mag, Double_t *ydist, Double_t *zdist){
1747   if (!esd_track)
1748     return -3;
1749
1750   AliExternalTrackParam* refParam = NULL;
1751   if (esd_track->GetOuterParam())
1752     refParam = new AliExternalTrackParam(*(esd_track->GetOuterParam()));
1753   if (!refParam)
1754     refParam = new AliExternalTrackParam(*(esd_track));
1755
1756   Int_t res = EstimateTrackDistance(refParam, stack, layerMask, trklLocalY, trklBinZ, mag, ydist, zdist);
1757
1758   if (refParam)
1759     delete refParam;
1760
1761   return res;
1762
1763 }
1764
1765 Int_t AliHLTTRDTriggerComponent::EstimateTrackDistance(AliExternalTrackParam *refParam,
1766                                                         const UShort_t stack,
1767                                                         const UShort_t layerMask,
1768                                                         const Float_t trklLocalY[6], const Int_t trklBinZ[6],
1769                                                         Double_t mag, Double_t *ydist, Double_t *zdist){
1770
1771   Float_t diff_y = 0;
1772   Float_t diff_z = 0;
1773   Int_t nLayers = 0;
1774   Double_t xtrkl[3];
1775   Double_t ptrkl[3];
1776   Double_t ptrkl2[3];
1777   UInt_t trklDet;
1778   UShort_t trklLayer;
1779   UInt_t stack_gtu;
1780   UShort_t stackInSector;
1781
1782   AliTRDpadPlane* padPlane;
1783
1784   for (UShort_t iLayer = 0; iLayer < 6; iLayer++){
1785     if ((layerMask >> iLayer) & 1){
1786       trklDet = stack*6 + iLayer;
1787       trklLayer = iLayer;
1788       stack_gtu = stack;
1789       stackInSector = stack % 5;
1790
1791       // local coordinates of the outer end point of the tracklet
1792       xtrkl[0] = AliTRDgeometry::AnodePos();
1793       xtrkl[1] = trklLocalY[iLayer];
1794
1795       padPlane = fTRDGeometry->GetPadPlane(trklLayer, stackInSector);
1796       if(stackInSector == 2){ // corrected version by Felix Muecke
1797         xtrkl[2] = padPlane->GetRowPos(trklBinZ[iLayer]) - (padPlane->GetRowSize(trklBinZ[iLayer]))/2. - padPlane->GetRowPos(6);
1798       } else {
1799         xtrkl[2] = padPlane->GetRowPos(trklBinZ[iLayer]) - (padPlane->GetRowSize(trklBinZ[iLayer]))/2. - padPlane->GetRowPos(8);
1800       }
1801
1802       // transform to global coordinates
1803       TGeoHMatrix *matrix = fTRDGeometry->GetClusterMatrix(trklDet);
1804       if (!matrix){
1805         LogError("invalid TRD cluster matrix in EstimateTrackDistance for detector %i", trklDet);
1806         return -5;
1807       }
1808       matrix->LocalToMaster(xtrkl, ptrkl);
1809       fTRDGeometry->RotateBack((stack/5) * 30, ptrkl, ptrkl2);  // ptrkl2 now contains the global position of the outer end point of the tracklet
1810
1811       // calculate parameterization of plane representing the tracklets layer
1812       Double_t layer_zero_local[3] = {0., 0.,   0.};
1813       Double_t layer_zero_global[3], layer_zero_global2[3];
1814
1815       matrix->LocalToMaster(layer_zero_local, layer_zero_global);
1816       fTRDGeometry->RotateBack(trklDet, layer_zero_global, layer_zero_global2); // layer_zero_global2 points to chamber origin in global coords
1817
1818       Double_t layer_ref_local[3] = {AliTRDgeometry::AnodePos(), 0.,  0.};
1819       Double_t layer_ref_global[3], layer_ref_global2[3];
1820
1821       matrix->LocalToMaster(layer_ref_local, layer_ref_global);
1822       fTRDGeometry->RotateBack(trklDet, layer_ref_global, layer_ref_global2); // layer_ref_global2 points to center anode pos within plane in global coords
1823
1824       Double_t n0[3] = {layer_ref_global2[0]-layer_zero_global2[0],
1825                         layer_ref_global2[1]-layer_zero_global2[1],
1826                         layer_ref_global2[2]-layer_zero_global2[2]};
1827
1828       Double_t n_len = TMath::Sqrt(n0[0]*n0[0] + n0[1]*n0[1] + n0[2]*n0[2]);
1829       if (n_len == 0.){ // This should never happen
1830         //printf("<ERROR> divison by zero in estimate_track_distance!");
1831         n_len = 1.;
1832       }
1833       Double_t n[3] = {n0[0]/n_len, n0[1]/n_len, n0[2]/n_len}; // normal vector of plane
1834
1835       Bool_t isects = TrackPlaneIntersect(refParam, layer_ref_global2, n, mag); // find intersection point between track and TRD layer
1836
1837       if (isects == kFALSE){ // extrapolation fails, because track never reaches the TRD radius
1838         return -1;
1839       }
1840
1841       Double_t m[2] = {ptrkl2[0] - layer_ref_global2[0], ptrkl2[1] - layer_ref_global2[1]};
1842       Double_t len_m = TMath::Sqrt(m[0]*m[0] + m[1]*m[1]);
1843       diff_y += len_m;
1844       diff_z += TMath::Abs(ptrkl2[2] - layer_ref_global2[2]);
1845       nLayers++;
1846     }
1847   }
1848
1849   if (nLayers > 0){
1850     *ydist = diff_y / nLayers;
1851     *zdist = diff_z / nLayers;
1852     return 0;
1853   } else {
1854     LogError("invalid number of contributing layers (%d) in EstimateTrackDistance()", nLayers);
1855     return -4;
1856   }
1857 }