20-apr-2005 NvE Id of owning device added to the output of AliSignal::Data().
[u/mrichter/AliRoot.git] / RALICE / icepack / iceconvert / IceF2k.cxx
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
70 ClassImp(IceF2k) // Class implementation to enable ROOT I/O
71
72 IceF2k::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 ///////////////////////////////////////////////////////////////////////////
109 IceF2k::~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 ///////////////////////////////////////////////////////////////////////////
131 TDatabasePDG* IceF2k::GetPDG()
132 {
133 // Provide pointer to the PDG database
134  return fPdg;
135 }
136 ///////////////////////////////////////////////////////////////////////////
137 AliObjMatrix* IceF2k::GetOMdbase()
138 {
139 // Provide pointer to the OM geometry, calib. etc... database
140  return fOmdb;
141 }
142 ///////////////////////////////////////////////////////////////////////////
143 AliDevice* IceF2k::GetFitdefs()
144 {
145 // Provide pointer to the fit definitions
146  return fFitdefs;
147 }
148 ///////////////////////////////////////////////////////////////////////////
149 void 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
196 /*************
197
198  // The LEDA specific output data
199  AliCalorimeter* ledaup=new AliCalorimeter(44,144);
200  AliCalorimeter* ledalw=new AliCalorimeter(40,144);
201
202  ledaup->SetName("LedaUp");
203  ledalw->SetName("LedaDown");
204
205  evt->InitLeda(ledaup);
206  evt->InitLeda(ledalw);
207
208  TDatime datim;
209  Float_t pos[3],costh;
210  AliSignal s;
211  s.SetName("CPV signal ADC");
212 ************************/
213
214  for (Int_t jentry=0; jentry<nentries; jentry++)
215  {
216   if (rdmc_revt(fInput,&fHeader,&fEvent) != 0) break;
217
218   // Reset the complete Event structure
219   evt->Reset();
220
221   evt->SetRunNumber(fEvent.nrun);
222   evt->SetEventNumber(fEvent.enr);
223   evt->SetMJD(fEvent.mjd,fEvent.secs,fEvent.nsecs);
224
225   PutMcTracks(evt);
226
227   PutRecoTracks(evt);
228
229   PutHits(evt);
230
231 /*********************************
232   datim.Set(Jdate,Jtime);
233   evt->SetDayTime(datim);
234   evt->SetProjectile(207,82,158);
235   evt->SetTarget(207,82,0);
236   evt->SetWeight(Jwscal);
237   evt->SetTrig(Itword);
238   evt->SetZdc(Zdc*1000.);
239   evt->SetMiracE(1000.*Emir,Emire,Emirh);
240   evt->SetMiracEt(Etm,Etme,Etmh);
241  
242   ledaup->Reset();
243   ledalw->Reset();
244   // Fill calorimeter with module data
245   for (Int_t i=0; i<Nmod; i++)
246   {
247    if (Adcl[i] > 3) // Adc cut of 3 to remove noise
248    {
249     if (Irowl[i] > 0) ledaup->SetSignal(Irowl[i],Icoll[i],Adcl[i]);
250     if (Irowl[i] < 0) ledalw->SetSignal(-Irowl[i],Icoll[i],Adcl[i]);
251    }
252   }
253
254   // Store associated CPV signals
255   for (Int_t j=0; j<Ncluv; j++)
256   {
257    s.Reset();
258    s.SetSignal(Iadccv[j]);
259    pos[1]=Thetacv[j]*pi/180.;
260    pos[2]=Phicv[j]*pi/180.;
261    costh=cos(pos[1]);
262    pos[0]=0;
263    if (costh) pos[0]=2103./costh;
264    s.SetPosition(pos,"sph");
265    pos[0]=0.4;
266    pos[1]=2.2;
267    pos[2]=0;
268    s.SetPositionErrors(pos,"car");
269    if (Phicv[j]>=0. && Phicv[j]<=180.)
270    {
271     ledaup->AddVetoSignal(s);
272    }
273    else
274    {
275     ledalw->AddVetoSignal(s);
276    }
277   }
278
279   evt->AddDevice(ledaup);
280   evt->AddDevice(ledalw);
281 ************************************/
282
283   if (!(jentry%printfreq))
284   {
285    evt->HeaderData();
286   }
287
288   // Write the complete structure to the output Tree
289   otree->Fill();
290  }
291
292  if (evt) delete evt;
293 }
294 ///////////////////////////////////////////////////////////////////////////
295 void IceF2k::FillOMdbase()
296 {
297 // Fill the database with geometry, calib. etc... parameters 
298 // for all the devices.
299
300  if (fHeader.nch<=0) return;
301
302  if (fOmdb)
303  {
304   fOmdb->Reset();
305  }
306  else
307  {
308   fOmdb=new AliObjMatrix();
309   fOmdb->SetNameTitle("OMDBASE","The OM geometry, calib. etc... database");
310   fOmdb->SetOwner();
311  }
312
313  IceAOM* dev=0;
314  Double_t pos[3]={0,0,0};
315  for (Int_t i=0; i<fHeader.nch; i++)
316  {
317   dev=new IceAOM();
318   dev->SetUniqueID(i+1);
319   dev->SetSlotName("TYPE",1);
320   dev->SetSlotName("ORIENT",2);
321   dev->SetSlotName("T0",3);
322   dev->SetSlotName("ALPHA",4);
323   dev->SetSlotName("KADC",5);
324   dev->SetSlotName("KTOT",6);
325   dev->SetSlotName("KTDC",7);
326
327   pos[0]=fHeader.x[i];
328   pos[1]=fHeader.y[i];
329   pos[2]=fHeader.z[i];
330   dev->SetPosition(pos,"car");
331   dev->SetSignal(fHeader.type[i],1);
332   dev->SetSignal((Float_t)fHeader.costh[i],2);
333   dev->SetSignal(fHeader.cal[i].t_0,3);
334   dev->SetSignal(fHeader.cal[i].alpha_t,4);
335   dev->SetSignal(fHeader.cal[i].beta_a,5);
336   dev->SetSignal(fHeader.cal[i].beta_tot,6);
337   dev->SetSignal(fHeader.cal[i].beta_t,7);
338   fOmdb->EnterObject(i+1,1,dev);
339  }
340 }
341 ///////////////////////////////////////////////////////////////////////////
342 void IceF2k::SetFitdefs()
343 {
344 // Obtain the names of the variables for each fit procedure from the
345 // f2000 header. Each different fit procedure is then stored as a separate
346 // hit of an AliDevice object and the various fit variables are stored
347 // as separate signal slots of the corresponding hit.
348 // Via the GetFitdefs() memberfunction this AliDevice object can be
349 // retrieved and stored in the ROOT output file if wanted.
350 // The name of the object is FitDefinitions and the stored data can be
351 // inspected via the AliDevice::Data() memberfunction and looks as follows :
352 //
353 //  *AliDevice::Data* Id :0 Name : FitDefinitions
354 //    Position Vector in car coordinates : 0 0 0
355 //    Err. in car coordinates : 0 0 0
356 //  The following 8 hits are registered :
357 //  *AliSignal::Data* Id :0
358 //    Position Vector in car coordinates : 0 0 0
359 //    Err. in car coordinates : 0 0 0
360 //    Owned by device : AliDevice Name : FitDefinitions
361 //    Slot : 1 Signal value : 1 name : id
362 //    Slot : 2 Signal value : 2 name : rchi2
363 //    Slot : 3 Signal value : 3 name : prob
364 //    Slot : 4 Signal value : 4 name : sigth
365 //    Slot : 5 Signal value : 5 name : covmin
366 //    Slot : 6 Signal value : 6 name : covmax
367 //    Slot : 7 Signal value : 7 name : cutflag
368 //    Slot : 8 Signal value : 8 name : chi2
369 //  *AliSignal::Data* Id :1
370 //    Position Vector in car coordinates : 0 0 0
371 //    Err. in car coordinates : 0 0 0
372 //    Owned by device : AliDevice Name : FitDefinitions
373 //    Slot : 1 Signal value : 1 name : id
374 //    Slot : 2 Signal value : 2 name : rchi2
375 //    Slot : 3 Signal value : 3 name : prob
376 // etc....  
377 //
378 // This memberfunction is based on the original idea/code by Adam Bouchta.
379
380  if (fHeader.n_fit<=0) return;
381
382  if (fFitdefs)
383  {
384   fFitdefs->Reset(1);
385  }
386  else
387  {
388   fFitdefs=new AliDevice();
389  }
390
391  fFitdefs->SetName("FitDefinitions");
392  fFitdefs->SetHitCopy (1);
393
394  AliSignal s;
395  s.Reset();
396
397  for (Int_t i=0; i<fHeader.n_fit; i++)
398  {
399   s.SetUniqueID(fHeader.def_fit[i].id);
400
401   for (Int_t j=0; j<fHeader.def_fit[i].nwords; j++)
402   {
403    s.SetSlotName(TString(fHeader.def_fit[i].words[j]),j+1);
404    s.SetSignal(j+1,j+1);
405   }
406
407   fFitdefs->AddHit(s);
408   s.Reset(1);
409  }
410 }
411 ///////////////////////////////////////////////////////////////////////////
412 void IceF2k::PutMcTracks(IceEvent* evt)
413 {
414 // Get the MC tracks from the F2000 file into the IcePack structure.
415 // Note : MC tracks are given negative track id's in the event structure.
416 // This memberfunction is based on the original code by Adam Bouchta.
417
418  if (!evt || fEvent.ntrack<=0) return;
419
420  // Loop over all the tracks and add them to the current event
421  AliTrack t;
422  Double_t vec[3];
423  AliPosition r;
424  Ali3Vector p;
425  Int_t tid=0;
426  Int_t idpdg=0;
427  Int_t idf2k=0;
428  for (Int_t i=0; i<fEvent.ntrack; i++)
429  {
430   t.Reset ();
431
432   // Beginpoint of the track
433   vec[0]=fEvent.gen[i].x;
434   vec[1]=fEvent.gen[i].y;
435   vec[2]=fEvent.gen[i].z;
436   r.SetPosition(vec,"car");
437   t.SetBeginPoint(r);
438
439   // Endpoint of the track
440   vec[0]+=fEvent.gen[i].length*fEvent.gen[i].px;
441   vec[1]+=fEvent.gen[i].length*fEvent.gen[i].py;
442   vec[2]+=fEvent.gen[i].length*fEvent.gen[i].pz;
443   r.SetPosition(vec,"car");
444   t.SetEndPoint(r);
445
446   // Momentum in GeV/c
447   vec[0]=fEvent.gen[i].e*fEvent.gen[i].px*1e-3;
448   vec[1]=fEvent.gen[i].e*fEvent.gen[i].py*1e-3;
449   vec[2]=fEvent.gen[i].e*fEvent.gen[i].pz*1e-3;
450   p.SetVector (vec,"car");
451   t.Set3Momentum(p);
452
453   // MC tracks are indicated by negative track id's
454   tid=fEvent.gen[i].tag;
455   t.SetId(-abs(tid));
456
457   idf2k=fEvent.gen[i].id;
458   idpdg=0;
459   if (idf2k>1000)
460   {
461    idpdg=idf2k+10000000;
462   }
463   else if (idf2k <= 48)
464   {
465    idpdg=fPdg->ConvertGeant3ToPdg(idf2k);
466   }
467   else
468   {
469    if (idf2k==201) idpdg=12;
470    if (idf2k==202) idpdg=14;
471    if (idf2k==203) idpdg=16;
472    if (idf2k==204) idpdg=-12;
473    if (idf2k==205) idpdg=-14;
474    if (idf2k==206) idpdg=-16;
475   }
476
477   t.SetParticleCode(idpdg);
478   t.SetName(fPdg->GetParticle(idpdg)->GetName());
479   t.SetTitle("MC track");
480   t.SetMass(fPdg->GetParticle(idpdg)->Mass());
481   t.SetCharge(fPdg->GetParticle(idpdg)->Charge()/3.);
482
483   evt->AddTrack(t);
484  }
485
486  // Create the pointers to each particle's parent particle.
487  Int_t txid=0;
488  Int_t parid=0;
489  for (Int_t itk=1; itk<=evt->GetNtracks (); itk++)
490  {
491   AliTrack* tx=evt->GetTrack(itk);
492
493   if (!tx) continue;
494
495   txid=tx->GetId();
496
497   parid=-1;
498   for (Int_t j=0; j<fEvent.ntrack; j++)
499   {
500    tid=fEvent.gen[j].tag;
501    if (-abs(tid) == txid) parid=fEvent.gen[j].parent;
502   }
503
504   if (parid<0) continue;
505
506   AliTrack* tpar=evt->GetIdTrack(-abs(parid));
507
508   if (!tpar) continue;
509
510   tx->SetParentTrack(tpar);
511  }
512 }
513 ///////////////////////////////////////////////////////////////////////////
514 void IceF2k::PutRecoTracks(IceEvent* evt)
515 {
516 // Get the reconstructed tracks from the F2000 file into the IcePack structure.
517 // Note : Reco tracks are given positive track id's in the event structure.
518 // This memberfunction is based on the original code by Adam Bouchta.
519
520  if (!evt || fEvent.nfit<=0) return;
521
522  // Loop over all the tracks and add them to the current event
523  AliTrack t;
524  Double_t vec[3];
525  AliPosition r;
526  Ali3Vector p;
527  Int_t tid=0;
528  Int_t idpdg=0;
529  Int_t idf2k=0;
530  for (Int_t i=0; i<fEvent.nfit; i++)
531  {
532   t.Reset ();
533
534   // Beginpoint of the track
535   vec[0]=fEvent.rec[i].x;
536   vec[1]=fEvent.rec[i].y;
537   vec[2]=fEvent.rec[i].z;
538   r.SetPosition(vec,"car");
539   t.SetBeginPoint(r);
540
541   // Endpoint of the track
542   vec[0]+=fEvent.rec[i].length*fEvent.rec[i].px;
543   vec[1]+=fEvent.rec[i].length*fEvent.rec[i].py;
544   vec[2]+=fEvent.rec[i].length*fEvent.rec[i].pz;
545   r.SetPosition(vec,"car");
546   t.SetEndPoint(r);
547
548   // Momentum in GeV/c
549   if (fEvent.rec[i].e > 0)
550   {
551    vec[0]=fEvent.rec[i].e*fEvent.rec[i].px*1e-3;
552    vec[1]=fEvent.rec[i].e*fEvent.rec[i].py*1e-3;
553    vec[2]=fEvent.rec[i].e*fEvent.rec[i].pz*1e-3;
554   }
555   else // Give the track a nominal momentum of 1 GeV/c
556   {
557    vec[0]=fEvent.rec[i].px;
558    vec[1]=fEvent.rec[i].py;
559    vec[2]=fEvent.rec[i].pz;
560   }
561   p.SetVector (vec,"car");
562   t.Set3Momentum(p);
563
564   // Use the fit number as track id
565   tid=fEvent.rec[i].tag;
566   t.SetId(abs(tid));
567
568   idf2k=fEvent.rec[i].id;
569   idpdg=0;
570   if (idf2k>1000)
571   {
572    idpdg=idf2k+10000000;
573   }
574   else if (idf2k <= 48)
575   {
576    idpdg=fPdg->ConvertGeant3ToPdg(idf2k);
577   }
578   else
579   {
580    if (idf2k==201) idpdg=12;
581    if (idf2k==202) idpdg=14;
582    if (idf2k==203) idpdg=16;
583    if (idf2k==204) idpdg=-12;
584    if (idf2k==205) idpdg=-14;
585    if (idf2k==206) idpdg=-16;
586   }
587
588   t.SetParticleCode(idpdg);
589   t.SetName(fPdg->GetParticle(idpdg)->GetName());
590   t.SetTitle("RECO track");
591   t.SetMass(fPdg->GetParticle(idpdg)->Mass());
592   t.SetCharge(fPdg->GetParticle(idpdg)->Charge()/3.);
593
594   // Retrieve the various fit parameters for this track
595   AliSignal* fitdata=fFitdefs->GetIdHit(i);
596   for (Int_t jval=0; jval<fEvent.fresult[i].nval; jval++)
597   {
598    fitdata->SetSignal(fEvent.fresult[i].val[jval],jval+1);
599   }
600
601   // Store the various fit parameters for this track
602   t.SetFitDetails(fitdata);
603
604   // Store the various reco tracks as track hypotheses.
605   // A copy of the first reco track is entered as a new track instance
606   // into the event and all reco tracks (incl. the first one) are
607   // stored as hypotheses linked to this new reco track.
608   if (i==0)
609   {
610    evt->AddTrack(t);
611    AliTrack* tx=evt->GetTrack(evt->GetNtracks());
612    Int_t nrec=evt->GetNtracks(1);
613    tx->SetId(nrec+1);
614   }
615   AliTrack* tx=evt->GetTrack(evt->GetNtracks());
616   if (tx) tx->AddTrackHypothesis(t);
617  }
618 }
619 ///////////////////////////////////////////////////////////////////////////
620 void IceF2k::PutHits(IceEvent* evt)
621 {
622 // Get the hit and waveform info from the F2000 file into the IcePack structure.
623 // This memberfunction is based on the original code by Adam Bouchta.
624
625  if (!evt) return;
626
627  // Loop over all the hits and add them to the current event
628  IceAOM om;
629  AliSignal s;
630  s.SetSlotName("ADC",1);
631  s.SetSlotName("LE",2);
632  s.SetSlotName("TOT",3);
633  Int_t chan=0;
634  Int_t maxchan=800;
635  if (fOmdb) maxchan=fHeader.nch;
636  IceAOM* omx=0;
637  AliSignal* sx=0;
638  Int_t tid=0;
639  AliTrack* tx=0;
640  for (Int_t i=0; i<fEvent.nhits; i++)
641  {
642   chan=fEvent.h[i].ch+1;
643   if (chan>maxchan) continue; // Channels 9001, 9002 etc are trigger channels
644
645   // Get corresponding device from the current event structure  
646   omx=(IceAOM*)evt->GetIdDevice(chan);
647   if (!omx)
648   {
649    if (fOmdb)
650    {
651     omx=(IceAOM*)fOmdb->GetObject(chan,1);
652     evt->AddDevice(omx);
653    }
654    else
655    {
656     om.Reset(1);
657     om.SetUniqueID(chan);
658     evt->AddDevice(om);
659    }
660    omx=(IceAOM*)evt->GetIdDevice(chan);
661   }
662
663   if (!omx) continue;
664
665   s.Reset();
666   s.SetUniqueID(fEvent.h[i].id);
667   s.SetSignal(fEvent.h[i].amp,1);
668   s.SetSignal(fEvent.h[i].t,2);
669   s.SetSignal(fEvent.h[i].tot,3);
670
671   omx->AddHit(s);
672
673   sx=omx->GetHit(omx->GetNhits());
674   if (!sx) continue;
675
676   // Bi-directional link between this hit and the track that caused the ADC value.
677   // This F2K info is probably only present for MC tracks.
678   tid=fEvent.h[i].ma;
679   if (tid > 0)
680   {
681    tx=evt->GetIdTrack(tid); // Reco tracks
682    if (!tx) tx=evt->GetIdTrack(-tid); // MC tracks
683    if (tx) sx->AddLink(tx);
684   }
685   else
686   {
687    if (tid == -2) s.SetNameTitle("N","Noise");
688    if (tid == -3) s.SetNameTitle("A","Afterpulse");
689   }
690  }
691
692  // Loop over all the waveforms and add the histo(s) to the corresponding OM's
693  TH1F histo;
694  Int_t nbins=0;
695  Float_t xlow=0;
696  Float_t xup=0;
697  TString hname;
698  for (Int_t iwf=0; iwf<fEvent.nwf; iwf++)
699  {
700   chan=fEvent.wf[iwf].om;
701   if (chan<=0 || chan>maxchan) continue; // Skip trigger channels
702
703   // Get corresponding device from the current event structure  
704   omx=(IceAOM*)evt->GetIdDevice(chan);
705   if (!omx)
706   {
707    if (fOmdb)
708    {
709     omx=(IceAOM*)fOmdb->GetObject(chan,1);
710     evt->AddDevice(omx);
711    }
712    else
713    {
714     om.Reset(1);
715     om.SetUniqueID(chan);
716     evt->AddDevice(om);
717    }
718    omx=(IceAOM*)evt->GetIdDevice(chan);
719   }
720
721   if (!omx) continue;
722
723   omx->SetSlotName("BASELINE",omx->GetNnames()+1);
724   omx->SetSignal(fEvent.wf[iwf].baseline,"BASELINE");
725
726   // Fill the waveform histogram
727   hname="OM";
728   hname+=chan;
729   hname+="-WF";
730   hname+=omx->GetNwaveforms()+1;
731
732   histo.Reset();
733   histo.SetName(hname.Data());
734   nbins=fEvent.wf[iwf].ndigi;
735   xlow=fEvent.wf[iwf].t_start;
736   xup=xlow+float(nbins)*fEvent.wf[iwf].t_bin;
737   histo.SetBins(nbins,xlow,xup);
738
739   for (Int_t jbin=1; jbin<=fEvent.wf[iwf].ndigi; jbin++)
740   {
741    histo.SetBinContent(jbin,fEvent.wf[iwf].digi[jbin]);
742   }
743
744   omx->SetWaveform(&histo,omx->GetNwaveforms()+1);
745  }
746
747  // Set bi-directional links between hits and reco track hypotheses.
748  // Note : Reco tracks are recognised by their positive id.
749  Int_t hid=0;
750  TObjArray* rectracks=evt->GetTracks(1);
751  for (Int_t jtk=0; jtk<rectracks->GetEntries(); jtk++)
752  {
753   tx=(AliTrack*)rectracks->At(jtk);
754   if (!tx) continue;
755   
756   for (Int_t jhyp=1; jhyp<=tx->GetNhypotheses(); jhyp++)
757   {
758    AliTrack* hypx=tx->GetTrackHypothesis(jhyp);
759    if (!hypx) continue;
760
761    // Loop over all combinations of F2K fits and used OM hits
762    for (Int_t k=0; k<fEvent.nfit_uses; k++)
763    {
764     if (fEvent.fit_uses[k].useid != hypx->GetId()) continue;
765     hid=fEvent.fit_uses[k].hitid;
766     sx=evt->GetIdHit(hid,"IceAOM");
767     if (sx) sx->AddLink(hypx);
768    }
769   }
770  }
771 }
772 ///////////////////////////////////////////////////////////////////////////