Another histos for lumi
[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 ////////////////////////////////////////////////////////////////
20 //  This class initializes the class AliITSgeom
21 //  The initialization is done starting from 
22 //  a geometry coded by means of the ROOT geometrical modeler
23 //  This initialization can be used both for simulation and reconstruction
24 ///////////////////////////////////////////////////////////////
25
26 #include <TArrayD.h>
27 #include <TArrayF.h>
28 #include <TStopwatch.h>
29 #include <TGeoManager.h>
30 #include <TGeoMatrix.h>
31 #include <TGeoVolume.h>
32 #include <TGeoShape.h>
33 #include <TGeoBBox.h>
34 #include <TGeoTrd1.h>
35 #include <TGeoTrd2.h>
36 #include <TGeoArb8.h>
37 #include <TGeoTube.h>
38 #include <TGeoCone.h>
39 #include <TGeoSphere.h>
40 #include <TGeoPara.h>
41 #include <TGeoPgon.h>
42 #include <TGeoPcon.h>
43 #include <TGeoEltu.h>
44 #include <TGeoHype.h>
45 #include <TMath.h>
46
47 #include "AliLog.h"
48 #include "AliITSsegmentationSPD.h"
49 #include "AliITSsegmentationSDD.h"
50 #include "AliITSsegmentationSSD.h"
51 #include "AliITSInitGeometry.h"
52 #include <TDatime.h>
53
54 ClassImp(AliITSInitGeometry)
55
56 const Bool_t AliITSInitGeometry::fgkOldSPDbarrel = kFALSE;
57 const Bool_t AliITSInitGeometry::fgkOldSDDbarrel = kFALSE;
58 const Bool_t AliITSInitGeometry::fgkOldSSDbarrel = kFALSE;
59 const Bool_t AliITSInitGeometry::fgkOldSDDcone   = kFALSE;
60 const Bool_t AliITSInitGeometry::fgkOldSSDcone   = kFALSE;
61 const Bool_t AliITSInitGeometry::fgkOldSPDshield = kFALSE;
62 const Bool_t AliITSInitGeometry::fgkOldSDDshield = kTRUE;
63 const Bool_t AliITSInitGeometry::fgkOldSSDshield = kTRUE;
64 const Bool_t AliITSInitGeometry::fgkOldServices  = kFALSE;
65 const Bool_t AliITSInitGeometry::fgkOldSupports  = kFALSE;
66 //______________________________________________________________________
67 AliITSInitGeometry::AliITSInitGeometry():
68 TObject(),                   // Base Class
69 fName(0),                    // Geometry name
70 fMinorVersion(-1),           // Minor version number/type
71 fMajorVersion(kvDefault),    // Major versin number
72 fTiming(kFALSE),             // Flag to start inilization timing
73 fSegGeom(kFALSE),            // Flag to switch between the old use of
74                              // AliITSgeomS?D class, or AliITSsegmentation
75                              // class in fShape of AliITSgeom class.
76 fDecode(kFALSE),             // Flag for new/old decoding
77 fDebug(0){                   // Debug flag
78     // Default Creator
79     // Inputs:
80     //   none.
81     // Outputs:
82     //   none.
83     // Return:
84     //   A default inilized AliITSInitGeometry object
85
86     fName = "Undefined";
87 }
88 //______________________________________________________________________
89 AliITSInitGeometry::AliITSInitGeometry(AliITSVersion_t version,
90                                        Int_t minorversion):
91 TObject(),                   // Base Class
92 fName(0),                    // Geometry name
93 fMinorVersion(minorversion), // Minor version number/type
94 fMajorVersion(version),      // Major versin number
95 fTiming(kFALSE),             // Flag to start inilization timing
96 fSegGeom(kFALSE),            // Flag to switch between the old use of
97                              // AliITSgeomS?D class, or AliITSsegmentation
98                              // class in fShape of AliITSgeom class.
99 fDecode(kFALSE),             // Flag for new/old decoding
100 fDebug(0){                   // Debug flag
101     // Default Creator
102     // Inputs:
103     //   none.
104     // Outputs:
105     //   none.
106     // Return:
107     //   A default inilized AliITSInitGeometry object
108
109   switch (version) {
110     case kv11:
111         fName="AliITSv11";
112         break;
113     case kvDefault:
114     default:
115         AliFatal(Form("Undefined geometry: fMajorVersion=%d, "
116                       "fMinorVersion= %d",(Int_t)fMajorVersion,fMinorVersion));
117         fName = "Undefined";
118         break;
119     } // switch
120 }
121 //______________________________________________________________________
122 AliITSgeom* AliITSInitGeometry::CreateAliITSgeom(){
123     // Creates and Initilizes the geometry transformation class AliITSgeom
124     // to values appropreate to this specific geometry. Now that
125     // the segmentation is part of AliITSgeom, the detector
126     // segmentations are also defined here.
127     // Inputs:
128     //   none.
129     // Outputs:
130     //   none.
131     // Return:
132     //   A pointer to a new properly inilized AliITSgeom class. If
133     //   pointer = 0 then failed to init.
134
135
136   AliITSVersion_t version = kvDefault;
137   Int_t minor = 0;
138   TDatime datetime;
139   TGeoVolume *itsV = gGeoManager->GetVolume("ITSV");
140   if(!itsV){
141     AliError("Can't find ITS volume ITSV, exiting - nothing done!");
142     return 0;
143   }// end if
144   const Char_t *title = itsV->GetTitle();
145   if(!ReadVersionString(title,(Int_t)strlen(title),version,minor,
146                         datetime))
147     Warning("UpdateInternalGeometry","Can't read title=%s\n",title);
148   SetTiming(kFALSE);
149   SetSegGeom(kFALSE);
150   SetDecoding(kFALSE);
151   AliITSgeom *geom = CreateAliITSgeom(version,minor);
152   AliDebug(1,"AliITSgeom object has been initialized from TGeo\n");
153   return geom;
154 }
155 //______________________________________________________________________
156 AliITSgeom* AliITSInitGeometry::CreateAliITSgeom(Int_t major,Int_t minor){
157     // Creates and Initilizes the geometry transformation class AliITSgeom
158     // to values appropreate to this specific geometry. Now that
159     // the segmentation is part of AliITSgeom, the detector
160     // segmentations are also defined here.
161     // Inputs:
162     //   Int_t major   major version, see AliITSVersion_t
163     //   Int_t minor   minor version
164     // Outputs:
165     //   none.
166     // Return:
167     //   A pointer to a new properly inilized AliITSgeom class. If
168     //   pointer = 0 then failed to init.
169
170     switch(major){
171     case kv11:
172         SetGeometryName("AliITSv11");
173         SetVersion(kv11,minor);
174         break;
175     case kvDefault:
176     default:
177         SetGeometryName("Undefined");
178         SetVersion(kvDefault,minor);
179         break;
180     } // end switch
181     AliITSgeom *geom = new AliITSgeom();
182     if(!InitAliITSgeom(geom)){ // Error initilization failed
183         delete geom;
184         geom = 0;
185     } // end if
186     return geom;
187 }
188 //______________________________________________________________________
189 Bool_t AliITSInitGeometry::InitAliITSgeom(AliITSgeom *geom){
190   // Initilizes the geometry transformation class AliITSgeom
191   // to values appropreate to this specific geometry. Now that
192   // the segmentation is part of AliITSgeom, the detector
193   // segmentations are also defined here.
194   // Inputs:
195   //   AliITSgeom *geom  A pointer to the AliITSgeom class
196   // Outputs:
197   //   AliITSgeom *geom  This pointer recreated and properly inilized.
198   // Return:
199   //   none.
200
201     if(!gGeoManager){
202         AliFatal("The geometry manager has not been initialized (e.g. "
203                  "TGeoManager::Import(\"geometry.root\")should be "
204                  "called in advance) - exit forced");
205         return kFALSE;
206     } // end if
207     switch(fMajorVersion) {
208     case kv11: {
209         return InitAliITSgeomV11(geom);
210     } break; // end case
211     case kvDefault: default: {
212         AliFatal("Undefined geometry");
213         return kFALSE;
214     } break; // end case
215     } // end switch
216     return kFALSE;
217 }
218 //______________________________________________________________________
219 void AliITSInitGeometry::TransposeTGeoHMatrix(TGeoHMatrix *m)const{
220     // Transpose the rotation matrix part of a TGeoHMatrix. This
221     // is needed because TGeo stores the transpose of the rotation
222     // matrix as compared to what AliITSgeomMatrix uses (and Geant3).
223     // Inputs:
224     //    TGeoHMatrix *m  The matrix to be transposed
225     // Outputs:
226     //    TGEoHMatrix *m  The transposed matrix
227     // Return:
228     //    none.
229     Int_t i;
230     Double_t r[9];
231
232     if(m==0) return; // no matrix to transpose.
233     for(i=0;i<9;i += 4) r[i] = m->GetRotationMatrix()[i]; // diagonals
234     r[1] = m->GetRotationMatrix()[3];
235     r[2] = m->GetRotationMatrix()[6];
236     r[3] = m->GetRotationMatrix()[1];
237     r[5] = m->GetRotationMatrix()[7];
238     r[6] = m->GetRotationMatrix()[2];
239     r[7] = m->GetRotationMatrix()[5];
240     m->SetRotation(r);
241     return;
242 }
243
244
245 //______________________________________________________________________
246 Bool_t AliITSInitGeometry::InitAliITSgeomV11(AliITSgeom *geom){
247   // Initilizes the geometry transformation class AliITSgeom
248   // Now that the segmentation is part of AliITSgeom, the detector
249   // segmentations are also defined here.
250   //
251   // Inputs:
252   //   AliITSgeom *geom  A pointer to the AliITSgeom class
253   // Outputs:
254   //   AliITSgeom *geom  This pointer recreated and properly inilized.
255   // LG
256
257   const Int_t kItype  = 0; // Type of transformation defined 0=> Geant
258   const Int_t klayers = 6; // number of layers in the ITS
259   const Int_t kladders[klayers]   = {20,40,14,22,34,38}; // Number of ladders
260   const Int_t kdetectors[klayers] = {4,4,6,8,22,25};// number of detector/lad
261   const AliITSDetector kIdet[6]   = {kSPD,kSPD,kSDD,kSDD,kSSD,kSSD};
262   const TString kPathbase = "/ALIC_1/ITSV_1/";
263
264   const char *pathSPDsens1, *pathSPDsens2;
265   pathSPDsens1="%sITSSPD_1/ITSSPDCarbonFiberSectorV_%d/ITSSPDSensitiveVirtualvolumeM0_1/ITSSPDlay1-Stave_%d/ITSSPDhalf-Stave%d_1/ITSSPDlay1-Ladder_%d/ITSSPDlay1-sensor_1";
266   pathSPDsens2="%sITSSPD_1/ITSSPDCarbonFiberSectorV_%d/ITSSPDSensitiveVirtualvolumeM0_1/ITSSPDlay2-Stave_%d/ITSSPDhalf-Stave%d_1/ITSSPDlay2-Ladder_%d/ITSSPDlay2-sensor_1";
267
268   const char *pathSDDsens1, *pathSDDsens2;
269   pathSDDsens1 = "%sITSsddLayer3_1/ITSsddLadd_%d/ITSsddSensor3_%d/ITSsddWafer3_%d/ITSsddSensitivL3_1";
270   pathSDDsens2 = "%sITSsddLayer4_1/ITSsddLadd_%d/ITSsddSensor4_%d/ITSsddWafer4_%d/ITSsddSensitivL4_1";
271
272   const char *pathSSDsens1, *pathSSDsens2;
273   pathSSDsens1 = "%sITSssdLayer5_1/ITSssdLay5Ladd_%d/ITSssdSensor5_%d/ITSssdSensitivL5_1";
274   pathSSDsens2 = "%sITSssdLayer6_1/ITSssdLay6Ladd_%d/ITSssdSensor6_%d/ITSssdSensitivL6_1";
275
276   const TString kNames[klayers] = {
277     pathSPDsens1, // lay=1
278     pathSPDsens2, // lay=2
279     pathSDDsens1, // lay=3
280     pathSDDsens2, // lay=4
281     pathSSDsens1, // lay=5
282     pathSSDsens2};// Lay=6
283   
284   Int_t mod,nmods=0, lay, lad, det, cpn0, cpn1, cpn2, cpnHS=1;
285   Double_t tran[3]={0.,0.,0.}, rot[10]={9*0.0,1.0};
286   TArrayD shapePar;
287   TString path, shapeName;
288   TGeoHMatrix matrix;
289   Bool_t initSeg[3]={kFALSE, kFALSE, kFALSE};
290   TStopwatch *time = 0x0;
291   if(fTiming) time = new TStopwatch();
292
293   if(fTiming) time->Start();
294   for(mod=0;mod<klayers;mod++) nmods += kladders[mod]*kdetectors[mod];
295   geom->Init(kItype,klayers,kladders,kdetectors,nmods);
296
297   for(mod=0; mod<nmods; mod++) {
298
299     DecodeDetectorLayers(mod,lay,lad,det);
300     geom->CreateMatrix(mod,lay,lad,det,kIdet[lay-1],tran,rot);
301     RecodeDetector(mod,cpn0,cpn1,cpn2);
302
303     if (kIdet[lay-1]==kSPD) { // we need 1 more copy number because of the half-stave
304       if (det<3) cpnHS = 0; else cpnHS = 1;
305       path.Form(kNames[lay-1].Data(),kPathbase.Data(),cpn0,cpn1,cpnHS,cpn2);
306     } else {
307       path.Form(kNames[lay-1].Data(),kPathbase.Data(),cpn0,cpn1,cpn2);
308     };
309
310     geom->GetGeomMatrix(mod)->SetPath(path);
311     GetTransformation(path.Data(),matrix);
312     geom->SetTrans(mod,matrix.GetTranslation());
313     TransposeTGeoHMatrix(&matrix); //Transpose TGeo's rotation matrixes
314     geom->SetRotMatrix(mod,matrix.GetRotationMatrix());
315     if(initSeg[kIdet[lay-1]]) continue;
316     GetShape(path,shapeName,shapePar);
317     if(shapeName.CompareTo("BOX")){
318       Error("InitITSgeom","Geometry changed without proper code update"
319             "or error in reading geometry. Shape is not BOX.");
320       return kFALSE;
321     } // end if
322   } // end for module
323
324   if(fTiming){
325     time->Stop();
326     time->Print();
327     delete time;
328   } // end if
329   return kTRUE;
330 }
331
332 //_______________________________________________________________________
333 Bool_t AliITSInitGeometry::GetTransformation(const TString &volumePath,
334                                              TGeoHMatrix &mat){
335     // Returns the Transformation matrix between the volume specified
336     // by the path volumePath and the Top or mater volume. The format
337     // of the path volumePath is as follows (assuming ALIC is the Top volume)
338     // "/ALIC_1/DDIP_1/S05I_2/S05H_1/S05G_3". Here ALIC is the top most
339     // or master volume which has only 1 instance of. Of all of the daughter
340     // volumes of ALICE, DDIP volume copy #1 is indicated. Similarly for
341     // the daughter volume of DDIP is S05I copy #2 and so on.
342     // Inputs:
343     //   TString& volumePath  The volume path to the specific volume
344     //                        for which you want the matrix. Volume name
345     //                        hierarchy is separated by "/" while the
346     //                        copy number is appended using a "_".
347     // Outputs:
348     //  TGeoHMatrix &mat      A matrix with its values set to those
349     //                        appropriate to the Local to Master transformation
350     // Return:
351     //   A logical value if kFALSE then an error occurred and no change to
352     //   mat was made.
353
354     // We have to preserve the modeler state
355
356     // Preserve the modeler state.
357     gGeoManager->PushPath();
358     if (!gGeoManager->cd(volumePath.Data())) {
359       gGeoManager->PopPath();
360       Error("GetTransformation","Error in cd-ing to %s",volumePath.Data());
361       return kFALSE;
362     } // end if !gGeoManager
363     mat = *gGeoManager->GetCurrentMatrix();
364     // Retstore the modeler state.
365     gGeoManager->PopPath();
366     return kTRUE;
367 }
368 //______________________________________________________________________
369 Bool_t AliITSInitGeometry::GetShape(const TString &volumePath,
370                                     TString &shapeType,TArrayD &par){
371     // Returns the shape and its parameters for the volume specified
372     // by volumeName.
373     // Inputs:
374     //   TString& volumeName  The volume name
375     // Outputs:
376     //   TString &shapeType   Shape type
377     //   TArrayD &par         A TArrayD of parameters with all of the
378     //                        parameters of the specified shape.
379     // Return:
380     //   A logical indicating whether there was an error in getting this
381     //   information
382     Int_t npar;
383     gGeoManager->PushPath();
384     if (!gGeoManager->cd(volumePath.Data())) {
385         gGeoManager->PopPath();
386         return kFALSE;
387     }
388     TGeoVolume * vol = gGeoManager->GetCurrentVolume();
389     gGeoManager->PopPath();
390     if (!vol) return kFALSE;
391     TGeoShape *shape = vol->GetShape();
392     TClass *classType = shape->IsA();
393     if (classType==TGeoBBox::Class()) {
394         shapeType = "BOX";
395         npar = 3;
396         par.Set(npar);
397         TGeoBBox *box = (TGeoBBox*)shape;
398         par.AddAt(box->GetDX(),0);
399         par.AddAt(box->GetDY(),1);
400         par.AddAt(box->GetDZ(),2);
401         return kTRUE;
402     } // end if
403     if (classType==TGeoTrd1::Class()) {
404         shapeType = "TRD1";
405         npar = 4;
406         par.Set(npar);
407         TGeoTrd1 *trd1 = (TGeoTrd1*)shape;
408         par.AddAt(trd1->GetDx1(),0);
409         par.AddAt(trd1->GetDx2(),1);
410         par.AddAt(trd1->GetDy(), 2);
411         par.AddAt(trd1->GetDz(), 3);
412         return kTRUE;
413     } // end if
414     if (classType==TGeoTrd2::Class()) {
415         shapeType = "TRD2";
416         npar = 5;
417         par.Set(npar);
418         TGeoTrd2 *trd2 = (TGeoTrd2*)shape;
419         par.AddAt(trd2->GetDx1(),0);
420         par.AddAt(trd2->GetDx2(),1);
421         par.AddAt(trd2->GetDy1(),2);
422         par.AddAt(trd2->GetDy2(),3);
423         par.AddAt(trd2->GetDz(), 4);
424         return kTRUE;
425     } // end if
426     if (classType==TGeoTrap::Class()) {
427         shapeType = "TRAP";
428         npar = 11;
429         par.Set(npar);
430         TGeoTrap *trap = (TGeoTrap*)shape;
431         Double_t tth = TMath::Tan(trap->GetTheta()*TMath::DegToRad());
432         par.AddAt(trap->GetDz(),0);
433         par.AddAt(tth*TMath::Cos(trap->GetPhi()*TMath::DegToRad()),1);
434         par.AddAt(tth*TMath::Sin(trap->GetPhi()*TMath::DegToRad()),2);
435         par.AddAt(trap->GetH1(),3);
436         par.AddAt(trap->GetBl1(),4);
437         par.AddAt(trap->GetTl1(),5);
438         par.AddAt(TMath::Tan(trap->GetAlpha1()*TMath::DegToRad()),6);
439         par.AddAt(trap->GetH2(),7);
440         par.AddAt(trap->GetBl2(),8);
441         par.AddAt(trap->GetTl2(),9);
442         par.AddAt(TMath::Tan(trap->GetAlpha2()*TMath::DegToRad()),10);
443         return kTRUE;
444     } // end if
445     if (classType==TGeoTube::Class()) {
446         shapeType = "TUBE";
447         npar = 3;
448         par.Set(npar);
449         TGeoTube *tube = (TGeoTube*)shape;
450         par.AddAt(tube->GetRmin(),0);
451         par.AddAt(tube->GetRmax(),1);
452         par.AddAt(tube->GetDz(),2);
453         return kTRUE;
454     } // end if
455     if (classType==TGeoTubeSeg::Class()) {
456         shapeType = "TUBS";
457         npar = 5;
458         par.Set(npar);
459         TGeoTubeSeg *tubs = (TGeoTubeSeg*)shape;
460         par.AddAt(tubs->GetRmin(),0);
461         par.AddAt(tubs->GetRmax(),1);
462         par.AddAt(tubs->GetDz(),2);
463         par.AddAt(tubs->GetPhi1(),3);
464         par.AddAt(tubs->GetPhi2(),4);
465         return kTRUE;
466     } // end if
467     if (classType==TGeoCone::Class()) {
468         shapeType = "CONE";
469         npar = 5;
470         par.Set(npar);
471         TGeoCone *cone = (TGeoCone*)shape;
472         par.AddAt(cone->GetDz(),0);
473         par.AddAt(cone->GetRmin1(),1);
474         par.AddAt(cone->GetRmax1(),2);
475         par.AddAt(cone->GetRmin2(),3);
476         par.AddAt(cone->GetRmax2(),4);
477         return kTRUE;
478     } // end if
479     if (classType==TGeoConeSeg::Class()) {
480         shapeType = "CONS";
481         npar = 7;
482         par.Set(npar);
483         TGeoConeSeg *cons = (TGeoConeSeg*)shape;
484         par.AddAt(cons->GetDz(),0);
485         par.AddAt(cons->GetRmin1(),1);
486         par.AddAt(cons->GetRmax1(),2);
487         par.AddAt(cons->GetRmin2(),3);
488         par.AddAt(cons->GetRmax2(),4);
489         par.AddAt(cons->GetPhi1(),5);
490         par.AddAt(cons->GetPhi2(),6);
491         return kTRUE;
492     } // end if
493     if (classType==TGeoSphere::Class()) {
494         shapeType = "SPHE";
495         npar = 6;
496         par.Set(npar);
497         
498         TGeoSphere *sphe = (TGeoSphere*)shape;
499         par.AddAt(sphe->GetRmin(),0);
500         par.AddAt(sphe->GetRmax(),1);
501         par.AddAt(sphe->GetTheta1(),2);
502         par.AddAt(sphe->GetTheta2(),3);
503         par.AddAt(sphe->GetPhi1(),4);
504         par.AddAt(sphe->GetPhi2(),5);
505         return kTRUE;
506     } // end if
507     if (classType==TGeoPara::Class()) {
508         shapeType = "PARA";
509         npar = 6;
510         par.Set(npar);
511         TGeoPara *para = (TGeoPara*)shape;
512         par.AddAt(para->GetX(),0);
513         par.AddAt(para->GetY(),1);
514         par.AddAt(para->GetZ(),2);
515         par.AddAt(para->GetTxy(),3);
516         par.AddAt(para->GetTxz(),4);
517         par.AddAt(para->GetTyz(),5);
518         return kTRUE;
519     } // end if
520     if (classType==TGeoPgon::Class()) {
521         shapeType = "PGON";
522         TGeoPgon *pgon = (TGeoPgon*)shape;
523         Int_t nz = pgon->GetNz();
524         const Double_t *rmin = pgon->GetRmin();
525         const Double_t *rmax = pgon->GetRmax();
526         const Double_t *z = pgon->GetZ();
527         npar = 4 + 3*nz;
528         par.Set(npar);
529         par.AddAt(pgon->GetPhi1(),0);
530         par.AddAt(pgon->GetDphi(),1);
531         par.AddAt(pgon->GetNedges(),2);
532         par.AddAt(pgon->GetNz(),3);
533         for (Int_t i=0; i<nz; i++) {
534             par.AddAt(z[i], 4+3*i);
535             par.AddAt(rmin[i], 4+3*i+1);
536             par.AddAt(rmax[i], 4+3*i+2);
537         }
538         return kTRUE;
539     } // end if
540     if (classType==TGeoPcon::Class()) {
541         shapeType = "PCON";
542         TGeoPcon *pcon = (TGeoPcon*)shape;
543         Int_t nz = pcon->GetNz();
544         const Double_t *rmin = pcon->GetRmin();
545         const Double_t *rmax = pcon->GetRmax();
546         const Double_t *z = pcon->GetZ();
547         npar = 3 + 3*nz;
548         par.Set(npar);
549         par.AddAt(pcon->GetPhi1(),0);
550         par.AddAt(pcon->GetDphi(),1);
551         par.AddAt(pcon->GetNz(),2);
552         for (Int_t i=0; i<nz; i++) {
553             par.AddAt(z[i], 3+3*i);
554             
555             par.AddAt(rmin[i], 3+3*i+1);
556             par.AddAt(rmax[i], 3+3*i+2);
557         }
558         return kTRUE;
559     } // end if
560     if (classType==TGeoEltu::Class()) {
561         shapeType = "ELTU";
562         npar = 3;
563         par.Set(npar);
564         TGeoEltu *eltu = (TGeoEltu*)shape;
565         par.AddAt(eltu->GetA(),0);
566         par.AddAt(eltu->GetB(),1);
567         par.AddAt(eltu->GetDz(),2);
568         return kTRUE;
569     } // end if
570     if (classType==TGeoHype::Class()) {
571         shapeType = "HYPE";
572         npar = 5;
573         par.Set(npar);
574         TGeoHype *hype = (TGeoHype*)shape;
575         par.AddAt(TMath::Sqrt(hype->RadiusHypeSq(0.,kTRUE)),0);
576         par.AddAt(TMath::Sqrt(hype->RadiusHypeSq(0.,kFALSE)),1);
577         par.AddAt(hype->GetDZ(),2);
578         par.AddAt(hype->GetStIn(),3);
579         par.AddAt(hype->GetStOut(),4);
580         return kTRUE;
581     } // end if
582     if (classType==TGeoGtra::Class()) {
583         shapeType = "GTRA";
584         npar = 12;
585         par.Set(npar);
586         TGeoGtra *trap = (TGeoGtra*)shape;
587         Double_t tth = TMath::Tan(trap->GetTheta()*TMath::DegToRad());
588         par.AddAt(trap->GetDz(),0);
589         par.AddAt(tth*TMath::Cos(trap->GetPhi()*TMath::DegToRad()),1);
590         par.AddAt(tth*TMath::Sin(trap->GetPhi()*TMath::DegToRad()),2);
591         par.AddAt(trap->GetH1(),3);
592         par.AddAt(trap->GetBl1(),4);
593         par.AddAt(trap->GetTl1(),5);
594         par.AddAt(TMath::Tan(trap->GetAlpha1()*TMath::DegToRad()),6);
595         par.AddAt(trap->GetH2(),7);
596         par.AddAt(trap->GetBl2(),8);
597         par.AddAt(trap->GetTl2(),9);
598         par.AddAt(TMath::Tan(trap->GetAlpha2()*TMath::DegToRad()),10);
599         par.AddAt(trap->GetTwistAngle(),11);
600         return kTRUE;
601     } // end if
602     if (classType==TGeoCtub::Class()) {
603         shapeType = "CTUB";
604         npar = 11;
605         par.Set(npar);
606         TGeoCtub *ctub = (TGeoCtub*)shape;
607         const Double_t *lx = ctub->GetNlow();
608         const Double_t *tx = ctub->GetNhigh();
609         par.AddAt(ctub->GetRmin(),0);
610         par.AddAt(ctub->GetRmax(),1);
611         par.AddAt(ctub->GetDz(),2);
612         par.AddAt(ctub->GetPhi1(),3);
613         par.AddAt(ctub->GetPhi2(),4);
614         par.AddAt(lx[0],5);
615         par.AddAt(lx[1],6);
616         par.AddAt(lx[2],7);
617         par.AddAt(tx[0],8);
618         par.AddAt(tx[1],9);
619         par.AddAt(tx[2],10);
620         return kTRUE;
621     } // end if
622     Error("GetShape","Getting shape parameters for shape %s not implemented",
623           shape->ClassName());
624     shapeType = "Unknown";
625     return kFALSE;
626 }
627 //______________________________________________________________________
628 void AliITSInitGeometry::DecodeDetector(
629     Int_t &mod,Int_t layer,Int_t cpn0,Int_t cpn1,Int_t cpn2) const {
630     // decode geometry into detector module number. There are two decoding
631     // Scheams. Old which does not follow the ALICE coordinate system
632     // requirements, and New which dose.
633     // Inputs:
634     //    Int_t layer    The ITS layer
635     //    Int_t cpn0     The lowest copy number
636     //    Int_t cpn1     The middle copy number
637     //    Int_t cpn2     the highest copy number
638     // Output:
639     //    Int_t &mod     The module number assoicated with this set
640     //                   of copy numbers.
641     // Return:
642     //    none.
643
644     // This is a FIXED switch yard function. I (Bjorn Nilsen) Don't 
645     // like them but I see not better way for the moment.
646     switch (fMajorVersion){
647     case kvDefault:{
648         Error("DecodeDetector","Major version = kvDefault, not supported");
649     }break;
650     case kv11:{
651         return DecodeDetectorv11(mod,layer,cpn0,cpn1,cpn2);
652     }break;
653     default:{
654         Error("DecodeDetector","Major version = %d, not supported",
655               (Int_t)fMajorVersion);
656         return;
657     }break;
658     } // end switch
659     return;
660 }
661 //______________________________________________________________________
662 void AliITSInitGeometry::RecodeDetector(Int_t mod,Int_t &cpn0,
663                                         Int_t &cpn1,Int_t &cpn2){
664     // decode geometry into detector module number. There are two decoding
665     // Scheams. Old which does not follow the ALICE coordinate system
666     // requirements, and New which dose.
667     // Inputs:
668     //    Int_t mod      The module number assoicated with this set
669     //                   of copy numbers.
670     // Output:
671     //    Int_t cpn0     The lowest copy number
672     //    Int_t cpn1     The middle copy number
673     //    Int_t cpn2     the highest copy number
674     // Return:
675     //    none.
676
677     // This is a FIXED switch yard function. I (Bjorn Nilsen) Don't 
678     // like them but I see not better way for the moment.
679     switch (fMajorVersion){
680     case kvDefault:{
681         Error("RecodeDetector","Major version = kvDefault, not supported");
682         return;
683     }
684     case kv11:{
685         return RecodeDetectorv11(mod,cpn0,cpn1,cpn2);
686     }break;
687     default:{
688         Error("RecodeDetector","Major version = %d, not supported",
689               (Int_t)fMajorVersion);
690         return;
691     }break;
692     } // end switch
693     return;
694 }
695 //______________________________________________________________________
696 void AliITSInitGeometry::DecodeDetectorLayers(Int_t mod,Int_t &layer,
697                                               Int_t &lad,Int_t &det){
698     // decode geometry into detector module number. There are two decoding
699     // Scheams. Old which does not follow the ALICE coordinate system
700     // requirements, and New which dose. Note, this use of layer ladder
701     // and detector numbers are strictly for internal use of this
702     // specific code. They do not represent the "standard" layer ladder
703     // or detector numbering except in a very old and obsoleate sence.
704     // Inputs:
705     //    Int_t mod      The module number assoicated with this set
706     //                   of copy numbers.
707     // Output:
708     //    Int_t lay     The layer number
709     //    Int_t lad     The ladder number
710     //    Int_t det     the dettector number
711     // Return:
712     //    none.
713
714     // This is a FIXED switch yard function. I (Bjorn Nilsen) Don't 
715     // like them but I see not better way for the moment.
716     switch (fMajorVersion) {
717     case kvDefault:{
718         Error("DecodeDetectorLayers",
719               "Major version = kvDefault, not supported");
720         return;
721     }break;
722     case kv11:{
723         return DecodeDetectorLayersv11(mod,layer,lad,det);
724     }break;
725     default:{
726         Error("DecodeDetectorLayers","Major version = %d, not supported",
727               (Int_t)fMajorVersion);
728         return;
729     }break;
730     } // end switch
731     return;
732 }
733
734 //______________________________________________________________________
735 void AliITSInitGeometry::DecodeDetectorv11(Int_t &mod,Int_t layer,
736                                  Int_t cpn0,Int_t cpn1,Int_t cpn2) const {
737     // decode geometry into detector module number
738     // Inputs:
739     //    Int_t layer    The ITS layer
740     //    Int_t cpn0     The lowest copy number
741     //    Int_t cpn1     The middle copy number
742     //    Int_t cpn2     the highest copy number
743     // Output:
744     //    Int_t &mod     The module number assoicated with this set
745     //                   of copy numbers.
746     // Return:
747     //    none.
748   const Int_t kDetPerLadderSPD[2]={2,4};
749   const Int_t kDetPerLadder[6]={4,4,6,8,22,25};
750   const Int_t kLadPerLayer[6]={20,40,14,22,34,38};
751   Int_t lad=-1,det=-1;
752   
753   switch(layer) {
754   case 1: case 2:{
755     lad = cpn1+kDetPerLadderSPD[layer-1]*(cpn0-1);
756     det = cpn2;
757   } break;
758   case 3: case 4:{
759     lad = cpn0+1;
760     det = cpn1+1;
761   } break;
762   case 5: case 6:{
763     lad = cpn0+1;
764     det = cpn1+1;
765   } break;
766   default:{
767   } break;
768   } // end switch
769   mod = 0;
770   for(Int_t i=0;i<layer-1;i++) mod += kLadPerLayer[i]*kDetPerLadder[i];
771   mod += kDetPerLadder[layer-1]*(lad-1)+det-1;// module start at zero.
772   return;
773 }
774
775
776 //______________________________________________________________________
777 void AliITSInitGeometry::RecodeDetectorv11(Int_t mod,Int_t &cpn0,
778                                            Int_t &cpn1,Int_t &cpn2) {
779     // decode geometry into detector module number using the new decoding
780     // Scheme.
781     // Inputs:
782     //    Int_t mod      The module number assoicated with this set
783     //                   of copy numbers.
784     // Output:
785     //    Int_t cpn0     The lowest copy number  (SPD sector or SDD/SSD ladder)
786     //    Int_t cpn1     The middle copy number  (SPD stave or SDD/SSD module)
787     //    Int_t cpn2     the highest copy number (SPD ladder or 1 for SDD/SSD)
788     // Return:
789     //    none.
790     const Int_t kDetPerLadderSPD[2]={2,4};
791     Int_t lay,lad,det;
792
793     DecodeDetectorLayersv11(mod,lay,lad,det);
794     if (lay<3) { // SPD
795         cpn2 = det;     // Detector 1-4
796         cpn0 = (lad+kDetPerLadderSPD[lay-1]-1)/kDetPerLadderSPD[lay-1];
797         cpn1 = (lad+kDetPerLadderSPD[lay-1]-1)%kDetPerLadderSPD[lay-1] + 1;
798     } else { // SDD and SSD
799         cpn2 = 1;
800         cpn1 = det;
801         cpn0 = lad;
802         if (lay<5) { // SDD
803           cpn1--;
804           cpn0--;
805         } else { //SSD
806           cpn1--;
807           cpn0--;
808         } // end if Lay<5/else
809     } // end if lay<3/else
810
811 }
812
813
814 //______________________________________________________________________
815 void AliITSInitGeometry::DecodeDetectorLayersv11(Int_t mod,Int_t &lay,
816                                                  Int_t &lad,Int_t &det) {
817
818   // decode module number into detector indices for v11
819   // mod starts from 0
820   // lay, lad, det start from 1
821
822   // Inputs:
823   //    Int_t mod      The module number associated with this set
824   //                   of copy numbers.
825   // Output:
826   //    Int_t lay     The layer number
827   //    Int_t lad     The ladder number
828   //    Int_t det     the dettector number
829
830   const Int_t kDetPerLadder[6] = {4,4,6,8,22,25};
831   const Int_t kLadPerLayer[6]  = {20,40,14,22,34,38};
832   
833   Int_t mod2 = 0;
834   lay  = 0;
835   
836   do {
837     mod2 += kLadPerLayer[lay]*kDetPerLadder[lay];
838     lay++;
839   } while(mod2<=mod); // end while
840   if(lay>6) Error("DecodeDetectorLayers","lay=%d>6",lay);
841
842   mod2 = kLadPerLayer[lay-1]*kDetPerLadder[lay-1] - mod2+mod;
843   lad = mod2/kDetPerLadder[lay-1];
844
845   if(lad>=kLadPerLayer[lay-1]||lad<0) Error("DecodeDetectorLayers",
846                                       "lad=%d not in the correct range",lad);
847   det = (mod2 - lad*kDetPerLadder[lay-1])+1;
848   if(det>kDetPerLadder[lay-1]||det<1) Error("DecodeDetectorLayers",
849                                       "det=%d not in the correct range",det);
850   lad++;
851 }
852
853 //______________________________________________________________________
854 Bool_t AliITSInitGeometry::WriteVersionString(Char_t *str,Int_t length,
855                         AliITSVersion_t maj,Int_t min,
856                         const Char_t *cvsDate,const Char_t *cvsRevision)const{
857     // fills the string str with the major and minor version number
858     // Inputs:
859     //   Char_t *str          The character string to hold the major 
860     //                        and minor version numbers in
861     //   Int_t  length        The maximum number of characters which 
862     //                        can be accomidated by this string. 
863     //                        str[length-1] must exist and will be set to zero
864     //   AliITSVersion_t maj  The major number
865     //   Int_t           min  The minor number
866     //   Char_t *cvsDate      The date string from cvs
867     //   Char_t *cvsRevision  The Revision string from cvs
868     // Outputs:
869     //   Char_t *str          The character string holding the major and minor
870     //                        version numbers. str[length-1] must exist
871     //                        and will be set to zero
872     // Return:
873     //   kTRUE if no errors
874     Char_t cvslikedate[30];
875     Int_t i,n,cvsDateLength,cvsRevisionLength;
876
877     cvsDateLength = (Int_t)strlen(cvsDate);
878     if(cvsDateLength>30){ // svn string, make a cvs like string
879         i=0;n=0;
880         do{
881             cvslikedate[i] = cvsDate[i];
882             if(cvsDate[i]=='+' || cvsDate[i++]=='-'){
883                 n++; // count number of -
884                 cvslikedate[i-1] = '/'; // replace -'s by /'s.
885             } // end if
886         } while(n<3&&i<30); // once additonal - of time zone reach exit
887         cvslikedate[i-1] = '$'; // put $ at end then zero.
888         for(;i<30;i++) cvslikedate[i]=0;// i starts wher do loop left off.
889     }else{
890         for(i=0;i<cvsDateLength&&i<30;i++) cvslikedate[i]=cvsDate[i];
891     }// end if
892     cvsDateLength = (Int_t)strlen(cvslikedate);
893     cvsRevisionLength = (Int_t)strlen(cvsRevision);
894     i = (Int_t)maj;
895     n = 50+(Int_t)(TMath::Log10(TMath::Abs((Double_t)i)))+1+
896         (Int_t)(TMath::Log10(TMath::Abs((Double_t)min)))+1
897         +cvsDateLength-6+cvsRevisionLength-10;
898     if(GetDebug()>1) printf("AliITSInitGeometry::WriteVersionString:"
899                         "length=%d major=%d minor=%d cvsDate=%s[%d] "
900                         "cvsRevision=%s[%d] n=%d\n",length,i,min,cvslikedate,
901                         cvsDateLength,cvsRevision,cvsRevisionLength,n);
902     if(i<0) n++;
903     if(min<0) n++;
904     if(length<n){// not enough space to write in output string.
905         Warning("WriteVersionString","Output string not long enough "
906                 "lenght=%d must be at least %d long\n",length,n);
907         return kFALSE;
908     } // end if length<n
909     char *cvsrevision = new char[cvsRevisionLength-10];
910     char *cvsdate = new char[cvsDateLength-6];
911     for(i=0;i<cvsRevisionLength-10;i++)
912         if(10+i<cvsRevisionLength-1)
913             cvsrevision[i] = cvsRevision[10+i]; else cvsrevision[i] = 0;
914     for(i=0;i<cvsDateLength-6;i++) if(6+i<cvsDateLength-1)
915         cvsdate[i] = cvslikedate[6+i]; else cvsdate[i] = 0;
916     for(i=0;i<length;i++) str[i] = 0; // zero it out for now.
917     i = (Int_t)maj;
918     snprintf(str,length-1,"Major Version= %d Minor Version= %d Revision: %s Date: %s",i,min,cvsrevision,cvsdate);
919     /* this gives compilation warnings on some compilers: descriptor zu
920     if(GetDebug()>1)printf("AliITSInitGeometry::WriteVersionString: "
921                        "n=%d str=%s revision[%zu] date[%zu]\n",
922                        n,str,strlen(cvsrevision),strlen(cvsdate));
923     */
924     delete[] cvsrevision;
925     delete[] cvsdate;
926     return kTRUE;
927 }
928 //______________________________________________________________________
929 Bool_t AliITSInitGeometry::ReadVersionString(const Char_t *str,Int_t length,
930                                              AliITSVersion_t &maj,Int_t &min,
931                                              TDatime &dt)const{
932     // fills the string str with the major and minor version number
933     // Inputs:
934     //   Char_t *str   The character string to holding the major and minor
935     //                 version numbers in
936     //   Int_t  length The maximum number of characters which can be
937     //                 accomidated by this string. str[length-1] must exist
938     // Outputs:
939     //   Char_t *str   The character string holding the major and minor
940     //                 version numbers unchanged. str[length-1] must exist.
941     //   AliITSVersion_t maj  The major number
942     //   Int_t           min  The minor number
943     //   TDatime         dt   The date and time of the cvs commit
944     // Return:
945     //   kTRUE if no errors
946     Bool_t ok;
947     Char_t cvsRevision[10],cvsDate[11],cvsTime[9];
948     Int_t i,m,n=strlen(str),year,month,day,hours,minuits,seconds;
949     memset(cvsRevision,0,10*sizeof(Char_t));
950     memset(cvsDate,0,11*sizeof(Char_t));    
951     memset(cvsTime,0,9*sizeof(Char_t));
952
953     if(GetDebug()>1)printf("AliITSInitGeometry::ReadVersionString:"
954                        "str=%s length=%d\n",
955                        str,length);
956     if(n<35) return kFALSE; // not enough space for numbers
957     m = sscanf(str,"Major Version= %d  Minor Version= %d Revision: %9s "
958                "Date: %10s %8s",&i,&min,cvsRevision,cvsDate,cvsTime);
959
960     // v11Hybrid geometry is treated as a v11 geometry
961     if(i == 110)i=11;
962
963     ok = m==5;
964     if(!ok) return !ok;
965     m = sscanf(cvsDate,"%d/%d/%d",&year,&month,&day);
966     ok = m==3;
967     if(!ok) return !ok;
968     m = sscanf(cvsTime,"%d:%d:%d",&hours,&minuits,&seconds);
969     ok = m==3;
970     if(!ok) return !ok;
971     dt.Set(year,month,day,hours,minuits,seconds);
972     if(GetDebug()>1)printf("AliITSInitGeometry::ReadVersionString: i=%d "
973                      "min=%d cvsRevision=%s cvsDate=%s cvsTime=%s m=%d\n",
974                        i,min,cvsRevision,cvsDate,cvsTime,m);
975     if(GetDebug()>1)printf("AliITSInitGeometry::ReadVersionString: year=%d"
976                        " month=%d day=%d hours=%d minuits=%d seconds=%d\n",
977                        year,month,day,hours,minuits,seconds);
978     switch (i){
979     case kv11:{
980         maj = kv11;
981     } break;
982     default:{
983         maj = kvDefault;
984     } break;
985     } // end switch
986     return ok;
987 }