PHOS/AliPHOS.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 /*
  17 $Log$
  18 Revision 1.9  1999/11/08 07:12:31  fca
  19 Minor corrections thanks to I.Hrivnacova
  20
  21 Revision 1.8  1999/09/29 09:24:23  fca
  22 Introduction of the Copyright and cvs Log
  23
  24 */
  25
  26 ////////////////////////////////////////////////
  27 //  Manager and hits classes for set:PHOS     //
  28 ////////////////////////////////////////////////
  29
  30 // --- ROOT system ---
  31 #include "TH1.h"
  32 #include "TRandom.h"
  33 #include "TFile.h"
  34 #include "TTree.h"
  35 #include "TBRIK.h"
  36 #include "TNode.h"
  37 #include "TMath.h"
  38
  39 // --- Standard library ---
  40 #include <stdio.h>
  41 #include <string.h>
  42 #include <stdlib.h>
  43 #include <iostream.h>
  44
  45 // --- galice header files ---
  46 #include "AliPHOS.h"
  47 #include "AliRun.h"
  48
  49 //______________________________________________________________________________
  50
  51
  52 ClassImp(AliPHOS)
  53
  54 //______________________________________________________________________________
  55
  56 AliPHOS::~AliPHOS(void)
  57 {
  58   delete fHits;                 // 28.12.1998
  59   delete fTreePHOS;             // 28.12.1998
  60   fCradles->Delete();
  61   delete fCradles;
  62 }
  63
  64 //______________________________________________________________________________
  65
  66 AliPHOS::AliPHOS() :
  67          fDebugLevel            (0),
  68          fTreePHOS              (NULL),
  69          fBranchNameOfCradles   ("AliPHOSCradles"),
  70          fTreeName              ("PHOS")
  71 {
  72    fIshunt   = 0;
  73
  74   if( NULL==(fCradles=new TObjArray) )
  75   {
  76     Error("AliPHOS","Can not create fCradles");
  77     exit(1);
  78   }
  79   DefPars();
  80 }
  81
  82 //______________________________________________________________________________
  83
  84 AliPHOS::AliPHOS(const char *name, const char *title)
  85        : AliDetector            (name,title),
  86          fDebugLevel            (0),
  87          fTreePHOS              (NULL),
  88          fBranchNameOfCradles   ("AliPHOSCradles"),
  89          fTreeName              ("PHOS")
  90 {
  91 //Begin_Html
  92 /*
  93 <img src="picts/aliphos.gif">
  94 */
  95 //End_Html
  96
  97    fHits   = new TClonesArray("AliPHOShit",  405);
  98
  99    fIshunt     =  0;
 100
 101    SetMarkerColor(kGreen);
 102    SetMarkerStyle(2);
 103    SetMarkerSize(0.4);
 104
 105   if( NULL==(fCradles=new TObjArray) ) {
 106      Error("AliPHOS","Can not create fCradles");
 107      exit(1);
 108   }
 109   DefPars();
 110 }
 111
 112 //______________________________________________________________________________
 113
 114 void AliPHOS::DefPars()
 115 {
 116       PHOSflags[0]=0;
 117       PHOSflags[1]=1;
 118       PHOSflags[2]=0;
 119       PHOSflags[3]=0;
 120       PHOSflags[4]=0;
 121       PHOSflags[5]=0;
 122       PHOSflags[6]=0;
 123       PHOSflags[7]=0;
 124       PHOSflags[8]=0;
 125       PHOScell[0]=2.2;
 126       PHOScell[1]=18.;
 127       PHOScell[2]=0.01;
 128       PHOScell[3]=0.01;
 129       PHOScell[4]=1.0;
 130       PHOScell[5]=0.1;
 131       PHOScell[6]=0.;
 132       PHOScell[7]=0.;
 133       PHOScell[8]=0.;
 134       PHOSradius=460.;
 135       PHOSsize[0]=104;
 136       PHOSsize[1]=88;
 137       PHOSsize[2]=4;
 138       PHOScradlesA=0.;
 139       PHOSextra[0]=0.001;
 140       PHOSextra[1]=6.95;
 141       PHOSextra[2]=4.;
 142       PHOSextra[3]=5.;
 143       PHOSextra[4]=2.;
 144       PHOSextra[5]=0.06;
 145       PHOSextra[6]=10.;
 146       PHOSextra[7]=3.;
 147       PHOSextra[8]=1.;
 148       PHOSTXW[0]=209.;
 149       PHOSTXW[1]=71.;
 150       PHOSTXW[2]=250.;
 151       PHOSAIR[0]=206.;
 152       PHOSAIR[1]=66.;
 153       PHOSAIR[2]=244.;
 154       PHOSFTI[0]=214.6;
 155       PHOSFTI[1]=80.;
 156       PHOSFTI[2]=260.;
 157       PHOSFTI[3]=467.;
 158 }
 159 //______________________________________________________________________________
 160
 161 void AliPHOS::AddHit(Int_t track, Int_t *vol, Float_t *hits)
 162 {
 163   TClonesArray &lhits = *fHits;
 164   new(lhits[fNhits++]) AliPHOShit(fIshunt,track,vol,hits);
 165 }
 166
 167 //___________________________________________
 168 void AliPHOS::BuildGeometry()
 169 {
 170
 171   TNode *Node, *Top;
 172
 173   const int kColorPHOS = kRed;
 174   //
 175   Top=gAlice->GetGeometry()->GetNode("alice");
 176
 177
 178   // PHOS
 179   Float_t pphi=12.9399462;
 180   new TRotMatrix("rot988","rot988",90,-3*pphi,90,90-3*pphi,0,0);
 181   new TRotMatrix("rot989","rot989",90,-  pphi,90,90-  pphi,0,0);
 182   new TRotMatrix("rot990","rot990",90,   pphi,90,90+  pphi,0,0);
 183   new TRotMatrix("rot991","rot991",90, 3*pphi,90,90+3*pphi,0,0);
 184   new TBRIK("S_PHOS","PHOS box","void",107.3,40,130);
 185   Top->cd();
 186   Node = new TNode("PHOS1","PHOS1","S_PHOS",-317.824921,-395.014343,0,"rot988");
 187   Node->SetLineColor(kColorPHOS);
 188   fNodes->Add(Node);
 189   Top->cd();
 190   Node = new TNode("PHOS2","PHOS2","S_PHOS",-113.532333,-494.124908,0,"rot989");
 191   fNodes->Add(Node);
 192   Node->SetLineColor(kColorPHOS);
 193   Top->cd();
 194   Node = new TNode("PHOS3","PHOS3","S_PHOS", 113.532333,-494.124908,0,"rot990");
 195   Node->SetLineColor(kColorPHOS);
 196   fNodes->Add(Node);
 197   Top->cd();
 198   Node = new TNode("PHOS4","PHOS4","S_PHOS", 317.824921,-395.014343,0,"rot991");
 199   Node->SetLineColor(kColorPHOS);
 200   fNodes->Add(Node);
 201 }
 202
 203 //___________________________________________
 204 void AliPHOS::CreateMaterials()
 205 {
 206 // *** DEFINITION OF AVAILABLE PHOS MATERIALS ***
 207
 208 // CALLED BY : PHOS_MEDIA
 209 // ORIGIN    : NICK VAN EIJNDHOVEN
 210
 211
 212
 213     Int_t   ISXFLD = gAlice->Field()->Integ();
 214     Float_t SXMGMX = gAlice->Field()->Max();
 215
 216 // --- The PbWO4 crystals ---
 217     Float_t ax[3] = { 207.19,183.85,16. };
 218     Float_t zx[3] = { 82.,74.,8. };
 219     Float_t wx[3] = { 1.,1.,4. };
 220     Float_t dx    = 8.28;
 221 // --- Stainless Steel ---
 222     Float_t as[5] = { 55.847,12.011,51.9961,58.69,28.0855 };
 223     Float_t zs[5] = { 26.,6.,24.,28.,14. };
 224     Float_t ws[5] = { .6392,8e-4,.2,.14,.02 };
 225     Float_t ds    = 8.;
 226 // --- The polysterene scintillator (CH) ---
 227     Float_t ap[2] = { 12.011,1.00794 };
 228     Float_t zp[2] = { 6.,1. };
 229     Float_t wp[2] = { 1.,1. };
 230     Float_t dp    = 1.032;
 231 // --- Tyvek (CnH2n)
 232     Float_t at[2] = { 12.011,1.00794 };
 233     Float_t zt[2] = { 6.,1. };
 234     Float_t wt[2] = { 1.,2. };
 235     Float_t dt    = .331;
 236 // --- Polystyrene foam ---
 237     Float_t af[2] = { 12.011,1.00794 };
 238     Float_t zf[2] = { 6.,1. };
 239     Float_t wf[2] = { 1.,1. };
 240     Float_t df    = .12;
 241 //--- Foam thermo insulation (actual chemical composition unknown yet!) ---
 242     Float_t ati[2] = { 12.011,1.00794 };
 243     Float_t zti[2] = { 6.,1. };
 244     Float_t wti[2] = { 1.,1. };
 245     Float_t dti    = .1;
 246 // --- Textolit (actual chemical composition unknown yet!) ---
 247     Float_t atx[2] = { 12.011,1.00794 };
 248     Float_t ztx[2] = { 6.,1. };
 249     Float_t wtx[2] = { 1.,1. };
 250     Float_t dtx    = 1.83;
 251
 252     Int_t *idtmed = fIdtmed->GetArray()-699;
 253
 254     AliMixture(  0, "PbWO4$",          ax, zx, dx, -3, wx);
 255     AliMixture(  1, "Polystyrene$",    ap, zp, dp, -2, wp);
 256     AliMaterial( 2, "Al$",             26.98, 13., 2.7, 8.9, 999);
 257 // ---                                Absorption length^ is ignored ---
 258     AliMixture(  3, "Tyvek$",           at, zt, dt, -2, wt);
 259     AliMixture(  4, "Foam$",            af, zf, df, -2, wf);
 260     AliMixture(  5, "Stainless Steel$", as, zs, ds, 5, ws);
 261     AliMaterial( 6, "Si$",              28.09, 14., 2.33, 9.36, 42.3);
 262     AliMixture(  7, "Thermo Insul.$",   ati, zti, dti, -2, wti);
 263     AliMixture(  8, "Textolit$",        atx, ztx, dtx, -2, wtx);
 264     AliMaterial(99, "Air$",             14.61, 7.3, .001205, 30420., 67500);
 265
 266     AliMedium(0, "PHOS Xtal    $", 0, 1, ISXFLD, SXMGMX, 10., .1, .1, .1, .1);
 267     AliMedium(2, "Al parts     $", 2, 0, ISXFLD, SXMGMX, 10., .1, .1, .001, .001);
 268     AliMedium(3, "Tyvek wrapper$", 3, 0, ISXFLD, SXMGMX, 10., .1, .1, .001, .001);
 269     AliMedium(4, "Polyst. foam $", 4, 0, ISXFLD, SXMGMX, 10., .1, .1, .1, .1);
 270     AliMedium(5, "Steel cover  $", 5, 0, ISXFLD, SXMGMX, 10., .1, .1, 1e-4, 1e-4);
 271     AliMedium(6, "Si PIN       $", 6, 0, ISXFLD, SXMGMX, 10., .1, .1, .01, .01);
 272     AliMedium(7, "Thermo Insul.$", 7, 0, ISXFLD, SXMGMX, 10., .1, .1, .1, .1);
 273     AliMedium(8, "Textolit     $", 8, 0, ISXFLD, SXMGMX, 10., .1, .1, .1, .1);
 274     AliMedium(99, "Air          $",99, 0, ISXFLD, SXMGMX, 10., 1., .1, .1, 10);
 275
 276 // --- Generate explicitly delta rays in the steel cover ---
 277     gMC->Gstpar(idtmed[704], "LOSS", 3.);
 278     gMC->Gstpar(idtmed[704], "DRAY", 1.);
 279 // --- and in aluminium parts ---
 280     gMC->Gstpar(idtmed[701], "LOSS", 3.);
 281     gMC->Gstpar(idtmed[701], "DRAY", 1.);
 282 }
 283
 284 //______________________________________________________________________________
 285
 286 void AliPHOS::AddPHOSCradles()
 287 {
 288   Int_t i;
 289   for(i=0;i<GetCradlesAmount();i++) {
 290
 291     int n = fCradles->GetEntries();
 292     fCradles->Add(new AliPHOSCradle( IsVersion(),            // geometry.
 293                                      GetCrystalSideSize    (),
 294                                      GetCrystalLength      (),
 295                                      GetWrapThickness      (),
 296                                      GetAirThickness       (),
 297                                      GetPIN_SideSize       (),
 298                                      GetPIN_Length         (),
 299                                      GetRadius             (),
 300                                      GetNz                 (),
 301                                      GetNphi               (),
 302                                      GetCradleAngle        (i)));
 303
 304     if( n+1 != fCradles->GetEntries() || NULL == fCradles->At(n) )
 305       {
 306         cout << "  Can not create or add AliPHOSCradle.\n";
 307         exit(1);
 308       }
 309   }
 310 }
 311
 312 //______________________________________________________________________________
 313
 314 Int_t AliPHOS::DistancetoPrimitive(Int_t , Int_t )
 315 {
 316    return 9999;
 317 }
 318
 319 //___________________________________________
 320 void AliPHOS::Init()
 321 {
 322   Int_t i;
 323   //
 324   printf("\n");
 325   for(i=0;i<35;i++) printf("*");
 326   printf(" PHOS_INIT ");
 327   for(i=0;i<35;i++) printf("*");
 328   printf("\n");
 329   //
 330   // Here the ABSO initialisation code (if any!)
 331   for(i=0;i<80;i++) printf("*");
 332   printf("\n");
 333 }
 334
 335 //______________________________________________________________________________
 336
 337 void AliPHOS::MakeBranch(Option_t *)
 338 {
 339 // ROOT output initialization to ROOT file.
 340 //
 341 // AliDetector::MakeBranch()  is always called.
 342 //
 343 // There will be also special tree "PHOS" with one branch "AliPHOSCradles"
 344 // if it was set next flag in the galice card file:
 345 //  * PHOSflags:    YES: X<>0   NO: X=0
 346 //  * PHOSflags(1) : -----X.  Create branch for TObjArray of AliPHOSCradle
 347 //     Examples:
 348 //     PHOSflags      1.
 349 //     PHOSflags 636301.
 350 // In that case special bit CradlesBranch_Bit will be set for AliPHOS
 351
 352   AliDetector::MakeBranch();
 353
 354   int i;
 355   float t = GetPHOS_flag(0)/10;
 356   i = (int) t;
 357   i = (int) ((t-i)*10);
 358   if( !i )
 359     return;
 360
 361   SetBit(CradlesBranch_Bit);
 362
 363   if( NULL==(fTreePHOS=new TTree(fTreeName.Data(),"PHOS events tree")) )
 364   {
 365     Error("MakeBranch","Can not create TTree");
 366     exit(1);
 367   }
 368
 369   if( NULL==fTreePHOS->GetCurrentFile() )
 370   {
 371     Error("MakeBranch","There is no opened ROOT file");
 372     exit(1);
 373   }
 374
 375   // Create a new branch in the current Root Tree.
 376
 377   if( NULL==fTreePHOS->Branch(fBranchNameOfCradles.Data(),"TObjArray",&fCradles,4000,0) )
 378   {
 379     Error("MakeBranch","Can not create branch");
 380     exit(1);
 381   }
 382
 383   printf("The branch %s has been created\n",fBranchNameOfCradles.Data());
 384 }
 385
 386 //______________________________________________________________________________
 387
 388 void AliPHOS::SetTreeAddress(void)
 389 {
 390 // ROOT input initialization.
 391 //
 392 // AliDetector::SetTreeAddress()  is always called.
 393 //
 394 // If CradlesBranch_Bit is set (see AliPHOS::MakeBranch) than fTreePHOS is
 395 // initilized.
 396
 397   AliDetector::SetTreeAddress();
 398
 399   if( !TestBit(CradlesBranch_Bit) )
 400     return;
 401
 402   if( NULL==(fTreePHOS=(TTree*)gDirectory->Get((char*)(fTreeName.Data()))  ) )
 403   {
 404     Error("SetTreeAddress","Can not find Tree \"%s\"\n",fTreeName.Data());
 405     exit(1);
 406   }
 407
 408   TBranch *branch = fTreePHOS->GetBranch(fBranchNameOfCradles.Data());
 409   if( NULL==branch )
 410   {
 411     Error("SetTreeAddress","Can not find branch %s in TTree:%s",fBranchNameOfCradles.Data(),fTreeName.Data());
 412     exit(1);
 413   }
 414
 415   branch->SetAddress(&fCradles);
 416 }
 417
 418 //______________________________________________________________________________
 419
 420 AliPHOSCradle *AliPHOS::GetCradleOfTheParticle(const TVector3 &p,const TVector3 &v) const
 421 {
 422 // For a given direction 'p' and source point 'v' returns pointer to AliPHOSCradle
 423 // in that direction or NULL if AliPHOSCradle was not found.
 424
 425   for( int m=0; m<fCradles->GetEntries(); m++ )
 426   {
 427     AliPHOS *PHOS = (AliPHOS *)this;     // Removing 'const'...
 428     AliPHOSCradle *cradle = (AliPHOSCradle *)PHOS->fCradles->operator[](m);
 429
 430     float x,y,l;
 431     const float d = cradle->GetRadius();
 432     cradle->GetXY(p,v,d,x,y,l);
 433
 434     if( l>0 && TMath::Abs(x)<cradle->GetNz  ()*cradle->GetCellSideSize()/2
 435             && TMath::Abs(y)<cradle->GetNphi()*cradle->GetCellSideSize()/2 )
 436       return cradle;
 437   }
 438
 439   return NULL;
 440 }
 441
 442 //______________________________________________________________________________
 443
 444 void AliPHOS::Reconstruction(Float_t signal_step, UInt_t min_signal_reject)
 445 {
 446 // Call AliPHOSCradle::Reconstruction(Float_t signal_step, UInt_t min_signal_reject)
 447 // for all AliPHOSCradles.
 448
 449   for( int i=0; i<fCradles->GetEntries(); i++ )
 450     GetCradle(i).Reconstruction(signal_step,min_signal_reject);
 451 }
 452
 453 //______________________________________________________________________________
 454
 455 void AliPHOS::ResetDigits(void)
 456 {
 457   AliDetector::ResetDigits();
 458
 459   for( int i=0; i<fCradles->GetEntries(); i++ )
 460     ((AliPHOSCradle*)(*fCradles)[i]) -> Clear();
 461 }
 462
 463 //______________________________________________________________________________
 464
 465 void AliPHOS::FinishEvent(void)
 466 {
 467 // Called at the end of each 'galice' event.
 468
 469   if( NULL!=fTreePHOS )
 470     fTreePHOS->Fill();
 471 }
 472
 473 //______________________________________________________________________________
 474
 475 void AliPHOS::FinishRun(void)
 476 {
 477 }
 478
 479 //______________________________________________________________________________
 480
 481 void AliPHOS::Print(Option_t *opt)
 482 {
 483 // Print PHOS information.
 484 // For each AliPHOSCradle the function AliPHOSCradle::Print(opt) is called.
 485
 486   AliPHOS &PHOS = *(AliPHOS *)this;     // Removing 'const'...
 487
 488   for( int i=0; i<fCradles->GetEntries(); i++ )
 489   {
 490     printf("PHOS cradle %d from %d\n",i+1, fCradles->GetEntries());
 491     PHOS.GetCradle(i).Print(opt);
 492     printf( "---------------------------------------------------\n");
 493   }
 494 }
 495
 496 //______________________________________________________________________________
 497 void AliPHOS::SetFlags(Float_t p1,Float_t p2,Float_t p3,Float_t p4,
 498                        Float_t p5,Float_t p6,Float_t p7,Float_t p8,Float_t p9)
 499 {
 500   PHOSflags[0]=p1;
 501   PHOSflags[1]=p2;
 502   PHOSflags[2]=p3;
 503   PHOSflags[3]=p4;
 504   PHOSflags[4]=p5;
 505   PHOSflags[5]=p6;
 506   PHOSflags[6]=p7;
 507   PHOSflags[7]=p8;
 508   PHOSflags[8]=p9;
 509 }
 510
 511 //______________________________________________________________________________
 512 void AliPHOS::SetCell(Float_t p1,Float_t p2,Float_t p3,Float_t p4,
 513                        Float_t p5,Float_t p6,Float_t p7,Float_t p8,Float_t p9)
 514 {
 515   PHOScell[0]=p1;
 516   PHOScell[1]=p2;
 517   PHOScell[2]=p3;
 518   PHOScell[3]=p4;
 519   PHOScell[4]=p5;
 520   PHOScell[5]=p6;
 521   PHOScell[6]=p7;
 522   PHOScell[7]=p8;
 523   PHOScell[8]=p9;
 524 }
 525
 526 //______________________________________________________________________________
 527 void AliPHOS::SetRadius(Float_t radius)
 528 {
 529    PHOSradius=radius;
 530 }
 531
 532 //______________________________________________________________________________
 533 void AliPHOS::SetCradleSize(Int_t nz, Int_t nphi, Int_t ncradles)
 534 {
 535    PHOSsize[0]=nz;
 536    PHOSsize[1]=nphi;
 537    PHOSsize[2]=ncradles;
 538 }
 539
 540 //______________________________________________________________________________
 541 void AliPHOS::SetCradleA(Float_t angle)
 542 {
 543    PHOScradlesA=angle;
 544 }
 545
 546 //______________________________________________________________________________
 547 void AliPHOS::SetExtra(Float_t p1,Float_t p2,Float_t p3,Float_t p4,
 548                        Float_t p5,Float_t p6,Float_t p7,Float_t p8,Float_t p9)
 549 {
 550    PHOSextra[0] = p1;
 551    PHOSextra[1] = p2;
 552    PHOSextra[2] = p3;
 553    PHOSextra[3] = p4;
 554    PHOSextra[4] = p5;
 555    PHOSextra[5] = p6;
 556    PHOSextra[6] = p7;
 557    PHOSextra[7] = p8;
 558    PHOSextra[8] = p9;
 559 }
 560
 561 //______________________________________________________________________________
 562 void AliPHOS::SetTextolitWall(Float_t dx, Float_t dy, Float_t dz)
 563 {
 564    PHOSTXW[0] = dx;
 565    PHOSTXW[1] = dy;
 566    PHOSTXW[2] = dz;
 567 }
 568
 569 //______________________________________________________________________________
 570 void AliPHOS::SetInnerAir(Float_t dx, Float_t dy, Float_t dz)
 571 {
 572    PHOSAIR[0] = dx;
 573    PHOSAIR[1] = dy;
 574    PHOSAIR[2] = dz;
 575 }
 576
 577 //______________________________________________________________________________
 578 void AliPHOS::SetFoam(Float_t dx, Float_t dy, Float_t dz, Float_t dr)
 579 {
 580    PHOSFTI[0] = dx;
 581    PHOSFTI[1] = dy;
 582    PHOSFTI[2] = dz;
 583    PHOSFTI[3] = dr;
 584 }
 585
 586 ClassImp(AliPHOSCradle)
 587
 588 //______________________________________________________________________________
 589
 590 AliPHOSCradle::AliPHOSCradle(void) {}
 591
 592 //______________________________________________________________________________
 593
 594 AliPHOSCradle::AliPHOSCradle( int   Geometry           ,
 595                               float CrystalSideSize    ,
 596                               float CrystalLength      ,
 597                               float WrapThickness      ,
 598                               float AirThickness       ,
 599                               float PIN_SideSize       ,
 600                               float PIN_Length         ,
 601                               float Radius             ,
 602                               int   Nz                 ,
 603                               int   Nphi               ,
 604                               float Angle              ) :
 605     fGeometry                   (Geometry),
 606 //  fCellEnergy                 (),
 607 //  fChargedTracksInPIN         (),
 608     fCrystalSideSize            (CrystalSideSize),
 609     fCrystalLength              (CrystalLength),
 610     fWrapThickness              (WrapThickness),
 611     fAirThickness               (AirThickness),
 612     fPIN_SideSize               (PIN_SideSize),
 613     fPIN_Length                 (PIN_Length),
 614     fRadius                     (Radius),
 615     fNz                         (Nz),
 616     fNphi                       (Nphi),
 617     fPhi                        (Angle)
 618 {
 619         fCellEnergy         = TH2F("CellE","Energy deposition in a cells",fNz,0,fNz,fNphi,0,fNphi);
 620         fCellEnergy           .SetDirectory(0);
 621         fChargedTracksInPIN = TH2S("PINCtracks","Amount of charged tracks in PIN",fNz,0,fNz,fNphi,0,fNphi);
 622         fChargedTracksInPIN   .SetDirectory(0);
 623 }
 624
 625 //______________________________________________________________________________
 626
 627 AliPHOSCradle::~AliPHOSCradle(void)        // 28.12.1998
 628 {
 629   fGammasReconstructed.Delete();
 630   fParticles          .Delete();
 631 }
 632
 633 //______________________________________________________________________________
 634
 635 void AliPHOSCradle::Clear(Option_t *)
 636 {
 637 // Clear digit. information.
 638
 639   fCellEnergy              .Reset();
 640   fChargedTracksInPIN      .Reset();
 641   GetParticles()           .Delete();
 642   GetParticles()           .Compress();
 643   GetGammasReconstructed() .Delete();
 644   GetGammasReconstructed() .Compress();
 645
 646 }
 647
 648 //______________________________________________________________________________
 649
 650 void AliPHOSCradle::GetXY(const TVector3 &p,const TVector3 &v,float R,float &x,float &y,float &l) const
 651 {
 652 // This function calculates hit position (x,y) in the CRADLE cells plain from particle in
 653 // the direction given by 'p' (not required to be normalized) and start point
 654 // given by 3-vector 'v'. So the particle trajectory is   t(l) = v + p*l
 655 // were 'l' is a number (distance from 'v' to CRADLE cells plain) and 't' is resulting
 656 // three-vector of trajectory point.
 657 //
 658 // After the call to this function user should test that l>=0 (the particle HITED the
 659 // plain) and (x,y) are in the region of CRADLE:
 660 //
 661 // Example:
 662 //   AliPHOSCradle cradle(......);
 663 //   TVector3 p(....), v(....);
 664 //   Float_t x,y,l;
 665 //   cradle.GetXY(p,v,x,y,l);
 666 //   if( l<0 || TMath::Abs(x)>cradle.GetNz()  *cradle.GetCellSideSize()/2
 667 //           || TMath::Abs(y)>cradle.GetNphi()*cradle.GetCellSideSize()/2 )
 668 //     cout << "Outside the CRADLE.\n";
 669
 670   // We have to create three vectors:
 671   //    s  - central point on the PHOS surface
 672   //    n1 - first vector in CRADLE plain
 673   //    n2 - second vector in CRADLE plain
 674   // This three vectors are orthonormalized.
 675
 676   double phi = fPhi/180*TMath::Pi();
 677   TVector3        n1(   0.0      ,   0.0      , 1.0 ),   // Z direction (X)
 678                   n2(  -sin(phi) ,   cos(phi) , 0 ),   // around beam (Y)
 679                   s ( R*cos(phi) , R*sin(phi) , 0 );   // central point
 680
 681   const double l1_min = 1e-2;
 682   double l1,
 683          p_n1 = p*n1,        // * - scalar product.
 684          p_n2 = p*n2,
 685          v_n1 = v*n1,
 686          v_n2 = v*n2,
 687          s_n1 = s*n1, // 0
 688          s_n2 = s*n2; // 0
 689
 690   if      ( TMath::Abs(l1=p.X()-n1.X()*p_n1-n2.X()*p_n2)>l1_min )
 691     { l = (-v.X()+s.X()+n1.X()*(v_n1-s_n1)+n2.X()*(v_n2-s_n2))/l1; }
 692   else if ( TMath::Abs(l1=p.Y()-n1.Y()*p_n1-n2.Y()*p_n2)>l1_min )
 693     { l = (-v.Y()+s.Y()+n1.Y()*(v_n1-s_n1)+n2.Y()*(v_n2-s_n2))/l1; }
 694   else if ( TMath::Abs(l1=p.Z()-n1.Z()*p_n1-n2.Z()*p_n2)>l1_min )
 695     { l = (-v.Z()+s.Z()+n1.Z()*(v_n1-s_n1)+n2.Z()*(v_n2-s_n2))/l1; }
 696
 697 //         double lx = (-v.X()+s.X()+n1.X()*(v.dot(n1)-s.dot(n1))+n2.X()*(v.dot(n2)-s.dot(n2)))/
 698 //                     (p.X()-n1.X()*p.dot(n1)-n2.X()*p.dot(n2)),
 699 //                ly = (-v.Y()+s.Y()+n1.Y()*(v.dot(n1)-s.dot(n1))+n2.Y()*(v.dot(n2)-s.dot(n2)))/
 700 //                     (p.Y()-n1.Y()*p.dot(n1)-n2.Y()*p.dot(n2)),
 701 //                lz = (-v.Z()+s.Z()+n1.Z()*(v.dot(n1)-s.dot(n1))+n2.Z()*(v.dot(n2)-s.dot(n2)))/
 702 //                     (p.Z()-n1.Z()*p.dot(n1)-n2.Z()*p.dot(n2));
 703 //         cout.form("x: %g %g %g %g\n",lx,-v.X()+s.X()+n1.X()*(v.dot(n1)-s.dot(n1))+n2.X()*(v.dot(n2)-s.dot(n2)),p.X()-n1.X()*p.dot(n1)-n2.X()*p.dot(n2));
 704 //         cout.form("y: %g %g %g %g\n",lx,-v.Y()+s.Y()+n1.Y()*(v.dot(n1)-s.dot(n1))+n2.Y()*(v.dot(n2)-s.dot(n2)),p.Y()-n1.Y()*p.dot(n1)-n2.Y()*p.dot(n2));
 705 //         cout.form("z: %g %g %g %g\n",lx,-v.Z()+s.Z()+n1.Z()*(v.dot(n1)-s.dot(n1))+n2.Z()*(v.dot(n2)-s.dot(n2)),p.Z()-n1.Z()*p.dot(n1)-n2.Z()*p.dot(n2));
 706 //         cout.form("lx,ly,lz =   %g,%g,%g\n",lx,ly,lz);
 707
 708   x = p_n1*l + v_n1 - s_n1;
 709   y = p_n2*l + v_n2 - s_n2;
 710 }
 711
 712 //______________________________________________________________________________
 713
 714 void AliPHOSCradle::Print(Option_t *opt)
 715 {
 716 // Print AliPHOSCradle information.
 717 //
 718 // options:  'd' - print energy deposition for EVERY cell
 719 //           'p' - print particles list that hit the cradle
 720 //           'r' - print list of reconstructed particles
 721
 722   AliPHOSCradle *cr = (AliPHOSCradle *)this;     // Removing 'const'...
 723
 724   printf("AliPHOSCradle:  Nz=%d  Nphi=%d, fPhi=%f, E=%g\n",fNz,fNphi,fPhi,
 725        cr->fCellEnergy.GetSumOfWeights());
 726
 727   if( NULL!=strchr(opt,'d') )
 728   {
 729     printf("\n\nCells Energy (in MeV):\n\n   |");
 730     for( int x=0; x<fNz; x++ )
 731       printf(" %4d|",x+1);
 732     printf("\n");
 733
 734     for( int y=fNphi-1; y>=0; y-- )
 735     {
 736       printf("%3d|",y+1);
 737       for( int x=0; x<fNz; x++ )
 738         printf("%6d",(int)(cr->fCellEnergy.GetBinContent(cr->fCellEnergy.GetBin(x,y))*1000));
 739       printf("\n");
 740     }
 741     printf("\n");
 742   }
 743
 744   if( NULL!=strchr(opt,'p') )
 745   {
 746     printf("This cradle was hit by %d particles\n",
 747          ((AliPHOSCradle*)this)->GetParticles().GetEntries());
 748     TObjArray &p=((AliPHOSCradle*)this)->GetParticles();
 749     for( int i=0; i<p.GetEntries(); i++ )
 750       ((AliPHOSgamma*)(p[i]))->Print();
 751   }
 752
 753   if( NULL!=strchr(opt,'p') )
 754   {
 755     printf("Amount of reconstructed gammas is %d\n",
 756          ((AliPHOSCradle*)this)->GetGammasReconstructed().GetEntries());
 757
 758     TObjArray &p=((AliPHOSCradle*)this)->GetGammasReconstructed();
 759     for( int i=0; i<p.GetEntries(); i++ )
 760       ((AliPHOSgamma*)(p[i]))->Print();
 761   }
 762 }
 763
 764 //______________________________________________________________________________
 765
 766 void AliPHOSCradle::Distortion(const TH2F *Noise, const TH2F *Stochastic, const TH2F *Calibration)
 767 {
 768 // This function changes histogram of cell energies fCellEnergy on the base of input
 769 // histograms Noise, Stochastic, Calibration. The histograms must have
 770 // size Nz x Nphi.
 771
 772   //////////////////////////////////
 773   // Testing the histograms size. //
 774   //////////////////////////////////
 775
 776   if( fNz!=fCellEnergy.GetNbinsX() || fNphi!=fCellEnergy.GetNbinsY() )
 777   {
 778     printf      ("Bad size of CellEnergy!   Must be:   Nz x Nphi = %d x %d\n"
 779                  "but size of CellEnergy is:  %d x %d\n",
 780                  fNz,fNphi,fCellEnergy.GetNbinsX(),fCellEnergy.GetNbinsY());
 781     exit(1);
 782   }
 783
 784   if( fNz!=fChargedTracksInPIN.GetNbinsX() || fNphi!=fChargedTracksInPIN.GetNbinsY() )
 785   {
 786     printf      ("Bad size of ChargedTracksInPIN!   Must be:   Nz x Nphi = %d x %d\n"
 787                  "but size of ChargedTracksInPIN is:  %d x %d\n",
 788                  fNz,fNphi,fChargedTracksInPIN.GetNbinsX(),fChargedTracksInPIN.GetNbinsY());
 789     exit(1);
 790   }
 791
 792   if( NULL!=Noise && (fNz!=Noise->GetNbinsX() || fNphi!=Noise->GetNbinsX()) )
 793   {
 794     printf      ("Bad size of Noise!   Must be:   Nz x Nphi = %d x %d\n"
 795                  "but size of Noise is:  %d x %d\n",
 796                  fNz,fNphi,fChargedTracksInPIN.GetNbinsX(),fChargedTracksInPIN.GetNbinsY());
 797     exit(1);
 798   }
 799
 800   if( NULL!=Stochastic && (fNz!=Stochastic->GetNbinsX() || fNphi!=Stochastic->GetNbinsX()) )
 801   {
 802     printf      ("Bad size of Stochastic!   Must be:   Nz x Nphi = %d x %d\n"
 803                  "but size of Stochastic is:  %d x %d\n",
 804                  fNz,fNphi,fChargedTracksInPIN.GetNbinsX(),fChargedTracksInPIN.GetNbinsY());
 805     exit(1);
 806   }
 807
 808   if( NULL!=Calibration && (fNz!=Calibration->GetNbinsX() || fNphi!=Calibration->GetNbinsX()) )
 809   {
 810     printf      ("Bad size of Calibration!   Must be:   Nz x Nphi = %d x %d\n"
 811                  "but size of Calibration is:  %d x %d\n",
 812                  fNz,fNphi,fChargedTracksInPIN.GetNbinsX(),fChargedTracksInPIN.GetNbinsY());
 813     exit(1);
 814   }
 815
 816   ////////////////////
 817   // Do distortion! //
 818   ////////////////////
 819
 820   for( int y=0; y<fNphi; y++ )
 821     for( int x=0; x<fNz; x++ )
 822     {
 823       const int n = fCellEnergy.GetBin(x,y);   // Bin number
 824       static TRandom r;
 825
 826       Float_t   E_old=fCellEnergy.GetBinContent(n),   E_new=E_old;
 827
 828       if( NULL!=Stochastic )
 829         E_new   = r.Gaus(E_old,sqrt(E_old)*GetDistortedValue(Stochastic,n));
 830
 831       if( NULL!=Calibration )
 832         E_new  *=  GetDistortedValue(Calibration,n);
 833
 834       if( NULL!=Noise )
 835         E_new  +=  GetDistortedValue(Noise,n);
 836
 837       fCellEnergy.SetBinContent(n,E_new);
 838     }
 839 }
 840
 841 ////////////////////////////////////////////////////////////////////////////////
 842
 843 TH2F* AliPHOSCradle::CreateHistForDistortion(const char *name, const char *title,
 844                                              Int_t Nx, Int_t Ny,
 845                                              Float_t MU_mu,    Float_t MU_sigma,
 846                                              Float_t SIGMA_mu, Float_t SIGMA_sigma)
 847 {
 848 // Create (new TH2F(...)) histogram with information (for every bin) that will
 849 // be used for VALUE creation.
 850 // Two values will be created for each bin:
 851 // MU    = TRandom::Gaus(MU_mu,MU_sigma)
 852 // and
 853 // SIGMA = TRandom::Gaus(SIGMA_mu,SIGMA_sigma)
 854 // The VALUE in a particluar bin will be equal
 855 // VALUE = TRandom::Gaus(MU,SIGMA)
 856 //
 857 // Do not forget to delete the histogram at the end of the work.
 858
 859   TH2F *h = new TH2F( name,title, Nx,1,Nx, Ny,1,Ny );
 860   if( h==NULL )
 861   {
 862     Error("CreateHistForDistortion","Can not create the histogram");
 863     exit(1);
 864   }
 865   h->SetDirectory(0);
 866
 867   for( int y=0; y<Ny; y++ )
 868     for( int x=0; x<Nx; x++ )
 869     {
 870       const int n = h->GetBin(x,y);
 871       h->SetBinContent(n,r.Gaus(   MU_mu,   MU_sigma));
 872       h->SetBinError  (n,r.Gaus(SIGMA_mu,SIGMA_sigma));
 873     }
 874
 875   return h;
 876 }
 877
 878 ////////////////////////////////////////////////////////////////////////////////
 879
 880 Float_t AliPHOSCradle::GetDistortedValue(const TH2F *h, UInt_t n)
 881 {
 882   return r.Gaus(((TH2F*)h)->GetBinContent(n),n);
 883 }
 884
 885 ////////////////////////////////////////////////////////////////////////////////
 886 //______________________________________________________________________________
 887
 888 #ifdef WIN32
 889   #define common_for_event_storing COMMON_FOR_EVENT_STORING
 890 #else
 891   #define common_for_event_storing common_for_event_storing_
 892 #endif
 893
 894 /* extern "C" */ struct
 895 {
 896   enum { crystals_matrix_amount_max=4, crystals_in_matrix_amount_max=40000 };
 897
 898   // Event-independent information
 899   UShort_t      crystals_matrix_amount_PHOS,
 900                 crystal_matrix_type,
 901                 amount_of_crystals_on_Z,
 902                 amount_of_crystals_on_PHI;
 903   Float_t       radius,
 904                 crystal_size,
 905                 crystal_length,
 906                 matrix_coordinate_Z             [crystals_matrix_amount_max],
 907                 matrix_coordinate_PHI           [crystals_matrix_amount_max];
 908   UInt_t        event_number;
 909   UShort_t      crystals_amount_with_amplitudes [crystals_matrix_amount_max],
 910                 crystals_amplitudes_Iad         [crystals_matrix_amount_max]
 911                                                 [crystals_in_matrix_amount_max][2];
 912 } common_for_event_storing;
 913
 914 //       integer*4 crystals_amount_max,crystals_in_matrix_amount_max,
 915 //      +          crystals_matrix_amount_max
 916 //       parameter (crystals_matrix_amount_max=4)
 917 //       parameter (crystals_in_matrix_amount_max=40000)
 918 //       parameter (crystals_amount_max =crystals_matrix_amount_max*
 919 //      +                                crystals_in_matrix_amount_max)
 920 //
 921 // * All units are in GeV, cm, radian
 922 //       real       crystal_amplitudes_unit, radius_unit,
 923 //      +           crystal_size_unit, crystal_length_unit,
 924 //      +           matrix_coordinate_Z_unit, matrix_coordinate_PHI_unit
 925 //       integer    crystal_amplitudes_in_units_min
 926 //       parameter (crystal_amplitudes_in_units_min        = 1)
 927 //       parameter (crystal_amplitudes_unit                = 0.001 ) ! 1.0  MeV
 928 //       parameter (radius_unit                            = 0.1   ) ! 0.1  cm
 929 //       parameter (crystal_size_unit                      = 0.01  ) ! 0.01 cm
 930 //       parameter (crystal_length_unit                    = 0.01  ) ! 0.01 cm
 931 //       parameter (matrix_coordinate_Z_unit               = 0.1   ) ! 0.1  cm
 932 //       parameter (matrix_coordinate_PHI_unit             = 1e-4  ) ! 1e-4 radian
 933 //
 934 //       integer*2 crystals_matrix_amount_PHOS, crystal_matrix_type,
 935 //      +          amount_of_crystals_on_Z, amount_of_crystals_on_PHI,
 936 //      +          crystals_amount_with_amplitudes, crystals_amplitudes_Iad
 937 //       integer*4 event_number
 938 //
 939 //       real      radius, crystal_size, crystal_length,
 940 //      +          matrix_coordinate_Z, matrix_coordinate_PHI
 941 //
 942 //       real      crystals_amplitudes, crystals_energy_total
 943 //       integer   event_file_unit_number
 944 //
 945 //       common /common_for_event_storing/
 946 //      + ! Event-independent information
 947 //      +        crystals_matrix_amount_PHOS,
 948 //      +        crystal_matrix_type,
 949 //      +        amount_of_crystals_on_Z,
 950 //      +        amount_of_crystals_on_PHI,
 951 //      +        radius,
 952 //      +        crystal_size,
 953 //      +        crystal_length,
 954 //      +        matrix_coordinate_Z     (crystals_matrix_amount_max),
 955 //      +        matrix_coordinate_PHI   (crystals_matrix_amount_max),
 956 //      +
 957 //      + ! Event-dependent information
 958 //      +        event_number,
 959 //      +        crystals_amount_with_amplitudes
 960 //      +                                (crystals_matrix_amount_max),
 961 //      +        crystals_amplitudes_Iad (2,crystals_in_matrix_amount_max,
 962 //      +                                 crystals_matrix_amount_max),
 963 //      +
 964 //      + ! These information don't store in data file
 965 //      +        crystals_amplitudes     (crystals_amount_max),
 966 //      +        crystals_energy_total,
 967 //      +        event_file_unit_number
 968
 969
 970 //      parameter (NGp=1000,nsps=10,nvertmax=1000)
 971 //         COMMON /GAMMA/KG,MW(ngp),ID(ngp),JD(ngp),E(ngp),E4(ngp),
 972 //      ,  XW(ngp),YW(ngp),ES(nsps,ngp),ET(nsps,ngp),ISsd(ngp),
 973 //      ,  IGDEV(ngp),ZGDEV(ngp),sigexy(3,ngp),Emimx(2,nsps,ngp),
 974 //      ,  kgfix,igfix(ngp),cgfix(3,ngp),sgfix(3,ngp),hiw(ngp),
 975 //      ,  wsw(nsps,ngp),h1w(ngp),h0w(ngp),raxay(5,ngp),
 976 //      ,  sigmaes0(nsps,ngp),dispeces(nsps,ngp),
 977 //      ,  igamvert(ngp)
 978
 979
 980 #ifdef WIN32
 981 #define rcgamma RCGAMMA
 982 #else
 983 #define rcgamma rcgamma_
 984 #endif
 985
 986 /* extern "C" */ struct
 987 {
 988   enum {NGP=1000, nsps=10, nvertmax=1000};
 989   int   recons_gammas_amount, mw[NGP],ID[NGP],JD[NGP];
 990   float E[NGP], E4[NGP], XW[NGP], YW[NGP], ES[NGP][nsps],ET[NGP][nsps],ISsd[NGP],
 991         igdev[NGP],Zgdev[NGP];
 992 //      sigexy(3,ngp),Emimx(2,nsps,ngp),
 993 //   ,  kgfix,igfix(ngp),cgfix(3,ngp),sgfix(3,ngp),hiw(ngp),
 994 //   ,  wsw(nsps,ngp),h1w(ngp),h0w(ngp),raxay(5,ngp),
 995 //   ,  sigmaes0(nsps,ngp),dispeces(nsps,ngp),
 996 //   ,  igamvert(ngp)
 997 } rcgamma;
 998
 999 #ifdef WIN32
1000 #define reconsfirst RECONSFIRST
1001 #define type_of_call _stdcall
1002 #else
1003 #define reconsfirst reconsfirst_
1004 #define type_of_call
1005 #endif
1006
1007 extern "C" void type_of_call reconsfirst(const float &,const float &);
1008
1009 void AliPHOSCradle::Reconstruction(Float_t signal_step, UInt_t min_signal_reject)
1010 {
1011 // Call of PHOS reconstruction program.
1012 // signal_step=0.001  GeV (1MeV)
1013 // min_signal_reject = 15 or 30 MeV
1014
1015
1016   common_for_event_storing.event_number                       = 0;  // We do not know event number?
1017   common_for_event_storing.crystals_matrix_amount_PHOS        = 1;
1018   common_for_event_storing.crystal_matrix_type                = 1; // 1 - rectangular
1019   common_for_event_storing.amount_of_crystals_on_Z            = fNz;
1020   common_for_event_storing.amount_of_crystals_on_PHI          = fNphi;
1021
1022   common_for_event_storing.radius                             = fRadius;
1023   common_for_event_storing.crystal_size                       = GetCellSideSize();
1024   common_for_event_storing.crystal_length                     = fCrystalLength;
1025
1026   common_for_event_storing.matrix_coordinate_Z            [0] = 0;
1027   common_for_event_storing.matrix_coordinate_PHI          [0] = fPhi;
1028
1029   #define  k    common_for_event_storing.crystals_amount_with_amplitudes[0]
1030   k=0;
1031
1032   for( int y=0; y<fNphi; y++ )
1033     for( int x=0; x<fNz; x++ )
1034     {
1035       UInt_t    n       = fCellEnergy.GetBin(x,y);
1036       UInt_t    signal  = (int) (fCellEnergy.GetBinContent(n)/signal_step);
1037       if( signal>=min_signal_reject )
1038       {
1039         common_for_event_storing.crystals_amplitudes_Iad[0][k][0] = signal;
1040         common_for_event_storing.crystals_amplitudes_Iad[0][k][1] = x + y*fNz;
1041         k++;
1042       }
1043     }
1044   #undef  k
1045
1046   GetGammasReconstructed().Delete();
1047   GetGammasReconstructed().Compress();
1048
1049   const float   stochastic_term   = 0.03,        // per cents over sqrt(E);  E in GeV
1050                 electronic_noise  = 0.01;        // GeV
1051   reconsfirst(stochastic_term,electronic_noise); // Call of reconstruction program.
1052
1053   for( int i=0; i<rcgamma.recons_gammas_amount; i++ )
1054   {
1055 //     new (GetGammasReconstructed().UncheckedAt(i) ) AliPHOSgamma;
1056 //     AliPHOSgamma &g = *(AliPHOSgamma*)(GetGammasReconstructed().UncheckedAt(i));
1057
1058     AliPHOSgamma *gggg = new AliPHOSgamma;
1059     if( NULL==gggg )
1060     {
1061       Error("Reconstruction","Can not create AliPHOSgamma");
1062       exit(1);
1063     }
1064
1065     GetGammasReconstructed().Add(gggg);
1066     AliPHOSgamma &g=*gggg;
1067
1068     Float_t thetta, alpha, betta, R=fRadius+rcgamma.Zgdev[i]/10;
1069
1070     g.fX      = rcgamma.YW[i]/10;
1071     g.fY      = rcgamma.XW[i]/10;
1072     g.fE      = rcgamma.E [i];
1073
1074     thetta      = atan(g.fX/R);
1075
1076     alpha = atan(g.fY/R);
1077     betta = fPhi/180*TMath::Pi() + alpha;
1078
1079     g.fPx = g.fE * cos(thetta) * cos(betta);
1080     g.fPy = g.fE * cos(thetta) * sin(betta);
1081     g.fPz = g.fE * sin(thetta);
1082   }
1083 }
1084
1085 //______________________________________________________________________________
1086 //______________________________________________________________________________
1087 //______________________________________________________________________________
1088 //______________________________________________________________________________
1089 //______________________________________________________________________________
1090
1091 ClassImp(AliPHOSgamma)
1092
1093 //______________________________________________________________________________
1094
1095 void AliPHOSgamma::Print(Option_t *)
1096 {
1097   float mass = fE*fE - fPx*fPx - fPy*fPy - fPz*fPz;
1098
1099   if( mass>=0 )
1100     mass =  sqrt( mass);
1101   else
1102     mass = -sqrt(-mass);
1103
1104   printf("XY=(%+7.2f,%+7.2f)  (%+7.2f,%+7.2f,%+7.2f;%7.2f)  mass=%8.4f  Ipart=%2d\n",
1105           fX,fY,fPx,fPy,fPz,fE,mass,fIpart);
1106 }
1107
1108 //______________________________________________________________________________
1109
1110 AliPHOSgamma &AliPHOSgamma::operator=(const AliPHOSgamma &g)
1111 {
1112   fX           = g.fX;
1113   fY           = g.fY;
1114   fE           = g.fE;
1115   fPx          = g.fPx;
1116   fPy          = g.fPy;
1117   fPz          = g.fPz;
1118   fIpart       = g.fIpart;
1119
1120   return *this;
1121 }
1122
1123 //______________________________________________________________________________
1124 //______________________________________________________________________________
1125 //______________________________________________________________________________
1126 //______________________________________________________________________________
1127 //______________________________________________________________________________
1128
1129 ClassImp(AliPHOShit)
1130
1131 //______________________________________________________________________________
1132
1133 AliPHOShit::AliPHOShit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits):
1134 AliHit(shunt, track)
1135 {
1136    Int_t i;
1137    for (i=0;i<5;i++) fVolume[i] = vol[i];
1138    fX       = hits[0];
1139    fY       = hits[1];
1140    fZ       = hits[2];
1141    fELOS    = hits[3];
1142 }
1143
1144 //______________________________________________________________________________