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