C++ TGeoManager::SetTopVolume方法代码示例

void OpNoviceGeom()
{
// Create root geometry corresponding to OpNovice example
  double fExpHall_x, fExpHall_y, fExpHall_z;
  double fTank_x, fTank_y, fTank_z;
  double fBubble_x, fBubble_y, fBubble_z;
  double a, z, density;
  int nelements;

  fExpHall_x = fExpHall_y = fExpHall_z = 10000.0;
  fTank_x    = fTank_y    = fTank_z    =  5000.0;
  fBubble_x  = fBubble_y  = fBubble_z  =  500.0;
  TGeoManager *geom = new TGeoManager("OpNoviceGeom", "GEANT4 OpNovice example");
// Air
//
  TGeoElement *N = new TGeoElement("Nitrogen", "N", z=7 , a=14.01);
  TGeoElement *O = new TGeoElement("Oxygen"  , "O", z=8 , a=16.00);
  TGeoMixture *air = new TGeoMixture("Air", 2, density=1.29E-3);
  air->AddElement(N, 0.7);
  air->AddElement(O, 0.3);
  TGeoMedium *medair = new TGeoMedium("Air", 1, air);
// Water
//
  TGeoElement *H = new TGeoElement("Hydrogen", "H", z=1 , a=1.01);
  TGeoMixture *water = new TGeoMixture("Water", 2, density=1.0); 
  water->AddElement(H,2);
  water->AddElement(O,1);
  TGeoMedium *medwater = new TGeoMedium("Water", 2, water);
// The experimental Hall
//
  TGeoBBox *expHall_box = new TGeoBBox("World",fExpHall_x,fExpHall_y,fExpHall_z);
  TGeoVolume *expHall_log = new TGeoVolume("World", expHall_box, medair);
  expHall_log->SetLineColor(1);
  expHall_log->SetVisContainers();
  geom->SetTopVolume(expHall_log);
// The Water Tank
//
  TGeoBBox *waterTank_box = new TGeoBBox("Tank",fTank_x,fTank_y,fTank_z);
  TGeoVolume *waterTank_log = new TGeoVolume("Tank", waterTank_box, medwater);
  waterTank_log->SetLineColor(kBlue);
  waterTank_log->SetTransparency(70);
  waterTank_log->SetVisContainers();
  expHall_log->AddNode(waterTank_log, 0);
// The Air Bubble
//
  TGeoBBox *bubbleAir_box = new TGeoBBox("Bubble",fBubble_x,fBubble_y,fBubble_z);
  TGeoVolume *bubbleAir_log = new TGeoVolume("Bubble", bubbleAir_box, medair);
  bubbleAir_log->SetLineColor(kCyan);
  bubbleAir_log->SetTransparency(70);
  waterTank_log->AddNode(bubbleAir_log, 0, new TGeoTranslation(0,2500,0));
  
  geom->CloseGeometry();
  geom->SetTopVisible(true);
  geom->Export("OpNoviceGeom.root");
}

void rootgeom()
{
//Definition of a simple geometry (the 4 ROOT characters)
//Author: Andrei Gheata
   
   // gStyle->SetCanvasPreferGL(true);
   gSystem->Load("libGeom");
   TGeoManager *geom = new TGeoManager("simple1", "Simple geometry");
   
   //--- define some materials
   TGeoMaterial *matVacuum = new TGeoMaterial("Vacuum", 0,0,0);
   TGeoMaterial *matAl = new TGeoMaterial("Al", 26.98,13,2.7);
//   //--- define some media
   TGeoMedium *Vacuum = new TGeoMedium("Vacuum",1, matVacuum);
   TGeoMedium *Al = new TGeoMedium("Root Material",2, matAl);

   //--- define the transformations
   TGeoTranslation *tr1 = new TGeoTranslation(20., 0, 0.);
   TGeoTranslation *tr2 = new TGeoTranslation(10., 0., 0.);
   TGeoTranslation *tr3 = new TGeoTranslation(10., 20., 0.);
   TGeoTranslation *tr4 = new TGeoTranslation(5., 10., 0.);
   TGeoTranslation *tr5 = new TGeoTranslation(20., 0., 0.);
   TGeoTranslation *tr6 = new TGeoTranslation(-5., 0., 0.);
   TGeoTranslation *tr7 = new TGeoTranslation(7.5, 7.5, 0.);
   TGeoRotation   *rot1 = new TGeoRotation("rot1", 90., 0., 90., 270., 0., 0.);
   TGeoCombiTrans *combi1 = new TGeoCombiTrans(7.5, -7.5, 0., rot1);
   TGeoTranslation *tr8 = new TGeoTranslation(7.5, -5., 0.);
   TGeoTranslation *tr9 = new TGeoTranslation(7.5, 20., 0.);
   TGeoTranslation *tr10 = new TGeoTranslation(85., 0., 0.);
   TGeoTranslation *tr11 = new TGeoTranslation(35., 0., 0.);
   TGeoTranslation *tr12 = new TGeoTranslation(-15., 0., 0.);
   TGeoTranslation *tr13 = new TGeoTranslation(-65., 0., 0.);
   
   TGeoTranslation  *tr14 = new TGeoTranslation(0,0,-100);
   TGeoCombiTrans *combi2 = new TGeoCombiTrans(0,0,100, 
                                   new TGeoRotation("rot2",90,180,90,90,180,0));
   TGeoCombiTrans *combi3 = new TGeoCombiTrans(100,0,0, 
                                   new TGeoRotation("rot3",90,270,0,0,90,180));
   TGeoCombiTrans *combi4 = new TGeoCombiTrans(-100,0,0, 
                                   new TGeoRotation("rot4",90,90,0,0,90,0));
   TGeoCombiTrans *combi5 = new TGeoCombiTrans(0,100,0, 
                                   new TGeoRotation("rot5",0,0,90,180,90,270));
   TGeoCombiTrans *combi6 = new TGeoCombiTrans(0,-100,0, 
                                   new TGeoRotation("rot6",180,0,90,180,90,90));
   
   //--- make the top container volume
   Double_t worldx = 110.;
   Double_t worldy = 50.;
   Double_t worldz = 5.;
   TGeoVolume *top = geom->MakeBox("TOP", Vacuum, 270., 270., 120.);
   geom->SetTopVolume(top);
   TGeoVolume *replica = geom->MakeBox("REPLICA", Vacuum,120,120,120);
   replica->SetVisibility(kFALSE);
   TGeoVolume *rootbox = geom->MakeBox("ROOT", Vacuum, 110., 50., 5.);
   rootbox->SetVisibility(kFALSE);
   
   //--- make letter 'R'
   TGeoVolume *R = geom->MakeBox("R", Vacuum, 25., 25., 5.);
   R->SetVisibility(kFALSE);
   TGeoVolume *bar1 = geom->MakeBox("bar1", Al, 5., 25, 5.);
   bar1->SetLineColor(kRed);
   R->AddNode(bar1, 1, tr1);
   TGeoVolume *bar2 = geom->MakeBox("bar2", Al, 5., 5., 5.);
   bar2->SetLineColor(kRed);
   R->AddNode(bar2, 1, tr2);
   R->AddNode(bar2, 2, tr3);
   TGeoVolume *tub1 = geom->MakeTubs("tub1", Al, 5., 15., 5., 90., 270.);
   tub1->SetLineColor(kRed);
   R->AddNode(tub1, 1, tr4);
   TGeoVolume *bar3 = geom->MakeArb8("bar3", Al, 5.);
   bar3->SetLineColor(kRed);
   TGeoArb8 *arb = (TGeoArb8*)bar3->GetShape();
   arb->SetVertex(0, 15., -5.);
   arb->SetVertex(1, 0., -25.);
   arb->SetVertex(2, -10., -25.);
   arb->SetVertex(3, 5., -5.);
   arb->SetVertex(4, 15., -5.);
   arb->SetVertex(5, 0., -25.);
   arb->SetVertex(6, -10., -25.);
   arb->SetVertex(7, 5., -5.);
   R->AddNode(bar3, 1, gGeoIdentity);
   
   //--- make letter 'O'
   TGeoVolume *O = geom->MakeBox("O", Vacuum, 25., 25., 5.);
   O->SetVisibility(kFALSE);
   TGeoVolume *bar4 = geom->MakeBox("bar4", Al, 5., 7.5, 5.);
   bar4->SetLineColor(kYellow);
   O->AddNode(bar4, 1, tr5);
   O->AddNode(bar4, 2, tr6);
   TGeoVolume *tub2 = geom->MakeTubs("tub1", Al, 7.5, 17.5, 5., 0., 180.);
   tub2->SetLineColor(kYellow);
   O->AddNode(tub2, 1, tr7);
   O->AddNode(tub2, 2, combi1);
   
   //--- make letter 'T'
   TGeoVolume *T = geom->MakeBox("T", Vacuum, 25., 25., 5.);
   T->SetVisibility(kFALSE);
   TGeoVolume *bar5 = geom->MakeBox("bar5", Al, 5., 20., 5.);
   bar5->SetLineColor(kBlue);
   T->AddNode(bar5, 1, tr8);
//.........这里部分代码省略.........

void robot()
{
  // Drawing a famous Korean robot, TaekwonV, using ROOT geometry class.
  // Name: robot.C
  // Author: Jin Hui Hwang, Dept. of Physics, Univ. of Seoul
  // Reviewed by Sunman Kim ([email protected])
  // Supervisor: Prof. Inkyu Park ([email protected])
  //
  // How to run: .x robot.C in ROOT terminal, then use OpenGL
  //
  // This macro was created for the evaluation of Computational Physics course in 2006.
  // We thank to Prof. Inkyu Park for his special lecture on ROOT and to all of ROOT team
  //

  TGeoManager *Robot = new TGeoManager("Robot", "This is Taegwon V");

  TGeoMaterial *vacuum = new TGeoMaterial("vacuum", 0, 0, 0);
  TGeoMaterial *Fe     = new TGeoMaterial("Fe", 55.845, 26, 7.87);

  TGeoMedium *Air  = new TGeoMedium("Vacuum", 0, vacuum);
  TGeoMedium *Iron = new TGeoMedium("Iron", 1, Fe);

  // create volume

  TGeoVolume *top = Robot->MakeBox("top", Air, 1000, 1000, 1000);
  Robot->SetTopVolume(top);
  Robot->SetTopVisible(0);
  // If you want to see the boundary, please input the number, 1 instead of 0.
  // Like this, geom->SetTopVisible(1);

  // head
  TGeoVolume *Band = Robot->MakeEltu("Band", Iron, 20, 20, 2.5);
  Band->SetLineColor(12);
  Band->SetFillColor(12);
  TGeoVolume *Band_b = Robot->MakeSphere("Band_b", Iron, 0, 2, 0, 180, 180, 360);
  Band_b->SetLineColor(2);
  Band_b->SetFillColor(2);
  TGeoVolume *Head = Robot->MakeSphere("Head", Iron, 0, 19, 0, 180, 180, 360);
  Head->SetLineColor(17);
  Head->SetFillColor(17);
  TGeoVolume *Horn = Robot->MakeSphere("Horn", Iron, 0, 10, 60, 180, 240, 300);

  // drawing head
  top->AddNodeOverlap(Band, 1, new TGeoTranslation(0, 0, 90));
  float Phi = 3.14;
  int N     = 10;

  for (int i = 0; i <= N; i++) {
    top->AddNodeOverlap(Band_b, 1, new TGeoCombiTrans(sin(2 * Phi / N * i) * 19, -cos(2 * Phi / N * i) * 19, 90,
                                                      new TGeoRotation("R1", -90 + (360 / N * i), -90, 90)));
  }
  top->AddNodeOverlap(Head, 1, new TGeoCombiTrans(0, 0, 87.5, new TGeoRotation("R2", 0, -90, 0)));

  char name[50];
  float pcs = 30;
  for (int i = 1; i < pcs; i++) {
    sprintf(name, "Horn%d", i);
    Horn = Robot->MakeSphere(name, Iron, 10 - 10 / pcs * i, 10, 180 - (120 / pcs) * i, 180 - ((120 / pcs) * (i - 1)),
                             240, 300);
    Horn->SetLineColor(2);
    Horn->SetFillColor(2);
    top->AddNodeOverlap(Horn, 1, new TGeoCombiTrans(0, 8, 102, new TGeoRotation("R2", 0, 140, 0)));
    top->AddNodeOverlap(Horn, 1, new TGeoCombiTrans(0, -8, 102, new TGeoRotation("R2", 180, 140, 0)));
  }

  // face
  TGeoVolume *Migan = Robot->MakeGtra("Migan", Iron, 3, 0, 0, 0, 3, 2, 11, 0, 3, 3, 11, 0);
  Migan->SetLineColor(17);
  Migan->SetFillColor(17);
  TGeoVolume *Ko = Robot->MakeGtra("Ko", Iron, 7, 0, 0, 0, 3, 1, 5, 0, 3, 2, 5, 0);
  Ko->SetLineColor(17);
  Ko->SetFillColor(17);
  TGeoVolume *Ko_m = Robot->MakeBox("Ko_m", Iron, 2, 8, 4);
  Ko_m->SetLineColor(17);
  Ko_m->SetFillColor(17);
  TGeoVolume *Bol_1 = Robot->MakeBox("Bol_1", Iron, 7, 5.5, 7);
  Bol_1->SetLineColor(17);
  Bol_1->SetFillColor(17);
  TGeoVolume *Bol_2 = Robot->MakeGtra("Bol_2", Iron, 1, 0, 0, 0, 7, 0, 9, 0, 7, 0, 9, 0);
  Bol_2->SetLineColor(17);
  Bol_2->SetFillColor(17);
  TGeoVolume *Noon = Robot->MakeBox("Noon", Iron, 1, 10, 5);
  Noon->SetLineColor(12);
  Noon->SetFillColor(12);
  TGeoVolume *Tuck = Robot->MakeBox("Tuck", Iron, 2, 10, 5.5);
  Tuck->SetLineColor(2);
  Tuck->SetFillColor(2);
  TGeoVolume *Tuck_1 = Robot->MakeBox("Tuck_1", Iron, 2, 9, 1);
  Tuck_1->SetLineColor(2);
  Tuck_1->SetFillColor(2);
  TGeoVolume *Tuck_2 = Robot->MakeBox("Tuck_2", Iron, 3, 1, 14);
  Tuck_2->SetLineColor(2);
  Tuck_2->SetFillColor(2);
  TGeoVolume *Tuck_j = Robot->MakeSphere("Tuck_j", Iron, 0, 3.5, 0, 180, 0, 360);
  Tuck_j->SetLineColor(5);
  Tuck_j->SetFillColor(5);
  TGeoVolume *Ear = Robot->MakeCons("Ear", Iron, 1, 0, 3, 0, 3, 0, 360);
  Ear->SetLineColor(5);
  Ear->SetFillColor(5);
  TGeoVolume *Ear_2 = Robot->MakeCone("Ear_2", Iron, 5, 0, 0, 0, 3);
//.........这里部分代码省略.........

void
show
(
 char* inpDir  = "",         // MuDST directory
 char* inpFile = "show.lis", // MuDST file(s);                      
 char* outFile = "show.root",// output tree file
 Int_t nFiles  = 50,         // # of MuDST file(s)
 Int_t nEvents = 100         // # of events
 )
  // remeber to adjust dbase timestamp below !!!! 
  // what a ... design
{ 
  //gErrorIgnoreLevel=1999;

  // load root/root4star libraries
  gROOT->LoadMacro("$STAR/StRoot/StMuDSTMaker/COMMON/macros/loadSharedLibraries.C");
  loadSharedLibraries();

  // load more libraries :)
  gSystem->Load("libmysqlclient");
  gSystem->Load("libGeom");
  gSystem->Load("StDbLib");
  gSystem->Load("StDbBroker");
  gSystem->Load("St_db_Maker");

  // load even more libraries (EEMC stuff) 
  gSystem->Load("StEEmcUtil");
  gSystem->Load("StEEmcDbMaker");
  gSystem->Load("StEEmcPoolTTM");

  // create the chain    
  chain = new StChain("StChain"); 
  //

  now = new TDatime;
  // for display
  TCanvas    *c1   = new TCanvas("eemc","eemc",10,10,1000,1000);
  TPaveLabel *tlab = new TPaveLabel(-0.99,+0.99,+0.99,+0.90,
				    "EEMC TOWERS & TPC TRACKS     Piotr A Zolnierczuk (IU)");

  eventInfo = new TPaveText (-0.99,-0.99,+0.0 ,-0.75);
  dateInfo  = new TPaveLabel(+0.60,-0.99,+0.99,-0.95,now->AsString());

  TGeoManager  *gm    = new TGeoManager("eemc", "eemc tower display");
  TGeoVolume   *top   = gm->MakeBox("star",0, 200., 200., 350.);
  TGeoVolume   *smbox = gm->MakeBox("smbox1",0, 2., 2., 2.);
  smbox->SetLineColor(kRed);
  // eemc 
  eemc  = new EEmcTTDisplay();
  eemc->SetMagneticField(0.5); // in Tesla
  eemc->SetShowExtrapolatedTracks(true);
 
  TGeoTranslation *etra = new TGeoTranslation(0.0,0.0,0.5*(eemc->getZ1()+eemc->getZ2()));
  top->AddNode(smbox, 1,NULL);
  top->AddNode(eemc(),1,etra);
  gm->SetTopVolume(top);
  gm->CloseGeometry();
  gm->SetVisLevel(4);
  gm->SetVisOption(0);

  c1->SetTheta(90);
  c1->SetPhi(0);

  top->Draw();
  tlab->Draw();

  gPad->Update();
  
  // now we add Makers to the chain...  some of that is black magic :) 
  muDstMk  = new StMuDstMaker(0,0,inpDir,inpFile,"",nFiles);       // muDST main chain
  StMuDbReader  *db       = StMuDbReader::instance();              // need the database
  St_db_Maker   *dbMk        = new St_db_Maker("StarDb", "MySQL:StarDb");   // need another db(?) 
  new StEEmcDbMaker("eemcDb");                     // need EEMC database  

  // now comment in/out/change the below if you want it your way
  dbMk->setTimeStampDay(20040331);   // format: yyyymmdd

  // finally after so many lines we arrive at the good stuff

  ttm = new  EEmcTTMMaker ("TTM",muDstMk,eemcDbMk);
  ttm->Summary(cout);    // 

  StMuDebug::setLevel(0);

  chain->Init();

  StEEmcDb *eemcDb = (StEEmcDb*)chain->GetDataSet("StEEmcDb");
  eemcDb->setSectors(1,12);            // request EEMC DB for sectors you need (dafault:1-12)
  eemcDb->setPreferedFlavor("onlped","eemcPMTped"); // request alternative flavor

  chain->ls(3);

  //---------------------------------------------------
  next();
}

void snoopy()
{
   TGeoManager *geom = new TGeoManager("snoopy", "Snoopy Detector");

   // define some materials
   TGeoMaterial *matVacuum = new TGeoMaterial("Vacuum", 0, 0, 0);
   TGeoMaterial *matAl     = new TGeoMaterial("Al", 26.98, 13, 2.7);
   TGeoMaterial *matFe     = new TGeoMaterial("Fe", 55.84, 26, 7.9);
      
   // define some media
   TGeoMedium *Vacuum = new TGeoMedium("Vacuum", 1, matVacuum);
   TGeoMedium *Al     = new TGeoMedium("Al", 2, matAl);
   TGeoMedium *Fe     = new TGeoMedium("Fe", 3, matFe);
      
   // create top volume
   TGeoVolume *top = geom->MakeBox("top", Vacuum, 600, 600, 2500);
   geom->SetTopVolume(top);
   geom->SetTopVisible(0);
   
   // first part of vacuum chamber up to veto station
   TGeoVolume *tub1 = geom->MakeTube("tub1", Al, 245, 250, 50);
   tub1->SetLineColor(18);  // silver/gray
   top->AddNode(tub1, 1, new TGeoTranslation(0, 0, -2450));
   
   // veto station
   TGeoBBox *detbox1 = new TGeoBBox("detbox1", 250, 250, 10);
   TGeoBBox *detbox2 = new TGeoBBox("detbox2", 245, 245, 10);
   
   TGeoCompositeShape *detcomp1 = new TGeoCompositeShape("detcomp1", "detbox1-detbox2");
   TGeoVolume *det1 = new TGeoVolume("det1", detcomp1,Vacuum);
   det1->SetLineColor(kRed);
   top->AddNode(det1, 1, new TGeoTranslation(0, 0, -2390));
   TGeoRotation r0;
   r0.SetAngles(15,0,0);
   TGeoTranslation t0(0, 0, -2370);
   TGeoCombiTrans c0(t0, r0);
   TGeoHMatrix *h0 = new TGeoHMatrix(c0);
   top->AddNode(det1, 11, h0);
   
   // second part of vacuum chamber up to first tracking station
   TGeoVolume *tub2 = geom->MakeTube("tub2", Al, 245, 250, 3880/2);  // 1890
   tub2->SetLineColor(18);
   top->AddNode(tub2, 1, new TGeoTranslation(0, 0, -440));
   
   // tracking station 1
   top->AddNode(det1, 2, new TGeoTranslation(0, 0, 1510));
   TGeoRotation r1;
   r1.SetAngles(15,0,0);
   TGeoTranslation t1(0, 0, 1530);
   TGeoCombiTrans c1(t1, r1);
   TGeoHMatrix *h1 = new TGeoHMatrix(c1);
   top->AddNode(det1, 3, h1);
   
   // third part of vacuum chamber up to second tracking station
   TGeoVolume *tub3 = geom->MakeTube("tub3", Al, 245, 250, 80);
   tub3->SetLineColor(18);
   top->AddNode(tub3, 1, new TGeoTranslation(0, 0, 1620));
   
   // tracking station 2
   top->AddNode(det1, 4, new TGeoTranslation(0, 0, 1710));
   TGeoRotation r2;
   r2.SetAngles(15,0,0);
   TGeoTranslation t2(0, 0, 1730);
   TGeoCombiTrans c2(t2, r2);
   TGeoHMatrix *h2 = new TGeoHMatrix(c2);
   top->AddNode(det1, 5, h2);
 
   // fourth part of vacuum chamber up to third tracking station and being covered by magnet
   TGeoVolume *tub4 = geom->MakeTube("tub4", Al, 245, 250, 200);
   tub4->SetLineColor(18);
   top->AddNode(tub4, 1, new TGeoTranslation(0, 0, 1940));

   // magnet yoke
   TGeoBBox *magyoke1 = new TGeoBBox("magyoke1", 350, 350, 125);
   TGeoBBox *magyoke2 = new TGeoBBox("magyoke2", 250, 250, 126);
   
   TGeoCompositeShape *magyokec = new TGeoCompositeShape("magyokec", "magyoke1-magyoke2");
   TGeoVolume *magyoke = new TGeoVolume("magyoke", magyokec, Fe);
   magyoke->SetLineColor(kBlue);
   //magyoke->SetTransparency(50);
   top->AddNode(magyoke, 1, new TGeoTranslation(0, 0, 1940));

   // magnet
   TGeoTubeSeg *magnet1a = new TGeoTubeSeg("magnet1a", 250, 300, 35, 45, 135);
   TGeoTubeSeg *magnet1b = new TGeoTubeSeg("magnet1b", 250, 300, 35, 45, 135);
   TGeoTubeSeg *magnet1c = new TGeoTubeSeg("magnet1c", 250, 270, 125, 45, 60);
   TGeoTubeSeg *magnet1d = new TGeoTubeSeg("magnet1d", 250, 270, 125, 120, 135);

   // magnet composite shape matrices   
   TGeoTranslation *m1 = new TGeoTranslation(0, 0, 160);
   m1->SetName("m1");
   m1->RegisterYourself();
   TGeoTranslation *m2 = new TGeoTranslation(0, 0, -160);
   m2->SetName("m2");
   m2->RegisterYourself();

   TGeoCompositeShape *magcomp1 = new TGeoCompositeShape("magcomp1", "magnet1a:m1+magnet1b:m2+magnet1c+magnet1d");
   TGeoVolume *magnet1 = new TGeoVolume("magnet1", magcomp1, Fe);
   magnet1->SetLineColor(kYellow);
   top->AddNode(magnet1, 1, new TGeoTranslation(0, 0, 1940));
//.........这里部分代码省略.........

void assembly()
{
//--- Definition of a simple geometry
    gSystem->Load("libGeom");
    TGeoManager *geom = new TGeoManager("Assemblies",
                                        "Geometry using assemblies");
    Int_t i;
    //--- define some materials
    TGeoMaterial *matVacuum = new TGeoMaterial("Vacuum", 0,0,0);
    TGeoMaterial *matAl = new TGeoMaterial("Al", 26.98,13,2.7);
//   //--- define some media
    TGeoMedium *Vacuum = new TGeoMedium("Vacuum",1, matVacuum);
    TGeoMedium *Al = new TGeoMedium("Aluminium",2, matAl);

    //--- make the top container volume
    TGeoVolume *top = geom->MakeBox("TOP", Vacuum, 1000., 1000., 100.);
    geom->SetTopVolume(top);

    // Make the elementary assembly of the whole structure
    TGeoVolume *tplate = new TGeoVolumeAssembly("TOOTHPLATE");

    Int_t ntooth = 5;
    Double_t xplate = 25;
    Double_t yplate = 50;
    Double_t xtooth = 10;
    Double_t ytooth = 0.5*yplate/ntooth;
    Double_t dshift = 2.*xplate + xtooth;
    Double_t xt,yt;

    TGeoVolume *plate = geom->MakeBox("PLATE", Al, xplate,yplate,1);
    plate->SetLineColor(kBlue);
    TGeoVolume *tooth = geom->MakeBox("TOOTH", Al, xtooth,ytooth,1);
    tooth->SetLineColor(kBlue);
    tplate->AddNode(plate,1);
    for (i=0; i<ntooth; i++) {
        xt = xplate+xtooth;
        yt = -yplate + (4*i+1)*ytooth;
        tplate->AddNode(tooth, i+1, new TGeoTranslation(xt,yt,0));
        xt = -xplate-xtooth;
        yt = -yplate + (4*i+3)*ytooth;
        tplate->AddNode(tooth, ntooth+i+1, new TGeoTranslation(xt,yt,0));
    }

    TGeoRotation *rot1 = new TGeoRotation();
    rot1->RotateX(90);
    TGeoRotation *rot;
    // Make a hexagone cell out of 6 tooth plates. These can zip together
    // without generating overlaps (they are self-contained)
    TGeoVolume *cell = new TGeoVolumeAssembly("CELL");
    for (i=0; i<6; i++) {
        Double_t phi =  60.*i;
        Double_t phirad = phi*TMath::DegToRad();
        Double_t xp = dshift*TMath::Sin(phirad);
        Double_t yp = -dshift*TMath::Cos(phirad);
        rot = new TGeoRotation(*rot1);
        rot->RotateZ(phi);
        cell->AddNode(tplate,i+1,new TGeoCombiTrans(xp,yp,0,rot));
    }

    // Make a row as an assembly of cells, then combine rows in a honeycomb
    // structure. This again works without any need to define rows as
    // "overlapping"
    TGeoVolume *row = new TGeoVolumeAssembly("ROW");
    Int_t ncells = 5;
    for (i=0; i<ncells; i++) {
        Double_t ycell = (2*i+1)*(dshift+10);
        row->AddNode(cell, ncells+i+1, new TGeoTranslation(0,ycell,0));
        row->AddNode(cell,ncells-i,new TGeoTranslation(0,-ycell,0));
    }

    Double_t dxrow = 3.*(dshift+10.)*TMath::Tan(30.*TMath::DegToRad());
    Double_t dyrow = dshift+10.;
    Int_t nrows = 5;
    for (i=0; i<nrows; i++) {
        Double_t xrow = 0.5*(2*i+1)*dxrow;
        Double_t yrow = 0.5*dyrow;
        if ((i%2)==0) yrow = -yrow;
        top->AddNode(row, nrows+i+1, new TGeoTranslation(xrow,yrow,0));
        top->AddNode(row, nrows-i, new TGeoTranslation(-xrow,-yrow,0));
    }

    //--- close the geometry
    geom->CloseGeometry();

    geom->SetVisLevel(4);
    geom->SetVisOption(0);
    top->Draw();
}

void event3D()
{
  gSystem->Load("libGeom");
   TGeoManager *geom = new TGeoManager("simple1", "Simple geometry");

   //--- define some materials
   TGeoMaterial *matVacuum = new TGeoMaterial("Vacuum", 0,0,0);
   TGeoMaterial *matAl = new TGeoMaterial("Al", 26.98,13,2.7);
//   //--- define some media
   TGeoMedium *Vacuum = new TGeoMedium("Vacuum",1, matVacuum);
   TGeoMedium *Al = new TGeoMedium("Root Material",2, matAl);

     //--- make the top container volume
   Double_t worldx = 110.;
   Double_t worldy = 50.;
   Double_t worldz = 5.;

   TGeoVolume *top = geom->MakeBox("TOP", Vacuum, 100., 100., 100.);
   geom->SetTopVolume(top);

   gGeoManager->GetVolume("TOP")->InvisibleAll();

   //Measurements (in cm)
   
   Double_t target_rad = 0.3;
   Double_t target_len = 48.58;

   Double_t RTPC_len = 41.;
   
   Double_t foil1_R1 = 2.;
   Double_t foil1_R2 = 2.0000018;

   Double_t foil2_R1 = 3.;
   Double_t foil2_R2 = 3.0000018;

   Double_t GEM1_R1 = 7.;
   Double_t GEM1_R2 = 7.0005;

   Double_t GEM2_R1 = 7.3;
   Double_t GEM2_R2 = 7.3005;

   Double_t GEM3_R1 = 7.6;
   Double_t GEM3_R2 = 7.6005;

   Double_t Z_offset = 0;

   Int_t trans_lvl = 70;
   
   TGeoTube *target = new TGeoTube("target_tube",0, target_rad, target_len/2);
   
   TGeoVolume *target_vol = new TGeoVolume("target_vol",target, Al); //(*)

   top->AddNode(target_vol,1,new TGeoTranslation(0,0, Z_offset));
   gGeoManager->GetVolume("target_vol")->SetTransparency(0);   

   TGeoTube *foil1 = new TGeoTube("foil1_tube", foil1_R1, foil1_R2, RTPC_len/2);
   
   TGeoVolume *foil1_vol = new TGeoVolume("foil1_vol", foil1, Al); //(*)

   top->AddNode(foil1_vol,1,new TGeoTranslation(0,0, Z_offset));


   TGeoTube *foil2 = new TGeoTube("foil2_tube", foil2_R1, foil2_R2, RTPC_len/2);
   
   TGeoVolume *foil2_vol = new TGeoVolume("foil2_vol", foil2, Al); //(*)

   top->AddNode(foil2_vol,1,new TGeoTranslation(0,0, Z_offset));


   TGeoTube *gem1 = new TGeoTube("gem1_tube", GEM1_R1, GEM1_R2, RTPC_len/2);
   TGeoVolume *gem1_vol = new TGeoVolume("gem1_vol", gem1, Al); //(*)
   gem1_vol->SetLineColor(kOrange);
   top->AddNode(gem1_vol,1,new TGeoTranslation(0,0, Z_offset));

   TGeoTube *gem2 = new TGeoTube("gem2_tube", GEM2_R1, GEM2_R2, RTPC_len/2);
   TGeoVolume *gem2_vol = new TGeoVolume("gem2_vol", gem2, Al); //(*)
   gem2_vol->SetLineColor(kGreen);
   top->AddNode(gem2_vol,1,new TGeoTranslation(0,0, Z_offset));

   TGeoTube *gem3 = new TGeoTube("gem3_tube", GEM3_R1, GEM3_R2, RTPC_len/2);
   TGeoVolume *gem3_vol = new TGeoVolume("gem3_vol", gem3, Al); //(*)
   top->AddNode(gem3_vol,1,new TGeoTranslation(0,0, Z_offset));
   gem3_vol->SetLineColor(kRed);
   gGeoManager->GetVolume("gem3_vol")->SetTransparency(trans_lvl);   

   
   //--- draw the ROOT box.
   // by default the picture will appear in the standard ROOT TPad.
   //if you have activated the following line in system.rootrc,
   //it will appear in the GL viewer
   //#Viewer3D.DefaultDrawOption:   ogl
   
   geom->CloseGeometry();
   
   geom->SetVisLevel(4);
  
   //   top->SetVisibility(kFALSE);
   top->Draw("ogl");
   gWorld = top;
}

//.........这里部分代码省略.........
  TGeoVolume* infinite_plate =
    geom->MakeTube( "Infinite_Plate", med_1, 0.0, plate_radius, plate_half_length );

  infinite_plate->SetUniqueID( 3 );
  infinite_plate->SetLineColor( 4 );

  //--------------------------------------------------------
  // Create the vacuum (R = "infinite" cm, h = "infinite")
  //--------------------------------------------------------

  // Set a radius of 7 cm ( 167 x the CSDA range = 0.4814;
  double vacuum_radius = 7.0;
  // Set a vacuum length large enough to cover the whole geometry
  double vacuum_half_length = 2.0*plate_half_length + 5.0 + calorimeter_thickness + 0.5;
  // Add region
  TGeoVolume* void_region =
    geom->MakeTube( "Void", void_med, 0.0, vacuum_radius, vacuum_half_length );

  void_region->SetUniqueID( 4 );
  void_region->SetLineColor( 6 );
  void_region->SetVisibility( true );

  //--------------------------------------------------------
  // Create the graveyard (encompasses the entire geometry)
  //--------------------------------------------------------

  // Add region
  TGeoVolume* graveyard_region =
    geom->MakeTube( "Graveyard", graveyard_med, 0.0, vacuum_radius + 0.5, vacuum_half_length + 0.5 );

  graveyard_region->SetUniqueID( 5 );
  graveyard_region->SetLineColor( 7 );
  graveyard_region->SetVisContainers( true );

//---------------------------------------------------------------------------//
// Heirarchy (Volume) Definitions
//---------------------------------------------------------------------------//

// Set graveyard as the top volume
geom->SetTopVolume(graveyard_region);

  //------------------------
  // Add vacuum to graveyard
  //------------------------
  // Get z position of the vacuum
  double vacuum_z = vacuum_half_length;
  // Add region
  graveyard_region->AddNode(void_region, 0 );

    //-----------------------------------
    // Add the calorimeter foil to vacuum
    //-----------------------------------

    // Set the front of the calorimeter at a z of 5 cm
    double calorimeter_z = 5.0 + half_cal_thickness;
    // Add the front calorimeter foil to vacuum
    void_region->AddNode( calorimeter_foil, 1, new TGeoTranslation(0,0,calorimeter_z) );

    //-----------------------------------
    // Add the front foil to vacuum
    //-----------------------------------

    if( front_thickness > 0.0 )
    {
      // Set the front foil 0.1 cm in front of the calorimeter
      double front_z = 4.9 - front_thickness/2.0;
      // Add the region
      void_region->AddNode( front_foil, 2, new TGeoTranslation(0,0,front_z) );
    }

    //-----------------------------------
    // Add the infinite plate to vacuum
    //-----------------------------------

    // Set the infinite plate 0.1 cm behind the  calorimeter
    double plate_z = 5.1 + calorimeter_thickness + plate_half_length;
    // Add region
    void_region->AddNode( infinite_plate, 3, new TGeoTranslation(0,0,plate_z) );


//---------------------------------------------------------------------------//
// Export and Drawing Capabilities
//---------------------------------------------------------------------------//

  // Close the geometry
  geom->CloseGeometry();

  // // Draw the geometry
  // geom->SetTopVisible();
  // graveyard_region->Draw();

  // TView *view = gPad->GetView();
  // view->ShowAxis();

  // Export the geometry
  geom->Export( "geom_Al.root" );

  // Finished
  exit(0);
}

void lego()
{
  // Drawing a figure, made of lego block, using ROOT geometry class.
  // Name: lego.C
  // Author: Soon Gi Kwon([email protected]), Dept. of Physics, Univ. of Seoul
  // Reviewed by Sunman Kim ([email protected])
  // Supervisor: Prof. Inkyu Park ([email protected])
  //
  // How to run: .x lego.C in ROOT terminal, then use OpenGL
  //
  // This macro was created for the evaluation of Computational Physics course in 2006.
  // We thank to Prof. Inkyu Park for his special lecture on ROOT and to all of ROOT team
  //

   TGeoManager *geom = new TGeoManager("geom","My first 3D geometry");


   TGeoMaterial *vacuum=new TGeoMaterial("vacuum",0,0,0);
   TGeoMaterial *Fe=new TGeoMaterial("Fe",55.845,26,7.87);



   TGeoMedium *Air=new TGeoMedium("Vacuum",0,vacuum);
   TGeoMedium *Iron=new TGeoMedium("Iron",1,Fe);


 // create volume
   TGeoVolume *top=geom->MakeBox("top",Air,100,100,100);
   geom->SetTopVolume(top);
   geom->SetTopVisible(0);
   // If you want to see the boundary, please input the number, 1 instead of 0.
   // Like this, geom->SetTopVisible(1);


//----------------------------------------------------------------------

TGeoVolume *ha1=geom->MakeSphere("ha1",Iron,0,10,80,90,0,360);
 ha1->SetLineColor(41);
top->AddNodeOverlap(ha1,1,new TGeoCombiTrans(0,0,4,new TGeoRotation("ha1",0,0,0)));

TGeoVolume *ha2=geom->MakeSphere("ha2",Iron,0,7,90,180,0,360);
 ha2->SetLineColor(41);
top->AddNodeOverlap(ha2,1,new TGeoCombiTrans(0,0,4,new TGeoRotation("ha2",0,180,0)));

TGeoVolume *ha3=geom->MakeSphere("ha3",Iron,0,7.3,80,90,0,360);
 ha3->SetLineColor(2);
top->AddNodeOverlap(ha3,1,new TGeoCombiTrans(0,0,4.8,new TGeoRotation("ha3",0,0,0)));


TGeoVolume *h1=geom->MakeTubs("h1",Iron,0,6,4.5,0,0);
 h1->SetLineColor(5);
top->AddNodeOverlap(h1,1,new TGeoCombiTrans(0,0,0,new TGeoRotation("h1",0,0,0)));

TGeoVolume *h2=geom->MakeSphere("h2",Iron,0,7.5,0,52.5,0,360);
 h2->SetLineColor(5);
top->AddNodeOverlap(h2,1,new TGeoCombiTrans(0,0,0,new TGeoRotation("h2",0,0,0)));

TGeoVolume *h3=geom->MakeSphere("h3",Iron,0,7.5,0,52.5,0,360);
 h3->SetLineColor(5);
top->AddNodeOverlap(h3,1,new TGeoCombiTrans(0,0,0,new TGeoRotation("h3",180,180,0)));

TGeoVolume *h4=geom->MakeTubs("h4",Iron,2.5,3.5,1.5,0,0);
 h4->SetLineColor(5);
top->AddNodeOverlap(h4,1,new TGeoCombiTrans(0,0,7.5,new TGeoRotation("h4",0,0,0)));



TGeoVolume *t1_1=geom->MakeTubs("t1_1",Iron,0,0.8,1,0,360);
 t1_1->SetLineColor(12);
top->AddNodeOverlap(t1_1,1,new TGeoCombiTrans(-5,2,1.5,new TGeoRotation("t1_1",-90,90,0)));

TGeoVolume *t2_1=geom->MakeTubs("t2_1",Iron,0,0.8,1,0,360);
 t2_1->SetLineColor(12);
top->AddNodeOverlap(t2_1,1,new TGeoCombiTrans(-5,-2,1.5,new TGeoRotation("t2_1",-90,90,0)));

TGeoVolume *fb1=geom->MakeTubs("fb1",Iron,2,2.3,1,100,260);
 fb1->SetLineColor(12);
top->AddNodeOverlap(fb1,1,new TGeoCombiTrans(-5,0,-1,new TGeoRotation("fb1",90,90,90)));



TGeoVolume *m1=geom->MakeBox("m1",Iron,7,8,4);
 m1->SetLineColor(2);
top->AddNodeOverlap(m1,1,new TGeoCombiTrans(0,0,-17,new TGeoRotation("m1",90,90,0)));

TGeoVolume *m2=geom->MakeTubs("m2",Iron,0,1,7,90,180);
 m2->SetLineColor(2);
top->AddNodeOverlap(m2,1,new TGeoCombiTrans(-3,0,-9,new TGeoRotation("m2",0,90,0)));

TGeoVolume *m3=geom->MakeTubs("m3",Iron,0,1,7,0,90);
 m3->SetLineColor(2);
top->AddNodeOverlap(m3,1,new TGeoCombiTrans(3,0,-9,new TGeoRotation("m3",0,90,0)));

TGeoVolume *m4=geom->MakeBox("m4",Iron,3,7,0.5);
 m4->SetLineColor(2);
top->AddNodeOverlap(m4,1,new TGeoCombiTrans(0,0,-8.5,new TGeoRotation("m4",90,0,90)));

TGeoVolume *m5=geom->MakeTubs("m5",Iron,0,1.5,1.2,0,0);
 m5->SetLineColor(5);
top->AddNodeOverlap(m5,1,new TGeoCombiTrans(0,0,-7.8,new TGeoRotation("m5",0,0,0)));
//.........这里部分代码省略.........

void glViewerLOD(Int_t reqNodes = 1000, Bool_t randomDist = kTRUE,
                 Bool_t reqSpheres = kTRUE, Bool_t reqTubes = kTRUE)
{
   TGeoManager * geom = new TGeoManager("LODTest", "GL viewer LOD test");
   geom->SetNsegments(4); // Doesn't matter keep low
   TGeoMaterial *matEmptySpace = new TGeoMaterial("EmptySpace", 0, 0, 0);
   TGeoMaterial *matSolid      = new TGeoMaterial("Solid"    , .938, 1., 10000.);

   TGeoMedium *medEmptySpace = new TGeoMedium("Empty", 1, matEmptySpace);
   TGeoMedium *medSolid      = new TGeoMedium("Solid", 1, matSolid);

   Double_t sizeBase = 20.0;
   Double_t worldRadius;
   if (randomDist) {
      worldRadius = pow(reqNodes,.5)*sizeBase;
   } else {
      worldRadius = pow(reqNodes,.3)*sizeBase;
   }

   TGeoVolume *top = geom->MakeBox
      ("WORLD", medEmptySpace, worldRadius, worldRadius, worldRadius);
   geom->SetTopVolume(top);

   gRandom->SetSeed();

   // Create random number of unique sphere shapes - up to 25% of
   // total placed sphere requested
   UInt_t volumeCount = gRandom->Integer(reqNodes/4)+1;
   TGeoVolume ** volumes = new TGeoVolume[volumeCount];
   TGeoVolume * volume;

   Double_t dummy;

   for (UInt_t i = 0; i < volumeCount; i++) {
      char name[128];
      sprintf(name, "Volume_%d", i);

      // Random volume shape
      Int_t type = -1;
      if (reqSpheres && reqTubes) {
         type = gRandom->Integer(2);
         if (type == 1)
            type += gRandom->Integer(3);
      }
      else if(reqSpheres)
         type = 0;
      else if(reqTubes)
         type = 1 + gRandom->Integer(3);

      // Random dimensions
      Double_t rMin = gRandom->Rndm() * sizeBase;
      Double_t rMax = rMin + gRandom->Rndm() * sizeBase * 2.0;
      Double_t dz   = pow(gRandom->Rndm(),2.0) * sizeBase * 15.0;
      Double_t phi1 = gRandom->Rndm() * 90.0;
      Double_t phi2 = phi1 + gRandom->Rndm() * 270.0;

      // Pick random color (not black)
      Int_t color = gRandom->Integer(50);
      if (color == kBlack) color += 1;

      switch (type) {
        case 0: {
            // GL viewer only supports solid spheres (0. inner radius)
            volumes[i] = geom->MakeSphere(name,  medSolid,  0., rMax);
            printf("Volume %d : Color %d, Sphere, Radius %f\n", i, color, rMax);
            break;
         }
         case 1: {
            volumes[i] = geom->MakeTube(name,  medSolid,  rMin, rMax, dz);
            printf("Volume %d : Color %d, Tube, Inner Radius %f, "
                   "Outer Radius %f, Length %f\n",
                   i, color, rMin, rMax, dz);
            break;
         }
         case 2: {
            volumes[i] = geom->MakeTubs(name,  medSolid,  rMin, rMax, dz,
                                        phi1, phi2);
            printf("Volume %d : Color %d, Tube Seg, Inner Radius %f, "
                   "Outer Radius %f, Length %f, Phi1 %f, Phi2 %f\n",
                   i, color, rMin, rMax, dz, phi1, phi2);
            break;
         }
         case 3: {
            Double_t n1[3], n2[3];
            n1[0] = gRandom->Rndm()*.5;
            n1[1] = gRandom->Rndm()*.5; n1[2] = -1.0 + gRandom->Rndm()*.5;
            n2[0] = gRandom->Rndm()*.5;
            n2[1] = gRandom->Rndm()*.5; n2[2] =  1.0 - gRandom->Rndm()*.5;

            volumes[i] = geom->MakeCtub(name,  medSolid,  rMin, rMax, dz,
                                        phi1, phi2, n1[0], n1[1], n1[2],
                                        n2[0], n2[1], n2[2]);
            printf("Volume %d : Color %d, Cut Tube, Inner Radius %f, "
                   "Outer Radius %f, Length %f, Phi1 %f, Phi2 %f, "
                   "n1 (%f,%f,%f), n2 (%f,%f,%f)\n",
                   i, color, rMin, rMax, dz, phi1, phi2,
                   n1[0], n1[1], n1[2], n2[0], n2[1], n2[2]);
            break;
         }
         default: {
//.........这里部分代码省略.........

TGeoVolume *VP_SimpleECal(int nphi = 4, int nz = 3, double density = 8.28)
{
  const double world_x = 9000.;
  const double world_y = 9000.;
  const double world_z = 16000;

  TGeoManager *geom = new TGeoManager("SimpleECal", "Simplification of the CMS ECal");

  TGeoMaterial *world_mat = new TGeoMaterial("world", 1, 2, 0);
  TGeoMedium *world_med   = new TGeoMedium("world_medium", 0, world_mat);
  TGeoVolume *world       = geom->MakeBox("top", world_med, world_x, world_y, world_z);

  geom->SetTopVolume(world);
  geom->SetTopVisible(1);

  int crystal_nphi    = nphi;
  int crystal_nz      = nz;
  double ecal_density = density;

  const int crystal_n    = crystal_nphi * crystal_nz;
  const double ecal_zmin = -3000.;
  const double ecal_zmax = 3000.;

  const double ecal_rmin = 10.;
  const double ecal_rmax = 5000.;
  const double ecal_dz   = 0.5 * (ecal_zmax - ecal_zmin) / crystal_nz;
  const double ecal_sphi = 0.;
  //    const G4double ecal_dphi  =  2.0*M_PI/crystal_nphi;
  // G4 seems to be in radian while TGeo seems to be in degree.
  // const double ecal_dphi  =  2.0*TMath::Pi()/crystal_nphi;
  const double ecal_dphi = 2.0 * 180 / crystal_nphi;

  int iptr = 0;

  TGeoElementTable *table = gGeoManager->GetElementTable();
  // TGeoElement* elPb = new TGeoElement( "Lead", "Pb", 82., 207.19*g/mole );
  // TGeoElement* elW = new TGeoElement( "Tungstenm", "W",74., 183.85*g/mole);
  // TGeoElement* elO = new TGeoElement( "Oxygen", "O2", 8., 16.*g/mole );
  TGeoElement *elPb = table->GetElement(82);
  TGeoElement *elW  = table->GetElement(74);
  TGeoElement *elO  = table->GetElement(8);
  // TGeoMaterial *ecal_mat = new TGeoMaterial("ecal",90,120,density);
  TGeoMixture *ecal_mat = new TGeoMixture("ecal_mat", 3, density);
  ecal_mat->AddElement(elPb, 1);
  ecal_mat->AddElement(elW, 1);
  ecal_mat->AddElement(elO, 4);
  TGeoMedium *ecal_med = new TGeoMedium("ecal_med", 0, ecal_mat);

  for (int j = 0; j < crystal_nz; ++j) {
    for (int i = 0; i < crystal_nphi; ++i) {

      iptr = i + j * crystal_nphi;

      TGeoVolume *ecal = geom->MakeTubs(TString::Format("ecal-%d-%d", j, i), ecal_med, ecal_rmin, ecal_rmax, ecal_dz,
                                        ecal_sphi + i * ecal_dphi, ecal_sphi + (i + 1) * ecal_dphi);
      ecal->SetLineColor(iptr);
      // top->AddNode(ecal,1,new TGeoCombiTrans(0,0,0,new TGeoRotation("ecal",0,0,0)));

      // GPThreeVector ecal_trans = GPThreeVector_create(0,0,ecal_zmin+(2.0*j+1.0)*ecal_dz);
      double dx                   = 0.0;
      double dy                   = 0.0;
      double dz                   = ecal_zmin + (2.0 * j + 1.0) * ecal_dz;
      TGeoTranslation *ecal_trans = new TGeoTranslation("", dx, dy, dz);

      //         GPLogicalVolume_Constructor(ecal_log+iptr, (GPVSolid*)ecal, ecal_mat);
      // GPVPhysicalVolume_Constructor(ecal_phy+iptr, idRot, ecal_trans, ecal_log+iptr);
      world->AddNode(ecal, iptr, ecal_trans);

      // Set mother
      //         GPVPhysicalVolume_SetMotherLogical(ecal_phy+iptr, world_log);

      //         addLogicalVolumePointers( ecal_log+iptr);
      //         addPhysicalVolumePointers( ecal_phy+iptr);
    }
  }

  // add daughter volume
  //   for ( int j=0; j < crystal_nz ; ++j ) {
  //      for ( int i=0; i < crystal_nphi ; ++i ) {
  //         iptr = i+j*crystal_nphi;
  //         addLogicalVolumeDaughter( world_log, ecal_phy+iptr);
  //      }
  //   }

  // Register world volume pointers for relocation
  //   addLogicalVolumePointers( world_log );
  //   addPhysicalVolumePointers( world_phy );

  geom->CloseGeometry();

  return world;
}

void nucleus(Int_t nProtons  = 40,Int_t  nNeutrons = 60) 
{
   Double_t NeutronRadius = 60,
            ProtonRadius = 60,
            NucleusRadius,
            distance = 60;
   Double_t vol = nProtons + nNeutrons;
   vol = 3 * vol / (4 * TMath::Pi());

   NucleusRadius = distance * TMath::Power(vol, 1./3.);
//   cout << "NucleusRadius: " << NucleusRadius << endl;

   TGeoManager * geom = new TGeoManager("nucleus", "Model of a nucleus");
   geom->SetNsegments(40);
   TGeoMaterial *matEmptySpace = new TGeoMaterial("EmptySpace", 0, 0, 0);
   TGeoMaterial *matProton     = new TGeoMaterial("Proton"    , .938, 1., 10000.);
   TGeoMaterial *matNeutron    = new TGeoMaterial("Neutron"   , .935, 0., 10000.);

   TGeoMedium *EmptySpace = new TGeoMedium("Empty", 1, matEmptySpace);
   TGeoMedium *Proton     = new TGeoMedium("Proton", 1, matProton);
   TGeoMedium *Neutron    = new TGeoMedium("Neutron",1, matNeutron);

//  the space where the nucleus lives (top container volume)

   Double_t worldx = 200.;
   Double_t worldy = 200.;
   Double_t worldz = 200.;
 
   TGeoVolume *top = geom->MakeBox("WORLD", EmptySpace, worldx, worldy, worldz); 
   geom->SetTopVolume(top);

   TGeoVolume * proton  = geom->MakeSphere("proton",  Proton,  0., ProtonRadius); 
   TGeoVolume * neutron = geom->MakeSphere("neutron", Neutron, 0., NeutronRadius); 
   proton->SetLineColor(kRed);
   neutron->SetLineColor(kBlue);

   Double_t x, y, z, dummy;
   Int_t i = 0; 
   while ( i<  nProtons) {
      gRandom->Rannor(x, y);
      gRandom->Rannor(z,dummy);
      if ( TMath::Sqrt(x*x + y*y + z*z) < 1) {
         x = (2 * x - 1) * NucleusRadius;
         y = (2 * y - 1) * NucleusRadius;
         z = (2 * z - 1) * NucleusRadius;
         top->AddNode(proton, i, new TGeoTranslation(x, y, z));
         i++;
      }
   }
   i = 0; 
   while ( i <  nNeutrons) {
      gRandom->Rannor(x, y);
      gRandom->Rannor(z,dummy);
      if ( TMath::Sqrt(x*x + y*y + z*z) < 1) {
         x = (2 * x - 1) * NucleusRadius;
         y = (2 * y - 1) * NucleusRadius;
         z = (2 * z - 1) * NucleusRadius;
         top->AddNode(neutron, i + nProtons, new TGeoTranslation(x, y, z));
         i++;
      }
   }
   geom->CloseGeometry();
   geom->SetVisLevel(4);
   top->Draw("ogl");
}

//.........这里部分代码省略.........
  geoFace->setMediaFile(medFile);
  geoFace->readMedia();
  TGeoManager* gGeoMan = gGeoManager;
  // --------------------------------------------------------------------------



  // -----------------   Get and create the required media    -----------------
  FairGeoMedia*   geoMedia = geoFace->getMedia();
  FairGeoBuilder* geoBuild = geoLoad->getGeoBuilder();
  
  // ---> pipe medium
  FairGeoMedium* fPipeMedium = geoMedia->getMedium(pipeMediumName.Data());
  TString fairError = "FairMedium " + pipeMediumName + " not found";
  if ( ! fPipeMedium ) Fatal("Main", fairError.Data());
  geoBuild->createMedium(fPipeMedium);
  TGeoMedium* pipeMedium = gGeoMan->GetMedium(pipeMediumName.Data());
  TString geoError = "Medium " + pipeMediumName + " not found";
  if ( ! pipeMedium ) Fatal("Main", geoError.Data());
  // ---> vacuum
  FairGeoMedium* mVacuum = geoMedia->getMedium("vacuum");
  if ( ! mVacuum ) Fatal("Main", "FairMedium vacuum not found");
  geoBuild->createMedium(mVacuum);
  TGeoMedium* vacuum = gGeoMan->GetMedium("vacuum");
  if ( ! vacuum ) Fatal("Main", "Medium vacuum not found");
  // --------------------------------------------------------------------------



  // --------------   Create geometry and top volume  -------------------------
  gGeoMan = (TGeoManager*)gROOT->FindObject("FAIRGeom");
  gGeoMan->SetName("PIPEgeom");
  TGeoVolume* top = new TGeoVolumeAssembly("TOP");
  gGeoMan->SetTopVolume(top);
  TGeoVolume* pipe = new TGeoVolumeAssembly(pipeName.Data());
  // --------------------------------------------------------------------------


  // -----   Create sections  -------------------------------------------------
  infoFile << endl << "Beam pipe section: " << pipe1name << endl;
  infoFile << setw(2) << "i" << setw(10) << "Z,mm" << setw(10) << "Rin,mm" << setw(10) << "Rout,mm" << setw(10) << "h,mm" << endl;
  //*
  TGeoVolume* pipe1    = MakePipe  (1, nSects1,  z1,  rin1,  rout1,  pipeMedium, &infoFile); 
  pipe1->SetLineColor(kGray);
  pipe->AddNode(pipe1, 0);
  
  infoFile << endl << "Beam pipe section: " << pipe2name << endl;
  infoFile << setw(2) << "i" << setw(10) << "Z,mm" << setw(10) << "Rin,mm" << setw(10) << "Rout,mm" << setw(10) << "h,mm" << endl;
  TGeoVolume* pipe2    = MakePipe  (2, nSects2,  z2,  rin2,  rout2,  pipeMedium, &infoFile); 
  pipe2->SetLineColor(kBlue);
  pipe->AddNode(pipe2, 0);
  TGeoVolume* pipevac1 = MakeVacuum(1, nSects01, z01, rin01, rout01, vacuum,     &infoFile); 
  pipevac1->SetLineColor(kCyan);
  pipe->AddNode(pipevac1, 0);
  
  infoFile << endl << "Beam pipe section: " << pipe3name << endl;
  infoFile << setw(2) << "i" << setw(10) << "Z,mm" << setw(10) << "Rin,mm" << setw(10) << "Rout,mm" << setw(10) << "h,mm" << endl;
  TGeoVolume* pipe3    = MakePipe  (3, nSects3,  z3,  rin3,  rout3,  pipeMedium, &infoFile); 
  pipe3->SetLineColor(kGreen);
  pipe->AddNode(pipe3, 0);
  TGeoVolume* pipevac2 = MakeVacuum(2, nSects02, z02, rin02, rout02, vacuum,     &infoFile); 
  pipevac2->SetLineColor(kCyan);
  pipe->AddNode(pipevac2, 0);
  //*/
  
  // define some rotation & translation for pipe4-pipe7 & pipevac3-pipevac5

int main(int argc, char * argv[])
{

   // generate benchmark cases
   TransformationMatrix const * identity = new TransformationMatrix(0,0,0,0,0,0);

   double L = 10.;
   double Lz = 10.;
   const double Sqrt2 = sqrt(2.);

   TGeoManager * geom = new TGeoManager("simple1", "ToyDetector");
   TGeoMaterial * matVacuum = new TGeoMaterial("Vacuum",0,0,0);

   TGeoMedium * vac = new TGeoMedium("Vacuum",1,matVacuum);
   TGeoVolume * world = geom->MakeBox("world", vac, L, L, Lz  );

   geom->SetTopVolume( world );

   TGeoVolume * boxlevel2 = geom->MakeBox( "boxlevel2", vac, Sqrt2*L/2./2., Sqrt2*L/2./2.,Lz);
   TGeoVolume * boxlevel3 = geom->MakeBox( "boxlevel3", vac, L/2./2.,L/2./2.,Lz);
   TGeoVolume * boxlevel1 = geom->MakeBox( "boxlevel1", vac, L/2.,L/2., Lz);

   boxlevel2->AddNode( boxlevel3, 1, new TGeoRotation("mat1",0,0,-45));
   boxlevel1->AddNode( boxlevel2, 1, new TGeoRotation("mat2",0,0,45));
   world->AddNode( boxlevel1, 1, new TGeoTranslation("trans1",-L/2., 0,0 ) );
   world->AddNode( boxlevel1, 2, new TGeoTranslation("trans2",L/2., 0,0 ) );
   geom->CloseGeometry();
   delete world->GetVoxels();
   world->SetVoxelFinder(0);

    // perform basic tests
   TGeoNavigator * nav = geom->GetCurrentNavigator(); // new TGeoNavigator(geom);
   StopWatch timer;

   timer.Start();
   int stepsdone=0;
   for(int n=0;n<1000;n++)
     {
   for(int i=0;i<100000;i++)
   // testing the NavigationAndStepInterface
   {
      int localstepsdone=0;
      double distancetravelled=0.;
      Vector3D p;
      PhysicalVolume::samplePoint( p, L, L, Lz, 1. );

      //std::cerr << p << std::endl;
      //setup point in world
      Vector3D d(1,0,0);

      TGeoNode const * vol;
      nav->SetCurrentPoint( p.x, p.y, p.z );
      nav->SetCurrentDirection( d.x, d.y, d.z);
      vol=nav->FindNode(p.x,p.y,p.z);

      while( vol!=NULL )
      {
         localstepsdone++;

         // do one step ( this will internally adjust the current point and so on )
         vol = nav->FindNextBoundaryAndStep(Utils::kInfinity);

         distancetravelled+=nav->GetStep();
         //double const * p = nav->GetCurrentPoint();
         //double const * pl = nav->GetLastPoint();
         //double const * cd = nav->GetCurrentDirection();

         //std::cerr << " proposed step: " << nav->GetStep();
         //std::cerr << " current point " << p[0] << " " << p[1] << " " << p[2] << std::endl;
         //std::cerr << " last point " << pl[0] << " " << pl[1] << " " << pl[2] << std::endl;
         //std::cerr << " current dir " << cd[0] << " " << cd[1] << " " << cd[2] << std::endl;
      }
   //   std::cerr << localstepsdone << " " << distancetravelled << std::endl;
      stepsdone+=localstepsdone;
   }
     }
   timer.Stop();
   std::cout << " time for 100000 particles " << timer.getDeltaSecs( ) << std::endl;
   std::cout << " average steps done " << stepsdone / 100000. << std::endl;
   std::cout << " time per step " << timer.getDeltaSecs()/stepsdone << std::endl;
}

/*
 * CreateLambdaGeometry.C
 *
 *  Created on: Mar 21, 2013
 *      Author: stockman
 */
startdetector1()
{
	  //-----------------------------
	Double_t positionOfDisk1inZ = 40.0;		//in cm
	//Double_t positionOfDisk2inZ = 53.0;		//in cm

	  //--------------------------------------------------------------------
	  gROOT->Macro("$VMCWORKDIR/gconfig/rootlogon.C");
	  
	  // Load this libraries
	  gSystem->Load("libGeoBase");
	  gSystem->Load("libParBase");
	  gSystem->Load("libBase");
	  gSystem->Load("libPndData");
	  gSystem->Load("libPassive");
	  
	  TString outfile= "../../geometry/startdetector.root";
	  TFile* fi = new TFile(outfile,"RECREATE");  
	  
	  FairGeoLoader* geoLoad = new FairGeoLoader("TGeo","FairGeoLoader");
	  FairGeoInterface *geoFace = geoLoad->getGeoInterface();
	  geoFace->setMediaFile("../../geometry/media_pnd.geo");	
	  geoFace->readMedia();
	  geoFace->print();

	  FairGeoMedia *Media =  geoFace->getMedia();
	  FairGeoBuilder *geobuild=geoLoad->getGeoBuilder();

	  FairGeoMedium *CbmMediumHYPdiamond  = Media->getMedium("HYPdiamond");
	  FairGeoMedium *CbmMediumCarbon  = Media->getMedium("carbon");
	  FairGeoMedium *CbmMediumCarbonFoam  = Media->getMedium("carbonfoam");
	  FairGeoMedium *CbmMediumAluminium = Media->getMedium("aluminium");
	  
	  Int_t nmed=geobuild->createMedium(CbmMediumHYPdiamond);
	  nmed=geobuild->createMedium(CbmMediumCarbon);
	  nmed=geobuild->createMedium(CbmMediumCarbonFoam);
	  nmed=geobuild->createMedium(CbmMediumAluminium);

	  TGeoManager* gGeoMan = (TGeoManager*)gROOT->FindObject("FAIRGeom");

	  TGeoVolume *top = new TGeoVolumeAssembly("top");
	                    
	  gGeoMan->SetTopVolume(top);
	  
	TGeoVolume *combinedLambdaDisks = new TGeoVolumeAssembly("CombinedLambdaDisks");
	TGeoVolumeAssembly* lambdaDisk = new TGeoVolumeAssembly("LambdaDisk");
	TGeoVolumeAssembly* largeRing = new TGeoVolumeAssembly("LargeRing");
	TGeoVolumeAssembly* smallRing = new TGeoVolumeAssembly("SmallRing");

//	lambdaDisk->AddNode(largeSensorVolume,0,trc1);
//	TGeoShape* largeSensorShape  = new TGeoTrd1("StripSensorActiveLargeShape", 3.4593/2, 1.9936/2,0.03/2, 5.5665/2);
	Double_t points [16];
	points[0] = -1.7297;		points[1] = 2.7833;
	points[2] = 1.7297; 		points[3] = 2.7833;
	points[4] = 0.9968; 		points[5] = -2.7833;
	points[6] = -0.9968; 		points[7] = -2.7833;

	points[8] = -1.7297; 		points[9] = 2.7833;
	points[10] = 1.7297; 		points[11] = 2.7833;
	points[12] = 0.9968; 		points[13] = -2.7833;
	points[14] = -0.9968; 		points[15] = -2.7833;

	TGeoShape* largeSensorShape  = new TGeoArb8("StripSensorActiveLargeShape", 0.03/2, points);
	TGeoVolume* largeSensorVolume = new TGeoVolume("StripSensorActiveLargeTrap",largeSensorShape,gGeoMan->GetMedium("HYPdiamond"));

//	TGeoRotation rotSensor("rotSensor",0,90,0);
	TGeoRotation rotSensor("rotSensor",0,0,0);
	TGeoTranslation transLargeSensor(0,(14.8669 + 5.5665)/2,0);
	TGeoCombiTrans combined1(transLargeSensor, rotSensor);

	for (int i = 0; i < 12; i++){
		TString rotName("rotLargeSens");
		rotName.Append(i+1);
		TGeoRotation rotSens(rotName.Data(), 0,0, i*2 * 360/24);
		TGeoHMatrix sens = rotSens * combined1;
		largeRing->AddNode(largeSensorVolume, i, new TGeoHMatrix(sens));
	}

	Double_t pointsSmall [16];
	pointsSmall[0] = -1.9579/2;		pointsSmall[1] = 5.3358/2;
	pointsSmall[2] = 1.9579/2; 		pointsSmall[3] = 5.3358/2;
	pointsSmall[4] = 0.5529/2; 		pointsSmall[5] = -5.3358/2;
	pointsSmall[6] = -0.5529/2; 		pointsSmall[7] = -5.3358/2;

	pointsSmall[8] = -1.9579/2;		pointsSmall[9] = 5.3358/2;
	pointsSmall[10] = 1.9579/2; 	pointsSmall[11] = 5.3358/2;
	pointsSmall[12] = 0.5529/2; 	pointsSmall[13] = -5.3358/2;
	pointsSmall[14] = -0.5529/2; 	pointsSmall[15] = -5.3358/2;

//	TGeoShape* smallSensorShape = new TGeoTrd1("StripSensorActiveSmallShape", 1.9579/2, 0.5529/2, 0.03/2, 5.3358/2);
	TGeoShape* smallSensorShape = new TGeoArb8("StripSensorActiveSmallShape", 0.03/2, pointsSmall);
	TGeoVolume* smallSensorVolume = new TGeoVolume("StripSensorActiveSmallTrap", smallSensorShape, gGeoMan->GetMedium("HYPdiamond"));

	TGeoTranslation transSmallSensor(0,5.3358/2+2.1,0);
//.........这里部分代码省略.........