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