Fixing warnings which are treated as errors
[u/mrichter/AliRoot.git] / ITS / AliITSInitGeometry.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 /*
17 $Id$
18 */
19 #include <TArrayD.h>
20 #include <TArrayF.h>
21 #include <TStopwatch.h>
22 #include <AliITSgeomSPD.h>
23 #include <AliITSgeomSDD.h>
24 #include <AliITSgeomSSD.h>
25 #include <AliITSsegmentationSPD.h>
26 #include <AliITSsegmentationSDD.h>
27 #include <AliITSsegmentationSSD.h>
28 #include <TGeoManager.h>
29 #include <TGeoVolume.h>
30 #include <TGeoShape.h>
31 #include <TGeoBBox.h>
32 #include <TGeoTrd1.h>
33 #include <TGeoTrd2.h>
34 #include <TGeoArb8.h>
35 #include <TGeoTube.h>
36 #include <TGeoCone.h>
37 #include <TGeoSphere.h>
38 #include <TGeoPara.h>
39 #include <TGeoPgon.h>
40 #include <TGeoPcon.h>
41 #include <TGeoEltu.h>
42 #include <TGeoHype.h>
43 #include <TClass.h>
44
45 #include <AliLog.h>
46 #include "AliITSgeom.h"
47 #include "AliITSInitGeometry.h"
48
49 ClassImp(AliITSInitGeometry)
50 //______________________________________________________________________
51 AliITSInitGeometry::AliITSInitGeometry():
52 TObject(),
53 fName(),
54 fMinorVersion(0),
55 fMajorVersion(0),
56 fTiming(kFALSE),
57 fSegGeom(kFALSE),
58 fDecode(kFALSE){
59     // Default Creator
60     // Inputs:
61     //   none.
62     // Outputs:
63     //   none.
64     // Return:
65     //   A default inilized AliITSInitGeometry object
66 }
67 //______________________________________________________________________
68 AliITSInitGeometry::AliITSInitGeometry(const Char_t *name,Int_t minorversion):
69 TObject(),
70 fName(name),
71 fMinorVersion(minorversion),
72 fMajorVersion(0),
73 fTiming(kFALSE),
74 fSegGeom(kFALSE),
75 fDecode(kFALSE){
76     // Default Creator
77     // Inputs:
78     //   none.
79     // Outputs:
80     //   none.
81     // Return:
82     //   A default inilized AliITSInitGeometry object
83
84     if(fName.CompareTo("AliITSvPPRasymmFMD")==0)if(fMinorVersion==1||
85                                                    fMinorVersion==2){
86         fMajorVersion=10;
87         return;
88     } // end if
89     // if not defined geometry error
90     Error("AliITSInitGeometry(name,version)"," Name must be AliITSvPPRasymmFMD"
91         " and version must be 1 or 2 for now.");
92     fMinorVersion = 0;
93     fName = "";
94     return;
95 }
96 //______________________________________________________________________
97 AliITSgeom* AliITSInitGeometry::CreateAliITSgeom(){
98     // Creates and Initilizes the geometry transformation class AliITSgeom
99     // to values appropreate to this specific geometry. Now that
100     // the segmentation is part of AliITSgeom, the detector
101     // segmentations are also defined here.
102     // Inputs:
103     //   none.
104     // Outputs:
105     //   none.
106     // Return:
107     //   A pointer to a new properly inilized AliITSgeom class. If
108     //   pointer = 0 then failed to init.
109
110     AliITSgeom *geom = new AliITSgeom();
111     if(!InitAliITSgeom(geom)){ // Error initilization failed
112         delete geom;
113         geom = 0;
114     } // end if
115     return geom;
116 }
117 //______________________________________________________________________
118 Bool_t AliITSInitGeometry::InitAliITSgeom(AliITSgeom *geom){
119     // Initilizes the geometry transformation class AliITSgeom
120     // to values appropreate to this specific geometry. Now that
121     // the segmentation is part of AliITSgeom, the detector
122     // segmentations are also defined here.
123     // Inputs:
124     //   AliITSgeom *geom  A pointer to the AliITSgeom class
125     // Outputs:
126     //   AliITSgeom *geom  This pointer recreated and properly inilized.
127     // Return:
128     //   none.
129
130     switch(fMajorVersion){
131     case 10:{ // only case defined so far
132         return InitAliITSgeomPPRasymmFMD(geom);
133     }break; // end case
134     default:{
135         Error("InitAliITSgeom","Undefine geomtery");
136         return kFALSE;
137     } break; // end case
138     } // end switch
139     return kFALSE;
140 }
141 //______________________________________________________________________
142 Bool_t AliITSInitGeometry::InitAliITSgeomPPRasymmFMD(AliITSgeom *geom){
143     // Initilizes the geometry transformation class AliITSgeom
144     // to values appropreate to this specific geometry. Now that
145     // the segmentation is part of AliITSgeom, the detector
146     // segmentations are also defined here.
147     // Inputs:
148     //   AliITSgeom *geom  A pointer to the AliITSgeom class
149     // Outputs:
150     //   AliITSgeom *geom  This pointer recreated and properly inilized.
151     // Return:
152     //   none.
153   //    const Double_t kcm2micron = 1.0E4;
154     const Int_t kItype=0; // Type of transormation defined 0=> Geant
155     const Int_t klayers = 6; // number of layers in the ITS
156     const Int_t kladders[klayers]   = {20,40,14,22,34,38}; // Number of ladders
157     const Int_t kdetectors[klayers] = {4,4,6,8,22,25};// number of detector/lad
158     const AliITSDetector idet[6]   = {kSPD,kSPD,kSDD,kSDD,kSSD,kSSD};
159     const TString pathbase = "/ALIC_1/ITSV_1/ITSD_1/";
160     const TString names[2][klayers] = {
161         {"%sIT12_1/I12A_%d/I10A_%d/I103_%d/I101_1/ITS1_1", // lay=1
162          "%sIT12_1/I12A_%d/I20A_%d/I1D3_%d/I1D1_1/ITS2_1", // lay=2
163          "%sIT34_1/I004_%d/I302_%d/ITS3_%d/", // lay=3
164          "%sIT34_1/I005_%d/I402_%d/ITS4_%d/", // lay=4
165          "%sIT56_1/I565_%d/I562_%d/ITS5_%d/", // lay=5
166          "%sIT56_1/I569_%d/I566_%d/ITS6_%d/"},// lay=6
167         {"%sIT12_1/I12B_%d/I10B_%d/I107_%d/I101_1/ITS1_1", // lay=1
168          "%sIT12_1/I12B_%d/I20B_%d/I1D7_%d/I1D1_1/ITS2_1", // lay=2
169          "%sIT34_1/I004_%d/I302_%d/ITS3_%d", // lay=3
170          "%sIT34_1/I005_%d/I402_%d/ITS4_%d", // lay=4
171          "%sIT56_1/I565_%d/I562_%d/ITS5_%d", // lay=5
172          "%sIT56_1/I569_%d/I566_%d/ITS6_%d"}};// Lay=6
173     /*
174       Int_t itsGeomTreeCopys[knlayers][3]= {{10, 2, 4},// lay=1
175       {10, 4, 4},// lay=2
176       {14, 6, 1},// lay=3
177       {22, 8, 1},// lay=4
178       {34,22, 1},// lay=5
179       {38,25, 1}};//lay=6
180     */
181     Int_t mod,nmods=0,lay,lad,det,cpn0,cpn1,cpn2;
182     Double_t tran[3]={0.0,0.0,0.0},rot[10]={9*0.0,1.0};
183     TArrayD shapePar;
184     TString path,shapeName;
185     TGeoHMatrix materix;
186     Bool_t initSeg[3]={kFALSE,kFALSE,kFALSE};
187     TStopwatch *time = 0x0;if(fTiming) time=new TStopwatch();
188
189     if(fTiming) time->Start();
190     for(mod=0;mod<klayers;mod++) nmods += kladders[mod]*kdetectors[mod];
191     geom->Init(kItype,klayers,kladders,kdetectors,nmods);
192     for(mod=0;mod<nmods;mod++){
193         DecodeDetectorLayers(mod,lay,lad,det); // Write
194         geom->CreateMatrix(mod,lay,lad,det,idet[lay-1],tran,rot);
195         RecodeDetector(mod,cpn0,cpn1,cpn2); // Write reusing lay,lad,det.
196         path.Form(names[fMinorVersion-1][lay-1].Data(),
197                   pathbase.Data(),cpn0,cpn1,cpn2);
198         geom->GetGeomMatrix(mod)->SetPath(path);
199         GetTransformation(path.Data(),materix);
200         geom->SetTrans(mod,materix.GetTranslation());
201         geom->SetRotMatrix(mod,materix.GetRotationMatrix());
202         if(initSeg[idet[lay-1]]) continue;
203         GetShape(path,shapeName,shapePar);
204         if(shapeName.CompareTo("BOX")){
205             Error("InitITSgeom","Geometry changed without proper code update"
206                   "or error in reading geometry. Shape is not BOX.");
207             return kFALSE;
208         } // end if
209         InitGeomShapePPRasymmFMD(idet[lay-1],initSeg,shapePar,geom);
210     } // end for module
211     if(fTiming){
212         time->Stop();
213         time->Print();
214         delete time;
215     } // end if
216     return kTRUE;
217 }
218 //______________________________________________________________________
219 Bool_t AliITSInitGeometry::InitGeomShapePPRasymmFMD(AliITSDetector idet,
220                                                        Bool_t *initSeg,
221                                                        TArrayD &shapePar,
222                                                        AliITSgeom *geom){
223     // Initilizes the geometry segmentation class AliITSgeomS?D, or
224     // AliITSsegmentationS?D depending on the vaule of fSegGeom,
225     // to values appropreate to this specific geometry. Now that
226     // the segmentation is part of AliITSgeom, the detector
227     // segmentations are also defined here.
228     // Inputs:
229     //   Int_t      lay    The layer number/name.
230     //   AliITSgeom *geom  A pointer to the AliITSgeom class
231     // Outputs:
232     //   AliITSgeom *geom  This pointer recreated and properly inilized.
233     // Return:
234     //   none.
235   //   const Double_t kcm2micron = 1.0E4;
236     const Double_t kmicron2cm = 1.0E-4;
237     Int_t i;
238     TArrayF shapeParF;
239
240     shapeParF.Set(shapePar.GetSize());
241     for(i=0;i<shapePar.GetSize();i++) shapeParF[i]=shapePar[i];
242     switch (idet){
243     case kSPD:{
244         initSeg[idet] = kTRUE;
245         AliITSgeomSPD *geomSPD = new AliITSgeomSPD425Short();
246         Float_t bx[256],bz[280];
247         for(i=000;i<256;i++) bx[i] =  50.0*kmicron2cm; // in x all are 50 microns.
248         for(i=000;i<160;i++) bz[i] = 425.0*kmicron2cm; // most are 425 microns
249         // except below
250         for(i=160;i<280;i++) bz[i] =   0.0*kmicron2cm; // Outside of detector.
251         bz[ 31] = bz[ 32] = 625.0*kmicron2cm; // first chip boundry
252         bz[ 63] = bz[ 64] = 625.0*kmicron2cm; // first chip boundry
253         bz[ 95] = bz[ 96] = 625.0*kmicron2cm; // first chip boundry
254         bz[127] = bz[128] = 625.0*kmicron2cm; // first chip boundry
255         bz[160] = 425.0*kmicron2cm;// Set so that there is no zero pixel size for fNz.
256         geomSPD->ReSetBins(shapeParF[1],256,bx,160,bz);
257         geom->ReSetShape(idet,geomSPD);
258     }break;
259     case kSDD:{
260         initSeg[idet] = kTRUE;
261         AliITSgeomSDD *geomSDD = new AliITSgeomSDD256(shapeParF.GetSize(),
262                                                       shapeParF.GetArray());
263         geom->ReSetShape(idet,geomSDD);
264     }break;
265     case kSSD:{
266         initSeg[idet] = kTRUE;
267         AliITSgeomSSD *geomSSD = new AliITSgeomSSD275and75(
268             shapeParF.GetSize(),shapeParF.GetArray());
269         geom->ReSetShape(idet,geomSSD);
270     }break;
271     default:{// Others, Note no kSDDp or kSSDp in this geometry.
272         geom->ReSetShape(idet,0);
273         Info("InitGeomShapePPRasymmFMD",
274              "default Dx=%f Dy=%f Dz=%f default=%d",
275              shapePar[0],shapePar[1],shapePar[2],idet);
276     }break;
277     } // end switch
278     return kTRUE;
279 }
280 //______________________________________________________________________
281 Bool_t AliITSInitGeometry::InitSegmentationPPRasymmFMD(AliITSDetector idet,
282                                                        Bool_t *initSeg,
283                                                        TArrayD &shapePar,
284                                                        AliITSgeom *geom){
285     // Initilizes the geometry segmentation class AliITSgeomS?D, or
286     // AliITSsegmentationS?D depending on the vaule of fSegGeom,
287     // to values appropreate to this specific geometry. Now that
288     // the segmentation is part of AliITSgeom, the detector
289     // segmentations are also defined here.
290     // Inputs:
291     //   Int_t      lay    The layer number/name.
292     //   AliITSgeom *geom  A pointer to the AliITSgeom class
293     // Outputs:
294     //   AliITSgeom *geom  This pointer recreated and properly inilized.
295     // Return:
296     //   none.
297     const Double_t kcm2micron = 1.0E4;
298     Int_t i;
299
300     switch (idet){
301     case kSPD:{
302         initSeg[idet] = kTRUE;
303         AliITSsegmentationSPD *segSPD = new AliITSsegmentationSPD();
304         segSPD->SetDetSize(2.*shapePar[0]*kcm2micron, // X
305                            2.*shapePar[2]*kcm2micron, // Z
306                            2.*shapePar[1]*kcm2micron);// Y  Microns
307         segSPD->SetNPads(256,160);// Number of Bins in x and z
308         Float_t bx[256],bz[280];
309         for(i=000;i<256;i++) bx[i] =  50.0; // in x all are 50 microns.
310         for(i=000;i<160;i++) bz[i] = 425.0; // most are 425 microns
311         // except below
312         for(i=160;i<280;i++) bz[i] =   0.0; // Outside of detector.
313         bz[ 31] = bz[ 32] = 625.0; // first chip boundry
314         bz[ 63] = bz[ 64] = 625.0; // first chip boundry
315         bz[ 95] = bz[ 96] = 625.0; // first chip boundry
316         bz[127] = bz[128] = 625.0; // first chip boundry
317         bz[160] = 425.0;// Set so that there is no zero pixel size for fNz.
318         segSPD->SetBinSize(bx,bz); // Based on AliITSgeomSPD for now.
319         geom->ReSetShape(idet,segSPD);
320     }break;
321     case kSDD:{
322         initSeg[idet] = kTRUE;
323         AliITSsegmentationSDD *segSDD = new AliITSsegmentationSDD();
324         segSDD->SetDetSize(shapePar[0]*kcm2micron, // X
325                            2.*shapePar[2]*kcm2micron, // Z
326                            2.*shapePar[1]*kcm2micron);// Y  Microns
327         segSDD->SetNPads(256,256);// Anodes, Samples
328         geom->ReSetShape(idet,segSDD);
329     }break;
330     case kSSD:{
331         initSeg[idet] = kTRUE;
332         AliITSsegmentationSSD *segSSD = new AliITSsegmentationSSD();
333         segSSD->SetDetSize(2.*shapePar[0]*kcm2micron, // X
334                            2.*shapePar[2]*kcm2micron, // Z
335                            2.*shapePar[1]*kcm2micron);// Y  Microns.
336         segSSD->SetPadSize(95.,0.); // strip x pitch in microns
337         segSSD->SetNPads(768,2); // number of strips on each side, sides.
338         segSSD->SetAngles(0.0075,0.0275); // strip angels rad P and N side.
339         segSSD->SetAnglesLay5(0.0075,0.0275);//strip angels rad P and N
340         segSSD->SetAnglesLay6(0.0275,0.0075);//strip angels rad P and N
341         geom->ReSetShape(idet,segSSD);
342     }break;
343     default:{// Others, Note no kSDDp or kSSDp in this geometry.
344         geom->ReSetShape(idet,0);
345         Info("InitSegmentationPPRasymmFMD",
346              "default segmentation Dx=%f Dy=%f Dz=%f default=%d",
347              shapePar[0],shapePar[1],shapePar[2],idet);
348     }break;
349     } // end switch
350     return kTRUE;
351 }
352 //______________________________________________________________________
353 Bool_t AliITSInitGeometry::GetTransformation(const TString &volumePath,
354                                              TGeoHMatrix &mat){
355     // Returns the Transformation matrix between the volume specified
356     // by the path volumePath and the Top or mater volume. The format
357     // of the path volumePath is as follows (assuming ALIC is the Top volume)
358     // "/ALIC_1/DDIP_1/S05I_2/S05H_1/S05G_3". Here ALIC is the top most
359     // or master volume which has only 1 instance of. Of all of the daughter
360     // volumes of ALICE, DDIP volume copy #1 is indicated. Similarly for
361     // the daughter volume of DDIP is S05I copy #2 and so on.
362     // Inputs:
363     //   TString& volumePath  The volume path to the specific volume
364     //                        for which you want the matrix. Volume name
365     //                        hierarchy is separated by "/" while the
366     //                        copy number is appended using a "_".
367     // Outputs:
368     //  TGeoHMatrix &mat      A matrix with its values set to those
369     //                        appropriate to the Local to Master transformation
370     // Return:
371     //   A logical value if kFALSE then an error occurred and no change to
372     //   mat was made.
373
374     // We have to preserve the modeler state
375
376     // Preserve the modeler state.
377     gGeoManager->PushPath();
378     if (!gGeoManager->cd(volumePath.Data())) {
379         gGeoManager->PopPath();
380         Error("GetTransformation","Error in cd-ing to ",volumePath.Data());
381         return kFALSE;
382     } // end if !gGeoManager
383     mat = *gGeoManager->GetCurrentMatrix();
384     // Retstore the modeler state.
385     gGeoManager->PopPath();
386     return kTRUE;
387 }
388 //______________________________________________________________________
389 Bool_t AliITSInitGeometry::GetShape(const TString &volumePath,
390                                     TString &shapeType,TArrayD &par){
391     // Returns the shape and its parameters for the volume specified
392     // by volumeName.
393     // Inputs:
394     //   TString& volumeName  The volume name
395     // Outputs:
396     //   TString &shapeType   Shape type
397     //   TArrayD &par         A TArrayD of parameters with all of the
398     //                        parameters of the specified shape.
399     // Return:
400     //   A logical indicating whether there was an error in getting this
401     //   information
402     Int_t npar;
403     gGeoManager->PushPath();
404     if (!gGeoManager->cd(volumePath.Data())) {
405         gGeoManager->PopPath();
406         return kFALSE;
407     }
408     TGeoVolume * vol = gGeoManager->GetCurrentVolume();
409     gGeoManager->PopPath();
410     if (!vol) return kFALSE;
411     TGeoShape *shape = vol->GetShape();
412     TClass *class_type = shape->IsA();
413     if (class_type==TGeoBBox::Class()) {
414         shapeType = "BOX";
415         npar = 3;
416         par.Set(npar);
417         TGeoBBox *box = (TGeoBBox*)shape;
418         par.AddAt(box->GetDX(),0);
419         par.AddAt(box->GetDY(),1);
420         par.AddAt(box->GetDZ(),2);
421         return kTRUE;
422     }
423     if (class_type==TGeoTrd1::Class()) {
424         shapeType = "TRD1";
425         npar = 4;
426         par.Set(npar);
427         TGeoTrd1 *trd1 = (TGeoTrd1*)shape;
428         par.AddAt(trd1->GetDx1(),0);
429         par.AddAt(trd1->GetDx2(),1);
430         par.AddAt(trd1->GetDy(), 2);
431         par.AddAt(trd1->GetDz(), 3);
432         return kTRUE;
433     }
434     if (class_type==TGeoTrd2::Class()) {
435         shapeType = "TRD2";
436         npar = 5;
437         par.Set(npar);
438         TGeoTrd2 *trd2 = (TGeoTrd2*)shape;
439         par.AddAt(trd2->GetDx1(),0);
440         par.AddAt(trd2->GetDx2(),1);
441         par.AddAt(trd2->GetDy1(),2);
442         par.AddAt(trd2->GetDy2(),3);
443         par.AddAt(trd2->GetDz(), 4);
444         return kTRUE;
445     }
446     if (class_type==TGeoTrap::Class()) {
447         shapeType = "TRAP";
448         npar = 11;
449         par.Set(npar);
450         TGeoTrap *trap = (TGeoTrap*)shape;
451         Double_t tth = TMath::Tan(trap->GetTheta()*TMath::DegToRad());
452         par.AddAt(trap->GetDz(),0);
453         par.AddAt(tth*TMath::Cos(trap->GetPhi()*TMath::DegToRad()),1);
454         par.AddAt(tth*TMath::Sin(trap->GetPhi()*TMath::DegToRad()),2);
455         par.AddAt(trap->GetH1(),3);
456         par.AddAt(trap->GetBl1(),4);
457         par.AddAt(trap->GetTl1(),5);
458         par.AddAt(TMath::Tan(trap->GetAlpha1()*TMath::DegToRad()),6);
459         par.AddAt(trap->GetH2(),7);
460         par.AddAt(trap->GetBl2(),8);
461         par.AddAt(trap->GetTl2(),9);
462         par.AddAt(TMath::Tan(trap->GetAlpha2()*TMath::DegToRad()),10);
463         return kTRUE;
464     }
465     if (class_type==TGeoTube::Class()) {
466         shapeType = "TUBE";
467         npar = 3;
468         par.Set(npar);
469         TGeoTube *tube = (TGeoTube*)shape;
470         par.AddAt(tube->GetRmin(),0);
471         par.AddAt(tube->GetRmax(),1);
472         par.AddAt(tube->GetDz(),2);
473         return kTRUE;
474     }
475     if (class_type==TGeoTubeSeg::Class()) {
476         shapeType = "TUBS";
477         npar = 5;
478         par.Set(npar);
479         TGeoTubeSeg *tubs = (TGeoTubeSeg*)shape;
480         par.AddAt(tubs->GetRmin(),0);
481         par.AddAt(tubs->GetRmax(),1);
482         par.AddAt(tubs->GetDz(),2);
483         par.AddAt(tubs->GetPhi1(),3);
484         par.AddAt(tubs->GetPhi2(),4);
485         return kTRUE;
486     }
487     if (class_type==TGeoCone::Class()) {
488         shapeType = "CONE";
489         npar = 5;
490         par.Set(npar);
491         TGeoCone *cone = (TGeoCone*)shape;
492         par.AddAt(cone->GetDz(),0);
493         par.AddAt(cone->GetRmin1(),1);
494         par.AddAt(cone->GetRmax1(),2);
495         par.AddAt(cone->GetRmin2(),3);
496         par.AddAt(cone->GetRmax2(),4);
497         return kTRUE;
498     }
499     if (class_type==TGeoConeSeg::Class()) {
500         shapeType = "CONS";
501         npar = 7;
502         par.Set(npar);
503         TGeoConeSeg *cons = (TGeoConeSeg*)shape;
504         par.AddAt(cons->GetDz(),0);
505         par.AddAt(cons->GetRmin1(),1);
506         par.AddAt(cons->GetRmax1(),2);
507         par.AddAt(cons->GetRmin2(),3);
508         par.AddAt(cons->GetRmax2(),4);
509         par.AddAt(cons->GetPhi1(),5);
510         par.AddAt(cons->GetPhi2(),6);
511         return kTRUE;
512     }
513     if (class_type==TGeoSphere::Class()) {
514         shapeType = "SPHE";
515         npar = 6;
516         par.Set(npar);
517         
518         TGeoSphere *sphe = (TGeoSphere*)shape;
519         par.AddAt(sphe->GetRmin(),0);
520         par.AddAt(sphe->GetRmax(),1);
521         par.AddAt(sphe->GetTheta1(),2);
522         par.AddAt(sphe->GetTheta2(),3);
523         par.AddAt(sphe->GetPhi1(),4);
524         par.AddAt(sphe->GetPhi2(),5);
525         return kTRUE;
526     }
527     if (class_type==TGeoPara::Class()) {
528         shapeType = "PARA";
529         npar = 6;
530         par.Set(npar);
531         TGeoPara *para = (TGeoPara*)shape;
532         par.AddAt(para->GetX(),0);
533         par.AddAt(para->GetY(),1);
534         par.AddAt(para->GetZ(),2);
535         par.AddAt(para->GetTxy(),3);
536         par.AddAt(para->GetTxz(),4);
537         par.AddAt(para->GetTyz(),5);
538         return kTRUE;
539     }
540     if (class_type==TGeoPgon::Class()) {
541         shapeType = "PGON";
542         TGeoPgon *pgon = (TGeoPgon*)shape;
543         Int_t nz = pgon->GetNz();
544         const Double_t *rmin = pgon->GetRmin();
545         const Double_t *rmax = pgon->GetRmax();
546         const Double_t *z = pgon->GetZ();
547         npar = 4 + 3*nz;
548         par.Set(npar);
549         par.AddAt(pgon->GetPhi1(),0);
550         par.AddAt(pgon->GetDphi(),1);
551         par.AddAt(pgon->GetNedges(),2);
552         par.AddAt(pgon->GetNz(),3);
553         for (Int_t i=0; i<nz; i++) {
554             par.AddAt(z[i], 4+3*i);
555             par.AddAt(rmin[i], 4+3*i+1);
556             par.AddAt(rmax[i], 4+3*i+2);
557         }
558         return kTRUE;
559     }
560     if (class_type==TGeoPcon::Class()) {
561         shapeType = "PCON";
562         TGeoPcon *pcon = (TGeoPcon*)shape;
563         Int_t nz = pcon->GetNz();
564         const Double_t *rmin = pcon->GetRmin();
565         const Double_t *rmax = pcon->GetRmax();
566         const Double_t *z = pcon->GetZ();
567         npar = 3 + 3*nz;
568         par.Set(npar);
569         par.AddAt(pcon->GetPhi1(),0);
570         par.AddAt(pcon->GetDphi(),1);
571         par.AddAt(pcon->GetNz(),2);
572         for (Int_t i=0; i<nz; i++) {
573             par.AddAt(z[i], 3+3*i);
574             
575             par.AddAt(rmin[i], 3+3*i+1);
576             par.AddAt(rmax[i], 3+3*i+2);
577         }
578         return kTRUE;
579     }
580     if (class_type==TGeoEltu::Class()) {
581         shapeType = "ELTU";
582         npar = 3;
583         par.Set(npar);
584         TGeoEltu *eltu = (TGeoEltu*)shape;
585         par.AddAt(eltu->GetA(),0);
586         par.AddAt(eltu->GetB(),1);
587         par.AddAt(eltu->GetDz(),2);
588         return kTRUE;
589     }
590     if (class_type==TGeoHype::Class()) {
591         shapeType = "HYPE";
592         npar = 5;
593         par.Set(npar);
594         TGeoHype *hype = (TGeoHype*)shape;
595         par.AddAt(TMath::Sqrt(hype->RadiusHypeSq(0.,kTRUE)),0);
596         par.AddAt(TMath::Sqrt(hype->RadiusHypeSq(0.,kFALSE)),1);
597         par.AddAt(hype->GetDZ(),2);
598         par.AddAt(hype->GetStIn(),3);
599         par.AddAt(hype->GetStOut(),4);
600         return kTRUE;
601     }
602     if (class_type==TGeoGtra::Class()) {
603         shapeType = "GTRA";
604         npar = 12;
605         par.Set(npar);
606         TGeoGtra *trap = (TGeoGtra*)shape;
607         Double_t tth = TMath::Tan(trap->GetTheta()*TMath::DegToRad());
608         par.AddAt(trap->GetDz(),0);
609         par.AddAt(tth*TMath::Cos(trap->GetPhi()*TMath::DegToRad()),1);
610         par.AddAt(tth*TMath::Sin(trap->GetPhi()*TMath::DegToRad()),2);
611         par.AddAt(trap->GetH1(),3);
612         par.AddAt(trap->GetBl1(),4);
613         par.AddAt(trap->GetTl1(),5);
614         par.AddAt(TMath::Tan(trap->GetAlpha1()*TMath::DegToRad()),6);
615         par.AddAt(trap->GetH2(),7);
616         par.AddAt(trap->GetBl2(),8);
617         par.AddAt(trap->GetTl2(),9);
618         par.AddAt(TMath::Tan(trap->GetAlpha2()*TMath::DegToRad()),10);
619         par.AddAt(trap->GetTwistAngle(),11);
620         return kTRUE;
621     }
622     if (class_type==TGeoCtub::Class()) {
623         shapeType = "CTUB";
624         npar = 11;
625         par.Set(npar);
626         TGeoCtub *ctub = (TGeoCtub*)shape;
627         const Double_t *lx = ctub->GetNlow();
628         const Double_t *tx = ctub->GetNhigh();
629         par.AddAt(ctub->GetRmin(),0);
630         par.AddAt(ctub->GetRmax(),1);
631         par.AddAt(ctub->GetDz(),2);
632         par.AddAt(ctub->GetPhi1(),3);
633         par.AddAt(ctub->GetPhi2(),4);
634         par.AddAt(lx[0],5);
635         par.AddAt(lx[1],6);
636         par.AddAt(lx[2],7);
637         par.AddAt(tx[0],8);
638         par.AddAt(tx[1],9);
639         par.AddAt(tx[2],10);
640         return kTRUE;
641     }
642     Error("GetShape","Getting shape parameters for shape %s not implemented",
643           shape->ClassName());
644     return kFALSE;
645 }
646 //______________________________________________________________________
647 void AliITSInitGeometry::DecodeDetector(Int_t &mod,Int_t layer,Int_t cpn0,
648                                         Int_t cpn1,Int_t cpn2){
649     // decode geometry into detector module number. There are two decoding
650     // Scheams. Old which does not follow the ALICE coordinate system
651     // requirements, and New which dose.
652     // Inputs:
653     //    Int_t layer    The ITS layer
654     //    Int_t cpn0     The lowest copy number
655     //    Int_t cpn1     The middle copy number
656     //    Int_t cpn2     the highest copy number
657     // Output:
658     //    Int_t &mod     The module number assoicated with this set
659     //                   of copy numbers.
660     // Return:
661     //    none.
662     const Int_t detPerLadderSPD[2]={2,4};
663     const Int_t detPerLadder[6]={4,4,6,8,22,25};
664     const Int_t ladPerLayer[6]={20,40,14,22,34,38};
665     Int_t lay=-1,lad=-1,det=-1,i;
666
667     if(fDecode){ // New decoding scheam
668         switch (layer){
669         case 1:{
670             lay = layer;
671             det = 5-cpn2;
672             if(cpn0==4&&cpn1==1) lad=1;
673             else if(cpn0==4&&cpn1==2) lad=20;
674             else if(cpn0<4){
675                 lad = 8-cpn1-detPerLadderSPD[layer-1]*(cpn0-1);
676             }else{ // cpn0>4
677                 lad = 28-cpn1-detPerLadderSPD[layer-1]*(cpn0-1);
678             } // end if
679         } break;
680         case 2:{
681             lay = layer;
682             det = 5-cpn2;
683             if(cpn0==4&&cpn1==1) lad=1;
684             else if(cpn0<4){
685                 lad = 14-cpn1-detPerLadderSPD[layer-1]*(cpn0-1);
686             }else{ // cpn0>4
687                 lad = 54-cpn1-detPerLadderSPD[layer-1]*(cpn0-1);
688             } // end if
689         } break;
690         case 3:{
691             lay = layer;
692             if(cpn0<5) lad = 5-cpn0;
693             else lad = 19-cpn0;
694             det = 7-cpn1;
695         } break;
696         case 4:{
697             lay = layer;
698             if(cpn0<7) lad = 7-cpn0;
699             else lad = 29-cpn0;
700             det = 9-cpn1;
701         } break;
702         case 5:{
703             lay = layer;
704             if(cpn0<10) lad = 10-cpn0;
705             else lad = 44-cpn0;
706             det = 23-cpn1;
707         } break;
708         case 6:{
709             lay = layer;
710             if(cpn0<9) lad = 9-cpn0;
711             else lad = 47-cpn0;
712             det = 26-cpn1;
713         } break;
714         } // end switch
715         mod = 0;
716         for(i=0;i<layer-1;i++) mod += ladPerLayer[i]*detPerLadder[i];
717         mod += detPerLadder[layer-1]*(lad-1)+det-1;// module start at zero.
718         return;
719     } // end if
720     // Old decoding scheam
721     switch(layer){
722     case 1: case 2:{
723         lay = layer;
724         lad = cpn1+detPerLadderSPD[layer-1]*(cpn0-1);
725         det = cpn2;
726         }break;
727     case 3: case 4:{
728         lay = layer;
729         lad = cpn0;
730         det = cpn1;
731         }break;
732     case 5: case 6:{
733         lay = layer;
734         lad = cpn0;
735         det = cpn1;
736         }break;
737     default:{
738         }break;
739     } // end switch
740     mod = 0;
741     for(i=0;i<layer-1;i++) mod += ladPerLayer[i]*detPerLadder[i];
742     mod += detPerLadder[layer-1]*(lad-1)+det-1;// module start at zero.
743     return;
744 }
745 //______________________________________________________________________
746 void AliITSInitGeometry::RecodeDetector(Int_t mod,Int_t &cpn0,
747                                         Int_t &cpn1,Int_t &cpn2){
748     // decode geometry into detector module number. There are two decoding
749     // Scheams. Old which does not follow the ALICE coordinate system
750     // requirements, and New which dose.
751     // Inputs:
752     //    Int_t mod      The module number assoicated with this set
753     //                   of copy numbers.
754     // Output:
755     //    Int_t cpn0     The lowest copy number
756     //    Int_t cpn1     The middle copy number
757     //    Int_t cpn2     the highest copy number
758     // Return:
759     //    none.
760     const Int_t itsGeomTreeCopys[6][3]= {{10, 2, 4},// lay=1
761                                          {10, 4, 4},// lay=2
762                                          {14, 6, 1},// lay=3
763                                          {22, 8, 1},// lay=4
764                                          {34,22, 1},// lay=5
765                                          {38,25, 1}};//lay=6
766     const Int_t detPerLadderSPD[2]={2,4};
767     //    const Int_t detPerLadder[6]={4,4,6,8,22,25};
768     //    const Int_t ladPerLayer[6]={20,40,14,22,34,38};
769     Int_t lay,lad,det;
770
771     cpn0 = cpn1 = cpn2 = 0;
772     DecodeDetectorLayers(mod,lay,lad,det);
773     if(fDecode){ // New decoding scheam
774         switch (lay){
775         case 1:{
776             cpn2 = 5-det;     // Detector 1-4
777             cpn1 = 1+(lad-1)%detPerLadderSPD[lay-1];
778             cpn0 = 5-(lad+detPerLadderSPD[lay-1])/detPerLadderSPD[lay-1];
779             if(mod>27) cpn0 = 15-(lad+detPerLadderSPD[lay-1])/
780                            detPerLadderSPD[lay-1];
781         } break;
782         case 2:{
783             cpn2 = 5-det;     // Detector 1-4
784             cpn1 = 4-(lad+2)%detPerLadderSPD[lay-1];
785             cpn0 = 1+(14-cpn1-lad)/detPerLadderSPD[lay-1];
786             if(mod>131) cpn0 = 1+(54-lad-cpn1)/detPerLadderSPD[lay-1];
787         } break;
788         case 3:{
789             cpn2 = 1;
790             if(lad<5) cpn0 = 5-lad;
791             else cpn0 = 19-lad;
792             cpn1 = 7-det;
793         } break;
794         case 4:{
795             cpn2 = 1;
796             if(lad<7) cpn0 = 7-lad;
797             else cpn0 = 29-lad;
798             cpn1 = 9-det;
799         } break;
800         case 5:{
801             cpn2 = 1;
802             if(lad<10) cpn0 = 10-lad;
803             else cpn0 = 44-lad;
804             cpn1 = 23-det;
805         } break;
806         case 6:{
807             cpn2 = 1;
808             if(lad<9) cpn0 = 9-lad;
809             else cpn0 = 47-lad;
810             cpn1 = 26-det;
811         } break;
812         default:{
813             Error("RecodeDetector","New: mod=%d lay=%d not 1-6.");
814             return;
815         } break;
816         } // end switch
817         if(cpn0<1||cpn1<1||cpn2<1||
818            cpn0>itsGeomTreeCopys[lay-1][0]||
819            cpn1>itsGeomTreeCopys[lay-1][1]||
820            cpn2>itsGeomTreeCopys[lay-1][2])
821             Error("RecodeDetector",
822                   "cpn0=%d cpn1=%d cpn2=%d mod=%d lay=%d lad=%d det=%d",
823                   cpn0,cpn1,cpn2,mod,lay,lad,det);
824         return;
825     } // end if
826     // Old encoding
827     switch (lay){
828     case 1: case 2:{
829         cpn2 = det;     // Detector 1-4
830         cpn0 = (lad+detPerLadderSPD[lay-1]-1)/detPerLadderSPD[lay-1];
831         cpn1 = (lad+detPerLadderSPD[lay-1]-1)%detPerLadderSPD[lay-1] + 1;
832     } break;
833     case 3: case 4: case 5 : case 6:{
834         cpn2 = 1;
835         cpn1 = det;
836         cpn0 = lad;
837     } break;
838     default:{
839         Error("RecodeDetector","Old: mod=%d lay=%d not 1-6.");
840         return;
841     } break;
842     } // end switch
843     if(cpn0<1||cpn1<1||cpn2<1||
844        cpn0>itsGeomTreeCopys[lay-1][0]||
845        cpn1>itsGeomTreeCopys[lay-1][1]||
846        cpn2>itsGeomTreeCopys[lay-1][2])
847         Error("RecodeDetector",
848               "cpn0=%d cpn1=%d cpn2=%d mod=%d lay=%d lad=%d det=%d",
849               cpn0,cpn1,cpn2,mod,lay,lad,det);
850     return;
851 }
852 //______________________________________________________________________
853 void AliITSInitGeometry::DecodeDetectorLayers(Int_t mod,Int_t &lay,
854                                               Int_t &lad,Int_t &det){
855     // decode geometry into detector module number. There are two decoding
856     // Scheams. Old which does not follow the ALICE coordinate system
857     // requirements, and New which dose. Note, this use of layer ladder
858     // and detector numbers are strictly for internal use of this
859     // specific code. They do not represent the "standard" layer ladder
860     // or detector numbering except in a very old and obsoleate sence.
861     // Inputs:
862     //    Int_t mod      The module number assoicated with this set
863     //                   of copy numbers.
864     // Output:
865     //    Int_t lay     The layer number
866     //    Int_t lad     The ladder number
867     //    Int_t det     the dettector number
868     // Return:
869     //    none.
870   //    const Int_t detPerLadderSPD[2]={2,4};
871     const Int_t detPerLadder[6]={4,4,6,8,22,25};
872     const Int_t ladPerLayer[6]={20,40,14,22,34,38};
873     Int_t mod2;
874
875     det  = 0;
876     lad  = 0;
877     lay  = 0;
878     mod2 = 0;
879     do{
880         mod2 += ladPerLayer[lay]*detPerLadder[lay];
881         lay++;
882     }while(mod2<=mod); // end while
883     if(lay>6||lay<1) Error("DecodeDetectorLayers","0<lay=%d>6",lay);
884     mod2 -= ladPerLayer[lay-1]*detPerLadder[lay-1];
885     do{
886         lad++;
887         mod2 += detPerLadder[lay-1];
888     }while(mod2<=mod); // end while
889     if(lad>ladPerLayer[lay-1]||lad<1) Error("DecodeDetectorLayera",
890             "lad=%d>ladPerLayer[lay-1=%d]=%d mod=%d mod2=%d",lad,lay-1,
891                                             ladPerLayer[lay-1],mod,mod2);
892     mod2 -= detPerLadder[lay-1];
893     det = mod-mod2+1;
894     if(det>detPerLadder[lay-1]||det<1) Error("DecodeDetectorLayers",
895            "det=%d>detPerLayer[lay-1=%d]=%d mod=%d mod2=%d lad=%d",det,
896                                   lay-1,detPerLadder[lay-1],mod,mod2,lad);
897     return;
898 }
899