X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=PHOS%2FAliPHOSGeometry.cxx;h=ecfd9403f2d710da0608cd839be72d8faea823d8;hb=34eaaa7ef77cf5c1e497d5a2218baa7163448bc9;hp=32b883a11ff4c9feebe91a073841296a2d9f9513;hpb=e5b7b511079ed49862daca7e4ceb7d7f8ec4d392;p=u%2Fmrichter%2FAliRoot.git diff --git a/PHOS/AliPHOSGeometry.cxx b/PHOS/AliPHOSGeometry.cxx index 32b883a11ff..ecfd9403f2d 100644 --- a/PHOS/AliPHOSGeometry.cxx +++ b/PHOS/AliPHOSGeometry.cxx @@ -139,7 +139,17 @@ void AliPHOSGeometry::Init(void) fPHOSParams[3] = emcParams[3] + fGeometryCPV->GetCPVBoxSize(1)/2. ; fIPtoUpperCPVsurface = fGeometryEMCA->GetIPtoOuterCoverDistance() - fGeometryCPV->GetCPVBoxSize(1) ; - + + //calculate offset to crystal surface + Float_t * inthermo = fGeometryEMCA->GetInnerThermoHalfSize() ; + Float_t * strip = fGeometryEMCA->GetStripHalfSize() ; + Float_t* splate = fGeometryEMCA->GetSupportPlateHalfSize(); + Float_t * crystal = fGeometryEMCA->GetCrystalHalfSize() ; + Float_t * pin = fGeometryEMCA->GetAPDHalfSize() ; + 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() @@ -228,7 +238,7 @@ void AliPHOSGeometry::SetPHOSAngles() } //____________________________________________________________________________ -Bool_t AliPHOSGeometry::AbsToRelNumbering(Int_t AbsId, Int_t * relid) const +Bool_t AliPHOSGeometry::AbsToRelNumbering(Int_t absId, Int_t * relid) const { // Converts the absolute numbering into the following array // relid[0] = PHOS Module number 1:fNModules @@ -238,7 +248,7 @@ Bool_t AliPHOSGeometry::AbsToRelNumbering(Int_t AbsId, Int_t * relid) const // relid[3] = Column number inside a PHOS module Bool_t rv = kTRUE ; - Float_t id = AbsId ; + Float_t id = absId ; Int_t phosmodulenumber = (Int_t)TMath:: Ceil( id / GetNCristalsInModule() ) ; @@ -262,13 +272,18 @@ Bool_t AliPHOSGeometry::AbsToRelNumbering(Int_t AbsId, Int_t * relid) const } return rv ; } - //____________________________________________________________________________ void AliPHOSGeometry::GetGlobal(const AliRecPoint* recPoint, TVector3 & gpos) const +{ + AliFatal(Form("Please use GetGlobalPHOS(recPoint,gpos) instead of GetGlobal!")); +} + +//____________________________________________________________________________ +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) ; @@ -281,34 +296,29 @@ void AliPHOSGeometry::GetGlobal(const AliRecPoint* recPoint, TVector3 & gpos) co Double_t dy ; if(tmpPHOS->IsEmc()){ sprintf(path,"/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1",tmpPHOS->GetPHOSMod()) ; - //calculate offset to crystal surface - Float_t * inthermo = fGeometryEMCA->GetInnerThermoHalfSize() ; - Float_t * strip = fGeometryEMCA->GetStripHalfSize() ; - Float_t* splate = fGeometryEMCA->GetSupportPlateHalfSize(); - Float_t * crystal = fGeometryEMCA->GetCrystalHalfSize() ; - Float_t * pin = fGeometryEMCA->GetAPDHalfSize() ; - Float_t * preamp = fGeometryEMCA->GetPreampHalfSize() ; - dy=-inthermo[1]+strip[1]+splate[1]+crystal[1]-fGeometryEMCA->GetAirGapLed()/2.+pin[1]+preamp[1] ; +// sprintf(path,"/ALIC_1/PHOS_%d",tmpPHOS->GetPHOSMod()) ; + dy=fCrystalShift ; } else{ sprintf(path,"/ALIC_1/PHOS_%d/PCPV_1",tmpPHOS->GetPHOSMod()) ; dy= GetCPVBoxSize(1)/2. ; //center of CPV module } 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); + if (m){ + m->LocalToMaster(pos,posC); + } else{ AliFatal("Geo matrixes are not loaded \n") ; } - if(tmpPHOS->IsEmc()) - gpos.SetXYZ(posC[0],posC[1],-posC[2]) ; - else - gpos.SetXYZ(posC[0],posC[1],posC[2]) ; + gpos.SetXYZ(posC[0],posC[1],posC[2]) ; } //____________________________________________________________________________ @@ -318,24 +328,15 @@ void AliPHOSGeometry::ImpactOnEmc(Double_t * vtx, Double_t theta, Double_t phi, // calculates the impact coordinates on PHOS of a neutral particle // emitted in the vertex vtx[3] with direction theta and phi in the ALICE global coordinate system TVector3 p(TMath::Sin(theta)*TMath::Cos(phi),TMath::Sin(theta)*TMath::Sin(phi),TMath::Cos(theta)) ; - TVector3 v(vtx[0],vtx[1],-vtx[2]) ; + TVector3 v(vtx[0],vtx[1],vtx[2]) ; - //calculate offset to crystal surface - Float_t * inthermo = fGeometryEMCA->GetInnerThermoHalfSize() ; - Float_t * strip = fGeometryEMCA->GetStripHalfSize() ; - Float_t* splate = fGeometryEMCA->GetSupportPlateHalfSize(); - Float_t * crystal = fGeometryEMCA->GetCrystalHalfSize() ; - Float_t * pin = fGeometryEMCA->GetAPDHalfSize() ; - Float_t * preamp = fGeometryEMCA->GetPreampHalfSize() ; - Float_t dy=-inthermo[1]+strip[1]+splate[1]+crystal[1]-fGeometryEMCA->GetAirGapLed()/2.+pin[1]+preamp[1] ; - if (!gGeoManager){ AliFatal("Geo manager not initialized\n"); } for(Int_t imod=1; imod<=GetNModules() ; imod++){ //create vector from (0,0,0) to center of crystal surface of imod module - Double_t tmp[3]={0.,-dy,0.} ; + Double_t tmp[3]={0.,-fCrystalShift,0.} ; char path[100] ; sprintf(path,"/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1",imod) ; @@ -382,25 +383,25 @@ Bool_t AliPHOSGeometry::Impact(const TParticle * particle) const } //____________________________________________________________________________ -Bool_t AliPHOSGeometry::RelToAbsNumbering(const Int_t * relid, Int_t & AbsId) const +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 + // absId = from 1 to fNModules * fNPhi * fNZ // CPV pad: - // AbsId = from N(total PHOS crystals) + 1 + // 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 = + 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 { // it is a CPV pad - AbsId = GetNPhi() * GetNZ() * GetNModules() // the offset to separate EMCA crystals from CPV pads + 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 @@ -426,14 +427,17 @@ void AliPHOSGeometry::RelPosInAlice(Int_t id, TVector3 & pos ) const char path[100] ; if(relid[1]==0){ //this is EMC - Float_t * xls = fGeometryEMCA->GetCrystalHalfSize() ; - Double_t ps[3]= {0.0,-xls[1],0.}; //Position incide the crystal + Double_t ps[3]= {0.0,-fCryCellShift,0.}; //Position incide the crystal Double_t psC[3]={0.0,0.0,0.}; //Global position //Shift and possibly apply misalignment corrections - Int_t nCellsInStrip=fGeometryEMCA->GetNCellsXInStrip() ; - Int_t strip=1+(relid[3]-1)/fGeometryEMCA->GetNCellsZInStrip()+((relid[2]-1)/nCellsInStrip)*fGeometryEMCA->GetNStripZ() ; - Int_t cell=fGeometryEMCA->GetNCellsZInStrip()*(1+(relid[2]-1)%nCellsInStrip)-(relid[3]-1)%fGeometryEMCA->GetNCellsZInStrip() ; + Int_t nCellsXInStrip=fGeometryEMCA->GetNCellsXInStrip() ; + Int_t nCellsZInStrip=fGeometryEMCA->GetNCellsZInStrip() ; + Int_t strip=1+((Int_t) TMath::Ceil((Double_t)relid[2]/nCellsXInStrip))*fGeometryEMCA->GetNStripZ()- + (Int_t) TMath::Ceil((Double_t)relid[3]/nCellsZInStrip) ; + Int_t cellraw= relid[3]%nCellsZInStrip ; + if(cellraw==0)cellraw=nCellsZInStrip ; + Int_t cell= ((relid[2]-1)%nCellsXInStrip)*nCellsZInStrip + cellraw ; sprintf(path,"/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1/PSTR_%d/PCEL_%d", relid[0],strip,cell) ; if (!gGeoManager->cd(path)){ @@ -444,7 +448,7 @@ void AliPHOSGeometry::RelPosInAlice(Int_t id, TVector3 & pos ) const else{ AliFatal("Geo matrixes are not loaded \n") ; } - pos.SetXYZ(psC[0],psC[1],-psC[2]) ; + pos.SetXYZ(psC[0],psC[1],psC[2]) ; } else{ //first calculate position with respect to CPV plain @@ -470,7 +474,7 @@ void AliPHOSGeometry::RelPosInAlice(Int_t id, TVector3 & pos ) const } //____________________________________________________________________________ -void AliPHOSGeometry::RelPosToAbsId(Int_t module, Double_t x, Double_t z, Int_t & AbsId) const +void AliPHOSGeometry::RelPosToAbsId(Int_t module, Double_t x, Double_t z, Int_t & absId) const { // converts local PHOS-module (x, z) coordinates to absId @@ -478,20 +482,23 @@ void AliPHOSGeometry::RelPosToAbsId(Int_t module, Double_t x, Double_t z, Int_t if (!gGeoManager){ AliFatal("Geo manager not initialized\n"); } - Int_t relid[4] ; + Double_t posL[3]={x,-fCrystalShift,-z} ; //Only for EMC!!! + Double_t posG[3] ; char path[100] ; sprintf(path,"/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1",module) ; - Double_t lpos[3]={x,0.0,z} ; - Double_t gpos[3]={0.,0.,0.} ; if (!gGeoManager->cd(path)){ - AliFatal("Geo manager can not find path \n"); + AliFatal("Geo manager can not find path \n"); + } + TGeoHMatrix *mPHOS = gGeoManager->GetCurrentMatrix(); + if (mPHOS){ + mPHOS->LocalToMaster(posL,posG); } - TGeoHMatrix *m = gGeoManager->GetCurrentMatrix(); - if (m) m->LocalToMaster(lpos,gpos); else{ AliFatal("Geo matrixes are not loaded \n") ; } - gGeoManager->FindNode(gpos[0],gpos[1],gpos[2]) ; + + Int_t relid[4] ; + gGeoManager->FindNode(posG[0],posG[1],posG[2]) ; //Check that path contains PSTR and extract strip number TString cpath(gGeoManager->GetPath()) ; Int_t indx = cpath.Index("PCEL") ; @@ -504,6 +511,7 @@ void AliPHOSGeometry::RelPosToAbsId(Int_t module, Double_t x, Double_t z, Int_t if(relid[3]<1)relid[3]=1 ; if(relid[2]>GetNPhi())relid[2]=GetNPhi() ; if(relid[3]>GetNZ())relid[3]=GetNZ() ; + RelToAbsNumbering(relid,absId) ; } else{ Int_t indx2 = cpath.Index("/",indx) ; @@ -515,18 +523,16 @@ void AliPHOSGeometry::RelPosToAbsId(Int_t module, Double_t x, Double_t z, Int_t indx2 = cpath.Index("/",indx) ; TString strip=cpath(indx+5,indx2-indx-5) ; Int_t iStrip = strip.Atoi() ; - relid[0] = module ; - relid[1] = 0 ; - Int_t raw = fGeometryEMCA->GetNCellsXInStrip()*((iStrip-1)/fGeometryEMCA->GetNStripZ()) + - 1 + (icell-1)/fGeometryEMCA->GetNCellsZInStrip() ; - Int_t col = fGeometryEMCA->GetNCellsZInStrip()*(1+(iStrip-1)%fGeometryEMCA->GetNStripZ()) - - (icell-1)%fGeometryEMCA->GetNCellsZInStrip() ; - if(col==0) col=GetNZ() ; - relid[2] = raw ; - relid[3] = col ; + + Int_t row = fGeometryEMCA->GetNStripZ() - (iStrip - 1) % (fGeometryEMCA->GetNStripZ()) ; + Int_t col = (Int_t) TMath::Ceil((Double_t) iStrip/(fGeometryEMCA->GetNStripZ())) -1 ; + + // Absid for 8x2-strips. Looks nice :) + absId = (module-1)*GetNCristalsInModule() + + row * 2 + (col*fGeometryEMCA->GetNCellsXInStrip() + (icell - 1) / 2)*GetNZ() - (icell & 1 ? 1 : 0); + } - RelToAbsNumbering(relid,AbsId) ; } //____________________________________________________________________________ @@ -547,14 +553,18 @@ void AliPHOSGeometry::RelPosInModule(const Int_t * relid, Float_t & x, Float_t & // x = - ( GetNPhi()/2. - relid[2] + 0.5 ) * GetCellStep() ; // position of Xtal with respect // z = - ( GetNZ() /2. - relid[3] + 0.5 ) * GetCellStep() ; // of center of PHOS module - Double_t pos[3]= {0.0,0.0,0.}; //Position incide the crystal (Y coordinalte is not important) + Double_t pos[3]= {0.0,-fCryCellShift,0.}; //Position incide the crystal Double_t posC[3]={0.0,0.0,0.}; //Global position //Shift and possibly apply misalignment corrections - Int_t nCellsInStrip=fGeometryEMCA->GetNCellsXInStrip() ; - Int_t strip=1+(relid[3]-1)/fGeometryEMCA->GetNCellsZInStrip()+((relid[2]-1)/nCellsInStrip)*fGeometryEMCA->GetNStripZ() ; - Int_t sellRaw = (relid[2]-1)%nCellsInStrip ; - Int_t cell=fGeometryEMCA->GetNCellsZInStrip()*(sellRaw+1) -(relid[3]-1)%fGeometryEMCA->GetNCellsZInStrip() ; + Int_t nCellsXInStrip=fGeometryEMCA->GetNCellsXInStrip() ; + Int_t nCellsZInStrip=fGeometryEMCA->GetNCellsZInStrip() ; +// Int_t strip=1+(relid[3]-1)/fGeometryEMCA->GetNCellsZInStrip()+((relid[2]-1)/nCellsInStrip)*fGeometryEMCA->GetNStripZ() ; + Int_t strip=1+((Int_t) TMath::Ceil((Double_t)relid[2]/nCellsXInStrip))*fGeometryEMCA->GetNStripZ()- + (Int_t) TMath::Ceil((Double_t)relid[3]/nCellsZInStrip) ; + Int_t cellraw= relid[3]%nCellsZInStrip ; + if(cellraw==0)cellraw=nCellsZInStrip ; + Int_t cell= ((relid[2]-1)%nCellsXInStrip)*nCellsZInStrip + cellraw ; sprintf(path,"/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1/PSTR_%d/PCEL_%d", relid[0],strip,cell) ; if (!gGeoManager->cd(path)){ @@ -576,10 +586,10 @@ void AliPHOSGeometry::RelPosInModule(const Int_t * relid, Float_t & x, Float_t & else{ AliFatal("Geo matrixes are not loaded \n") ; } +//printf("RelPosInMod: posL=[%f,%f,%f]\n",posL[0],posL[1],posL[2]) ; //printf("old: x=%f, z=%f \n",x,z); x=posL[0] ; - z=-posL[1]; -//printf("Geo: x=%f, z=%f \n",x,z); + z=posL[1]; return ; } else{//CPV @@ -648,7 +658,7 @@ void AliPHOSGeometry::Global2Local(TVector3& localPosition, { // Transforms a global position of the rec.point to the local coordinate system //Return to PHOS local system - Double_t posG[3]={globalPosition.X(),globalPosition.Y(),-globalPosition.Z()} ; + Double_t posG[3]={globalPosition.X(),globalPosition.Y(),globalPosition.Z()} ; Double_t posL[3]={0.,0.,0.} ; char path[100] ; sprintf(path,"/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1",module) ; @@ -661,7 +671,7 @@ void AliPHOSGeometry::Global2Local(TVector3& localPosition, else{ AliFatal("Geo matrixes are not loaded \n") ; } - localPosition.SetXYZ(posL[0],posL[1],posL[2]) ; + localPosition.SetXYZ(posL[0],posL[1]+fCrystalShift,-posL[2]) ; /* Float_t angle = GetPHOSAngle(module); // (40,20,0,-20,-40) degrees @@ -677,31 +687,25 @@ void AliPHOSGeometry::Local2Global(Int_t mod, Float_t x, Float_t z, { char path[100] ; sprintf(path,"/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1",mod) ; +// sprintf(path,"/ALIC_1/PHOS_%d",mod) ; if (!gGeoManager->cd(path)){ AliFatal("Geo manager can not find path \n"); } - //calculate offset to crystal surface - Float_t * inthermo = fGeometryEMCA->GetInnerThermoHalfSize() ; - Float_t * strip = fGeometryEMCA->GetStripHalfSize() ; - Float_t* splate = fGeometryEMCA->GetSupportPlateHalfSize(); - Float_t * crystal = fGeometryEMCA->GetCrystalHalfSize() ; - Float_t * pin = fGeometryEMCA->GetAPDHalfSize() ; - Float_t * preamp = fGeometryEMCA->GetPreampHalfSize() ; - Float_t dy=-inthermo[1]+strip[1]+splate[1]+crystal[1]-fGeometryEMCA->GetAirGapLed()/2.+pin[1]+preamp[1] ; - Double_t posL[3]={x,-dy,z} ; + Double_t posL[3]={x,-fCrystalShift,-z} ; //Only for EMC!!! Double_t posG[3] ; TGeoHMatrix *mPHOS = gGeoManager->GetCurrentMatrix(); - if (mPHOS) mPHOS->LocalToMaster(posL,posG); + if (mPHOS){ + mPHOS->LocalToMaster(posL,posG); + } else{ AliFatal("Geo matrixes are not loaded \n") ; } - globalPosition.SetXYZ(posG[0],posG[1],-posG[2]) ; - - + globalPosition.SetXYZ(posG[0],posG[1],posG[2]) ; } //____________________________________________________________________________ void AliPHOSGeometry::GetIncidentVector(const TVector3 &vtx, Int_t module, Float_t x,Float_t z, TVector3 &vInc) const { //Calculates vector pointing from vertex to current poisition in module local frame + //Note that PHOS local system and ALICE global have opposite z directions Global2Local(vInc,vtx,module) ; vInc.SetXYZ(vInc.X()+x,vInc.Y(),vInc.Z()+z) ;