Some leftover obsolete statements removed from IceF2k.cxx
[u/mrichter/AliRoot.git] / RALICE / icepack / iceconvert / IceF2k.cxx
CommitLineData
f67e2651 1/*******************************************************************************
2 * Copyright(c) 2003, IceCube Experiment at the South Pole. All rights reserved.
3 *
4 * Author: The IceCube RALICE-based Offline 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.
12 * The authors make no claims about the suitability of this software for
13 * any purpose. It is provided "as is" without express or implied warranty.
14 *******************************************************************************/
15
16// $Id$
17
18///////////////////////////////////////////////////////////////////////////
19// Class IceF2k
20// Conversion of Amanda F2K data into IceEvent physics event structures.
21//
22// Usage example :
23// ---------------
24//
25// gSystem->Load("ralice");
26// gSystem->Load("icepack");
27// gSystem->Load("iceconvert");
28//
29// // Output file for the event structures
30// TFile* ofile=new TFile("events.root","RECREATE","F2K data in IceEvent structure");
31// TTree* otree=new TTree("T","Data of an Amanda run");
32//
33// // Limit the number of entries for testing
34// Int_t nentries=300;
35//
36// // Print frequency to produce a short summary print every printfreq events
37// Int_t printfreq=10;
38//
39// // Split level for the output structures
40// Int_t split=2;
41//
42// // Buffer size for the output structures
43// Int_t bsize=32000;
44//
45// IceF2k q("run8000.f2k",split,bsize);
46// q.Loop(otree,nentries,printfreq);
47//
48// // Select various objects to be added to the output file
49//
50// AliObjMatrix* omdb=q.GetOMdbase();
51// if (omdb) omdb->Write();
52//
53// AliDevice* fitdefs=q.GetFitdefs();
54// if (fitdefs) fitdefs->Write();
55//
56// TDatabasePDG* pdg=q.GetPDG();
57// if (pdg) pdg->Write();
58//
59// // Close output file
60// ofile->Write();
61// ofile->Close();
62//
63//--- Author: Nick van Eijndhoven 11-mar-2005 Utrecht University
64//- Modified: NvE $Date$ Utrecht University
65///////////////////////////////////////////////////////////////////////////
66
67#include "IceF2k.h"
68#include "Riostream.h"
69
70ClassImp(IceF2k) // Class implementation to enable ROOT I/O
71
72IceF2k::IceF2k(char* fname,Int_t split,Int_t bsize)
73{
74// Default constructor.
75// Initialise the input file and data structres to be converted.
76// Also the required split level and buffer size of the output tree
77// can be specified in this constructor.
78// By default tree=0, split=0 and bsize=32000.
79
80 fSplit=split;
81 fBsize=bsize;
82
83 fPdg=0;
84 fOmdb=0;
85 fFitdefs=0;
86
87 if (!fname)
88 {
89 cout << " *IceF2k ctor* No data input file specified." << endl;
90 return;
91 }
92
93 // Open the input file in the default ascii format (autodetection) for reading
94 fInput=rdmc_mcopen(fname,"r",RDMC_DEFAULT_ASCII_F);
95
96 if (!fInput)
97 {
98 cout << " *IceF2k ctor* No input file found with name : " << fname << endl;
99 return;
100 }
101
102 // Initialise the event structure
103 rdmc_init_mevt(&fEvent);
104
105 // Read the file header information
106 rdmc_rarr(fInput,&fHeader);
107}
108///////////////////////////////////////////////////////////////////////////
109IceF2k::~IceF2k()
110{
111// Default destructor.
112 if (fPdg)
113 {
114 delete fPdg;
115 fPdg=0;
116 }
117
118 if (fOmdb)
119 {
120 delete fOmdb;
121 fOmdb=0;
122 }
123
124 if (fFitdefs)
125 {
126 delete fFitdefs;
127 fFitdefs=0;
128 }
129}
130///////////////////////////////////////////////////////////////////////////
131TDatabasePDG* IceF2k::GetPDG()
132{
133// Provide pointer to the PDG database
134 return fPdg;
135}
136///////////////////////////////////////////////////////////////////////////
137AliObjMatrix* IceF2k::GetOMdbase()
138{
139// Provide pointer to the OM geometry, calib. etc... database
140 return fOmdb;
141}
142///////////////////////////////////////////////////////////////////////////
143AliDevice* IceF2k::GetFitdefs()
144{
145// Provide pointer to the fit definitions
146 return fFitdefs;
147}
148///////////////////////////////////////////////////////////////////////////
149void IceF2k::Loop(TTree* otree,Int_t nentries,Int_t printfreq)
150{
151// Loop over the specified number of entries and convert the
152// F2K data into the IceEvent structure.
153// The output will be written on the output tree specified as "otree".
154// If otree=0, a default standard output tree will be created.
155// If nentries<0 (default) all the entries of the input file
156// will be processed.
157// Every "printfreq" events a short event summary will be printed.
158// The default value is printfreq=1.
159
160 if (!fInput || fSplit<0) return;
161
162 if (!otree) otree=new TTree("T","F2K Data");
163
164 Double_t pi=acos(-1.);
165
166 IceEvent* evt=new IceEvent();
167
168 evt->SetTrackCopy(1);
169 evt->SetDevCopy(1);
170
171 // Branch in the tree for the event structure
172 otree->Branch("IceEvent","IceEvent",&evt,fBsize,fSplit);
173
174 // Create the particle database and extend it with some F2000 specific definitions
175 if (!fPdg) fPdg=new TDatabasePDG();
176 Double_t me=fPdg->GetParticle(11)->Mass();
177 fPdg->AddParticle("brems" ,"brems" ,0,1,0,0,"none",10001001,0,0);
178 fPdg->AddParticle("deltae" ,"deltae" ,me,1,0,-3,"Lepton",10001002,0,0);
179 fPdg->AddParticle("pairprod","pairprod",0,1,0,0,"none",10001003,0,0);
180 fPdg->AddParticle("nucl_int","nucl_Int",0,1,0,0,"none",10001004,0,0);
181 fPdg->AddParticle("mu_pair" ,"mu_pair" ,0,1,0,0,"none",10001005,0,0);
182 fPdg->AddParticle("hadrons" ,"hadrons" ,0,1,0,0,"none",10001006,0,0);
183 fPdg->AddParticle("fiberlaser","fiberlaser",0,1,0,0,"none",10002100,0,0);
184 fPdg->AddParticle("n2laser" ,"n2laser" ,0,1,0,0,"none",10002101,0,0);
185 fPdg->AddParticle("yaglaser" ,"yaglaser" ,0,1,0,0,"none",10002201,0,0);
186 fPdg->AddParticle("z_primary","z_primary",0,1,0,0,"none",10003000,0,0);
187 fPdg->AddParticle("a_primary","a_primary",0,1,0,0,"none",10003500,0,0);
188
189 // Fill the database with geometry, calib. etc... parameters
190 // for all the devices
191 FillOMdbase();
192
193 // Set the fit definitions according to the F2000 header info
194 SetFitdefs();
195
f67e2651 196 for (Int_t jentry=0; jentry<nentries; jentry++)
197 {
198 if (rdmc_revt(fInput,&fHeader,&fEvent) != 0) break;
199
200 // Reset the complete Event structure
201 evt->Reset();
202
203 evt->SetRunNumber(fEvent.nrun);
204 evt->SetEventNumber(fEvent.enr);
205 evt->SetMJD(fEvent.mjd,fEvent.secs,fEvent.nsecs);
206
207 PutMcTracks(evt);
208
209 PutRecoTracks(evt);
210
211 PutHits(evt);
212
f67e2651 213 if (!(jentry%printfreq))
214 {
215 evt->HeaderData();
216 }
217
218 // Write the complete structure to the output Tree
219 otree->Fill();
220 }
221
222 if (evt) delete evt;
223}
224///////////////////////////////////////////////////////////////////////////
225void IceF2k::FillOMdbase()
226{
227// Fill the database with geometry, calib. etc... parameters
228// for all the devices.
229
230 if (fHeader.nch<=0) return;
231
232 if (fOmdb)
233 {
234 fOmdb->Reset();
235 }
236 else
237 {
238 fOmdb=new AliObjMatrix();
239 fOmdb->SetNameTitle("OMDBASE","The OM geometry, calib. etc... database");
240 fOmdb->SetOwner();
241 }
242
243 IceAOM* dev=0;
244 Double_t pos[3]={0,0,0};
245 for (Int_t i=0; i<fHeader.nch; i++)
246 {
247 dev=new IceAOM();
248 dev->SetUniqueID(i+1);
249 dev->SetSlotName("TYPE",1);
250 dev->SetSlotName("ORIENT",2);
251 dev->SetSlotName("T0",3);
252 dev->SetSlotName("ALPHA",4);
253 dev->SetSlotName("KADC",5);
254 dev->SetSlotName("KTOT",6);
255 dev->SetSlotName("KTDC",7);
256
257 pos[0]=fHeader.x[i];
258 pos[1]=fHeader.y[i];
259 pos[2]=fHeader.z[i];
260 dev->SetPosition(pos,"car");
261 dev->SetSignal(fHeader.type[i],1);
262 dev->SetSignal((Float_t)fHeader.costh[i],2);
263 dev->SetSignal(fHeader.cal[i].t_0,3);
264 dev->SetSignal(fHeader.cal[i].alpha_t,4);
265 dev->SetSignal(fHeader.cal[i].beta_a,5);
266 dev->SetSignal(fHeader.cal[i].beta_tot,6);
267 dev->SetSignal(fHeader.cal[i].beta_t,7);
268 fOmdb->EnterObject(i+1,1,dev);
269 }
270}
271///////////////////////////////////////////////////////////////////////////
272void IceF2k::SetFitdefs()
273{
274// Obtain the names of the variables for each fit procedure from the
275// f2000 header. Each different fit procedure is then stored as a separate
276// hit of an AliDevice object and the various fit variables are stored
277// as separate signal slots of the corresponding hit.
278// Via the GetFitdefs() memberfunction this AliDevice object can be
279// retrieved and stored in the ROOT output file if wanted.
280// The name of the object is FitDefinitions and the stored data can be
281// inspected via the AliDevice::Data() memberfunction and looks as follows :
282//
283// *AliDevice::Data* Id :0 Name : FitDefinitions
284// Position Vector in car coordinates : 0 0 0
285// Err. in car coordinates : 0 0 0
286// The following 8 hits are registered :
287// *AliSignal::Data* Id :0
288// Position Vector in car coordinates : 0 0 0
289// Err. in car coordinates : 0 0 0
290// Owned by device : AliDevice Name : FitDefinitions
291// Slot : 1 Signal value : 1 name : id
292// Slot : 2 Signal value : 2 name : rchi2
293// Slot : 3 Signal value : 3 name : prob
294// Slot : 4 Signal value : 4 name : sigth
295// Slot : 5 Signal value : 5 name : covmin
296// Slot : 6 Signal value : 6 name : covmax
297// Slot : 7 Signal value : 7 name : cutflag
298// Slot : 8 Signal value : 8 name : chi2
299// *AliSignal::Data* Id :1
300// Position Vector in car coordinates : 0 0 0
301// Err. in car coordinates : 0 0 0
302// Owned by device : AliDevice Name : FitDefinitions
303// Slot : 1 Signal value : 1 name : id
304// Slot : 2 Signal value : 2 name : rchi2
305// Slot : 3 Signal value : 3 name : prob
306// etc....
307//
308// This memberfunction is based on the original idea/code by Adam Bouchta.
309
310 if (fHeader.n_fit<=0) return;
311
312 if (fFitdefs)
313 {
314 fFitdefs->Reset(1);
315 }
316 else
317 {
318 fFitdefs=new AliDevice();
319 }
320
321 fFitdefs->SetName("FitDefinitions");
322 fFitdefs->SetHitCopy (1);
323
324 AliSignal s;
325 s.Reset();
326
327 for (Int_t i=0; i<fHeader.n_fit; i++)
328 {
329 s.SetUniqueID(fHeader.def_fit[i].id);
330
331 for (Int_t j=0; j<fHeader.def_fit[i].nwords; j++)
332 {
333 s.SetSlotName(TString(fHeader.def_fit[i].words[j]),j+1);
334 s.SetSignal(j+1,j+1);
335 }
336
337 fFitdefs->AddHit(s);
338 s.Reset(1);
339 }
340}
341///////////////////////////////////////////////////////////////////////////
342void IceF2k::PutMcTracks(IceEvent* evt)
343{
344// Get the MC tracks from the F2000 file into the IcePack structure.
345// Note : MC tracks are given negative track id's in the event structure.
346// This memberfunction is based on the original code by Adam Bouchta.
347
348 if (!evt || fEvent.ntrack<=0) return;
349
350 // Loop over all the tracks and add them to the current event
351 AliTrack t;
352 Double_t vec[3];
353 AliPosition r;
354 Ali3Vector p;
355 Int_t tid=0;
356 Int_t idpdg=0;
357 Int_t idf2k=0;
358 for (Int_t i=0; i<fEvent.ntrack; i++)
359 {
360 t.Reset ();
361
362 // Beginpoint of the track
363 vec[0]=fEvent.gen[i].x;
364 vec[1]=fEvent.gen[i].y;
365 vec[2]=fEvent.gen[i].z;
366 r.SetPosition(vec,"car");
367 t.SetBeginPoint(r);
368
369 // Endpoint of the track
370 vec[0]+=fEvent.gen[i].length*fEvent.gen[i].px;
371 vec[1]+=fEvent.gen[i].length*fEvent.gen[i].py;
372 vec[2]+=fEvent.gen[i].length*fEvent.gen[i].pz;
373 r.SetPosition(vec,"car");
374 t.SetEndPoint(r);
375
376 // Momentum in GeV/c
377 vec[0]=fEvent.gen[i].e*fEvent.gen[i].px*1e-3;
378 vec[1]=fEvent.gen[i].e*fEvent.gen[i].py*1e-3;
379 vec[2]=fEvent.gen[i].e*fEvent.gen[i].pz*1e-3;
380 p.SetVector (vec,"car");
381 t.Set3Momentum(p);
382
383 // MC tracks are indicated by negative track id's
384 tid=fEvent.gen[i].tag;
385 t.SetId(-abs(tid));
386
387 idf2k=fEvent.gen[i].id;
388 idpdg=0;
389 if (idf2k>1000)
390 {
391 idpdg=idf2k+10000000;
392 }
393 else if (idf2k <= 48)
394 {
395 idpdg=fPdg->ConvertGeant3ToPdg(idf2k);
396 }
397 else
398 {
399 if (idf2k==201) idpdg=12;
400 if (idf2k==202) idpdg=14;
401 if (idf2k==203) idpdg=16;
402 if (idf2k==204) idpdg=-12;
403 if (idf2k==205) idpdg=-14;
404 if (idf2k==206) idpdg=-16;
405 }
406
407 t.SetParticleCode(idpdg);
408 t.SetName(fPdg->GetParticle(idpdg)->GetName());
409 t.SetTitle("MC track");
410 t.SetMass(fPdg->GetParticle(idpdg)->Mass());
411 t.SetCharge(fPdg->GetParticle(idpdg)->Charge()/3.);
412
413 evt->AddTrack(t);
414 }
415
416 // Create the pointers to each particle's parent particle.
417 Int_t txid=0;
418 Int_t parid=0;
419 for (Int_t itk=1; itk<=evt->GetNtracks (); itk++)
420 {
421 AliTrack* tx=evt->GetTrack(itk);
422
423 if (!tx) continue;
424
425 txid=tx->GetId();
426
427 parid=-1;
428 for (Int_t j=0; j<fEvent.ntrack; j++)
429 {
430 tid=fEvent.gen[j].tag;
431 if (-abs(tid) == txid) parid=fEvent.gen[j].parent;
432 }
433
434 if (parid<0) continue;
435
436 AliTrack* tpar=evt->GetIdTrack(-abs(parid));
437
438 if (!tpar) continue;
439
440 tx->SetParentTrack(tpar);
441 }
442}
443///////////////////////////////////////////////////////////////////////////
444void IceF2k::PutRecoTracks(IceEvent* evt)
445{
446// Get the reconstructed tracks from the F2000 file into the IcePack structure.
447// Note : Reco tracks are given positive track id's in the event structure.
448// This memberfunction is based on the original code by Adam Bouchta.
449
450 if (!evt || fEvent.nfit<=0) return;
451
452 // Loop over all the tracks and add them to the current event
453 AliTrack t;
454 Double_t vec[3];
455 AliPosition r;
456 Ali3Vector p;
457 Int_t tid=0;
458 Int_t idpdg=0;
459 Int_t idf2k=0;
460 for (Int_t i=0; i<fEvent.nfit; i++)
461 {
462 t.Reset ();
463
464 // Beginpoint of the track
465 vec[0]=fEvent.rec[i].x;
466 vec[1]=fEvent.rec[i].y;
467 vec[2]=fEvent.rec[i].z;
468 r.SetPosition(vec,"car");
469 t.SetBeginPoint(r);
470
471 // Endpoint of the track
472 vec[0]+=fEvent.rec[i].length*fEvent.rec[i].px;
473 vec[1]+=fEvent.rec[i].length*fEvent.rec[i].py;
474 vec[2]+=fEvent.rec[i].length*fEvent.rec[i].pz;
475 r.SetPosition(vec,"car");
476 t.SetEndPoint(r);
477
478 // Momentum in GeV/c
479 if (fEvent.rec[i].e > 0)
480 {
481 vec[0]=fEvent.rec[i].e*fEvent.rec[i].px*1e-3;
482 vec[1]=fEvent.rec[i].e*fEvent.rec[i].py*1e-3;
483 vec[2]=fEvent.rec[i].e*fEvent.rec[i].pz*1e-3;
484 }
485 else // Give the track a nominal momentum of 1 GeV/c
486 {
487 vec[0]=fEvent.rec[i].px;
488 vec[1]=fEvent.rec[i].py;
489 vec[2]=fEvent.rec[i].pz;
490 }
491 p.SetVector (vec,"car");
492 t.Set3Momentum(p);
493
494 // Use the fit number as track id
495 tid=fEvent.rec[i].tag;
496 t.SetId(abs(tid));
497
498 idf2k=fEvent.rec[i].id;
499 idpdg=0;
500 if (idf2k>1000)
501 {
502 idpdg=idf2k+10000000;
503 }
504 else if (idf2k <= 48)
505 {
506 idpdg=fPdg->ConvertGeant3ToPdg(idf2k);
507 }
508 else
509 {
510 if (idf2k==201) idpdg=12;
511 if (idf2k==202) idpdg=14;
512 if (idf2k==203) idpdg=16;
513 if (idf2k==204) idpdg=-12;
514 if (idf2k==205) idpdg=-14;
515 if (idf2k==206) idpdg=-16;
516 }
517
518 t.SetParticleCode(idpdg);
519 t.SetName(fPdg->GetParticle(idpdg)->GetName());
520 t.SetTitle("RECO track");
521 t.SetMass(fPdg->GetParticle(idpdg)->Mass());
522 t.SetCharge(fPdg->GetParticle(idpdg)->Charge()/3.);
523
524 // Retrieve the various fit parameters for this track
525 AliSignal* fitdata=fFitdefs->GetIdHit(i);
526 for (Int_t jval=0; jval<fEvent.fresult[i].nval; jval++)
527 {
528 fitdata->SetSignal(fEvent.fresult[i].val[jval],jval+1);
529 }
530
531 // Store the various fit parameters for this track
532 t.SetFitDetails(fitdata);
533
534 // Store the various reco tracks as track hypotheses.
535 // A copy of the first reco track is entered as a new track instance
536 // into the event and all reco tracks (incl. the first one) are
537 // stored as hypotheses linked to this new reco track.
538 if (i==0)
539 {
540 evt->AddTrack(t);
541 AliTrack* tx=evt->GetTrack(evt->GetNtracks());
542 Int_t nrec=evt->GetNtracks(1);
543 tx->SetId(nrec+1);
544 }
545 AliTrack* tx=evt->GetTrack(evt->GetNtracks());
546 if (tx) tx->AddTrackHypothesis(t);
547 }
548}
549///////////////////////////////////////////////////////////////////////////
550void IceF2k::PutHits(IceEvent* evt)
551{
552// Get the hit and waveform info from the F2000 file into the IcePack structure.
553// This memberfunction is based on the original code by Adam Bouchta.
554
555 if (!evt) return;
556
557 // Loop over all the hits and add them to the current event
558 IceAOM om;
559 AliSignal s;
560 s.SetSlotName("ADC",1);
561 s.SetSlotName("LE",2);
562 s.SetSlotName("TOT",3);
563 Int_t chan=0;
564 Int_t maxchan=800;
565 if (fOmdb) maxchan=fHeader.nch;
566 IceAOM* omx=0;
567 AliSignal* sx=0;
568 Int_t tid=0;
569 AliTrack* tx=0;
570 for (Int_t i=0; i<fEvent.nhits; i++)
571 {
572 chan=fEvent.h[i].ch+1;
573 if (chan>maxchan) continue; // Channels 9001, 9002 etc are trigger channels
574
575 // Get corresponding device from the current event structure
576 omx=(IceAOM*)evt->GetIdDevice(chan);
577 if (!omx)
578 {
579 if (fOmdb)
580 {
581 omx=(IceAOM*)fOmdb->GetObject(chan,1);
582 evt->AddDevice(omx);
583 }
584 else
585 {
586 om.Reset(1);
587 om.SetUniqueID(chan);
588 evt->AddDevice(om);
589 }
590 omx=(IceAOM*)evt->GetIdDevice(chan);
591 }
592
593 if (!omx) continue;
594
595 s.Reset();
596 s.SetUniqueID(fEvent.h[i].id);
597 s.SetSignal(fEvent.h[i].amp,1);
598 s.SetSignal(fEvent.h[i].t,2);
599 s.SetSignal(fEvent.h[i].tot,3);
600
601 omx->AddHit(s);
602
603 sx=omx->GetHit(omx->GetNhits());
604 if (!sx) continue;
605
606 // Bi-directional link between this hit and the track that caused the ADC value.
607 // This F2K info is probably only present for MC tracks.
608 tid=fEvent.h[i].ma;
609 if (tid > 0)
610 {
611 tx=evt->GetIdTrack(tid); // Reco tracks
612 if (!tx) tx=evt->GetIdTrack(-tid); // MC tracks
613 if (tx) sx->AddLink(tx);
614 }
615 else
616 {
617 if (tid == -2) s.SetNameTitle("N","Noise");
618 if (tid == -3) s.SetNameTitle("A","Afterpulse");
619 }
620 }
621
622 // Loop over all the waveforms and add the histo(s) to the corresponding OM's
623 TH1F histo;
624 Int_t nbins=0;
625 Float_t xlow=0;
626 Float_t xup=0;
627 TString hname;
628 for (Int_t iwf=0; iwf<fEvent.nwf; iwf++)
629 {
630 chan=fEvent.wf[iwf].om;
631 if (chan<=0 || chan>maxchan) continue; // Skip trigger channels
632
633 // Get corresponding device from the current event structure
634 omx=(IceAOM*)evt->GetIdDevice(chan);
635 if (!omx)
636 {
637 if (fOmdb)
638 {
639 omx=(IceAOM*)fOmdb->GetObject(chan,1);
640 evt->AddDevice(omx);
641 }
642 else
643 {
644 om.Reset(1);
645 om.SetUniqueID(chan);
646 evt->AddDevice(om);
647 }
648 omx=(IceAOM*)evt->GetIdDevice(chan);
649 }
650
651 if (!omx) continue;
652
653 omx->SetSlotName("BASELINE",omx->GetNnames()+1);
654 omx->SetSignal(fEvent.wf[iwf].baseline,"BASELINE");
655
656 // Fill the waveform histogram
657 hname="OM";
658 hname+=chan;
659 hname+="-WF";
660 hname+=omx->GetNwaveforms()+1;
661
662 histo.Reset();
663 histo.SetName(hname.Data());
664 nbins=fEvent.wf[iwf].ndigi;
665 xlow=fEvent.wf[iwf].t_start;
666 xup=xlow+float(nbins)*fEvent.wf[iwf].t_bin;
667 histo.SetBins(nbins,xlow,xup);
668
669 for (Int_t jbin=1; jbin<=fEvent.wf[iwf].ndigi; jbin++)
670 {
671 histo.SetBinContent(jbin,fEvent.wf[iwf].digi[jbin]);
672 }
673
674 omx->SetWaveform(&histo,omx->GetNwaveforms()+1);
675 }
676
677 // Set bi-directional links between hits and reco track hypotheses.
678 // Note : Reco tracks are recognised by their positive id.
679 Int_t hid=0;
680 TObjArray* rectracks=evt->GetTracks(1);
681 for (Int_t jtk=0; jtk<rectracks->GetEntries(); jtk++)
682 {
683 tx=(AliTrack*)rectracks->At(jtk);
684 if (!tx) continue;
685
686 for (Int_t jhyp=1; jhyp<=tx->GetNhypotheses(); jhyp++)
687 {
688 AliTrack* hypx=tx->GetTrackHypothesis(jhyp);
689 if (!hypx) continue;
690
691 // Loop over all combinations of F2K fits and used OM hits
692 for (Int_t k=0; k<fEvent.nfit_uses; k++)
693 {
694 if (fEvent.fit_uses[k].useid != hypx->GetId()) continue;
695 hid=fEvent.fit_uses[k].hitid;
696 sx=evt->GetIdHit(hid,"IceAOM");
697 if (sx) sx->AddLink(hypx);
698 }
699 }
700 }
701}
702///////////////////////////////////////////////////////////////////////////