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