TRD/AliTRDgtuSim.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /* $Id: AliTRDgtuSim.cxx 28397 2008-09-02 09:33:00Z cblume $ */
  17
  18 ////////////////////////////////////////////////////////////////////////////
  19 //                                                                        //
  20 //  GTU simulation                                                        //
  21 //                                                                        //
  22 //  Authors: J. Klein (Jochen.Klein@cern.ch)                              //
  23 //                                                                        //
  24 ////////////////////////////////////////////////////////////////////////////
  25
  26 #include <stdio.h>
  27 #include <fstream>
  28 #include <string>
  29
  30 #include "TFile.h"
  31 #include "TROOT.h"
  32 #include "TObject.h"
  33 #include "TClonesArray.h"
  34
  35 #include "AliRun.h"
  36 #include "AliRunLoader.h"
  37 #include "AliLoader.h"
  38 #include "AliTreeLoader.h"
  39 #include "AliLog.h"
  40 #include "AliESDTrdTrack.h"
  41 #include "AliESDTrdTracklet.h"
  42
  43 #include "AliTRDgtuSim.h"
  44 #include "AliTRDfeeParam.h"
  45 #include "AliTRDgtuTMU.h"
  46 #include "AliTRDtrackGTU.h"
  47 #include "AliTRDtrackletWord.h"
  48 #include "AliTRDtrackletMCM.h"
  49 #include "AliESDEvent.h"
  50
  51 ClassImp(AliTRDgtuSim)
  52
  53 AliTRDgtuSim::AliTRDgtuSim(AliRunLoader *rl)
  54   : TObject(),
  55   fRunLoader(rl),
  56   fFeeParam(AliTRDfeeParam::Instance()),
  57   fTMU(0x0),
  58   fTrackletArray(0x0),
  59   fTrackTree(0x0),
  60   fTrackletTree(0x0)
  61 {
  62 //   fTrackletTree = new TTree("gtutracklets", "Tree with GTU tracklets");
  63 //   fTrackletTree->SetDirectory(0);
  64 }
  65
  66 AliTRDgtuSim::~AliTRDgtuSim()
  67 {
  68   // destructor
  69
  70   if (fTrackletArray)
  71     fTrackletArray->Clear();
  72   delete fTrackletArray;
  73   delete fTrackletTree;
  74 }
  75
  76 Bool_t AliTRDgtuSim::RunGTUFromTrackletFile(TString filename, Int_t event, Int_t noev)
  77 {
  78   // run the GTU from a file of tracklets
  79   // used for comparison to VHDL simulation
  80
  81   ifstream input(filename.Data());
  82
  83   std::string str;
  84   TString string;
  85   int lineno = 0;
  86
  87   Int_t iEventPrev = -1;
  88   Int_t iStackPrev = -1;
  89   Int_t iSecPrev = -1;
  90   Int_t iSec = -1;
  91   Int_t iStack = -1;
  92   Int_t iLink = -1;
  93   Int_t iEvent = -1;
  94   Int_t evcnt = -1;
  95
  96   fTMU = 0x0;
  97
  98   TClonesArray trklArray("AliTRDtrackletWord", 100);
  99   TClonesArray trklArrayGTU("AliTRDtrackletGTU", 100);
 100
 101   AliDebug(1, Form("--------- Reading from %s ----------", filename.Data()));
 102   while (getline(input, str)) {
 103     lineno++;
 104     string = str;
 105
 106     TObjArray *tokens = string.Tokenize(" ");
 107     if (tokens->GetEntriesFast() < 7) {
 108       AliWarning(Form("Invalid input in line %i, too few parameters", lineno));
 109       continue;
 110     }
 111
 112     if ( ((TObjString*) tokens->At(0))->GetString().Atoi() < event)
 113       continue;
 114
 115     iEvent = ((TObjString*) tokens->At(0))->GetString().Atoi();
 116     iSec = ((TObjString*) tokens->At(1))->GetString().Atoi();
 117     iStack = ((TObjString*) tokens->At(2))->GetString().Atoi();
 118     iLink = 2 * ((TObjString*) tokens->At(3))->GetString().Atoi() + ((TObjString*) tokens->At(4))->GetString().Atoi();
 119
 120     if ((iEvent != iEventPrev) ||
 121         (iStack != iStackPrev) ||
 122         (iSec != iSecPrev)) {
 123       if(fTMU) {
 124         TList *listOfTracks = new TList();
 125         fTMU->SetStack(iStackPrev);
 126         fTMU->SetSector(iSecPrev);
 127         fTMU->RunTMU(listOfTracks);
 128         AliDebug(1,Form("--- There are %i tracks. Writing ...", listOfTracks->GetEntries()));
 129         WriteTracksToTree(listOfTracks);
 130         // fTMU->WriteTrackletsToTree(fTrackletTree);
 131         // WriteTracksToLoader(listOfTracks);
 132         WriteTracksToDataFile(listOfTracks, iEventPrev);
 133         if (listOfTracks->GetEntries() > 0)
 134           AliDebug(2,Form("   %4.1f GeV/c", ((AliTRDtrackGTU*) listOfTracks->At(0))->GetPt() ));
 135         delete fTMU;
 136         fTMU = new AliTRDgtuTMU();
 137         delete listOfTracks;
 138         listOfTracks = 0x0;
 139       } else {
 140         fTMU = new AliTRDgtuTMU();
 141       }
 142       iStackPrev = iStack;
 143       iSecPrev = iSec;
 144       iEventPrev = iEvent;
 145       evcnt++;
 146       if (evcnt == noev)
 147         break;
 148     }
 149     for (Int_t i = 5; i < tokens->GetEntriesFast(); i++) {
 150       UInt_t trackletWord = 0;
 151       sscanf(((TObjString*) tokens->At(i))->GetString().Data(), "%i", &trackletWord);
 152       if (trackletWord == 0x10001000)
 153         break;
 154       AliDebug(2, Form("link: %2i trkl: %2i - %s -> 0x%08x",
 155                        iLink, i-4, ((TObjString*) tokens->At(i))->GetString().Data(), trackletWord));
 156       AliTRDtrackletWord *tracklet = new (trklArray[trklArray.GetEntriesFast()])       AliTRDtrackletWord(trackletWord);
 157       AliTRDtrackletGTU   *trkl    = new (trklArrayGTU[trklArrayGTU.GetEntriesFast()]) AliTRDtrackletGTU(tracklet);
 158       if (fTMU)
 159         fTMU->AddTracklet(trkl, iLink);
 160     }
 161   }
 162
 163   if (fTMU && evcnt < noev) {
 164     TList *listOfTracks = new TList();
 165     fTMU->SetStack(iStackPrev);
 166     fTMU->SetSector(iSecPrev);
 167     fTMU->RunTMU(listOfTracks);
 168     WriteTracksToTree(listOfTracks);
 169     // fTMU->WriteTrackletsToTree(fTrackletTree);
 170     // WriteTracksToLoader(listOfTracks);
 171     WriteTracksToDataFile(listOfTracks, iEventPrev);
 172     delete fTMU;
 173     delete listOfTracks;
 174     fTMU = 0x0;
 175   }
 176
 177   AliInfo(Form("Analyzed %i events", evcnt));
 178   return kTRUE;
 179 }
 180
 181 Bool_t AliTRDgtuSim::RunGTU(AliLoader *loader, AliESDEvent *esd)
 182 {
 183   // run the GTU on tracklets taken from the loader
 184   // if specified the GTU tracks are written to the ESD event
 185
 186   if (!fFeeParam->GetTracklet())
 187     return kFALSE;
 188
 189   if (fTrackletArray)
 190     fTrackletArray->Clear();
 191
 192   if (loader) {
 193     if (!LoadTracklets(loader)) {
 194         AliError("Could not load the tracklets. Nothing done ...");
 195         return kFALSE;
 196     }
 197   }
 198   else {
 199     LoadTracklets(esd);
 200   }
 201
 202     AliDebug(1, Form("running on %i tracklets", fTrackletArray->GetEntriesFast()));
 203
 204     Int_t iStackPrev = -1;
 205     Int_t iSecPrev = -1;
 206     Int_t iSec = -1;
 207     Int_t iStack = -1;
 208     Int_t iLink = -1;
 209
 210     if (fTMU) {
 211         delete fTMU;
 212         fTMU = 0x0;
 213     }
 214
 215     TList *listOfTracks = new TList();
 216
 217     TIter next(fTrackletArray);
 218
 219     while (AliTRDtrackletGTU *trkl = (AliTRDtrackletGTU*) next()) {
 220         iSec = trkl->GetDetector() / 30;
 221         iStack = (trkl->GetDetector() % 30) / 6;
 222         iLink = trkl->GetHCId() % 12;
 223
 224         if (iStack != iStackPrev || iSec != iSecPrev) {
 225             if(fTMU) {
 226                 fTMU->SetStack(iStackPrev);
 227                 fTMU->SetSector(iSecPrev);
 228                 fTMU->RunTMU(listOfTracks);
 229                 WriteTracksToTree(listOfTracks);
 230 //              fTMU->WriteTrackletsToTree(fTrackletTree);
 231                 WriteTracksToLoader(listOfTracks);
 232                 WriteTracksToESD(listOfTracks, esd);
 233                 fTMU->Reset();
 234                 listOfTracks->Delete();
 235             } else {
 236                 fTMU = new AliTRDgtuTMU();
 237             }
 238             iStackPrev = iStack;
 239             iSecPrev = iSec;
 240             AliDebug(1, Form("now in sec %i, stack %i", iSec, iStack));
 241         }
 242         AliDebug(1, Form("adding tracklet: 0x%08x in sec %i stack %i link %i",
 243                          trkl->GetTrackletWord(), trkl->GetDetector() / 30, (trkl->GetDetector() % 30) / 6, trkl->GetHCId() % 12));
 244         if (fTMU) {
 245           fTMU->AddTracklet(trkl, iLink);
 246         }
 247     }
 248
 249     if (fTMU) {
 250         fTMU->SetStack(iStackPrev);
 251         fTMU->SetSector(iSecPrev);
 252         fTMU->RunTMU(listOfTracks);
 253         WriteTracksToTree(listOfTracks);
 254 //      fTMU->WriteTrackletsToTree(fTrackletTree);
 255         WriteTracksToLoader(listOfTracks);
 256         WriteTracksToESD(listOfTracks, esd);
 257         delete fTMU;
 258         fTMU = 0x0;
 259         listOfTracks->Delete();
 260     }
 261
 262     delete listOfTracks;
 263
 264     return kTRUE;
 265 }
 266
 267 Bool_t AliTRDgtuSim::LoadTracklets(const AliESDEvent *const esd)
 268 {
 269   AliDebug(1,"Loading tracklets from ESD event ...");
 270
 271   if (!fTrackletArray)
 272     fTrackletArray = new TClonesArray("AliTRDtrackletGTU", 1000);
 273
 274   for (Int_t iTracklet = 0; iTracklet < esd->GetNumberOfTrdTracklets(); iTracklet++) {
 275     AliESDTrdTracklet *trkl = esd->GetTrdTracklet(iTracklet);
 276     new ((*fTrackletArray)[fTrackletArray->GetEntries()]) AliTRDtrackletGTU(trkl);
 277   }
 278
 279   return kTRUE;
 280 }
 281
 282 Bool_t AliTRDgtuSim::LoadTracklets(AliLoader *const loader)
 283 {
 284   // load the tracklets using the given loader
 285
 286   AliDebug(1,"Loading tracklets ...");
 287
 288   if (!fFeeParam->GetTracklet())
 289     return kFALSE;
 290
 291   if (!loader) {
 292     AliError("No loader given!");
 293     return kFALSE;
 294   }
 295
 296   AliDataLoader *trackletLoader = loader->GetDataLoader("tracklets");
 297   if (!trackletLoader) {
 298       AliError("No tracklet loader found!");
 299       return kFALSE;
 300   }
 301
 302   trackletLoader->Load();
 303   TTree *trackletTree = 0x0;
 304
 305   // simulated tracklets
 306   trackletTree = trackletLoader->Tree();
 307   if (trackletTree) {
 308     TBranch *trklbranch = trackletTree->GetBranch("mcmtrklbranch");
 309     if (trklbranch) {
 310       if (!fTrackletArray)
 311         fTrackletArray = new TClonesArray("AliTRDtrackletGTU", 1000);
 312
 313       AliTRDtrackletMCM *trkl = 0x0;
 314       trklbranch->SetAddress(&trkl);
 315       for (Int_t iTracklet = 0; iTracklet < trklbranch->GetEntries(); iTracklet++) {
 316         trklbranch->GetEntry(iTracklet);
 317         new ((*fTrackletArray)[fTrackletArray->GetEntries()]) AliTRDtrackletGTU(new AliTRDtrackletMCM(*trkl));
 318       }
 319       return kTRUE;
 320     }
 321   }
 322
 323   // raw tracklets
 324   AliTreeLoader *tl = (AliTreeLoader*) trackletLoader->GetBaseLoader("tracklets-raw");
 325   trackletTree = tl ? tl->Load(), tl->Tree() : 0x0;
 326
 327   if (trackletTree) {
 328     if (!fTrackletArray)
 329       fTrackletArray = new TClonesArray("AliTRDtrackletGTU", 1000);
 330
 331     Int_t hc;
 332     TClonesArray *ar = 0x0;
 333     trackletTree->SetBranchAddress("hc", &hc);
 334     trackletTree->SetBranchAddress("trkl", &ar);
 335
 336     for (Int_t iEntry = 0; iEntry < trackletTree->GetEntries(); iEntry++) {
 337       trackletTree->GetEntry(iEntry);
 338       AliDebug(2, Form("%i tracklets in HC %i", ar->GetEntriesFast(), hc));
 339       for (Int_t iTracklet = 0; iTracklet < ar->GetEntriesFast(); iTracklet++) {
 340         AliTRDtrackletWord *trklWord = (AliTRDtrackletWord*) (*ar)[iTracklet];
 341         new((*fTrackletArray)[fTrackletArray->GetEntriesFast()]) AliTRDtrackletGTU(new AliTRDtrackletWord(trklWord->GetTrackletWord(), hc));
 342       }
 343     }
 344     return kTRUE;
 345   }
 346
 347   AliError("No raw tracklet tree found\n");
 348
 349   return kFALSE;
 350 }
 351
 352 Bool_t AliTRDgtuSim::WriteTracksToDataFile(TList *listOfTracks, Int_t event)
 353 {
 354   // write the found tracks to a data file
 355   // used for comparison to VHDL simulation
 356
 357     Int_t sm = 0;
 358     Int_t stack = 0;
 359
 360     FILE *out;
 361     out = fopen("test.data", "a");
 362
 363     AliDebug(1,Form("%i tracks found in event %i", listOfTracks->GetSize(), event));
 364     // fprintf(out, "0 %5i %2i %i  00000000\n", event, sm, stack);
 365     for (Int_t i = 0; i < listOfTracks->GetSize(); i++) {
 366         AliTRDtrackGTU *trk = (AliTRDtrackGTU*) listOfTracks->At(i);
 367         sm = trk->GetSector();
 368         stack = trk->GetStack();
 369
 370         ULong64_t trackWord = 1;
 371         AppendBits(trackWord,   1, 0);
 372         AppendBits(trackWord,   6, trk->GetTrackletMask());
 373         AppendBits(trackWord,  18, (Int_t) trk->GetA());
 374         AppendBits(trackWord,  18, (Int_t) trk->GetB());
 375         AppendBits(trackWord,  12, (Int_t) trk->GetC());
 376         AppendBits(trackWord,   8, trk->GetPID());
 377         fprintf(out, "ev. %i sec. %i stack %i - track word: 0x%016llx, ",
 378                 event, sm, stack, trackWord);
 379
 380         trackWord = 0;
 381         AppendBits(trackWord, 11, 0); // flags
 382         AppendBits(trackWord,  3, 0);
 383         AppendBits(trackWord, 13, trk->GetYapprox());
 384         AppendBits(trackWord,  6, trk->GetTrackletIndex(5));
 385         AppendBits(trackWord,  6, trk->GetTrackletIndex(4));
 386         AppendBits(trackWord,  6, trk->GetTrackletIndex(3));
 387         AppendBits(trackWord,  6, trk->GetTrackletIndex(2));
 388         AppendBits(trackWord,  6, trk->GetTrackletIndex(1));
 389         AppendBits(trackWord,  6, trk->GetTrackletIndex(0));
 390         fprintf(out, "extended track word: 0x%016llx\n", trackWord);
 391
 392         fprintf(out, "1 %5i %2i %2i %3i %3i %3i %3i %3i %3i %3i %4i %f\n", event, sm, stack, trk->GetTrackletMask(),
 393                trk->GetTrackletIndex(5),
 394                trk->GetTrackletIndex(4),
 395                trk->GetTrackletIndex(3),
 396                trk->GetTrackletIndex(2),
 397                trk->GetTrackletIndex(1),
 398                trk->GetTrackletIndex(0),
 399                trk->GetPtInt(),
 400                trk->GetPt());
 401     }
 402     fclose(out);
 403     return kTRUE;
 404 }
 405
 406 Bool_t AliTRDgtuSim::WriteTracksToTree(TList *listOfTracks, Int_t /*event*/)
 407 {
 408   // write the tracks to the tree for intermediate storage
 409
 410   AliDebug(1,Form("Writing %i tracks to the tree...", listOfTracks->GetEntries()));
 411
 412   if (!listOfTracks)
 413     return kFALSE;
 414
 415   if (listOfTracks->GetEntries() <= 0)
 416     return kTRUE;
 417
 418   if (!fTrackTree) {
 419     fTrackTree = new TTree("gtutracks", "GTU tracks");
 420     fTrackTree->SetDirectory(0);
 421   }
 422
 423   AliTRDtrackGTU *trk = 0x0;
 424   TBranch *branch = fTrackTree->GetBranch("TRDgtuTrack");
 425   if (!branch) {
 426       branch = fTrackTree->Branch("TRDgtuTrack", "AliTRDtrackGTU", &trk, 32000, 99);
 427   }
 428
 429   TIter next(listOfTracks);
 430   while ((trk = (AliTRDtrackGTU*) next())) {
 431       trk->CookLabel();
 432       branch->SetAddress(&trk);
 433       fTrackTree->Fill();
 434   }
 435   fTrackTree->ResetBranchAddress(branch);
 436
 437   return kTRUE;
 438 }
 439
 440 Bool_t AliTRDgtuSim::WriteTreesToFile() const {
 441   // write the trees holding tracklets and tracks to file
 442
 443   TFile *f = TFile::Open("TRD.GtuTracking.root", "RECREATE");
 444   f->cd();
 445   if (fTrackTree)
 446     f->WriteTObject(fTrackTree);
 447   if (fTrackletTree)
 448     f->WriteTObject(fTrackletTree);
 449   f->Close();
 450   return kTRUE;
 451 }
 452
 453 Bool_t AliTRDgtuSim::WriteTracksToESD(const TList * const listOfTracks, AliESDEvent *esd)
 454 {
 455   // fill the found tracks to the given ESD event
 456
 457     if (esd) {
 458         TIter next(listOfTracks);
 459         while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
 460             AliESDTrdTrack *trdtrack = trk->CreateTrdTrack();
 461             if (trdtrack->GetLabel() < 0)
 462               trdtrack->SetLabel(-2);
 463             esd->AddTrdTrack(trdtrack);
 464             delete trdtrack;
 465         }
 466     }
 467     return kTRUE;
 468 }
 469
 470 Bool_t AliTRDgtuSim::WriteTracksToLoader(const TList * const listOfTracks)
 471 {
 472   // write the GTU tracks to the dedicated loader
 473   // these tracks contain more information than the ones in the ESD
 474
 475   if (!fTrackTree) {
 476     AliDebug(1, "No track tree found!");
 477     return kFALSE;
 478   }
 479
 480   AliRunLoader *rl = AliRunLoader::Instance();
 481   AliDataLoader *dl = 0x0;
 482   if (rl)
 483     dl = rl->GetLoader("TRDLoader")->GetDataLoader("gtutracks");
 484   if (!dl) {
 485     AliError("Could not get the GTU-track data loader!");
 486     return kFALSE;
 487   }
 488
 489   TTree *trackTree = dl->Tree();
 490   if (!trackTree) {
 491     dl->MakeTree();
 492     trackTree = dl->Tree();
 493   }
 494
 495   AliTRDtrackGTU *trk = 0x0;
 496
 497   if (!trackTree->GetBranch("TRDtrackGTU"))
 498     trackTree->Branch("TRDtrackGTU", "AliTRDtrackGTU", &trk, 32000);
 499
 500   TIter next(listOfTracks);
 501   while ((trk = (AliTRDtrackGTU*) next())) {
 502     trackTree->SetBranchAddress("TRDtrackGTU", &trk);
 503     trackTree->Fill();
 504   }
 505
 506   dl->WriteData("OVERWRITE");
 507
 508   return kTRUE;
 509 }