X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PHOS%2FAliPHOSGeometry.cxx;h=0fae4bfa94c6e4aa879b42b97183d5f64df1f89a;hb=f79e2d8626030694970cbe186947a6b409261f4b;hp=088fc60d3a33d573e646ac58f27ae20c683b2239;hpb=6ad0bfa0999be42354e6caa187cf84fabe5cfa39;p=u%2Fmrichter%2FAliRoot.git diff --git a/PHOS/AliPHOSGeometry.cxx b/PHOS/AliPHOSGeometry.cxx index 088fc60d3a3..0fae4bfa94c 100644 --- a/PHOS/AliPHOSGeometry.cxx +++ b/PHOS/AliPHOSGeometry.cxx @@ -13,511 +13,287 @@ * provided "as is" without express or implied warranty. * **************************************************************************/ +/* $Id$ */ + //_________________________________________________________________________ -// Geometry class for PHOS version SUBATECH -//*-- Author : Y. Schutz SUBATECH -////////////////////////////////////////////////////////////////////////////// +// Geometry class for PHOS : singleton +// PHOS consists of the electromagnetic calorimeter (EMCA) +// and a charged particle veto either in the Subatech's version (PPSD) +// or in the IHEP's one (CPV). +// The EMCA/PPSD/CPV modules are parametrized so that any configuration +// can be easily implemented +// The title is used to identify the version of CPV used. +// +// -- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (RRC "KI" & SUBATECH) // --- ROOT system --- #include "TVector3.h" #include "TRotation.h" +#include "TParticle.h" +#include +#include // --- Standard library --- -#include -#include - // --- AliRoot header files --- - +#include "AliLog.h" #include "AliPHOSGeometry.h" -#include "AliPHOSPpsdRecPoint.h" -#include "AliConst.h" +#include "AliPHOSEMCAGeometry.h" +#include "AliPHOSRecPoint.h" ClassImp(AliPHOSGeometry) - AliPHOSGeometry * AliPHOSGeometry::fGeom = 0 ; - -//____________________________________________________________________________ -AliPHOSGeometry::~AliPHOSGeometry(void) -{ - fRotMatrixArray->Delete() ; - delete fRotMatrixArray ; -} +// these initialisations are needed for a singleton +AliPHOSGeometry * AliPHOSGeometry::fgGeom = 0 ; +Bool_t AliPHOSGeometry::fgInit = kFALSE ; //____________________________________________________________________________ -Bool_t AliPHOSGeometry::AbsToRelNumbering(const Int_t AbsId, Int_t * RelId) -{ - // RelId[0] = PHOS Module number 1:fNModules - // RelId[1] = 0 if PbW04 - // = PPSD Module number 1:fNumberOfModulesPhi*fNumberOfModulesZ*2 (2->up and bottom level) - // RelId[2] = Row number inside a PHOS or PPSD module - // RelId[3] = Column number inside a PHOS or PPSD module - - Bool_t rv = kTRUE ; - Float_t Id = AbsId ; - - Int_t PHOSModuleNumber = (Int_t)TMath:: Ceil( Id / ( GetNPhi() * GetNZ() ) ) ; - - if ( PHOSModuleNumber > GetNModules() ) { // its a PPSD pad - - Id -= GetNPhi() * GetNZ() * GetNModules() ; - Float_t tempo = 2 * GetNumberOfModulesPhi() * GetNumberOfModulesZ() * GetNumberOfPadsPhi() * GetNumberOfPadsZ() ; - RelId[0] = (Int_t)TMath::Ceil( Id / tempo ) ; - Id -= ( RelId[0] - 1 ) * tempo ; - RelId[1] = (Int_t)TMath::Ceil( Id / ( GetNumberOfPadsPhi() * GetNumberOfPadsZ() ) ) ; - Id -= ( RelId[1] - 1 ) * GetNumberOfPadsPhi() * GetNumberOfPadsZ() ; - RelId[2] = (Int_t)TMath::Ceil( Id / GetNumberOfPadsPhi() ) ; - RelId[3] = (Int_t) ( Id - ( RelId[2] - 1 ) * GetNumberOfPadsPhi() ) ; - } - else { // its a PW04 crystal - - RelId[0] = PHOSModuleNumber ; - RelId[1] = 0 ; - Id -= ( PHOSModuleNumber - 1 ) * GetNPhi() * GetNZ() ; - RelId[2] = (Int_t)TMath::Ceil( Id / GetNPhi() ) ; - RelId[3] = (Int_t)( Id - ( RelId[2] - 1 ) * GetNPhi() ) ; - } - return rv ; -} -//____________________________________________________________________________ -void AliPHOSGeometry::EmcModuleCoverage(const Int_t mod, Double_t & tm, Double_t & tM, Double_t & pm, Double_t & pM, Option_t * opt) -{ - // calculates the angular coverage in theta and phi of a EMC module - - Double_t conv ; - if ( opt == kRadian ) - conv = 1. ; - else if ( opt == kDegre ) - conv = 180. / TMath::Pi() ; - else { - cout << " AliPHOSGeometry::EmcXtalCoverage : " << opt << " unknown option; result in radian " << endl ; - conv = 1. ; - } - - Float_t phi = GetPHOSAngle(mod) * (TMath::Pi() / 180.) ; - Float_t Y0 = GetIPtoOuterCoverDistance() + GetUpperPlateThickness() - + GetSecondUpperPlateThickness() + GetUpperCoolingPlateThickness() ; - - Double_t angle = TMath::ATan( GetCrystalSize(0)*GetNPhi() / (2 * Y0) ) ; - phi = phi + 1.5 * TMath::Pi() ; // to follow the convention of the particle generator(PHOS is between 230 and 310 deg.) - Double_t m = phi - angle ; - Double_t M = phi + angle ; - pM = TMath::Max(M, m) * conv ; - pm = TMath::Min(M, m) * conv ; - - angle = TMath::ATan( GetCrystalSize(2)*GetNZ() / (2 * Y0) ) ; - M = TMath::Pi() / 2. + angle ; // to follow the convention of the particle generator(PHOS is at 90 deg.) - m = TMath::Pi() / 2. - angle ; - tM = TMath::Max(M, m) * conv ; - tm = TMath::Min(M, m) * conv ; - -} - -//____________________________________________________________________________ -void AliPHOSGeometry::EmcXtalCoverage(Double_t & theta, Double_t & phi, Option_t * opt) +AliPHOSGeometry::AliPHOSGeometry() : + AliPHOSGeoUtils(), + fAngle(0.f), + fPHOSAngle(0), + fIPtoUpperCPVsurface(0), + fCrystalShift(0), + fCryCellShift(0), + fRotMatrixArray(0) { - // calculates the angular coverage in theta and phi of a single crystal in a EMC module - - Double_t conv ; - if ( opt == kRadian ) - conv = 1. ; - else if ( opt == kDegre ) - conv = 180. / TMath::Pi() ; - else { - cout << " AliPHOSGeometry::EmcXtalCoverage : " << opt << " unknown option; result in radian " << endl ; - conv = 1. ; - } + // default ctor + // must be kept public for root persistency purposes, but should never be called by the outside world + fgGeom = 0 ; - Float_t Y0 = GetIPtoOuterCoverDistance() + GetUpperPlateThickness() - + GetSecondUpperPlateThickness() + GetUpperCoolingPlateThickness() ; - theta = 2 * TMath::ATan( GetCrystalSize(2) / (2 * Y0) ) * conv ; - phi = 2 * TMath::ATan( GetCrystalSize(0) / (2 * Y0) ) * conv ; -} - + fPHOSParams[0] = 0.; + fPHOSParams[1] = 0.; + fPHOSParams[2] = 0.; + fPHOSParams[3] = 0.; +} //____________________________________________________________________________ -void AliPHOSGeometry::ImpactOnEmc(const Double_t theta, const Double_t phi, Int_t & ModuleNumber, Double_t & z, Double_t & x) +AliPHOSGeometry::AliPHOSGeometry(const AliPHOSGeometry & rhs) + : AliPHOSGeoUtils(rhs), + fAngle(rhs.fAngle), + fPHOSAngle(0), + fIPtoUpperCPVsurface(rhs.fIPtoUpperCPVsurface), + fCrystalShift(rhs.fCrystalShift), + fCryCellShift(rhs.fCryCellShift), + fRotMatrixArray(0) { - // calculates the impact coordinates of a neutral particle - // emitted in direction theta and phi in ALICE - - // searches for the PHOS EMC module - ModuleNumber = 0 ; - Double_t tm, tM, pm, pM ; - Int_t index = 1 ; - while ( ModuleNumber == 0 && index <= GetNModules() ) { - EmcModuleCoverage(index, tm, tM, pm, pM) ; - if ( (theta >= tm && theta <= tM) && (phi >= pm && phi <= pM ) ) - ModuleNumber = index ; - index++ ; - } - if ( ModuleNumber != 0 ) { - Float_t phi0 = GetPHOSAngle(ModuleNumber) * (TMath::Pi() / 180.) + 1.5 * TMath::Pi() ; - Float_t Y0 = GetIPtoOuterCoverDistance() + GetUpperPlateThickness() - + GetSecondUpperPlateThickness() + GetUpperCoolingPlateThickness() ; - Double_t angle = phi - phi0; - x = Y0 * TMath::Tan(angle) ; - angle = theta - TMath::Pi() / 2 ; - z = Y0 * TMath::Tan(angle) ; - } + Fatal("cpy ctor", "not implemented") ; } //____________________________________________________________________________ -void AliPHOSGeometry::GetGlobal(const AliRecPoint* RecPoint, TVector3 & gpos, TMatrix & gmat) -{ - - AliPHOSRecPoint * tmpPHOS = (AliPHOSRecPoint *) RecPoint ; - TVector3 LocalPosition ; - - tmpPHOS->GetLocalPosition(gpos) ; - - - if ( tmpPHOS->IsEmc() ) // it is a EMC crystal - { gpos.SetY( -(GetIPtoOuterCoverDistance() + GetUpperPlateThickness() + - GetSecondUpperPlateThickness() + GetUpperCoolingPlateThickness()) ) ; - - } - else - { // it is a PPSD pad - AliPHOSPpsdRecPoint * tmpPpsd = (AliPHOSPpsdRecPoint *) RecPoint ; - if (tmpPpsd->GetUp() ) // it is an upper module - { - gpos.SetY(-( GetIPtoOuterCoverDistance() - GetMicromegas2Thickness() - - GetLeadToMicro2Gap() - GetLeadConverterThickness() - - GetMicro1ToLeadGap() - GetMicromegas1Thickness() / 2.0 ) ) ; - } - else // it is a lower module - gpos.SetY(-( GetIPtoOuterCoverDistance() - GetMicromegas2Thickness() / 2.0) ) ; - } - - Float_t Phi = GetPHOSAngle( tmpPHOS->GetPHOSMod()) ; - Double_t const RADDEG = 180.0 / kPI ; - Float_t rPhi = Phi / RADDEG ; - - TRotation Rot ; - Rot.RotateZ(-rPhi) ; // a rotation around Z by angle - - TRotation dummy = Rot.Invert() ; // to transform from original frame to rotate frame - gpos.Transform(Rot) ; // rotate the baby - +AliPHOSGeometry::AliPHOSGeometry(const Text_t* name, const Text_t* title) + : AliPHOSGeoUtils(name, title), + fAngle(0.f), + fPHOSAngle(0), + fIPtoUpperCPVsurface(0), + fCrystalShift(0), + fCryCellShift(0), + fRotMatrixArray(0) +{ + // ctor only for internal usage (singleton) + Init() ; + fgGeom = this; } //____________________________________________________________________________ -void AliPHOSGeometry::GetGlobal(const AliRecPoint* RecPoint, TVector3 & gpos) +AliPHOSGeometry::~AliPHOSGeometry(void) { - AliPHOSRecPoint * tmpPHOS = (AliPHOSRecPoint *) RecPoint ; - TVector3 LocalPosition ; - tmpPHOS->GetLocalPosition(gpos) ; + // dtor - - if ( tmpPHOS->IsEmc() ) // it is a EMC crystal - { gpos.SetY( -(GetIPtoOuterCoverDistance() + GetUpperPlateThickness() + - GetSecondUpperPlateThickness() + GetUpperCoolingPlateThickness()) ) ; - } - else - { // it is a PPSD pad - AliPHOSPpsdRecPoint * tmpPpsd = (AliPHOSPpsdRecPoint *) RecPoint ; - if (tmpPpsd->GetUp() ) // it is an upper module - { - gpos.SetY(-( GetIPtoOuterCoverDistance() - GetMicromegas2Thickness() - - GetLeadToMicro2Gap() - GetLeadConverterThickness() - - GetMicro1ToLeadGap() - GetMicromegas1Thickness() / 2.0 ) ) ; - } - else // it is a lower module - gpos.SetY(-( GetIPtoOuterCoverDistance() - GetMicromegas2Thickness() / 2.0) ) ; - } - - Float_t Phi = GetPHOSAngle( tmpPHOS->GetPHOSMod()) ; - Double_t const RADDEG = 180.0 / kPI ; - Float_t rPhi = Phi / RADDEG ; - - TRotation Rot ; - Rot.RotateZ(-rPhi) ; // a rotation around Z by angle - - TRotation dummy = Rot.Invert() ; // to transform from original frame to rotate frame - gpos.Transform(Rot) ; // rotate the baby + if (fRotMatrixArray) fRotMatrixArray->Delete() ; + if (fRotMatrixArray) delete fRotMatrixArray ; + if (fPHOSAngle ) delete[] fPHOSAngle ; } //____________________________________________________________________________ void AliPHOSGeometry::Init(void) { - fRotMatrixArray = new TObjArray(fNModules) ; - - cout << "PHOS geometry setup: parameters for option " << fName << " " << fTitle << endl ; - if ( ((strcmp( fName, "default" )) == 0) || ((strcmp( fName, "GPS2" )) == 0) ) { - fInit = kTRUE ; - this->InitPHOS() ; - this->InitPPSD() ; - this->SetPHOSAngles() ; - } - else { - fInit = kFALSE ; - cout << "PHOS Geometry setup: option not defined " << fName << endl ; - } -} - -//____________________________________________________________________________ -void AliPHOSGeometry::InitPHOS(void) -{ - // PHOS - - fNPhi = 64 ; - fNZ = 64 ; - fNModules = 5 ; + // Initializes the PHOS parameters : + // IHEP is the Protvino CPV (cathode pad chambers) - fPHOSAngle[0] = 0.0 ; // Module position angles are set in CreateGeometry() - fPHOSAngle[1] = 0.0 ; - fPHOSAngle[2] = 0.0 ; - fPHOSAngle[3] = 0.0 ; - - fXtlSize[0] = 2.2 ; - fXtlSize[1] = 18.0 ; - fXtlSize[2] = 2.2 ; + fgInit = kTRUE ; - // all these numbers coming next are subject to changes + fAngle = 20; - fOuterBoxThickness[0] = 2.8 ; - fOuterBoxThickness[1] = 5.0 ; - fOuterBoxThickness[2] = 5.0 ; - - fUpperPlateThickness = 4.0 ; - - fSecondUpperPlateThickness = 5.0 ; - fCrystalSupportHeight = 6.95 ; - fCrystalWrapThickness = 0.01 ; - fCrystalHolderThickness = 0.005 ; - fModuleBoxThickness = 2.0 ; - fIPtoOuterCoverDistance = 447.0 ; - fIPtoCrystalSurface = 460.0 ; + fPHOSAngle = new Float_t[fNModules] ; - fPinDiodeSize[0] = 1.0 ; - fPinDiodeSize[1] = 0.1 ; - fPinDiodeSize[2] = 1.0 ; + const Float_t * emcParams = fGeometryEMCA->GetEMCParams() ; - fUpperCoolingPlateThickness = 0.06 ; - fSupportPlateThickness = 10.0 ; - fLowerThermoPlateThickness = 3.0 ; - fLowerTextolitPlateThickness = 1.0 ; - fGapBetweenCrystals = 0.03 ; + fPHOSParams[0] = TMath::Max((Double_t)fGeometryCPV->GetCPVBoxSize(0)/2., + (Double_t)(emcParams[0] - (emcParams[1]-emcParams[0])* + fGeometryCPV->GetCPVBoxSize(1)/2/emcParams[3])); + fPHOSParams[1] = emcParams[1] ; + fPHOSParams[2] = TMath::Max((Double_t)emcParams[2], (Double_t)fGeometryCPV->GetCPVBoxSize(2)/2.); + fPHOSParams[3] = emcParams[3] + fGeometryCPV->GetCPVBoxSize(1)/2. ; - fTextolitBoxThickness[0] = 1.5 ; - fTextolitBoxThickness[1] = 0.0 ; - fTextolitBoxThickness[2] = 3.0 ; - - fAirThickness[0] = 1.56 ; - fAirThickness[1] = 20.5175 ; - fAirThickness[2] = 2.48 ; - - Float_t XtalModulePhiSize = fNPhi * ( fXtlSize[0] + 2 * fGapBetweenCrystals ) ; - Float_t XtalModuleZSize = fNZ * ( fXtlSize[2] + 2 * fGapBetweenCrystals ) ; - - // The next dimensions are calculated from the above parameters - - fOuterBoxSize[0] = XtalModulePhiSize + 2 * ( fAirThickness[0] + fModuleBoxThickness - + fTextolitBoxThickness[0] + fOuterBoxThickness[0] ) ; - fOuterBoxSize[1] = ( fXtlSize[1] + fCrystalSupportHeight + fCrystalWrapThickness + fCrystalHolderThickness ) - + 2 * (fAirThickness[1] + fModuleBoxThickness + fTextolitBoxThickness[1] + fOuterBoxThickness[1] ) ; - fOuterBoxSize[2] = XtalModuleZSize + 2 * ( fAirThickness[2] + fModuleBoxThickness - + fTextolitBoxThickness[2] + fOuterBoxThickness[2] ) ; - - fTextolitBoxSize[0] = fOuterBoxSize[0] - 2 * fOuterBoxThickness[0] ; - fTextolitBoxSize[1] = fOuterBoxSize[1] - fOuterBoxThickness[1] - fUpperPlateThickness ; - fTextolitBoxSize[2] = fOuterBoxSize[2] - 2 * fOuterBoxThickness[2] ; - - fAirFilledBoxSize[0] = fTextolitBoxSize[0] - 2 * fTextolitBoxThickness[0] ; - fAirFilledBoxSize[1] = fTextolitBoxSize[1] - fSecondUpperPlateThickness ; - fAirFilledBoxSize[2] = fTextolitBoxSize[2] - 2 * fTextolitBoxThickness[2] ; + fIPtoUpperCPVsurface = fGeometryEMCA->GetIPtoOuterCoverDistance() - fGeometryCPV->GetCPVBoxSize(1) ; + + //calculate offset to crystal surface + const Float_t * inthermo = fGeometryEMCA->GetInnerThermoHalfSize() ; + const Float_t * strip = fGeometryEMCA->GetStripHalfSize() ; + const Float_t * splate = fGeometryEMCA->GetSupportPlateHalfSize(); + const Float_t * crystal = fGeometryEMCA->GetCrystalHalfSize() ; + const Float_t * pin = fGeometryEMCA->GetAPDHalfSize() ; + const Float_t * preamp = fGeometryEMCA->GetPreampHalfSize() ; + fCrystalShift=-inthermo[1]+strip[1]+splate[1]+crystal[1]-fGeometryEMCA->GetAirGapLed()/2.+pin[1]+preamp[1] ; + fCryCellShift=crystal[1]-(fGeometryEMCA->GetAirGapLed()-2*pin[1]-2*preamp[1])/2; + + Int_t index ; + for ( index = 0; index < fNModules; index++ ) + fPHOSAngle[index] = 0.0 ; // Module position angles are set in CreateGeometry() -} + fRotMatrixArray = new TObjArray(fNModules) ; -//____________________________________________________________________________ -void AliPHOSGeometry::InitPPSD(void) -{ - // PPSD + // Geometry parameters are calculated + + SetPHOSAngles(); + Double_t const kRADDEG = 180.0 / TMath::Pi() ; + Float_t r = GetIPtoOuterCoverDistance() + fPHOSParams[3] - GetCPVBoxSize(1) ; + for (Int_t iModule=0; iModule( fgGeom ) ; } //____________________________________________________________________________ AliPHOSGeometry * AliPHOSGeometry::GetInstance(const Text_t* name, const Text_t* title) { + // Returns the pointer of the unique instance + // Creates it with the specified options (name, title) if it does not exist yet + AliPHOSGeometry * rv = 0 ; - if ( fGeom == 0 ) { - fGeom = new AliPHOSGeometry(name, title) ; - rv = (AliPHOSGeometry * ) fGeom ; + if ( fgGeom == 0 ) { + if ( strcmp(name,"") == 0 ) + rv = 0 ; + else { + fgGeom = new AliPHOSGeometry(name, title) ; + if ( fgInit ) + rv = (AliPHOSGeometry * ) fgGeom ; + else { + rv = 0 ; + delete fgGeom ; + fgGeom = 0 ; + } + } } else { - if ( strcmp(fGeom->GetName(), name) != 0 ) { - cout << "AliPHOSGeometry : current geometry is " << fGeom->GetName() << endl - << " you cannot call " << name << endl ; - } + if ( strcmp(fgGeom->GetName(), name) != 0 ) + ::Error("GetInstance", "Current geometry is %s. You cannot call %s", + fgGeom->GetName(), name) ; else - rv = (AliPHOSGeometry *) fGeom ; + rv = (AliPHOSGeometry *) fgGeom ; } return rv ; } //____________________________________________________________________________ -Bool_t AliPHOSGeometry::RelToAbsNumbering(const Int_t * RelId, Int_t & AbsId) -{ - - // AbsId = 1:fNModules * fNPhi * fNZ -> PbWO4 - // AbsId = 1:fNModules * 2 * (fNumberOfModulesPhi * fNumberOfModulesZ) * fNumberOfPadsPhi * fNumberOfPadsZ -> PPSD +void AliPHOSGeometry::SetPHOSAngles() +{ + // Calculates the position of the PHOS modules in ALICE global coordinate system + // in ideal geometry + + Double_t const kRADDEG = 180.0 / TMath::Pi() ; + Float_t pphi = 2 * TMath::ATan( GetOuterBoxSize(0) / ( 2.0 * GetIPtoUpperCPVsurface() ) ) ; + pphi *= kRADDEG ; + if (pphi > fAngle){ + AliError(Form("PHOS modules overlap!\n pphi = %f fAngle = %f", + pphi, fAngle)); - Bool_t rv = kTRUE ; - - if ( RelId[1] > 0 ) { // its a PPSD pad - - AbsId = GetNPhi() * GetNZ() * GetNModules() // the offset to separate emcal crystals from PPSD pads - + ( RelId[0] - 1 ) * GetNumberOfModulesPhi() * GetNumberOfModulesZ() // the pads offset of PHOS modules - * GetNumberOfPadsPhi() * GetNumberOfPadsZ() * 2 - + ( RelId[1] - 1 ) * GetNumberOfPadsPhi() * GetNumberOfPadsZ() // the pads offset of PPSD modules - + ( RelId[2] - 1 ) * GetNumberOfPadsPhi() // the pads offset of a PPSD row - + RelId[3] ; // the column number - } - else { - if ( RelId[1] == 0 ) { // its a Phos crystal - AbsId = ( RelId[0] - 1 ) * GetNPhi() * GetNZ() // the offset of PHOS modules - + ( RelId[2] - 1 ) * GetNPhi() // the offset of a xtal row - + RelId[3] ; // the column number - } } - - return rv ; + pphi = fAngle; + + for( Int_t i = 1; i <= fNModules ; i++ ) { + Float_t angle = pphi * ( i - fNModules / 2.0 - 0.5 ) ; + fPHOSAngle[i-1] = - angle ; + } } - //____________________________________________________________________________ - -void AliPHOSGeometry::RelPosInAlice(const Int_t Id, TVector3 & pos ) +void AliPHOSGeometry::GetGlobal(const AliRecPoint* , TVector3 & ) const { - if (Id > 0) { + AliFatal(Form("Please use GetGlobalPHOS(recPoint,gpos) instead of GetGlobal!")); +} - Int_t RelId[4] ; +//____________________________________________________________________________ +void AliPHOSGeometry::GetGlobalPHOS(const AliPHOSRecPoint* recPoint, TVector3 & gpos) const +{ + // Calculates the coordinates of a RecPoint and the error matrix in the ALICE global coordinate system - AbsToRelNumbering(Id , RelId) ; + const AliPHOSRecPoint * tmpPHOS = recPoint ; + TVector3 localposition ; - Int_t PHOSModule = RelId[0] ; - - Float_t Y0 ; + tmpPHOS->GetLocalPosition(gpos) ; - if ( RelId[1] == 0 ) // it is a PbW04 crystal - { Y0 = -(GetIPtoOuterCoverDistance() + GetUpperPlateThickness() - + GetSecondUpperPlateThickness() + GetUpperCoolingPlateThickness()) ; + if (!gGeoManager){ + AliFatal("Geo manager not initialized\n"); } - if ( RelId[1] > 0 ) { // its a PPSD pad - if ( RelId[1] > GetNumberOfModulesPhi() * GetNumberOfModulesZ() ) // its an bottom module - { - Y0 = -( GetIPtoOuterCoverDistance() - GetMicromegas2Thickness() / 2.0) ; - } - else // its an upper module - Y0 = -( GetIPtoOuterCoverDistance() - GetMicromegas2Thickness() - GetLeadToMicro2Gap() - - GetLeadConverterThickness() - GetMicro1ToLeadGap() - GetMicromegas1Thickness() / 2.0) ; + //construct module name + char path[100] ; + Double_t dy ; + if(tmpPHOS->IsEmc()){ + TString spath="/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1"; + snprintf(path,spath.Length(),spath.Data(),tmpPHOS->GetPHOSMod()) ; +// sprintf(path,"/ALIC_1/PHOS_%d",tmpPHOS->GetPHOSMod()) ; + dy=fCrystalShift ; } - - Float_t x, z ; - RelPosInModule(RelId, x, z) ; - - pos.SetX(x) ; - pos.SetZ(z) ; - pos.SetY( TMath::Sqrt(x*x + z*z + Y0*Y0) ) ; - - - - Float_t Phi = GetPHOSAngle( PHOSModule) ; - Double_t const RADDEG = 180.0 / kPI ; - Float_t rPhi = Phi / RADDEG ; - - TRotation Rot ; - Rot.RotateZ(-rPhi) ; // a rotation around Z by angle - - TRotation dummy = Rot.Invert() ; // to transform from original frame to rotate frame - - pos.Transform(Rot) ; // rotate the baby + else{ + TString spath="/ALIC_1/PHOS_%d/PCPV_1"; + snprintf(path,spath.Length(),spath.Data(),tmpPHOS->GetPHOSMod()) ; + dy= GetCPVBoxSize(1)/2. ; //center of CPV module } - else { - pos.SetX(0.); - pos.SetY(0.); - pos.SetZ(0.); - } -} + Double_t pos[3]={gpos.X(),gpos.Y()-dy,gpos.Z()} ; + if(tmpPHOS->IsEmc()) + pos[2]=-pos[2] ; //Opposite z directions in EMC matrix and local frame!!! + Double_t posC[3] = {}; + //now apply possible shifts and rotations + if (!gGeoManager->cd(path)){ + AliFatal("Geo manager can not find path \n"); + } + TGeoHMatrix *m = gGeoManager->GetCurrentMatrix(); + if (m){ + m->LocalToMaster(pos,posC); + } + else{ + AliFatal("Geo matrixes are not loaded \n") ; + } + gpos.SetXYZ(posC[0],posC[1],posC[2]) ; -//____________________________________________________________________________ -void AliPHOSGeometry::RelPosInModule(const Int_t * RelId, Float_t & x, Float_t & z) -{ - Int_t PPSDModule ; - Int_t Row = RelId[2] ; //offset along z axiz - Int_t Column = RelId[3] ; //offset along x axiz - - Float_t PadSizeZ = GetPPSDModuleSize(2)/ GetNumberOfPadsZ(); - Float_t PadSizeX = GetPPSDModuleSize(0)/ GetNumberOfPadsPhi(); - - if ( RelId[1] == 0 ) { // its a PbW04 crystal - x = -( GetNPhi()/2. - Row + 0.5 ) * GetCrystalSize(0) ; // position ox Xtal with respect - z = ( GetNZ() /2. - Column + 0.5 ) * GetCrystalSize(2) ; // of center of PHOS module - } - else { - if ( RelId[1] > GetNumberOfModulesPhi() * GetNumberOfModulesZ() ) - PPSDModule = RelId[1]-GetNumberOfModulesPhi() * GetNumberOfModulesZ(); - else PPSDModule = RelId[1] ; - Int_t ModRow = 1+(Int_t)TMath::Ceil( (Float_t)PPSDModule / GetNumberOfModulesPhi()-1. ) ; - Int_t ModCol = PPSDModule - ( ModRow-1 ) * GetNumberOfModulesPhi() ; - Float_t x0 = ( GetNumberOfModulesPhi() / 2. - ModRow + 0.5 ) * GetPPSDModuleSize(0) ; - Float_t z0 = ( GetNumberOfModulesZ() / 2. - ModCol + 0.5 ) * GetPPSDModuleSize(2) ; - x = - ( GetNumberOfPadsPhi()/2. - Row - 0.5 ) * PadSizeX + x0 ; // position of pad with respect - z = ( GetNumberOfPadsZ()/2. - Column - 0.5 ) * PadSizeZ - z0 ; // of center of PHOS module - } } - //____________________________________________________________________________ -void AliPHOSGeometry:: SetPHOSAngles() -{ - Double_t const RADDEG = 180.0 / kPI ; - Float_t PPHI = TMath::ATan( fOuterBoxSize[0] / ( 2.0 * fIPtoOuterCoverDistance ) ) ; - PPHI *= RADDEG ; - - for( Int_t i = 1; i <= fNModules ; i++ ) { - Float_t angle = PPHI * 2 * ( i - fNModules / 2.0 - 0.5 ) ; - fPHOSAngle[i-1] = - angle ; - } + +void AliPHOSGeometry::GetModuleCenter(TVector3& center, + const char *det, + Int_t module) const +{ + // Returns a position of the center of the CPV or EMC module + // in ideal (not misaligned) geometry + Float_t rDet = 0.; + if (strcmp(det,"CPV") == 0) rDet = GetIPtoCPVDistance (); + else if (strcmp(det,"EMC") == 0) rDet = GetIPtoCrystalSurface(); + else + AliFatal(Form("Wrong detector name %s",det)); + + Float_t angle = GetPHOSAngle(module); // (40,20,0,-20,-40) degrees + angle *= TMath::Pi()/180; + angle += 3*TMath::Pi()/2.; + center.SetXYZ(rDet*TMath::Cos(angle), rDet*TMath::Sin(angle), 0.); }