f521af05ddadecbe53773e2b0132204ad9a8255f
[u/mrichter/AliRoot.git] / ITS / AliITSvtest.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 $Log$
18 Revision 1.7  2001/02/09 20:06:26  nilsen
19 Fixed bug in distructor. Can't distroy fixxed length arrays. Thanks Peter.
20
21 Revision 1.6  2001/02/09 00:05:31  nilsen
22 Added fMajor/MinorVersion variables and made other changes to better make
23 use of the new code changes in AliITSgeom related classes.
24
25 Revision 1.5  2001/01/30 09:23:14  hristov
26 Streamers removed (R.Brun)
27
28 Revision 1.4  2001/01/18 06:25:09  barbera
29 ITS geometry using test Euclid files
30
31 Revision 1.1.2.8  2000/10/05 20:28:18  nilsen
32 Now using root generated streamer function.
33
34 Revision 1.1.2.7  2000/07/31 13:51:22  barbera
35 Updated from the release
36
37 Revision 1.2  2000/07/10 16:07:19  fca
38 Release version of ITS code
39
40 Revision 1.1.2.2  2000/03/02 21:53:36  nilsen
41 to make it compatable with the changes in AliRun/AliModule.
42
43 Revision 1.1.2.1  2000/01/12 20:19:03  nilsen
44         The changes made with this latest inclusion of code is very large.
45 Many of the new files were added just in December when P. Cerello added his
46 SDD simulations to the distrobutions. Also added are some file of P. Skowronski
47 for SSD cluster finding and ghost RecPoints. None of this "new" code has been
48 proporly tested. Other code new to this cvs repository is explained in the
49 ITS Off-line web page. In general the changes are too large to give a resonable
50 discription of them but probably should be taken as the starting point for
51 the developement branch (ITS-working).
52     B. S. Nilsen
53
54 Revision 1.13  1999/10/16 19:49:00  BSN
55 $Name$
56 $Author$
57 $Id$
58 */
59
60 ///////////////////////////////////////////////////////////////////////////////
61 //                                                                           //
62 //  Inner Traking System version Test                                        //
63 //  This class contains the base procedures for the Inner Tracking System    //
64 //                                                                           //
65 // Authors: R. Barbera, B. S. Nilsen.                                        //
66 // version  Test                                                             //
67 // Created October 16 1999.                                                  //
68 //                                                                           //
69 ///////////////////////////////////////////////////////////////////////////////
70 #include <stdio.h>
71 #include <stdlib.h>
72 #include <TMath.h>
73 #include <TGeometry.h>
74 #include <TNode.h>
75 #include <TTUBE.h>
76 #include <TFile.h>    // only required for Tracking function?
77 #include <TCanvas.h>
78 #include <TObjArray.h>
79 #include <TObjString.h>
80 #include <TClonesArray.h>
81 #include <TLorentzVector.h>
82 #include <TBRIK.h>
83 #include <TSystem.h>
84
85 #include "AliMC.h"
86 #include "AliRun.h"
87 #include "../TGeant3/TGeant3.h"
88 #include "AliITSGeant3Geometry.h"
89 #include "AliITShit.h"
90 #include "AliITS.h"
91 #include "AliITSvtest.h"
92 #include "AliITSgeom.h"
93 #include "AliITSgeomSPD.h"
94 #include "AliITSgeomSDD.h"
95 #include "AliITSgeomSSD.h"
96
97 ClassImp(AliITSvtest)
98  
99 //_____________________________________________________________________________
100 AliITSvtest::AliITSvtest() {
101     // Standard constructor for the ITS
102     Int_t i;
103
104     fIdN    = 0;
105     fIdName = 0;
106     fIdSens = 0;
107     fEuclidOut    = kFALSE; // Don't write Euclide file
108     fGeomDetOut   = kFALSE; // Don't write .det file
109     fGeomDetIn    = kTRUE; // Read .det file
110     fMajorVersion = IsVersion();
111     fMinorVersion = -1;
112     for(i=0;i<60;i++) fRead[i] = '\0';
113     for(i=0;i<60;i++) fWrite[i] = '\0';
114     for(i=0;i<60;i++) fEuclidGeomDet[i] = '\0';
115 }
116 //____________________________________________________________________________
117 AliITSvtest::AliITSvtest(const AliITSvtest &source){
118 ////////////////////////////////////////////////////////////////////////
119 //     Copy Constructor for ITS test version.
120 ////////////////////////////////////////////////////////////////////////
121     if(&source == this) return;
122     Warning("Copy Constructor","Not allowed to copy AliITSvtest");
123     return;
124 }
125 //_____________________________________________________________________________
126 AliITSvtest& AliITSvtest::operator=(const AliITSvtest &source){
127 ////////////////////////////////////////////////////////////////////////
128 //    Assignment operator for the ITS version 1.
129 ////////////////////////////////////////////////////////////////////////
130         if(&source == this) return *this;
131         Warning("= operator","Not allowed to copy AliITSvtest");
132         return *this;
133 }
134 //_____________________________________________________________________________
135 AliITSvtest::~AliITSvtest() {
136     // Standard destructor for the ITS
137 }
138 //_____________________________________________________________________________
139 AliITSvtest::AliITSvtest(const char *fileeuc,const char *filetme,
140                          const char *name, const char *title) 
141     : AliITS(name, title){
142     //
143     // Standard constructor for the ITS
144     //
145     fIdN    = 6;
146     fIdName    = new TString[fIdN];
147     fIdName[0] = "ITS1";
148     fIdName[1] = "ITS2";
149     fIdName[2] = "ITS3";
150     fIdName[3] = "ITS4";
151     fIdName[4] = "ITS5";
152     fIdName[5] = "ITS6";
153     fIdSens    = new Int_t[fIdN];
154     for (Int_t i=0;i<fIdN;i++) fIdSens[i] = 0;
155     fMajorVersion = IsVersion();
156     fMinorVersion = 1;
157     fEuclidOut    = kFALSE; // Don't write Euclide file
158     fGeomDetOut   = kFALSE; // Don't write .det file
159     fGeomDetIn    = kTRUE; // Read .det file
160
161     fEuclidMaterial = filetme;
162     fEuclidGeometry = fileeuc;
163     strncpy(fEuclidGeomDet,"$ALICE_ROOT/ITS/ITSgeometry_PPR.det",60);
164     strncpy(fRead,fEuclidGeomDet,60);
165     strncpy(fWrite,fEuclidGeomDet,60);
166 //  The .det file for the geometry must have the same name as fileeuc with
167 //  .euc replaced by .det.
168 }
169
170  
171 //_____________________________________________________________________________
172 void AliITSvtest::CreateMaterials(){
173   //
174   // Read materials for the ITS
175   //
176     char *filtmp;
177 //
178   filtmp = gSystem->ExpandPathName(fEuclidMaterial.Data());
179 //  FILE *file = fopen(fEuclidMaterial.Data(),"r");
180   FILE *file = fopen(filtmp,"r");
181   if(file) {
182     fclose(file);
183 //    ReadEuclidMedia(fEuclidMaterial.Data(),this);
184     ReadEuclidMedia(filtmp);
185   } else {
186     Error("CreateMaterials"," THE MEDIA FILE %s DOES NOT EXIST !",
187 //        fEuclidMaterial.Data());
188           filtmp);
189     exit(1);
190   } // end if(file)
191 }
192
193 //_____________________________________________________________________________
194 void AliITSvtest::CreateGeometry(){
195 //////////////////////////////////////////////////////////////////////
196 ////////////////////////////////////////////////////////////////////////
197 // Read geometry for the ITS
198 //
199     char topvol[5];
200     char *filtmp;
201 //
202   filtmp = gSystem->ExpandPathName(fEuclidGeometry.Data());
203   FILE *file = fopen(filtmp,"r");
204   delete [] filtmp;
205   if(file) {
206     fclose(file);
207     printf("Ready to read Euclid geometry file\n");
208     ReadEuclid(fEuclidGeometry.Data(),topvol);
209     printf("Read in euclid geometries\n");
210   } else {
211     Error("CreateGeometry"," THE GEOM FILE %s DOES NOT EXIST !",
212           fEuclidGeometry.Data());
213     exit(1);
214   } // end if(file)
215   //
216   //---Place the ITS ghost volume ITSV in its mother volume (ALIC) and make it
217   //     invisible
218   //
219   gMC->Gspos("ITSV",1,"ALIC",0,0,0,0,"ONLY");
220   //
221   //---Outputs the geometry tree in the EUCLID/CAD format
222   
223     if (fEuclidOut) {
224       gMC->WriteEuclid("ITSgeometry", "ITSV", 1, 5);
225     } // end if (fEuclidOut)
226     cout <<"finished with euclid geometrys"<< endl;
227 }
228 //______________________________________________________________________
229 void AliITSvtest::InitAliITSgeom(){
230 //     Based on the geometry tree defined in Geant 3.21, this
231 // routine initilizes the Class AliITSgeom from the Geant 3.21 ITS geometry
232 // sturture.
233     if(!(dynamic_cast<TGeant3*>(gMC))) {
234         Error("InitAliITSgeom",
235                 "Wrong Monte Carlo. InitAliITSgeom uses TGeant3 calls");
236         return;
237     } // end if
238     cout << "Reading Geometry transformation directly from Geant 3." << endl;
239     const Int_t nlayers = 6;
240     const Int_t ndeep = 9;
241     Int_t itsGeomTreeNames[nlayers][ndeep],lnam[20],lnum[20];
242     Int_t nlad[nlayers],ndet[nlayers];
243     Double_t t[3],r[10];
244     Float_t  par[20],att[20];
245     Int_t    npar,natt,idshape,imat,imed;
246     AliITSGeant3Geometry *ig = new AliITSGeant3Geometry();
247     Int_t mod,lay,lad,det,i,j,k;
248     char *names[nlayers][ndeep] = {
249      {"ALIC","ITSV","ITSD","IT12","I12B","I10B","I107","I101","ITS1"}, // lay=1
250      {"ALIC","ITSV","ITSD","IT12","I12B","I20B","I1D7","I1D1","ITS2"}, // lay=2
251      {"ALIC","ITSV","ITSD","IT34","I004","I302","ITS3","    ","    "}, // lay=3
252      {"ALIC","ITSV","ITSD","IT34","I005","I402","ITS4","    ","    "}, // lay=4
253      {"ALIC","ITSV","ITSD","IT56","I565","I562","ITS5","    ","    "}, // lay=5
254      {"ALIC","ITSV","ITSD","IT56","I569","I566","ITS6","    ","    "}};// lay=6
255     Int_t itsGeomTreeCopys[nlayers][ndeep] = {{1,1,1,1,10, 2, 4, 1, 1},// lay=1
256                                               {1,1,1,1,10, 4, 4, 1, 1},// lay=2
257                                               {1,1,1,1,14, 6, 1, 0, 0},// lay=3
258                                               {1,1,1,1,22, 8, 1, 0, 0},// lay=4
259                                               {1,1,1,1,34,22, 1, 0, 0},// lay=5
260                                               {1,1,1,1,38,25, 1, 0, 0}};//lay=6
261
262     // Sorry, but this is not very pritty code. It should be replaced
263     // at some point with a version that can search through the geometry
264     // tree its self.
265     cout << "Reading Geometry informaton from Geant3 common blocks" << endl;
266     for(i=0;i<20;i++) lnam[i] = lnum[i] = 0;
267     for(i=0;i<nlayers;i++)for(j=0;j<ndeep;j++) 
268         itsGeomTreeNames[i][j] = ig->StringToInt(names[i][j]);
269     mod = 0;
270     for(i=0;i<nlayers;i++){
271         k = 1;
272         for(j=0;j<ndeep;j++) if(itsGeomTreeCopys[i][j]!=0)
273             k *= TMath::Abs(itsGeomTreeCopys[i][j]);
274         mod += k;
275     } // end for i
276
277     if(fITSgeom!=0) delete fITSgeom;
278     nlad[0]=20;nlad[1]=40;nlad[2]=14;nlad[3]=22;nlad[4]=34;nlad[5]=38;
279     ndet[0]=4;ndet[1]=4;ndet[2]=6;ndet[3]=8;ndet[4]=22;ndet[5]=25;
280     fITSgeom = new AliITSgeom(0,6,nlad,ndet,mod);
281     mod = -1;
282     for(lay=1;lay<=nlayers;lay++){
283         for(j=0;j<ndeep;j++) lnam[j] = itsGeomTreeNames[lay-1][j];
284         for(j=0;j<ndeep;j++) lnum[j] = itsGeomTreeCopys[lay-1][j];
285         switch (lay){
286         case 1: case 2: // layers 1 and 2 are a bit special
287             lad = 0;
288             for(j=1;j<=itsGeomTreeCopys[lay-1][4];j++){
289                 lnum[4] = j;
290                 for(k=1;k<=itsGeomTreeCopys[lay-1][5];k++){
291                     lad++;
292                     lnum[5] = k;
293                     for(det=1;det<=itsGeomTreeCopys[lay-1][6];det++){
294                         lnum[6] = det;
295                         mod++;
296                         ig->GetGeometry(ndeep,lnam,lnum,t,r,idshape,npar,natt,
297                                         par,att,imat,imed);
298                         fITSgeom->CreatMatrix(mod,lay,lad,det,kSPD,t,r);
299                         if(!(fITSgeom->IsShapeDefined((Int_t)kSPD)))
300                             if(fMinorVersion==1){
301                              fITSgeom->ReSetShape(kSPD,
302                                                   new AliITSgeomSPD425Short());
303                             } else if(fMinorVersion==2)
304                              fITSgeom->ReSetShape(kSPD,
305                                                   new AliITSgeomSPD425Short());
306                     } // end for det
307                 } // end for k
308             } // end for j
309             break;
310         case 3: case 4: case 5: case 6: // layers 3-6
311             lnum[6] = 1;
312             for(lad=1;lad<=itsGeomTreeCopys[lay-1][4];lad++){
313                 lnum[4] = lad;
314                 for(det=1;det<=itsGeomTreeCopys[lay-1][5];det++){
315                     lnum[5] = det;
316                     mod++;
317                     ig->GetGeometry(7,lnam,lnum,t,r,idshape,npar,natt,
318                                     par,att,imat,imed);
319                     switch (lay){
320                     case 3: case 4:
321                         fITSgeom->CreatMatrix(mod,lay,lad,det,kSDD,t,r);
322                         if(!(fITSgeom->IsShapeDefined(kSDD))) 
323                             fITSgeom->ReSetShape(kSDD,new AliITSgeomSDD256());
324                             break;
325                         case 5:
326                             fITSgeom->CreatMatrix(mod,lay,lad,det,kSSD,t,r);
327                             if(!(fITSgeom->IsShapeDefined(kSSD))) 
328                                 fITSgeom->ReSetShape(kSSD,new AliITSgeomSSD275and75());
329                             break;
330                         case 6:
331                             fITSgeom->CreatMatrix(mod,lay,lad,det,kSSDp,t,r);
332                             if(!(fITSgeom->IsShapeDefined(kSSDp))) 
333                                 fITSgeom->ReSetShape(kSSDp,new AliITSgeomSSD75and275());
334                             break;
335                         } // end switch
336                 } // end for det
337             } // end for lad
338             break;
339         } // end switch
340     } // end for lay
341     return;
342 }
343 //_____________________________________________________________________________
344 void AliITSvtest::Init(){
345 ////////////////////////////////////////////////////////////////////////
346 //     Initialise the ITS after it has been created.
347 ////////////////////////////////////////////////////////////////////////
348     Int_t i;
349
350     cout << endl;
351     for(i=0;i<29;i++) cout << "*";cout << " ITSvtest_Init ";
352     for(i=0;i<28;i++) cout << "*";cout << endl;
353 //
354     if(fRead[0]=='\0') strncpy(fRead,fEuclidGeomDet,60);
355     if(fWrite[0]=='\0') strncpy(fWrite,fEuclidGeomDet,60);
356     if(fITSgeom!=0) delete fITSgeom;
357     fITSgeom = new AliITSgeom();
358     if(fGeomDetIn) fITSgeom->ReadNewFile(fRead);
359     if(!fGeomDetIn) this->InitAliITSgeom();
360     if(fGeomDetOut) fITSgeom->WriteNewFile(fWrite);
361     AliITS::Init();
362 //
363     for(i=0;i<72;i++) cout << "*";
364     cout << endl;
365 }
366 //_____________________________________________________________________________
367 void AliITSvtest::StepManager(){
368   //
369   // Called for every step in the ITS
370   //
371   Int_t          copy, id;
372   Int_t          copy1,copy2;
373   Float_t        hits[8];
374   Int_t          vol[4];
375   TLorentzVector position, momentum;
376   TClonesArray   &lhits = *fHits;
377   //
378   // Track status
379   vol[3] = 0;
380   if(gMC->IsTrackInside())      vol[3] +=  1;
381   if(gMC->IsTrackEntering())    vol[3] +=  2;
382   if(gMC->IsTrackExiting())     vol[3] +=  4;
383   if(gMC->IsTrackOut())         vol[3] +=  8;
384   if(gMC->IsTrackDisappeared()) vol[3] += 16;
385   if(gMC->IsTrackStop())        vol[3] += 32;
386   if(gMC->IsTrackAlive())       vol[3] += 64;
387   //
388   // Fill hit structure.
389   if(!(gMC->TrackCharge())) return;
390   //
391   // Only entering charged tracks
392   if((id = gMC->CurrentVolID(copy)) == fIdSens[0]) {
393       vol[0] = 1;
394       id = gMC->CurrentVolOffID(0,copy);
395       //detector copy in the ladder = 1<->4  (ITS1)
396       vol[1] = copy;
397       gMC->CurrentVolOffID(1,copy1);
398       //ladder copy in the module   = 1<->2  (I186)
399       gMC->CurrentVolOffID(2,copy2);
400       //module copy in the layer    = 1<->10 (I132)
401       vol[2] = copy1+(copy2-1)*2;//# of ladders in one module  = 2
402   } else if(id == fIdSens[1]){
403       vol[0] = 2;
404       id = gMC->CurrentVolOffID(0,copy);
405       //detector copy in the ladder = 1<->4  (ITS2)
406       vol[1] = copy;
407       gMC->CurrentVolOffID(1,copy1);
408       //ladder copy in the module   = 1<->4  (I131)
409       gMC->CurrentVolOffID(2,copy2);
410       //module copy in the layer    = 1<->10 (I132)
411       vol[2] = copy1+(copy2-1)*4;//# of ladders in one module  = 4
412   } else if(id == fIdSens[2]){
413       vol[0] = 3;
414       id = gMC->CurrentVolOffID(1,copy);
415       //detector copy in the ladder = 1<->5  (ITS3 is inside I314)
416       vol[1] = copy;
417       id = gMC->CurrentVolOffID(2,copy);
418       //ladder copy in the layer    = 1<->12 (I316)
419       vol[2] = copy;
420   } else if(id == fIdSens[3]){
421       vol[0] = 4;
422       id = gMC->CurrentVolOffID(1,copy);
423       //detector copy in the ladder = 1<->8  (ITS4 is inside I414)
424       vol[1] = copy;
425       id = gMC->CurrentVolOffID(2,copy);
426       //ladder copy in the layer    = 1<->22 (I417)
427       vol[2] = copy;
428   }else if(id == fIdSens[4]){
429       vol[0] = 5;
430       id = gMC->CurrentVolOffID(1,copy);
431       //detector copy in the ladder = 1<->23  (ITS5 is inside I562)
432       vol[1] = copy;
433       id = gMC->CurrentVolOffID(2,copy);
434      //ladder copy in the layer    = 1<->34 (I565)
435       vol[2] = copy;
436   }else if(id == fIdSens[5]){
437       vol[0] = 6;
438       id = gMC->CurrentVolOffID(1,copy);
439       //detector copy in the ladder = 1<->26  (ITS6 is inside I566)
440       vol[1] = copy;
441       id = gMC->CurrentVolOffID(2,copy);
442       //ladder copy in the layer = 1<->38 (I569)
443       vol[2] = copy;
444   } else {
445       return; // not an ITS volume?
446   } // end if/else if (gMC->CurentVolID(copy) == fIdSens[i])
447 //
448   gMC->TrackPosition(position);
449   gMC->TrackMomentum(momentum);
450   hits[0]=position[0];
451   hits[1]=position[1];
452   hits[2]=position[2];
453   hits[3]=momentum[0];
454   hits[4]=momentum[1];
455   hits[5]=momentum[2];
456   hits[6]=gMC->Edep();
457   hits[7]=gMC->TrackTime();
458   // Fill hit structure with this new hit.
459   new(lhits[fNhits++]) AliITShit(fIshunt,gAlice->CurrentTrack(),vol,hits);
460   return;
461 }
462