TGeant3/AliG3toRoot.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.1  2001/07/09 11:43:01  morsch
  19 Class Responsible for G3 -> Root geometry conversion.
  20
  21 */
  22
  23
  24 #include "AliG3toRoot.h"
  25 #include "AliG3Volume.h"
  26 #include "AliG3Medium.h"
  27 #include "AliG3Material.h"
  28 #include "AliNode.h"
  29 #include "AliRun.h"
  30
  31 #include <TClonesArray.h>
  32 #include <TList.h>
  33 #include <TObjArray.h>
  34 #include <TFolder.h>
  35 #include <TRotMatrix.h>
  36 #include <TFile.h>
  37 #include <TGeometry.h>
  38 #include <TBRIK.h>
  39
  40 ClassImp(AliG3toRoot)
  41     AliG3toRoot::AliG3toRoot()
  42 {
  43 //  Constructor
  44     TFolder* top = gROOT->GetRootFolder ()->AddFolder ("G3toRoot", "Geometry Tree");
  45     fTopFolder=top->AddFolder("ALIC", "Top Volume");
  46     gROOT->GetListOfBrowsables ()->Add (top, "G3toRoot Folders");
  47
  48 //
  49     fVolumes   = new TClonesArray("AliG3Volume",1000);
  50     fMaterials = new TClonesArray("AliG3Material",1000);
  51     fMedia     = new TClonesArray("AliG3Medium",1000);
  52     fRotations = 0;
  53     fExpand    = 0;
  54
  55 //  Store local copy of zebra bank entries
  56     TGeant3 *geant3 = (TGeant3*) gMC;
  57     if (geant3) {
  58         fZlq=geant3->Lq();
  59         fZq=geant3->Q();
  60         fZiq=geant3->Iq();
  61         fGclink=geant3->Gclink();
  62         fGcnum=geant3->Gcnum();
  63     }
  64 // Geometry
  65     fGeometry = new TGeometry("AliceGeom","Detailed Geometry");
  66 }
  67 void AliG3toRoot::G3toRoot()
  68 {
  69     ReadGeometryTree();
  70     ReadMaterials();
  71     ReadRotations();
  72     ConvertToRootShapes();
  73
  74 //
  75     TFolder* topFolder = gROOT->GetRootFolder ()->AddFolder ("G3", "G3");
  76     TFolder* geoFolder = topFolder->AddFolder("G3toRoot Shapes", "  ");
  77     geoFolder->Add(fGeometry);
  78     gROOT->GetListOfBrowsables ()->Add(geoFolder, "G3toRoot Shapes");
  79
  80     TFile* file = new TFile("g3toRootShapes.root", "recreate");
  81     geoFolder->Write();
  82     file->Close();
  83
  84     file = new TFile("g3toRootFolders.root", "recreate");
  85     fTopFolder->Write();
  86     file->Close();
  87
  88 }
  89
  90 void AliG3toRoot::ReadGeometryTree()
  91 {
  92 //
  93 // Copy zebra volume tree into ROOT LisTree
  94 //
  95     char *vname;
  96     char /* *namec, */ *tmp;
  97     char namec[30];
  98
  99     AliG3Volume  *volume;
 100
 101     Int_t nst=0;
 102     Int_t nlevel=1;
 103     Int_t newlevel=nlevel;
 104     Float_t mpar[3] = {2000., 2000., 2000.};
 105
 106
 107     volume = new AliG3Volume("ALIC");
 108     volume->SetIdVolume(((TGeant3*)gMC)->VolId("ALIC"));
 109     volume->SetIdCopy(0);
 110     volume->SetShape(1);
 111     volume->SetPosition(0.,0.,0.);
 112     volume->SetParameters(3,mpar);
 113     volume->SetItem(fTopFolder);
 114     (*fVolumes)[0]=volume;
 115 //
 116 //  Loop over volumes for which information has been collected
 117     while(nlevel>nst) {
 118         for (Int_t i=nst; i<nlevel; i++)
 119         {
 120             TFolder *itemi, *item2;
 121 // GEANT3 volume number
 122             Int_t ivol=TMath::Abs(Volume(i)->GetIdVolume());
 123 // Copy number
 124 // icopy=1 normal positioning
 125 // icopy>1 positioning with parposp
 126 // icopy<0 division
 127             Int_t icopy = Volume(i)->GetIdCopy();
 128 // Medium and material number, handle special case of divisions
 129             Int_t imat, imed;
 130
 131             if (icopy <0) {
 132                 imed=Medium(-ivol);
 133                 imat=Material(-ivol);
 134             } else {
 135                 imed=Medium(ivol);
 136                 imat=Material(ivol);
 137             }
 138
 139             TArrayF pos;
 140             pos = Volume(i)->Position(0);
 141
 142 //
 143 // Number of children
 144             Int_t nch = NChildren(ivol);
 145             strcpy(namec,((TGeant3*)gMC)->VolName(ivol));
 146             if (nch >= 0) {
 147                 printf("\n %s has %d children  \n ", namec,  nch);
 148             } else {
 149                 printf("\n %s has  divisions \n ", namec);
 150             }
 151 //
 152 // Name to be used in ListTree
 153             vname = new char[5];
 154             strncpy(vname,namec, 4);
 155             vname[4]='\0';
 156
 157             if (icopy >1) {
 158                 sprintf(namec,"%s*%3dPos",namec,icopy);
 159             } else if (icopy <0) {
 160                 sprintf(namec,"%s*%3dDiv",namec,-icopy);
 161             }
 162             if (i>0) {
 163                 itemi=(TFolder*) Volume(i)->GetItem();
 164             } else {
 165                 itemi=fTopFolder;
 166             }
 167             Volume(i)->SetName(namec);
 168             Volume(i)->SetIdMedium(imed);
 169             Volume(i)->SetIdMaterial(imat);
 170
 171 //
 172 // Add volume to list tree
 173
 174             if (i>0) {
 175                 item2 = AddVolume(new AliG3Volume(*Volume(i)), itemi, namec);
 176             } else {
 177                 item2 = fTopFolder;
 178                 fTopFolder->Add(new AliG3Volume(*Volume(i)));
 179             }
 180 //
 181 // Collect children information
 182 //
 183 // nch < 0: Children by division
 184             if (nch < 0) {
 185 //
 186 // Geant volume number
 187                 Int_t icvol=Child(ivol,1);
 188 // Name
 189                 strcpy(namec,((TGeant3*)gMC)->VolName(-icvol));
 190                 tmp = new char[4];
 191                 strncpy(tmp,(char *) &namec, 4);
 192                 volume = new AliG3Volume(namec);
 193                 volume->SetIdVolume(-icvol);
 194
 195                 Int_t jvo  = fZlq[fGclink->jvolum-ivol];
 196                 Int_t jdiv = fZlq[jvo-1];
 197 // Axis
 198                 Int_t iaxe = Int_t (fZq[jdiv+1]);
 199 // System volume number
 200 //              Int_t ivin = Int_t (fZq[jdiv+2]);
 201 // Number of divisions
 202                 Int_t ndiv = Int_t (fZq[jdiv+3]);
 203 // Start value
 204                 Float_t startc =    fZq[jdiv+4];
 205 // Step
 206                 Float_t step   =    fZq[jdiv+5];
 207
 208                 Int_t ish   = Volume(i)->Shape();
 209                 Int_t npar  = Volume(i)->NParam();
 210 // Copy here shape parameters from mother
 211 // Expand divisions later (when needed)
 212                 TArrayF apar;
 213                 apar.Set(npar);
 214                 Volume(i)->Parameters(0, apar);
 215                 Float_t* par = new Float_t[npar];
 216                 for (Int_t jj=0; jj<npar; jj++) par[jj]=apar.At(jj);
 217                 volume->SetIdCopy(-ndiv);
 218                 volume->SetDivision(ndiv, iaxe, startc, step);
 219                 volume->SetParameters(npar, par);
 220                 volume->SetPosition(0.,0.,0.);
 221                 volume->SetShape(ish);
 222                 volume->SetRotMatrix(0);
 223 // Link to mother
 224                 volume->SetItem(item2);
 225                 (*fVolumes)[newlevel]=volume;
 226                 newlevel++;
 227 //
 228 // Children by positioning
 229             } else if (nch > 0) {
 230                 Int_t nnew=0;
 231 // Loop over children
 232                 for (Int_t j=0; j<nch; j++)
 233                 {
 234                     Int_t jvo  = fZlq[fGclink->jvolum-ivol];
 235                     Int_t jin  = fZlq[jvo-j-1];
 236                     Int_t ivin = Int_t (fZq[jin + 2 ]);
 237                     Int_t jvin = fZlq[fGclink->jvolum - ivin];
 238                     Int_t ish  = Int_t (fZq[jvin + 2]);
 239                     Int_t irot = Int_t(fZq[jin +4]);
 240                     Float_t x    =  fZq[jin +5];
 241                     Float_t y    =  fZq[jin +6];
 242                     Float_t z    =  fZq[jin +7];
 243                     Int_t ndata = Int_t(fZiq[jin-1]);
 244                     Int_t npar;
 245                     Float_t* par;
 246 //
 247 // Get shape parameters
 248 //
 249 // parp or normal positioning ?
 250                     if (ndata == 8) {
 251                         npar  = Int_t (fZq[jvin+5]);
 252                         par = new Float_t[npar];
 253                         if (npar>0) {
 254                             for(Int_t jj=0; jj<npar; jj++)
 255                                 par[jj] = fZq[jvin+7+jj];
 256                         }
 257                     } else {
 258                         npar =  Int_t (fZq[jin+9]);
 259                         par = new Float_t[npar];
 260                         for(Int_t jj=0;jj<npar;jj++)
 261                             par[jj]=fZq[jin+10+jj];
 262                     }
 263 //
 264 // Find out if this volume was already positioned and count copies
 265                     Int_t icvol=Child(ivol,j+1);
 266                     icvol = TMath::Abs(icvol);
 267                     Bool_t inList=kFALSE;
 268                     AliG3Volume* copy0=0;
 269
 270                     for (Int_t k=0; k<nnew; k++) {
 271                         if (icvol==
 272                             Volume(newlevel-k-1)->GetIdVolume())
 273                         {
 274                             copy0  = Volume(newlevel-k-1);
 275                             Volume(newlevel-k-1)->AddCopy();
 276                             inList=kTRUE;
 277                         }
 278                     }
 279 //
 280 // New child
 281 //
 282 // Name
 283                     strcpy(namec,((TGeant3*)gMC)->VolName(icvol));
 284                     tmp = new char[4];
 285                     strncpy(tmp,(char *) &namec, 4);
 286                     volume = new AliG3Volume(namec);
 287                     volume->SetPosition(x,y,z);
 288                     volume->SetShape(ish);
 289                     volume->SetParameters(npar, par);
 290                     volume->SetRotMatrix(irot);
 291                     if (ndata != 8) volume->SetPosp(kTRUE);
 292                     volume->SetIdVolume(icvol);
 293                     volume->SetIdCopy(1);
 294 // Link to mother
 295                     volume->SetItem(item2);
 296                     if (!inList) {
 297                         (*fVolumes)[newlevel]=volume;
 298                         newlevel++;
 299                         nnew++;
 300                     } else {
 301                         copy0->AddCopy(volume);
 302                     }
 303                 }
 304             }
 305         }
 306 // Move one level deaper
 307         nst=nlevel;
 308         nlevel=newlevel;
 309     }
 310 }
 311
 312 void AliG3toRoot::ReadMaterials()
 313 {
 314 //
 315 // Puts media and material names into ComboBox and
 316 // collects material information
 317 //
 318     Float_t a, z, dens, radl, absl;
 319     a=z=dens=radl=absl=-1;
 320     Int_t npar=0;
 321     Int_t imat, isvol, ifield;
 322     Float_t fieldm, tmaxfd, stemax, deemax, epsil, stmin;
 323     Float_t par[50];
 324     Float_t ubuf[50];
 325     Int_t nwbuf;
 326
 327     char natmed[21], namate[21];
 328 //
 329 // Loop over media
 330     Int_t nEntries=0;
 331
 332     for(Int_t itm=1;itm<=fGcnum->ntmed;itm++) {
 333         Int_t jtm  = fZlq[fGclink->jtmed-itm];
 334         if (jtm > 0) {
 335             nEntries++;
 336 //
 337 // Get medium parameters
 338             ((TGeant3*)(gMC))->
 339                 Gftmed(itm, natmed, imat, isvol, ifield, fieldm,
 340                        tmaxfd, stemax, deemax, epsil, stmin, ubuf, &nwbuf);
 341             strncpy(natmed,(char*)&fZiq[jtm+1],20);
 342             natmed[20]='\0';
 343 //
 344 // Create new medium object
 345             AliG3Medium * medium =
 346                 new AliG3Medium(itm, imat, natmed, isvol, ifield, fieldm,
 347                                       tmaxfd, stemax, deemax, epsil, stmin);
 348             (*fMedia)[nEntries-1]=medium;
 349         { //Begin local scope for j
 350           for (Int_t j=1; j<=22; j++) {
 351             medium->SetPar(j,Cut(itm,j));
 352           }
 353         } //End local scope for j
 354         { //Begin local scope for j
 355           for (Int_t j=23; j<=26; j++) {
 356             medium->SetPar(j,Cut(itm,j+3));
 357           }
 358         } //End local scope for j
 359         { //Begin local scope for j
 360           for (Int_t j=27; j<=29; j++) {
 361             medium->SetPar(j,Cut(itm,j+4));
 362           }
 363         } //End local scope for j
 364 //
 365 // Associated material
 366             imat =  Int_t (fZq[jtm+6]);
 367             Int_t jma  =  Int_t (fZlq[fGclink->jmate-imat]);
 368 //
 369 // Get material parameters
 370             ((TGeant3*)(gMC))->Gfmate (imat,namate,a,z,dens,radl,absl,par,npar);
 371             strncpy(namate,(char *)&fZiq[jma+1],20);
 372             namate[20]='\0';
 373 //
 374 // Create new material object
 375             AliG3Material* material = new AliG3Material(namate, " ", a,z,dens,radl,absl);
 376             material->SetId(imat);
 377
 378             (*fMaterials)[nEntries-1] = material;
 379             material->Dump();
 380         }
 381     }
 382 }
 383
 384 void AliG3toRoot::ReadRotations()
 385 {
 386 // Read rotation matrices
 387     Int_t nrot = 0;
 388     Int_t irm, jrm;
 389     char name[9];
 390     if(fGclink->jrotm) nrot = Int_t(fZiq[fGclink->jrotm-2]);
 391     if (nrot) {
 392          fRotations = new TObjArray(nrot);
 393          for(irm=1;irm<=nrot;irm++) {
 394              jrm  = fZlq[fGclink->jrotm-irm];
 395              sprintf(name, "R%d", irm);
 396          //PH        (*fRotations)[irm-1] =
 397            //PH          new TRotMatrix(name, "Rotation", fZq[jrm+11], fZq[jrm+12], fZq[jrm+13],
 398            //PH                         fZq[jrm+14], fZq[jrm+15], fZq[jrm+16]);
 399              fRotations->AddAt(
 400                  new TRotMatrix(name, "Rotation", fZq[jrm+11], fZq[jrm+12], fZq[jrm+13],
 401                                 fZq[jrm+14], fZq[jrm+15], fZq[jrm+16]), irm-1);
 402          }
 403     }
 404 }
 405
 406 Int_t AliG3toRoot::NChildren(Int_t idvol)
 407 {
 408 //
 409 // Return number of children for volume idvol
 410     Int_t jvo = fZlq[fGclink->jvolum-idvol];
 411     Int_t nin = Int_t(fZq[jvo+3]);
 412     return nin;
 413 }
 414
 415 Int_t AliG3toRoot::Child(Int_t idvol, Int_t idc)
 416 {
 417 //
 418 // Return GEANT id of child number idc of volume idvol
 419     Int_t jvo = fZlq[fGclink->jvolum-idvol];
 420     Int_t nin=idc;
 421     Int_t jin = fZlq[jvo-nin];
 422     Int_t numb =  Int_t (fZq[jin +3]);
 423     if (numb > 1) {
 424         return -Int_t(fZq[jin+2]);
 425     } else {
 426         return Int_t(fZq[jin+2]);
 427     }
 428 }
 429
 430 Int_t AliG3toRoot::Medium(Int_t idvol)
 431 {
 432 //
 433 // Return medium number for volume idvol.
 434 // If idvol is negative the volume results from a division.
 435     Int_t imed;
 436     if (idvol > 0) {
 437         Int_t jvo = fZlq[fGclink->jvolum-idvol];
 438         imed = Int_t(fZq[jvo+4]);
 439     } else {
 440         idvol=-idvol;
 441         Int_t jdiv = fZlq[fGclink->jvolum-idvol];
 442         Int_t ivin = Int_t ( fZq[jdiv+2]);
 443         Int_t jvin = fZlq[fGclink->jvolum-ivin];
 444         imed = Int_t (fZq[jvin+4]);
 445     }
 446     return imed;
 447 }
 448
 449 Int_t AliG3toRoot::Material(Int_t idvol)
 450 {
 451 // Return material number for volume idvol.
 452 // If idvol is negative the volume results from a division.
 453
 454     Int_t imed=Medium(idvol);
 455     Int_t jtm  = fZlq[fGclink->jtmed-imed];
 456     return Int_t (fZq[jtm+6]);
 457 }
 458
 459
 460 Float_t AliG3toRoot::Cut(Int_t imed, Int_t icut)
 461 {
 462 // Return cuts icut for medium idmed
 463 //
 464     Int_t jtm  = fZlq[fGclink->jtmed-imed];
 465 //
 466 //  Have the defaults been modified ??
 467     Int_t jtmn = fZlq[jtm];
 468     if (jtmn >0) {
 469         jtm=jtmn;
 470
 471     } else {
 472         jtm=fGclink->jtmed;
 473     }
 474
 475     return Float_t (fZq[jtm+icut]);
 476 }
 477
 478 TFolder* AliG3toRoot::AddVolume(TObject * obj, TFolder *parent, const char* name)
 479 {
 480 // Add item to the list tree
 481     TFolder* newFolder = parent->AddFolder(name, "volume");
 482     newFolder->Add(obj);
 483     return newFolder;
 484 }
 485
 486 AliG3Volume* AliG3toRoot::Volume(Int_t id)
 487
 488 {
 489 // Volume at id
 490     return (AliG3Volume *) (fVolumes->UncheckedAt(id));
 491 }
 492
 493
 494 void AliG3toRoot::ConvertToRootShapes(TFolder *item, AliNode** node, Int_t nNodes)
 495 {
 496 // Convert the geometry represented by TFolders into root shapes
 497 //
 498     AliG3Volume* volume;
 499     TArrayF pos;
 500     Int_t irot, imat;
 501     Int_t npos = 0;
 502     Bool_t top=kFALSE;
 503
 504     char nameV[20];
 505     char nameN[20];
 506
 507     if (!item) {
 508         item = fTopFolder;
 509         top  = kTRUE;
 510     }
 511     if (nNodes == 0) nNodes = 1;
 512
 513     TList* folders = (TList*) item->GetListOfFolders();
 514     TIter next(folders);
 515
 516     npos = 0;
 517     volume = ((AliG3Volume *) item->FindObject(item->GetName()));
 518     Int_t ncopy = volume->NCopies();
 519     AliNode** newnode = new AliNode*[(ncopy+1)*nNodes];
 520
 521     for (Int_t ino=0; ino<nNodes; ino++) {
 522
 523         pos    = volume->Position(0);
 524         irot   = volume->RotMatrix();
 525         imat   = volume->Material();
 526
 527         sprintf(nameV,"v%s%d", volume->GetName(), npos);
 528         sprintf(nameN,"%s%d" , volume->GetName(), npos);
 529
 530
 531         imat   = volume->Material();
 532         Int_t nmat = fMaterials->GetEntriesFast();
 533         AliG3Material* mat=0;
 534
 535         for (Int_t mati=0; mati<nmat; mati++) {
 536             mat = (AliG3Material*)
 537                 (fMaterials->UncheckedAt(mati));
 538             if ((mat->Id())==imat) break;
 539         }
 540
 541
 542         volume->CreateTShape(nameV, mat);
 543
 544         if (!top) {
 545             node[ino]->cd();
 546             newnode[npos] = new AliNode(nameN," ",nameV, pos[0], pos[1], pos[2]);
 547             newnode[npos]->SetDivision(volume->Ndiv(), volume->Axis(),
 548                                        volume->Step(), volume->StartC());
 549             if (irot > 0) newnode[npos]->SetMatrix((TRotMatrix*) (*fRotations)[irot-1]);
 550             npos++;
 551
 552         } else {
 553             new TMaterial("void","Vacuum",0,0,0);  //Everything is void
 554             TBRIK *brik   = new TBRIK("D_alice","alice volume",
 555                                       "void",2000,2000,3000);
 556             brik->SetVisibility(0);
 557             newnode[npos] = new AliNode("alice","alice","D_alice");
 558             npos++;
 559         }
 560
 561         if (ncopy) {
 562             for (Int_t icop = 0; icop < ncopy; icop++) {
 563                 AliG3Volume* copy = volume->Copy(icop);
 564
 565                 sprintf(nameN,"%s%d" , volume->GetName(), icop+1);
 566                 if (copy->Posp()) {
 567                     sprintf(nameV,"v%s%d", volume->GetName(), icop+1);
 568                     volume->CreateTShape(nameV, mat);
 569                 }
 570                 node[ino]->cd();
 571                 pos    = copy->Position(0);
 572                 irot   = copy->RotMatrix();
 573                 newnode[npos] = new AliNode(nameN," ",nameV, pos[0], pos[1], pos[2]);
 574                 newnode[npos]->SetDivision(volume->Ndiv(), volume->Axis(),
 575                                            volume->Step(), volume->StartC());
 576                 if (irot >0) newnode[npos]->SetMatrix((TRotMatrix*) (*fRotations)[irot-1]);
 577                 npos++;
 578             } // copy loop
 579         } // copies exist
 580     } // node loop
 581     TObject* obj;
 582
 583     while ((obj = next()))
 584     {
 585         if ((AliG3Volume*) obj == volume) continue;
 586         item = (TFolder*) obj;
 587
 588         ConvertToRootShapes(item, newnode, npos);
 589     }
 590
 591 //  Expand divisions of demand
 592     if (fExpand) ExpandDivisions();
 593 }
 594
 595 void AliG3toRoot::ExpandDivisions(AliNode* node)
 596 {
 597 // Expand volume divisions
 598 //
 599
 600     Bool_t top = kFALSE;
 601
 602     if (!node) {
 603         node = (AliNode*) fGeometry->GetNode("alice");
 604         top = kTRUE;
 605     }
 606     if (!top && node->Ndiv() > 0 && node->Axis()>0) {
 607
 608         node->ExpandDivisions();
 609         node->GetParent()->RecursiveRemove(node);
 610         ExpandDivisions();
 611
 612     } else {
 613         TList* sons = node->GetListOfNodes();
 614         TIter next(sons);
 615         AliNode* son=0;
 616
 617         while((son = (AliNode*)next())) {
 618             ExpandDivisions(son);
 619         }
 620     }
 621 }
 622
 623
 624
 625
 626
 627
 628