1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
16 /* $Id: AliTRDgtuSim.cxx 28397 2008-09-02 09:33:00Z cblume $ */
18 ////////////////////////////////////////////////////////////////////////////
22 // Authors: J. Klein (Jochen.Klein@cern.ch) //
24 ////////////////////////////////////////////////////////////////////////////
32 #include "TClonesArray.h"
35 #include "AliRunLoader.h"
36 #include "AliLoader.h"
38 #include "AliESDTrdTrack.h"
40 #include "AliTRDgtuSim.h"
41 #include "AliTRDgtuTMU.h"
42 #include "AliTRDtrackGTU.h"
43 #include "AliTRDtrackletWord.h"
44 #include "AliTRDtrackletMCM.h"
45 #include "AliESDEvent.h"
47 ClassImp(AliTRDgtuSim)
49 AliTRDgtuSim::AliTRDgtuSim(AliRunLoader *rl)
57 fTrackletTree = new TTree("gtutracklets", "Tree with GTU tracklets");
58 fTrackletTree->SetDirectory(0);
61 AliTRDgtuSim::~AliTRDgtuSim()
66 fTrackletArray->Delete();
67 delete fTrackletArray;
71 Bool_t AliTRDgtuSim::RunGTUFromTrackletFile(TString filename, Int_t event, Int_t noev)
73 // run the GTU from a file of tracklets
74 // used for comparison to VHDL simulation
76 AliInfo(Form("Running the GTU simulation on file: %s", filename.Data()));
77 ifstream input(filename.Data());
83 Int_t iEventPrev = -1;
84 Int_t iStackPrev = -1;
94 AliDebug(5,"--------- Reading from file ----------");
95 while (getline(input, str)) {
98 AliDebug(5,Form("Line %i : %s", lineno, string.Data()));
100 TObjArray *tokens = string.Tokenize(" ");
101 if (tokens->GetEntriesFast() < 7) {
102 AliWarning(Form("Invalid input in line %i, too few parameters", lineno));
106 if ( ((TObjString*) tokens->At(0))->GetString().Atoi() < event)
109 iEvent = ((TObjString*) tokens->At(0))->GetString().Atoi();
110 iSec = ((TObjString*) tokens->At(1))->GetString().Atoi();
111 iStack = ((TObjString*) tokens->At(2))->GetString().Atoi();
112 iLink = 2 * ((TObjString*) tokens->At(3))->GetString().Atoi() + ((TObjString*) tokens->At(4))->GetString().Atoi();
114 if (iEvent != iEventPrev || iStack != iStackPrev || iSec != iSecPrev) {
116 TList *listOfTracks = new TList();
117 fTMU->SetStack(iStackPrev);
118 fTMU->SetSector(iSecPrev);
119 fTMU->RunTMU(listOfTracks);
120 AliDebug(1,Form("--- There are %i tracks. Writing ...", listOfTracks->GetEntries()));
121 WriteTracksToTree(listOfTracks);
122 fTMU->WriteTrackletsToTree(fTrackletTree);
123 WriteTracksToDataFile(listOfTracks, iEventPrev);
124 if (listOfTracks->GetEntries() > 0)
125 AliDebug(2,Form(" %d GeV/c", ((AliTRDtrackGTU*) listOfTracks->At(0))->GetPt() ));
127 fTMU = new AliTRDgtuTMU();
131 fTMU = new AliTRDgtuTMU();
140 for (Int_t i = 5; i < tokens->GetEntriesFast(); i++) {
141 UInt_t trackletWord = 0;
142 sscanf(((TObjString*) tokens->At(i))->GetString().Data(), "%x", &trackletWord);
143 if (trackletWord == 0x10001000)
145 AliDebug(2,Form("%i. tracklet: %s -> 0x%08x", i-4, ((TObjString*) tokens->At(i))->GetString().Data(), trackletWord));
146 AliTRDtrackletWord *trkl = new AliTRDtrackletWord(trackletWord);
147 fTMU->AddTracklet(trkl, iLink);
151 if (fTMU && evcnt < noev) {
152 TList *listOfTracks = new TList();
153 fTMU->SetStack(iStackPrev);
154 fTMU->SetSector(iSecPrev);
155 fTMU->RunTMU(listOfTracks);
156 WriteTracksToTree(listOfTracks);
157 fTMU->WriteTrackletsToTree(fTrackletTree);
158 WriteTracksToDataFile(listOfTracks, iEventPrev);
164 AliInfo(Form("Analyzed %i events", evcnt));
168 Bool_t AliTRDgtuSim::RunGTU(AliLoader *loader, AliESDEvent *esd)
170 // run the GTU on tracklets taken from the loader
171 // if specified the GTU tracks are written to the ESD event
173 if (!LoadTracklets(loader)) {
174 AliError("Could not load the tracklets. Nothing done ...");
178 Int_t iStackPrev = -1;
189 TList *listOfTracks = new TList();
191 TIter next(fTrackletArray);
192 AliTRDtrackletBase *trkl;
194 while ((trkl = (AliTRDtrackletBase*) next())) {
195 iSec = trkl->GetDetector() / 30;
196 iStack = (trkl->GetDetector() % 30) / 6;
197 iLink = 2 * (trkl->GetDetector() % 6) + (trkl->GetYbin() < 0 ? 0 : 1);
199 if (iStack != iStackPrev || iSec != iSecPrev) {
201 fTMU->SetStack(iStackPrev);
202 fTMU->SetSector(iSecPrev);
203 fTMU->RunTMU(listOfTracks);
204 WriteTracksToTree(listOfTracks);
205 fTMU->WriteTrackletsToTree(fTrackletTree);
206 WriteTracksToESD(listOfTracks, esd);
208 listOfTracks->Delete();
210 fTMU = new AliTRDgtuTMU();
215 fTMU->AddTracklet(trkl, iLink);
219 fTMU->SetStack(iStackPrev);
220 fTMU->SetSector(iSecPrev);
221 fTMU->RunTMU(listOfTracks);
222 WriteTracksToTree(listOfTracks);
223 fTMU->WriteTrackletsToTree(fTrackletTree);
224 WriteTracksToESD(listOfTracks, esd);
227 listOfTracks->Delete();
235 Bool_t AliTRDgtuSim::LoadTracklets(AliLoader *const loader)
237 // load the tracklets using the given loader
239 AliDebug(1,"Loading tracklets ...");
242 AliError("No loader given!");
246 AliDataLoader *trackletLoader = loader->GetDataLoader("tracklets");
247 if (!trackletLoader) {
248 AliError("No tracklet loader found!");
252 trackletLoader->Load();
254 TTree *trackletTree = trackletLoader->Tree();
256 AliError("No tracklet tree found");
261 TBranch *trklbranch = trackletTree->GetBranch("mcmtrklbranch");
264 fTrackletArray = new TClonesArray("AliTRDtrackletMCM", 1000);
265 else if ((TClass::GetClass("AliTRDtrackletMCM"))->InheritsFrom(fTrackletArray->Class()))
266 fTrackletArray->Delete();
268 fTrackletArray->Delete();
269 delete fTrackletArray;
270 fTrackletArray = new TClonesArray("AliTRDtrackletMCM", 1000);
273 AliTRDtrackletMCM *trkl = new AliTRDtrackletMCM;
274 trklbranch->SetAddress(&trkl);
275 for (Int_t iTracklet = 0; iTracklet < trklbranch->GetEntries(); iTracklet++) {
276 trklbranch->GetEntry(iTracklet);
277 new ((*fTrackletArray)[fTrackletArray->GetEntries()]) AliTRDtrackletMCM(*trkl);
282 trklbranch = trackletTree->GetBranch("trkbranch");
285 AliError("Could not get trkbranch");
290 fTrackletArray = new TClonesArray("AliTRDtrackletWord", 1000);
291 else if ((TClass::GetClass("AliTRDtrackletWord"))->InheritsFrom(fTrackletArray->Class()))
292 fTrackletArray->Delete();
294 fTrackletArray->Delete();
295 delete fTrackletArray;
296 fTrackletArray = new TClonesArray("AliTRDtrackletWord", 1000);
300 UInt_t *leaves = new UInt_t[258];
301 AliDebug(1,Form("No. of entries: %i", trklbranch->GetEntries()));
303 for (Int_t iEntry = 0; iEntry < trklbranch->GetEntries(); iEntry++) {
304 trklbranch->SetAddress(leaves);
305 trklbranch->GetEntry(iEntry);
306 for (Int_t iTracklet = 0; iTracklet < 256; iTracklet++) {
307 if (leaves[2 + iTracklet] == 0)
309 new((*fTrackletArray)[notrkl]) AliTRDtrackletWord(leaves[2 + iTracklet], 2*leaves[0] + leaves[1]);
312 AliDebug(2,Form("Entry: %3i: Det: %3i, side: %i, 1st tracklet: 0x%08x, no: %i", iEntry, leaves[0], leaves[1], leaves[2], notrkl));
318 Bool_t AliTRDgtuSim::WriteTracksToDataFile(TList *listOfTracks, Int_t event)
320 // write the found tracks to a data file
321 // used for comparison to VHDL simulation
327 out = fopen("test.data", "a");
329 AliDebug(1,Form("%i tracks found in event %i", listOfTracks->GetSize(), event));
330 fprintf(out, "0 %5i %2i %i 00000000\n", event, sm, stack);
331 for (Int_t i = 0; i < listOfTracks->GetSize(); i++) {
332 AliTRDtrackGTU *trk = (AliTRDtrackGTU*) listOfTracks->At(i);
333 sm = trk->GetSector();
334 stack = trk->GetStack();
335 fprintf(out, "1 %5i %2i %2i %3i %3i %3i %3i %3i %3i %3i %4i %f\n", event, sm, stack, trk->GetTrackletMask(),
336 trk->GetTrackletIndex(5),
337 trk->GetTrackletIndex(4),
338 trk->GetTrackletIndex(3),
339 trk->GetTrackletIndex(2),
340 trk->GetTrackletIndex(1),
341 trk->GetTrackletIndex(0),
349 Bool_t AliTRDgtuSim::WriteTracksToTree(TList *listOfTracks, Int_t /*event*/)
351 // write the tracks to the tree for intermediate storage
353 AliDebug(1,Form("Writing %i tracks to the tree...", listOfTracks->GetEntries()));
358 if (listOfTracks->GetEntries() <= 0)
362 fTrackTree = new TTree("gtutracks", "GTU tracks");
363 fTrackTree->SetDirectory(0);
366 AliTRDtrackGTU *trk = 0x0;
367 TBranch *branch = fTrackTree->GetBranch("TRDgtuTrack");
369 branch = fTrackTree->Branch("TRDgtuTrack", "AliTRDtrackGTU", &trk, 32000, 99);
372 TIter next(listOfTracks);
373 while ((trk = (AliTRDtrackGTU*) next())) {
375 branch->SetAddress(&trk);
378 fTrackTree->ResetBranchAddress(branch);
383 Bool_t AliTRDgtuSim::WriteTreesToFile() const {
384 // write the trees holding tracklets and tracks to file
386 TFile *f = TFile::Open("TRD.GtuTracking.root", "RECREATE");
389 f->WriteTObject(fTrackTree);
391 f->WriteTObject(fTrackletTree);
396 Bool_t AliTRDgtuSim::WriteTracksToESD(const TList * const listOfTracks, AliESDEvent *esd)
398 // fill the found tracks to the given ESD event
401 TIter next(listOfTracks);
402 while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
403 AliESDTrdTrack *trdtrack = trk->CreateTrdTrack();
404 esd->AddTrdTrack(trdtrack);
411 Bool_t AliTRDgtuSim::WriteTracksToLoader()
413 // write the GTU tracks to the dedicated loader
414 // these tracks contain more information than the ones in the ESD
417 AliError("No track tree found!");
421 AliRunLoader *rl = AliRunLoader::Instance();
422 AliDataLoader *dl = 0x0;
424 dl = rl->GetLoader("TRDLoader")->GetDataLoader("gtutracks");
426 AliError("Could not get the GTU-track data loader!");
430 TTree *trackTree = dl->Tree();
433 trackTree = dl->Tree();
436 AliTRDtrackGTU *trk = 0x0;
437 TBranch *trkbranch = trackTree->GetBranch("TRDtrackGTU");
439 trkbranch = trackTree->Branch("TRDtrackGTU", "AliTRDtrackGTU", &trk, 32000);
441 TBranch *branch = fTrackTree->GetBranch("TRDgtuTrack");
445 AliDebug(1,Form("Found %i tracks", branch->GetEntries()));
447 for (Int_t iTrack = 0; iTrack < branch->GetEntries(); iTrack++) {
448 branch->SetAddress(&trk);
449 branch->GetEntry(iTrack);
450 trkbranch->SetAddress(&trk);
454 dl->WriteData("OVERWRITE");