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"
37 #include "AliTreeLoader.h"
39 #include "AliESDTrdTrack.h"
41 #include "AliTRDgtuSim.h"
42 #include "AliTRDfeeParam.h"
43 #include "AliTRDgtuTMU.h"
44 #include "AliTRDtrackGTU.h"
45 #include "AliTRDtrackletWord.h"
46 #include "AliTRDtrackletMCM.h"
47 #include "AliESDEvent.h"
49 ClassImp(AliTRDgtuSim)
51 AliTRDgtuSim::AliTRDgtuSim(AliRunLoader *rl)
54 fFeeParam(AliTRDfeeParam::Instance()),
60 fTrackletTree = new TTree("gtutracklets", "Tree with GTU tracklets");
61 fTrackletTree->SetDirectory(0);
64 AliTRDgtuSim::~AliTRDgtuSim()
69 fTrackletArray->Delete();
70 delete fTrackletArray;
74 Bool_t AliTRDgtuSim::RunGTUFromTrackletFile(TString filename, Int_t event, Int_t noev)
76 // run the GTU from a file of tracklets
77 // used for comparison to VHDL simulation
79 AliInfo(Form("Running the GTU simulation on file: %s", filename.Data()));
80 ifstream input(filename.Data());
86 Int_t iEventPrev = -1;
87 Int_t iStackPrev = -1;
97 AliDebug(5,"--------- Reading from file ----------");
98 while (getline(input, str)) {
101 AliDebug(5,Form("Line %i : %s", lineno, string.Data()));
103 TObjArray *tokens = string.Tokenize(" ");
104 if (tokens->GetEntriesFast() < 7) {
105 AliWarning(Form("Invalid input in line %i, too few parameters", lineno));
109 if ( ((TObjString*) tokens->At(0))->GetString().Atoi() < event)
112 iEvent = ((TObjString*) tokens->At(0))->GetString().Atoi();
113 iSec = ((TObjString*) tokens->At(1))->GetString().Atoi();
114 iStack = ((TObjString*) tokens->At(2))->GetString().Atoi();
115 iLink = 2 * ((TObjString*) tokens->At(3))->GetString().Atoi() + ((TObjString*) tokens->At(4))->GetString().Atoi();
117 if (iEvent != iEventPrev || iStack != iStackPrev || iSec != iSecPrev) {
119 TList *listOfTracks = new TList();
120 fTMU->SetStack(iStackPrev);
121 fTMU->SetSector(iSecPrev);
122 fTMU->RunTMU(listOfTracks);
123 AliDebug(1,Form("--- There are %i tracks. Writing ...", listOfTracks->GetEntries()));
124 WriteTracksToTree(listOfTracks);
125 fTMU->WriteTrackletsToTree(fTrackletTree);
126 WriteTracksToDataFile(listOfTracks, iEventPrev);
127 if (listOfTracks->GetEntries() > 0)
128 AliDebug(2,Form(" %4.1f GeV/c", ((AliTRDtrackGTU*) listOfTracks->At(0))->GetPt() ));
130 fTMU = new AliTRDgtuTMU();
134 fTMU = new AliTRDgtuTMU();
143 for (Int_t i = 5; i < tokens->GetEntriesFast(); i++) {
144 UInt_t trackletWord = 0;
145 sscanf(((TObjString*) tokens->At(i))->GetString().Data(), "%i", &trackletWord);
146 if (trackletWord == 0x10001000)
148 AliDebug(2,Form("%i. tracklet: %s -> 0x%08x", i-4, ((TObjString*) tokens->At(i))->GetString().Data(), trackletWord));
149 AliTRDtrackletWord *trkl = new AliTRDtrackletWord(trackletWord);
151 fTMU->AddTracklet(trkl, iLink);
155 if (fTMU && evcnt < noev) {
156 TList *listOfTracks = new TList();
157 fTMU->SetStack(iStackPrev);
158 fTMU->SetSector(iSecPrev);
159 fTMU->RunTMU(listOfTracks);
160 WriteTracksToTree(listOfTracks);
161 fTMU->WriteTrackletsToTree(fTrackletTree);
162 WriteTracksToDataFile(listOfTracks, iEventPrev);
168 AliInfo(Form("Analyzed %i events", evcnt));
172 Bool_t AliTRDgtuSim::RunGTU(AliLoader *loader, AliESDEvent *esd)
174 // run the GTU on tracklets taken from the loader
175 // if specified the GTU tracks are written to the ESD event
177 if (!fFeeParam->GetTracklet())
180 if (!LoadTracklets(loader)) {
181 AliError("Could not load the tracklets. Nothing done ...");
185 AliDebug(1, Form("running on %i tracklets", fTrackletArray->GetEntriesFast()));
187 Int_t iStackPrev = -1;
198 TList *listOfTracks = new TList();
200 TIter next(fTrackletArray);
201 AliTRDtrackletBase *trkl;
203 while ((trkl = (AliTRDtrackletBase*) next())) {
204 iSec = trkl->GetDetector() / 30;
205 iStack = (trkl->GetDetector() % 30) / 6;
206 iLink = 2 * (trkl->GetDetector() % 6) + (trkl->GetYbin() < 0 ? 0 : 1);
208 if (iStack != iStackPrev || iSec != iSecPrev) {
210 fTMU->SetStack(iStackPrev);
211 fTMU->SetSector(iSecPrev);
212 fTMU->RunTMU(listOfTracks);
213 WriteTracksToTree(listOfTracks);
214 fTMU->WriteTrackletsToTree(fTrackletTree);
215 WriteTracksToESD(listOfTracks, esd);
217 listOfTracks->Delete();
219 fTMU = new AliTRDgtuTMU();
224 AliDebug(1, Form("adding tracklet: 0x%08x", trkl->GetTrackletWord()));
226 fTMU->AddTracklet(trkl, iLink);
230 fTMU->SetStack(iStackPrev);
231 fTMU->SetSector(iSecPrev);
232 fTMU->RunTMU(listOfTracks);
233 WriteTracksToTree(listOfTracks);
234 fTMU->WriteTrackletsToTree(fTrackletTree);
235 WriteTracksToESD(listOfTracks, esd);
238 listOfTracks->Delete();
246 Bool_t AliTRDgtuSim::LoadTracklets(AliLoader *const loader)
248 // load the tracklets using the given loader
250 AliDebug(1,"Loading tracklets ...");
252 if (!fFeeParam->GetTracklet())
256 AliError("No loader given!");
260 AliDataLoader *trackletLoader = loader->GetDataLoader("tracklets");
261 if (!trackletLoader) {
262 AliError("No tracklet loader found!");
266 trackletLoader->Load();
267 TTree *trackletTree = 0x0;
269 // simulated tracklets
270 trackletTree = trackletLoader->Tree();
272 TBranch *trklbranch = trackletTree->GetBranch("mcmtrklbranch");
275 fTrackletArray = new TClonesArray("AliTRDtrackletMCM", 1000);
276 else if ((TClass::GetClass("AliTRDtrackletMCM"))->InheritsFrom(fTrackletArray->Class()))
277 fTrackletArray->Delete();
279 fTrackletArray->Delete();
280 delete fTrackletArray;
281 fTrackletArray = new TClonesArray("AliTRDtrackletMCM", 1000);
284 AliTRDtrackletMCM *trkl = new AliTRDtrackletMCM;
285 trklbranch->SetAddress(&trkl);
286 for (Int_t iTracklet = 0; iTracklet < trklbranch->GetEntries(); iTracklet++) {
287 trklbranch->GetEntry(iTracklet);
288 new ((*fTrackletArray)[fTrackletArray->GetEntries()]) AliTRDtrackletMCM(*trkl);
295 AliTreeLoader *tl = (AliTreeLoader*) trackletLoader->GetBaseLoader("tracklets-raw");
296 trackletTree = tl ? tl->Load(), tl->Tree() : 0x0;
300 fTrackletArray = new TClonesArray("AliTRDtrackletWord", 1000);
301 else if ((TClass::GetClass("AliTRDtrackletWord"))->InheritsFrom(fTrackletArray->Class()))
302 fTrackletArray->Delete();
304 fTrackletArray->Delete();
305 delete fTrackletArray;
306 fTrackletArray = new TClonesArray("AliTRDtrackletWord", 1000);
310 TClonesArray *ar = 0x0;
311 trackletTree->SetBranchAddress("hc", &hc);
312 trackletTree->SetBranchAddress("trkl", &ar);
314 for (Int_t iEntry = 0; iEntry < trackletTree->GetEntries(); iEntry++) {
315 trackletTree->GetEntry(iEntry);
316 printf("%i tracklets in HC %i\n", ar->GetEntriesFast(), hc);
317 for (Int_t iTracklet = 0; iTracklet < ar->GetEntriesFast(); iTracklet++) {
318 AliTRDtrackletWord *trklWord = (AliTRDtrackletWord*) (*ar)[iTracklet];
319 new((*fTrackletArray)[fTrackletArray->GetEntriesFast()]) AliTRDtrackletWord(trklWord->GetTrackletWord(), hc);
325 AliError("No raw tracklet tree found\n");
330 Bool_t AliTRDgtuSim::WriteTracksToDataFile(TList *listOfTracks, Int_t event)
332 // write the found tracks to a data file
333 // used for comparison to VHDL simulation
339 out = fopen("test.data", "a");
341 AliDebug(1,Form("%i tracks found in event %i", listOfTracks->GetSize(), event));
342 fprintf(out, "0 %5i %2i %i 00000000\n", event, sm, stack);
343 for (Int_t i = 0; i < listOfTracks->GetSize(); i++) {
344 AliTRDtrackGTU *trk = (AliTRDtrackGTU*) listOfTracks->At(i);
345 sm = trk->GetSector();
346 stack = trk->GetStack();
347 fprintf(out, "1 %5i %2i %2i %3i %3i %3i %3i %3i %3i %3i %4i %f\n", event, sm, stack, trk->GetTrackletMask(),
348 trk->GetTrackletIndex(5),
349 trk->GetTrackletIndex(4),
350 trk->GetTrackletIndex(3),
351 trk->GetTrackletIndex(2),
352 trk->GetTrackletIndex(1),
353 trk->GetTrackletIndex(0),
361 Bool_t AliTRDgtuSim::WriteTracksToTree(TList *listOfTracks, Int_t /*event*/)
363 // write the tracks to the tree for intermediate storage
365 AliDebug(1,Form("Writing %i tracks to the tree...", listOfTracks->GetEntries()));
370 if (listOfTracks->GetEntries() <= 0)
374 fTrackTree = new TTree("gtutracks", "GTU tracks");
375 fTrackTree->SetDirectory(0);
378 AliTRDtrackGTU *trk = 0x0;
379 TBranch *branch = fTrackTree->GetBranch("TRDgtuTrack");
381 branch = fTrackTree->Branch("TRDgtuTrack", "AliTRDtrackGTU", &trk, 32000, 99);
384 TIter next(listOfTracks);
385 while ((trk = (AliTRDtrackGTU*) next())) {
387 branch->SetAddress(&trk);
390 fTrackTree->ResetBranchAddress(branch);
395 Bool_t AliTRDgtuSim::WriteTreesToFile() const {
396 // write the trees holding tracklets and tracks to file
398 TFile *f = TFile::Open("TRD.GtuTracking.root", "RECREATE");
401 f->WriteTObject(fTrackTree);
403 f->WriteTObject(fTrackletTree);
408 Bool_t AliTRDgtuSim::WriteTracksToESD(const TList * const listOfTracks, AliESDEvent *esd)
410 // fill the found tracks to the given ESD event
413 TIter next(listOfTracks);
414 while (AliTRDtrackGTU *trk = (AliTRDtrackGTU*) next()) {
415 AliESDTrdTrack *trdtrack = trk->CreateTrdTrack();
416 esd->AddTrdTrack(trdtrack);
423 Bool_t AliTRDgtuSim::WriteTracksToLoader()
425 // write the GTU tracks to the dedicated loader
426 // these tracks contain more information than the ones in the ESD
429 AliDebug(1, "No track tree found!");
433 AliRunLoader *rl = AliRunLoader::Instance();
434 AliDataLoader *dl = 0x0;
436 dl = rl->GetLoader("TRDLoader")->GetDataLoader("gtutracks");
438 AliError("Could not get the GTU-track data loader!");
442 TTree *trackTree = dl->Tree();
445 trackTree = dl->Tree();
448 AliTRDtrackGTU *trk = 0x0;
449 if (!trackTree->GetBranch("TRDtrackGTU"))
450 trackTree->Branch("TRDtrackGTU", "AliTRDtrackGTU", &trk, 32000);
452 AliDebug(1,Form("Found %lld tracks", fTrackTree->GetEntries()));
454 for (Int_t iTrack = 0; iTrack < fTrackTree->GetEntries(); iTrack++) {
455 fTrackTree->SetBranchAddress("TRDgtuTrack", &trk);
456 fTrackTree->GetEntry(iTrack);
457 trackTree->SetBranchAddress("TRDtrackGTU", &trk);
461 dl->WriteData("OVERWRITE");