0fae4bfa94c6e4aa879b42b97183d5f64df1f89a
[u/mrichter/AliRoot.git] / PHOS / AliPHOSGeometry.cxx
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 /* $Id$ */
17
18 //_________________________________________________________________________
19 // Geometry class  for PHOS : singleton  
20 // PHOS consists of the electromagnetic calorimeter (EMCA)
21 // and a charged particle veto either in the Subatech's version (PPSD)
22 // or in the IHEP's one (CPV).
23 // The EMCA/PPSD/CPV modules are parametrized so that any configuration
24 // can be easily implemented 
25 // The title is used to identify the version of CPV used.
26 //                  
27 // -- Author: Yves Schutz (SUBATECH) & Dmitri Peressounko (RRC "KI" & SUBATECH)
28
29 // --- ROOT system ---
30
31 #include "TVector3.h"
32 #include "TRotation.h" 
33 #include "TParticle.h"
34 #include <TGeoManager.h>
35 #include <TGeoMatrix.h>
36
37 // --- Standard library ---
38
39 // --- AliRoot header files ---
40 #include "AliLog.h"
41 #include "AliPHOSGeometry.h"
42 #include "AliPHOSEMCAGeometry.h" 
43 #include "AliPHOSRecPoint.h"
44
45 ClassImp(AliPHOSGeometry)
46
47 // these initialisations are needed for a singleton
48 AliPHOSGeometry  * AliPHOSGeometry::fgGeom = 0 ;
49 Bool_t             AliPHOSGeometry::fgInit = kFALSE ;
50
51 //____________________________________________________________________________
52 AliPHOSGeometry::AliPHOSGeometry() : 
53                     AliPHOSGeoUtils(),
54                     fAngle(0.f),
55                     fPHOSAngle(0),
56                     fIPtoUpperCPVsurface(0),
57                     fCrystalShift(0),
58                     fCryCellShift(0),
59                     fRotMatrixArray(0)
60 {
61     // default ctor 
62     // must be kept public for root persistency purposes, but should never be called by the outside world
63     fgGeom          = 0 ;
64
65     fPHOSParams[0] = 0.;
66     fPHOSParams[1] = 0.;
67     fPHOSParams[2] = 0.;
68     fPHOSParams[3] = 0.;
69 }  
70
71 //____________________________________________________________________________
72 AliPHOSGeometry::AliPHOSGeometry(const AliPHOSGeometry & rhs)
73                     : AliPHOSGeoUtils(rhs),
74                       fAngle(rhs.fAngle),
75                       fPHOSAngle(0),
76                       fIPtoUpperCPVsurface(rhs.fIPtoUpperCPVsurface),
77                       fCrystalShift(rhs.fCrystalShift),
78                       fCryCellShift(rhs.fCryCellShift),
79                       fRotMatrixArray(0)
80 {
81   Fatal("cpy ctor", "not implemented") ; 
82 }
83
84 //____________________________________________________________________________
85 AliPHOSGeometry::AliPHOSGeometry(const Text_t* name, const Text_t* title) 
86                   : AliPHOSGeoUtils(name, title),
87                     fAngle(0.f),
88                     fPHOSAngle(0),
89                     fIPtoUpperCPVsurface(0),
90                     fCrystalShift(0),
91                     fCryCellShift(0),
92                     fRotMatrixArray(0)
93
94   // ctor only for internal usage (singleton)
95   Init() ; 
96   fgGeom = this;
97 }
98
99 //____________________________________________________________________________
100 AliPHOSGeometry::~AliPHOSGeometry(void)
101 {
102   // dtor
103
104   if (fRotMatrixArray) fRotMatrixArray->Delete() ; 
105   if (fRotMatrixArray) delete fRotMatrixArray ; 
106   if (fPHOSAngle     ) delete[] fPHOSAngle ; 
107 }
108
109 //____________________________________________________________________________
110 void AliPHOSGeometry::Init(void)
111 {
112   // Initializes the PHOS parameters :
113   //  IHEP is the Protvino CPV (cathode pad chambers)
114   
115   fgInit     = kTRUE ; 
116
117   fAngle        = 20;
118
119   
120   fPHOSAngle = new Float_t[fNModules] ;
121   
122   const Float_t * emcParams = fGeometryEMCA->GetEMCParams() ;
123   
124   fPHOSParams[0] =  TMath::Max((Double_t)fGeometryCPV->GetCPVBoxSize(0)/2., 
125                                (Double_t)(emcParams[0] - (emcParams[1]-emcParams[0])*
126                                           fGeometryCPV->GetCPVBoxSize(1)/2/emcParams[3]));
127   fPHOSParams[1] = emcParams[1] ;
128   fPHOSParams[2] = TMath::Max((Double_t)emcParams[2], (Double_t)fGeometryCPV->GetCPVBoxSize(2)/2.);
129   fPHOSParams[3] = emcParams[3] + fGeometryCPV->GetCPVBoxSize(1)/2. ;
130   
131   fIPtoUpperCPVsurface = fGeometryEMCA->GetIPtoOuterCoverDistance() - fGeometryCPV->GetCPVBoxSize(1) ;
132
133   //calculate offset to crystal surface
134   const Float_t * inthermo = fGeometryEMCA->GetInnerThermoHalfSize() ;
135   const Float_t * strip = fGeometryEMCA->GetStripHalfSize() ;
136   const Float_t * splate = fGeometryEMCA->GetSupportPlateHalfSize();
137   const Float_t * crystal = fGeometryEMCA->GetCrystalHalfSize() ;
138   const Float_t * pin = fGeometryEMCA->GetAPDHalfSize() ;
139   const Float_t * preamp = fGeometryEMCA->GetPreampHalfSize() ;
140   fCrystalShift=-inthermo[1]+strip[1]+splate[1]+crystal[1]-fGeometryEMCA->GetAirGapLed()/2.+pin[1]+preamp[1] ;
141   fCryCellShift=crystal[1]-(fGeometryEMCA->GetAirGapLed()-2*pin[1]-2*preamp[1])/2;
142  
143   Int_t index ;
144   for ( index = 0; index < fNModules; index++ )
145     fPHOSAngle[index] = 0.0 ; // Module position angles are set in CreateGeometry()
146   
147   fRotMatrixArray = new TObjArray(fNModules) ; 
148
149   // Geometry parameters are calculated
150
151   SetPHOSAngles();
152   Double_t const kRADDEG = 180.0 / TMath::Pi() ;
153   Float_t r = GetIPtoOuterCoverDistance() + fPHOSParams[3] - GetCPVBoxSize(1) ;
154   for (Int_t iModule=0; iModule<fNModules; iModule++) {
155     fModuleCenter[iModule][0] = r * TMath::Sin(fPHOSAngle[iModule] / kRADDEG );
156     fModuleCenter[iModule][1] =-r * TMath::Cos(fPHOSAngle[iModule] / kRADDEG );
157     fModuleCenter[iModule][2] = 0.;
158     
159     fModuleAngle[iModule][0][0] =  90;
160     fModuleAngle[iModule][0][1] =   fPHOSAngle[iModule];
161     fModuleAngle[iModule][1][0] =   0;
162     fModuleAngle[iModule][1][1] =   0;
163     fModuleAngle[iModule][2][0] =  90;
164     fModuleAngle[iModule][2][1] = 270 + fPHOSAngle[iModule];
165   }
166
167 }
168
169 //____________________________________________________________________________
170 AliPHOSGeometry *  AliPHOSGeometry::GetInstance() 
171
172   // Returns the pointer of the unique instance; singleton specific
173   
174   return static_cast<AliPHOSGeometry *>( fgGeom ) ; 
175 }
176
177 //____________________________________________________________________________
178 AliPHOSGeometry *  AliPHOSGeometry::GetInstance(const Text_t* name, const Text_t* title) 
179 {
180   // Returns the pointer of the unique instance
181   // Creates it with the specified options (name, title) if it does not exist yet
182
183   AliPHOSGeometry * rv = 0  ; 
184   if ( fgGeom == 0 ) {
185     if ( strcmp(name,"") == 0 ) 
186       rv = 0 ;
187     else {    
188       fgGeom = new AliPHOSGeometry(name, title) ;
189       if ( fgInit )
190         rv = (AliPHOSGeometry * ) fgGeom ;
191       else {
192         rv = 0 ; 
193         delete fgGeom ; 
194         fgGeom = 0 ; 
195       }
196     }
197   }
198   else {
199     if ( strcmp(fgGeom->GetName(), name) != 0 ) 
200       ::Error("GetInstance", "Current geometry is %s. You cannot call %s", 
201                       fgGeom->GetName(), name) ; 
202     else
203       rv = (AliPHOSGeometry *) fgGeom ; 
204   } 
205   return rv ; 
206 }
207
208 //____________________________________________________________________________
209 void AliPHOSGeometry::SetPHOSAngles() 
210
211   // Calculates the position of the PHOS modules in ALICE global coordinate system
212   // in ideal geometry
213   
214   Double_t const kRADDEG = 180.0 / TMath::Pi() ;
215   Float_t pphi =  2 * TMath::ATan( GetOuterBoxSize(0)  / ( 2.0 * GetIPtoUpperCPVsurface() ) ) ;
216   pphi *= kRADDEG ;
217   if (pphi > fAngle){ 
218     AliError(Form("PHOS modules overlap!\n pphi = %f fAngle = %f", 
219                   pphi, fAngle));
220
221   }
222   pphi = fAngle;
223   
224   for( Int_t i = 1; i <= fNModules ; i++ ) {
225     Float_t angle = pphi * ( i - fNModules / 2.0 - 0.5 ) ;
226     fPHOSAngle[i-1] = -  angle ;
227   } 
228 }
229 //____________________________________________________________________________
230 void AliPHOSGeometry::GetGlobal(const AliRecPoint* , TVector3 & ) const
231 {
232   AliFatal(Form("Please use GetGlobalPHOS(recPoint,gpos) instead of GetGlobal!"));
233 }
234
235 //____________________________________________________________________________
236 void AliPHOSGeometry::GetGlobalPHOS(const AliPHOSRecPoint* recPoint, TVector3 & gpos) const
237 {
238   // Calculates the coordinates of a RecPoint and the error matrix in the ALICE global coordinate system
239  
240   const AliPHOSRecPoint * tmpPHOS = recPoint ;  
241   TVector3 localposition ;
242
243   tmpPHOS->GetLocalPosition(gpos) ;
244
245   if (!gGeoManager){
246     AliFatal("Geo manager not initialized\n");
247   }
248   //construct module name
249   char path[100] ; 
250   Double_t dy ;
251   if(tmpPHOS->IsEmc()){
252     TString spath="/ALIC_1/PHOS_%d/PEMC_1/PCOL_1/PTIO_1/PCOR_1/PAGA_1/PTII_1";
253     snprintf(path,spath.Length(),spath.Data(),tmpPHOS->GetPHOSMod()) ;
254 //    sprintf(path,"/ALIC_1/PHOS_%d",tmpPHOS->GetPHOSMod()) ;
255     dy=fCrystalShift ;
256   }
257   else{
258     TString spath="/ALIC_1/PHOS_%d/PCPV_1";
259     snprintf(path,spath.Length(),spath.Data(),tmpPHOS->GetPHOSMod()) ;
260     dy= GetCPVBoxSize(1)/2. ; //center of CPV module 
261   }
262   Double_t pos[3]={gpos.X(),gpos.Y()-dy,gpos.Z()} ;
263   if(tmpPHOS->IsEmc())
264     pos[2]=-pos[2] ; //Opposite z directions in EMC matrix and local frame!!!
265   Double_t posC[3] = {};
266   //now apply possible shifts and rotations
267   if (!gGeoManager->cd(path)){
268     AliFatal("Geo manager can not find path \n");
269   }
270   TGeoHMatrix *m = gGeoManager->GetCurrentMatrix();
271   if (m){
272      m->LocalToMaster(pos,posC);
273   }
274   else{
275     AliFatal("Geo matrixes are not loaded \n") ;
276   }
277   gpos.SetXYZ(posC[0],posC[1],posC[2]) ;
278
279 }
280 //____________________________________________________________________________
281
282 void AliPHOSGeometry::GetModuleCenter(TVector3& center, 
283                                       const char *det,
284                                       Int_t module) const
285 {
286   // Returns a position of the center of the CPV or EMC module
287   // in ideal (not misaligned) geometry
288   Float_t rDet = 0.;
289   if      (strcmp(det,"CPV") == 0) rDet  = GetIPtoCPVDistance   ();
290   else if (strcmp(det,"EMC") == 0) rDet  = GetIPtoCrystalSurface();
291   else 
292     AliFatal(Form("Wrong detector name %s",det));
293
294   Float_t angle = GetPHOSAngle(module); // (40,20,0,-20,-40) degrees
295   angle *= TMath::Pi()/180;
296   angle += 3*TMath::Pi()/2.;
297   center.SetXYZ(rDet*TMath::Cos(angle), rDet*TMath::Sin(angle), 0.);
298 }
299