TRD/AliTRDCalibTask.cxx

   1 \r
   2 /**************************************************************************\r
   3  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
   4  *                                                                        *\r
   5  * Author: The ALICE Off-line Project.                                    *\r
   6  * Contributors are mentioned in the code where appropriate.              *\r
   7  *                                                                        *\r
   8  * Permission to use, copy, modify and distribute this software and its   *\r
   9  * documentation strictly for non-commercial purposes is hereby granted   *\r
  10  * without fee, provided that the above copyright notice appears in all   *\r
  11  * copies and that both the copyright notice and this permission notice   *\r
  12  * appear in the supporting documentation. The authors make no claims     *\r
  13  * about the suitability of this software for any purpose. It is          *\r
  14  * provided "as is" without express or implied warranty.                  *\r
  15  **************************************************************************/\r
  16 \r
  17 /////////////////////////////////////////////////////////////////////////////////\r
  18 //                                                                             \r
  19 // AliTRDCalibTask                                                               \r
  20 //                                                                             \r
  21 // Offline TRD calibration task                                \r
  22 //                        \r
  23 // Author:\r
  24 //   R. Bailhache (R.Bailhache@gsi.de)\r
  25 //                            \r
  26 //////////////////////////////////////////////////////////////////////////////////////\r
  27 \r
  28 \r
  29 \r
  30 #include "Riostream.h"\r
  31 #include "TChain.h"\r
  32 #include "TTree.h"\r
  33 #include "TFile.h"\r
  34 #include "TProfile2D.h"\r
  35 #include "TH2I.h"\r
  36 #include "TH1F.h"\r
  37 #include "TH2F.h"\r
  38 #include "TList.h"\r
  39 #include "TMath.h"\r
  40 #include "TObject.h"\r
  41 #include "TObjArray.h"\r
  42 #include "TString.h"\r
  43 #include "TCanvas.h"\r
  44 #include "TLegend.h"\r
  45 #include "TStyle.h"\r
  46 #include "TLine.h"\r
  47 \r
  48 #include "AliAnalysisTask.h"\r
  49 #include "AliAnalysisManager.h"\r
  50 \r
  51 #include "AliESDVertex.h"\r
  52 #include "AliESDEvent.h"\r
  53 #include "AliESDfriend.h"\r
  54 #include "AliESDInputHandler.h"\r
  55 #include "AliESDtrack.h"\r
  56 #include "AliESDfriendTrack.h"\r
  57 #include "AliTRDtrackV1.h"\r
  58 #include "AliTRDseedV1.h"\r
  59 #include "AliTRDcluster.h"\r
  60 #include "AliTRDgeometry.h"\r
  61 #include "AliESDtrackCuts.h"\r
  62 #include "AliESDVertex.h"\r
  63 #include "AliTRDCalDet.h"\r
  64 \r
  65 #include "AliTRDCalibraFillHisto.h"\r
  66 #include "AliTRDCalibraVdriftLinearFit.h" \r
  67 \r
  68 #include "AliTRDcalibDB.h"\r
  69 \r
  70 \r
  71 #include "AliTRDCalibTask.h"\r
  72 \r
  73 \r
  74 ClassImp(AliTRDCalibTask)\r
  75 \r
  76 //________________________________________________________________________\r
  77   AliTRDCalibTask::AliTRDCalibTask(const char *name) \r
  78     : AliAnalysisTask(name, ""), fESD(0), fESDfriend(0),\r
  79       fkEsdTrack(0),\r
  80       fFriendTrack(0),\r
  81       fCalibObject(0),\r
  82       fTrdTrack(0),\r
  83       fCl(0),\r
  84       fListHist(0),\r
  85       fTRDCalibraFillHisto(0),\r
  86       fNEvents(0),\r
  87       fNbTRDTrack(0),\r
  88       fNbTRDTrackOffline(0),\r
  89       fNbTRDTrackStandalone(0),\r
  90       fNbTPCTRDtrack(0),\r
  91       fNbTimeBin(0),\r
  92       fNbTimeBinOffline(0),\r
  93       fNbTimeBinStandalone(0),\r
  94       fNbClusters(0),\r
  95       fNbClustersOffline(0),\r
  96       fNbClustersStandalone(0),\r
  97       fNbTracklets(0),\r
  98       fNbTrackletsOffline(0),\r
  99       fNbTrackletsStandalone(0),\r
 100       fCH2dSum(0),\r
 101       fPH2dSum(0),\r
 102       fCH2dSM(0),\r
 103       fPH2dSM(0),\r
 104       fHisto2d(kTRUE),\r
 105       fVector2d(kFALSE),\r
 106       fVdriftLinear(kTRUE),\r
 107       fNbTimeBins(0),\r
 108       fSelectedTrigger(new TObjArray()),\r
 109       fRejected(kTRUE),\r
 110       fEsdTrackCuts(0),\r
 111       fRequirePrimaryVertex(kFALSE),\r
 112       fMinNbContributors(0),\r
 113       fLow(0),\r
 114       fHigh(30),\r
 115       fFillZero(kFALSE),\r
 116       fNormalizeNbOfCluster(kFALSE),\r
 117       fRelativeScale(0.0),\r
 118       fMaxCluster(100.0),\r
 119       fNbMaxCluster(2),\r
 120       fOfflineTracks(kFALSE),\r
 121       fStandaloneTracks(kFALSE),\r
 122       fCalDetGain(0x0),\r
 123       fMaxEvent(0),\r
 124       fCounter(0),\r
 125       fDebug(0)\r
 126 {\r
 127   //\r
 128   // Default constructor\r
 129   //\r
 130 \r
 131   fNz[0] = 0;\r
 132   fNz[1] = 0;\r
 133   fNz[2] = 0;\r
 134   \r
 135   fNrphi[0] = 0;\r
 136   fNrphi[1] = 0;\r
 137   fNrphi[2] = 0;\r
 138 \r
 139   // Define input and output slots here\r
 140   // Input slot #0 works with a TChain\r
 141   DefineInput(0, TChain::Class());\r
 142     \r
 143   // Output slot #0 writes into a TList container\r
 144   DefineOutput(0, TList::Class());\r
 145   \r
 146    \r
 147 }\r
 148 //____________________________________________________________________________________\r
 149 AliTRDCalibTask::~AliTRDCalibTask()\r
 150 {\r
 151   //\r
 152   // AliTRDCalibTask destructor\r
 153   //\r
 154 \r
 155   // Pointeur\r
 156   if(fNEvents) delete fNEvents;\r
 157   if(fNbTRDTrack) delete fNbTRDTrack;\r
 158   if(fNbTRDTrackOffline) delete fNbTRDTrackOffline;\r
 159   if(fNbTRDTrackStandalone) delete fNbTRDTrackStandalone;\r
 160   if(fNbTPCTRDtrack) delete fNbTPCTRDtrack;\r
 161   if(fNbTimeBin) delete fNbTimeBin;\r
 162   if(fNbTimeBinOffline) delete fNbTimeBinOffline;\r
 163   if(fNbTimeBinStandalone) delete fNbTimeBinStandalone;\r
 164   if(fNbClusters) delete fNbClusters;\r
 165   if(fNbClustersOffline) delete fNbClustersOffline;\r
 166   if(fNbClustersStandalone) delete fNbClustersStandalone;\r
 167   if(fNbTracklets) delete fNbTracklets;\r
 168   if(fNbTrackletsOffline) delete fNbTrackletsOffline;\r
 169   if(fNbTrackletsStandalone) delete fNbTrackletsStandalone;\r
 170   if(fCH2dSum) delete fCH2dSum;\r
 171   if(fPH2dSum) delete fPH2dSum;\r
 172   if(fCH2dSM) delete fCH2dSM;\r
 173   if(fPH2dSM) delete fPH2dSM;\r
 174   if(fCalDetGain) delete fCalDetGain;\r
 175   \r
 176   if(fSelectedTrigger) {\r
 177     fSelectedTrigger->Delete();\r
 178     delete fSelectedTrigger;\r
 179   }\r
 180   if(fEsdTrackCuts) {\r
 181     delete fEsdTrackCuts;\r
 182   }\r
 183   \r
 184 }\r
 185 //________________________________________________________________________\r
 186 void AliTRDCalibTask::ConnectInputData(Option_t *) \r
 187 {\r
 188   // Connect ESD or AOD here\r
 189   // Called once\r
 190   \r
 191   TTree* tree = dynamic_cast<TTree*> (GetInputData(0)); //pointer wird "umgecastet" auf anderen Variablentyp\r
 192   if (!tree) {\r
 193     Printf("ERROR: Could not read chain from input slot 0");\r
 194   } else {\r
 195     \r
 196     AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
 197     \r
 198     if (!esdH) {\r
 199       Printf("ERROR: Could not get ESDInputHandler");\r
 200     } else {\r
 201       fESD = esdH->GetEvent();\r
 202       esdH->SetActiveBranches("ESDfriend*");\r
 203       Printf("*** CONNECTED NEW EVENT ****");\r
 204     }\r
 205     \r
 206   }\r
 207 }\r
 208 //________________________________________________________________________\r
 209 void AliTRDCalibTask::CreateOutputObjects() \r
 210 {\r
 211   //\r
 212   // Create the histos\r
 213   //\r
 214 \r
 215   // Number of time bins\r
 216   if(fNbTimeBins==0) {\r
 217     AliTRDcalibDB *cal = AliTRDcalibDB::Instance();\r
 218     fNbTimeBins = cal->GetNumberOfTimeBinsDCS();\r
 219     if(fNbTimeBins <= 0){ \r
 220       AliWarning(Form("No of TimeBins from DB [%d] use default [30]", fNbTimeBins));\r
 221       fNbTimeBins = 30;\r
 222     }\r
 223   }\r
 224   \r
 225   // instance calibration \r
 226   fTRDCalibraFillHisto = AliTRDCalibraFillHisto::Instance();\r
 227   fTRDCalibraFillHisto->SetHisto2d(fHisto2d); // choose to use histograms\r
 228   fTRDCalibraFillHisto->SetVector2d(fVector2d); // choose to use vectors\r
 229   fTRDCalibraFillHisto->SetCH2dOn();  // choose to calibrate the gain\r
 230   fTRDCalibraFillHisto->SetPH2dOn();  // choose to calibrate the drift velocity\r
 231   fTRDCalibraFillHisto->SetPRF2dOn(); // choose to look at the PRF\r
 232   fTRDCalibraFillHisto->SetLinearFitterOn(fVdriftLinear); // Other possibility vdrift VDRIFT\r
 233   fTRDCalibraFillHisto->SetLinearFitterDebugOn(fVdriftLinear); // Other possibility vdrift\r
 234   for(Int_t k = 0; k < 3; k++){\r
 235     if(((fNz[k] != 10) && (fNrphi[k] != 10)) && ((fNz[k] != 100) && (fNrphi[k] != 100))) {\r
 236       fTRDCalibraFillHisto->SetNz(k,fNz[k]);                                    // Mode calibration\r
 237       fTRDCalibraFillHisto->SetNrphi(k,fNrphi[k]);                             // Mode calibration\r
 238     }\r
 239     else {\r
 240       if((fNz[k] == 100) && (fNrphi[k] == 100))  {\r
 241         if(fVector2d) AliInfo("The mode all together is not supported by the vector method");\r
 242         fTRDCalibraFillHisto->SetAllTogether(k);\r
 243       }\r
 244       if((fNz[k] == 10) && (fNrphi[k] == 10))  {\r
 245         if(fVector2d) AliInfo("The mode per supermodule is not supported by the vector method");\r
 246         fTRDCalibraFillHisto->SetPerSuperModule(k);\r
 247       }\r
 248     }\r
 249   }\r
 250   // Variables for how to fill\r
 251   fTRDCalibraFillHisto->SetFillWithZero(fFillZero);\r
 252   fTRDCalibraFillHisto->SetNormalizeNbOfCluster(fNormalizeNbOfCluster); \r
 253   fTRDCalibraFillHisto->SetMaxCluster(fMaxCluster);\r
 254   fTRDCalibraFillHisto->SetNbMaxCluster(fNbMaxCluster);\r
 255   \r
 256   // Init with 30 timebins\r
 257   fTRDCalibraFillHisto->Init2Dhistos(fNbTimeBins); // initialise the histos\r
 258   fTRDCalibraFillHisto->SetNumberClusters(fLow); // At least 11 clusters\r
 259   fTRDCalibraFillHisto->SetNumberClustersf(fHigh); // At least 11 clusters\r
 260   fRelativeScale = fTRDCalibraFillHisto->GetRelativeScale(); // Get the relative scale for the gain\r
 261   \r
 262   // For testing only\r
 263   if(fDebug > 2) fTRDCalibraFillHisto->SetDebugLevel(1); //debug stuff\r
 264   \r
 265   // output list\r
 266   fListHist = new TList();\r
 267   if(fHisto2d) {  \r
 268     fListHist->Add(fTRDCalibraFillHisto->GetCH2d());\r
 269     fListHist->Add(fTRDCalibraFillHisto->GetPH2d()); \r
 270     fListHist->Add(fTRDCalibraFillHisto->GetPRF2d());\r
 271   } \r
 272   if(fVdriftLinear) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetVdriftLinearFit());\r
 273   if(fVector2d) fListHist->Add((TObject *) fTRDCalibraFillHisto->GetCalibraVector()); //calibra vector  \r
 274   fNEvents = new TH1I("NEvents","NEvents", 2, 0, 2);\r
 275   fListHist->Add(fNEvents);\r
 276   \r
 277   /////////////////////////////////////////\r
 278   // First debug level\r
 279   ///////////////////////////////////////\r
 280   if(fDebug > 0) {\r
 281     \r
 282     fPH2dSM = new TProfile2D("PH2dSM","Nz10Nrphi10"\r
 283                             ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)\r
 284                            ,18,0,18);\r
 285     fPH2dSM->SetYTitle("Det/pad groups");\r
 286     fPH2dSM->SetXTitle("time [#mus]");\r
 287     fPH2dSM->SetZTitle("<PH> [a.u.]");\r
 288     fPH2dSM->SetStats(0);\r
 289     //\r
 290     fCH2dSM = new TH2I("CH2dSM","Nz10Nrphi10",50,0,300,18,0,18);\r
 291     fCH2dSM->SetYTitle("Det/pad groups");\r
 292     fCH2dSM->SetXTitle("charge deposit [a.u]");\r
 293     fCH2dSM->SetZTitle("counts");\r
 294     fCH2dSM->SetStats(0);\r
 295     fCH2dSM->Sumw2();\r
 296     //\r
 297     fPH2dSum = new TProfile2D("PH2dSum","Nz100Nrphi100"\r
 298                             ,fNbTimeBins,-0.05,(Double_t)((fNbTimeBins-0.5)/10.0)\r
 299                             ,1,0,1);\r
 300     fPH2dSum->SetYTitle("Det/pad groups");\r
 301     fPH2dSum->SetXTitle("time [#mus]");\r
 302     fPH2dSum->SetZTitle("<PH> [a.u.]");\r
 303     fPH2dSum->SetStats(0);\r
 304     //\r
 305     fCH2dSum = new TH2I("CH2dSum","Nz100Nrphi100",50,0,300,1,0,1);\r
 306     fCH2dSum->SetYTitle("Det/pad groups");\r
 307     fCH2dSum->SetXTitle("charge deposit [a.u]");\r
 308     fCH2dSum->SetZTitle("counts");\r
 309     fCH2dSum->SetStats(0);\r
 310     fCH2dSum->Sumw2();\r
 311     \r
 312     // Add them\r
 313     fListHist->Add(fPH2dSM);\r
 314     fListHist->Add(fCH2dSM);\r
 315     fListHist->Add(fPH2dSum);\r
 316     fListHist->Add(fCH2dSum);\r
 317   }\r
 318 \r
 319   /////////////////////////////////////////\r
 320   // Second debug level\r
 321   ///////////////////////////////////////\r
 322   if(fDebug > 1) {\r
 323 \r
 324     fNbTRDTrack = new TH1F("TRDTrack","TRDTrack",50,0,50);\r
 325     fNbTRDTrack->Sumw2();\r
 326     fNbTRDTrackOffline = new TH1F("TRDTrackOffline","TRDTrackOffline",50,0,50);\r
 327     fNbTRDTrackOffline->Sumw2();\r
 328     fNbTRDTrackStandalone = new TH1F("TRDTrackStandalone","TRDTrackStandalone",50,0,50);\r
 329     fNbTRDTrackStandalone->Sumw2();\r
 330     fNbTPCTRDtrack = new TH2F("NbTPCTRDtrack","NbTPCTRDtrack",100,0,100,100,0,100);\r
 331     fNbTPCTRDtrack->Sumw2();\r
 332     //\r
 333     fNbTimeBin = new TH1F("NbTimeBin","NbTimeBin",35,0,35);\r
 334     fNbTimeBin->Sumw2();\r
 335     fNbTimeBinOffline = new TH1F("NbTimeBinOffline","NbTimeBinOffline",35,0,35);\r
 336     fNbTimeBinOffline->Sumw2();\r
 337     fNbTimeBinStandalone = new TH1F("NbTimeBinStandalone","NbTimeBinStandalone",35,0,35);\r
 338     fNbTimeBinStandalone->Sumw2();\r
 339     //\r
 340     fNbClusters = new TH1F("NbClusters","",35,0,35);\r
 341     fNbClusters->Sumw2();\r
 342     fNbClustersOffline = new TH1F("NbClustersOffline","",35,0,35);\r
 343     fNbClustersOffline->Sumw2();\r
 344     fNbClustersStandalone = new TH1F("NbClustersStandalone","",35,0,35);\r
 345     fNbClustersStandalone->Sumw2();\r
 346     //\r
 347     fNbTracklets = new TH1F("NbTracklets","NbTracklets",540,0.,540.);\r
 348     fNbTracklets->Sumw2();\r
 349     fNbTrackletsOffline = new TH1F("NbTrackletsOffline","NbTrackletsOffline",540,0.,540.);\r
 350     fNbTrackletsOffline->Sumw2();\r
 351     fNbTrackletsStandalone = new TH1F("NbTrackletsStandalone","NbTrackletsStandalone",540,0.,540.);\r
 352     fNbTrackletsStandalone->Sumw2();\r
 353    \r
 354     fListHist->Add(fNbTRDTrack);\r
 355     fListHist->Add(fNbTRDTrackOffline);\r
 356     fListHist->Add(fNbTRDTrackStandalone);\r
 357     fListHist->Add(fNbTPCTRDtrack);\r
 358     \r
 359     fListHist->Add(fNbTimeBin);\r
 360     fListHist->Add(fNbTimeBinOffline);\r
 361     fListHist->Add(fNbTimeBinStandalone);\r
 362     fListHist->Add(fNbClusters);\r
 363     fListHist->Add(fNbClustersOffline);\r
 364     fListHist->Add(fNbClustersStandalone);\r
 365     fListHist->Add(fNbTracklets);\r
 366     fListHist->Add(fNbTrackletsOffline);\r
 367     fListHist->Add(fNbTrackletsStandalone);\r
 368     \r
 369   }\r
 370   \r
 371  }\r
 372  //________________________________________________________________________\r
 373  void AliTRDCalibTask::Exec(Option_t *) \r
 374  {\r
 375    //\r
 376    // Filling of the histos\r
 377    //\r
 378 \r
 379    AliLog::SetGlobalLogLevel(AliLog::kError);\r
 380 \r
 381    if (!fESD) {\r
 382      Printf("ERROR: fESD not available");\r
 383      return;\r
 384    }\r
 385 \r
 386    fCounter++;\r
 387    if((fMaxEvent != 0) && (fMaxEvent < fCounter)) return;\r
 388 \r
 389    ///////////////////////////////\r
 390    // Require a primary vertex\r
 391    //////////////////////////////\r
 392    if(fRequirePrimaryVertex) {\r
 393      //printf("test\n");\r
 394      const AliESDVertex* vtxESD = 0x0;\r
 395      vtxESD = fESD->GetPrimaryVertexTPC() ;\r
 396      if(!vtxESD){\r
 397        //printf("No primary at all vertex\n");\r
 398        return;\r
 399      }\r
 400      Int_t nCtrb = vtxESD->GetNContributors();\r
 401      //printf("Number of contributors %d\n",nCtrb);\r
 402      if(nCtrb < fMinNbContributors) return;\r
 403    }\r
 404 \r
 405    //printf("Primary vertex passed\n");\r
 406 \r
 407    ///////////////////\r
 408    // Check trigger\r
 409    ///////////////////\r
 410    //printf("Class Fired %s\n",(const char*)fESD->GetFiredTriggerClasses());\r
 411    Bool_t pass = kTRUE;\r
 412    Int_t numberOfTriggerSelected = fSelectedTrigger->GetEntriesFast();\r
 413    //printf("numberofTriggerSelected %d\n",numberOfTriggerSelected);\r
 414    if(fRejected) {\r
 415      pass = kTRUE;\r
 416      for(Int_t k = 0; k < numberOfTriggerSelected; k++){\r
 417        const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);\r
 418        const TString tString=obString->GetString();\r
 419        if(fESD->IsTriggerClassFired((const char*)tString)) {\r
 420          pass = kFALSE;\r
 421        }\r
 422      }\r
 423    }\r
 424    else {\r
 425      pass = kFALSE;\r
 426      for(Int_t k = 0; k < numberOfTriggerSelected; k++){\r
 427        const TObjString *const obString=(TObjString*)fSelectedTrigger->At(k);\r
 428        const TString tString=obString->GetString();\r
 429        if(fESD->IsTriggerClassFired((const char*)tString)) {\r
 430          pass = kTRUE;\r
 431        }\r
 432      }\r
 433    }\r
 434    if(!pass) return;   \r
 435 \r
 436    //printf("Trigger passed\n");\r
 437   \r
 438    fNEvents->Fill(1);\r
 439 \r
 440    // In total\r
 441    Int_t nbTrdTracks = 0;\r
 442    // standalone\r
 443    Int_t nbTrdTracksStandalone = 0;\r
 444    // offline\r
 445    Int_t nbTrdTracksOffline = 0;\r
 446    // TPC\r
 447    Int_t nbtrackTPC = 0;\r
 448 \r
 449    Double_t nbTracks = fESD->GetNumberOfTracks();\r
 450    //printf("Number of tracks %f\n",nbTracks);  \r
 451 \r
 452    if (nbTracks <= 0.0) {\r
 453      \r
 454      if(fDebug > 1) {\r
 455        fNbTRDTrack->Fill(nbTrdTracks);\r
 456        fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);\r
 457        fNbTRDTrackOffline->Fill(nbTrdTracksOffline);\r
 458      }\r
 459      PostData(0, fListHist);\r
 460      return;\r
 461    }\r
 462 \r
 463    fESDfriend = (AliESDfriend *)fESD->FindListObject("AliESDfriend");\r
 464    fESD->SetESDfriend(fESDfriend);\r
 465    if(!fESDfriend) return;   \r
 466    \r
 467    //printf("has friends\n");\r
 468 \r
 469    ////////////////////////////////////\r
 470    // Check the number of TPC tracks\r
 471    ///////////////////////////////////\r
 472    for(Int_t itrk = 0; itrk < nbTracks; itrk++){\r
 473      // Get ESD track\r
 474      fkEsdTrack = fESD->GetTrack(itrk);\r
 475      ULong_t status = fkEsdTrack->GetStatus(); \r
 476      if(status&(AliESDtrack::kTPCout)) nbtrackTPC++;\r
 477      // Check that the calibration object is here\r
 478      if(fESDfriend && (fESDfriend->GetTrack(itrk))) {\r
 479        fFriendTrack = new AliESDfriendTrack(*(fESDfriend->GetTrack(itrk)));\r
 480        Int_t counteer = 0;\r
 481        Int_t icalib=0;\r
 482        while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){\r
 483          if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;\r
 484          counteer++;\r
 485        }\r
 486        if(counteer > 0) {\r
 487          nbTrdTracks++;\r
 488          if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {\r
 489            nbTrdTracksStandalone++;\r
 490          }\r
 491          if((status&(AliESDtrack::kTRDin))) {\r
 492            nbTrdTracksOffline++;\r
 493          }\r
 494        }\r
 495        if(fFriendTrack) delete fFriendTrack;\r
 496      }\r
 497    }\r
 498    //printf("Number of TPC tracks %d, TRD %d\n",nbtrackTPC,nbTrdTracks);\r
 499    if((nbtrackTPC>0) && (nbTrdTracks > (3.0*nbtrackTPC))) pass = kFALSE;\r
 500    \r
 501    if(fDebug > 1) {\r
 502      \r
 503      fNbTRDTrack->Fill(nbTrdTracks);\r
 504      fNbTRDTrackStandalone->Fill(nbTrdTracksStandalone);\r
 505      fNbTRDTrackOffline->Fill(nbTrdTracksOffline);\r
 506      fNbTPCTRDtrack->Fill(nbTrdTracks,nbtrackTPC);\r
 507    \r
 508    }\r
 509 \r
 510    if(!pass) {\r
 511      PostData(0, fListHist);\r
 512      return;\r
 513    }\r
 514   \r
 515    /////////////////////////////////////\r
 516    // Loop on AliESDtrack\r
 517    ////////////////////////////////////\r
 518       \r
 519    for(int itrk=0; itrk < nbTracks; itrk++){\r
 520 \r
 521      // Get ESD track\r
 522      fkEsdTrack = fESD->GetTrack(itrk);\r
 523      if(!fkEsdTrack) continue;\r
 524 \r
 525      // Quality cuts on the AliESDtrack\r
 526      if((fEsdTrackCuts) && (!fEsdTrackCuts->IsSelected((AliVParticle *)fkEsdTrack))) {\r
 527        //printf("Not a good track\n");\r
 528        continue;\r
 529      }\r
 530      \r
 531      // Other cuts\r
 532      Bool_t good = kTRUE;\r
 533      Bool_t standalonetrack = kFALSE;\r
 534      Bool_t offlinetrack = kFALSE;\r
 535      ULong_t status = fkEsdTrack->GetStatus();\r
 536      \r
 537      if(!(fESDfriend->GetTrack(itrk)))  continue;   \r
 538      \r
 539      fFriendTrack = new AliESDfriendTrack(*(fESDfriend->GetTrack(itrk)));\r
 540      \r
 541      //////////////////////////////////////\r
 542      // Loop on calibration objects\r
 543      //////////////////////////////////////\r
 544      Int_t icalib=0;\r
 545      while((fCalibObject = (TObject *)(fFriendTrack->GetCalibObject(icalib++)))){\r
 546        if(strcmp(fCalibObject->IsA()->GetName(), "AliTRDtrackV1") != 0) continue;\r
 547 \r
 548        if((status&(AliESDtrack::kTRDout)) && (!(status&(AliESDtrack::kTRDin)))) {\r
 549          standalonetrack = kTRUE;\r
 550        }\r
 551        if((status&(AliESDtrack::kTRDin))) {\r
 552          offlinetrack = kTRUE;\r
 553        }\r
 554        if(fOfflineTracks){\r
 555          if(!offlinetrack){\r
 556            good = kFALSE;\r
 557          }\r
 558        }\r
 559        else if(fStandaloneTracks){\r
 560          if(!standalonetrack){\r
 561            good = kFALSE;\r
 562          }\r
 563        }\r
 564        \r
 565        fTrdTrack = (AliTRDtrackV1 *)fCalibObject;\r
 566        if(good) fTRDCalibraFillHisto->UpdateHistogramsV1(fTrdTrack);\r
 567        \r
 568        //////////////////////////////////\r
 569        // Debug \r
 570        ////////////////////////////////\r
 571 \r
 572        if(fDebug > 0) {\r
 573          \r
 574          \r
 575          Int_t nbtracklets = 0;\r
 576          \r
 577          // Check some stuff\r
 578          Bool_t standalonetracklet = kFALSE;  \r
 579          const AliTRDseedV1 *tracklet = 0x0;\r
 580          //////////////////////////////////////\r
 581          // Loop tracklets\r
 582          ///////////////////////////////////// \r
 583          for(Int_t itr = 0; itr < 6; itr++){\r
 584            \r
 585            if(!(tracklet = fTrdTrack->GetTracklet(itr))) continue;\r
 586            if(!tracklet->IsOK()) continue;\r
 587            nbtracklets++;\r
 588            standalonetracklet = kFALSE; \r
 589            if(tracklet->IsStandAlone()) standalonetracklet = kTRUE;\r
 590 \r
 591            Int_t nbclusters = 0;\r
 592            Double_t phtb[AliTRDseedV1::kNtb];\r
 593            memset(phtb, 0, AliTRDseedV1::kNtb*sizeof(Double_t));\r
 594            Double_t sum = 0.0;\r
 595            Float_t normalisation = 6.67;\r
 596            Int_t detector = 0;\r
 597            Int_t sector = 0;\r
 598            //Int_t crossrow = 0;\r
 599            \r
 600            // Check no shared clusters\r
 601            //for(int icc=AliTRDseedV1::kNtb; icc<AliTRDseedV1::kNclusters; icc++){\r
 602            //  if((fcl = tracklet->GetClusters(icc)))  crossrow = 1;\r
 603            // }\r
 604            \r
 605            // Loop on clusters\r
 606            for(int ic=0; ic<AliTRDseedV1::kNtb; ic++){\r
 607              \r
 608              if(!(fCl = tracklet->GetClusters(ic))) continue;\r
 609              nbclusters++;\r
 610              Int_t time = fCl->GetPadTime();\r
 611              Float_t ch =  tracklet->GetdQdl(ic);\r
 612              Float_t qcl = TMath::Abs(fCl->GetQ());\r
 613              detector = fCl->GetDetector();       \r
 614              if((time>-1) && (time<fNbTimeBins)) phtb[time]=qcl;\r
 615              if((fCalDetGain) && (fCalDetGain->GetValue(detector) > 0.0)) sum += ch*fCalDetGain->GetValue(detector)/normalisation;      \r
 616              else sum += ch/normalisation;\r
 617              \r
 618              if(fDebug > 1) {\r
 619                fNbTimeBin->Fill(time);\r
 620                if(tracklet->IsStandAlone()) fNbTimeBinStandalone->Fill(time);\r
 621                else fNbTimeBinOffline->Fill(time);\r
 622              }\r
 623            }\r
 624            sector = AliTRDgeometry::GetSector(detector);\r
 625            \r
 626            if(fDebug > 1) {\r
 627              fNbTracklets->Fill(detector);\r
 628              if(tracklet->IsStandAlone()) fNbTrackletsStandalone->Fill(detector);\r
 629              else fNbTrackletsOffline->Fill(detector);\r
 630            \r
 631              fNbClusters->Fill(nbclusters);\r
 632              if(tracklet->IsStandAlone())  fNbClustersStandalone->Fill(nbclusters);\r
 633              else  fNbClustersOffline->Fill(nbclusters);\r
 634            }       \r
 635 \r
 636            if(fDebug > 0) {\r
 637              if((nbclusters > fLow) && (nbclusters < fHigh)){\r
 638                if(fRelativeScale > 0.0) sum = sum/fRelativeScale;              \r
 639                fCH2dSM->Fill(sum,sector+0.5);\r
 640                fCH2dSum->Fill(sum,0.5);\r
 641                Bool_t checknoise = kTRUE;\r
 642                if(fMaxCluster > 0) {\r
 643                  if(phtb[0] > fMaxCluster) checknoise = kFALSE;\r
 644                  if(fNbTimeBins > fNbMaxCluster) {\r
 645                    for(Int_t k = (fNbTimeBins-fNbMaxCluster); k < fNbTimeBins; k++){\r
 646                      if(phtb[k] > fMaxCluster) checknoise = kFALSE;\r
 647                    }\r
 648                  }\r
 649                }\r
 650                if(checknoise) {        \r
 651                  for(int ic=0; ic<fNbTimeBins; ic++){\r
 652                    if(fFillZero) {\r
 653                      fPH2dSum->Fill((Double_t)(ic/10.0),0.5,(Double_t)phtb[ic]);\r
 654                      fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);\r
 655                    }\r
 656                    else {\r
 657                      if(phtb[ic] > 0.0) {\r
 658                        fPH2dSum->Fill((Double_t)(ic/10.0),0.0,(Double_t)phtb[ic]);\r
 659                        fPH2dSM->Fill((Double_t)(ic/10.0),sector+0.5,(Double_t)phtb[ic]);\r
 660                      }\r
 661                    }\r
 662                  }\r
 663                }\r
 664              }\r
 665            }\r
 666          } // loop on tracklets\r
 667     \r
 668        } // debug\r
 669        \r
 670      }// while calibration objects\r
 671      \r
 672      delete fFriendTrack;\r
 673      \r
 674    } // loop ESD track\r
 675    \r
 676    // Post output data\r
 677    PostData(0, fListHist);\r
 678  }     \r
 679 //________________________________________________________________________\r
 680 void AliTRDCalibTask::Terminate(Option_t *) \r
 681 {\r
 682   //\r
 683   // Terminate\r
 684   //\r
 685   \r
 686   if(fTRDCalibraFillHisto) fTRDCalibraFillHisto->DestroyDebugStreamer();\r
 687 \r
 688  \r
 689 }\r
 690 //_______________________________________________________\r
 691 Bool_t AliTRDCalibTask::Load(const Char_t *filename)\r
 692 {\r
 693   //\r
 694   // Generic container loader\r
 695   //\r
 696 \r
 697   if(!TFile::Open(filename)){\r
 698     //AliWarning(Form("Couldn't open file %s.", filename));\r
 699     return kFALSE;\r
 700   }\r
 701   TList *o = 0x0;\r
 702   if(!(o = (TList*)gFile->Get(GetName()))){\r
 703     //AliWarning("Missing histogram container.");\r
 704     return kFALSE;\r
 705   }\r
 706   fListHist = (TList*)o->Clone(GetName());\r
 707   gFile->Close();\r
 708   return kTRUE;\r
 709 }\r
 710 //________________________________________________________________________\r
 711 void AliTRDCalibTask::Plot() \r
 712 {\r
 713   //\r
 714   // Plot the histos stored in the TList\r
 715   //\r
 716  \r
 717   if(!fListHist) return;\r
 718 \r
 719   /////////////////////////////////////\r
 720   // Take the debug stuff\r
 721   /////////////////////////////////////\r
 722 \r
 723   TH1I *nEvents  = (TH1I *) fListHist->FindObject("NEvents");\r
 724 \r
 725   TH1F *trdTrack = (TH1F *) fListHist->FindObject("TRDTrack");\r
 726   TH1F *trdTrackOffline = (TH1F *) fListHist->FindObject("TRDTrackOffline");\r
 727   TH1F *trdTrackStandalone = (TH1F *) fListHist->FindObject("TRDTrackStandalone");\r
 728 \r
 729   TH2F *tpctrdTrack = (TH2F *) fListHist->FindObject("NbTPCTRDtrack");\r
 730 \r
 731   TH1F *nbTimeBin = (TH1F *) fListHist->FindObject("NbTimeBin");\r
 732   TH1F *nbTimeBinOffline = (TH1F *) fListHist->FindObject("NbTimeBinOffline");\r
 733   TH1F *nbTimeBinStandalone = (TH1F *) fListHist->FindObject("NbTimeBinStandalone");\r
 734 \r
 735   TH1F *nbClusters = (TH1F *) fListHist->FindObject("NbClusters");\r
 736   TH1F *nbClustersOffline = (TH1F *) fListHist->FindObject("NbClustersOffline");\r
 737   TH1F *nbClustersStandalone = (TH1F *) fListHist->FindObject("NbClustersStandalone");\r
 738 \r
 739   TH1F *nbTracklets = (TH1F *) fListHist->FindObject("NbTracklets");\r
 740   TH1F *nbTrackletsOffline = (TH1F *) fListHist->FindObject("NbTrackletsOffline");\r
 741   TH1F *nbTrackletsStandalone = (TH1F *) fListHist->FindObject("NbTrackletsStandalone");\r
 742   \r
 743   /////////////////////////////////////\r
 744   // Take the calibration objects\r
 745   /////////////////////////////////////\r
 746 \r
 747   TH2I *ch2d = (TH2I *) fListHist->FindObject("CH2d");\r
 748   TProfile2D *ph2d = (TProfile2D *) fListHist->FindObject("PH2d");\r
 749 \r
 750   TH2I *ch2dSum = (TH2I *) fListHist->FindObject("CH2dSum");\r
 751   TProfile2D *ph2dSum = (TProfile2D *) fListHist->FindObject("PH2dSum");\r
 752 \r
 753   TH2I *ch2dSM = (TH2I *) fListHist->FindObject("CH2dSM");\r
 754   TProfile2D *ph2dSM = (TProfile2D *) fListHist->FindObject("PH2dSM");\r
 755   \r
 756   AliTRDCalibraVdriftLinearFit *linearfit = (AliTRDCalibraVdriftLinearFit *) fListHist->FindObject("AliTRDCalibraVdriftLinearFit");\r
 757   \r
 758   ////////////////////////////////////////////////\r
 759   // Add the AliTRDCalibraVdriftLinearFit\r
 760   ///////////////////////////////////////////////\r
 761   \r
 762   Int_t first = 0;\r
 763   TH2S *histolinearfitsum = 0x0;\r
 764   \r
 765   if(linearfit) {\r
 766     for(Int_t det = 0; det < 540; det++) {\r
 767       if(linearfit->GetLinearFitterHisto(det)){\r
 768         if(TMath::Abs(first)<0.0001){\r
 769           histolinearfitsum = linearfit->GetLinearFitterHisto(det);\r
 770           first += 1;\r
 771         }\r
 772         else {\r
 773           histolinearfitsum ->Add(linearfit->GetLinearFitterHisto(det));\r
 774         }\r
 775       }\r
 776     }\r
 777   }\r
 778 \r
 779   ///////////////////////////\r
 780   // Style\r
 781   //////////////////////////\r
 782 \r
 783   gStyle->SetPalette(1);\r
 784   gStyle->SetOptStat(1111);\r
 785   gStyle->SetOptFit(1111);\r
 786   gStyle->SetPadBorderMode(0);\r
 787   gStyle->SetCanvasColor(10);\r
 788   gStyle->SetPadLeftMargin(0.13);\r
 789   gStyle->SetPadRightMargin(0.13);\r
 790 \r
 791   /////////////////////////\r
 792   // Plot\r
 793   ////////////////////////\r
 794 \r
 795  if(nEvents) {\r
 796    \r
 797     TCanvas *debugEvents = new TCanvas("cNEvents","cNEvents",10,10,510,510);\r
 798     debugEvents->cd(1);\r
 799     if(nEvents) nEvents->Draw();\r
 800       \r
 801   }\r
 802 \r
 803   if(trdTrack || tpctrdTrack) {\r
 804     \r
 805     TCanvas *debugtrdtpcTrack = new TCanvas("TRDtracktpctrdtrack","TRDtracktpctrdtrack",10,10,510,510);\r
 806     debugtrdtpcTrack->Divide(2,1);\r
 807     debugtrdtpcTrack->cd(1);\r
 808     if(trdTrack) trdTrack->Draw();\r
 809     if(trdTrackOffline) trdTrackOffline->Draw("same");\r
 810     if(trdTrackStandalone) trdTrackStandalone->Draw("same");\r
 811     TLegend *leg = new TLegend(0.4,0.6,0.89,0.89);\r
 812     if(trdTrack) leg->AddEntry(trdTrack,"All","p");\r
 813     if(trdTrackOffline) leg->AddEntry(trdTrackOffline,"Offline","p");\r
 814     if(trdTrackStandalone) leg->AddEntry(trdTrackStandalone,"Standalone","p");\r
 815     leg->Draw("same");\r
 816     debugtrdtpcTrack->cd(2);\r
 817     if(tpctrdTrack) tpctrdTrack->Draw();\r
 818     TLine *line = new TLine(0.0,0.0,100.0,100.0);\r
 819     line->Draw("same");\r
 820     \r
 821   }\r
 822  \r
 823   if(nbTimeBin || nbTracklets || nbClusters) {\r
 824     \r
 825     TCanvas *debugTracklets = new TCanvas("TRDtimebintrackletcluster","TRDtimebintrackletcluster",10,10,510,510);\r
 826     debugTracklets->Divide(3,1);\r
 827     debugTracklets->cd(1);\r
 828     if(nbTimeBin) nbTimeBin->Draw();\r
 829     if(nbTimeBinOffline) nbTimeBinOffline->Draw("same");\r
 830     if(nbTimeBinStandalone) nbTimeBinStandalone->Draw("same");\r
 831     TLegend *lega = new TLegend(0.4,0.6,0.89,0.89);\r
 832     if(nbTimeBin) lega->AddEntry(nbTimeBin,"All","p");\r
 833     if(nbTimeBinOffline) lega->AddEntry(nbTimeBinOffline,"Offline","p");\r
 834     if(nbTimeBinStandalone) lega->AddEntry(nbTimeBinStandalone,"Standalone","p");\r
 835     lega->Draw("same");\r
 836     debugTracklets->cd(2);\r
 837     if(nbTracklets) nbTracklets->Draw();\r
 838     if(nbTrackletsOffline) nbTrackletsOffline->Draw("same");\r
 839     if(nbTrackletsStandalone) nbTrackletsStandalone->Draw("same");\r
 840     TLegend *legb = new TLegend(0.4,0.6,0.89,0.89);\r
 841     if(nbTracklets) legb->AddEntry(nbTracklets,"All","p");\r
 842     if(nbTrackletsOffline) legb->AddEntry(nbTrackletsOffline,"Offline","p");\r
 843     if(nbTrackletsStandalone) legb->AddEntry(nbTrackletsStandalone,"Standalone","p");\r
 844     legb->Draw("same");\r
 845     debugTracklets->cd(3);\r
 846     if(nbClusters) nbClusters->Draw();\r
 847     if(nbClustersOffline) nbClustersOffline->Draw("same");\r
 848     if(nbClustersStandalone) nbClustersStandalone->Draw("same");\r
 849     TLegend *legc = new TLegend(0.4,0.6,0.89,0.89);\r
 850     if(nbClusters) legc->AddEntry(nbClusters,"All","p");\r
 851     if(nbClustersOffline) legc->AddEntry(nbClustersOffline,"Offline","p");\r
 852     if(nbClustersStandalone) legc->AddEntry(nbClustersStandalone,"Standalone","p");\r
 853     legc->Draw("same");\r
 854   \r
 855   }\r
 856 \r
 857   if(ch2dSum || ph2dSum || histolinearfitsum) {\r
 858     \r
 859     TCanvas *debugSum = new TCanvas("SumCalibrationObjects","SumCalibrationObjects",10,10,510,510);\r
 860     debugSum->Divide(3,1);\r
 861     debugSum->cd(1);\r
 862     if(ch2dSum) ch2dSum->Draw("lego");\r
 863     debugSum->cd(2);\r
 864     if(ph2dSum) ph2dSum->Draw("lego");\r
 865     debugSum->cd(3);\r
 866     if(histolinearfitsum) histolinearfitsum->Draw();\r
 867   \r
 868   }\r
 869 \r
 870   if(ch2dSM || ph2dSM) {\r
 871     \r
 872     TCanvas *debugSM = new TCanvas("SMCalibrationObjects","SMCalibrationObjects",10,10,510,510);\r
 873     debugSM->Divide(2,1);\r
 874     debugSM->cd(1);\r
 875     if(ch2dSM) ch2dSM->Draw("lego");\r
 876     debugSM->cd(2);\r
 877     if(ph2dSM) ph2dSM->Draw("lego");\r
 878      \r
 879   }\r
 880 \r
 881   if(ch2d || ph2d) {\r
 882     \r
 883     TCanvas *debug = new TCanvas("CalibrationObjects","CalibrationObjects",10,10,510,510);\r
 884     debug->Divide(2,1);\r
 885     debug->cd(1);\r
 886     if(ch2d) ch2d->Draw("lego");\r
 887     debug->cd(2);\r
 888     if(ph2d) ph2d->Draw("lego");\r
 889      \r
 890   }\r
 891  \r
 892 }\r