JETAN/AliJetDummyGeo.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 // Temporarily added to define part of the EMCal geometry
  18 // necessary for the jet finder
  19 //
  20 // M. Estienne
  21 // Magali.Estienne@cern.ch
  22 //
  23
  24 #include <Riostream.h>
  25 #include <assert.h>
  26 #include <TList.h>
  27
  28 // --- Root header files ---
  29 #include <TVector3.h>
  30 #include <TGeoManager.h>
  31 #include <TGeoMatrix.h>
  32 #include <TGeoNode.h>
  33
  34 // --- AliRoot header files ---
  35 #include "AliLog.h"
  36 #include "AliJetDummyGeo.h"
  37 #include "AliJetDummyShishKebabTrd1Module.h"
  38
  39 using std::cout;
  40 using std::endl;
  41 ClassImp(AliJetDummyGeo)
  42
  43 AliJetDummyGeo::AliJetDummyGeo():
  44   fArm1EtaMin(-0.7),
  45   fArm1EtaMax(+0.7),
  46   fArm1PhiMin(80.0),
  47   fArm1PhiMax(200.0),
  48   fNumberOfSuperModules(12),
  49   fSteelFrontThick(0.0),
  50   fLateralSteelStrip(0.01),
  51   fIPDistance(428.0),
  52   fPhiGapForSM(2.),
  53   fNPhi(12),
  54   fNZ(24),
  55   fPhiModuleSize(12.26 - fPhiGapForSM / Float_t(fNPhi)),
  56   fEtaModuleSize(fPhiModuleSize),
  57   fNPHIdiv(2),
  58   fNETAdiv(2),
  59   fPhiTileSize(fPhiModuleSize/Double_t(fNPHIdiv) - fLateralSteelStrip),
  60   fEtaTileSize(fEtaModuleSize/Double_t(fNETAdiv) - fLateralSteelStrip),
  61   fNECLayers(77),
  62   fECScintThick(0.16),
  63   fECPbRadThickness(0.16),
  64   fSampling(12.327),
  65   fTrd1Angle(1.5),
  66   fNCellsInModule(fNPHIdiv*fNETAdiv),
  67   fNCellsInSupMod(fNCellsInModule*fNPhi*fNZ),
  68   fNCells(fNCellsInSupMod*fNumberOfSuperModules-fNCellsInSupMod),
  69   fLongModuleSize(fNECLayers*(fECScintThick + fECPbRadThickness)),
  70   f2Trd1Dx2(fEtaModuleSize + 2.*fLongModuleSize*TMath::Tan(fTrd1Angle*TMath::DegToRad()/2.)),
  71   fShellThickness(TMath::Sqrt(fLongModuleSize*fLongModuleSize + f2Trd1Dx2*f2Trd1Dx2)+fSteelFrontThick),
  72   fZLength(2.*ZFromEtaR(fIPDistance+fShellThickness,fArm1EtaMax)),
  73   fEtaMaxOfTRD1(0.),
  74   fPhiBoundariesOfSM(0),
  75   fPhiCentersOfSM(0),
  76   fCentersOfCellsEtaDir(0),
  77   fCentersOfCellsXDir(0),
  78   fCentersOfCellsPhiDir(0),
  79   fEtaCentersOfCells(0),
  80   fPhiCentersOfCells(0),
  81   fShishKebabTrd1Modules(0),
  82   fDebug(0)
  83 {
  84   // Constructor
  85
  86   // Local coordinates
  87   fParSM[0] = GetShellThickness()/2.;
  88   fParSM[1] = GetPhiModuleSize() * GetNPhi()/2.;
  89   fParSM[2] = 350./2.;
  90
  91   fEnvelop[0] = fIPDistance; // mother volume inner radius
  92   fEnvelop[1] = fIPDistance + fShellThickness; // mother volume outer r.
  93   fEnvelop[2] = 1.00001*fZLength; // add some padding for mother volume.
  94
  95   // SM phi boundaries - (0,1),(2,3) .. (10,11) - has the same boundaries; Nov 7, 2006
  96   fPhiBoundariesOfSM.Set(fNumberOfSuperModules);
  97   fPhiCentersOfSM.Set(fNumberOfSuperModules/2);
  98   fPhiBoundariesOfSM[0] = TMath::PiOver2() - TMath::ATan2(fParSM[1] , fIPDistance); // 1th and 2th modules)
  99   fPhiBoundariesOfSM[1] = TMath::PiOver2() + TMath::ATan2(fParSM[1] , fIPDistance);
 100   fPhiCentersOfSM[0]    = TMath::PiOver2();
 101   for(Int_t i=1; i<=4; i++) { // from 2th ro 9th
 102     fPhiBoundariesOfSM[2*i]   = fPhiBoundariesOfSM[0] + 20.*TMath::DegToRad()*i;
 103     fPhiBoundariesOfSM[2*i+1] = fPhiBoundariesOfSM[1] + 20.*TMath::DegToRad()*i;
 104     fPhiCentersOfSM[i]        = fPhiCentersOfSM[0]    + 20.*TMath::DegToRad()*i;
 105   }
 106   fPhiBoundariesOfSM[11] = 190.*TMath::DegToRad();
 107   fPhiBoundariesOfSM[10] = fPhiBoundariesOfSM[11] - TMath::ATan2((fParSM[1]) , fIPDistance);
 108   fPhiCentersOfSM[5]     = (fPhiBoundariesOfSM[10]+fPhiBoundariesOfSM[11])/2.;
 109
 110   Int_t nphism  = GetNumberOfSuperModules()/2;
 111   Double_t dphi = (GetArm1PhiMax() - GetArm1PhiMin())/nphism;
 112   Double_t rpos = (GetEnvelop(0) + GetEnvelop(1))/2.;
 113   Double_t phi, phiRad, xpos, ypos, zpos;
 114   for(Int_t i=0; i<nphism; i++){
 115     phi    = GetArm1PhiMin() + dphi*(2*i+1)/2.; // phi= 90, 110, 130, 150, 170, 190
 116     phiRad = phi*TMath::Pi()/180.;
 117     xpos = rpos * TMath::Cos(phiRad);
 118     ypos = rpos * TMath::Sin(phiRad);
 119     zpos = fParSM[2];
 120     if(i==5) {
 121       xpos += (fParSM[1]/2. * TMath::Sin(phiRad));
 122       ypos -= (fParSM[1]/2. * TMath::Cos(phiRad));
 123     }
 124     // pozitive z
 125     Int_t ind = 2*i;
 126     TGeoRotation *geoRot0 = new TGeoRotation("geoRot0", 90.0, phi, 90.0, 90.0+phi, 0.0, 0.0);
 127     fMatrixOfSM[ind] = new TGeoCombiTrans(Form("EmcalSM%2.2i",ind),
 128                                           xpos,ypos, zpos, geoRot0);
 129     // negaive z
 130     ind++;
 131     Double_t phiy = 90. + phi + 180.;
 132     if(phiy>=360.) phiy -= 360.;
 133     TGeoRotation *geoRot1 = new TGeoRotation("geoRot1", 90.0, phi, 90.0, phiy, 180.0, 0.0);
 134     fMatrixOfSM[ind] = new TGeoCombiTrans(Form("EmcalSM%2.2i",ind),
 135                                           xpos,ypos,-zpos, geoRot1);
 136   } // for
 137
 138   CreateListOfTrd1Modules();
 139
 140   if(fDebug > 0){
 141     for(Int_t i=0; i<6; i++){
 142       cout << "fMatrixOfSM[" << i << "]: " << fMatrixOfSM[i] << endl;
 143     }
 144     cout << "fArm1EtaMin: " << fArm1EtaMin << endl;
 145     cout << "fArm1EtaMax: " << fArm1EtaMax << endl;
 146     cout << "fArm1PhiMin: " << fArm1PhiMin << endl;
 147     cout << "fArm1PhiMax: " << fArm1PhiMax << endl;
 148     cout << "fNumberOfSuperModules: " << fNumberOfSuperModules << endl;
 149     cout << "fSteelFrontThick: " << fSteelFrontThick << endl;
 150     cout << "fIPDistance: " << fIPDistance  << endl;
 151     cout << "fZLength: " << fZLength << endl;
 152     cout << "fPhiGapForSM: " << fPhiGapForSM << endl;
 153     cout << "fNPhi: " << fNPhi << endl;
 154     cout << "fNZ: " << fNZ << endl;
 155     cout << "fPhiModuleSize: " << fPhiModuleSize << endl;
 156     cout << "fEtaModuleSize: " << fEtaModuleSize << endl;
 157     cout << "fNPHIdiv: " << fNPHIdiv << endl;
 158     cout << "fNETAdiv: " << fNETAdiv << endl;
 159     cout << "fNECLayers: " << fNECLayers << endl;
 160     cout << "fECScintThick: " <<  fECScintThick << endl;
 161     cout << "fECPbRadThickness: " << fECPbRadThickness << endl;
 162     cout << "fSampling: " << fSampling << endl;
 163     cout << "fTrd1Angle: " << fTrd1Angle << endl;
 164     cout << "fNCellsInModule: " << fNCellsInModule << endl;
 165     cout << "fNCellsInSupMod: " << fNCellsInSupMod << endl;
 166     cout << "fNCells: " << fNCells << endl;
 167     cout << "fLongModuleSize: " <<  fLongModuleSize << endl;
 168     cout << "f2Trd1Dx2: " << f2Trd1Dx2 << endl;
 169     cout << "fShellThickness: " << fShellThickness << endl;
 170     cout << "fEtaMaxOfTRD1: " << fEtaMaxOfTRD1 << endl;
 171   }
 172 }
 173
 174 AliJetDummyGeo::AliJetDummyGeo(const AliJetDummyGeo& geom):
 175   TObject(geom),
 176   fArm1EtaMin(geom.fArm1EtaMin),
 177   fArm1EtaMax(geom.fArm1EtaMax),
 178   fArm1PhiMin(geom.fArm1PhiMin),
 179   fArm1PhiMax(geom.fArm1PhiMax),
 180   fNumberOfSuperModules(geom.fNumberOfSuperModules),
 181   fSteelFrontThick(geom.fSteelFrontThick),
 182   fLateralSteelStrip(geom.fLateralSteelStrip),
 183   fIPDistance(geom.fIPDistance),
 184   fPhiGapForSM(geom.fPhiGapForSM),
 185   fNPhi(geom.fNPhi),
 186   fNZ(geom.fNZ),
 187   fPhiModuleSize(geom.fPhiModuleSize),
 188   fEtaModuleSize(geom.fEtaModuleSize),
 189   fNPHIdiv(geom.fNPHIdiv),
 190   fNETAdiv(geom.fNETAdiv),
 191   fPhiTileSize(geom.fPhiTileSize),
 192   fEtaTileSize(geom.fEtaTileSize),
 193   fNECLayers(geom.fNECLayers),
 194   fECScintThick(geom.fECScintThick),
 195   fECPbRadThickness(geom.fECPbRadThickness),
 196   fSampling(geom.fSampling),
 197   fTrd1Angle(geom.fTrd1Angle),
 198   fNCellsInModule(geom.fNCellsInModule),
 199   fNCellsInSupMod(geom.fNCellsInSupMod),
 200   fNCells(geom.fNCells),
 201   fLongModuleSize(geom.fLongModuleSize),
 202   f2Trd1Dx2(geom.f2Trd1Dx2),
 203   fShellThickness(geom.fShellThickness),
 204   fZLength(geom.fZLength),
 205   fEtaMaxOfTRD1(geom.fEtaMaxOfTRD1),
 206   fPhiBoundariesOfSM(geom.fPhiBoundariesOfSM),
 207   fPhiCentersOfSM(geom.fPhiCentersOfSM),
 208   fCentersOfCellsEtaDir(geom.fCentersOfCellsEtaDir),
 209   fCentersOfCellsXDir(geom.fCentersOfCellsXDir),
 210   fCentersOfCellsPhiDir(geom.fCentersOfCellsPhiDir),
 211   fEtaCentersOfCells(geom.fEtaCentersOfCells),
 212   fPhiCentersOfCells(geom.fPhiCentersOfCells),
 213   fShishKebabTrd1Modules(geom.fShishKebabTrd1Modules),
 214   fDebug(geom.fDebug)
 215 {
 216   // Constructor
 217   // Local coordinates
 218   fParSM[0] = GetShellThickness()/2.;
 219   fParSM[1] = GetPhiModuleSize() * GetNPhi()/2.;
 220   fParSM[2] = 350./2.;
 221
 222   fEnvelop[0] = fIPDistance; // mother volume inner radius
 223   fEnvelop[1] = fIPDistance + fShellThickness; // mother volume outer r.
 224   fEnvelop[2] = 1.00001*fZLength; // add some padding for mother volume.
 225
 226   // SM phi boundaries - (0,1),(2,3) .. (10,11) - has the same boundaries; Nov 7, 2006
 227   fPhiBoundariesOfSM.Set(fNumberOfSuperModules);
 228   fPhiCentersOfSM.Set(fNumberOfSuperModules/2);
 229   fPhiBoundariesOfSM[0] = TMath::PiOver2() - TMath::ATan2(fParSM[1] , fIPDistance); // 1th and 2th modules)
 230   fPhiBoundariesOfSM[1] = TMath::PiOver2() + TMath::ATan2(fParSM[1] , fIPDistance);
 231   fPhiCentersOfSM[0]    = TMath::PiOver2();
 232   for(Int_t i=1; i<=4; i++) { // from 2th ro 9th
 233     fPhiBoundariesOfSM[2*i]   = fPhiBoundariesOfSM[0] + 20.*TMath::DegToRad()*i;
 234     fPhiBoundariesOfSM[2*i+1] = fPhiBoundariesOfSM[1] + 20.*TMath::DegToRad()*i;
 235     fPhiCentersOfSM[i]        = fPhiCentersOfSM[0]    + 20.*TMath::DegToRad()*i;
 236   }
 237   fPhiBoundariesOfSM[11] = 190.*TMath::DegToRad();
 238   fPhiBoundariesOfSM[10] = fPhiBoundariesOfSM[11] - TMath::ATan2((fParSM[1]) , fIPDistance);
 239   fPhiCentersOfSM[5]     = (fPhiBoundariesOfSM[10]+fPhiBoundariesOfSM[11])/2.;
 240
 241   Int_t nphism  = GetNumberOfSuperModules()/2;
 242   Double_t dphi = (GetArm1PhiMax() - GetArm1PhiMin())/nphism;
 243   Double_t rpos = (GetEnvelop(0) + GetEnvelop(1))/2.;
 244   Double_t phi, phiRad, xpos, ypos, zpos;
 245   for(Int_t i=0; i<nphism; i++){
 246     phi    = GetArm1PhiMin() + dphi*(2*i+1)/2.; // phi= 90, 110, 130, 150, 170, 190
 247     phiRad = phi*TMath::Pi()/180.;
 248     xpos = rpos * TMath::Cos(phiRad);
 249     ypos = rpos * TMath::Sin(phiRad);
 250     zpos = fParSM[2];
 251     if(i==5) {
 252       xpos += (fParSM[1]/2. * TMath::Sin(phiRad));
 253       ypos -= (fParSM[1]/2. * TMath::Cos(phiRad));
 254     }
 255     // pozitive z
 256     Int_t ind = 2*i;
 257     TGeoRotation *geoRot0 = new TGeoRotation("geoRot0", 90.0, phi, 90.0, 90.0+phi, 0.0, 0.0);
 258     fMatrixOfSM[ind] = new TGeoCombiTrans(Form("EmcalSM%2.2i",ind),
 259                                           xpos,ypos, zpos, geoRot0);
 260     // negaive z
 261     ind++;
 262     Double_t phiy = 90. + phi + 180.;
 263     if(phiy>=360.) phiy -= 360.;
 264     TGeoRotation *geoRot1 = new TGeoRotation("geoRot1", 90.0, phi, 90.0, phiy, 180.0, 0.0);
 265     fMatrixOfSM[ind] = new TGeoCombiTrans(Form("EmcalSM%2.2i",ind),
 266                                           xpos,ypos,-zpos, geoRot1);
 267
 268     delete geoRot0;
 269     delete geoRot1;
 270
 271   } // for
 272
 273   CreateListOfTrd1Modules();
 274
 275 }
 276
 277 //------------------------------------------------------------------------------------
 278 AliJetDummyGeo::~AliJetDummyGeo()
 279 {
 280   // Destructor
 281   // delete [] fMatrixOfSM;
 282 }
 283
 284 //------------------------------------------------------------------------------------
 285 void AliJetDummyGeo::EtaPhiFromIndex(Int_t absId, Float_t& eta, Float_t& phi) const
 286 {
 287   // Nov 16, 2006- float to double
 288   // version for TRD1 only
 289   static TVector3 vglob;
 290   GetGlobal(absId, vglob);
 291   eta = vglob.Eta();
 292   phi = vglob.Phi();
 293
 294 }
 295
 296 //------------------------------------------------------------------------------------
 297 void AliJetDummyGeo::GetGlobal(const Double_t *loc, Double_t *glob, Int_t ind) const
 298 {
 299   // Figure out the global numbering of a given supermodule from the
 300   // local numbering
 301   // Alice numbering - Jun 03,2006
 302   //  if(fMatrixOfSM[0] == 0) GetTransformationForSM();
 303
 304   if(ind>=0 && ind < GetNumberOfSuperModules()) {
 305     fMatrixOfSM[ind]->LocalToMaster(loc, glob);
 306   }
 307 }
 308
 309 //------------------------------------------------------------------------------------
 310 void AliJetDummyGeo::GetGlobal(Int_t absId , Double_t glob[3]) const
 311 {
 312   // Alice numbering scheme - Jun 03, 2006
 313   static Int_t nSupMod, nModule, nIphi, nIeta;
 314   static Double_t loc[3];
 315
 316   glob[0]=glob[1]=glob[2]=0.0; // bad case
 317   if(RelPosCellInSModule(absId, loc)) {
 318     GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
 319     fMatrixOfSM[nSupMod]->LocalToMaster(loc, glob);
 320   }
 321 }
 322
 323 //------------------------------------------------------------------------------------
 324 void AliJetDummyGeo::GetGlobal(Int_t absId , TVector3 &vglob) const
 325 {
 326   // Alice numbering scheme - Jun 03, 2006
 327   static Double_t glob[3];
 328
 329   GetGlobal(absId, glob);
 330   vglob.SetXYZ(glob[0], glob[1], glob[2]);
 331
 332 }
 333
 334 //------------------------------------------------------------------------------------
 335 Bool_t AliJetDummyGeo::RelPosCellInSModule(Int_t absId, Double_t loc[3]) const
 336 {
 337   // Alice numbering scheme - Jun 03, 2006
 338   loc[0] = loc[1] = loc[2]=0.0;
 339   if(RelPosCellInSModule(absId, loc[0],loc[1],loc[2])) {
 340     return kTRUE;
 341   }
 342   return kFALSE;
 343 }
 344
 345 //------------------------------------------------------------------------------------
 346 Bool_t AliJetDummyGeo::RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t &yr, Double_t &zr) const
 347 {
 348   // Look to see what the relative position inside a given cell is
 349   // for a recpoint.
 350   // Alice numbering scheme - Jun 08, 2006
 351   // In:
 352   // absId   - cell is as in Geant,     0<= absId   < fNCells;
 353   // OUT:
 354   // xr,yr,zr - x,y,z coordinates of cell with absId inside SM
 355
 356   // Shift index taking into account the difference between standard SM
 357   // and SM of half size in phi direction
 358   const Int_t phiIndexShift = fCentersOfCellsPhiDir.GetSize()/4; // Nov 22, 2006; was 6 for cas 2X2
 359   static Int_t nSupMod, nModule, nIphi, nIeta, iphi, ieta;
 360   if(!CheckAbsCellId(absId)) return kFALSE;
 361
 362   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
 363   GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta);
 364
 365   xr = fCentersOfCellsXDir.At(ieta);
 366   zr = fCentersOfCellsEtaDir.At(ieta);
 367
 368   if(nSupMod<10) {
 369     yr = fCentersOfCellsPhiDir.At(iphi);
 370   } else {
 371     yr = fCentersOfCellsPhiDir.At(iphi + phiIndexShift);
 372   }
 373
 374   return kTRUE;
 375 }
 376
 377 //------------------------------------------------------------------------------------
 378 Bool_t  AliJetDummyGeo::CheckAbsCellId(Int_t absId) const
 379 {
 380   // May 31, 2006; only trd1 now
 381   if(absId<0 || absId >= fNCells) return kFALSE;
 382   else                            return kTRUE;
 383 }
 384
 385 //------------------------------------------------------------------------------------
 386 Bool_t AliJetDummyGeo::GetCellIndex(Int_t absId,Int_t &nSupMod,Int_t &nModule,Int_t &nIphi,Int_t &nIeta) const
 387 {
 388   // 21-sep-04; 19-oct-05;
 389   // May 31, 2006; ALICE numbering scheme:
 390   //
 391   // In:
 392   // absId   - cell is as in Geant,     0<= absId   < fNCells;
 393   // Out:
 394   // nSupMod - super module(SM) number, 0<= nSupMod < fNumberOfSuperModules;
 395   // nModule  - module number in SM,    0<= nModule  < fNCellsInSupMod/fNCellsInSupMod or(/2) for tow last SM (10th and 11th);
 396   // nIphi   - cell number in phi driection inside module; 0<= nIphi < fNPHIdiv;
 397   // nIeta   - cell number in eta driection inside module; 0<= nIeta < fNETAdiv;
 398   //
 399   static Int_t tmp=0, sm10=0;
 400   if(!CheckAbsCellId(absId)) return kFALSE;
 401
 402   sm10 = fNCellsInSupMod*10;
 403   if(absId >= sm10) { // 110 degree case; last two supermodules
 404     nSupMod = (absId-sm10) / (fNCellsInSupMod/2) + 10;
 405     tmp     = (absId-sm10) % (fNCellsInSupMod/2);
 406   } else {
 407     nSupMod = absId / fNCellsInSupMod;
 408     tmp     = absId % fNCellsInSupMod;
 409   }
 410
 411   nModule  = tmp / fNCellsInModule;
 412   tmp     = tmp % fNCellsInModule;
 413   nIphi   = tmp / fNPHIdiv;
 414   nIeta   = tmp % fNPHIdiv;
 415
 416   return kTRUE;
 417 }
 418
 419 //------------------------------------------------------------------------------------
 420 void AliJetDummyGeo::GetCellPhiEtaIndexInSModule(Int_t nSupMod, Int_t nModule, Int_t nIphi, Int_t nIeta, Int_t &iphi, Int_t &ieta) const
 421 {
 422   //
 423   // Added nSupMod; Nov 25, 05
 424   // Alice numbering scheme  - Jun 01,2006
 425   // IN:
 426   // nSupMod - super module(SM) number, 0<= nSupMod < fNumberOfSuperModules;
 427   // nModule  - module number in SM,     0<= nModule  < fNCellsInSupMod/fNCellsInSupMod or(/2) for tow last SM (10th and 11th);
 428   // nIphi   - cell number in phi driection inside module; 0<= nIphi < fNPHIdiv;
 429   // nIeta   - cell number in eta driection inside module; 0<= nIeta < fNETAdiv;
 430   //
 431  // OUT:
 432   // ieta, iphi - indexes of cell(tower) in two dimensional grid of SM
 433   // ieta - have to change from 0 to (fNZ*fNETAdiv-1)
 434   // iphi - have to change from 0 to (fNPhi*fNPHIdiv-1 or fNPhi*fNPHIdiv/2-1)
 435   //
 436   static Int_t iphim, ietam;
 437
 438   GetModulePhiEtaIndexInSModule(nSupMod,nModule, iphim, ietam);
 439   //  ieta  = ietam*fNETAdiv + (1-nIeta); // x(module) = -z(SM)
 440   ieta  = ietam*fNETAdiv + (fNETAdiv - 1 - nIeta); // x(module) = -z(SM)
 441   iphi  = iphim*fNPHIdiv + nIphi;     // y(module) =  y(SM)
 442
 443 }
 444
 445
 446 //------------------------------------------------------------------------------------
 447 void AliJetDummyGeo::GetModulePhiEtaIndexInSModule(Int_t nSupMod, Int_t nModule,  Int_t &iphim, Int_t &ietam) const
 448 {
 449   // added nSupMod; - 19-oct-05 !
 450   // Alice numbering scheme        - Jun 01,2006
 451   // ietam, iphi - indexes of module in two dimensional grid of SM
 452   // ietam - have to change from 0 to fNZ-1
 453   // iphim - have to change from 0 to nphi-1 (fNPhi-1 or fNPhi/2-1)
 454   static Int_t nphi;
 455
 456   if(nSupMod>=10) nphi = fNPhi/2;
 457   else            nphi = fNPhi;
 458
 459   ietam = nModule/nphi;
 460   iphim = nModule%nphi;
 461 }
 462
 463 //------------------------------------------------------------------------------------
 464 Bool_t AliJetDummyGeo::GetAbsCellIdFromEtaPhi(Double_t eta, Double_t phi, Int_t &absId) const
 465 {
 466   // Nov 17,2006
 467   // stay here - phi problem as usual
 468   static Int_t nSupMod, i, ieta, iphi, etaShift, nphi;
 469   static Double_t absEta=0.0, d=0.0, dmin=0.0, phiLoc;
 470   absId = nSupMod = - 1;
 471   if(SuperModuleNumberFromEtaPhi(eta, phi, nSupMod)) {
 472     // phi index first
 473     phi    = TVector2::Phi_0_2pi(phi);
 474     phiLoc = phi - fPhiCentersOfSM[nSupMod/2];
 475     nphi   = fPhiCentersOfCells.GetSize();
 476     if(nSupMod>=10) {
 477       phiLoc = phi - 190.*TMath::DegToRad();
 478       nphi /= 2;
 479     }
 480
 481     dmin   = TMath::Abs(fPhiCentersOfCells[0]-phiLoc);
 482     iphi   = 0;
 483     for(i=1; i<nphi; i++) {
 484       d = TMath::Abs(fPhiCentersOfCells[i] - phiLoc);
 485       if(d < dmin) {
 486         dmin = d;
 487         iphi = i;
 488       }
 489     }
 490     // odd SM are turned with respect of even SM - reverse indexes
 491
 492     // eta index
 493     absEta   = TMath::Abs(eta);
 494     etaShift = iphi*fCentersOfCellsEtaDir.GetSize();
 495     dmin     = TMath::Abs(fEtaCentersOfCells[etaShift]-absEta);
 496     ieta     = 0;
 497     for(i=1; i<fCentersOfCellsEtaDir.GetSize(); i++) {
 498       d = TMath::Abs(fEtaCentersOfCells[i+etaShift] - absEta);
 499       if(d < dmin) {
 500         dmin = d;
 501         ieta = i;
 502       }
 503     }
 504
 505     if(eta<0) iphi = (nphi-1) - iphi;
 506     absId = GetAbsCellIdFromCellIndexes(nSupMod, iphi, ieta);
 507
 508     return kTRUE;
 509   }
 510   return kFALSE;
 511 }
 512
 513 //------------------------------------------------------------------------------------
 514 Bool_t AliJetDummyGeo::SuperModuleNumberFromEtaPhi(Double_t eta, Double_t phi, Int_t &nSupMod) const
 515 {
 516   // Return false if phi belongs a phi cracks between SM
 517
 518   static Int_t i;
 519
 520   if(TMath::Abs(eta) > fEtaMaxOfTRD1) return kFALSE;
 521
 522   phi = TVector2::Phi_0_2pi(phi); // move phi to (0,2pi) boundaries
 523   for(i=0; i<6; i++) {
 524     if(phi>=fPhiBoundariesOfSM[2*i] && phi<=fPhiBoundariesOfSM[2*i+1]) {
 525       nSupMod = 2*i;
 526       if(eta < 0.0) nSupMod++;
 527       return kTRUE;
 528     }
 529   }
 530   return kFALSE;
 531 }
 532
 533 //------------------------------------------------------------------------------------
 534 Int_t  AliJetDummyGeo::GetAbsCellIdFromCellIndexes(Int_t nSupMod, Int_t iphi, Int_t ieta) const
 535 {
 536   // Transition from super module number(nSupMod) and cell indexes (ieta,iphi) to absId
 537   static Int_t ietam, iphim, nModule;
 538   static Int_t nIeta, nIphi; // cell indexes in module
 539
 540   GetModuleIndexesFromCellIndexesInSModule(nSupMod, iphi, ieta, ietam, iphim, nModule);
 541
 542   nIeta = ieta%fNETAdiv;
 543   nIeta = fNETAdiv - 1 - nIeta;
 544   nIphi = iphi%fNPHIdiv;
 545
 546   return GetAbsCellId(nSupMod, nModule, nIphi, nIeta);
 547 }
 548
 549 //------------------------------------------------------------------------------------
 550 void  AliJetDummyGeo::GetModuleIndexesFromCellIndexesInSModule(Int_t nSupMod, Int_t iphi, Int_t ieta,
 551                         Int_t &iphim, Int_t &ietam, Int_t &nModule) const
 552 {
 553   // Transition from cell indexes (ieta,iphi) to module indexes (ietam,iphim, nModule)
 554   static Int_t nphi;
 555   nphi  = GetNumberOfModuleInPhiDirection(nSupMod);
 556
 557   ietam  = ieta/fNETAdiv;
 558   iphim  = iphi/fNPHIdiv;
 559   nModule = ietam * nphi + iphim;
 560 }
 561
 562 //------------------------------------------------------------------------------------
 563 Int_t AliJetDummyGeo::GetAbsCellId(Int_t nSupMod, Int_t nModule, Int_t nIphi, Int_t nIeta) const
 564 {
 565   // 27-aug-04;
 566   // corr. 21-sep-04;
 567   //       13-oct-05; 110 degree case
 568   // May 31, 2006; ALICE numbering scheme:
 569   // 0 <= nSupMod < fNumberOfSuperModules
 570   // 0 <= nModule  < fNPHI * fNZ ( fNPHI * fNZ/2 for fKey110DEG=1)
 571   // 0 <= nIphi   < fNPHIdiv
 572   // 0 <= nIeta   < fNETAdiv
 573   // 0 <= absid   < fNCells
 574   static Int_t id=0; // have to change from 0 to fNCells-1
 575   if(nSupMod >= 10) { // 110 degree case; last two supermodules
 576     id  = fNCellsInSupMod*10 + (fNCellsInSupMod/2)*(nSupMod-10);
 577   } else {
 578     id  = fNCellsInSupMod*nSupMod;
 579   }
 580   id += fNCellsInModule *nModule;
 581   id += fNPHIdiv *nIphi;
 582   id += nIeta;
 583   if(id<0 || id >= fNCells) {
 584     id = -TMath::Abs(id); // if negative something wrong
 585   }
 586   return id;
 587 }
 588
 589 //------------------------------------------------------------------------------------
 590 Bool_t AliJetDummyGeo::GetPhiBoundariesOfSMGap(Int_t nPhiSec, Double_t &phiMin, Double_t &phiMax) const
 591 {
 592   // 0<= nPhiSec <=4; phi in rad
 593   // 0;  gap boundaries between  0th&2th  | 1th&3th SM
 594   // 1;  gap boundaries between  2th&4th  | 3th&5th SM
 595   // 2;  gap boundaries between  4th&6th  | 5th&7th SM
 596   // 3;  gap boundaries between  6th&8th  | 7th&9th SM
 597   // 4;  gap boundaries between  8th&10th | 9th&11th SM
 598   if(nPhiSec<0 || nPhiSec >4) return kFALSE;
 599   phiMin = fPhiBoundariesOfSM[2*nPhiSec+1];
 600   phiMax = fPhiBoundariesOfSM[2*nPhiSec+2];
 601   return kTRUE;
 602 }
 603
 604 //------------------------------------------------------------------------------------
 605 void AliJetDummyGeo::CreateListOfTrd1Modules()
 606 {
 607   // Generate the list of Trd1 modules
 608   // which will make up the EMCAL
 609   // geometry
 610     printf("CreateListOfTrd1Modules() \n");
 611
 612   AliJetDummyShishKebabTrd1Module *mod=0, *mTmp=0; // current module
 613   if(fShishKebabTrd1Modules == 0) {
 614     fShishKebabTrd1Modules = new TList;
 615     fShishKebabTrd1Modules->SetName("ListOfTRD1");
 616     for(Int_t iz=0; iz< GetNZ(); iz++) {
 617       if(iz==0) {
 618         mod  = new AliJetDummyShishKebabTrd1Module(TMath::Pi()/2.,this);
 619       } else {
 620         mTmp  = new AliJetDummyShishKebabTrd1Module(*mod);
 621         mod   = mTmp;
 622       }
 623       fShishKebabTrd1Modules->Add(mod);
 624     }
 625   } else {
 626     printf(" Already exits : ");
 627   }
 628   mod = (AliJetDummyShishKebabTrd1Module*)fShishKebabTrd1Modules->At(fShishKebabTrd1Modules->GetSize()-1);
 629   fEtaMaxOfTRD1 = mod->GetMaxEtaOfModule(0);
 630
 631   //  printf(" fShishKebabTrd1Modules has %i modules : max eta %5.4f \n",
 632   //              fShishKebabTrd1Modules->GetSize(),fEtaMaxOfTRD1);
 633
 634   // Feb 20,2006;
 635   // Jun 01, 2006 - ALICE numbering scheme
 636   // define grid for cells in eta(z) and x directions in local coordinates system of SM
 637   // Works just for 2x2 case only -- ?? start here
 638   //
 639   //
 640   // Define grid for cells in phi(y) direction in local coordinates system of SM
 641   // as for 2X2 as for 3X3 - Nov 8,2006
 642   //
 643
 644   //  printf(" Cells grid in phi directions : size %i\n", fCentersOfCellsPhiDir.GetSize());
 645
 646   Int_t ind=0; // this is phi index
 647   Int_t ieta=0, nModule=0, iphiTemp;
 648   Double_t xr, zr, theta, phi, eta, r, x,y;
 649   xr = 0.;
 650   zr = 0.;
 651
 652   TVector3 vglob;
 653   Double_t ytCenterModule=0.0, ytCenterCell=0.0;
 654
 655   fCentersOfCellsPhiDir.Set(fNPhi*fNPHIdiv);
 656   fPhiCentersOfCells.Set(fNPhi*fNPHIdiv);
 657
 658   Double_t R0 = GetIPDistance() + GetLongModuleSize()/2.;
 659   for(Int_t it=0; it<fNPhi; it++) { // cycle on modules
 660     ytCenterModule = -fParSM[1] + fPhiModuleSize*(2*it+1)/2;  // center of module
 661     for(Int_t ic=0; ic<fNPHIdiv; ic++) { // cycle on cells in module
 662       if(fNPHIdiv==2) {
 663         ytCenterCell = ytCenterModule + fPhiTileSize *(2*ic-1)/2.;
 664       } else if(fNPHIdiv==3){
 665         ytCenterCell = ytCenterModule + fPhiTileSize *(ic-1);
 666       } else if(fNPHIdiv==1){
 667         ytCenterCell = ytCenterModule;
 668       }
 669       fCentersOfCellsPhiDir.AddAt(ytCenterCell,ind);
 670       // Define grid on phi direction
 671       // Grid is not the same for different eta bin;
 672       // Effect is small but is still here
 673       phi = TMath::ATan2(ytCenterCell, R0);
 674       fPhiCentersOfCells.AddAt(phi, ind);
 675
 676       //      printf(" ind %2.2i : y %8.3f ", ind, fCentersOfCellsPhiDir.At(ind));
 677       ind++;
 678     }
 679   }
 680
 681   fCentersOfCellsEtaDir.Set(fNZ *fNETAdiv);
 682   fCentersOfCellsXDir.Set(fNZ *fNETAdiv);
 683   fEtaCentersOfCells.Set(fNZ *fNETAdiv * fNPhi*fNPHIdiv);
 684   if(fDebug>1) AliInfo(Form(" Cells grid in eta directions : size %i\n", fCentersOfCellsEtaDir.GetSize()));
 685   for(Int_t it=0; it<fNZ; it++) {
 686     AliJetDummyShishKebabTrd1Module *trd1 = GetShishKebabModule(it);
 687     nModule = fNPhi*it;
 688     for(Int_t ic=0; ic<fNETAdiv; ic++) {
 689       if(fNPHIdiv==2) {
 690         trd1->GetCenterOfCellInLocalCoordinateofSM(ic, xr, zr);      // case of 2X2
 691         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
 692       } if(fNPHIdiv==3) {
 693         trd1->GetCenterOfCellInLocalCoordinateofSM3X3(ic, xr, zr);  // case of 3X3
 694         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
 695       } if(fNPHIdiv==1) {
 696         trd1->GetCenterOfCellInLocalCoordinateofSM1X1(xr, zr);      // case of 1X1
 697         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
 698       }
 699       fCentersOfCellsXDir.AddAt(Float_t(xr) - fParSM[0],ieta);
 700       fCentersOfCellsEtaDir.AddAt(Float_t(zr) - fParSM[2],ieta);
 701       // Define grid on eta direction for each bin in phi
 702       for(Int_t iphi=0; iphi<fCentersOfCellsPhiDir.GetSize(); iphi++) {
 703         x = xr + trd1->GetRadius();
 704         y = fCentersOfCellsPhiDir[iphi];
 705         r = TMath::Sqrt(x*x + y*y + zr*zr);
 706         theta = TMath::ACos(zr/r);
 707         eta   = AliJetDummyShishKebabTrd1Module::ThetaToEta(theta);
 708         //        ind   = ieta*fCentersOfCellsPhiDir.GetSize() + iphi;
 709         ind   = iphi*fCentersOfCellsEtaDir.GetSize() + ieta;
 710         fEtaCentersOfCells.AddAt(eta, ind);
 711       }
 712       //printf(" ieta %i : xr + trd1->GetRadius() %f : zr %f : eta %f \n", ieta, xr + trd1->GetRadius(), zr, eta);
 713     }
 714   }
 715
 716   if(fDebug>1){
 717     for(Int_t i=0; i<fCentersOfCellsEtaDir.GetSize(); i++) {
 718       AliInfo(Form(" ind %2.2i : z %8.3f : x %8.3f", i+1,
 719                    fCentersOfCellsEtaDir.At(i),fCentersOfCellsXDir.At(i)));
 720     }
 721   }
 722 }
 723
 724 //------------------------------------------------------------------------------------
 725 AliJetDummyShishKebabTrd1Module* AliJetDummyGeo::GetShishKebabModule(Int_t neta)
 726 {
 727   // This method was too long to be included in the header file - the
 728   // rule checker complained about it's length, so we move it here. It
 729   // returns the shishkebabmodule at a given eta index point.
 730
 731   static AliJetDummyShishKebabTrd1Module* trd1=0;
 732   if(fShishKebabTrd1Modules && neta>=0 && neta<fShishKebabTrd1Modules->GetSize()) {
 733     trd1 = (AliJetDummyShishKebabTrd1Module*)fShishKebabTrd1Modules->At(neta);
 734   } else trd1 = 0;
 735   return trd1;
 736 }
 737
 738 //------------------------------------------------------------------------------------
 739 void  AliJetDummyGeo::GetTransformationForSM()
 740 {
 741   // Uses the geometry manager to load the transformation matrix
 742   // for the supermodules
 743   // Unused after 19 Jan, 2007 - keep for compatibility;
 744
 745   return;
 746   static Bool_t transInit=kFALSE;
 747   if(transInit) return;
 748
 749   Int_t i=0;
 750   if(gGeoManager == 0) {
 751     Info("CreateTransformationForSM() "," Load geometry : TGeoManager::Import()");
 752     assert(0);
 753   }
 754
 755   TGeoNode *tn = gGeoManager->GetTopNode();
 756   TGeoNode *node=0, *xen1 = 0;
 757   for(i=0; i<tn->GetNdaughters(); i++) {
 758     node = tn->GetDaughter(i);
 759     TString ns(node->GetName());
 760     if(ns.Contains(GetNameOfEMCALEnvelope())) {
 761       xen1 = node;
 762       break;
 763     }
 764   }
 765
 766   if(!xen1) {
 767     Info("CreateTransformationForSM() "," geometry has not EMCAL envelope with name %s",
 768     GetNameOfEMCALEnvelope());
 769     assert(0);
 770   }
 771   AliInfo(Form(" i %i : EMCAL Envelope is %s : #SM %i \n", i, xen1->GetName(), xen1->GetNdaughters()));
 772   for(i=0; i<xen1->GetNdaughters(); i++) {
 773     TGeoNodeMatrix *sm = (TGeoNodeMatrix*)xen1->GetDaughter(i);
 774     fMatrixOfSM[i] = sm->GetMatrix();
 775   }
 776   AliInfo(Form("transInit %d: ", transInit));
 777   transInit = kTRUE;
 778 }
 779