11-oct-2005 NvE SetMJD invoked at the end of AliTimestamp::Add instead of just updati...
[u/mrichter/AliRoot.git] / RALICE / AliDevice.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$
17
18 ///////////////////////////////////////////////////////////////////////////
19 // Class AliDevice
20 // Signal (Hit) handling of a generic device.
21 // Basically this class provides a user interface to group and handle
22 // various instances of AliSignal objects, called generically "hits".
23 // An AliDevice object itself has (in addition to hit storage) also the
24 // complete functionality of the class AliSignal.
25 //
26 // Example :
27 // =========
28 //
29 // AliDevice m;
30 // // Set user defined status word to indicate e.g. readout electronics version 
31 // m.SetStatus(100201);
32 // m.SetHitCopy(1);
33 // m.SetName("OM123");
34 //
35 // Float_t pos[3]={1,2,3};
36 // m.SetPosition(pos,"car");
37 //
38 // AliSignal s;
39 //
40 // s.Reset(1);
41 // s.SetName("OM123 Hit 1");
42 // s.SetSlotName("ADC");
43 // s.SetSignal(10);
44 // s.SetSlotName("LE",2);
45 // s.SetSignal(-100,2);
46 // s.SetSlotName("TOT",3);
47 // s.SetSignal(-1000,3);
48 // m.AddHit(s);
49 //
50 // s.Reset(1);
51 // s.SetName("OM123 Hit 2");
52 // s.SetSlotName("ADC");
53 // s.SetSignal(11);
54 // s.SetSlotName("LE",2);
55 // s.SetSignal(-101,2);
56 // s.SetSlotName("TOT",3);
57 // s.SetSignal(1001,3);
58 // m.AddHit(s);
59 //
60 // s.Reset(1);
61 // s.SetName("OM123 Hit 3");
62 // s.SetSlotName("ADC");
63 // s.SetSignal(12);
64 // s.SetSlotName("LE",2);
65 // s.SetSignal(-102,2);
66 // s.SetSlotName("TOT",3);
67 // s.SetSignal(-1002,3);
68 // m.AddHit(s);
69 //
70 // TObjArray* ordered=m.SortHits("TOT");
71 // nhits=ordered->GetEntries();
72 // for (Int_t i=0; i<nhits; i++)
73 // {
74 //  AliSignal* sx=(AliSignal*)ordered->At(i);
75 //  if (sx) sx->Data();
76 // }
77 //
78 //--- Author: Nick van Eijndhoven 23-jun-2004 Utrecht University
79 //- Modified: NvE $Date$ Utrecht University
80 ///////////////////////////////////////////////////////////////////////////
81
82 #include "AliDevice.h"
83 #include "Riostream.h"
84  
85 ClassImp(AliDevice) // Class implementation to enable ROOT I/O
86  
87 AliDevice::AliDevice() : AliSignal()
88 {
89 // Default constructor.
90 // The user definable status word is set to zero.
91 // By default private copies of the recorded hits will be made.
92 // This implies that by default the device will own the registered hits.
93 // See the SetHitCopy() memberfunction for further details.
94  fStatus=0;
95  fHitCopy=1;
96  fHits=0;
97  fOrdered=0;
98  fMarkers=0;
99 }
100 ///////////////////////////////////////////////////////////////////////////
101 AliDevice::~AliDevice()
102 {
103 // Default destructor.
104
105  // Remove backward links to this device from the hits
106  // which were not owned by it.
107  if (!fHitCopy)
108  {
109   for (Int_t ih=1; ih<=GetNhits(); ih++)
110   {
111    AliSignal* sx=GetHit(ih);
112    if (sx) sx->ResetLinks(this);
113   }
114  }
115
116  if (fHits)
117  {
118   delete fHits;
119   fHits=0;
120  }
121
122  if (fOrdered)
123  {
124   delete fOrdered;
125   fOrdered=0;
126  }
127
128  if (fMarkers)
129  {
130   delete fMarkers;
131   fMarkers=0;
132  }
133 }
134 ///////////////////////////////////////////////////////////////////////////
135 AliDevice::AliDevice(const AliDevice& dev) : AliSignal(dev)
136 {
137 // Copy constructor.
138
139  fHits=0;
140  fOrdered=0;
141  fMarkers=0;
142
143  fStatus=dev.GetStatus();
144  fHitCopy=dev.GetHitCopy();
145
146  Int_t nhits=dev.GetNhits();
147  if (nhits)
148  {
149   fHits=new TObjArray(nhits);
150   if (fHitCopy) fHits->SetOwner();
151   for (Int_t ih=1; ih<=nhits; ih++)
152   {
153    AliSignal* sx=dev.GetHit(ih);
154    if (fHitCopy)
155    {
156     fHits->Add(sx->Clone());
157     AliSignal* s=(AliSignal*)fHits->Last();
158     s->ResetLinks((AliDevice*)&dev);
159     s->SetDevice(this);
160    }
161    else
162    {
163     sx->AddLink(this);
164     fHits->Add(sx);
165    }
166   }
167  }
168 }
169 ///////////////////////////////////////////////////////////////////////////
170 void AliDevice::Reset(Int_t mode)
171 {
172 // Reset registered hits and AliSignal attributes.
173 // Note : The status word and HitCopy flag are NOT modified.
174 //        Use SetStatus() and SetHitCopy() to modify these parameters. 
175 // See AliSignal::Reset() for further details.
176  RemoveHits();
177  AliSignal::Reset(mode);
178 }
179 ///////////////////////////////////////////////////////////////////////////
180 void AliDevice::SetHitCopy(Int_t j)
181 {
182 // (De)activate the creation of private copies of the AliSignals added as hits.
183 // j=0 ==> No private copies are made; pointers of original hits are stored.
184 // j=1 ==> Private copies of the hits are made and these pointers are stored.
185 //
186 // Note : Once the storage contains pointer(s) to hit(s) one cannot
187 //        change the HitCopy mode anymore.
188 //        To change the HitCopy mode for an existing AliDevice containing
189 //        hits one first has to invoke either RemoveHits() or Reset().
190  if (!fHits)
191  {
192   if (j==0 || j==1)
193   {
194    fHitCopy=j;
195   }
196   else
197   {
198    cout << "*AliDevice::SetHitCopy* Invalid argument : " << j << endl;
199   }
200  }
201  else
202  {
203   cout << "*AliDevice::SetHitCopy* Storage already contained hits."
204        << "  ==> HitCopy mode not changed." << endl; 
205  }
206 }
207 ///////////////////////////////////////////////////////////////////////////
208 Int_t AliDevice::GetHitCopy() const
209 {
210 // Provide value of the HitCopy mode.
211 // 0 ==> No private copies are made; pointers of original hits are stored.
212 // 1 ==> Private copies of the hits are made and these pointers are stored.
213  return fHitCopy;
214 }
215 ///////////////////////////////////////////////////////////////////////////
216 void AliDevice::SetStatus(Int_t word)
217 {
218 // Set a user defined status word for this device.
219  fStatus=word;
220 }
221 ///////////////////////////////////////////////////////////////////////////
222 Int_t AliDevice::GetStatus() const
223 {
224 // Provide the user defined status word for this device.
225  return fStatus;
226 }
227 ///////////////////////////////////////////////////////////////////////////
228 void AliDevice::AddHit(AliSignal& s)
229 {
230 // Register an AliSignal object as a hit to this device.
231 // Note : In case this device owns the AliSignal object, the pointer to
232 //        this device will be stored in the special owning device
233 //        pointer of the AliSignal object.
234 //        In case this device does not own the AliSignal object, a (backward)
235 //        link to this device is added to the first slot of the AliSignal
236 //        if there was no link to this device already present.
237
238  if (!fHits)
239  {
240   fHits=new TObjArray(1);
241   if (fHitCopy) fHits->SetOwner();
242  }
243
244  // Check if this signal is already stored for this device.
245  Int_t nhits=GetNhits();
246  for (Int_t i=0; i<nhits; i++)
247  {
248   if (&s==fHits->At(i)) return; 
249  }
250
251  // Check for existing (backward) link to this device.
252  Int_t nlinks=s.GetNlinks(this);
253
254  if (fHitCopy)
255  {
256   fHits->Add(s.Clone());
257   // Remove unnecessary backward link(s) from the various slots
258   // and set the owning link to this device
259   AliSignal* sx=(AliSignal*)fHits->Last();
260   if (nlinks) sx->ResetLinks(this);
261   sx->SetDevice(this);
262  }
263  else
264  {
265   fHits->Add(&s);
266   // Set (backward) link to the this device
267   if (!nlinks) s.AddLink(this);
268  }
269 }
270 ///////////////////////////////////////////////////////////////////////////
271 void AliDevice::RemoveHit(AliSignal& s)
272 {
273 // Remove AliSignal object registered as a hit from this device.
274  if (fHits)
275  {
276   AliSignal* test=(AliSignal*)fHits->Remove(&s);
277   if (test)
278   {
279    fHits->Compress();
280    if (fHitCopy) delete test;
281   }
282  }
283  if (fOrdered)
284  {
285   AliSignal* test=(AliSignal*)fOrdered->Remove(&s);
286   if (test) fOrdered->Compress();
287  }
288 }
289 ///////////////////////////////////////////////////////////////////////////
290 void AliDevice::RemoveHits()
291 {
292 // Remove all AliSignal objects registered as hits from this device.
293  if (fHits)
294  {
295   delete fHits;
296   fHits=0;
297  }
298  if (fOrdered)
299  {
300   delete fOrdered;
301   fOrdered=0;
302  }
303  if (fMarkers)
304  {
305   delete fMarkers;
306   fMarkers=0;
307  }
308 }
309 ///////////////////////////////////////////////////////////////////////////
310 Int_t AliDevice::GetNhits() const
311 {
312 // Provide the number of registered hits for this device.
313  Int_t nhits=0;
314  if (fHits) nhits=fHits->GetEntries();
315  return nhits;
316 }
317 ///////////////////////////////////////////////////////////////////////////
318 AliSignal* AliDevice::GetHit(Int_t j) const
319 {
320 // Provide the AliSignal object registered as hit number j.
321 // Note : j=1 denotes the first hit.
322  if (!fHits) return 0;
323
324  if ((j >= 1) && (j <= GetNhits()))
325  {
326   return (AliSignal*)fHits->At(j-1);
327  }
328  else
329  {
330   return 0;
331  }
332 }
333 ///////////////////////////////////////////////////////////////////////////
334 AliSignal* AliDevice::GetIdHit(Int_t id) const
335 {
336 // Return the hit with unique identifier "id".
337  if (!fHits || id<0) return 0;
338
339  AliSignal* sx=0;
340  Int_t sid=0;
341  for (Int_t i=0; i<GetNhits(); i++)
342  {
343   sx=(AliSignal*)fHits->At(i);
344   if (sx)
345   {
346    sid=sx->GetUniqueID();
347    if (id==sid) return sx;
348   }
349  }
350  return 0; // No matching id found
351 }
352 ///////////////////////////////////////////////////////////////////////////
353 TObjArray* AliDevice::GetHits()
354 {
355 // Provide the references to all the registered hits.
356  return fHits;
357 }
358 ///////////////////////////////////////////////////////////////////////////
359 void AliDevice::ShowHit(Int_t j) const
360 {
361 // Show data of the registered j-th hit.
362 // If j=0 all associated hits will be shown.
363 // The default is j=0.
364  if (!j)
365  {
366   Int_t nhits=GetNhits();
367   for (Int_t ih=1; ih<=nhits; ih++)
368   {
369    AliSignal* sx=GetHit(ih);
370    if (sx) sx->Data();
371   }
372  }
373  else
374  {
375   AliSignal* s=GetHit(j);
376   if (s) s->Data();
377  }
378 }
379 ///////////////////////////////////////////////////////////////////////////
380 void AliDevice::Data(TString f) const
381 {
382 // Print the device and all registered hit info according to the specified
383 // coordinate frame.
384  AliSignal::Data(f);
385  Int_t nhits=GetNhits();
386  if (nhits)
387  {
388   cout << " The following " << nhits << " hits are registered : " << endl;
389   ShowHit();
390  }
391  else
392  {
393   cout << " No hits have been registered for this device." << endl;
394  }
395 }
396 ///////////////////////////////////////////////////////////////////////////
397 void AliDevice::GetExtremes(Float_t& vmin,Float_t& vmax,Int_t idx,TObjArray* hits,Int_t mode) const
398 {
399 // Provide the min. and max. signal values of an array of hits.
400 // The input argument "idx" denotes the index of the signal slots to be investigated.
401 // The default is idx=1;
402 // In case hits=0 (default), the registered hits of the current device are used. 
403 // Signals which were declared as "Dead" will be rejected.
404 // The gain etc... corrected signals will be used in the process as specified
405 // by the  "mode" argument. The definition of this "mode" parameter corresponds to
406 // the description provided in the GetSignal memberfunction of class AliSignal.
407 // The default is mode=1 (for backward compatibility reasons).
408
409  vmin=0;
410  vmax=0;
411
412  if (!hits) hits=fHits;
413  
414  if (idx<=0 || !hits) return;
415
416  Int_t nhits=hits->GetEntries();
417
418  Float_t sig=0;
419  for (Int_t i=0; i<nhits; i++)
420  {
421   AliSignal* sx=(AliSignal*)hits->At(i);
422
423   if (!sx) continue;
424   if (idx > sx->GetNvalues()) continue; // User specified slotindex out of range for this signal
425   if (sx->GetDeadValue(idx)) continue;  // Only take alive signals
426
427   sig=sx->GetSignal(idx,mode);
428   if (i==0)
429   {
430    vmin=sig;
431    vmax=sig;
432   }
433   else
434   {
435    if (sig<vmin) vmin=sig;
436    if (sig>vmax) vmax=sig;
437   }
438  }
439 }
440 ///////////////////////////////////////////////////////////////////////////
441 void AliDevice::GetExtremes(Float_t& vmin,Float_t& vmax,TString name,TObjArray* hits,Int_t mode) const
442 {
443 // Provide the min. and max. signal values of an array of hits.
444 // The input argument "name" denotes the name of the signal slots to be investigated.
445 // In case hits=0 (default), the registered hits of the current device are used. 
446 // Signals which were declared as "Dead" will be rejected.
447 // The gain etc... corrected signals will be used in the process as specified
448 // by the  "mode" argument. The definition of this "mode" parameter corresponds to
449 // the description provided in the GetSignal memberfunction of class AliSignal.
450 // The default is mode=1 (for backward compatibility reasons).
451
452  vmin=0;
453  vmax=0;
454
455  if (!hits) hits=fHits;
456  
457  if (!hits) return;
458
459  Int_t nhits=hits->GetEntries();
460
461  Int_t idx=0; // The signal slotindex to perform the sorting on
462
463  Float_t sig=0;
464  for (Int_t i=0; i<nhits; i++)
465  {
466   AliSignal* sx=(AliSignal*)hits->At(i);
467
468   if (!sx) continue;
469
470   // Obtain the slotindex corresponding to the user selection
471   idx=sx->GetSlotIndex(name);
472   if (!idx) continue;
473
474   if (sx->GetDeadValue(idx)) continue; // Only take alive signals
475
476   sig=sx->GetSignal(idx,mode);
477   if (i==0)
478   {
479    vmin=sig;
480    vmax=sig;
481   }
482   else
483   {
484    if (sig<vmin) vmin=sig;
485    if (sig>vmax) vmax=sig;
486   }
487  }
488 }
489 ///////////////////////////////////////////////////////////////////////////
490 TObjArray* AliDevice::SortHits(Int_t idx,Int_t mode,TObjArray* hits,Int_t mcal)
491 {
492 // Order the references to an array of hits by looping over the input array "hits"
493 // and checking the signal value. The ordered array is returned as a TObjArray.
494 // In case hits=0 (default), the registered hits of the current device are used. 
495 // Note that the original hit array is not modified.
496 // A "hit" represents an abstract object which is derived from AliSignal.
497 // The user can specify the index of the signal slot to perform the sorting on.
498 // By default the slotindex will be 1.
499 // Via the "mode" argument the user can specify ordering in decreasing
500 // order (mode=-1) or ordering in increasing order (mode=1).
501 // The default is mode=-1.
502 // Signals which were declared as "Dead" will be rejected.
503 // The gain etc... corrected signals will be used in the ordering process as
504 // specified by the "mcal" argument. The definition of this "mcal" parameter
505 // corresponds to the signal correction mode described in the GetSignal
506 // memberfunction of class AliSignal.
507 // The default is mcal=1 (for backward compatibility reasons).
508
509  if (fOrdered)
510  {
511   delete fOrdered;
512   fOrdered=0;
513  }
514
515  if (!hits) hits=fHits;
516  
517  if (idx<=0 || abs(mode)!=1 || !hits) return fOrdered;
518
519  Int_t nhits=hits->GetEntries();
520  if (!nhits)
521  {
522   return fOrdered;
523  }
524  else
525  {
526   fOrdered=new TObjArray(nhits);
527  }
528
529  Int_t nord=0;
530  for (Int_t i=0; i<nhits; i++) // Loop over all hits of the array
531  {
532   AliSignal* s=(AliSignal*)hits->At(i);
533
534   if (!s) continue;
535
536   if (idx > s->GetNvalues()) continue; // User specified slotindex out of range for this signal
537   if (s->GetDeadValue(idx)) continue;  // Only take alive signals
538  
539   if (nord == 0) // store the first hit with a signal at the first ordered position
540   {
541    nord++;
542    fOrdered->AddAt(s,nord-1);
543    continue;
544   }
545  
546   for (Int_t j=0; j<=nord; j++) // put hit in the right ordered position
547   {
548    if (j == nord) // module has smallest (mode=-1) or largest (mode=1) signal seen so far
549    {
550     nord++;
551     fOrdered->AddAt(s,j); // add hit at the end
552     break; // go for next hit
553    }
554  
555    if (mode==-1 && s->GetSignal(idx,mcal) < ((AliSignal*)fOrdered->At(j))->GetSignal(idx,mcal)) continue;
556    if (mode==1 && s->GetSignal(idx,mcal) > ((AliSignal*)fOrdered->At(j))->GetSignal(idx,mcal)) continue;
557  
558    nord++;
559    for (Int_t k=nord-1; k>j; k--) // create empty position
560    {
561     fOrdered->AddAt(fOrdered->At(k-1),k);
562    }
563    fOrdered->AddAt(s,j); // put hit at empty position
564    break; // go for next matrix module
565   }
566  }
567  return fOrdered;
568 }
569 ///////////////////////////////////////////////////////////////////////////
570 TObjArray* AliDevice::SortHits(TString name,Int_t mode,TObjArray* hits,Int_t mcal)
571 {
572 // Order the references to an array of hits by looping over the input array "hits"
573 // and checking the signal value. The ordered array is returned as a TObjArray.
574 // In case hits=0 (default), the registered hits of the current device are used. 
575 // Note that the input array is not modified.
576 // A "hit" represents an abstract object which is derived from AliSignal.
577 // The user can specify the name of the signal slot to perform the sorting on.
578 // In case no matching slotname is found, the signal will be skipped.
579 // Via the "mode" argument the user can specify ordering in decreasing
580 // order (mode=-1) or ordering in increasing order (mode=1).
581 // The default is mode=-1.
582 // Signals which were declared as "Dead" will be rejected.
583 // The gain etc... corrected signals will be used in the ordering process as
584 // specified by the "mcal" argument. The definition of this "mcal" parameter
585 // corresponds to the signal correction mode described in the GetSignal
586 // memberfunction of class AliSignal.
587 // The default is mcal=1 (for backward compatibility reasons).
588
589  if (fOrdered)
590  {
591   delete fOrdered;
592   fOrdered=0;
593  }
594
595  if (!hits) hits=fHits;
596  
597  if (abs(mode)!=1 || !hits) return fOrdered;
598
599  Int_t nhits=hits->GetEntries();
600  if (!nhits)
601  {
602   return fOrdered;
603  }
604  else
605  {
606   fOrdered=new TObjArray(nhits);
607  }
608
609  Int_t idx=0; // The signal slotindex to perform the sorting on
610
611  Int_t nord=0;
612  for (Int_t i=0; i<nhits; i++) // loop over all hits of the array
613  {
614   AliSignal* s=(AliSignal*)hits->At(i);
615
616   if (!s) continue;
617
618   // Obtain the slotindex corresponding to the user selection
619   idx=s->GetSlotIndex(name);
620   if (!idx) continue;
621
622   if (s->GetDeadValue(idx)) continue; // only take alive signals
623  
624   if (nord == 0) // store the first hit with a signal at the first ordered position
625   {
626    nord++;
627    fOrdered->AddAt(s,nord-1);
628    continue;
629   }
630  
631   for (Int_t j=0; j<=nord; j++) // put hit in the right ordered position
632   {
633    if (j == nord) // module has smallest (mode=-1) or largest (mode=1) signal seen so far
634    {
635     nord++;
636     fOrdered->AddAt(s,j); // add hit at the end
637     break; // go for next hit
638    }
639  
640    if (mode==-1 && s->GetSignal(idx,mcal) < ((AliSignal*)fOrdered->At(j))->GetSignal(idx,mcal)) continue;
641    if (mode==1 && s->GetSignal(idx,mcal) > ((AliSignal*)fOrdered->At(j))->GetSignal(idx,mcal)) continue;
642  
643    nord++;
644    for (Int_t k=nord-1; k>j; k--) // create empty position
645    {
646     fOrdered->AddAt(fOrdered->At(k-1),k);
647    }
648    fOrdered->AddAt(s,j); // put hit at empty position
649    break; // go for next matrix module
650   }
651  }
652  return fOrdered;
653 }
654 ///////////////////////////////////////////////////////////////////////////
655 void AliDevice::DisplayHits(Int_t idx,Float_t scale,TObjArray* hits,Int_t dp,Int_t mode,Int_t mcol)
656 {
657 // 3D color display of an array hits.
658 // The user can specify the index (default=1) of the signal slot to perform the display for.
659 // The marker size will indicate the absolute value of the signal (specified by the slotindex)
660 // as a percentage of the input argument "scale".
661 // In case scale<0 the maximum absolute signal value encountered in the hit array will be used
662 // to define the 100% scale. The default is scale=-1.
663 // In case hits=0 (default), the registered hits of the current device are used. 
664 // Note that the input array is not modified.
665 // In case dp=1 the device position will be used, otherwise the hit position will
666 // be used in the display. The default is dp=0.
667 // Via the "mcol" argument the user can specify the marker color (see TPolyMarker3D).
668 // The default is mcol=blue.
669 // Signals which were declared as "Dead" will not be displayed.
670 // The gain etc... corrected signals will be used to determine the marker size.
671 // The gain correction is performed according to "mode" argument. The definition of this
672 // "mode" parameter corresponds to the description provided in the GetSignal
673 // memberfunction of class AliSignal.
674 // The default is mode=1 (for backward compatibility reasons).
675 //
676 // Note :
677 // ------
678 // Before any display activity, a TCanvas and a TView have to be initiated
679 // first by the user like for instance
680 // 
681 // TCanvas* c1=new TCanvas("c1","c1");
682 // TView* view=new TView(1);
683 // view->SetRange(-1000,-1000,-1000,1000,1000,1000);
684 // view->ShowAxis();
685
686  Int_t thisdev=0; // Indicate whether this is the owning device or not 
687  if (!hits)
688  {
689   hits=fHits;
690   thisdev=1;
691  }
692  
693  if (idx<=0 || !hits) return;
694
695  Int_t nhits=hits->GetEntries();
696  if (!nhits) return;
697
698  Float_t sigmax=fabs(scale);
699  if (scale<0)
700  {
701   Float_t vmin,vmax;
702   GetExtremes(vmin,vmax,idx,hits);
703   sigmax=fabs(vmax);
704   if (fabs(vmin)>sigmax) sigmax=fabs(vmin);
705  }
706
707  if (sigmax <=0) return;
708
709  if (fMarkers)
710  {
711   delete fMarkers;
712   fMarkers=0;
713  }
714  fMarkers=new TObjArray(nhits);
715  fMarkers->SetOwner();
716
717  Float_t pos[3];
718  GetPosition(pos,"car");
719
720  Float_t sig=0;
721  for (Int_t ih=0; ih<nhits; ih++)
722  {
723   AliSignal* sx=(AliSignal*)hits->At(ih);
724   if (!sx) continue;
725   if (!dp)
726   {
727    sx->GetPosition(pos,"car");
728   }
729   else
730   {
731    if (!thisdev)
732    {
733     AliDevice* dev=sx->GetDevice();
734     if (dev) dev->GetPosition(pos,"car");
735    }
736   }
737   sig=sx->GetSignal(idx,mode);
738   TPolyMarker3D* m=new TPolyMarker3D();
739   m->SetMarkerStyle(8);
740   m->SetMarkerColor(mcol);
741   m->SetMarkerSize(100.*fabs(sig)/sigmax);
742   m->SetPoint(0,pos[0],pos[1],pos[2]);
743   fMarkers->Add(m);
744   m->Draw();
745  }
746 }
747 ///////////////////////////////////////////////////////////////////////////
748 void AliDevice::DisplayHits(TString name,Float_t scale,TObjArray* hits,Int_t dp,Int_t mode,Int_t mcol)
749 {
750 // 3D color display of an array hits.
751 // The user can specify the name of the signal slot to perform the display for.
752 // The marker size will indicate the absolute value of the signal (specified by the slotname)
753 // as a percentage of the input argument "scale".
754 // In case scale<0 the maximum absolute signal value encountered in the hit array will be used
755 // to define the 100% scale. The default is scale=-1.
756 // In case hits=0 (default), the registered hits of the current device are used. 
757 // Note that the input array is not modified.
758 // In case dp=1 the device position will be used, otherwise the hit position will
759 // be used in the display. The default is dp=0.
760 // The marker size will indicate the percentage of the maximum encountered value
761 // of the absolute value of the name-specified input signal slots.
762 // Via the "mcol" argument the user can specify the marker color (see TPolyMarker3D).
763 // The default is mcol=blue.
764 // Signals which were declared as "Dead" will not be displayed.
765 // The gain etc... corrected signals will be used to determine the marker size.
766 // The gain correction is performed according to "mode" argument. The definition of this
767 // "mode" parameter corresponds to the description provided in the GetSignal
768 // memberfunction of class AliSignal.
769 // The default is mode=1 (for backward compatibility reasons).
770 //
771 // Note :
772 // ------
773 // Before any display activity, a TCanvas and a TView have to be initiated
774 // first by the user like for instance
775 // 
776 // TCanvas* c1=new TCanvas("c1","c1");
777 // TView* view=new TView(1);
778 // view->SetRange(-1000,-1000,-1000,1000,1000,1000);
779 // view->ShowAxis();
780
781  Int_t thisdev=0; // Indicate whether this is the owning device or not 
782  if (!hits)
783  {
784   hits=fHits;
785   thisdev=1;
786  }
787  
788  if (!hits) return;
789
790  Int_t nhits=hits->GetEntries();
791
792  if (!nhits) return;
793
794  Float_t sigmax=fabs(scale);
795  if (scale<0)
796  {
797   Float_t vmin,vmax;
798   GetExtremes(vmin,vmax,name,hits);
799   sigmax=fabs(vmax);
800   if (fabs(vmin)>sigmax) sigmax=fabs(vmin);
801  }
802
803  if (sigmax <=0) return;
804
805  if (fMarkers)
806  {
807   delete fMarkers;
808   fMarkers=0;
809  }
810  fMarkers=new TObjArray(nhits);
811  fMarkers->SetOwner();
812
813  Float_t pos[3];
814  GetPosition(pos,"car");
815
816  Int_t idx=0; // The slot index corresponding to the user specified name
817  Float_t sig=0;
818  for (Int_t ih=0; ih<nhits; ih++)
819  {
820   AliSignal* sx=(AliSignal*)hits->At(ih);
821   if (!sx) continue;
822   idx=sx->GetSlotIndex(name);
823   if (!idx) continue;
824   if (!dp)
825   {
826    sx->GetPosition(pos,"car");
827   }
828   else
829   {
830    if (!thisdev)
831    {
832     AliDevice* dev=sx->GetDevice();
833     if (dev) dev->GetPosition(pos,"car");
834    }
835   }
836   sig=sx->GetSignal(idx,mode);
837   TPolyMarker3D* m=new TPolyMarker3D();
838   m->SetMarkerStyle(8);
839   m->SetMarkerColor(mcol);
840   m->SetMarkerSize(100.*fabs(sig)/sigmax);
841   m->SetPoint(0,pos[0],pos[1],pos[2]);
842   fMarkers->Add(m);
843   m->Draw();
844  }
845 }
846 ///////////////////////////////////////////////////////////////////////////
847 TObject* AliDevice::Clone(const char* name) const
848 {
849 // Make a deep copy of the current object and provide the pointer to the copy.
850 // This memberfunction enables automatic creation of new objects of the
851 // correct type depending on the object type, a feature which may be very useful
852 // for containers like AliEvent when adding objects in case the
853 // container owns the objects. This feature allows e.g. AliEvent
854 // to store either AliDevice objects or objects derived from AliDevice
855 // via tha AddDevice memberfunction, provided these derived classes also have
856 // a proper Clone memberfunction. 
857
858  AliDevice* dev=new AliDevice(*this);
859  if (name)
860  {
861   if (strlen(name)) dev->SetName(name);
862  }
863  return dev;
864 }
865 ///////////////////////////////////////////////////////////////////////////