X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=EMCAL%2FAliEMCALGeometry.cxx;h=264b8ff62dbbf7eb898a3fbfdadd8da8b350434b;hb=b97c26da4b2d527d08a1b07186508293b98604c3;hp=3bb75983a05b8a28513b84f1171593dd85663626;hpb=a58900484fe5e4bec3b1e6a85c9344574d75c867;p=u%2Fmrichter%2FAliRoot.git diff --git a/EMCAL/AliEMCALGeometry.cxx b/EMCAL/AliEMCALGeometry.cxx index 3bb75983a05..264b8ff62db 100644 --- a/EMCAL/AliEMCALGeometry.cxx +++ b/EMCAL/AliEMCALGeometry.cxx @@ -1,4 +1,4 @@ - /************************************************************************** +/************************************************************************** * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * * * * Author: The ALICE Off-line Project. * @@ -35,7 +35,7 @@ #include #include -#include + //#include #include #include #include @@ -46,11 +46,11 @@ // -- ALICE Headers. //#include "AliConst.h" +#include "AliLog.h" // --- EMCAL headers #include "AliEMCALGeometry.h" #include "AliEMCALShishKebabTrd1Module.h" -//#include "AliRecPoint.h" #include "AliEMCALRecPoint.h" #include "AliEMCALDigit.h" #include "AliEMCALHistoUtilities.h" @@ -63,14 +63,22 @@ AliEMCALGeometry *AliEMCALGeometry::fgGeom = 0; Bool_t AliEMCALGeometry::fgInit = kFALSE; AliEMCALAlignData *AliEMCALGeometry::fgAlignData = 0; -TString name; // contains name of geometry -char *additionalOpts[]={"nl=", // number of sampling layers - "pbTh=", // cm, Thickness of the Pb - "scTh=" // cm, Thickness of the Sc -}; -int nAdditionalOpts = sizeof(additionalOpts) / sizeof(char*); +AliEMCALGeometry::AliEMCALGeometry() : AliGeometry() +{ + // default ctor only for internal usage (singleton) + // must be kept public for root persistency purposes, but should never be called by the outside world + // CreateListOfTrd1Modules(); + AliDebug(2, "AliEMCALGeometry : default ctor "); +} +//______________________________________________________________________ +AliEMCALGeometry::AliEMCALGeometry(const Text_t* name, const Text_t* title) : +AliGeometry(name, title) {// ctor only for internal usage (singleton) + AliDebug(2, Form("AliEMCALGeometry(%s,%s) ", name,title)); + Init(); + CreateListOfTrd1Modules(); +} //______________________________________________________________________ AliEMCALGeometry::~AliEMCALGeometry(void){ // dtor @@ -86,11 +94,19 @@ void AliEMCALGeometry::Init(void){ // SHISH_25 or SHISH_62 // 11-oct-05 - correction for pre final design // Feb 06,2006 - decrease the weight of EMCAL + + fAdditionalOpts[0] = "nl="; // number of sampling layers (fNECLayers) + fAdditionalOpts[1] = "pbTh="; // cm, Thickness of the Pb (fECPbRadThick) + fAdditionalOpts[2] = "scTh="; // cm, Thickness of the Sc (fECScintThick) + fAdditionalOpts[3] = "latSS="; // cm, Thickness of lateral steel strip (fLateralSteelStrip) + + fNAdditionalOpts = sizeof(fAdditionalOpts) / sizeof(char*); + fgInit = kFALSE; // Assume failed until proven otherwise. - name = GetName(); - name.ToUpper(); + fGeoName = GetName(); + fGeoName.ToUpper(); fKey110DEG = 0; - if(name.Contains("110DEG")) fKey110DEG = 1; // for GetAbsCellId + if(fGeoName.Contains("110DEG")) fKey110DEG = 1; // for GetAbsCellId fShishKebabTrd1Modules = 0; fTrd2AngleY = f2Trd2Dy2 = fEmptySpace = fTubsR = fTubsTurnAngle = 0; @@ -105,7 +121,7 @@ void AliEMCALGeometry::Init(void){ for(int i=0; i<12; i++) fMatrixOfSM[i] = 0; // geometry - if(name.Contains("SHISH")){ // Only shahslyk now + if(fGeoName.Contains("SHISH")){ // Only shahslyk now // 7-sep-05; integration issue fArm1PhiMin = 80.0; // 60 -> 80 fArm1PhiMax = 180.0; // 180 -> 190 @@ -125,9 +141,9 @@ void AliEMCALGeometry::Init(void){ fNECLayers = 62; fECScintThick = fECPbRadThickness = 0.2; fSampling = 1.; // 30-aug-04 - should be calculated - if(name.Contains("TWIST")) { // all about EMCAL module + if(fGeoName.Contains("TWIST")) { // all about EMCAL module fNZ = 27; // 16-sep-04 - } else if(name.Contains("TRD")) { + } else if(fGeoName.Contains("TRD")) { fIPDistance = 428.0; // 11-may-05 fSteelFrontThick = 0.0; // 3.17 -> 0.0; 28-mar-05 : no stell plate fNPhi = 12; @@ -137,12 +153,12 @@ void AliEMCALGeometry::Init(void){ fTrd1Angle = 1.3; // in degree // 18-nov-04; 1./0.08112=12.327 // http://pdsfweb01.nersc.gov/~pavlinov/ALICE/SHISHKEBAB/RES/linearityAndResolutionForTRD1.html - if(name.Contains("TRD1")) { // 30-jan-05 + if(fGeoName.Contains("TRD1")) { // 30-jan-05 // for final design fPhiGapForSM = 2.; // cm, only for final TRD1 geometry - if(name.Contains("MAY05") || name.Contains("WSUC") || name.Contains("FINAL")){ + if(fGeoName.Contains("MAY05") || fGeoName.Contains("WSUC") || fGeoName.Contains("FINAL")){ fNumberOfSuperModules = 12; // 20-may-05 - if(name.Contains("WSUC")) fNumberOfSuperModules = 1; // 27-may-05 + if(fGeoName.Contains("WSUC")) fNumberOfSuperModules = 1; // 27-may-05 fNECLayers = 77; // (13-may-05 from V.Petrov) fPhiModuleSize = 12.5; // 20-may-05 - rectangular shape fEtaModuleSize = 11.9; @@ -153,18 +169,18 @@ void AliEMCALGeometry::Init(void){ fNZ = 24; fTrd1Angle = 1.5; // 1.3 or 1.5 - if(name.Contains("FINAL")) { // 9-sep-05 + if(fGeoName.Contains("FINAL")) { // 9-sep-05 fNumberOfSuperModules = 10; - if(name.Contains("110DEG")) { + if(fGeoName.Contains("110DEG")) { fNumberOfSuperModules = 12;// last two modules have size 10 degree in phi (1801.64 @@ -178,54 +194,55 @@ void AliEMCALGeometry::Init(void){ fTubsTurnAngle = 3.; } fNPHIdiv = fNETAdiv = 2; // 13-oct-04 - division again - if(name.Contains("3X3")) { // 23-nov-04 + if(fGeoName.Contains("3X3")) { // 23-nov-04 fNPHIdiv = fNETAdiv = 3; - } else if(name.Contains("4X4")) { + } else if(fGeoName.Contains("4X4")) { fNPHIdiv = fNETAdiv = 4; } } - fPhiTileSize = fPhiModuleSize/2. - fLateralSteelStrip; // 13-may-05 - fEtaTileSize = fEtaModuleSize/2. - fLateralSteelStrip; // 13-may-05 - - if(name.Contains("25")){ + if(fGeoName.Contains("25")){ fNECLayers = 25; fECScintThick = fECPbRadThickness = 0.5; } - if(name.Contains("WSUC")){ // 18-may-05 - about common structure + if(fGeoName.Contains("WSUC")){ // 18-may-05 - about common structure fShellThickness = 30.; // should be change fNPhi = fNZ = 4; } - CheckAditionalOptions(); + CheckAdditionalOptions(); + DefineSamplingFraction(); + + fPhiTileSize = fPhiModuleSize/2. - fLateralSteelStrip; // 13-may-05 + fEtaTileSize = fEtaModuleSize/2. - fLateralSteelStrip; // 13-may-05 // constant for transition absid <--> indexes fNCellsInTower = fNPHIdiv*fNETAdiv; fNCellsInSupMod = fNCellsInTower*fNPhi*fNZ; fNCells = fNCellsInSupMod*fNumberOfSuperModules; - if(name.Contains("110DEG")) fNCells -= fNCellsInSupMod; + if(fGeoName.Contains("110DEG")) fNCells -= fNCellsInSupMod; fLongModuleSize = fNECLayers*(fECScintThick + fECPbRadThickness); - if(name.Contains("MAY05")) fLongModuleSize += (fFrontSteelStrip + fPassiveScintThick); + if(fGeoName.Contains("MAY05")) fLongModuleSize += (fFrontSteelStrip + fPassiveScintThick); // 30-sep-04 - if(name.Contains("TRD")) { + if(fGeoName.Contains("TRD")) { f2Trd1Dx2 = fEtaModuleSize + 2.*fLongModuleSize*TMath::Tan(fTrd1Angle*TMath::DegToRad()/2.); - if(name.Contains("TRD2")) { // 27-jan-05 + if(fGeoName.Contains("TRD2")) { // 27-jan-05 f2Trd2Dy2 = fPhiModuleSize + 2.*fLongModuleSize*TMath::Tan(fTrd2AngleY*TMath::DegToRad()/2.); } } - } else Fatal("Init", "%s is an undefined geometry!", name.Data()) ; + } else Fatal("Init", "%s is an undefined geometry!", fGeoName.Data()) ; fNPhiSuperModule = fNumberOfSuperModules/2; if(fNPhiSuperModule<1) fNPhiSuperModule = 1; //There is always one more scintillator than radiator layer because of the first block of aluminium fShellThickness = fAlFrontThick + fGap2Active + fNECLayers*GetECScintThick()+(fNECLayers-1)*GetECPbRadThick(); - if(name.Contains("SHISH")) { + if(fGeoName.Contains("SHISH")) { fShellThickness = fSteelFrontThick + fLongModuleSize; - if(name.Contains("TWIST")) { // 13-sep-04 + if(fGeoName.Contains("TWIST")) { // 13-sep-04 fShellThickness = TMath::Sqrt(fLongModuleSize*fLongModuleSize + fPhiModuleSize*fEtaModuleSize); fShellThickness += fSteelFrontThick; - } else if(name.Contains("TRD")) { // 1-oct-04 + } else if(fGeoName.Contains("TRD")) { // 1-oct-04 fShellThickness = TMath::Sqrt(fLongModuleSize*fLongModuleSize + f2Trd1Dx2*f2Trd1Dx2); fShellThickness += fSteelFrontThick; // Local coordinates @@ -247,18 +264,18 @@ void AliEMCALGeometry::Init(void){ fgInit = kTRUE; - if (kTRUE) { - printf("Init: geometry of EMCAL named %s is as follows:\n", name.Data()); + if (AliDebugLevel()>=2) { + printf("Init: geometry of EMCAL named %s is as follows:\n", fGeoName.Data()); printf( " ECAL : %d x (%f cm Pb, %f cm Sc) \n", GetNECLayers(), GetECPbRadThick(), GetECScintThick() ) ; printf(" fSampling %5.2f \n", fSampling ); - if(name.Contains("SHISH")){ + if(fGeoName.Contains("SHISH")){ printf(" fIPDistance %6.3f cm \n", fIPDistance); if(fSteelFrontThick>0.) printf(" fSteelFrontThick %6.3f cm \n", fSteelFrontThick); printf(" fNPhi %i | fNZ %i \n", fNPhi, fNZ); printf(" fNCellsInTower %i : fNCellsInSupMod %i : fNCells %i\n",fNCellsInTower, fNCellsInSupMod, fNCells); - if(name.Contains("MAY05")){ + if(fGeoName.Contains("MAY05")){ printf(" fFrontSteelStrip %6.4f cm (thickness of front steel strip)\n", fFrontSteelStrip); printf(" fLateralSteelStrip %6.4f cm (thickness of lateral steel strip)\n", @@ -272,20 +289,20 @@ void AliEMCALGeometry::Init(void){ printf(" fLongModuleSize %6.3f cm \n", fLongModuleSize); printf(" #supermodule in phi direction %i \n", fNPhiSuperModule ); } - if(name.Contains("TRD")) { + if(fGeoName.Contains("TRD")) { printf(" fTrd1Angle %7.4f\n", fTrd1Angle); printf(" f2Trd1Dx2 %7.4f\n", f2Trd1Dx2); - if(name.Contains("TRD2")) { + if(fGeoName.Contains("TRD2")) { printf(" fTrd2AngleY %7.4f\n", fTrd2AngleY); printf(" f2Trd2Dy2 %7.4f\n", f2Trd2Dy2); printf(" fTubsR %7.2f cm\n", fTubsR); printf(" fTubsTurnAngle %7.4f\n", fTubsTurnAngle); printf(" fEmptySpace %7.4f cm\n", fEmptySpace); - } else if(name.Contains("TRD1") && name.Contains("FINAL")){ + } else if(fGeoName.Contains("TRD1") && fGeoName.Contains("FINAL")){ printf("SM dimensions(TRD1) : dx %7.2f dy %7.2f dz %7.2f (SMOD, BOX)\n", fParSM[0],fParSM[1],fParSM[2]); printf(" fPhiGapForSM %7.4f cm \n", fPhiGapForSM); - if(name.Contains("110DEG"))printf(" Last two modules have size 10 degree in phi (180Delete(); delete fArrayOpts; @@ -315,32 +338,54 @@ void AliEMCALGeometry::CheckAditionalOptions() TString addOpt = o->String(); Int_t indj=-1; - for(Int_t j=0; j option |%s| unavailable : ** look to the file AliEMCALGeometry.h **\n", - addOpt.Data()); + AliDebug(2,Form(" option |%s| unavailable : ** look to the file AliEMCALGeometry.h **\n", + addOpt.Data())); assert(0); } else { - printf(" option |%s| is valid : number %i : |%s|\n", - addOpt.Data(), indj, additionalOpts[indj]); + AliDebug(2,Form(" option |%s| is valid : number %i : |%s|\n", + addOpt.Data(), indj, fAdditionalOpts[indj])); if (addOpt.Contains("NL=",TString::kIgnoreCase)) {// number of sampling layers sscanf(addOpt.Data(),"NL=%i", &fNECLayers); - printf(" fNECLayers %i (new) \n", fNECLayers); - } else if(addOpt.Contains("PBTH=",TString::kIgnoreCase)) {//Thickness of the Pb + AliDebug(2,Form(" fNECLayers %i (new) \n", fNECLayers)); + } else if(addOpt.Contains("PBTH=",TString::kIgnoreCase)) {//Thickness of the Pb(fECPbRadThicknes) sscanf(addOpt.Data(),"PBTH=%f", &fECPbRadThickness); - } else if(addOpt.Contains("SCTH=",TString::kIgnoreCase)) {//Thickness of the Sc + } else if(addOpt.Contains("SCTH=",TString::kIgnoreCase)) {//Thickness of the Sc(fECScintThick) sscanf(addOpt.Data(),"SCTH=%f", &fECScintThick); + } else if(addOpt.Contains("LATSS=",TString::kIgnoreCase)) {// Thickness of lateral steel strip (fLateralSteelStrip) + sscanf(addOpt.Data(),"LATSS=%f", &fLateralSteelStrip); + AliDebug(2,Form(" fLateralSteelStrip %f (new) \n", fLateralSteelStrip)); } } } } +void AliEMCALGeometry::DefineSamplingFraction() +{ + // Jun 05,2006 + // Look http://rhic.physics.wayne.edu/~pavlinov/ALICE/SHISHKEBAB/RES/linearityAndResolutionForTRD1.html + // Keep for compatibilty + // + if(fNECLayers == 69) { // 10% layer reduction + fSampling = 12.55; + } else if(fNECLayers == 61) { // 20% layer reduction + fSampling = 12.80; + } else if(fNECLayers == 77) { + if (fECScintThick>0.175 && fECScintThick<0.177) { // 10% Pb thicknes reduction + fSampling = 10.5; // fECScintThick = 0.176, fECPbRadThickness=0.144; + } else if(fECScintThick>0.191 && fECScintThick<0.193) { // 20% Pb thicknes reduction + fSampling = 8.93; // fECScintThick = 0.192, fECPbRadThickness=0.128; + } + } +} + //____________________________________________________________________________ void AliEMCALGeometry::FillTRU(const TClonesArray * digits, TClonesArray * ampmatrix, TClonesArray * timeRmatrix) { @@ -471,7 +516,7 @@ AliEMCALGeometry* AliEMCALGeometry::GetInstance(const Text_t* name, AliEMCALGeometry * rv = 0; if ( fgGeom == 0 ) { if ( strcmp(name,"") == 0 ) rv = 0; - else { + else { fgGeom = new AliEMCALGeometry(name, title); if ( fgInit ) rv = (AliEMCALGeometry * ) fgGeom; else { @@ -481,11 +526,9 @@ AliEMCALGeometry* AliEMCALGeometry::GetInstance(const Text_t* name, } // end if fgInit } // end if strcmp(name,"") }else{ - if ( strcmp(fgGeom->GetName(), name) != 0 ) { - printf("\ncurrent geometry is ") ; - printf(fgGeom->GetName()); - printf("\n you cannot call "); - printf(name); + if ( strcmp(fgGeom->GetName(), name) != 0) { + printf("\ncurrent geometry is %s : ", fgGeom->GetName()); + printf(" you cannot call %s ", name); }else{ rv = (AliEMCALGeometry *) fgGeom; } // end @@ -493,262 +536,8 @@ AliEMCALGeometry* AliEMCALGeometry::GetInstance(const Text_t* name, return rv; } -// These methods are obsolete but use in AliEMCALRecPoint - keep it now -//______________________________________________________________________ -Int_t AliEMCALGeometry::TowerIndex(Int_t ieta,Int_t iphi) const { - // Returns the tower index number from the based on the Z and Phi - // index numbers. - // Inputs: - // Int_t ieta // index along z axis [1-fNZ] - // Int_t iphi // index along phi axis [1-fNPhi] - // Outputs: - // none. - // Returned - // Int_t index // Tower index number - - if ( (ieta <= 0 || ieta>GetNEta()) || - (iphi <= 0 || iphi>GetNPhi())) { - Error("TowerIndex", "Unexpected parameters eta = %d phi = %d!", ieta, iphi) ; - return -1; - } - return ( (iphi - 1)*GetNEta() + ieta ); -} - -//______________________________________________________________________ -void AliEMCALGeometry::TowerIndexes(Int_t index,Int_t &ieta,Int_t &iphi) const { - // Inputs: - // Int_t index // Tower index number [1-fNZ*fNPhi] - // Outputs: - // Int_t ieta // index allong z axis [1-fNZ] - // Int_t iphi // index allong phi axis [1-fNPhi] - // Returned - // none. - - Int_t nindex = 0; - - if ( IsInECA(index) ) { // ECAL index - nindex = index ; - } - else { - Error("TowerIndexes", "Unexpected Id number!") ; - ieta = -1; - iphi = -1; - return; - } - - if (nindex%GetNZ()) - iphi = nindex / GetNZ() + 1 ; - else - iphi = nindex / GetNZ() ; - ieta = nindex - (iphi - 1) * GetNZ() ; - - if (gDebug==2) - printf("TowerIndexes: index=%d,%d, ieta=%d, iphi = %d", index, nindex,ieta, iphi) ; - return; - -} - -//______________________________________________________________________ -void AliEMCALGeometry::EtaPhiFromIndex(Int_t index,Float_t &eta,Float_t &phi) const { - // given the tower index number it returns the based on the eta and phi - // of the tower. - // Inputs: - // Int_t index // Tower index number [1-fNZ*fNPhi] - // Outputs: - // Float_t eta // eta of center of tower in pseudorapidity - // Float_t phi // phi of center of tower in degrees - // Returned - // none. - Int_t ieta, iphi; - Float_t deta, dphi ; - - TowerIndexes(index,ieta,iphi); - - if (gDebug == 2) - printf("EtaPhiFromIndex: index = %d, ieta = %d, iphi = %d", index, ieta, iphi) ; - - deta = (GetArm1EtaMax()-GetArm1EtaMin())/(static_cast(GetNEta())); - eta = GetArm1EtaMin() + ((static_cast(ieta) - 0.5 ))*deta; - - dphi = (GetArm1PhiMax() - GetArm1PhiMin())/(static_cast(GetNPhi())); // in degrees. - phi = GetArm1PhiMin() + dphi*(static_cast(iphi) - 0.5);//iphi range [1-fNphi]. -} - -//______________________________________________________________________ -Int_t AliEMCALGeometry::TowerIndexFromEtaPhi(Float_t eta,Float_t phi) const { - // returns the tower index number based on the eta and phi of the tower. - // Inputs: - // Float_t eta // eta of center of tower in pseudorapidity - // Float_t phi // phi of center of tower in degrees - // Outputs: - // none. - // Returned - // Int_t index // Tower index number [1-fNZ*fNPhi] - - Int_t ieta,iphi; - - ieta = static_cast ( 1 + (static_cast(GetNEta()) * (eta - GetArm1EtaMin()) / (GetArm1EtaMax() - GetArm1EtaMin())) ) ; - - if( ieta <= 0 || ieta > GetNEta() ) { - Error("TowerIndexFromEtaPhi", "Unexpected (eta, phi) = (%f, %f) value, outside of EMCAL!", eta, phi) ; - return -1 ; - } - - iphi = static_cast ( 1 + (static_cast(GetNPhi()) * (phi - GetArm1PhiMin()) / (GetArm1PhiMax() - GetArm1PhiMin())) ) ; - - if( iphi <= 0 || iphi > GetNPhi() ) { - Error("TowerIndexFromEtaPhi", "Unexpected (eta, phi) = (%f, %f) value, outside of EMCAL!", eta, phi) ; - return -1 ; - } - - return TowerIndex(ieta,iphi); -} - -//______________________________________________________________________ -Bool_t AliEMCALGeometry::AbsToRelNumbering(Int_t AbsId, Int_t *relid) const { - // Converts the absolute numbering into the following array/ - // relid[0] = Row number inside EMCAL - // relid[1] = Column number inside EMCAL - // Input: - // Int_t AbsId // Tower index number [1-2*fNZ*fNPhi] - // Outputs: - // Int_t *relid // array of 2. Described above. - Bool_t rv = kTRUE ; - Int_t ieta=0,iphi=0,index=AbsId; - - TowerIndexes(index,ieta,iphi); - relid[0] = ieta; - relid[1] = iphi; - - return rv; -} - -//______________________________________________________________________ -void AliEMCALGeometry::PosInAlice(const Int_t *relid, Float_t &theta, Float_t &phi) const -{ - // Converts the relative numbering into the local EMCAL-module (x, z) - // coordinates - Int_t ieta = relid[0]; // offset along x axis - Int_t iphi = relid[1]; // offset along z axis - Int_t index; - Float_t eta; - - index = TowerIndex(ieta,iphi); - EtaPhiFromIndex(index,eta,phi); - //theta = 180.*(2.0*TMath::ATan(TMath::Exp(-eta)))/TMath::Pi(); - theta = 2.0*TMath::ATan(TMath::Exp(-eta)); - - // correct for distance to IP - Float_t d = GetIP2ECASection() - GetIPDistance() ; - - Float_t correction = 1 + d/GetIPDistance() ; - Float_t tantheta = TMath::Tan(theta) * correction ; - theta = TMath::ATan(tantheta) * TMath::RadToDeg() ; - if (theta < 0 ) - theta += 180. ; - - return; -} - -//______________________________________________________________________ -void AliEMCALGeometry::PosInAlice(Int_t absid, Float_t &theta, Float_t &phi) const -{ - // Converts the relative numbering into the local EMCAL-module (x, z) - // coordinates - Int_t relid[2] ; - AbsToRelNumbering(absid, relid) ; - Int_t ieta = relid[0]; // offset along x axis - Int_t iphi = relid[1]; // offset along z axis - Int_t index; - Float_t eta; - - index = TowerIndex(ieta,iphi); - EtaPhiFromIndex(index,eta,phi); - theta = 2.0*TMath::ATan(TMath::Exp(-eta)) ; - - // correct for distance to IP - Float_t d = 0. ; - if (IsInECA(absid)) - d = GetIP2ECASection() - GetIPDistance() ; - else { - Error("PosInAlice", "Unexpected id # %d!", absid) ; - return; - } - - Float_t correction = 1 + d/GetIPDistance() ; - Float_t tantheta = TMath::Tan(theta) * correction ; - theta = TMath::ATan(tantheta) * TMath::RadToDeg() ; - if (theta < 0 ) - theta += 180. ; - - return; -} - -//______________________________________________________________________ -void AliEMCALGeometry::XYZFromIndex(const Int_t *relid,Float_t &x,Float_t &y, Float_t &z) const { - // given the tower relative number it returns the X, Y and Z - // of the tower. - - // Outputs: - // Float_t x // x of center of tower in cm - // Float_t y // y of center of tower in cm - // Float_t z // z of centre of tower in cm - // Returned - // none. - - Float_t eta,theta, phi,cylradius=0. ; - - Int_t ieta = relid[0]; // offset along x axis - Int_t iphi = relid[1]; // offset along z axis. - Int_t index; - - index = TowerIndex(ieta,iphi); - EtaPhiFromIndex(index,eta,phi); - theta = 180.*(2.0*TMath::ATan(TMath::Exp(-eta)))/TMath::Pi(); - - cylradius = GetIP2ECASection() ; - - Double_t kDeg2Rad = TMath::DegToRad() ; - x = cylradius * TMath::Cos(phi * kDeg2Rad ) ; - y = cylradius * TMath::Sin(phi * kDeg2Rad ) ; - z = cylradius / TMath::Tan(theta * kDeg2Rad ) ; - - return; -} - -//______________________________________________________________________ -void AliEMCALGeometry::XYZFromIndex(Int_t absid, TVector3 &v) const { - // given the tower relative number it returns the X, Y and Z - // of the tower. - - // Outputs: - // Float_t x // x of center of tower in cm - // Float_t y // y of center of tower in cm - // Float_t z // z of centre of tower in cm - // Returned - // none. - - Float_t theta, phi,cylradius=0. ; - - PosInAlice(absid, theta, phi) ; - - if ( IsInECA(absid) ) - cylradius = GetIP2ECASection() ; - else { - Error("XYZFromIndex", "Unexpected Tower section") ; - return; - } - - Double_t kDeg2Rad = TMath::DegToRad() ; - v.SetX(cylradius * TMath::Cos(phi * kDeg2Rad ) ); - v.SetY(cylradius * TMath::Sin(phi * kDeg2Rad ) ); - v.SetZ(cylradius / TMath::Tan(theta * kDeg2Rad ) ) ; - - return; -} - Bool_t AliEMCALGeometry::IsInEMCAL(Double_t x, Double_t y, Double_t z) const { - // Checks whether point is inside the EMCal volume + // Checks whether point is inside the EMCal volume, used in AliEMCALv*.cxx // // Code uses cylindrical approximation made of inner radius (for speed) // @@ -780,24 +569,26 @@ Bool_t AliEMCALGeometry::IsInEMCAL(Double_t x, Double_t y, Double_t z) const { // == Shish-kebab cases == // Int_t AliEMCALGeometry::GetAbsCellId(Int_t nSupMod, Int_t nTower, Int_t nIphi, Int_t nIeta) const -{ // 27-aug-04; +{ + // 27-aug-04; // corr. 21-sep-04; // 13-oct-05; 110 degree case - // 1 <= nSupMod <= fNumberOfSuperModules - // 1 <= nTower <= fNPHI * fNZ ( fNPHI * fNZ/2 for fKey110DEG=1) - // 1 <= nIphi <= fNPHIdiv - // 1 <= nIeta <= fNETAdiv - // 1 <= absid <= fNCells - static Int_t id=0; // have to change from 1 to fNCells - if(fKey110DEG == 1 && nSupMod > 10) { // 110 degree case; last two supermodules - id = fNCellsInSupMod*10 + (fNCellsInSupMod/2)*(nSupMod-11); + // May 31, 2006; ALICE numbering scheme: + // 0 <= nSupMod < fNumberOfSuperModules + // 0 <= nTower < fNPHI * fNZ ( fNPHI * fNZ/2 for fKey110DEG=1) + // 0 <= nIphi < fNPHIdiv + // 0 <= nIeta < fNETAdiv + // 0 <= absid < fNCells + static Int_t id=0; // have to change from 0 to fNCells-1 + if(fKey110DEG == 1 && nSupMod >= 10) { // 110 degree case; last two supermodules + id = fNCellsInSupMod*10 + (fNCellsInSupMod/2)*(nSupMod-10); } else { - id = fNCellsInSupMod*(nSupMod-1); + id = fNCellsInSupMod*nSupMod; } - id += fNCellsInTower *(nTower-1); - id += fNPHIdiv *(nIphi-1); + id += fNCellsInTower *nTower; + id += fNPHIdiv *nIphi; id += nIeta; - if(id<=0 || id > fNCells) { + if(id<0 || id >= fNCells) { // printf(" wrong numerations !!\n"); // printf(" id %6i(will be force to -1)\n", id); // printf(" fNCells %6i\n", fNCells); @@ -805,96 +596,140 @@ Int_t AliEMCALGeometry::GetAbsCellId(Int_t nSupMod, Int_t nTower, Int_t nIphi, I // printf(" nTower %6i\n", nTower); // printf(" nIphi %6i\n", nIphi); // printf(" nIeta %6i\n", nIeta); - id = -TMath::Abs(id); + id = -TMath::Abs(id); // if negative something wrong } return id; } -Bool_t AliEMCALGeometry::CheckAbsCellId(Int_t ind) const -{ // 17-niv-04 - analog of IsInECA - if(name.Contains("TRD")) { - if(ind<=0 || ind > fNCells) return kFALSE; - else return kTRUE; - } else return IsInECA(ind); +Bool_t AliEMCALGeometry::CheckAbsCellId(Int_t absId) const +{ + // May 31, 2006; only trd1 now + if(absId<0 || absId >= fNCells) return kFALSE; + else return kTRUE; } Bool_t AliEMCALGeometry::GetCellIndex(Int_t absId,Int_t &nSupMod,Int_t &nTower,Int_t &nIphi,Int_t &nIeta) const -{ // 21-sep-04 - // 19-oct-05; +{ + // 21-sep-04; 19-oct-05; + // May 31, 2006; ALICE numbering scheme: static Int_t tmp=0, sm10=0; - if(absId<=0 || absId>fNCells) { -// Info("GetCellIndex"," wrong abs Id %i !! \n", absId); - return kFALSE; - } + if(!CheckAbsCellId(absId)) return kFALSE; + sm10 = fNCellsInSupMod*10; - if(fKey110DEG == 1 && absId > sm10) { // 110 degree case; last two supermodules - nSupMod = (absId-1-sm10) / (fNCellsInSupMod/2) + 11; - tmp = (absId-1-sm10) % (fNCellsInSupMod/2); + if(fKey110DEG == 1 && absId >= sm10) { // 110 degree case; last two supermodules + nSupMod = (absId-sm10) / (fNCellsInSupMod/2) + 10; + tmp = (absId-sm10) % (fNCellsInSupMod/2); } else { - nSupMod = (absId-1) / fNCellsInSupMod + 1; - tmp = (absId-1) % fNCellsInSupMod; + nSupMod = absId / fNCellsInSupMod; + tmp = absId % fNCellsInSupMod; } - nTower = tmp / fNCellsInTower + 1; + nTower = tmp / fNCellsInTower; tmp = tmp % fNCellsInTower; - nIphi = tmp / fNPHIdiv + 1; - nIeta = tmp % fNPHIdiv + 1; + nIphi = tmp / fNPHIdiv; + nIeta = tmp % fNPHIdiv; return kTRUE; } -void AliEMCALGeometry::GetTowerPhiEtaIndexInSModule(Int_t nSupMod, Int_t nTower, int &iphit, int &ietat) const -{ // added nSupMod; have to check - 19-oct-05 ! +void AliEMCALGeometry::GetModulePhiEtaIndexInSModule(Int_t nSupMod, Int_t nTower, int &iphim, int &ietam) const +{ + // added nSupMod; have to check - 19-oct-05 ! + // Alice numbering scheme - Jun 01,2006 static Int_t nphi; - if(fKey110DEG == 1 && nSupMod>=11) nphi = fNPhi/2; + if(fKey110DEG == 1 && nSupMod>=10) nphi = fNPhi/2; else nphi = fNPhi; - ietat = (nTower-1)/nphi + 1; // have to change from 1 to fNZ - iphit = (nTower-1)%nphi + 1; // have to change from 1 to fNPhi + ietam = nTower/nphi; // have to change from 0 to fNZ-1 + iphim = nTower%nphi; // have to change from 0 to fNPhi-1 } void AliEMCALGeometry::GetCellPhiEtaIndexInSModule(Int_t nSupMod, Int_t nTower, Int_t nIphi, Int_t nIeta, int &iphi, int &ieta) const -{ // added nSupMod; Nov 25, 05 - static Int_t iphit, ietat; - - GetTowerPhiEtaIndexInSModule(nSupMod,nTower, iphit, ietat); - // have to change from 1 to fNZ*fNETAdiv - ieta = (ietat-1)*fNETAdiv + (3-nIeta); // x(module) = -z(SM) - // iphi - have to change from 1 to fNPhi*fNPHIdiv - iphi = (iphit-1)*fNPHIdiv + nIphi; // y(module) = y(SM) +{ + // added nSupMod; Nov 25, 05 + // Alice numbering scheme - Jun 01,2006 + static Int_t iphim, ietam; + + GetModulePhiEtaIndexInSModule(nSupMod,nTower, iphim, ietam); + // have to change from 0 to (fNZ*fNETAdiv-1) + ieta = ietam*fNETAdiv + (1-nIeta); // x(module) = -z(SM) + // iphi - have to change from 0 to (fNPhi*fNPHIdiv-1) + iphi = iphim*fNPHIdiv + nIphi; // y(module) = y(SM) } Int_t AliEMCALGeometry::GetSuperModuleNumber(Int_t absId) const { + //return the number of the + //supermodule given the absolute + //ALICE numbering + static Int_t nSupMod, nTower, nIphi, nIeta; GetCellIndex(absId, nSupMod, nTower, nIphi, nIeta); return nSupMod; } // Methods for AliEMCALRecPoint - Feb 19, 2006 -Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t &yr, Double_t &zr) +Bool_t AliEMCALGeometry::RelPosCellInSModule(Int_t absId, Double_t &xr, Double_t &yr, Double_t &zr) const { + // Look to see what the relative + // position inside a given cell is + // for a recpoint. + // Alice numbering scheme - Jun 08, 2006 + static Int_t nSupMod, nTower, nIphi, nIeta, iphi, ieta; + static Int_t phiIndexShift=6; if(!CheckAbsCellId(absId)) return kFALSE; GetCellIndex(absId, nSupMod, nTower, nIphi, nIeta); GetCellPhiEtaIndexInSModule(nSupMod,nTower,nIphi,nIeta, iphi, ieta); - xr = fXCentersOfCells->At(ieta-1); - zr = fEtaCentersOfCells->At(ieta-1); + xr = fXCentersOfCells.At(ieta); + zr = fEtaCentersOfCells.At(ieta); - yr = fPhiCentersOfCells->At(iphi-1); + if(nSupMod<10) { + yr = fPhiCentersOfCells.At(iphi); + } else { + yr = fPhiCentersOfCells.At(iphi + phiIndexShift); + cout<<" absId "<Add(mod); } } else { - cout<<" Already exits : "; + AliDebug(2,Form(" Already exits : ")); } - cout<<" fShishKebabTrd1Modules "<< fShishKebabTrd1Modules << " has " - << fShishKebabTrd1Modules->GetSize() << " modules" <GetSize())); // Feb 20,2006; + // Jun 01, 2006 - ALICE numbering scheme // define grid for cells in eta(z) and x directions in local coordinates system of SM - fEtaCentersOfCells = new TArrayD(fNZ *fNETAdiv); - fXCentersOfCells = new TArrayD(fNZ *fNETAdiv); - printf(" Cells grid in eta directions : size %i\n", fEtaCentersOfCells->GetSize()); + // fEtaCentersOfCells = new TArrayD(fNZ *fNETAdiv); + // fXCentersOfCells = new TArrayD(fNZ *fNETAdiv); + fEtaCentersOfCells.Set(fNZ *fNETAdiv); + fXCentersOfCells.Set(fNZ *fNETAdiv); + AliDebug(2,Form(" Cells grid in eta directions : size %i\n", fEtaCentersOfCells.GetSize())); Int_t iphi=0, ieta=0, nTower=0; Double_t xr, zr; for(Int_t it=0; itGetCenterOfCellInLocalCoordinateofSM(ic+1, xr, zr); - GetCellPhiEtaIndexInSModule(1, nTower, 1, ic+1, iphi, ieta); // don't depend from phi - fXCentersOfCells->AddAt(float(xr) - fParSM[0],ieta-1); - fEtaCentersOfCells->AddAt(float(zr) - fParSM[2],ieta-1); + trd1->GetCenterOfCellInLocalCoordinateofSM(ic, xr, zr); + GetCellPhiEtaIndexInSModule(0, nTower, 0, ic, iphi, ieta); // don't depend from phi - ieta in action + fXCentersOfCells.AddAt(float(xr) - fParSM[0],ieta); + fEtaCentersOfCells.AddAt(float(zr) - fParSM[2],ieta); } } - for(Int_t i=0; iGetSize(); i++) { - printf(" ind %2.2i : z %8.3f : x %8.3f", i+1, fEtaCentersOfCells->At(i),fXCentersOfCells->At(i)); - if(i%2 != 0) printf("\n"); + for(Int_t i=0; iGetSize()); + // fPhiCentersOfCells = new TArrayD(fNPhi*fNPHIdiv); + fPhiCentersOfCells.Set(fNPhi*fNPHIdiv); + AliDebug(2,Form(" Cells grid in phi directions : size %i\n", fPhiCentersOfCells.GetSize())); Int_t ind=0; for(Int_t it=0; itAddAt(ytLeftCenterCell,ind); - printf(" ind %2.2i : y %8.3f ", ind, fPhiCentersOfCells->At(ind)); ind++; - if(ic == fNPHIdiv-1) printf("\n"); + fPhiCentersOfCells.AddAt(ytLeftCenterCell,ind); + AliDebug(2,Form(" ind %2.2i : y %8.3f ", ind, fPhiCentersOfCells.At(ind))); + ind++; } } - printf("\n"); } void AliEMCALGeometry::GetTransformationForSM() { + //Uses the geometry manager to + //load the transformation matrix + //for the supermodules + static Bool_t transInit=kFALSE; if(transInit) return; @@ -961,23 +803,23 @@ void AliEMCALGeometry::GetTransformationForSM() assert(0); } TGeoNode *tn = gGeoManager->GetTopNode(); - TGeoNode *node=0, *XEN1 = 0; + TGeoNode *node=0, *xen1 = 0; for(i=0; iGetNdaughters(); i++) { node = tn->GetDaughter(i); TString ns(node->GetName()); if(ns.Contains(GetNameOfEMCALEnvelope())) { - XEN1 = node; + xen1 = node; break; } } - if(!XEN1) { + if(!xen1) { Info("CreateTransformationForSM() "," geometry has not EMCAL envelope with name %s", GetNameOfEMCALEnvelope()); assert(0); } - printf(" i %i : EMCAL Envelope is %s : #SM %i \n", i, XEN1->GetName(), XEN1->GetNdaughters()); - for(i=0; iGetNdaughters(); i++) { - TGeoNodeMatrix *sm = (TGeoNodeMatrix*)XEN1->GetDaughter(i); + printf(" i %i : EMCAL Envelope is %s : #SM %i \n", i, xen1->GetName(), xen1->GetNdaughters()); + for(i=0; iGetNdaughters(); i++) { + TGeoNodeMatrix *sm = (TGeoNodeMatrix*)xen1->GetDaughter(i); fMatrixOfSM[i] = sm->GetMatrix(); //Compiler doesn't like this syntax... // printf(" %i : matrix %x \n", i, fMatrixOfSM[i]); @@ -985,39 +827,92 @@ void AliEMCALGeometry::GetTransformationForSM() transInit = kTRUE; } -void AliEMCALGeometry::GetGlobal(const Double_t *loc, Double_t *glob, int nsm) const +void AliEMCALGeometry::GetGlobal(const Double_t *loc, Double_t *glob, int ind) const { + // Figure out the global numbering + // of a given supermodule from the + // local numbering + // Alice numbering - Jun 03,2006 // if(fMatrixOfSM[0] == 0) GetTransformationForSM(); - static int ind; - ind = nsm-1; + if(ind>=0 && ind < GetNumberOfSuperModules()) { fMatrixOfSM[ind]->LocalToMaster(loc, glob); } } -void AliEMCALGeometry::GetGlobal(Int_t /* absId */, TVector3 & /* vglob */) const -{ // have to be defined -} - -void AliEMCALGeometry::GetGlobal(const TVector3 &vloc, TVector3 &vglob, int nsm) const +void AliEMCALGeometry::GetGlobal(const TVector3 &vloc, TVector3 &vglob, int ind) const { + //Figure out the global numbering + //of a given supermodule from the + //local numbering given a 3-vector location + static Double_t tglob[3], tloc[3]; vloc.GetXYZ(tloc); - GetGlobal(tloc, tglob, nsm); + GetGlobal(tloc, tglob, ind); vglob.SetXYZ(tglob[0], tglob[1], tglob[2]); } +void AliEMCALGeometry::GetGlobal(Int_t absId , double glob[3]) const +{ + // Alice numbering scheme - Jun 03, 2006 + static Int_t nSupMod, nModule, nIphi, nIeta; + static double loc[3]; + + glob[0]=glob[1]=glob[2]=0.0; // bad case + if(RelPosCellInSModule(absId, loc)) { + GetCellIndex(absId, nSupMod, nModule, nIphi, nIeta); + fMatrixOfSM[nSupMod]->LocalToMaster(loc, glob); + } +} + +void AliEMCALGeometry::GetGlobal(Int_t absId , TVector3 &vglob) const +{ + // Alice numbering scheme - Jun 03, 2006 + static Double_t glob[3]; + + GetGlobal(absId, glob); + vglob.SetXYZ(glob[0], glob[1], glob[2]); + +} + void AliEMCALGeometry::GetGlobal(const AliRecPoint *rp, TVector3 &vglob) const { + // Figure out the global numbering + // of a given supermodule from the + // local numbering for RecPoints + static TVector3 vloc; - static Int_t nSupMod, nTower, nIphi, nIeta; + static Int_t nSupMod, nModule, nIphi, nIeta; AliRecPoint *rpTmp = (AliRecPoint*)rp; // const_cast ?? if(!rpTmp) return; AliEMCALRecPoint *rpEmc = (AliEMCALRecPoint*)rpTmp; - GetCellIndex(rpEmc->GetAbsId(0), nSupMod, nTower, nIphi, nIeta); + GetCellIndex(rpEmc->GetAbsId(0), nSupMod, nModule, nIphi, nIeta); rpTmp->GetLocalPosition(vloc); GetGlobal(vloc, vglob, nSupMod); } +void AliEMCALGeometry::EtaPhiFromIndex(Int_t absId,Float_t &eta,Float_t &phi) const +{ + // Jun 03, 2006 - version for TRD1 + static TVector3 vglob; + GetGlobal(absId, vglob); + eta = vglob.Eta(); + phi = vglob.Phi(); +} + +AliEMCALShishKebabTrd1Module* AliEMCALGeometry::GetShishKebabModule(Int_t neta=0) +{ + //This method was too long to be + //included in the header file - the + //rule checker complained about it's + //length, so we move it here. It returns the + //shishkebabmodule at a given eta index point. + + static AliEMCALShishKebabTrd1Module* trd1=0; + if(fShishKebabTrd1Modules && neta>=0 && netaGetSize()) { + trd1 = (AliEMCALShishKebabTrd1Module*)fShishKebabTrd1Modules->At(neta); + } else trd1 = 0; + return trd1; +}