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