X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PHOS%2FAliPHOSGeometry.cxx;h=0fae4bfa94c6e4aa879b42b97183d5f64df1f89a;hb=f1731579adf7c07dc59b41b159d7c7857a96c1c2;hp=bde9dd2639d4e7d4be4fe05792bb76442bfd54b2;hpb=22b8277f7422dfce93e0bb1f6d9c6b1e51c47cfb;p=u%2Fmrichter%2FAliRoot.git diff --git a/PHOS/AliPHOSGeometry.cxx b/PHOS/AliPHOSGeometry.cxx index bde9dd2639d..0fae4bfa94c 100644 --- a/PHOS/AliPHOSGeometry.cxx +++ b/PHOS/AliPHOSGeometry.cxx @@ -24,31 +24,77 @@ // can be easily implemented // The title is used to identify the version of CPV used. // -//*-- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (RRC "KI" & SUBATECH) +// -- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (RRC "KI" & SUBATECH) // --- ROOT system --- #include "TVector3.h" #include "TRotation.h" -#include "TFolder.h" -#include "TROOT.h" +#include "TParticle.h" +#include +#include // --- Standard library --- -#include - // --- AliRoot header files --- - +#include "AliLog.h" #include "AliPHOSGeometry.h" #include "AliPHOSEMCAGeometry.h" #include "AliPHOSRecPoint.h" -#include "AliConst.h" -ClassImp(AliPHOSGeometry) ; +ClassImp(AliPHOSGeometry) // these initialisations are needed for a singleton -AliPHOSGeometry * AliPHOSGeometry::fgGeom = 0 ; -Bool_t AliPHOSGeometry::fgInit = kFALSE ; +AliPHOSGeometry * AliPHOSGeometry::fgGeom = 0 ; +Bool_t AliPHOSGeometry::fgInit = kFALSE ; + +//____________________________________________________________________________ +AliPHOSGeometry::AliPHOSGeometry() : + AliPHOSGeoUtils(), + fAngle(0.f), + fPHOSAngle(0), + fIPtoUpperCPVsurface(0), + fCrystalShift(0), + fCryCellShift(0), + fRotMatrixArray(0) +{ + // default ctor + // must be kept public for root persistency purposes, but should never be called by the outside world + fgGeom = 0 ; + + fPHOSParams[0] = 0.; + fPHOSParams[1] = 0.; + fPHOSParams[2] = 0.; + fPHOSParams[3] = 0.; +} + +//____________________________________________________________________________ +AliPHOSGeometry::AliPHOSGeometry(const AliPHOSGeometry & rhs) + : AliPHOSGeoUtils(rhs), + fAngle(rhs.fAngle), + fPHOSAngle(0), + fIPtoUpperCPVsurface(rhs.fIPtoUpperCPVsurface), + fCrystalShift(rhs.fCrystalShift), + fCryCellShift(rhs.fCryCellShift), + fRotMatrixArray(0) +{ + Fatal("cpy ctor", "not implemented") ; +} + +//____________________________________________________________________________ +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; +} //____________________________________________________________________________ AliPHOSGeometry::~AliPHOSGeometry(void) @@ -59,49 +105,65 @@ AliPHOSGeometry::~AliPHOSGeometry(void) if (fRotMatrixArray) delete fRotMatrixArray ; if (fPHOSAngle ) delete[] fPHOSAngle ; } -//____________________________________________________________________________ +//____________________________________________________________________________ void AliPHOSGeometry::Init(void) { // Initializes the PHOS parameters : // IHEP is the Protvino CPV (cathode pad chambers) - TString test(GetName()) ; - if (test != "IHEP" ) { - Fatal("Init", "%s is not a known geometry (choose among IHEP)", test.Data() ) ; - } - fgInit = kTRUE ; - - fNModules = 5; + fAngle = 20; - - fGeometryEMCA = new AliPHOSEMCAGeometry(); - - fGeometryCPV = new AliPHOSCPVGeometry (); - - fGeometrySUPP = new AliPHOSSupportGeometry(); + fPHOSAngle = new Float_t[fNModules] ; - Float_t * emcParams = fGeometryEMCA->GetEMCParams() ; + const Float_t * emcParams = fGeometryEMCA->GetEMCParams() ; fPHOSParams[0] = TMath::Max((Double_t)fGeometryCPV->GetCPVBoxSize(0)/2., - (Double_t)(emcParams[0]*(fGeometryCPV->GetCPVBoxSize(1)+emcParams[3]) - - emcParams[1]* fGeometryCPV->GetCPVBoxSize(1))/emcParams[3] ) ; + (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. ; 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() - this->SetPHOSAngles() ; fRotMatrixArray = new TObjArray(fNModules) ; - + + // 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; iModuleGetName(), name) != 0 ) - ::Error("GetInstance", "Current geometry is %s. You cannot call %s", fgGeom->GetName(), name) ; + ::Error("GetInstance", "Current geometry is %s. You cannot call %s", + fgGeom->GetName(), name) ; else rv = (AliPHOSGeometry *) fgGeom ; } @@ -146,12 +209,14 @@ AliPHOSGeometry * AliPHOSGeometry::GetInstance(const Text_t* name, const Text_t void AliPHOSGeometry::SetPHOSAngles() { // Calculates the position of the PHOS modules in ALICE global coordinate system + // in ideal geometry - Double_t const kRADDEG = 180.0 / kPI ; + Double_t const kRADDEG = 180.0 / TMath::Pi() ; Float_t pphi = 2 * TMath::ATan( GetOuterBoxSize(0) / ( 2.0 * GetIPtoUpperCPVsurface() ) ) ; pphi *= kRADDEG ; if (pphi > fAngle){ - Error("SetPHOSAngles", "PHOS modules overlap!\n pphi = %f fAngle = %f", pphi, fAngle); + AliError(Form("PHOS modules overlap!\n pphi = %f fAngle = %f", + pphi, fAngle)); } pphi = fAngle; @@ -161,280 +226,74 @@ void AliPHOSGeometry::SetPHOSAngles() fPHOSAngle[i-1] = - angle ; } } - //____________________________________________________________________________ -Bool_t AliPHOSGeometry::AbsToRelNumbering(const Int_t AbsId, Int_t * relid) const +void AliPHOSGeometry::GetGlobal(const AliRecPoint* , TVector3 & ) const { - // Converts the absolute numbering into the following array/ - // relid[0] = PHOS Module number 1:fNModules - // relid[1] = 0 if PbW04 - // = -1 if CPV - // relid[2] = Row number inside a PHOS module - // relid[3] = Column number inside a PHOS module - - Bool_t rv = kTRUE ; - Float_t id = AbsId ; - - Int_t phosmodulenumber = (Int_t)TMath:: Ceil( id / GetNCristalsInModule() ) ; - - if ( phosmodulenumber > GetNModules() ) { // it is a CPV pad - - id -= GetNPhi() * GetNZ() * GetNModules() ; - Float_t nCPV = GetNumberOfCPVPadsPhi() * GetNumberOfCPVPadsZ() ; - relid[0] = (Int_t) TMath::Ceil( id / nCPV ) ; - relid[1] = -1 ; - id -= ( relid[0] - 1 ) * nCPV ; - relid[2] = (Int_t) TMath::Ceil( id / GetNumberOfCPVPadsZ() ) ; - relid[3] = (Int_t) ( id - ( relid[2] - 1 ) * GetNumberOfCPVPadsZ() ) ; - } - else { // it is a PW04 crystal - - relid[0] = phosmodulenumber ; - relid[1] = 0 ; - id -= ( phosmodulenumber - 1 ) * GetNPhi() * GetNZ() ; - relid[2] = (Int_t)TMath::Ceil( id / GetNZ() ) ; - relid[3] = (Int_t)( id - ( relid[2] - 1 ) * GetNZ() ) ; - } - return rv ; -} - -//____________________________________________________________________________ -void AliPHOSGeometry::EmcModuleCoverage(const Int_t mod, Double_t & tm, Double_t & tM, Double_t & pm, Double_t & pM, Option_t * opt) const -{ - // calculates the angular coverage in theta and phi of one EMC (=PHOS) module - - Double_t conv ; - if ( opt == Radian() ) - conv = 1. ; - else if ( opt == Degre() ) - conv = 180. / TMath::Pi() ; - else { - Warning("EmcModuleCoverage", "%s unknown option; result in radian", opt) ; - conv = 1. ; - } - - Float_t phi = GetPHOSAngle(mod) * (TMath::Pi() / 180.) ; - Float_t y0 = GetIPtoCrystalSurface() ; - Float_t * strip = fGeometryEMCA->GetStripHalfSize() ; - Float_t x0 = fGeometryEMCA->GetNStripX()*strip[0] ; - Float_t z0 = fGeometryEMCA->GetNStripZ()*strip[2] ; - Double_t angle = TMath::ATan( x0 / y0 ) ; - phi = phi + 1.5 * TMath::Pi() ; // to follow the convention of the particle generator(PHOS is between 220 and 320 deg.) - Double_t max = phi - angle ; - Double_t min = phi + angle ; - pM = TMath::Max(max, min) * conv ; - pm = TMath::Min(max, min) * conv ; - - angle = TMath::ATan( z0 / y0 ) ; - max = TMath::Pi() / 2. + angle ; // to follow the convention of the particle generator(PHOS is at 90 deg.) - min = TMath::Pi() / 2. - angle ; - tM = TMath::Max(max, min) * conv ; - tm = TMath::Min(max, min) * conv ; - -} - -//____________________________________________________________________________ -void AliPHOSGeometry::EmcXtalCoverage(Double_t & theta, Double_t & phi, Option_t * opt) const -{ - // calculates the angular coverage in theta and phi of a single crystal in a EMC(=PHOS) module - - Double_t conv ; - if ( opt == Radian() ) - conv = 1. ; - else if ( opt == Degre() ) - conv = 180. / TMath::Pi() ; - else { - Warning("EmcXtalCoverage", "%s unknown option; result in radian", opt) ; - conv = 1. ; - } - - Float_t y0 = GetIPtoCrystalSurface() ; - theta = 2 * TMath::ATan( GetCrystalSize(2) / (2 * y0) ) * conv ; - phi = 2 * TMath::ATan( GetCrystalSize(0) / (2 * y0) ) * conv ; + AliFatal(Form("Please use GetGlobalPHOS(recPoint,gpos) instead of GetGlobal!")); } - //____________________________________________________________________________ -void AliPHOSGeometry::GetGlobal(const AliRecPoint* RecPoint, TVector3 & gpos, TMatrix & gmat) const +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 - AliPHOSRecPoint * tmpPHOS = (AliPHOSRecPoint *) RecPoint ; + const AliPHOSRecPoint * tmpPHOS = recPoint ; TVector3 localposition ; tmpPHOS->GetLocalPosition(gpos) ; - - if ( tmpPHOS->IsEmc() ) // it is a EMC crystal - { gpos.SetY( - GetIPtoCrystalSurface()) ; - - } - else - { // it is a CPV - gpos.SetY(- GetIPtoUpperCPVsurface() ) ; - } - - Float_t phi = GetPHOSAngle( tmpPHOS->GetPHOSMod()) ; - Double_t const kRADDEG = 180.0 / kPI ; - Float_t rphi = phi / kRADDEG ; - - 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 - -} - -//____________________________________________________________________________ -void AliPHOSGeometry::GetGlobal(const AliRecPoint* RecPoint, TVector3 & gpos) const -{ - // Calculates the coordinates of a RecPoint in the ALICE global coordinate system - - AliPHOSRecPoint * tmpPHOS = (AliPHOSRecPoint *) RecPoint ; - TVector3 localposition ; - tmpPHOS->GetLocalPosition(gpos) ; - - - if ( tmpPHOS->IsEmc() ) // it is a EMC crystal - { gpos.SetY( - GetIPtoCrystalSurface() ) ; - } - else - { // it is a CPV - gpos.SetY(- GetIPtoUpperCPVsurface() ) ; - } - - Float_t phi = GetPHOSAngle( tmpPHOS->GetPHOSMod()) ; - Double_t const kRADDEG = 180.0 / kPI ; - Float_t rphi = phi / kRADDEG ; - - 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 -} - -//____________________________________________________________________________ -void AliPHOSGeometry::ImpactOnEmc(const Double_t theta, const Double_t phi, Int_t & ModuleNumber, Double_t & z, Double_t & x) const -{ - // calculates the impact coordinates on PHOS of a neutral particle - // emitted in the direction theta and phi in the ALICE global coordinate system - - // 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 (!gGeoManager){ + AliFatal("Geo manager not initialized\n"); } - if ( ModuleNumber != 0 ) { - Float_t phi0 = GetPHOSAngle(ModuleNumber) * (TMath::Pi() / 180.) + 1.5 * TMath::Pi() ; - Float_t y0 = GetIPtoCrystalSurface() ; - Double_t angle = phi - phi0; - x = y0 * TMath::Tan(angle) ; - angle = theta - TMath::Pi() / 2 ; - z = y0 * TMath::Tan(angle) ; + //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 ; } -} - -Bool_t AliPHOSGeometry::Impact(const TParticle * particle) const -{ - Bool_t In=kFALSE; - Int_t ModuleNumber=0; - Double_t z,x; - ImpactOnEmc(particle->Theta(),particle->Phi(),ModuleNumber,z,x); - if(ModuleNumber) In=kTRUE; - else In=kFALSE; - return In; -} - -//____________________________________________________________________________ -Bool_t AliPHOSGeometry::RelToAbsNumbering(const Int_t * relid, Int_t & AbsId) const -{ - // Converts the relative numbering into the absolute numbering - // EMCA crystals: - // AbsId = from 1 to fNModules * fNPhi * fNZ - // CPV pad: - // AbsId = from N(total PHOS crystals) + 1 - // to NCPVModules * fNumberOfCPVPadsPhi * fNumberOfCPVPadsZ - - Bool_t rv = kTRUE ; - - if ( relid[1] == 0 ) { // it is a Phos crystal - AbsId = - ( relid[0] - 1 ) * GetNPhi() * GetNZ() // the offset of PHOS modules - + ( relid[2] - 1 ) * GetNZ() // the offset along phi - + relid[3] ; // the offset along z + 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 { // it is a CPV pad - AbsId = GetNPhi() * GetNZ() * GetNModules() // the offset to separate EMCA crystals from CPV pads - + ( relid[0] - 1 ) * GetNumberOfCPVPadsPhi() * GetNumberOfCPVPadsZ() // the pads offset of PHOS modules - + ( relid[2] - 1 ) * GetNumberOfCPVPadsZ() // the pads offset of a CPV row - + relid[3] ; // the column number + 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"); } - - return rv ; -} + 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::RelPosInAlice(const Int_t id, TVector3 & pos ) const +void AliPHOSGeometry::GetModuleCenter(TVector3& center, + const char *det, + Int_t module) const { - // Converts the absolute numbering into the global ALICE coordinate system - - - Int_t relid[4] ; - - AbsToRelNumbering(id , relid) ; - - Int_t phosmodule = relid[0] ; - - Float_t y0 = 0 ; - - if ( relid[1] == 0 ) // it is a PbW04 crystal - y0 = - GetIPtoCrystalSurface() ; - else - y0 = - GetIPtoUpperCPVsurface() ; - - Float_t x, z ; - RelPosInModule(relid, x, z) ; - - pos.SetX(x) ; - pos.SetZ(z) ; - pos.SetY(y0) ; - - Float_t phi = GetPHOSAngle( phosmodule) ; - Double_t const kRADDEG = 180.0 / kPI ; - Float_t rphi = phi / kRADDEG ; - - 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 -} - -//____________________________________________________________________________ -void AliPHOSGeometry::RelPosInModule(const Int_t * relid, Float_t & x, Float_t & z) const -{ - // Converts the relative numbering into the local PHOS-module (x, z) coordinates - // Note: sign of z differs from that in the previous version (Yu.Kharlov, 12 Oct 2000) - - Int_t row = relid[2] ; //offset along x axis - Int_t column = relid[3] ; //offset along z axis - - - if ( relid[1] == 0 ) { // its a PbW04 crystal - x = - ( GetNPhi()/2. - row + 0.5 ) * GetCellStep() ; // position of Xtal with respect - z = ( GetNZ() /2. - column + 0.5 ) * GetCellStep() ; // of center of PHOS module - } - else { - x = - ( GetNumberOfCPVPadsPhi()/2. - row - 0.5 ) * GetPadSizePhi() ; // position of pad with respect - z = ( GetNumberOfCPVPadsZ() /2. - column - 0.5 ) * GetPadSizeZ() ; // of center of PHOS module - } + // 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.); } +