PHOS/AliPHOSv1.cxx

   1 /////////////////////////////////////////////////////////
   2 //  Manager and hits classes for set:PHOS version 1    //
   3 /////////////////////////////////////////////////////////
   4
   5 // --- ROOT system ---
   6 #include "TH1.h"
   7 #include "TRandom.h"
   8 #include "TFile.h"
   9 #include "TTree.h"
  10 #include "TBRIK.h"
  11 #include "TNode.h"
  12 #include "TMath.h"
  13
  14 // --- galice header files ---
  15 #include "AliPHOSv1.h"
  16 #include "AliRun.h"
  17 #include "TGeant3.h"
  18
  19 ClassImp(AliPHOSv1)
  20
  21 //______________________________________________________________________________
  22
  23
  24 AliPHOSv1::AliPHOSv1()
  25 {
  26 }
  27
  28 //______________________________________________________________________________
  29
  30 AliPHOSv1::AliPHOSv1(const char *name, const char *title)
  31           : AliPHOS(name, title)
  32 {
  33 }
  34
  35 //___________________________________________
  36 void AliPHOSv1::CreateGeometry()
  37 {
  38
  39   AliMC* pMC = AliMC::GetMC();
  40
  41   AliPHOS *PHOS_tmp = (AliPHOS*)gAlice->GetModule("PHOS");
  42   if( NULL==PHOS_tmp )
  43   {
  44     printf("There isn't PHOS detector!\n");
  45     return;
  46   }
  47 //  AliPHOS &PHOS = *PHOS_tmp;
  48
  49   //////////////////////////////////////////////////////////////////////////////
  50
  51   Int_t                 rotation_matrix_number=0;
  52   Float_t               par[11],
  53                         x,y,z;
  54
  55   const float           cell_length             = GetCrystalLength()+GetAirThickness()+GetWrapThickness()+GetPIN_Length(),
  56                         cell_side_size          = GetCrystalSideSize()+2*GetAirThickness()+2*GetWrapThickness(),
  57 //                        cell_angle              = 180/kPI * 2 * atan(cell_side_size/2 / GetRadius()),        // radians
  58                         cradle_thikness         = cell_length + GetCPV_Thickness() + GetCPV_PHOS_Distance(),
  59                         distance_to_CPV         = GetRadius() - GetCPV_Thickness() - GetCPV_PHOS_Distance();
  60
  61   //////////////////////////////////////////////////////////////////////////////
  62   // CELL volume and subvolumes creation
  63   //////////////////////////////////////////////////////////////////////////////
  64
  65   par[0] = GetCrystalSideSize()/2 + GetWrapThickness();
  66   par[1] = GetCrystalSideSize()/2 + GetWrapThickness();
  67   par[2] = GetCrystalLength()  /2 + GetWrapThickness()/2;
  68   pMC->Gsvolu("WRAP","BOX ",GetPHOS_IDTMED_Tyvek(),par,3);
  69
  70   par[0] = GetCrystalSideSize()/2;
  71   par[1] = GetCrystalSideSize()/2;
  72   par[2] = GetCrystalLength()/2;
  73   pMC->Gsvolu("CRST","BOX ",GetPHOS_IDTMED_PbWO4(),par,3);
  74
  75   // PIN
  76   par[0] = GetPIN_SideSize()/2;
  77   par[1] = GetPIN_SideSize()/2;
  78   par[2] = GetPIN_Length()/2;
  79   pMC->Gsvolu("PIN ","BOX ",GetPHOS_IDTMED_PIN(),par,3);
  80
  81   //////////////////////////////////////////////////////////////////////////////
  82   // CRADLE,CPV creation.
  83   //////////////////////////////////////////////////////////////////////////////
  84
  85   par[0] = cell_side_size/2 * GetNz();
  86   par[1] = cell_side_size/2 * GetNphi();
  87   par[2] = cradle_thikness/2;
  88   pMC->Gsvolu("PHOS","BOX ",GetPHOS_IDTMED_AIR(),par,3);
  89
  90 //par[0] : the same as above
  91 //par[1] : the same as above
  92   par[2] = GetCPV_Thickness()/2;
  93   pMC->Gsvolu("CPV ","BOX ",GetPHOS_IDTMED_CPV(),par,3);
  94
  95   x = 0;
  96   y = 0;
  97   z = (cell_length+GetCPV_PHOS_Distance())/2;
  98   pMC->Gspos("CPV ",1,"PHOS",x,y,z,0,"ONLY");
  99
 100   par[0] = cell_side_size/2 * GetNz();
 101   par[1] = cell_side_size/2 * GetNphi();
 102   par[2] = cell_length/2;
 103   pMC->Gsvolu("CRS0","BOX ",GetPHOS_IDTMED_AIR(),par,3);
 104
 105   x = 0;
 106   y = 0;
 107   z = -(cradle_thikness-cell_length)/2;
 108   pMC->Gspos("CRS0",1,"PHOS",x,y,z,0,"ONLY");
 109
 110   pMC->Gsdvn("CRS1","CRS0",GetNphi(),2);
 111   pMC->Gsdvn("CELL","CRS1",GetNz()  ,1);
 112
 113   //////////////////////////////////////////////////////////////////////////////
 114   // CELL creation
 115   //////////////////////////////////////////////////////////////////////////////
 116
 117   x = 0;
 118   y = 0;
 119   z = -GetWrapThickness()/2;
 120   pMC->Gspos("CRST",1,"WRAP",x,y,z,0,"ONLY");
 121
 122   x = 0;
 123   y = 0;
 124   z = GetPIN_Length()/2;
 125   pMC->Gspos("WRAP",1,"CELL",x,y,z,0,"ONLY");
 126
 127   x = 0;
 128   y = 0;
 129   z = -GetCrystalLength()/2-GetWrapThickness()/2;
 130   pMC->Gspos("PIN ",1,"CELL",x,y,z,0,"ONLY");
 131
 132   //////////////////////////////////////////////////////////////////////////////
 133   // CELL has been created.
 134   //////////////////////////////////////////////////////////////////////////////
 135
 136 //   int n=0;
 137 //   z = -(GetCPV_Thickness()+GetCPV_PHOS_Distance())/2;
 138 //
 139 //   for( int iy=0; iy<GetNphi(); iy++ )
 140 //   {
 141 //     y = (iy-(GetNphi()-1)/2.)*cell_side_size;
 142 //     for( int ix=0; ix<GetNz(); ix++ )
 143 //     {
 144 //       x = (ix-(GetNz()-1)/2.)*cell_side_size;
 145 //       pMC->Gspos("CELL",++n,"PHOS",x,y,z,0,"ONLY");
 146 //     }
 147 //   }
 148
 149   //////////////////////////////////////////////////////////////////////////////
 150   // End of CRADLE creation.
 151   //////////////////////////////////////////////////////////////////////////////
 152
 153
 154   //////////////////////////////////////////////////////////////////////////////
 155   // PHOS creation
 156   //////////////////////////////////////////////////////////////////////////////
 157
 158   for( int i=0; i<GetCradlesAmount(); i++ )
 159   {
 160     float c                = distance_to_CPV,           // Distance to CPV
 161           l                = cell_side_size*GetNphi()/2,      // Cradle half size around beam (for rect. geom.)
 162           cradle_angle     = 360/kPI*atan(l/c),
 163           cradle_angle_pos = -90+(i-(GetCradlesAmount()-1)/2.) * (cradle_angle+GetAngleBetweenCradles());
 164     // Cradles are numerated in clock reversed order. (general way of angle increment)
 165
 166     float   r       = GetRadius() + cradle_thikness/2;
 167     x = r*cos(cradle_angle_pos*kPI/180);
 168     y = r*sin(cradle_angle_pos*kPI/180);
 169     z = 0;
 170     AliMatrix(rotation_matrix_number, 0,0 , 90,90+cradle_angle_pos , 90,180+cradle_angle_pos);
 171     pMC->Gspos("PHOS",i+1,"ALIC",x,y,z,rotation_matrix_number,"ONLY");
 172
 173     GetCradleAngle(i) = cradle_angle_pos;
 174 //
 175 //    int n = PHOS.fCradles->GetEntries();
 176 //    PHOS.fCradles->Add(new AliPHOSCradle( 1,            // geometry.
 177 //                                          GetCrystalSideSize    (),
 178 //                                          GetCrystalLength      (),
 179 //                                          GetWrapThickness      (),
 180 //                                          GetAirThickness       (),
 181 //                                          GetPIN_SideSize       (),
 182 //                                          GetPIN_Length         (),
 183 //                                          GetRadius             (),
 184 //                                          GetCPV_Thickness      (),
 185 //                                          GetCPV_PHOS_Distance  (),
 186 //                                          GetNz                 (),
 187 //                                          GetNphi               (),
 188 //                                          cradle_angle_pos      ));
 189 //
 190 //    if( n+1 != PHOS.fCradles->GetEntries() ||
 191 //        NULL == PHOS.fCradles->At(n) )
 192 //    {
 193 //      cout << "  Can not create or add AliPHOSCradle.\n";
 194 //      exit(1);
 195 //    }
 196   }
 197   AddPHOSCradles();
 198
 199   //////////////////////////////////////////////////////////////////////////////
 200   // All is done.
 201   // Print some information.
 202   //////////////////////////////////////////////////////////////////////////////
 203 }
 204
 205 void AliPHOSv1::StepManager()
 206 {
 207   static Bool_t inwold=0;   // Status of previous ctrak->inwvol
 208   AliMC *MC = AliMC::GetMC();
 209   Int_t copy;
 210
 211   int cradle_number, cell_Z, cell_Phi;  // Variables that describe cell position.
 212
 213   if( MC->GetMedium() == GetPHOS_IDTMED_PIN() && (MC->TrackInside() || MC->TrackExiting()==2) && inwold && MC->TrackCharge()!=0 )
 214   {
 215     // GEANT particle just have entered into PIN diode.
 216
 217     AliPHOS &PHOS = *(AliPHOS*)gAlice->GetModule("PHOS");
 218
 219     MC->CurrentVolOff(4,0,copy);
 220     cradle_number  = copy-1;
 221     MC->CurrentVolOff(1,0,copy);
 222     cell_Z         = copy-1;
 223     MC->CurrentVolOff(2,0,copy);
 224     cell_Phi       = copy-1;
 225 /*
 226         cradle_number  = cvolu->number[cvolu->nlevel-5]-1;
 227         cell_Z         = cvolu->number[cvolu->nlevel-2]-1;
 228         cell_Phi       = cvolu->number[cvolu->nlevel-3]-1;
 229 */
 230
 231     TH2S &h = PHOS.GetCradle(cradle_number).fChargedTracksInPIN;
 232     h.AddBinContent(h.GetBin(cell_Z,cell_Phi));
 233   }
 234
 235   //////////////////////////////////////////////////////////////////////////////
 236
 237   if( MC->GetMedium() == GetPHOS_IDTMED_PbWO4() )
 238   {
 239     // GEANT particle into crystal.
 240
 241     AliPHOS &PHOS = *(AliPHOS*)gAlice->GetModule("PHOS");
 242
 243     MC->CurrentVolOff(5,0,copy);
 244     cradle_number  = copy-1;
 245     MC->CurrentVolOff(2,0,copy);
 246     cell_Z         = copy-1;
 247     MC->CurrentVolOff(3,0,copy);
 248     cell_Phi       = copy-1;
 249 /*
 250         cradle_number  = cvolu->number[cvolu->nlevel-6]-1;
 251         cell_Z         = cvolu->number[cvolu->nlevel-3]-1;
 252         cell_Phi       = cvolu->number[cvolu->nlevel-4]-1;
 253 */
 254     TH2F &h = PHOS.GetCradle(cradle_number).fCellEnergy;
 255     h.AddBinContent(h.GetBin(cell_Z,cell_Phi),MC->Edep());
 256   }
 257
 258   //////////////////////////////////////////////////////////////////////////////
 259
 260   if( MC->GetMedium()==GetPHOS_IDTMED_CPV() && (MC->TrackInside() || MC->TrackExiting()) && inwold )
 261   {
 262     // GEANT particle just have entered into CPV detector.
 263
 264     AliPHOS &PHOS = *(AliPHOS*)gAlice->GetModule("PHOS");
 265
 266     MC->CurrentVolOff(1,0,cradle_number);
 267     cradle_number--;
 268 //        cradle_number  = cvolu->number[cvolu->nlevel-2]-1;
 269
 270     // Save CPV x,y hits position of charged particles.
 271
 272     AliPHOSCradle  &cradle = PHOS.GetCradle(cradle_number);
 273
 274     Float_t   xyz[3];
 275     MC->TrackPosition(xyz);
 276     TVector3          p(xyz[0],xyz[1],xyz[2]),v;
 277
 278     float x,y,l;
 279     float R = cradle.GetRadius() - cradle.GetCPV_PHOS_Distance() - cradle.GetCPV_Thikness();
 280     cradle.GetXY(p,v,R,x,y,l);
 281     if( PHOS.fDebugLevel>0 )
 282       if( l<0 )
 283         printf("PHOS_STEP:  warning: negative distance to CPV!! %f\n", l);
 284
 285     // Store current particle in the list of Cradle particles.
 286     Float_t  pmom[4];
 287     MC->TrackMomentum(pmom);
 288     float     Px      =       pmom[0] * pmom[3],
 289               Py      =       pmom[1] * pmom[3],
 290               Pz      =       pmom[2] * pmom[3];
 291     Int_t     Ipart   =       MC->TrackPid();
 292
 293 //     TClonesArray &P=cradle.GetParticles();
 294 //     new( P[P.GetEntries()] ) AliPHOSgamma(x,0,y,0,ctrak->getot,0,Px,Py,Pz);
 295     cradle.GetParticles().Add(new AliPHOSgamma(x,y,MC->Etot(),Px,Py,Pz,Ipart));
 296
 297     if( MC->TrackCharge()!=0 )
 298       cradle.AddCPVHit(x,y);
 299   }
 300
 301   inwold=MC->TrackEntering();         // Save current status of GEANT variable.
 302 }
 303