add set and getter for neutral energy fraction
[u/mrichter/AliRoot.git] / EMCAL / AliEMCALShishKebabTrd1Module.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-2010, 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 // Main class for TRD1 geometry of Shish-Kebab case.
20 // Author: Alexei Pavlinov(WSU).
21 // Sep 20004 - Nov 2006; Apr 2010
22 // See web page with description of Shish-Kebab geometries:
23 // http://pdsfweb01.nersc.gov/~pavlinov/ALICE/SHISHKEBAB/RES/shishkebabALICE.html
24 // Nov 9,2006 - added case of 3X3
25 //_________________________________________________________________________
26
27 #include "AliLog.h"
28 #include "AliEMCALShishKebabTrd1Module.h"
29 #include "AliEMCALEMCGeometry.h"
30
31 #include <cassert>
32
33 ClassImp(AliEMCALShishKebabTrd1Module)
34
35   //AliEMCALEMCGeometry *AliEMCALShishKebabTrd1Module::fgGeometry=0; 
36   Double_t AliEMCALShishKebabTrd1Module::fga=0.; 
37   Double_t AliEMCALShishKebabTrd1Module::fga2=0.; 
38   Double_t AliEMCALShishKebabTrd1Module::fgb=0.; 
39   Double_t AliEMCALShishKebabTrd1Module::fgr=0.; 
40   Double_t AliEMCALShishKebabTrd1Module::fgangle=0.;   // around one degree 
41   Double_t AliEMCALShishKebabTrd1Module::fgtanBetta=0; //
42
43 //_____________________________________________________________________________
44 AliEMCALShishKebabTrd1Module::AliEMCALShishKebabTrd1Module(Double_t theta, AliEMCALEMCGeometry *g) 
45   : TNamed(),
46     fGeometry(g),
47     fOK(),
48     fA(0.),
49     fB(0.),
50     fThetaA(0.),
51     fTheta(theta),
52     fOK1(),
53     fOK2(),
54     fOB(),
55     fOB1(),
56     fOB2(),
57     fThetaOB1(0.),
58     fThetaOB2(0.),
59     fOK3X3(),
60     fORB(), 
61     fORT()
62
63   TString snam(g->GetName());
64   Int_t key=0;
65   if(snam.Contains("FIRSTYEARv1",TString::kIgnoreCase)) key=1;
66   if(snam.Contains("COMPLETEv1",TString::kIgnoreCase)) key=1;
67   if(GetParameters()) {
68     DefineFirstModule(key);
69   }
70   DefineName(fTheta);
71   AliDebug(10,Form("AliEMCALShishKebabTrd1Module - first module key=%i:  theta %1.4f geometry %s\n",key,fTheta, g->GetName()));
72 }
73
74 //_____________________________________________________________________________
75 AliEMCALShishKebabTrd1Module::AliEMCALShishKebabTrd1Module(AliEMCALShishKebabTrd1Module &leftNeighbor) 
76   : TNamed(),
77     fGeometry(leftNeighbor.fGeometry),
78     fOK(),
79     fA(0.),
80     fB(0.),
81     fThetaA(0.),
82     fTheta(0.),
83     fOK1(),
84     fOK2(),
85     fOB(),
86     fOB1(),
87     fOB2(),
88     fThetaOB1(0.),
89     fThetaOB2(0.),
90     fOK3X3(),
91     fORB(), 
92     fORT()
93
94   //  printf("** Left Neighbor : %s **\n", leftNeighbor.GetName());
95   fTheta  = leftNeighbor.GetTheta() - fgangle; 
96
97   TObject::SetUniqueID(leftNeighbor.GetUniqueID()+1);
98
99   Init(leftNeighbor.GetA(),leftNeighbor.GetB());
100 }
101
102 //________________________________________________________________
103 AliEMCALShishKebabTrd1Module::AliEMCALShishKebabTrd1Module(const AliEMCALShishKebabTrd1Module& mod) 
104   : TNamed(mod.GetName(),mod.GetTitle()),
105     fGeometry(mod.fGeometry),
106     fOK(mod.fOK),
107     fA(mod.fA),
108     fB(mod.fB),
109     fThetaA(mod.fThetaA),
110     fTheta(mod.fTheta),
111     fOK1(mod.fOK1),
112     fOK2(mod.fOK2),
113     fOB(mod.fOB),
114     fOB1(mod.fOB1),
115     fOB2(mod.fOB2),
116     fThetaOB1(mod.fThetaOB1),
117     fThetaOB2(mod.fThetaOB2),
118     fORB(mod.fORB),
119     fORT(mod.fORT)
120 {
121   //copy ctor
122   for (Int_t i=0; i<3; i++) fOK3X3[i] = mod.fOK3X3[i];
123 }
124
125 //________________________________________________________________
126 void AliEMCALShishKebabTrd1Module::Init(Double_t A, Double_t B)
127
128   // Define parameter module from parameters A,B from previos.
129   Double_t yl = (fgb/2)*TMath::Sin(fTheta) + (fga/2)*TMath::Cos(fTheta) + fgr, y = yl;
130   Double_t xl = (yl - B) / A;     // y=A*x+B
131
132   //  Double_t xp1 = (fga/2. + fgb/2.*fgtanBetta)/(TMath::Sin(fTheta) + fgtanBetta*TMath::Cos(fTheta));
133   //  printf(" xp1 %9.3f \n ", xp1);
134   // xp1 == xp => both methods give the same results - 3-feb-05
135   Double_t alpha = TMath::Pi()/2. + fgangle/2;
136   Double_t xt = (fga+fga2)*TMath::Tan(fTheta)*TMath::Tan(alpha)/(4.*(1.-TMath::Tan(fTheta)*TMath::Tan(alpha)));
137   Double_t yt = xt / TMath::Tan(fTheta), xp = TMath::Sqrt(xt*xt + yt*yt);
138   Double_t x  = xl + xp;
139   fOK.Set(x, y);
140   //  printf(" yl %9.3f | xl %9.3f | xp %9.3f \n", yl, xl, xp);
141
142   // have to define A and B; 
143   Double_t yCprev = fgr + fga*TMath::Cos(fTheta);
144   Double_t xCprev = (yCprev - B) / A;
145   Double_t xA     = xCprev + fga*TMath::Sin(fTheta), yA = fgr;
146
147   fThetaA = fTheta - fgangle/2.;
148   fA = TMath::Tan(fThetaA); // !!
149   fB = yA - fA*xA;
150
151   DefineAllStaff();
152 }
153
154 void AliEMCALShishKebabTrd1Module::DefineAllStaff()
155
156   //Define some parameters      
157   DefineName(fTheta);
158   // Centers of cells - 2X2 case
159   Double_t kk1 = (fga+fga2)/(2.*4.); // kk1=kk2 
160
161   Double_t xk1 = fOK.X() - kk1*TMath::Sin(fTheta);
162   Double_t yk1 = fOK.Y() + kk1*TMath::Cos(fTheta) - fgr;
163   fOK1.Set(xk1,yk1);
164
165   Double_t xk2 = fOK.X() + kk1*TMath::Sin(fTheta);
166   Double_t yk2 = fOK.Y() - kk1*TMath::Cos(fTheta) - fgr;
167   fOK2.Set(xk2,yk2);
168
169   // Centers of cells - 3X3 case; Nov 9,2006
170   fOK3X3[1].Set(fOK.X(), fOK.Y()-fgr); // coincide with module center
171
172   kk1 = ((fga+fga2)/4. + fga/6.)/2.; 
173
174   xk1 = fOK.X() - kk1*TMath::Sin(fTheta);
175   yk1 = fOK.Y() + kk1*TMath::Cos(fTheta) - fgr;
176   fOK3X3[0].Set(xk1,yk1);
177
178   xk2 = fOK.X() + kk1*TMath::Sin(fTheta);
179   yk2 = fOK.Y() - kk1*TMath::Cos(fTheta) - fgr;
180   fOK3X3[2].Set(xk2,yk2);
181
182   // May 15, 2006; position of module(cells) center face 
183   fOB.Set(fOK.X()-fgb/2.*TMath::Cos(fTheta),  fOK.Y()-fgb/2.*TMath::Sin(fTheta)-fgr);
184   fOB1.Set(fOB.X()-fga/4.*TMath::Sin(fTheta), fOB.Y()+fga/4.*TMath::Cos(fTheta));
185   fOB2.Set(fOB.X()+fga/4.*TMath::Sin(fTheta), fOB.Y()-fga/4.*TMath::Cos(fTheta));
186   // Jul 30, 2007 - for taking into account a position of shower maximum
187   fThetaOB1 = fTheta - fgangle/4.; // ??
188   fThetaOB2 = fTheta + fgangle/4.;
189
190   // Position of right/top point of module
191   // Gives the posibility to estimate SM size in z direction
192   Double_t xBottom     = (fgr - fB) / fA;
193   Double_t yBottom     = fgr;
194   fORB.Set(xBottom, yBottom);
195
196   Double_t l = fgb / TMath::Cos(fgangle/2.); // length of lateral module side
197   Double_t xTop = xBottom + l*TMath::Cos(TMath::ATan(fA));
198   Double_t yTop = fA*xTop + fB;
199   fORT.Set(xTop, yTop);
200 }
201
202 //_____________________________________________________________________________
203 void AliEMCALShishKebabTrd1Module::DefineFirstModule(const Int_t key)
204 {
205   //    Oct 23-25, 2010
206   //  key=0 - zero tilt of first module;
207   //  key=1 - angle=fgangle/2 = 0.75 degree.
208   //  This is what we have in produced SM. 
209
210   // Define first module
211   if(key==0) {
212     // theta in radians ; first object theta=pi/2.
213     fTheta = TMath::PiOver2();
214     fOK.Set(fga2/2., fgr + fgb/2.); // position the center of module vs o
215
216     // parameters of right line : y = A*z + B in system where zero point is IP.
217     fThetaA = fTheta - fgangle/2.;
218     fA      = TMath::Tan(fThetaA);
219     Double_t xA = fga/2. + fga2/2.;
220     Double_t yA = fgr;
221     fB = yA - fA*xA;
222
223   } else if(key==1) {
224     // theta in radians ; first object theta = 90-0.75 = 89.25 degree
225     fTheta = 89.25*TMath::DegToRad();
226     Double_t al1 = fgangle/2.;
227     Double_t x = 0.5*(fga*TMath::Cos(al1) + fgb*TMath::Sin(al1));
228     Double_t y = 0.5*(fgb + fga*TMath::Sin(al1))*TMath::Cos(al1);
229     fOK.Set(x, fgr + y);
230     // parameters of right line : y = A*z + B in system where zero point is IP.
231     fThetaA = fTheta - fgangle/2.;
232     fA      = TMath::Tan(fThetaA);
233     Double_t xA = fga*TMath::Cos(al1);
234     Double_t yA = fgr;
235     fB = yA - fA*xA;
236
237   } else {
238     printf("<E> key=%i : wrong case \n",key);
239     assert(0);
240   }
241   TObject::SetUniqueID(1); //
242
243   DefineAllStaff();
244 }
245
246 //_____________________________________________________________________________
247 void AliEMCALShishKebabTrd1Module::DefineName(Double_t theta)
248 {
249   // Define name of object
250   SetName(Form("%2i(%5.2f)", TObject::GetUniqueID(), theta*TMath::RadToDeg()));
251 }
252
253 //_____________________________________________________________________________
254 Bool_t AliEMCALShishKebabTrd1Module::GetParameters()
255 {
256   
257   // Get needing module parameters from EMCAL geometry
258   
259   if(!fGeometry) {
260     Warning("GetParameters()"," No geometry ");
261     return kFALSE; 
262   }  
263   
264   TString sn(fGeometry->GetName()); // 2-Feb-05
265   sn.ToUpper();
266   
267   fga        = (Double_t)fGeometry->GetEtaModuleSize();
268   fgb        = (Double_t)fGeometry->GetLongModuleSize();
269   fgangle    = Double_t(fGeometry->GetTrd1Angle())*TMath::DegToRad();
270   fgtanBetta = TMath::Tan(fgangle/2.);
271   fgr        = (Double_t)fGeometry->GetIPDistance();
272
273   fgr += fGeometry->GetSteelFrontThickness();
274
275   fga2       = Double_t(fGeometry->Get2Trd1Dx2());
276   //PH  PrintShish(0);
277   return kTRUE;
278 }
279
280 // service methods
281 //_____________________________________________________________________________
282 void AliEMCALShishKebabTrd1Module::PrintShish(int pri) const
283 {
284   // service method
285   if(pri>=0) { 
286     if(pri >= 1) {
287       printf("PrintShish() \n a %7.3f:%7.3f | b %7.2f | r %7.2f \n TRD1 angle %7.6f(%5.2f) | tanBetta %7.6f", 
288       fga, fga2, fgb, fgr, fgangle, fgangle*TMath::RadToDeg(), fgtanBetta);
289       printf(" fTheta %f : %5.2f : cos(theta) %f\n", 
290       fTheta, GetThetaInDegree(),TMath::Cos(fTheta)); 
291       printf(" OK : %i |%s| theta %f :  phi = %f(%5.2f) \n", 
292       GetUniqueID(), GetName(), fTheta, fOK.Phi(),fOK.Phi()*TMath::RadToDeg());
293     }
294     printf(" y %9.3f x %9.3f xrb %9.3f (right bottom on r=%9.3f ) \n", 
295            fOK.X(), fOK.Y(), fORB.X(),fORB.Y());
296     if(pri>=2) {
297       printf(" A %f B %f | fThetaA %7.6f(%5.2f)\n", fA,fB, fThetaA,fThetaA*TMath::RadToDeg());
298       printf(" fOK  : X %9.4f: Y %9.4f : eta  %5.3f\n",  fOK.X(), fOK.Y(), GetEtaOfCenterOfModule());
299       printf(" fOK1 : X %9.4f: Y %9.4f :   (local, ieta=2)\n", fOK1.X(), fOK1.Y());
300       printf(" fOK2 : X %9.4f: Y %9.4f :   (local, ieta=1)\n\n", fOK2.X(), fOK2.Y());
301       printf(" fOB  : X %9.4f: Y %9.4f \n", fOB.X(), fOB.Y());
302       printf(" fOB1 : X %9.4f: Y %9.4f (local, ieta=2)\n", fOB1.X(), fOB1.Y());
303       printf(" fOB2 : X %9.4f: Y %9.4f (local, ieta=1)\n", fOB2.X(), fOB2.Y());
304       // 3X3 
305       printf(" 3X3 \n");
306       for(int ieta=0; ieta<3; ieta++) {
307         printf(" fOK3X3[%i] : X %9.4f: Y %9.4f (local) \n", ieta, fOK3X3[ieta].X(), fOK3X3[ieta].Y());
308       }
309       //      fOK.Dump();
310       GetMaxEtaOfModule(pri);
311     }
312   }
313 }
314
315 //_____________________________________________________________________________
316 Double_t  AliEMCALShishKebabTrd1Module::GetThetaInDegree() const 
317 {
318   return fTheta*TMath::RadToDeg();
319 }
320
321 //_____________________________________________________________________________
322 Double_t  AliEMCALShishKebabTrd1Module::GetEtaOfCenterOfModule() const 
323
324   return -TMath::Log(TMath::Tan(fOK.Phi()/2.));
325 }
326
327 //_____________________________________________________________________________
328 void AliEMCALShishKebabTrd1Module::GetPositionAtCenterCellLine(Int_t ieta, Double_t dist, TVector2 &v)
329 {
330   // Jul 30, 2007
331   static Double_t theta=0., x=0., y=0.;
332   if(ieta==0) {    
333     v     = fOB2;
334     theta = fTheta;
335   } else if(ieta==1) {
336     v     = fOB1;
337     theta = fTheta;
338   } else {
339     assert(0);
340   }
341   
342   x = v.X() + TMath::Cos(theta) * dist;
343   y = v.Y() + TMath::Sin(theta) * dist;
344   //  printf(" GetPositionAtCenterCellLine() %s : dist %f : ieta %i : x %f %f v.X() | y %f %f v.Y() : cos %f sin %f \n", 
345   //GetName(), dist, ieta, v.X(),x, y,v.Y(),TMath::Cos(theta),TMath::Sin(theta));
346   v.Set(x,y);
347 }  
348
349
350 //_____________________________________________________________________________
351 Double_t  AliEMCALShishKebabTrd1Module::GetMaxEtaOfModule(int pri) const 
352
353   // Right bottom point of module
354   Double_t thetaBottom = TMath::ATan2(fORB.Y(),fORB.X());
355   Double_t etaBottom   = ThetaToEta(thetaBottom);
356   // Right top point of module
357   Double_t thetaTop = TMath::ATan2(fORT.Y(),fORT.X());
358   Double_t etaTop   = ThetaToEta(thetaTop);
359
360   if(pri) { 
361     printf(" Right bottom point of module : eta %5.4f : theta %6.4f (%6.2f) : x(zglob) %7.2f y(phi) %5.2f \n", 
362            etaBottom, thetaBottom, thetaBottom * TMath::RadToDeg(), fORB.X(), fORB.Y());
363     printf(" Right    top point of module : eta %5.4f : theta %6.4f (%6.2f) : x(zglob) %7.2f y(phi) %5.2f \n", 
364            etaTop, thetaTop, thetaTop * TMath::RadToDeg(), fORT.X(), fORT.Y());
365   }
366   return etaBottom>etaTop ? etaBottom : etaTop;
367 }