LCOV - code coverage report
Current view: top level - EMCAL/EMCALUtils - AliEMCALGeometry.cxx (source / functions) Hit Total Coverage
Test: coverage.info Lines: 295 724 40.7 %
Date: 2016-06-14 17:26:59 Functions: 23 55 41.8 %

          Line data    Source code
       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             : // --- ROOT system ---
      17             : #include <TParticle.h>
      18             : #include <TGeoManager.h>
      19             : #include <TGeoMatrix.h>
      20             : #include <TGeoBBox.h>
      21             : #include <TList.h>
      22             : #include <TBrowser.h>
      23             : 
      24             : // --- Standard library ---
      25             : //#include <Riostream.h>
      26             : 
      27             : // --- AliRoot header files ---
      28             : #include "AliLog.h"
      29             : #include "AliEMCALGeometry.h"
      30             : #include "AliEMCALShishKebabTrd1Module.h"
      31             : #include "AliEMCALTriggerMappingV1.h"
      32             : #include "AliEMCALTriggerMappingV2.h"
      33             : 
      34          72 : ClassImp(AliEMCALGeometry)
      35             : 
      36             : // these initialisations are needed for a singleton
      37             : AliEMCALGeometry  *AliEMCALGeometry::fgGeom      = 0;
      38             : const Char_t*      AliEMCALGeometry::fgkDefaultGeometryName = "EMCAL_COMPLETE12SMV1_DCAL_8SM";
      39             : 
      40             : //____________________________________________________________________________
      41           0 : AliEMCALGeometry::AliEMCALGeometry():
      42           0 :   fEMCGeometry(0x0),fTriggerMapping(0x0),fGeoName(0),
      43           0 :   fKey110DEG(0),fnSupModInDCAL(0),fNCellsInSupMod(0),fNETAdiv(0),fNPHIdiv(0),
      44           0 :   fNCellsInModule(0),fPhiBoundariesOfSM(0x0),fPhiCentersOfSM(0x0),
      45           0 :   fPhiCentersOfSMSec(0x0),fPhiCentersOfCells(0x0),fCentersOfCellsEtaDir(0x0),
      46           0 :   fCentersOfCellsPhiDir(0x0),fEtaCentersOfCells(0x0),
      47           0 :   fNCells(0),fNPhi(0),fCentersOfCellsXDir(0x0),fArm1EtaMin(0),
      48           0 :   fArm1EtaMax(0),fArm1PhiMin(0),fArm1PhiMax(0),fEtaMaxOfTRD1(0),
      49           0 :   fDCALPhiMin(0),fDCALPhiMax(0),fEMCALPhiMax(0),fDCALStandardPhiMax(0),
      50           0 :   fDCALInnerExtandedEta(0),fShishKebabTrd1Modules(0),fPhiModuleSize(0.),
      51           0 :   fEtaModuleSize(0.),fPhiTileSize(0.),fEtaTileSize(0.),fNZ(0),
      52           0 :   fIPDistance(0.),fLongModuleSize(0.),fShellThickness(0.),
      53           0 :   fZLength(0.),fSampling(0.),fUseExternalMatrices(kFALSE)
      54           0 : {
      55             :   // default ctor 
      56             :   // must be kept public for root persistency purposes, but should never be called by the outside world
      57           0 :   fEnvelop[0] = 0.;
      58           0 :   fEnvelop[1] = 0.;
      59           0 :   fEnvelop[2] = 0.;
      60           0 :   fParSM[0]   = 0.;
      61           0 :   fParSM[1]   = 0.;
      62           0 :   fParSM[2]   = 0.;
      63           0 :   for (Int_t i=0;i<AliEMCALGeoParams::fgkEMCALModules;i++)
      64           0 :     fkSModuleMatrix[i]=0 ;
      65           0 : }  
      66             : 
      67             : //____________________________________________________________________________
      68             : AliEMCALGeometry::AliEMCALGeometry(const AliEMCALGeometry & geo)
      69           0 :   : TNamed(geo),
      70           0 :     fEMCGeometry(geo.fEMCGeometry),fTriggerMapping(geo.fTriggerMapping),fGeoName(geo.fGeoName),
      71           0 :     fKey110DEG(geo.fKey110DEG),fnSupModInDCAL(geo.fnSupModInDCAL),fNCellsInSupMod(geo.fNCellsInSupMod),fNETAdiv(geo.fNETAdiv),fNPHIdiv(geo.fNPHIdiv),
      72           0 :     fNCellsInModule(geo.fNCellsInModule),fPhiBoundariesOfSM(geo.fPhiBoundariesOfSM),fPhiCentersOfSM(geo.fPhiCentersOfSM),
      73           0 :     fPhiCentersOfSMSec(geo.fPhiCentersOfSMSec),fPhiCentersOfCells(geo.fPhiCentersOfCells),fCentersOfCellsEtaDir(geo.fCentersOfCellsEtaDir),
      74           0 :     fCentersOfCellsPhiDir(geo.fCentersOfCellsPhiDir),fEtaCentersOfCells(geo.fEtaCentersOfCells),
      75           0 :     fNCells(geo.fNCells),fNPhi(geo.fNPhi),fCentersOfCellsXDir(geo.fCentersOfCellsXDir),fArm1EtaMin(geo.fArm1EtaMin),
      76           0 :     fArm1EtaMax(geo.fArm1EtaMax),fArm1PhiMin(geo.fArm1PhiMin),fArm1PhiMax(geo.fArm1PhiMax),fEtaMaxOfTRD1(geo.fEtaMaxOfTRD1),
      77           0 :     fDCALPhiMin(geo.fDCALPhiMin),fDCALPhiMax(geo.fDCALPhiMax),fEMCALPhiMax(geo.fEMCALPhiMax),fDCALStandardPhiMax(geo.fDCALStandardPhiMax),
      78           0 :     fDCALInnerExtandedEta(geo.fDCALInnerExtandedEta),fShishKebabTrd1Modules(geo.fShishKebabTrd1Modules),fPhiModuleSize(geo.fPhiModuleSize),
      79           0 :     fEtaModuleSize(geo.fEtaModuleSize),fPhiTileSize(geo.fPhiTileSize),fEtaTileSize(geo.fEtaTileSize),fNZ(geo.fNZ),
      80           0 :     fIPDistance(geo.fIPDistance),fLongModuleSize(geo.fLongModuleSize),fShellThickness(geo.fShellThickness),
      81           0 :     fZLength(geo.fZLength),fSampling(geo.fSampling),fUseExternalMatrices(geo.fUseExternalMatrices)
      82           0 : {
      83             :   // Copy constarctor
      84           0 :   fEnvelop[0] = geo.fEnvelop[0];
      85           0 :   fEnvelop[1] = geo.fEnvelop[1];
      86           0 :   fEnvelop[2] = geo.fEnvelop[2];
      87           0 :   fParSM[0]   = geo.fParSM[0];
      88           0 :   fParSM[1]   = geo.fParSM[1];
      89           0 :   fParSM[2]   = geo.fParSM[2];
      90           0 :   for (Int_t i=0;i<AliEMCALGeoParams::fgkEMCALModules;i++)
      91           0 :     fkSModuleMatrix[i]=0 ;
      92           0 : }
      93             : 
      94             : //____________________________________________________________________________
      95             : AliEMCALGeometry::AliEMCALGeometry(const Text_t* name,   const Text_t* title,
      96             :                                    const Text_t* mcname, const Text_t* mctitle) 
      97           5 :   : TNamed(name, title),
      98           5 :     fEMCGeometry(0x0),fTriggerMapping(0x0),fGeoName(0),
      99           5 :     fKey110DEG(0),fnSupModInDCAL(0),fNCellsInSupMod(0),fNETAdiv(0),fNPHIdiv(0),
     100          10 :     fNCellsInModule(0),fPhiBoundariesOfSM(0x0),fPhiCentersOfSM(0x0),
     101          15 :     fPhiCentersOfSMSec(0x0),fPhiCentersOfCells(0x0),fCentersOfCellsEtaDir(0x0),
     102          10 :     fCentersOfCellsPhiDir(0x0),fEtaCentersOfCells(0x0),
     103          10 :     fNCells(0),fNPhi(0),fCentersOfCellsXDir(0x0),fArm1EtaMin(0),
     104           5 :     fArm1EtaMax(0),fArm1PhiMin(0),fArm1PhiMax(0),fEtaMaxOfTRD1(0),
     105           5 :     fDCALPhiMin(0),fDCALPhiMax(0),fEMCALPhiMax(0),fDCALStandardPhiMax(0),
     106           5 :     fDCALInnerExtandedEta(0),fShishKebabTrd1Modules(0),fPhiModuleSize(0.),
     107           5 :     fEtaModuleSize(0.),fPhiTileSize(0.),fEtaTileSize(0.),fNZ(0),
     108           5 :     fIPDistance(0.),fLongModuleSize(0.),fShellThickness(0.),
     109           5 :     fZLength(0.),fSampling(0.), fUseExternalMatrices(kFALSE)
     110          25 : { 
     111             :   // ctor only for normal usage 
     112             :   
     113          15 :   fEMCGeometry = new AliEMCALEMCGeometry(name,title,mcname,mctitle);
     114          15 :   fGeoName = fEMCGeometry->GetGeoName();
     115           5 :   fKey110DEG = fEMCGeometry->GetKey110DEG();
     116           5 :   fnSupModInDCAL = fEMCGeometry->GetnSupModInDCAL();
     117           5 :   fNCellsInSupMod = fEMCGeometry->GetNCellsInSupMod();
     118           5 :   fNETAdiv = fEMCGeometry->GetNETAdiv();
     119           5 :   fNPHIdiv = fEMCGeometry->GetNPHIdiv();
     120           5 :   fNCellsInModule = fNPHIdiv*fNETAdiv;
     121             :   static int i=0;
     122           5 :   Int_t nSMod = fEMCGeometry->GetNumberOfSuperModules();
     123           5 :   fPhiBoundariesOfSM.Set(nSMod);
     124           5 :   fPhiCentersOfSM.Set(nSMod/2);
     125           5 :   fPhiCentersOfSMSec.Set(nSMod/2);
     126         150 :   for(Int_t sm=0; sm<nSMod; sm++) {
     127          70 :     i = sm/2;
     128         210 :     fEMCGeometry->GetPhiBoundariesOfSM(sm,fPhiBoundariesOfSM[2*i],fPhiBoundariesOfSM[2*i+1]);
     129             :   }
     130             : 
     131           5 :   Double_t phiMin =  0.;
     132           5 :   Double_t phiMax =  0.;
     133         150 :   for(Int_t sm=0; sm<nSMod; sm++) {
     134          70 :     fEMCGeometry->GetPhiBoundariesOfSM(sm,phiMin,phiMax);
     135          70 :     i=sm/2;
     136         210 :     fPhiCentersOfSM[i] = fEMCGeometry->GetPhiCenterOfSM(sm);
     137         210 :     fPhiCentersOfSMSec[i] = fEMCGeometry->GetPhiCenterOfSMSec(sm);
     138             :   }
     139           5 :   fNCells = fEMCGeometry->GetNCells();
     140           5 :   fNPhi = fEMCGeometry->GetNPhi();
     141           5 :   fEnvelop[0] = fEMCGeometry->GetEnvelop(0);
     142           5 :   fEnvelop[1] = fEMCGeometry->GetEnvelop(1);
     143           5 :   fEnvelop[2] = fEMCGeometry->GetEnvelop(2);
     144           5 :   fParSM[0]   = fEMCGeometry->GetSuperModulesPar(0);
     145           5 :   fParSM[1]   = fEMCGeometry->GetSuperModulesPar(1);
     146           5 :   fParSM[2]   = fEMCGeometry->GetSuperModulesPar(2);
     147           5 :   fArm1EtaMin = fEMCGeometry->GetArm1EtaMin();
     148           5 :   fArm1EtaMax = fEMCGeometry->GetArm1EtaMax();
     149           5 :   fArm1PhiMin = fEMCGeometry->GetArm1PhiMin();
     150           5 :   fArm1PhiMax = fEMCGeometry->GetArm1PhiMax();
     151           5 :   fDCALPhiMin = fEMCGeometry->GetDCALPhiMin();
     152           5 :   fDCALPhiMax = fEMCGeometry->GetDCALPhiMax();
     153           5 :   fEMCALPhiMax = fEMCGeometry->GetEMCALPhiMax();
     154           5 :   fDCALStandardPhiMax = fEMCGeometry->GetDCALStandardPhiMax();
     155           5 :   fDCALInnerExtandedEta = fEMCGeometry->GetDCALInnerExtandedEta();
     156           5 :   fShellThickness = fEMCGeometry->GetShellThickness();
     157           5 :   fZLength    = fEMCGeometry->GetZLength();
     158           5 :   fSampling   = fEMCGeometry->GetSampling();
     159           5 :   fEtaModuleSize = fEMCGeometry->GetEtaModuleSize();
     160           5 :   fPhiModuleSize = fEMCGeometry->GetPhiModuleSize();
     161           5 :   fEtaTileSize = fEMCGeometry->GetEtaTileSize();
     162           5 :   fPhiTileSize = fEMCGeometry->GetPhiTileSize();
     163           5 :   fNZ          = fEMCGeometry->GetNZ();
     164           5 :   fIPDistance  = fEMCGeometry->GetIPDistance();
     165           5 :   fLongModuleSize = fEMCGeometry->GetLongModuleSize();
     166             : 
     167           5 :   CreateListOfTrd1Modules();
     168             : 
     169         230 :   for(Int_t smod=0; smod < AliEMCALGeoParams::fgkEMCALModules; smod++)
     170         110 :     fkSModuleMatrix[smod]=0 ;   
     171             :         
     172          20 :   if (AliDebugLevel()>=2) {
     173           0 :     fEMCGeometry->Print();
     174           0 :     PrintGeometryGeoUtils();
     175             :   }
     176             :   
     177          10 :   AliLog::Message(AliLog::kInfo, Form("Name <<%s>>",name),
     178             :                   MODULENAME(), "AliEMCALGeometry", FUNCTIONNAME(), __FILE__, __LINE__);  
     179             :   
     180          15 :   if ((fEMCGeometry->GetGeoName()).Contains("DCAL")) {
     181           6 :     fTriggerMapping = new AliEMCALTriggerMappingV2(52, this);
     182           2 :     AliLog::Message(AliLog::kInfo, "EMCAL Trigger Mapping Version V2 Enabled",
     183             :                     MODULENAME(), "AliEMCALGeometry", FUNCTIONNAME(), __FILE__, __LINE__);
     184             : 
     185             :   } else { 
     186           9 :     fTriggerMapping = new AliEMCALTriggerMappingV1(32, this);
     187           3 :     AliLog::Message(AliLog::kInfo, "EMCAL Trigger Mapping Version V1 Enabled",
     188             :                     MODULENAME(), "AliEMCALGeometry", FUNCTIONNAME(), __FILE__, __LINE__);
     189             :   }
     190          10 : }
     191             : 
     192             : //____________________________________________________________________________
     193             : AliEMCALGeometry & AliEMCALGeometry::operator = (const AliEMCALGeometry  & /*rvalue*/) 
     194             : { 
     195             :   //assing operator
     196           0 :   Fatal("assignment operator", "not implemented") ; 
     197           0 :   return *this ;
     198             : }
     199             : 
     200             : //____________________________________________________________________________
     201             : AliEMCALGeometry::~AliEMCALGeometry(void)
     202           0 : {
     203             :   // dtor
     204           0 :   if (this==fgGeom)
     205             :   {
     206           0 :     AliError("Do not call delete on me");
     207           0 :     return;
     208             :   }
     209             :   
     210           0 :   if (fEMCGeometry)
     211             :   {
     212           0 :     for(Int_t smod = 0 ; smod < fEMCGeometry->GetNumberOfSuperModules(); smod++)
     213             :     {
     214           0 :       if(fkSModuleMatrix[smod])
     215           0 :         delete fkSModuleMatrix[smod] ;
     216             :       
     217           0 :       fkSModuleMatrix[smod]=0 ;
     218             :     }
     219             :     
     220           0 :     delete fEMCGeometry; // fEMCGeometry = 0 ;
     221             :   }
     222             :   
     223           0 :   if (fTriggerMapping) delete fTriggerMapping;
     224           0 : }
     225             : 
     226             : //______________________________________________________________________
     227             : AliEMCALGeometry *  AliEMCALGeometry::GetInstance()
     228             : { 
     229             :   // Returns the pointer of the unique instance
     230             :   
     231      143214 :   AliEMCALGeometry * rv = static_cast<AliEMCALGeometry *>( fgGeom );
     232       71607 :   return rv; 
     233             : }
     234             : 
     235             : //
     236             : /// \return the pointer of the unique instance
     237             : ///
     238             : //______________________________________________________________________
     239             : AliEMCALGeometry* AliEMCALGeometry::GetInstance(const Text_t* name,   const Text_t* title,
     240             :                                                 const Text_t* mcname, const Text_t* mctitle )
     241             : {
     242             :   AliEMCALGeometry * rv = 0; 
     243             :   
     244          18 :   if ( fgGeom == 0 ) 
     245             :   {
     246          10 :     if ( strcmp(name,"") == 0 ) 
     247             :     { // get default geometry
     248           5 :       fgGeom = new AliEMCALGeometry(fgkDefaultGeometryName, title, mcname, mctitle);
     249           0 :     } 
     250             :     else 
     251             :     {
     252          10 :       fgGeom = new AliEMCALGeometry(name, title,mcname,mctitle);
     253             :     }  // end if strcmp(name,"")
     254             :     
     255           5 :     if ( AliEMCALEMCGeometry::fgInit ) 
     256             :     {
     257           5 :       rv = (AliEMCALGeometry * ) fgGeom;
     258           5 :     }
     259             :     else 
     260             :     {
     261             :       rv = 0; 
     262           0 :       delete fgGeom; 
     263           0 :       fgGeom = 0; 
     264             :     } // end if fgInit
     265             :   }
     266             :   else
     267             :   {
     268           4 :     if ( strcmp(fgGeom->GetName(), name) != 0)
     269             :     {
     270           4 :       printf("\n current geometry is %s : ", fgGeom->GetName());
     271           4 :       printf(" you should not call %s \n",name);
     272           4 :     } // end 
     273             :     
     274           4 :     rv = (AliEMCALGeometry *) fgGeom; 
     275             :   }  // end if fgGeom
     276             :   
     277           9 :   return rv; 
     278           0 : }
     279             : 
     280             : ///
     281             : /// Instanciate geometry depending on the run number
     282             : ///
     283             : /// \return the pointer of the unique instance
     284             : ///
     285             : //___________________________________________________________________________
     286             : AliEMCALGeometry* AliEMCALGeometry::GetInstanceFromRunNumber(Int_t runNumber,       TString geoName,
     287             :                                                              const Text_t* mcname,  const Text_t* mctitle )
     288             : {
     289             :   //printf("AliEMCALGeometry::GetInstanceFromRunNumber() - run %d, geoName <<%s>> \n",runNumber,geoName.Data());
     290             :   
     291          18 :   static bool showInfo = !(getenv("HLT_ONLINE_MODE") && strcmp(getenv("HLT_ONLINE_MODE"), "on") == 0);
     292             :   
     293           3 :   if      ( runNumber >= 104064 && runNumber < 140000  ) 
     294             :   {
     295             :     // 2009-2010 runs
     296             :     // First year geometry, 4 SM.
     297             :     
     298           0 :     if (showInfo)
     299             :     {
     300           0 :       if(!geoName.Contains("FIRSTYEARV1") && geoName!="")
     301             :       { 
     302           0 :         printf("AliEMCALGeometry::GetInstanceFromRunNumber() *** ATTENTION *** \n");
     303           0 :         printf("\t Specified geometry name <<%s>> for run %d is not considered! \n",geoName.Data(),runNumber);
     304           0 :         printf("\t In use <<EMCAL_FIRSTYEARV1>>, check run number and year\n");
     305           0 :       }
     306             :       else 
     307             :       {
     308           0 :         printf("AliEMCALGeometry::GetInstanceFromRunNumber() - Initialized geometry with name <<EMCAL_FIRSTYEARV1>>\n");
     309             :       }
     310             :     }
     311             :     
     312           0 :     return AliEMCALGeometry::GetInstance("EMCAL_FIRSTYEARV1","EMCAL",mcname,mctitle) ;
     313             :   }
     314           3 :   else if ( runNumber >= 140000 && runNumber <= 170593 )
     315             :   {
     316             :     // Almost complete EMCAL geometry, 10 SM. Year 2011 configuration
     317             :     
     318           0 :     if (showInfo)
     319             :     {
     320           0 :       if(!geoName.Contains("COMPLETEV1") && geoName!="")
     321             :       {
     322           0 :         printf("AliEMCALGeometry::GetInstanceFromRunNumber() *** ATTENTION *** \n");
     323           0 :         printf("\t Specified geometry name <<%s>> for run %d is not considered! \n",geoName.Data(),runNumber);
     324           0 :         printf("\t In use <<EMCAL_COMPLETEV1>>, check run number and year\n");
     325           0 :       }
     326             :       else 
     327             :       {
     328           0 :         printf("AliEMCALGeometry::GetInstanceFromRunNumber() - Initialized geometry with name <<EMCAL_COMPLETEV1>>\n");
     329             :       }
     330             :     }
     331           0 :     return AliEMCALGeometry::GetInstance("EMCAL_COMPLETEV1","EMCAL",mcname,mctitle) ;
     332             :   }
     333           6 :   else if ( runNumber >  176000 && runNumber <= 197692 )
     334             :   {
     335             :     // Complete EMCAL geometry, 12 SM. Year 2012 and on
     336             :     // The last 2 SM were not active, anyway they were there.
     337             :     
     338           3 :     if (showInfo)
     339             :     {
     340           0 :       if(!geoName.Contains("COMPLETE12SMV1") && geoName!="")
     341             :       {
     342           0 :         printf("AliEMCALGeometry::GetInstanceFromRunNumber() *** ATTENTION *** \n");
     343           0 :         printf("\t Specified geometry name <<%s>> for run %d is not considered! \n",geoName.Data(),runNumber);
     344           0 :         printf("\t In use <<EMCAL_COMPLETE12SMV1>>, check run number and year\n");
     345           0 :       }
     346             :       else 
     347             :       {
     348           0 :         printf("AliEMCALGeometry::GetInstanceFromRunNumber() - Initialized geometry with name <<EMCAL_COMPLETE12SMV1>>\n");
     349             :       }
     350             :     }
     351           0 :     return AliEMCALGeometry::GetInstance("EMCAL_COMPLETE12SMV1","EMCAL",mcname,mctitle) ;
     352             :   }
     353             :   else // Run 2
     354             :   {
     355             :     // EMCAL + DCAL geometry, 20 SM. Year 2015 and on
     356             :     
     357           3 :     if (showInfo)
     358             :     {
     359           5 :       if(!geoName.Contains("DCAL_8SM") && geoName!="")
     360             :       {
     361           1 :         printf("AliEMCALGeometry::GetInstanceFromRunNumber() *** ATTENTION *** \n");
     362           1 :         printf("\t Specified geometry name <<%s>> for run %d is not considered! \n",geoName.Data(),runNumber);
     363           1 :         printf("\t In use <<EMCAL_COMPLETE12SMV1_DCAL_8SM>>, check run number and year\n");
     364           1 :       }
     365             :       else 
     366             :       {
     367           2 :         printf("AliEMCALGeometry::GetInstanceFromRunNumber() - Initialized geometry with name <<EMCAL_COMPLETE12SMV1_DCAL_8SM>>\n");
     368             :       }
     369             :     }
     370           3 :     return AliEMCALGeometry::GetInstance("EMCAL_COMPLETE12SMV1_DCAL_8SM","EMCAL",mcname,mctitle) ;
     371             :   }  
     372           3 : }
     373             : 
     374             : //________________________________________________________________________________________________
     375             : void AliEMCALGeometry::Browse(TBrowser* b)
     376             : {
     377             :   //Browse the modules
     378           0 :   if(fShishKebabTrd1Modules) b->Add(fShishKebabTrd1Modules);
     379           0 : }
     380             : 
     381             : //________________________________________________________________________________________________
     382             : Bool_t AliEMCALGeometry::IsFolder() const
     383             : {
     384             :   //Check if fShishKebabTrd1Modules is in folder
     385           0 :   if(fShishKebabTrd1Modules) return kTRUE;
     386           0 :   else                       return kFALSE;
     387           0 : }
     388             : 
     389             : //________________________________________________________________________________________________
     390             : void AliEMCALGeometry::GetGlobal(const Double_t *loc, Double_t *glob, int ind) const
     391             : {
     392             :   // Figure out the global numbering
     393             :   // of a given supermodule from the
     394             :   // local numbering and the transformation
     395             :   // matrix stored by the geometry manager (allows for misaligned
     396             :   // geometry)
     397             :         
     398         194 :   const TGeoHMatrix* m = GetMatrixForSuperModule(ind);
     399          97 :   if(m) {
     400          97 :     m->LocalToMaster(loc, glob);
     401          97 :   } else {
     402           0 :     AliFatal("Geo matrixes are not loaded \n") ;
     403             :   }
     404          97 : }
     405             : 
     406             : //________________________________________________________________________________________________
     407             : void AliEMCALGeometry::GetGlobal(const TVector3 &vloc, TVector3 &vglob, int ind) const
     408             : {
     409             :   //Figure out the global numbering
     410             :   //of a given supermodule from the
     411             :   //local numbering given a 3-vector location
     412             : 
     413             :   static Double_t tglob[3], tloc[3];
     414           0 :   vloc.GetXYZ(tloc);
     415           0 :   GetGlobal(tloc, tglob, ind);
     416           0 :   vglob.SetXYZ(tglob[0], tglob[1], tglob[2]);
     417           0 : }
     418             : 
     419             : //________________________________________________________________________________________________
     420             : void AliEMCALGeometry::GetGlobal(Int_t absId , double glob[3]) const
     421             : {
     422             :   // Alice numbering scheme - Jun 03, 2006
     423             :   static Int_t nSupMod=-1, nModule=-1, nIphi=-1, nIeta=-1;
     424             :   static double loc[3];
     425             : 
     426           0 :   glob[0]=glob[1]=glob[2]=0.0; // bad case
     427           0 :   if(RelPosCellInSModule(absId, loc)) {
     428           0 :     GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
     429           0 :     const TGeoHMatrix* m = GetMatrixForSuperModule(nSupMod);
     430           0 :     if(m) {
     431           0 :       m->LocalToMaster(loc, glob);
     432           0 :     } else {
     433           0 :       AliFatal("Geo matrixes are not loaded \n") ;
     434             :     }
     435           0 :   }
     436           0 : }
     437             : 
     438             : //___________________________________________________________________
     439             : void AliEMCALGeometry::GetGlobal(Int_t absId , TVector3 &vglob) const
     440             : {
     441             :   // Alice numbering scheme - Jun 03, 2006
     442             :   static Double_t glob[3];
     443             : 
     444           0 :   GetGlobal(absId, glob);
     445           0 :   vglob.SetXYZ(glob[0], glob[1], glob[2]);
     446           0 : }
     447             : 
     448             : //______________________________________________________________________
     449             : void AliEMCALGeometry::PrintCellIndexes(Int_t absId, int pri, const char *tit) const
     450             : {
     451             :   // Service methods
     452           0 :   Int_t nSupMod, nModule, nIphi, nIeta;
     453           0 :   Int_t iphi, ieta;
     454           0 :   TVector3 vg;
     455             : 
     456           0 :   GetCellIndex(absId,  nSupMod, nModule, nIphi, nIeta);
     457           0 :   printf(" %s | absId : %i -> nSupMod %i nModule %i nIphi %i nIeta %i \n", tit, absId,  nSupMod, nModule, nIphi, nIeta);
     458           0 :   if(pri>0) {
     459           0 :     GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi,ieta);
     460           0 :     printf(" local SM index : iphi %i : ieta %i \n", iphi,ieta);
     461           0 :     GetGlobal(absId, vg);
     462           0 :     printf(" vglob : mag %7.2f : perp %7.2f : z %7.2f : eta %6.4f : phi %6.4f(%6.2f) \n", 
     463           0 :            vg.Mag(), vg.Perp(), vg.Z(), vg.Eta(), vg.Phi(), vg.Phi()*TMath::RadToDeg());
     464             :   }
     465           0 : }
     466             : 
     467             : void AliEMCALGeometry::PrintLocalTrd1(Int_t pri) const
     468             : {
     469             :   // For comparing with numbers from drawing
     470           0 :   for(Int_t i=0; i<GetShishKebabTrd1Modules()->GetSize(); i++){
     471           0 :     printf(" %s | ", GetShishKebabModule(i)->GetName());
     472           0 :     if(i==0 && pri<1) GetShishKebabModule(i)->PrintShish(1);
     473           0 :     else     GetShishKebabModule(i)->PrintShish(pri);
     474             :   }
     475           0 : }
     476             : 
     477             : //________________________________________________________________________________________________
     478             : void AliEMCALGeometry::EtaPhiFromIndex(Int_t absId,Double_t &eta,Double_t &phi) const
     479             : {
     480             :   // Nov 16, 2006- float to double
     481             :   // version for TRD1 only
     482           0 :   static TVector3 vglob;
     483           0 :   GetGlobal(absId, vglob);
     484           0 :   eta = vglob.Eta();
     485           0 :   phi = vglob.Phi();
     486           0 : }
     487             : 
     488             : //________________________________________________________________________________________________
     489             : void AliEMCALGeometry::EtaPhiFromIndex(Int_t absId,Float_t &eta,Float_t &phi) const
     490             : {
     491             :   // Nov 16,2006 - should be discard in future
     492           0 :   static TVector3 vglob;
     493           0 :   GetGlobal(absId, vglob);
     494           0 :   eta = float(vglob.Eta());
     495           0 :   phi = float(vglob.Phi());
     496           0 : }
     497             : 
     498             : //
     499             : // == Shish-kebab cases ==
     500             : //
     501             : //________________________________________________________________________________________________
     502             : Int_t AliEMCALGeometry::GetAbsCellId(Int_t nSupMod, Int_t nModule, Int_t nIphi, Int_t nIeta) const
     503             : { 
     504             :   // 27-aug-04; 
     505             :   // corr. 21-sep-04; 
     506             :   //       13-oct-05; 110 degree case
     507             :   // May 31, 2006; ALICE numbering scheme:
     508             :   // 0 <= nSupMod < fNumberOfSuperModules
     509             :   // 0 <= nModule  < fNPHI * fNZ ( fNPHI * fNZ/2 for fKey110DEG=1)
     510             :   // 0 <= nIphi   < fNPHIdiv
     511             :   // 0 <= nIeta   < fNETAdiv
     512             :   // 0 <= absid   < fNCells
     513             :   Int_t id=0; // have to change from 0 to fNCells-1
     514      894622 :   for( int i = 0 ; i < nSupMod; i++) {
     515      698181 :     if(      GetSMType(i) == kEMCAL_Standard) id += fNCellsInSupMod;
     516        1747 :     else if( GetSMType(i) == kEMCAL_Half)     id += fNCellsInSupMod/2;
     517        2960 :     else if( GetSMType(i) == kEMCAL_3rd)      id += fNCellsInSupMod/3;
     518        1068 :     else if( GetSMType(i) == kDCAL_Standard)  id += 2*fNCellsInSupMod/3;
     519           0 :     else if( GetSMType(i) == kDCAL_Ext)       id += fNCellsInSupMod/3;
     520             :     else {
     521           0 :       AliError(Form("Uknown SuperModule Type !!"));
     522             :     }
     523             :   }
     524             :   
     525       64898 :   id += fNCellsInModule *nModule;
     526       64898 :   id += fNPHIdiv *nIphi;
     527       64898 :   id += nIeta;
     528       64898 :   if( !CheckAbsCellId(id) ) {
     529           0 :     id = -TMath::Abs(id); // if negative something wrong
     530           0 :   }
     531       64898 :   return id;
     532             : }
     533             : 
     534             : //________________________________________________________________________________________________
     535             : void  AliEMCALGeometry::GetModuleIndexesFromCellIndexesInSModule(Int_t nSupMod, Int_t iphi, Int_t ieta, 
     536             :                         Int_t &iphim, Int_t &ietam, Int_t &nModule) const
     537             : {
     538             :   // Transition from cell indexes (ieta,iphi) to module indexes (ietam,iphim, nModule)
     539             :   static Int_t nphi=-1;
     540      129796 :   nphi  = GetNumberOfModuleInPhiDirection(nSupMod);  
     541             : 
     542       64898 :   ietam  = ieta/fNETAdiv;
     543       64898 :   iphim  = iphi/fNPHIdiv;
     544       64898 :   nModule = ietam * nphi + iphim; 
     545       64898 : }
     546             : 
     547             : //________________________________________________________________________________________________
     548             : Int_t  AliEMCALGeometry::GetAbsCellIdFromCellIndexes(Int_t nSupMod, Int_t iphi, Int_t ieta) const
     549             : {
     550             :   // Transition from super module number(nSupMod) and cell indexes (ieta,iphi) to absId
     551             :   
     552             :   // Check if the indeces correspond to existing SM or tower indeces
     553      389406 :   if(iphi    < 0 || iphi    >= AliEMCALGeoParams::fgkEMCALRows || 
     554      129802 :      ieta    < 0 || ieta    >= AliEMCALGeoParams::fgkEMCALCols ||
     555      129802 :      nSupMod < 0 || nSupMod >= GetNumberOfSuperModules()         )
     556             :   {
     557           9 :     AliDebug(1,Form("Wrong cell indexes : SM %d, column (eta) %d, row (phi) %d", nSupMod,ieta,iphi));
     558           3 :     return -1 ;
     559             :   }
     560             :   
     561             :   static Int_t ietam=-1, iphim=-1, nModule=-1;
     562             :   static Int_t nIeta=-1, nIphi=-1; // cell indexes in module
     563             : 
     564       64898 :   GetModuleIndexesFromCellIndexesInSModule(nSupMod, iphi, ieta, ietam, iphim, nModule);
     565             : 
     566       64898 :   nIeta = ieta%fNETAdiv;
     567       64898 :   nIeta = fNETAdiv - 1 - nIeta;
     568       64898 :   nIphi = iphi%fNPHIdiv;
     569             :   
     570       64898 :   return GetAbsCellId(nSupMod, nModule, nIphi, nIeta);
     571       64901 : }
     572             : 
     573             : //________________________________________________________________________________________________
     574             : Bool_t AliEMCALGeometry::SuperModuleNumberFromEtaPhi(Double_t eta, Double_t phi, Int_t &nSupMod) const
     575             : { 
     576             :   // Return false if phi belongs a phi cracks between SM
     577             :  
     578             :   static Int_t i=0;
     579             : 
     580           0 :   if(TMath::Abs(eta) > fEtaMaxOfTRD1) return kFALSE;
     581             : 
     582           0 :   phi = TVector2::Phi_0_2pi(phi); // move phi to (0,2pi) boundaries
     583           0 :   Int_t nphism = fEMCGeometry->GetNumberOfSuperModules()/2;
     584           0 :   for(i=0; i<nphism; i++) {
     585           0 :     if(phi>=fPhiBoundariesOfSM[2*i] && phi<=fPhiBoundariesOfSM[2*i+1]) {
     586           0 :       nSupMod = 2*i;
     587           0 :       if(eta < 0.0) nSupMod++;
     588           0 :       if( GetSMType(nSupMod) == kDCAL_Standard) {// Gap between DCAL
     589           0 :         if(TMath::Abs(eta) < GetNEta()/3*(GetEMCGeometry()->GetTrd1Angle())*TMath::DegToRad()) return kFALSE;
     590             :       }
     591           0 :       AliDebug(1,Form("eta %f phi %f(%5.2f) : nSupMod %i : #bound %i", eta,phi,phi*TMath::RadToDeg(), nSupMod,i));
     592           0 :       return kTRUE;
     593             :     }
     594             :   }
     595           0 :   return kFALSE;
     596           0 : }
     597             : 
     598             : ///
     599             : /// Online mapping and numbering is the same for EMCal and DCal SMs but:
     600             : ///  - DCal odd SM (13,15,17) has online cols: 16-47; offline cols 0-31.
     601             : ///  - Even DCal SMs have the same numbering online and offline 0-31.
     602             : ///  - DCal 1/3 SM (18,19), online rows 16-23; offline rows 0-7
     603             : ///
     604             : /// Here shift the online cols or rows depending on the
     605             : /// super-module number to match the offline mapping.
     606             : ///
     607             : /// \param sm: super module number of the channel/cell
     608             : /// \param iphi: row/phi cell index, modified for DCal
     609             : /// \param ieta: column/eta index, modified for DCal
     610             : ///
     611             : //________________________________________________________________________________________________
     612             : void AliEMCALGeometry::ShiftOnlineToOfflineCellIndexes(Int_t sm, Int_t & iphi, Int_t & ieta) const 
     613             : {
     614         116 :   if ( sm == 13 || sm == 15 || sm == 17 )
     615             :   {
     616             :     // DCal odd SMs
     617           0 :     ieta -= 16; // Same cabling mapping as for EMCal, not considered offline.
     618           0 :   }
     619          58 :   else if ( sm == 18 || sm == 19 )
     620             :   {
     621             :     // DCal 1/3 SMs
     622           0 :     iphi -= 16; // Needed due to cabling mistake.
     623           0 :   }
     624          58 : }
     625             : 
     626             : ///
     627             : /// Here shift the DCal online cols or rows depending on the
     628             : /// super-module number to match the online mapping. 
     629             : ///
     630             : /// Reverse procedure to the one in the method above
     631             : /// ShiftOnlineToOfflineCellIndexes().
     632             : ///
     633             : /// \param sm: super module number of the channel/cell
     634             : /// \param iphi: row/phi cell index, modified for DCal
     635             : /// \param ieta: column/eta index, modified for DCal
     636             : ///
     637             : //________________________________________________________________________________________________
     638             : void AliEMCALGeometry::ShiftOfflineToOnlineCellIndexes(Int_t sm, Int_t & iphi, Int_t & ieta) const 
     639             : {
     640         116 :   if ( sm == 13 || sm == 15 || sm == 17 )
     641             :   {
     642             :     // DCal odd SMs
     643           0 :     ieta += 16; // Same cabling mapping as for EMCal, not considered offline.
     644           0 :   }
     645          58 :   else if ( sm == 18 || sm == 19 )
     646             :   {
     647             :     // DCal 1/3 SMs
     648           0 :     iphi += 16; // Needed due to cabling mistake.
     649           0 :   }
     650          58 : }
     651             : 
     652             : //________________________________________________________________________________________________
     653             : Bool_t AliEMCALGeometry::GetAbsCellIdFromEtaPhi(Double_t eta, Double_t phi, Int_t &absId) const
     654             : {
     655             : 
     656             :   // Nov 17,2006
     657             :   // stay here - phi problem as usual 
     658             :   static Int_t nSupMod=-1, i=0, ieta=-1, iphi=-1, etaShift=0, neta=-1, nphi=-1;
     659             :   static Double_t absEta=0.0, d=0.0, dmin=0.0, phiLoc=0;
     660           0 :   absId = nSupMod = - 1;
     661           0 :   if(SuperModuleNumberFromEtaPhi(eta, phi, nSupMod)) {
     662             :     // phi index first
     663           0 :     phi    = TVector2::Phi_0_2pi(phi);
     664           0 :     phiLoc = phi - fPhiCentersOfSMSec[nSupMod/2];
     665           0 :     nphi   = fPhiCentersOfCells.GetSize();
     666           0 :     if (     GetSMType(nSupMod) == kEMCAL_Half ) nphi  /= 2;
     667           0 :     else if( GetSMType(nSupMod) == kEMCAL_3rd )  nphi  /= 3;
     668           0 :     else if( GetSMType(nSupMod) == kDCAL_Ext )   nphi  /= 3;
     669             :     
     670           0 :     dmin   = TMath::Abs(fPhiCentersOfCells[0]-phiLoc);
     671           0 :     iphi   = 0;
     672           0 :     for(i=1; i<nphi; i++) {
     673           0 :       d = TMath::Abs(fPhiCentersOfCells[i] - phiLoc);
     674           0 :       if(d < dmin) {
     675           0 :         dmin = d;
     676           0 :         iphi = i;
     677           0 :       }
     678             :       //printf(" i %i : d %f : dmin %f : fPhiCentersOfCells[i] %f \n", i, d, dmin, fPhiCentersOfCells[i]);
     679             :     }
     680             :     // odd SM are turned with respect of even SM - reverse indexes
     681           0 :     AliDebug(2,Form(" iphi %i : dmin %f (phi %f, phiLoc %f ) ", iphi, dmin, phi, phiLoc));
     682             : 
     683             :     // eta index
     684           0 :     absEta   = TMath::Abs(eta);
     685           0 :     neta     = fCentersOfCellsEtaDir.GetSize();
     686           0 :     etaShift = iphi*neta;
     687           0 :     ieta     = 0;
     688           0 :     if( GetSMType(nSupMod) == kDCAL_Standard) ieta += 16; //jump 16 cells for DCSM
     689           0 :     dmin     = TMath::Abs(fEtaCentersOfCells[etaShift + ieta]-absEta);
     690           0 :     for(i= ieta+1 ; i<neta; i++) {
     691           0 :       d = TMath::Abs(fEtaCentersOfCells[i+etaShift] - absEta);
     692           0 :       if(d < dmin) {
     693           0 :         dmin = d;
     694           0 :         ieta = i;
     695           0 :       }
     696             :     }
     697           0 :     if( GetSMType(nSupMod) == kDCAL_Standard) ieta -= 16; //jump 16 cells for DCSM
     698             : 
     699           0 :     AliDebug(2,Form(" ieta %i : dmin %f (eta=%f) : nSupMod %i ", ieta, dmin, eta, nSupMod));
     700             :      
     701             :    //patch for mapping following alice convention  
     702           0 :     if(nSupMod%2 == 0) {// 47 + 16 -ieta for DCSM, 47 - ieta for others, revert the ordering on A side in order to keep convention.
     703           0 :       ieta = (neta -1)-ieta;
     704           0 :       if( GetSMType(nSupMod) == kDCAL_Standard) ieta -= 16; //recover cells for DCSM
     705             :     }
     706             : 
     707           0 :     absId = GetAbsCellIdFromCellIndexes(nSupMod, iphi, ieta);
     708           0 :     return kTRUE;
     709             :   }
     710           0 :   return kFALSE;
     711           0 : }
     712             : 
     713             : //________________________________________________________________________________________________
     714             : Bool_t  AliEMCALGeometry::CheckAbsCellId(Int_t absId) const
     715             : { 
     716             :   // May 31, 2006; only trd1 now
     717      420711 :   if(absId<0 || absId >= fNCells) return kFALSE;
     718      140233 :   else                            return kTRUE;
     719      140237 : }
     720             : 
     721             : //________________________________________________________________________________________________
     722             : Bool_t AliEMCALGeometry::GetCellIndex(Int_t absId,Int_t &nSupMod,Int_t &nModule,Int_t &nIphi,Int_t &nIeta) const
     723             : { 
     724             :   // 21-sep-04; 19-oct-05;
     725             :   // May 31, 2006; ALICE numbering scheme:
     726             :   // 
     727             :   // In:
     728             :   // absId   - cell is as in Geant,     0<= absId   < fNCells;
     729             :   // Out:
     730             :   // nSupMod - super module(SM) number, 0<= nSupMod < fNumberOfSuperModules;
     731             :   // nModule  - module number in SM,     0<= nModule  < fNCellsInSupMod/fNCellsInSupMod or(/2) for tow last SM (10th and 11th);
     732             :   // nIphi   - cell number in phi driection inside module; 0<= nIphi < fNPHIdiv; 
     733             :   // nIeta   - cell number in eta driection inside module; 0<= nIeta < fNETAdiv; 
     734             :   // 
     735      149364 :   if(!CheckAbsCellId(absId)) return kFALSE;
     736             : 
     737       74682 :   static Int_t tmp = absId;
     738             :   Int_t test = absId;
     739             :  
     740      811798 :   for(nSupMod = -1; test >= 0; ) {
     741      662446 :     nSupMod++;
     742      662446 :     tmp = test;
     743     1189360 :     if(      GetSMType(nSupMod) == kEMCAL_Standard) test -= fNCellsInSupMod;
     744      135532 :     else if( GetSMType(nSupMod) == kEMCAL_Half)     test -= fNCellsInSupMod/2;
     745      183389 :     else if( GetSMType(nSupMod) == kEMCAL_3rd)      test -= fNCellsInSupMod/3;
     746      170742 :     else if( GetSMType(nSupMod) == kDCAL_Standard)  test -= 2*fNCellsInSupMod/3;
     747        9216 :     else if( GetSMType(nSupMod) == kDCAL_Ext)       test -= fNCellsInSupMod/3;
     748             :     else {
     749           0 :       AliError(Form("Uknown SuperModule Type !!"));
     750           0 :       return kFALSE;
     751             :     }
     752             :   }
     753       74676 :   nModule = tmp / fNCellsInModule;
     754       74676 :   tmp     = tmp % fNCellsInModule;
     755       74676 :   nIphi   = tmp / fNPHIdiv;
     756       74676 :   nIeta   = tmp % fNPHIdiv;
     757             : 
     758       74676 :   return kTRUE;
     759       74680 : }
     760             : 
     761             : //________________________________________________________________________________________________
     762             : Int_t  AliEMCALGeometry::GetSuperModuleNumber(Int_t absId)  const
     763             : {
     764             :   // Return the number of the  supermodule given the absolute
     765             :   // ALICE numbering id
     766             : 
     767             :   static Int_t nSupMod=-1, nModule=-1, nIphi=-1, nIeta=-1;
     768         382 :   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
     769         191 :   return nSupMod;
     770             : } 
     771             : 
     772             : //________________________________________________________________________________________________
     773             : void AliEMCALGeometry::GetModulePhiEtaIndexInSModule(Int_t nSupMod, Int_t nModule,  int &iphim, int &ietam) const
     774             : { 
     775             :   // added nSupMod; - 19-oct-05 !
     776             :   // Alice numbering scheme        - Jun 01,2006 
     777             :   // ietam, iphi - indexes of module in two dimensional grid of SM
     778             :   // ietam - have to change from 0 to fNZ-1
     779             :   // iphim - have to change from 0 to nphi-1 (fNPhi-1 or fNPhi/2-1)
     780             :   static Int_t nphi=-1;
     781      279524 :   if(      GetSMType(nSupMod) == kEMCAL_Half )  nphi = fNPhi/2; // halfSM
     782      143653 :   else if( GetSMType(nSupMod) == kEMCAL_3rd  )  nphi = fNPhi/3; // 1/3 SM
     783      138943 :   else if( GetSMType(nSupMod) == kDCAL_Ext   )  nphi = fNPhi/3; // 1/3 SM
     784      132799 :   else                                          nphi = fNPhi;   // full SM
     785             :   
     786      139762 :   ietam = nModule/nphi;
     787      139762 :   iphim = nModule%nphi;
     788      139762 : }
     789             : 
     790             : //________________________________________________________________________________________________
     791             : void AliEMCALGeometry::GetCellPhiEtaIndexInSModule(Int_t nSupMod, Int_t nModule, Int_t nIphi, Int_t nIeta, 
     792             : int &iphi, int &ieta) const
     793             : { 
     794             :   // 
     795             :   // Added nSupMod; Nov 25, 05
     796             :   // Alice numbering scheme  - Jun 01,2006 
     797             :   // IN:
     798             :   // nSupMod - super module(SM) number, 0<= nSupMod < fNumberOfSuperModules;
     799             :   // nModule  - module number in SM,     0<= nModule  < fNCellsInSupMod/fNCellsInSupMod or(/2) for tow last SM (10th and 11th);
     800             :   // nIphi   - cell number in phi driection inside module; 0<= nIphi < fNPHIdiv; 
     801             :   // nIeta   - cell number in eta driection inside module; 0<= nIeta < fNETAdiv; 
     802             :   // 
     803             :   // OUT:
     804             :   // ieta, iphi - indexes of cell(tower) in two dimensional grid of SM
     805             :   // ieta - have to change from 0 to (fNZ*fNETAdiv-1)
     806             :   // iphi - have to change from 0 to (fNPhi*fNPHIdiv-1 or fNPhi*fNPHIdiv/2-1)
     807             :   //
     808             :   static Int_t iphim=-1, ietam=-1;
     809             : 
     810      279136 :   GetModulePhiEtaIndexInSModule(nSupMod,nModule, iphim, ietam); 
     811             :   //  ieta  = ietam*fNETAdiv + (1-nIeta); // x(module) = -z(SM) 
     812      139568 :   ieta  = ietam*fNETAdiv + (fNETAdiv - 1 - nIeta); // x(module) = -z(SM) 
     813      139568 :   iphi  = iphim*fNPHIdiv + nIphi;     // y(module) =  y(SM) 
     814             : 
     815      279136 :   if(iphi<0 || ieta<0)
     816           0 :   AliDebug(1,Form(" nSupMod %i nModule %i nIphi %i nIeta %i => ieta %i iphi %i\n", 
     817             :   nSupMod, nModule, nIphi, nIeta, ieta, iphi));
     818      139568 : }
     819             : 
     820             : // Methods for AliEMCALRecPoint - Feb 19, 2006
     821             : //________________________________________________________________________________________________
     822             : Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t &yr, Double_t &zr) const
     823             : {
     824             :   // Look to see what the relative
     825             :   // position inside a given cell is
     826             :   // for a recpoint.
     827             :   // Alice numbering scheme - Jun 08, 2006
     828             :   // In:
     829             :   // absId   - cell is as in Geant,     0<= absId   < fNCells;
     830             :   // OUT:
     831             :   // xr,yr,zr - x,y,z coordinates of cell with absId inside SM 
     832             : 
     833             :   // Shift index taking into account the difference between standard SM 
     834             :   // and SM of half (or one third) size in phi direction
     835             :  
     836           0 :   const Int_t kNphiIndex = fCentersOfCellsPhiDir.GetSize(); 
     837           0 :   Double_t  zshift = 0.5*GetDCALInnerEdge();
     838             :       
     839             :   static Int_t nSupMod=-1, nModule=-1, nIphi=-1, nIeta=-1, iphi=-1, ieta=-1;
     840           0 :   if(!CheckAbsCellId(absId)) return kFALSE;
     841             : 
     842           0 :   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
     843           0 :   GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta); 
     844             :         
     845             :   //Get eta position. Careful with ALICE conventions (increase index decrease eta)      
     846           0 :   Int_t ieta2 = ieta;
     847           0 :   if(nSupMod%2 == 0) {
     848           0 :     ieta2 = (fCentersOfCellsEtaDir.GetSize()-1)-ieta;// 47-ieta, revert the ordering on A side in order to keep convention.
     849           0 :   }
     850           0 :   if( GetSMType(nSupMod) == kDCAL_Standard && nSupMod%2 ) ieta2 += 16; // DCAL revert the ordering on C side ...
     851           0 :   zr = fCentersOfCellsEtaDir.At(ieta2); 
     852           0 :   if( GetSMType(nSupMod) == kDCAL_Standard ) zr -= zshift; // DCAL shift (SMALLER SM)
     853           0 :   xr = fCentersOfCellsXDir.At(ieta2);
     854             : 
     855             :   //Get phi position. Careful with ALICE conventions (increase index increase phi)
     856           0 :   Int_t iphi2 = iphi;
     857           0 :   if( GetSMType(nSupMod) == kDCAL_Ext ) {
     858           0 :   if(nSupMod%2 != 0)   iphi2 = (kNphiIndex/3 -1)-iphi;  // 7-iphi [1/3SM], revert the ordering on C side in order to keep convention.
     859           0 :     yr = fCentersOfCellsPhiDir.At(iphi2 + kNphiIndex/3);
     860           0 :   } else if( GetSMType(nSupMod) == kEMCAL_Half ){
     861           0 :     if(nSupMod%2 != 0)   iphi2 = (kNphiIndex/2 -1)-iphi;  //11-iphi [1/2SM], revert the ordering on C side in order to keep convention.
     862           0 :     yr = fCentersOfCellsPhiDir.At(iphi2 + kNphiIndex/4);
     863           0 :   } else if( GetSMType(nSupMod) == kEMCAL_3rd ){
     864           0 :     if(nSupMod%2 != 0)   iphi2 = (kNphiIndex/3 -1)-iphi;  // 7-iphi [1/3SM], revert the ordering on C side in order to keep convention.
     865           0 :     yr = fCentersOfCellsPhiDir.At(iphi2 + kNphiIndex/3);
     866           0 :   } else {
     867           0 :     if(nSupMod%2 != 0)   iphi2 = (kNphiIndex   -1)-iphi;// 23-iphi, revert the ordering on C side in order to keep conventi
     868           0 :     yr = fCentersOfCellsPhiDir.At(iphi2);
     869             :   } 
     870           0 :   AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
     871             : 
     872             :   return kTRUE;
     873           0 : }
     874             : 
     875             : //________________________________________________________________________________________________
     876             : Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t loc[3]) const
     877             : {
     878             :   // Look to see what the relative
     879             :   // position inside a given cell is
     880             :   // for a recpoint.    // Alice numbering scheme - Jun 03, 2006
     881           0 :   loc[0] = loc[1] = loc[2]=0.0;
     882           0 :   if(RelPosCellInSModule(absId, loc[0],loc[1],loc[2])) {
     883           0 :     return kTRUE;
     884             :   }
     885           0 :   return kFALSE;
     886           0 : }
     887             : 
     888             : //________________________________________________________________________________________________
     889             : Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, TVector3 &vloc) const
     890             : {
     891             :   // Look to see what the relative
     892             :   // position inside a given cell is
     893             :   // for a recpoint.  
     894             :   // Alice numbering scheme - Jun 03, 2006
     895             :   static Double_t loc[3];
     896           0 :   if(RelPosCellInSModule(absId,loc)) {
     897           0 :     vloc.SetXYZ(loc[0], loc[1], loc[2]);
     898           0 :     return kTRUE;
     899             :   } else {
     900           0 :     vloc.SetXYZ(0,0,0);
     901           0 :     return kFALSE;
     902             :   }
     903           0 : }
     904             : 
     905             : //________________________________________________________________________________________________
     906             : Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t distEff, Double_t &xr, Double_t &yr, Double_t &zr) const
     907             : {
     908             :   // Jul 30, 2007 - taking into account position of shower max
     909             :   // Look to see what the relative
     910             :   // position inside a given cell is
     911             :   // for a recpoint.
     912             :   // In:
     913             :   // absId   - cell is as in Geant,     0<= absId   < fNCells;
     914             :   // e       - cluster energy
     915             :   // OUT:
     916             :   // xr,yr,zr - x,y,z coordinates of cell with absId inside SM 
     917             :   
     918             :   // Shift index taking into account the difference between standard SM 
     919             :   // and SM of half (or one third) size in phi direction
     920             :   
     921         388 :   const Int_t kNphiIndex = fCentersOfCellsPhiDir.GetSize();
     922         194 :   Double_t  zshift = 0.5*GetDCALInnerEdge();
     923             :   Int_t kDCalshift = 8;//wangml DCal cut first 8 modules(16 cells)
     924             :    
     925             :   static Int_t nSupMod=0, nModule=-1, nIphi=-1, nIeta=-1, iphi=-1, ieta=-1;
     926             :   static Int_t iphim=-1, ietam=-1;
     927             :   static AliEMCALShishKebabTrd1Module *mod = 0;
     928         200 :   static TVector2 v;
     929         194 :   if(!CheckAbsCellId(absId)) return kFALSE;
     930             :   
     931         194 :   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
     932         194 :   GetModulePhiEtaIndexInSModule(nSupMod, nModule, iphim, ietam);
     933         194 :   GetCellPhiEtaIndexInSModule(nSupMod,nModule,nIphi,nIeta, iphi, ieta); 
     934             :   
     935             :   //Get eta position. Careful with ALICE conventions (increase index decrease eta)      
     936         194 :   if(nSupMod%2 == 0) {             
     937          32 :     ietam = (fCentersOfCellsEtaDir.GetSize()/2-1)-ietam;// 24-ietam, revert the ordering on A side in order to keep convention.
     938          64 :     if(nIeta == 0) nIeta = 1;
     939           0 :     else           nIeta = 0;
     940             :   }
     941         194 :   if( GetSMType(nSupMod) == kDCAL_Standard && nSupMod%2) ietam += kDCalshift; // DCAL revert the ordering on C side ....
     942         194 :   mod = GetShishKebabModule(ietam);
     943         194 :   mod ->GetPositionAtCenterCellLine(nIeta, distEff, v); 
     944         194 :   xr = v.Y() - fParSM[0];
     945         194 :   zr = v.X() - fParSM[2];
     946         194 :   if( GetSMType(nSupMod) == kDCAL_Standard ) zr -= zshift; // DCAL shift (SMALLER SM)
     947             :  
     948             :   //Get phi position. Careful with ALICE conventions (increase index increase phi)
     949         194 :   Int_t iphi2 = iphi;
     950         388 :   if( GetSMType(nSupMod) == kDCAL_Ext ) {
     951         194 :      if(nSupMod%2 != 0)  iphi2 = (kNphiIndex/3 -1)-iphi;  // 7-iphi [1/3SM], revert the ordering on C side in order to keep convention.
     952           0 :      yr = fCentersOfCellsPhiDir.At(iphi2 + kNphiIndex/3);
     953         388 :    } else if( GetSMType(nSupMod) == kEMCAL_Half ){
     954         194 :      if(nSupMod%2 != 0)  iphi2 = (kNphiIndex/2 -1)-iphi;  //11-iphi [1/2SM], revert the ordering on C side in order to keep convention.
     955           0 :      yr = fCentersOfCellsPhiDir.At(iphi2 + kNphiIndex/2);
     956         388 :    } else if( GetSMType(nSupMod) == kEMCAL_3rd ){
     957         202 :      if(nSupMod%2 != 0)  iphi2 = (kNphiIndex/3 -1)-iphi;  // 7-iphi [1/3SM], revert the ordering on C side in order to keep convention.
     958           4 :      yr = fCentersOfCellsPhiDir.At(iphi2 + kNphiIndex/3);
     959           4 :    } else {
     960         348 :      if(nSupMod%2 != 0)  iphi2 = (kNphiIndex   -1)-iphi;// 23-iphi, revert the ordering on C side in order to keep convention.
     961         190 :      yr = fCentersOfCellsPhiDir.At(iphi2);
     962             :    }
     963             :   
     964         582 :   AliDebug(1,Form("absId %i nSupMod %i iphi %i ieta %i xr %f yr %f zr %f ",absId,nSupMod,iphi,ieta,xr,yr,zr));
     965             :   
     966             :   return kTRUE;
     967         194 : }
     968             : 
     969             : //________________________________________________________________________________________________
     970             : void AliEMCALGeometry::CreateListOfTrd1Modules()
     971             : {
     972             :   // Generate the list of Trd1 modules
     973             :   // which will make up the EMCAL
     974             :   // geometry
     975             :   // key: look to the AliEMCALShishKebabTrd1Module::
     976             : 
     977          20 :   AliDebug(2,Form(" AliEMCALGeometry::CreateListOfTrd1Modules() started "));
     978             : 
     979             :   AliEMCALShishKebabTrd1Module *mod=0, *mTmp=0; // current module
     980           5 :   if(fShishKebabTrd1Modules == 0) {
     981          10 :     fShishKebabTrd1Modules = new TList;
     982           5 :     fShishKebabTrd1Modules->SetName("ListOfTRD1");
     983         250 :     for(int iz=0; iz< fEMCGeometry->GetNZ(); iz++) {
     984         240 :       if(iz==0) {
     985             :         //        mod  = new AliEMCALShishKebabTrd1Module(TMath::Pi()/2.,this);
     986         125 :         mod  = new AliEMCALShishKebabTrd1Module(TMath::Pi()/2.,fEMCGeometry);
     987           5 :       } else {
     988         115 :         mTmp  = new AliEMCALShishKebabTrd1Module(*mod);
     989             :         mod   = mTmp;
     990             :       }
     991         120 :       fShishKebabTrd1Modules->Add(mod);
     992             :     }
     993           5 :   } else {
     994           0 :     AliDebug(2,Form(" Already exits : "));
     995             :   }
     996           5 :   mod = (AliEMCALShishKebabTrd1Module*)fShishKebabTrd1Modules->At(fShishKebabTrd1Modules->GetSize()-1);
     997           5 :   fEtaMaxOfTRD1 = mod->GetMaxEtaOfModule(0);
     998             : 
     999          15 :   AliDebug(2,Form(" fShishKebabTrd1Modules has %i modules : max eta %5.4f \n",
    1000             :                   fShishKebabTrd1Modules->GetSize(),fEtaMaxOfTRD1));
    1001             :   // Feb 20,2006;
    1002             :   // Jun 01, 2006 - ALICE numbering scheme
    1003             :   // define grid for cells in eta(z) and x directions in local coordinates system of SM
    1004             :   // Works just for 2x2 case only -- ?? start here
    1005             :   //
    1006             :   //
    1007             :   // Define grid for cells in phi(y) direction in local coordinates system of SM
    1008             :   // as for 2X2 as for 3X3 - Nov 8,2006
    1009             :   //
    1010          15 :   AliDebug(2,Form(" Cells grid in phi directions : size %i\n", fCentersOfCellsPhiDir.GetSize()));
    1011             :   Int_t ind=0; // this is phi index
    1012           5 :   Int_t ieta=0, nModule=0, iphiTemp;
    1013           5 :   Double_t xr=0., zr=0., theta=0., phi=0., eta=0., r=0., x=0.,y=0.;
    1014           5 :   TVector3 vglob;
    1015             :   Double_t ytCenterModule=0.0, ytCenterCell=0.0;
    1016             : 
    1017           5 :   fCentersOfCellsPhiDir.Set(fNPhi*fNPHIdiv);
    1018           5 :   fPhiCentersOfCells.Set(fNPhi*fNPHIdiv);
    1019             : 
    1020           5 :   Double_t r0 = fIPDistance + fLongModuleSize/2.;
    1021         130 :   for(Int_t it=0; it<fNPhi; it++) { // cycle on modules
    1022          60 :     ytCenterModule = -fParSM[1] + fPhiModuleSize*(2*it+1)/2;  // center of module
    1023         360 :     for(Int_t ic=0; ic<fNPHIdiv; ic++) { // cycle on cells in module
    1024         120 :       if(fNPHIdiv==2) {
    1025         120 :         ytCenterCell = ytCenterModule + fPhiTileSize *(2*ic-1)/2.;
    1026         120 :       } else if(fNPHIdiv==3){
    1027           0 :         ytCenterCell = ytCenterModule + fPhiTileSize *(ic-1);
    1028           0 :       } else if(fNPHIdiv==1){
    1029             :         ytCenterCell = ytCenterModule;
    1030           0 :       }
    1031         120 :       fCentersOfCellsPhiDir.AddAt(ytCenterCell,ind);
    1032             :       // Define grid on phi direction
    1033             :       // Grid is not the same for different eta bin;
    1034             :       // Effect is small but is still here
    1035         120 :       phi = TMath::ATan2(ytCenterCell, r0);
    1036         120 :       fPhiCentersOfCells.AddAt(phi, ind);
    1037             : 
    1038         600 :       AliDebug(2,Form(" ind %2.2i : y %8.3f ", ind, fCentersOfCellsPhiDir.At(ind)));
    1039         120 :       ind++;
    1040             :     }
    1041             :   }
    1042             : 
    1043           5 :   fCentersOfCellsEtaDir.Set(fNZ *fNETAdiv);
    1044           5 :   fCentersOfCellsXDir.Set(fNZ *fNETAdiv);
    1045           5 :   fEtaCentersOfCells.Set(fNZ *fNETAdiv * fNPhi*fNPHIdiv);
    1046          25 :   AliDebug(2,Form(" Cells grid in eta directions : size %i\n", fCentersOfCellsEtaDir.GetSize()));
    1047         250 :   for(Int_t it=0; it<fNZ; it++) {
    1048         120 :     AliEMCALShishKebabTrd1Module *trd1 = GetShishKebabModule(it);
    1049         120 :     nModule = fNPhi*it;
    1050         720 :     for(Int_t ic=0; ic<fNETAdiv; ic++) {
    1051         240 :       if(fNPHIdiv==2) {
    1052         240 :         trd1->GetCenterOfCellInLocalCoordinateofSM(ic, xr, zr);      // case of 2X2
    1053         240 :         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
    1054         240 :       } if(fNPHIdiv==3) {
    1055           0 :         trd1->GetCenterOfCellInLocalCoordinateofSM3X3(ic, xr, zr);  // case of 3X3
    1056           0 :         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
    1057         240 :       } if(fNPHIdiv==1) {
    1058           0 :         trd1->GetCenterOfCellInLocalCoordinateofSM1X1(xr, zr);      // case of 1X1
    1059           0 :         GetCellPhiEtaIndexInSModule(0, nModule, 0, ic, iphiTemp, ieta);
    1060             :       }
    1061         240 :       fCentersOfCellsXDir.AddAt(float(xr) - fParSM[0],ieta);
    1062         240 :       fCentersOfCellsEtaDir.AddAt(float(zr) - fParSM[2],ieta);
    1063             :       // Define grid on eta direction for each bin in phi
    1064       12000 :       for(int iphi=0; iphi<fCentersOfCellsPhiDir.GetSize(); iphi++) {
    1065        5760 :         x = xr + trd1->GetRadius();
    1066       11520 :         y = fCentersOfCellsPhiDir[iphi];
    1067        5760 :         r = TMath::Sqrt(x*x + y*y + zr*zr);
    1068        5760 :         theta = TMath::ACos(zr/r);
    1069        5760 :         eta   = AliEMCALShishKebabTrd1Module::ThetaToEta(theta);
    1070             :         //        ind   = ieta*fCentersOfCellsPhiDir.GetSize() + iphi;
    1071        5760 :         ind   = iphi*fCentersOfCellsEtaDir.GetSize() + ieta;
    1072        5760 :         fEtaCentersOfCells.AddAt(eta, ind);
    1073             :       }
    1074             :       //printf(" ieta %i : xr + trd1->GetRadius() %f : zr %f : eta %f \n", ieta, xr + trd1->GetRadius(), zr, eta);
    1075             :     }
    1076             :   }
    1077         490 :   for(Int_t i=0; i<fCentersOfCellsEtaDir.GetSize(); i++) {
    1078        1200 :     AliDebug(2,Form(" ind %2.2i : z %8.3f : x %8.3f", i+1,
    1079             :                     fCentersOfCellsEtaDir.At(i),fCentersOfCellsXDir.At(i)));
    1080             :   }
    1081             : 
    1082           5 : }
    1083             : 
    1084             : //________________________________________________________________________________________________
    1085             : AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta) const
    1086             : {
    1087             :   //This method was too long to be
    1088             :   //included in the header file - the
    1089             :   //rule checker complained about it's
    1090             :   //length, so we move it here.  It returns the
    1091             :   //shishkebabmodule at a given eta index point.
    1092             : 
    1093             :   static AliEMCALShishKebabTrd1Module* trd1=0;
    1094         942 :   if(fShishKebabTrd1Modules && neta>=0 && neta<fShishKebabTrd1Modules->GetSize()) {
    1095         314 :     trd1 = (AliEMCALShishKebabTrd1Module*)fShishKebabTrd1Modules->At(neta);
    1096         314 :   } else trd1 = 0;
    1097         314 :   return trd1;
    1098             : }
    1099             : 
    1100             : //___________________________________________________________________
    1101             : void AliEMCALGeometry::PrintGeometryGeoUtils()
    1102             : {
    1103             :   //Print information from geometry
    1104           0 :   fEMCGeometry->PrintGeometry();
    1105             : 
    1106           0 :   printf(" fShishKebabTrd1Modules has %i modules : max eta %5.4f \n", 
    1107           0 :          fShishKebabTrd1Modules->GetSize(),fEtaMaxOfTRD1);
    1108             :   
    1109           0 :   printf("\n Cells grid in eta directions : size %i\n", fCentersOfCellsEtaDir.GetSize());
    1110           0 :   for(Int_t i=0; i<fCentersOfCellsEtaDir.GetSize(); i++) {
    1111           0 :     printf(" ind %2.2i : z %8.3f : x %8.3f \n", i, 
    1112           0 :            fCentersOfCellsEtaDir.At(i),fCentersOfCellsXDir.At(i));
    1113             :     int ind=0; // Nov 21,2006
    1114           0 :     for(Int_t iphi=0; iphi<fCentersOfCellsPhiDir.GetSize(); iphi++) {
    1115           0 :       ind = iphi*fCentersOfCellsEtaDir.GetSize() + i;
    1116           0 :       printf("%6.4f ", fEtaCentersOfCells[ind]);
    1117           0 :       if((iphi+1)%12 == 0) printf("\n");
    1118             :     }
    1119           0 :     printf("\n");
    1120             :     
    1121             :   }
    1122             : 
    1123           0 :   printf("\n Cells grid in phi directions : size %i\n", fCentersOfCellsPhiDir.GetSize());
    1124           0 :   for(Int_t i=0; i<fCentersOfCellsPhiDir.GetSize(); i++) {
    1125           0 :     double phi=fPhiCentersOfCells.At(i);
    1126           0 :     printf(" ind %2.2i : y %8.3f : phi %7.5f(%6.2f) \n", i, fCentersOfCellsPhiDir.At(i), 
    1127           0 :            phi, phi*TMath::RadToDeg());
    1128             :   }
    1129           0 : }
    1130             : 
    1131             : //____________________________________________________________________________
    1132             : Bool_t  AliEMCALGeometry::Impact(const TParticle * particle) const 
    1133             : {
    1134             :   // Tells if a particle enters EMCAL
    1135             :   Bool_t in=kFALSE;
    1136           0 :   Int_t absID=0;
    1137           0 :   TVector3 vtx(particle->Vx(),particle->Vy(),particle->Vz());
    1138           0 :   TVector3 vimpact(0,0,0);
    1139           0 :   ImpactOnEmcal(vtx,particle->Theta(),particle->Phi(),absID,vimpact);
    1140           0 :   if(absID>=0) 
    1141           0 :     in=kTRUE;
    1142           0 :   return in;
    1143           0 : }
    1144             : //____________________________________________________________________________
    1145             : void AliEMCALGeometry::ImpactOnEmcal(TVector3 vtx, Double_t theta, Double_t phi, 
    1146             :                                      Int_t & absId, TVector3 & vimpact) const
    1147             : {
    1148             :   // calculates the impact coordinates on EMCAL (centre of a tower/not on EMCAL surface) 
    1149             :   // of a neutral particle  
    1150             :   // emitted in the vertex vtx[3] with direction theta and phi in the ALICE global coordinate system
    1151             : 
    1152           0 :   TVector3 p(TMath::Sin(theta)*TMath::Cos(phi),TMath::Sin(theta)*TMath::Sin(phi),TMath::Cos(theta)) ;
    1153             : 
    1154           0 :   vimpact.SetXYZ(0,0,0);
    1155           0 :   absId=-1;
    1156           0 :   if(phi==0 || theta==0) return;
    1157             : 
    1158           0 :   TVector3 direction;
    1159           0 :   Double_t factor = (fIPDistance-vtx[1])/p[1];
    1160           0 :   direction = vtx + factor*p;
    1161             : 
    1162             :   //from particle direction -> tower hitted
    1163           0 :   GetAbsCellIdFromEtaPhi(direction.Eta(),direction.Phi(),absId);
    1164             :   
    1165             :   //tower absID hitted -> tower/module plane (evaluated at the center of the tower)
    1166           0 :   Int_t nSupMod=-1, nModule=-1, nIphi=-1, nIeta=-1;
    1167           0 :   Double_t loc[3],loc2[3],loc3[3];
    1168           0 :   Double_t glob[3]={},glob2[3]={},glob3[3]={};
    1169             :   
    1170           0 :   if(!RelPosCellInSModule(absId,loc)) return;
    1171             :   
    1172             :   //loc is cell center of tower
    1173           0 :   GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta);
    1174             : 
    1175             :   //look at 2 neighbours-s cell using nIphi={0,1} and nIeta={0,1}
    1176             :   Int_t nIphi2=-1,nIeta2=-1,absId2=-1,absId3=-1;
    1177           0 :   if(nIeta==0) nIeta2=1;
    1178             :   else nIeta2=0;
    1179           0 :   absId2=GetAbsCellId(nSupMod,nModule,nIphi,nIeta2);  
    1180           0 :   if(nIphi==0) nIphi2=1;
    1181             :   else nIphi2=0;
    1182           0 :   absId3=GetAbsCellId(nSupMod,nModule,nIphi2,nIeta);
    1183             : 
    1184             :   //2nd point on emcal cell plane
    1185           0 :   if(!RelPosCellInSModule(absId2,loc2)) return;
    1186             :     
    1187             :   //3rd point on emcal cell plane
    1188           0 :   if(!RelPosCellInSModule(absId3,loc3)) return;
    1189             :     
    1190             :   // Get Matrix
    1191           0 :   const TGeoHMatrix* m = GetMatrixForSuperModule(nSupMod);
    1192           0 :   if(m) {
    1193           0 :     m->LocalToMaster(loc, glob);
    1194           0 :     m->LocalToMaster(loc2, glob2);
    1195           0 :     m->LocalToMaster(loc3, glob3);
    1196             :   } else {
    1197           0 :     AliFatal("Geo matrixes are not loaded \n") ;
    1198             :   }
    1199             : 
    1200             :   //Equation of Plane from glob,glob2,glob3 (Ax+By+Cz+D=0)
    1201           0 :   Double_t a = glob[1]*(glob2[2]-glob3[2]) + glob2[1]*(glob3[2]-glob[2]) + glob3[1]*(glob[2]-glob2[2]);
    1202           0 :   Double_t b = glob[2]*(glob2[0]-glob3[0]) + glob2[2]*(glob3[0]-glob[0]) + glob3[2]*(glob[0]-glob2[0]);
    1203           0 :   Double_t c = glob[0]*(glob2[1]-glob3[1]) + glob2[0]*(glob3[1]-glob[1]) + glob3[0]*(glob[1]-glob2[1]);
    1204           0 :   Double_t d = glob[0]*(glob2[1]*glob3[2]-glob3[1]*glob2[2]) + glob2[0]*(glob3[1]*glob[2]-glob[1]*glob3[2]) + glob3[0]*(glob[1]*glob2[2]-glob2[1]*glob[2]);
    1205             :   d=-d;
    1206             :   
    1207             :   //shift equation of plane from tower/module center to surface along vector (A,B,C) normal to tower/module plane
    1208           0 :   Double_t dist = fLongModuleSize/2.;
    1209           0 :   Double_t norm = TMath::Sqrt(a*a+b*b+c*c);
    1210             :   Double_t glob4[3]={};
    1211           0 :   TVector3 dir(a,b,c);
    1212           0 :   TVector3 point(glob[0],glob[1],glob[2]); 
    1213           0 :   if(point.Dot(dir)<0) dist*=-1;
    1214           0 :   glob4[0]=glob[0]-dist*a/norm;
    1215           0 :   glob4[1]=glob[1]-dist*b/norm;
    1216           0 :   glob4[2]=glob[2]-dist*c/norm;
    1217           0 :   d = glob4[0]*a +  glob4[1]*b +  glob4[2]*c ;
    1218           0 :   d = -d;
    1219             : 
    1220             :   //Line determination (2 points for equation of line : vtx and direction)
    1221             :   //impact between line (particle) and plane (module/tower plane)
    1222           0 :   Double_t den = a*(vtx(0)-direction(0)) + b*(vtx(1)-direction(1)) + c*(vtx(2)-direction(2));
    1223           0 :   if(den==0){
    1224           0 :     printf("ImpactOnEmcal() No solution :\n");
    1225           0 :     return;
    1226             :   }
    1227             :   
    1228           0 :   Double_t length = a*vtx(0)+b*vtx(1)+c*vtx(2)+d;
    1229           0 :   length /=den;
    1230             :   
    1231           0 :   vimpact.SetXYZ(vtx(0)+length*(direction(0)-vtx(0)),vtx(1)+length*(direction(1)-vtx(1)),vtx(2)+length*(direction(2)-vtx(2)));
    1232             :   
    1233             :   //shift vimpact from tower/module surface to center along vector (A,B,C) normal to tower/module plane
    1234           0 :   vimpact.SetXYZ(vimpact(0)+dist*a/norm,vimpact(1)+dist*b/norm,vimpact(2)+dist*c/norm);
    1235             :   
    1236             :   return;
    1237           0 : }
    1238             : 
    1239             : //_____________________________________________________________________________
    1240             : Bool_t AliEMCALGeometry::IsInEMCAL(Double_t x, Double_t y, Double_t z) const 
    1241             : {
    1242             :   // Checks whether point is inside the EMCal volume 
    1243           0 :   if( IsInEMCALOrDCAL(x,y,z) == 1 ) return kTRUE;
    1244           0 :   else return kFALSE;
    1245           0 : }
    1246             : 
    1247             : //_____________________________________________________________________________
    1248             : Bool_t AliEMCALGeometry::IsInDCAL(Double_t x, Double_t y, Double_t z) const 
    1249             : {
    1250             :   // Checks whether point is inside the DCal volume
    1251           0 :   if( IsInEMCALOrDCAL(x,y,z) == 2 ) return kTRUE;
    1252           0 :   else return kFALSE;
    1253           0 : }
    1254             : 
    1255             : //_____________________________________________________________________________
    1256             : Int_t AliEMCALGeometry::IsInEMCALOrDCAL(Double_t x, Double_t y, Double_t z) const 
    1257             : {
    1258             :   // Checks whether point is inside the EMCal volume (included DCal), used in AliEMCALv*.cxx
    1259             :   //
    1260             :   // Code uses cylindrical approximation made of inner radius (for speed)
    1261             :   //
    1262             :   // Points behind EMCAl/DCal, i.e. R > outer radius, but eta, phi in acceptance 
    1263             :   // are considered to inside
    1264             : 
    1265         290 :   Double_t r=sqrt(x*x+y*y);
    1266             : 
    1267         221 :   if ( r <= fEnvelop[0] ) return 0;
    1268             :   else {
    1269          69 :     Double_t theta = TMath::ATan2(r,z);
    1270             :     Double_t eta;
    1271          69 :     if(theta == 0)  eta = 9999;
    1272          69 :     else            eta = -TMath::Log(TMath::Tan(theta/2.));
    1273         138 :     if (eta < fArm1EtaMin || eta > fArm1EtaMax) return 0;
    1274             : 
    1275          69 :     Double_t phi = TMath::ATan2(y,x) * 180./TMath::Pi();
    1276          78 :     if (phi < 0) phi += 360;  // phi should go from 0 to 360 in this case
    1277             : 
    1278         203 :     if (      phi >= fArm1PhiMin         && phi <= fEMCALPhiMax ) return 1;
    1279          16 :     else if ( phi >= fDCALPhiMin         && phi <= fDCALStandardPhiMax && TMath::Abs(eta) > fDCALInnerExtandedEta ) return 2;
    1280           0 :     else if ( phi > fDCALStandardPhiMax  && phi <= fDCALPhiMax  ) return 2;
    1281           0 :     else return 0;
    1282             :   } 
    1283         145 : }
    1284             : 
    1285             : ///
    1286             : /// Provides shift-rotation matrix for EMCAL from externally set matrix or 
    1287             : /// from TGeoManager
    1288             : /// \return alignment matrix for a super module number
    1289             : /// \param smod: super module number
    1290             : ///
    1291             : //____________________________________________________________________________
    1292             : const TGeoHMatrix * AliEMCALGeometry::GetMatrixForSuperModule(Int_t smod) const 
    1293             : {       
    1294         291 :   if(smod < 0 || smod > fEMCGeometry->GetNumberOfSuperModules()) 
    1295           0 :     AliFatal(Form("Wrong supermodule index -> %d",smod));
    1296             :                 
    1297             :   // Use matrices set externally
    1298         291 :   if(!gGeoManager || (gGeoManager && fUseExternalMatrices))
    1299             :   {
    1300           0 :     if(fkSModuleMatrix[smod])
    1301             :     {
    1302           0 :       return fkSModuleMatrix[smod] ;
    1303             :     }
    1304             :     else
    1305             :     {
    1306           0 :       AliInfo("Stop:");
    1307           0 :       printf("\t Can not find EMCAL misalignment matrixes\n") ;
    1308           0 :       printf("\t Either import TGeoManager from geometry.root or \n");
    1309           0 :       printf("\t read stored matrixes from AliESD Header:  \n") ;   
    1310           0 :       printf("\t AliEMCALGeometry::SetMisalMatrixes(header->GetEMCALMisalMatrix()) \n") ;
    1311           0 :       AliFatal("") ;
    1312             :     }  
    1313           0 :   }//external matrices
    1314             :   
    1315             :   // If gGeoManager exists, take matrix from it
    1316         194 :   if(gGeoManager) return GetMatrixForSuperModuleFromGeoManager(smod);
    1317             :   
    1318           0 :   return 0 ;
    1319          97 : }
    1320             : 
    1321             : 
    1322             : ///
    1323             : /// Provides shift-rotation matrix for EMCAL from fkSModuleMatrix[smod]
    1324             : /// Unsafe method, not to be used in reconstruction, just check there is 
    1325             : /// something in the array of matrices without crashing, for EVE checks.
    1326             : ///
    1327             : /// \return alignment matrix for a super module number
    1328             : /// \param smod: super module number
    1329             : ///
    1330             : //______________________________________________________________________________________
    1331             : const TGeoHMatrix * AliEMCALGeometry::GetMatrixForSuperModuleFromArray(Int_t smod) const 
    1332             : {       
    1333           0 :   if(smod < 0 || smod > fEMCGeometry->GetNumberOfSuperModules()) 
    1334           0 :     AliFatal(Form("Wrong supermodule index -> %d",smod));
    1335             :   
    1336           0 :   return fkSModuleMatrix[smod] ;
    1337             : }
    1338             : 
    1339             : 
    1340             : ///
    1341             : /// Provides shift-rotation matrix for EMCAL from the TGeoManager.
    1342             : /// \return alignment matrix for a super module number
    1343             : /// \param smod: super module number
    1344             : ///
    1345             : //____________________________________________________________________________
    1346             : const TGeoHMatrix * AliEMCALGeometry::GetMatrixForSuperModuleFromGeoManager(Int_t smod) const 
    1347             : {  
    1348             :   const Int_t buffersize = 255;
    1349         194 :   char  path[buffersize] ;
    1350             :   Int_t tmpType = -1;
    1351             :   Int_t smOrder = 0;
    1352             :   
    1353             :   //Get the order for SM
    1354        1378 :   for( Int_t i = 0; i < smod+1; i++)
    1355             :   {
    1356         592 :     if(GetSMType(i) == tmpType) 
    1357             :     {
    1358         493 :       smOrder++;
    1359         493 :     } 
    1360             :     else 
    1361             :     {
    1362          99 :       tmpType = GetSMType(i);
    1363             :       smOrder = 1;
    1364             :     }
    1365             :   } 
    1366             :   
    1367          97 :   Int_t   smType = GetSMType(smod);
    1368          97 :   TString smName = "";
    1369             :   
    1370         192 :   if      ( smType == kEMCAL_Standard ) smName = "SMOD";
    1371           2 :   else if ( smType == kEMCAL_Half )     smName = "SM10";
    1372           4 :   else if ( smType == kEMCAL_3rd )      smName = "SM3rd";
    1373           0 :   else if ( smType == kDCAL_Standard )  smName = "DCSM";
    1374           0 :   else if ( smType == kDCAL_Ext )       smName = "DCEXT";
    1375           0 :   else AliError("Unkown SM Type!!");
    1376             :   
    1377         194 :   snprintf(path,buffersize,"/ALIC_1/XEN1_1/%s_%d", smName.Data(), smOrder) ;
    1378             :   
    1379         194 :   if (!gGeoManager->cd(path))
    1380           0 :     AliFatal(Form("Geo manager can not find path %s!\n",path));
    1381             :   
    1382          97 :   return gGeoManager->GetCurrentMatrix();
    1383          97 : }
    1384             : 
    1385             : 
    1386             : //__________________________________________________________________________________________________________________
    1387             : void AliEMCALGeometry::RecalculateTowerPosition(Float_t drow, Float_t dcol, const Int_t sm, const Float_t depth,
    1388             :                                                 const Float_t misaligTransShifts[15], const Float_t misaligRotShifts[15], Float_t global[3]) const
    1389             : { 
    1390             :   //Transform clusters cell position into global with alternative method, taking into account the depth calculation.
    1391             :   //Input are: the tower indeces, 
    1392             :   //           supermodule, 
    1393             :   //           particle type (photon 0, electron 1, hadron 2 )
    1394             :   //           misalignment shifts to global position in case of need.
    1395             :   // Federico.Ronchetti@cern.ch
    1396             :     
    1397             :   // To use in a print later
    1398             :   Float_t droworg = drow;
    1399             :   Float_t dcolorg = dcol;
    1400             :   
    1401           0 :   if(gGeoManager){
    1402             :     //Recover some stuff
    1403             : 
    1404           0 :     const Int_t nSMod = fEMCGeometry->GetNumberOfSuperModules();
    1405             :  
    1406           0 :     gGeoManager->cd("ALIC_1/XEN1_1");
    1407           0 :     TGeoNode        *geoXEn1 = gGeoManager->GetCurrentNode();
    1408           0 :     TGeoNodeMatrix  *geoSM[nSMod];        
    1409           0 :     TGeoVolume      *geoSMVol[nSMod];     
    1410           0 :     TGeoShape       *geoSMShape[nSMod];    
    1411           0 :     TGeoBBox        *geoBox[nSMod];        
    1412           0 :     TGeoMatrix      *geoSMMatrix[nSMod];       
    1413             :     
    1414           0 :     for(int iSM = 0; iSM < nSMod; iSM++) {  
    1415           0 :       geoSM[iSM]       = dynamic_cast<TGeoNodeMatrix *>(geoXEn1->GetDaughter(iSM));
    1416           0 :       geoSMVol[iSM]    = geoSM[iSM]->GetVolume(); 
    1417           0 :       geoSMShape[iSM]  = geoSMVol[iSM]->GetShape();
    1418           0 :       geoBox[iSM]      = dynamic_cast<TGeoBBox *>(geoSMShape[iSM]);
    1419           0 :       geoSMMatrix[iSM] = geoSM[iSM]->GetMatrix();
    1420             :     }
    1421             :     
    1422           0 :     if(sm % 2 == 0) {
    1423           0 :       dcol = 47. - dcol;
    1424           0 :       drow = 23. - drow;
    1425           0 :     }
    1426             :     
    1427             :     Int_t istrip = 0;
    1428             :     Float_t z0   = 0;
    1429             :     Float_t zb   = 0;
    1430             :     Float_t zIs = 0;
    1431             :     
    1432             :     Float_t x,y,z; // return variables in terry's RF
    1433             :     
    1434             :     //***********************************************************
    1435             :     //Do not like this: too many hardcoded values, is it not already stored somewhere else?
    1436             :     //                : need more comments in the code 
    1437             :     //***********************************************************
    1438             :     
    1439             :     Float_t dz = 6.0;   // base cell width in eta
    1440             :     Float_t dx = 6.004; // base cell width in phi
    1441             :     
    1442             :     
    1443             :     //Float_t L = 26.04; // active tower length for hadron (lead+scint+paper)
    1444             :     // we use the geant numbers 13.87*2=27.74
    1445             :     Float_t teta1 = 0.;
    1446             :       
    1447             :     //Do some basic checks
    1448           0 :     if (dcol >= 47.5 || dcol<-0.5) {
    1449           0 :       AliError(Form("Bad tower coordinate dcol=%f, where dcol >= 47.5 || dcol<-0.5; org: %f", dcol, dcolorg));
    1450           0 :       return;
    1451             :     }
    1452           0 :     if (drow >= 23.5 || drow<-0.5) {
    1453           0 :       AliError(Form("Bad tower coordinate drow=%f, where drow >= 23.5 || drow<-0.5; org: %f", drow, droworg));
    1454           0 :       return;
    1455             :     }
    1456           0 :     if (sm >= nSMod || sm < 0) {
    1457           0 :       AliError(Form("Bad SM number sm=%d, where sm >= %d || sm < 0", nSMod, sm));
    1458           0 :       return;
    1459             :     }    
    1460             :     
    1461           0 :     istrip = int ((dcol+0.5)/2);
    1462             :     
    1463             :     // tapering angle
    1464           0 :     teta1 = TMath::DegToRad() * istrip * 1.5;
    1465             :     
    1466             :     // calculation of module corner along z 
    1467             :     // as a function of strip
    1468             :     
    1469           0 :     for (int is=0; is<= istrip; is++) {
    1470             :       
    1471           0 :       teta1 = TMath::DegToRad() * (is*1.5 + 0.75);
    1472           0 :       if(is==0)
    1473           0 :         zIs = zIs + 2*dz*TMath::Cos(teta1);
    1474             :       else
    1475           0 :         zIs = zIs + 2*dz*TMath::Cos(teta1) + 2*dz*TMath::Sin(teta1)*TMath::Tan(teta1-0.75*TMath::DegToRad());
    1476             :       
    1477             :     }
    1478             :     
    1479           0 :     z0 = dz*(dcol-2*istrip+0.5);
    1480           0 :     zb = (2*dz-z0-depth*TMath::Tan(teta1));
    1481             :     
    1482           0 :     z = zIs - zb*TMath::Cos(teta1);
    1483           0 :     y = depth/TMath::Cos(teta1) + zb*TMath::Sin(teta1);
    1484             :     
    1485           0 :     x = (drow + 0.5)*dx;
    1486             :     
    1487             :     // moving the origin from terry's RF
    1488             :     // to the GEANT one
    1489             :     
    1490           0 :     double xx =  y - geoBox[sm]->GetDX();
    1491           0 :     double yy = -x + geoBox[sm]->GetDY(); 
    1492           0 :     double zz =  z - geoBox[sm]->GetDZ(); 
    1493           0 :     const double localIn[3] = {xx, yy, zz};
    1494           0 :     double dglobal[3];
    1495             :     //geoSMMatrix[sm]->Print();
    1496             :     //printf("TFF Local    (row = %d, col = %d, x = %3.2f,  y = %3.2f, z = %3.2f)\n", iroworg, icolorg, localIn[0], localIn[1], localIn[2]);
    1497           0 :     geoSMMatrix[sm]->LocalToMaster(localIn, dglobal);
    1498             :     //printf("TFF Global   (row = %2.0f, col = %2.0f, x = %3.2f,  y = %3.2f, z = %3.2f)\n", drow, dcol, dglobal[0], dglobal[1], dglobal[2]);
    1499             :     
    1500             :     //apply global shifts
    1501           0 :     if(sm == 2 || sm == 3) {//sector 1
    1502           0 :       global[0] = dglobal[0] + misaligTransShifts[3] + misaligRotShifts[3]*TMath::Sin(TMath::DegToRad()*20) ; 
    1503           0 :       global[1] = dglobal[1] + misaligTransShifts[4] + misaligRotShifts[4]*TMath::Cos(TMath::DegToRad()*20) ; 
    1504           0 :       global[2] = dglobal[2] + misaligTransShifts[5];
    1505           0 :     }
    1506           0 :     else if(sm == 0 || sm == 1){//sector 0
    1507           0 :       global[0] = dglobal[0] + misaligTransShifts[0]; 
    1508           0 :       global[1] = dglobal[1] + misaligTransShifts[1]; 
    1509           0 :       global[2] = dglobal[2] + misaligTransShifts[2];
    1510           0 :     }
    1511             :     else {
    1512           0 :       AliInfo("Careful, correction not implemented yet!");
    1513           0 :       global[0] = dglobal[0] ;
    1514           0 :       global[1] = dglobal[1] ;
    1515           0 :       global[2] = dglobal[2] ;
    1516             :     }
    1517           0 :   }
    1518             :   else{
    1519           0 :     AliFatal("Geometry boxes information, check that geometry.root is loaded\n");
    1520             :   }
    1521           0 : }
    1522             : 
    1523             : //__________________________________________________________________________________________________________________
    1524             : void AliEMCALGeometry::SetMisalMatrix(const TGeoHMatrix * m, Int_t smod) 
    1525             : {
    1526             :   // Method to set shift-rotational matrixes from ESDHeader
    1527             :   // Move from header due to coding violations : Dec 2,2011 by PAI
    1528           0 :   fUseExternalMatrices = kTRUE;
    1529             : 
    1530           0 :   if (smod >= 0 && smod < fEMCGeometry->GetNumberOfSuperModules()){
    1531           0 :     if(!fkSModuleMatrix[smod]) fkSModuleMatrix[smod] = new TGeoHMatrix(*m) ; //Set only if not set yet
    1532           0 :   } else AliFatal(Form("Wrong supermodule index -> %d",smod));
    1533           0 : }
    1534             : 
    1535             : //__________________________________________________________________________________________________________________
    1536             : Bool_t AliEMCALGeometry::IsDCALSM(Int_t iSupMod) const
    1537             : {
    1538           0 :   if( fEMCGeometry->GetEMCSystem()[iSupMod] == kDCAL_Standard || fEMCGeometry->GetEMCSystem()[iSupMod] == kDCAL_Ext ) return kTRUE;
    1539           0 :   return kFALSE;
    1540           0 : }
    1541             : 
    1542             : //__________________________________________________________________________________________________________________
    1543             : Bool_t AliEMCALGeometry::IsDCALExtSM(Int_t iSupMod) const
    1544             : {
    1545           0 :   if( fEMCGeometry->GetEMCSystem()[iSupMod] == kDCAL_Ext ) return kTRUE;
    1546           0 :   return kFALSE;
    1547           0 : }

Generated by: LCOV version 1.11