source: trunk/examples/extended/geometry/olap/src/OlapDetConstr.cc @ 1346

Last change on this file since 1346 was 1337, checked in by garnier, 14 years ago

tag geant4.9.4 beta 1 + modifs locales

File size: 20.6 KB
Line 
1//
2// ********************************************************************
3// * License and Disclaimer                                           *
4// *                                                                  *
5// * The  Geant4 software  is  copyright of the Copyright Holders  of *
6// * the Geant4 Collaboration.  It is provided  under  the terms  and *
7// * conditions of the Geant4 Software License,  included in the file *
8// * LICENSE and available at  http://cern.ch/geant4/license .  These *
9// * include a list of copyright holders.                             *
10// *                                                                  *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work  make  any representation or  warranty, express or implied, *
14// * regarding  this  software system or assume any liability for its *
15// * use.  Please see the license in the file  LICENSE  and URL above *
16// * for the full disclaimer and the limitation of liability.         *
17// *                                                                  *
18// * This  code  implementation is the result of  the  scientific and *
19// * technical work of the GEANT4 collaboration.                      *
20// * By using,  copying,  modifying or  distributing the software (or *
21// * any work based  on the software)  you  agree  to acknowledge its *
22// * use  in  resulting  scientific  publications,  and indicate your *
23// * acceptance of all terms of the Geant4 Software license.          *
24// ********************************************************************
25//
26//
27// $Id: OlapDetConstr.cc,v 1.4 2006/06/29 17:22:52 gunter Exp $
28// GEANT4 tag $Name: geant4-09-04-beta-01 $
29//
30//
31// --------------------------------------------------------------
32// OlapDetConstr
33//
34// Author: Martin Liendl - Martin.Liendl@cern.ch
35//
36// --------------------------------------------------------------
37//
38#include "OlapDetConstr.hh"
39#include "OlapGenerator.hh"
40
41#include "SolidAnalyser.hh"
42#include "OlapManager.hh"
43
44#include "G4UImanager.hh"
45#include "G4VVisManager.hh"
46#include "G4VPhysicalVolume.hh"
47#include "G4PVPlacement.hh"
48#include "G4PVReplica.hh"
49#include "G4PVParameterised.hh"
50#include "G4LogicalVolume.hh"
51#include "G4RunManager.hh"
52#include "G4VisExtent.hh"
53#include "G4VisAttributes.hh"
54#include "G4LogicalVolumeStore.hh"
55#include "G4Box.hh"
56#include "G4Polyline.hh"
57#include "G4Polycone.hh"
58#include "G4Polyhedra.hh"
59#include "G4ThreeVector.hh"
60
61// debugging:
62//#define OLAPDEBUG1
63//#define OlapDetConstr_debug
64
65OlapDetConstr::OlapDetConstr(G4VUserDetectorConstruction * aGeometry,
66                             G4VPhysicalVolume *aWorld) :
67  theTheta(0),
68  thePhi(0),
69  theAlpha(0),
70  theNewWorldRot(0),
71  theFullGeometry(aGeometry), 
72  theWorld(aWorld),
73  theNewWorld(0),
74  theNewLV(0),
75  theNewWorldBox(new G4Box("world",1.,1.,1.)),
76  syncVis(true),
77  nrLV(0)
78{ 
79   //theNewWorldLV = new G4LogicalVolume
80   visMother =   new G4VisAttributes(G4Colour(7.0, 0.0, 0.0, 0.4));   // red
81   visMother->SetForceWireframe(true);
82   visDaughterA = new G4VisAttributes(G4Colour(.0, 0.0, 1.0, 0.3)); // blue
83   visDaughterA->SetForceWireframe(false);
84   visDaughterB = new G4VisAttributes(G4Colour(0., 0., 1.0, 0.3)); // light-blue
85   visDaughterB->SetForceWireframe(false);
86   visWorld =    new G4VisAttributes(G4Colour(1., 0, 0, 1.)); // red
87   visWorld->SetForceWireframe(false);
88   visWorld->SetVisibility(true);
89   
90   visFullWorld = new G4VisAttributes(G4Colour(1.,1.,1.,0.9));
91   visFullWorld->SetForceWireframe(true);
92   visWorld->SetVisibility(false);
93   
94   visInvisible = new G4VisAttributes(G4Colour(1.,1.,1.,0.9));
95   visInvisible->SetForceWireframe(true);
96   visInvisible->SetVisibility(false);
97   
98   visFirstLevel = new G4VisAttributes(G4Colour(1.,1.,1,1));
99   visFirstLevel->SetForceWireframe(true);
100   visFirstLevel->SetVisibility(true);
101}
102
103
104OlapDetConstr::~OlapDetConstr()
105{
106   //delete theNewWorld;
107   //delete theNewLV; not the owner!
108   //delete theFullGeometry;
109   //delete visMother;
110   //delete visDaughterA;
111   //delete visDaughterB;
112   //delete visWorld;
113   
114   //delete visFullWorld;
115   //delete visInvisible;
116   //delete visFirstLevel;
117}
118
119
120void OlapDetConstr::SetRotation(G4double t, G4double p, G4double a)
121{
122   theTheta = t;
123   thePhi = p;
124   theAlpha = a;
125}
126
127
128G4VPhysicalVolume * OlapDetConstr::Construct()
129{
130   if (!theWorld)
131        theWorld = theFullGeometry->Construct();
132   if (!theWorld)
133   {
134     G4cerr << "OlapDetConstr::Construct(): could not create full world." << G4endl;
135     G4cerr << "                            exiting ..." << G4endl;
136     G4Exception("ERROR - OlapDetConstr::Construct()");
137   } 
138   
139   // logical volume which will have a box with its dimensions derived from
140   // the extent of the 'new-mother' volume positioned inside
141   // the size of this NewWorld will be set in 'SetNewWorld(..)'
142   theNewWorldLV = 
143      new G4LogicalVolume(theNewWorldBox,
144                          theWorld->GetLogicalVolume()->GetMaterial(),
145                          "NewWorld"
146                          );
147
148   // physical volume which will serve as the 'small' NewWorld to be
149   // used for overlap checking
150   theNewWorld =
151      new G4PVPlacement( 0,               // rotation
152                         G4ThreeVector(), // translation
153                            "NewWorld", // name of phys = name of logical!!!
154                         theNewWorldLV,     // own logical vol
155                         0, false, 0      // no mother, not MANY, cpNr=0
156                       ); 
157   
158   theNewWorldLV->SetVisAttributes(visWorld);                       
159
160   //ResetColors();
161   
162   
163   theNewLV = theWorld->GetLogicalVolume();
164   nrLV = G4LogicalVolumeStore::GetInstance()->size();
165   G4cout << "OlapDetConstr::Construct(): nr of lvs=" << nrLV << G4endl;
166   return theWorld;
167}
168
169
170G4VPhysicalVolume * OlapDetConstr::SetNewWorld(G4LogicalVolume * aMotherLV,
171                                               G4bool debugFlag)
172{
173   if (debugFlag)
174   {
175      G4cout << "mother: " << aMotherLV->GetName() << " : " 
176             << aMotherLV->GetSolid()->GetName() << G4endl;
177      return 0;
178   }
179   
180   //ML: ResetColors(theNewLV);
181   theNewLV = aMotherLV;
182   
183   ConstructNewWorld();
184   
185   // set new world to Geant4-Kernel
186   G4RunManager::GetRunManager()->DefineWorldVolume(theNewWorld);
187   
188   // try to set the OlapGenerator   
189   const OlapGenerator * aGen = dynamic_cast<const OlapGenerator *>
190         (G4RunManager::GetRunManager()->GetUserPrimaryGeneratorAction());
191   if ( aGen )
192   {
193     OlapGenerator * aNonConstGen = const_cast<OlapGenerator*>(aGen);
194       
195     // 'optimal' extent
196     // aNonConstGen->SetExtent(theNewLV->GetSolid()->GetExtent());
197       
198     // extent which allows arbitrary rotations ...
199     G4Box * bx = dynamic_cast<G4Box*>
200                (theNewWorld->GetLogicalVolume()->GetSolid());
201     aNonConstGen->SetExtent(2.*bx->GetXHalfLength());
202   }
203   else
204   {
205      G4cerr << "Warning: Primary generator is not OlapGenerator!" << G4endl
206                   << "Overlap Detection will not work!" << G4endl;
207   }
208   return theNewWorld;
209}
210
211
212G4VPhysicalVolume * OlapDetConstr::GetNewWorld()
213{
214   return theNewWorld; 
215}
216
217
218G4VPhysicalVolume * OlapDetConstr::GetFullWorld()
219{
220   return theWorld;
221}
222
223
224void OlapDetConstr::DrawPolyOutline()
225{
226   SolidAnalyser::GetSolidAnalyser();
227   G4VSolid * solid = theNewLV->GetSolid();
228   G4Polyline line, endline;
229   G4double ex(0), first_r(0), first_z(0);
230   G4double maxz(0), maxr(0);
231   G4int nr(0); 
232   G4Polycone  * s  = dynamic_cast<G4Polycone*>(solid);
233   G4Polyhedra * sh = dynamic_cast<G4Polyhedra*>(solid);
234   if (s) 
235   { 
236      nr = s->GetNumRZCorner();
237     
238      first_r = s->GetCorner(0).r;
239      first_z = s->GetCorner(0).z;
240      for (G4int i=0; i<nr; i++)
241      {
242          G4double r = s->GetCorner(i).r;
243          G4double z = s->GetCorner(i).z;
244          //line.push_back(G4ThreeVector(r,0.,z)); 
245          line.push_back(G4ThreeVector(r,z,0)); 
246          if (r>maxr)
247            maxr=r;
248          if (z>maxz)
249            maxz=z; 
250          G4cout << G4ThreeVector(r,0,z) << G4endl;
251          //endline.push_back(G4ThreeVector(s->GetCorner(nr-1).r, 0.,
252          //                                s->GetCorner(nr-1).z));
253          endline.push_back(G4ThreeVector(s->GetCorner(nr-1).r,
254                                          s->GetCorner(nr-1).z,0));
255       }
256    }
257    else if (sh) 
258    { 
259      nr = sh->GetNumRZCorner();
260
261      first_r = sh->GetCorner(0).r;
262      first_z = sh->GetCorner(0).z;
263      for (G4int i=0; i<nr; i++)
264      {
265          G4double r = sh->GetCorner(i).r;
266          G4double z = sh->GetCorner(i).z;
267          line.push_back(G4ThreeVector(r,0.,z)); 
268          if (r>maxr)
269            maxr=r;
270          if (z>maxz)
271            maxz=z; 
272          G4cout << G4ThreeVector(r,0,z) << G4endl;
273          endline.push_back(G4ThreeVector(sh->GetCorner(nr-1).r, 0.,
274                                          sh->GetCorner(nr-1).z));
275      }
276    }
277    else
278      return;
279
280     //endline.push_back(G4ThreeVector(first_r,0.,first_z));
281     endline.push_back(G4ThreeVector(first_r,first_z,0.));
282     if (maxz>=maxr)
283       ex = maxz;
284     else
285       ex = maxr;
286     
287     G4cout << "ex=" << ex << G4endl; 
288     theNewWorldBox->SetXHalfLength(ex + ex/30.);
289     theNewWorldBox->SetYHalfLength(ex + ex/30.);
290     theNewWorldBox->SetZHalfLength(ex + ex/30.);
291   
292   G4Colour colour;
293   colour = G4Colour(1.,0.,0.);
294   line.SetVisAttributes(theNewLV->GetVisAttributes()->GetColor());
295   G4VVisManager::GetConcreteInstance()->Draw(line);
296   colour = G4Colour(0.,1.,0.);
297   endline.SetVisAttributes(theNewLV->GetVisAttributes()->GetColor());
298   //endline.SetVisAttributes(attribs3);
299   G4VVisManager::GetConcreteInstance()->Draw(endline);
300   colour = G4Colour(0.,0.,1.);
301   G4Polyline ax;
302   //ax.push_back(G4ThreeVector(0.,0.,-ex));
303   //ax.push_back(G4ThreeVector(0.,0.,ex));
304   ax.push_back(G4ThreeVector(0.,-ex,0.));
305   ax.push_back(G4ThreeVector(0.,ex,0.));
306   G4VisAttributes attribs2(colour);
307   ax.SetVisAttributes(attribs2);
308   G4VVisManager::GetConcreteInstance()->Draw(ax);
309}
310
311void OlapDetConstr::ConstructNewWorld()
312{ 
313   // delete the current NewWorld
314   DeleteNewWorld();
315
316   #ifdef OlapDetConstr_debug
317   G4int checkb = G4LogicalVolumeStore::GetInstance()->size(); 
318     G4cerr << "OlapDetConstr::ConstructNewWorld():"
319            << " nr of lvs after DeleteNewWorld is: " << checkb << G4endl;
320   #endif
321   
322   #ifdef somewhen_but_not_now
323   OlapManager * m = OlapManager::GetOlapManager();
324   // If polymode is set, don't create a overlap-detection world
325   // but draw the outlines of the mother polycone/-hedra (if it is one)
326   if (m->IsPolyMode())
327   {
328      DrawPolyOutline();
329      return;
330   }
331   #endif
332   G4int nr = theNewLV->GetNoDaughters();
333   
334   // create a NewWorld with a 2% larger extent than the
335   // solid of the NewMother
336   
337   const G4VisExtent & ext = theNewLV->GetSolid()->GetExtent();
338   G4double extX = (ext.GetXmax()-ext.GetXmin())/2.;
339   G4double extY = (ext.GetYmax()-ext.GetYmin())/2.;
340   G4double extZ = (ext.GetZmax()-ext.GetZmin())/2.;
341   G4ThreeVector pos( (ext.GetXmax()+ext.GetXmin())/2. ,
342                      (ext.GetYmax()+ext.GetYmin())/2. ,
343                      (ext.GetZmax()+ext.GetZmin())/2.
344                     ); 
345   pos=-pos;                     
346 
347   /*
348   extX += extX/1000.;
349   extY += extY/1000.;
350   extZ += extZ/1000.;
351   */
352   G4double worldDim = std::sqrt(extX*extX + extY*extY + extZ*extZ);
353   worldDim += worldDim/100.;
354   // this automatically sets the dimensions of the
355   // solid in theNewWorldLV, because it's a ptr to the same solid
356   theNewWorldBox->SetXHalfLength(worldDim);
357   theNewWorldBox->SetYHalfLength(worldDim);
358   theNewWorldBox->SetZHalfLength(worldDim);
359     
360   //G4LogicalVolume * aInbetweenLV =
361   // new G4LogicalVolume(theNewWorldBox, theNewLV->GetMaterial(), "INBETWEEN");
362   
363   //aInbetweenLV->SetVisAttributes(visWorld);
364                           
365   // create the mother with one layer of daughters
366   G4LogicalVolume * aNewMotherLV =
367     new G4LogicalVolume(theNewLV->GetSolid(),
368                         theNewLV->GetMaterial(),
369                         theNewLV->GetName()
370                         );
371   
372   //LM
373   // below the mother is fixed colored instead as taken from vis-atts of the
374   // full detector
375   //aNewMotherLV->SetVisAttributes(visMother);
376   aNewMotherLV->SetVisAttributes(theNewLV->GetVisAttributes());
377
378#ifdef OLAPDEBUG2
379   G4cout << " mother: " << theNewLV->GetName() << G4endl;
380#endif         
381 
382   delete theNewWorldRot;
383   G4ThreeVector rotAxis(std::cos(thePhi)*std::sin(theTheta),
384                         std::sin(thePhi)*std::sin(theTheta),
385                         std::cos(theTheta));
386   
387   theNewWorldRot = new G4RotationMatrix(rotAxis,theAlpha);
388   
389   /*
390   new G4PVPlacement( theNewWorldRot, G4ThreeVector(), aInbetweenLV,
391                      theNewLV->GetName(),
392                      theNewWorldLV, false, 0);
393   */
394                                                                                                     
395   // G4VPhysicalVolume * aNewMotherPV =
396   new G4PVPlacement( theNewWorldRot,   // rotation
397                      pos, // translation
398                      aNewMotherLV,     // own logical vol
399                      theNewLV->GetName(), // name of phys = name of logical!!!
400                     
401                      theNewWorldLV,    // mother logical vol
402                      //aInbetweenLV,
403                     
404                      false, 0      // not MANY, cpNr=0
405                      ); 
406
407 
408  // loop over daughters an copy them into the new world ...
409  //      - but don't iterate recursively over deeper layers of daughters!
410  for (G4int i=0; i<nr; i++) {
411   
412    G4VPhysicalVolume * pv = theNewLV->GetDaughter(i);
413    #ifdef OLAPDEBUG2
414      G4cout << i << " ";
415    #endif 
416
417    // ---------------------------------------------------------------
418    // only do the stuff, if OlapManager.NoOlapMap[lv] is set to false!
419    // otherwise the volume should not be considered for olap detection
420   
421    OlapManager * m = OlapManager::GetOlapManager();
422    if (m->NoOlapMap[pv->GetLogicalVolume()])
423      continue; 
424    // ---------------------------------------------------------------               
425   
426    G4LogicalVolume * myLogical = 
427      new G4LogicalVolume(pv->GetLogicalVolume()->GetSolid(),
428                          pv->GetLogicalVolume()->GetMaterial(),
429                          pv->GetLogicalVolume()->GetName()
430                          );
431                         
432  //LM                         
433  // take vis atts from whole detector instead of fixed values!
434  myLogical->SetVisAttributes(pv->GetLogicalVolume()->GetVisAttributes());
435     
436    // select between placements, replicas and parameterizations
437    G4PVPlacement * pvPlace; 
438    if (dynamic_cast<G4PVPlacement*>(pv)) 
439    {
440       pvPlace = dynamic_cast<G4PVPlacement*>(pv);
441       new G4PVPlacement( pvPlace->GetRotation(),
442                          pvPlace->GetTranslation(),
443                          myLogical,     // my own logical vol
444                          pvPlace->GetName(),
445                          aNewMotherLV,  // my mother logical vol
446                          false,         // not many
447                          pvPlace->GetCopyNo()
448                        ); 
449       #ifdef OLAPDEBUG2
450         G4cout << "  daughter: " << i << " " << pvPlace->GetName()
451                << " lv=" << pv->GetLogicalVolume()->GetName() 
452                << " solid=" << pv->GetLogicalVolume()->GetSolid()->GetEntityType() << G4endl;
453       #endif
454       
455    } 
456    else if ( pv->IsReplicated() && ! pv->GetParameterisation() ) // Replicas
457    {
458       // retrieve replication data
459       EAxis aAxis;
460       G4int nReplicas;
461       G4double aWidth, anOffset;
462       G4bool aDummy;
463       pv->GetReplicationData(aAxis, nReplicas, aWidth, anOffset, aDummy);
464       
465       new G4PVReplica( pv->GetName(),
466                        myLogical,
467                        aNewMotherLV,
468                        aAxis,
469                        nReplicas,
470                        aWidth,
471                        anOffset
472                      );       
473                       
474       #ifdef OLAPDEBUG2
475         G4cout << "  repl.daughter: " << i << " " << pv->GetName()
476                << " lv=" << pv->GetLogicalVolume()->GetName() 
477                << " solid=" << pv->GetLogicalVolume()->GetSolid()->GetEntityType() << G4endl;         
478       #endif                 
479           
480    }
481    else if ( pv->IsReplicated() && pv->GetParameterisation() )
482    {
483       // retrieve replication/parameterisation data
484       EAxis aAxis;
485       G4int nReplicas;
486       G4double aWidth, anOffset;
487       G4bool aDummy;
488       pv->GetReplicationData(aAxis, nReplicas, aWidth, anOffset, aDummy);
489       G4VPVParameterisation *pParam = pv->GetParameterisation();
490       
491       new G4PVParameterised( pv->GetName(),
492                              myLogical,
493                              aNewMotherLV,
494                              aAxis,
495                              nReplicas,
496                              pParam
497                            );
498       #ifdef OLAPDEBUG2
499         G4cout << "  param.daughter: " << i << " " << pv->GetName() 
500                << " lv=" << pv->GetLogicalVolume()->GetName() 
501                << " solid="
502                << pv->GetLogicalVolume()->GetSolid()->GetEntityType()
503                << G4endl;
504       #endif                 
505   
506    }// end:placement?replica?parameterisation
507   
508  } // end:Loop over Daughters                                   
509
510   #ifdef OlapDetConstr_debug
511   G4int check = G4LogicalVolumeStore::GetInstance()->size(); 
512   G4cerr << "OlapDetConstr::ConstructNewWorld():"
513          << " nr of lvs after construction is: " << check
514          << G4endl << G4endl;           
515   #endif
516   
517}
518
519
520void OlapDetConstr::DeleteNewWorld()
521{
522   #ifdef OlapDetConstr_debug
523   G4int checkb = G4LogicalVolumeStore::GetInstance()->size(); 
524   G4cerr << "OlapDetConstr::DeleteNewWorld():"
525            << " nr of lvs before delete is: " << checkb << G4endl;
526   #endif
527
528   G4int nr = theNewWorldLV->GetNoDaughters();
529   
530   if (nr==0)
531   {
532     G4cerr << "OlapDetConstr::DeleteNewWorld(): no daughter in NewWorld!"
533            << G4endl << G4endl;
534     return;
535   } 
536     
537   if (nr>1 || nr<0)  // only one daughter is allowed!!!!
538   { 
539      G4cerr << "OlapDetConstr::DeleteNewWorld(): "
540             << "too many daughters in NewWorldLV!" << nr << G4endl;
541      G4cerr << "                                 exiting ..." << G4endl;
542      G4Exception("ERROR - OlapDetConstr::DeleteNewWorld()");
543   }
544     
545   G4LogicalVolume * aMother;
546   G4VPhysicalVolume * aDaughter;
547   
548   aMother = theNewWorldLV->GetDaughter(0)->GetLogicalVolume();
549   
550   G4int nrd = aMother->GetNoDaughters();
551   for (G4int i=0; i<nrd; i++)
552   {
553     aDaughter = aMother->GetDaughter(i);
554     G4LogicalVolume * tmp = aDaughter->GetLogicalVolume();
555     //aMother->RemoveDaughter(aDaughter);
556     #ifdef OlapDetConstr_debug
557       G4cerr << "    deleting: i=" << i << " of " << nrd
558              << " : vol=" << tmp->GetName() << G4endl;
559     #endif
560     delete tmp;
561     delete aDaughter;
562   } 
563   
564   aDaughter = theNewWorldLV->GetDaughter(0);
565   theNewWorldLV->RemoveDaughter(aDaughter);
566   delete aDaughter;
567   delete aMother;
568
569   #ifdef OlapDetConstr_debug
570   G4int check = G4LogicalVolumeStore::GetInstance()->size(); 
571   if (check!=nrLV) {
572     G4cerr << "OlapDetConstr::DeleteNewWorld():"
573            << " nr of lvs should be: " << nrLV << " but is: "
574            << check << G4endl;
575   }           
576   #endif
577}
578
579
580// set everything in the full geometry to invisible
581void OlapDetConstr::ResetColors()
582{
583/*
584   G4ColorMap * aColMap = G4ColorMap::GetColorMap();
585   aColMap->SetAllInvisible();   
586*/   
587   //ColorFirstLevel();   
588}
589
590
591// colors the first layer of daughters of the full geometry
592void OlapDetConstr::ColorFirstLevel()
593{
594/* 
595   G4ColorMap * aColMap = G4ColorMap::GetColorMap();
596
597   G4LogicalVolume * lvWorld = theWorld->GetLogicalVolume();
598   aColMap->SetVisible(lvWorld,true);
599   G4int nr = lvWorld->GetNoDaughters();
600   for (G4int i = 0; i<nr; i++) {
601     aColMap->SetVisible(lvWorld->GetDaughter(i)->GetLogicalVolume(),true);
602   } 
603*/   
604}
605
606
607// sets NewWorld invisible
608void OlapDetConstr::ResetColors(G4LogicalVolume*)
609{
610/*
611   G4ColorMap * aColMap = G4ColorMap::GetColorMap();
612   G4int nr = lv->GetNoDaughters();
613   for (G4int i =0; i<nr; i++) {
614     aColMap->SetVisible(lv->GetDaughter(i)->GetLogicalVolume(),false);
615   }
616   aColMap->SetVisible(lv,false);
617*/   
618   //ColorFirstLevel();
619
620}
621
622
623G4VPhysicalVolume * OlapDetConstr::SetFullWorld()
624{
625   G4RunManager::GetRunManager()->DefineWorldVolume(theWorld);
626   
627   return theWorld;
628}
629
630
631G4LogicalVolume * OlapDetConstr::GetOriginalWorld()
632{
633   return theNewLV;
634}   
Note: See TracBrowser for help on using the repository browser.