00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033 #include "G4VSceneHandler.hh"
00034
00035 #include "G4ios.hh"
00036 #include <sstream>
00037
00038 #include "G4VisManager.hh"
00039 #include "G4VGraphicsSystem.hh"
00040 #include "G4VViewer.hh"
00041 #include "G4VSolid.hh"
00042 #include "G4RotationMatrix.hh"
00043 #include "G4ThreeVector.hh"
00044 #include "G4VPhysicalVolume.hh"
00045 #include "G4Material.hh"
00046 #include "G4Polyline.hh"
00047 #include "G4Scale.hh"
00048 #include "G4Text.hh"
00049 #include "G4Circle.hh"
00050 #include "G4Square.hh"
00051 #include "G4Polymarker.hh"
00052 #include "G4Polyhedron.hh"
00053 #include "G4NURBS.hh"
00054 #include "G4Visible.hh"
00055 #include "G4VisAttributes.hh"
00056 #include "G4VModel.hh"
00057 #include "G4TrajectoriesModel.hh"
00058 #include "G4Box.hh"
00059 #include "G4Cons.hh"
00060 #include "G4Tubs.hh"
00061 #include "G4Trd.hh"
00062 #include "G4Trap.hh"
00063 #include "G4Sphere.hh"
00064 #include "G4Para.hh"
00065 #include "G4Torus.hh"
00066 #include "G4Polycone.hh"
00067 #include "G4Polyhedra.hh"
00068 #include "G4DisplacedSolid.hh"
00069 #include "G4LogicalVolume.hh"
00070 #include "G4PhysicalVolumeModel.hh"
00071 #include "G4ModelingParameters.hh"
00072 #include "G4VTrajectory.hh"
00073 #include "G4VTrajectoryPoint.hh"
00074 #include "G4HitsModel.hh"
00075 #include "G4VHit.hh"
00076 #include "G4VDigi.hh"
00077 #include "G4ScoringManager.hh"
00078 #include "G4DefaultLinearColorMap.hh"
00079 #include "Randomize.hh"
00080 #include "G4StateManager.hh"
00081 #include "G4RunManager.hh"
00082 #include "G4Run.hh"
00083 #include "G4Transform3D.hh"
00084 #include "G4AttHolder.hh"
00085 #include "G4AttDef.hh"
00086 #include "G4PhysicalConstants.hh"
00087
00088 G4VSceneHandler::G4VSceneHandler (G4VGraphicsSystem& system, G4int id, const G4String& name):
00089 fSystem (system),
00090 fSceneHandlerId (id),
00091 fViewCount (0),
00092 fpViewer (0),
00093 fpScene (0),
00094 fMarkForClearingTransientStore (true),
00095
00096
00097
00098 fReadyForTransients (true),
00099 fProcessingSolid (false),
00100 fProcessing2D (false),
00101 fpModel (0),
00102 fNestingDepth (0),
00103 fpVisAttribs (0)
00104 {
00105 G4VisManager* pVMan = G4VisManager::GetInstance ();
00106 fpScene = pVMan -> GetCurrentScene ();
00107 if (name == "") {
00108 std::ostringstream ost;
00109 ost << fSystem.GetName () << '-' << fSceneHandlerId;
00110 fName = ost.str();
00111 }
00112 else {
00113 fName = name;
00114 }
00115 fTransientsDrawnThisEvent = pVMan->GetTransientsDrawnThisEvent();
00116 fTransientsDrawnThisRun = pVMan->GetTransientsDrawnThisRun();
00117 }
00118
00119 G4VSceneHandler::~G4VSceneHandler () {
00120 G4VViewer* last;
00121 while( ! fViewerList.empty() ) {
00122 last = fViewerList.back();
00123 fViewerList.pop_back();
00124 delete last;
00125 }
00126 }
00127
00128 void G4VSceneHandler::PreAddSolid (const G4Transform3D& objectTransformation,
00129 const G4VisAttributes& visAttribs) {
00130 fObjectTransformation = objectTransformation;
00131 fpVisAttribs = &visAttribs;
00132 fProcessingSolid = true;
00133 }
00134
00135 void G4VSceneHandler::PostAddSolid () {
00136 fpVisAttribs = 0;
00137 fProcessingSolid = false;
00138 if (fReadyForTransients) {
00139 fTransientsDrawnThisEvent = true;
00140 fTransientsDrawnThisRun = true;
00141 }
00142 }
00143
00144 void G4VSceneHandler::BeginPrimitives
00145 (const G4Transform3D& objectTransformation) {
00146
00147
00148 fNestingDepth++;
00149 if (fNestingDepth > 1)
00150 G4Exception
00151 ("G4VSceneHandler::BeginPrimitives",
00152 "visman0101", FatalException,
00153 "Nesting detected. It is illegal to nest Begin/EndPrimitives.");
00154 fObjectTransformation = objectTransformation;
00155 }
00156
00157 void G4VSceneHandler::EndPrimitives () {
00158 if (fNestingDepth <= 0)
00159 G4Exception("G4VSceneHandler::EndPrimitives",
00160 "visman0102", FatalException, "Nesting error.");
00161 fNestingDepth--;
00162 if (fReadyForTransients) {
00163 fTransientsDrawnThisEvent = true;
00164 fTransientsDrawnThisRun = true;
00165 }
00166 }
00167
00168 void G4VSceneHandler::BeginPrimitives2D
00169 (const G4Transform3D& objectTransformation) {
00170 fNestingDepth++;
00171 if (fNestingDepth > 1)
00172 G4Exception
00173 ("G4VSceneHandler::BeginPrimitives2D",
00174 "visman0103", FatalException,
00175 "Nesting detected. It is illegal to nest Begin/EndPrimitives.");
00176 fObjectTransformation = objectTransformation;
00177 fProcessing2D = true;
00178 }
00179
00180 void G4VSceneHandler::EndPrimitives2D () {
00181 if (fNestingDepth <= 0)
00182 G4Exception("G4VSceneHandler::EndPrimitives2D",
00183 "visman0104", FatalException, "Nesting error.");
00184 fNestingDepth--;
00185 if (fReadyForTransients) {
00186 fTransientsDrawnThisEvent = true;
00187 fTransientsDrawnThisRun = true;
00188 }
00189 fProcessing2D = false;
00190 }
00191
00192 void G4VSceneHandler::BeginModeling () {
00193 }
00194
00195 void G4VSceneHandler::EndModeling ()
00196 {
00197 fpModel = 0;
00198 }
00199
00200 void G4VSceneHandler::ClearStore () {}
00201
00202 void G4VSceneHandler::ClearTransientStore () {}
00203
00204 void G4VSceneHandler::AddSolid (const G4Box& box) {
00205 RequestPrimitives (box);
00206
00207
00208
00209
00210
00211
00212
00213
00214
00215
00216
00217 }
00218
00219 void G4VSceneHandler::AddSolid (const G4Tubs& tubs) {
00220 RequestPrimitives (tubs);
00221 }
00222
00223 void G4VSceneHandler::AddSolid (const G4Cons& cons) {
00224 RequestPrimitives (cons);
00225 }
00226
00227 void G4VSceneHandler::AddSolid (const G4Trd& trd) {
00228 RequestPrimitives (trd);
00229 }
00230
00231 void G4VSceneHandler::AddSolid (const G4Trap& trap) {
00232 RequestPrimitives (trap);
00233 }
00234
00235 void G4VSceneHandler::AddSolid (const G4Sphere& sphere) {
00236 RequestPrimitives (sphere );
00237 }
00238
00239 void G4VSceneHandler::AddSolid (const G4Para& para) {
00240 RequestPrimitives (para);
00241 }
00242
00243 void G4VSceneHandler::AddSolid (const G4Torus& torus) {
00244 RequestPrimitives (torus);
00245 }
00246
00247 void G4VSceneHandler::AddSolid (const G4Polycone& polycone) {
00248 RequestPrimitives (polycone);
00249 }
00250
00251 void G4VSceneHandler::AddSolid (const G4Polyhedra& polyhedra) {
00252 RequestPrimitives (polyhedra);
00253 }
00254
00255 void G4VSceneHandler::AddSolid (const G4VSolid& solid) {
00256 RequestPrimitives (solid);
00257 }
00258
00259 void G4VSceneHandler::AddCompound (const G4VTrajectory& traj) {
00260 G4TrajectoriesModel* trajectoriesModel =
00261 dynamic_cast<G4TrajectoriesModel*>(fpModel);
00262 if (!trajectoriesModel) G4Exception
00263 ("G4VSceneHandler::AddCompound(const G4VTrajectory&)",
00264 "visman0105", FatalException, "Not a G4TrajectoriesModel.");
00265 else {
00266 if (trajectoriesModel->IsDrawingModeSet()) {
00267 traj.DrawTrajectory(trajectoriesModel->GetDrawingMode());
00268 } else {
00269 traj.DrawTrajectory();
00270 }
00271 }
00272 }
00273
00274 void G4VSceneHandler::AddCompound (const G4VHit& hit) {
00275
00276 const_cast<G4VHit&>(hit).Draw();
00277 }
00278
00279 void G4VSceneHandler::AddCompound (const G4VDigi& digi) {
00280
00281 const_cast<G4VDigi&>(digi).Draw();
00282 }
00283
00284 void G4VSceneHandler::AddCompound (const G4THitsMap<G4double>& hits) {
00285
00286 G4bool scoreMapHits = false;
00287 G4ScoringManager* scoringManager = G4ScoringManager::GetScoringManagerIfExist();
00288 if (scoringManager) {
00289 size_t nMeshes = scoringManager->GetNumberOfMesh();
00290 for (size_t iMesh = 0; iMesh < nMeshes; ++iMesh) {
00291 G4VScoringMesh* mesh = scoringManager->GetMesh(iMesh);
00292 if (mesh && mesh->IsActive()) {
00293 MeshScoreMap scoreMap = mesh->GetScoreMap();
00294 for(MeshScoreMap::const_iterator i = scoreMap.begin();
00295 i != scoreMap.end(); ++i) {
00296 const G4String& scoreMapName = i->first;
00297 const G4THitsMap<G4double>* foundHits = i->second;
00298 if (foundHits == &hits) {
00299 G4DefaultLinearColorMap colorMap("G4VSceneHandlerColorMap");
00300 scoreMapHits = true;
00301 mesh->DrawMesh(scoreMapName, &colorMap);
00302 }
00303 }
00304 }
00305 }
00306 }
00307 if (scoreMapHits) {
00308 static G4bool first = true;
00309 if (first) {
00310 first = false;
00311 G4cout <<
00312 "Scoring map drawn with default parameters."
00313 "\n To get gMocren file for gMocren browser:"
00314 "\n /vis/open gMocrenFile"
00315 "\n /vis/viewer/flush"
00316 "\n Many other options available with /score/draw... commands."
00317 "\n You might want to \"/vis/viewer/set/autoRefresh false\"."
00318 << G4endl;
00319 }
00320 } else {
00321
00322 const_cast<G4THitsMap<G4double>&>(hits).DrawAllHits();
00323 }
00324 }
00325
00326 void G4VSceneHandler::AddViewerToList (G4VViewer* pViewer) {
00327 fViewerList.push_back (pViewer);
00328 }
00329
00330 void G4VSceneHandler::AddPrimitive (const G4Scale& scale) {
00331
00332 const G4double margin(0.01);
00333
00334
00335 const G4double oneMinusMargin (1. - margin);
00336
00337 const G4VisExtent& sceneExtent = fpScene->GetExtent();
00338
00339
00340 const G4double length(scale.GetLength());
00341 const G4double halfLength(length / 2.);
00342 const G4double tickLength(length / 20.);
00343 const G4double piBy2(halfpi);
00344
00345
00346 const G4double xmin = sceneExtent.GetXmin();
00347 const G4double xmax = sceneExtent.GetXmax();
00348 const G4double ymin = sceneExtent.GetYmin();
00349 const G4double ymax = sceneExtent.GetYmax();
00350 const G4double zmin = sceneExtent.GetZmin();
00351 const G4double zmax = sceneExtent.GetZmax();
00352
00353
00354 G4Polyline scaleLine, tick11, tick12, tick21, tick22;
00355 G4VisAttributes visAtts(*scale.GetVisAttributes());
00356 scaleLine.SetVisAttributes(&visAtts);
00357 tick11.SetVisAttributes(&visAtts);
00358 tick12.SetVisAttributes(&visAtts);
00359 tick21.SetVisAttributes(&visAtts);
00360 tick22.SetVisAttributes(&visAtts);
00361
00362
00363
00364 G4Point3D r1(G4Point3D(-halfLength, 0., 0.));
00365 G4Point3D r2(G4Point3D( halfLength, 0., 0.));
00366 scaleLine.push_back(r1);
00367 scaleLine.push_back(r2);
00368 G4Point3D ticky(0., tickLength, 0.);
00369 G4Point3D tickz(0., 0., tickLength);
00370 tick11.push_back(r1 + ticky);
00371 tick11.push_back(r1 - ticky);
00372 tick12.push_back(r1 + tickz);
00373 tick12.push_back(r1 - tickz);
00374 tick21.push_back(r2 + ticky);
00375 tick21.push_back(r2 - ticky);
00376 tick22.push_back(r2 + tickz);
00377 tick22.push_back(r2 - tickz);
00378 G4Point3D textPosition(0., tickLength, 0.);
00379
00380
00381
00382 G4Transform3D transformation;
00383 if (scale.GetAutoPlacing()) {
00384 G4Transform3D rotation;
00385 switch (scale.GetDirection()) {
00386 case G4Scale::x:
00387 break;
00388 case G4Scale::y:
00389 rotation = G4RotateZ3D(piBy2);
00390 break;
00391 case G4Scale::z:
00392 rotation = G4RotateY3D(piBy2);
00393 break;
00394 }
00395 G4double sxmid(scale.GetXmid());
00396 G4double symid(scale.GetYmid());
00397 G4double szmid(scale.GetZmid());
00398 sxmid = xmin + oneMinusMargin * (xmax - xmin);
00399 symid = ymin + margin * (ymax - ymin);
00400 szmid = zmin + oneMinusMargin * (zmax - zmin);
00401 switch (scale.GetDirection()) {
00402 case G4Scale::x:
00403 sxmid -= halfLength;
00404 break;
00405 case G4Scale::y:
00406 symid += halfLength;
00407 break;
00408 case G4Scale::z:
00409 szmid -= halfLength;
00410 break;
00411 }
00412 G4Translate3D translation(sxmid, symid, szmid);
00413 transformation = translation * rotation;
00414 } else {
00415 if (fpModel) transformation = fpModel->GetTransformation();
00416 }
00417
00418
00419
00420
00421
00422 AddPrimitive(scaleLine.transform(transformation));
00423 AddPrimitive(tick11.transform(transformation));
00424 AddPrimitive(tick12.transform(transformation));
00425 AddPrimitive(tick21.transform(transformation));
00426 AddPrimitive(tick22.transform(transformation));
00427 G4Text text(scale.GetAnnotation(),textPosition.transform(transformation));
00428 text.SetScreenSize(12.);
00429 AddPrimitive(text);
00430 }
00431
00432 void G4VSceneHandler::AddPrimitive (const G4Polymarker& polymarker) {
00433 switch (polymarker.GetMarkerType()) {
00434 default:
00435 case G4Polymarker::dots:
00436 {
00437 for (size_t iPoint = 0; iPoint < polymarker.size (); iPoint++) {
00438 G4Circle dot (polymarker);
00439 dot.SetPosition (polymarker[iPoint]);
00440 dot.SetWorldSize (0.);
00441 dot.SetScreenSize (0.1);
00442 AddPrimitive (dot);
00443 }
00444 }
00445 break;
00446 case G4Polymarker::circles:
00447 {
00448 for (size_t iPoint = 0; iPoint < polymarker.size (); iPoint++) {
00449 G4Circle circle (polymarker);
00450 circle.SetPosition (polymarker[iPoint]);
00451 AddPrimitive (circle);
00452 }
00453 }
00454 break;
00455 case G4Polymarker::squares:
00456 {
00457 for (size_t iPoint = 0; iPoint < polymarker.size (); iPoint++) {
00458 G4Square square (polymarker);
00459 square.SetPosition (polymarker[iPoint]);
00460 AddPrimitive (square);
00461 }
00462 }
00463 break;
00464 }
00465 }
00466
00467 void G4VSceneHandler::RemoveViewerFromList (G4VViewer* pViewer) {
00468 fViewerList.remove(pViewer);
00469 }
00470
00471 void G4VSceneHandler::SetScene (G4Scene* pScene) {
00472 fpScene = pScene;
00473
00474 G4ViewerListIterator i;
00475 for (i = fViewerList.begin(); i != fViewerList.end(); i++) {
00476 (*i) -> SetNeedKernelVisit (true);
00477 }
00478 }
00479
00480 void G4VSceneHandler::RequestPrimitives (const G4VSolid& solid) {
00481 BeginPrimitives (fObjectTransformation);
00482 G4NURBS* pNURBS = 0;
00483 G4Polyhedron* pPolyhedron = 0;
00484 switch (fpViewer -> GetViewParameters () . GetRepStyle ()) {
00485 case G4ViewParameters::nurbs:
00486 pNURBS = solid.CreateNURBS ();
00487 if (pNURBS) {
00488 static G4bool warned = false;
00489 if (!warned) {
00490 warned = true;
00491 G4cout <<
00492 "!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!"
00493 "!!!!! NURBS are deprecated and will be removed in the next major release."
00494 "\n!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!"
00495 << G4endl;
00496 }
00497 pNURBS -> SetVisAttributes (fpVisAttribs);
00498 AddPrimitive (*pNURBS);
00499 delete pNURBS;
00500 break;
00501 }
00502 else {
00503 G4VisManager::Verbosity verbosity = G4VisManager::GetVerbosity();
00504 if (verbosity >= G4VisManager::errors) {
00505 G4cout <<
00506 "ERROR: G4VSceneHandler::RequestPrimitives"
00507 "\n NURBS not available for "
00508 << solid.GetName () << G4endl;
00509 G4cout << "Trying polyhedron." << G4endl;
00510 }
00511 }
00512
00513 case G4ViewParameters::polyhedron:
00514 default:
00515 G4Polyhedron::SetNumberOfRotationSteps (GetNoOfSides (fpVisAttribs));
00516 pPolyhedron = solid.GetPolyhedron ();
00517 G4Polyhedron::ResetNumberOfRotationSteps ();
00518 if (pPolyhedron) {
00519 pPolyhedron -> SetVisAttributes (fpVisAttribs);
00520 AddPrimitive (*pPolyhedron);
00521 }
00522 else {
00523 G4VisManager::Verbosity verbosity = G4VisManager::GetVerbosity();
00524 if (verbosity >= G4VisManager::errors) {
00525 G4cout <<
00526 "ERROR: G4VSceneHandler::RequestPrimitives"
00527 "\n Polyhedron not available for " << solid.GetName () <<
00528 ".\n This means it cannot be visualized on most systems."
00529 "\n Contact the Visualization Coordinator." << G4endl;
00530 }
00531 }
00532 break;
00533 }
00534 EndPrimitives ();
00535 }
00536
00537 void G4VSceneHandler::ProcessScene () {
00538
00539
00540
00541
00542 if (!fpScene) return;
00543
00544 G4VisManager* visManager = G4VisManager::GetInstance();
00545
00546 if (!visManager->GetConcreteInstance()) return;
00547
00548 G4VisManager::Verbosity verbosity = visManager->GetVerbosity();
00549
00550 fReadyForTransients = false;
00551
00552
00553
00554
00555 G4bool tmpMarkForClearingTransientStore = fMarkForClearingTransientStore;
00556 fMarkForClearingTransientStore = false;
00557
00558
00559
00560 const std::vector<G4Scene::Model>& runDurationModelList =
00561 fpScene -> GetRunDurationModelList ();
00562
00563 if (runDurationModelList.size ()) {
00564 if (verbosity >= G4VisManager::confirmations) {
00565 G4cout << "Traversing scene data..." << G4endl;
00566 }
00567
00568 BeginModeling ();
00569
00570
00571 G4ModelingParameters* pMP = CreateModelingParameters ();
00572
00573 for (size_t i = 0; i < runDurationModelList.size (); i++) {
00574 if (runDurationModelList[i].fActive) {
00575 G4VModel* pModel = runDurationModelList[i].fpModel;
00576
00577
00578
00579
00580 pModel -> SetModelingParameters (pMP);
00581 SetModel (pModel);
00582 pModel -> DescribeYourselfTo (*this);
00583 pModel -> SetModelingParameters (0);
00584 }
00585 }
00586
00587 delete pMP;
00588 EndModeling ();
00589 }
00590
00591 fReadyForTransients = true;
00592
00593
00594
00595 G4StateManager* stateManager = G4StateManager::GetStateManager();
00596 G4ApplicationState state = stateManager->GetCurrentState();
00597 if (state == G4State_Idle || state == G4State_GeomClosed) {
00598
00599 visManager->SetEventRefreshing(true);
00600
00601 if (visManager->GetRequestedEvent()) {
00602 DrawEvent(visManager->GetRequestedEvent());
00603
00604 } else {
00605
00606 G4RunManager* runManager = G4RunManager::GetRunManager();
00607 if (runManager) {
00608 const G4Run* run = runManager->GetCurrentRun();
00609 const std::vector<const G4Event*>* events =
00610 run? run->GetEventVector(): 0;
00611 size_t nKeptEvents = 0;
00612 if (events) nKeptEvents = events->size();
00613 if (nKeptEvents) {
00614
00615 if (fpScene->GetRefreshAtEndOfEvent()) {
00616
00617 if (verbosity >= G4VisManager::confirmations) {
00618 G4cout << "Refreshing event..." << G4endl;
00619 }
00620 const G4Event* event = 0;
00621 if (events && events->size()) event = events->back();
00622 if (event) DrawEvent(event);
00623
00624 } else {
00625
00626 if (verbosity >= G4VisManager::confirmations) {
00627 G4cout << "Refreshing events in run..." << G4endl;
00628 }
00629 for (size_t i = 0; i < nKeptEvents; ++i) {
00630 const G4Event* event = (*events)[i];
00631 if (event) DrawEvent(event);
00632 }
00633
00634 if (!fpScene->GetRefreshAtEndOfRun()) {
00635 if (verbosity >= G4VisManager::warnings) {
00636 G4cout <<
00637 "WARNING: Cannot refresh events accumulated over more"
00638 "\n than one runs. Refreshed just the last run."
00639 << G4endl;
00640 }
00641 }
00642 }
00643 }
00644 }
00645 }
00646 visManager->SetEventRefreshing(false);
00647 }
00648
00649
00650
00651 if (state == G4State_Idle || state == G4State_GeomClosed) {
00652 DrawEndOfRunModels();
00653 }
00654
00655 fMarkForClearingTransientStore = tmpMarkForClearingTransientStore;
00656 }
00657
00658 void G4VSceneHandler::DrawEvent(const G4Event* event)
00659 {
00660 const std::vector<G4Scene::Model>& EOEModelList =
00661 fpScene -> GetEndOfEventModelList ();
00662 size_t nModels = EOEModelList.size();
00663 if (nModels) {
00664 G4ModelingParameters* pMP = CreateModelingParameters();
00665 pMP->SetEvent(event);
00666 for (size_t i = 0; i < nModels; i++) {
00667 if (EOEModelList[i].fActive) {
00668 G4VModel* pModel = EOEModelList[i].fpModel;
00669 pModel -> SetModelingParameters(pMP);
00670 SetModel (pModel);
00671 pModel -> DescribeYourselfTo (*this);
00672 pModel -> SetModelingParameters(0);
00673 }
00674 }
00675 delete pMP;
00676 SetModel (0);
00677 }
00678 }
00679
00680 void G4VSceneHandler::DrawEndOfRunModels()
00681 {
00682 const std::vector<G4Scene::Model>& EORModelList =
00683 fpScene -> GetEndOfRunModelList ();
00684 size_t nModels = EORModelList.size();
00685 if (nModels) {
00686 G4ModelingParameters* pMP = CreateModelingParameters();
00687 pMP->SetEvent(0);
00688 for (size_t i = 0; i < nModels; i++) {
00689 if (EORModelList[i].fActive) {
00690 G4VModel* pModel = EORModelList[i].fpModel;
00691 pModel -> SetModelingParameters(pMP);
00692 SetModel (pModel);
00693 pModel -> DescribeYourselfTo (*this);
00694 pModel -> SetModelingParameters(0);
00695 }
00696 }
00697 delete pMP;
00698 SetModel (0);
00699 }
00700 }
00701
00702 G4ModelingParameters* G4VSceneHandler::CreateModelingParameters ()
00703 {
00704
00705 const G4ViewParameters& vp = fpViewer -> GetViewParameters ();
00706
00707
00708 G4ModelingParameters::DrawingStyle modelDrawingStyle =
00709 G4ModelingParameters::wf;
00710 switch (vp.GetDrawingStyle ()) {
00711 default:
00712 case G4ViewParameters::wireframe:
00713 modelDrawingStyle = G4ModelingParameters::wf;
00714 break;
00715 case G4ViewParameters::hlr:
00716 modelDrawingStyle = G4ModelingParameters::hlr;
00717 break;
00718 case G4ViewParameters::hsr:
00719 modelDrawingStyle = G4ModelingParameters::hsr;
00720 break;
00721 case G4ViewParameters::hlhsr:
00722 modelDrawingStyle = G4ModelingParameters::hlhsr;
00723 break;
00724 }
00725
00726
00727 G4bool reallyCullCovered =
00728 vp.IsCullingCovered()
00729 && !vp.IsSection ()
00730 && !vp.IsCutaway ()
00731 ;
00732
00733 G4ModelingParameters* pModelingParams = new G4ModelingParameters
00734 (vp.GetDefaultVisAttributes (),
00735 modelDrawingStyle,
00736 vp.IsCulling (),
00737 vp.IsCullingInvisible (),
00738 vp.IsDensityCulling (),
00739 vp.GetVisibleDensity (),
00740 reallyCullCovered,
00741 vp.GetNoOfSides ()
00742 );
00743
00744 pModelingParams->SetWarning
00745 (G4VisManager::GetVerbosity() >= G4VisManager::warnings);
00746
00747 pModelingParams->SetExplodeFactor(vp.GetExplodeFactor());
00748 pModelingParams->SetExplodeCentre(vp.GetExplodeCentre());
00749
00750 pModelingParams->SetSectionSolid(CreateSectionSolid());
00751 pModelingParams->SetCutawaySolid(CreateCutawaySolid());
00752
00753
00754 pModelingParams->SetVisAttributesModifiers(vp.GetVisAttributesModifiers());
00755
00756 return pModelingParams;
00757 }
00758
00759 G4VSolid* G4VSceneHandler::CreateSectionSolid()
00760 {
00761 G4VSolid* sectioner = 0;
00762 const G4ViewParameters& vp = fpViewer->GetViewParameters();
00763 if (vp.IsSection () ) {
00764 G4double radius = fpScene->GetExtent().GetExtentRadius();
00765 G4double safe = radius + fpScene->GetExtent().GetExtentCentre().mag();
00766 G4VSolid* sectionBox =
00767 new G4Box("_sectioner", safe, safe, 1.e-5 * radius);
00768 const G4Plane3D& sp = vp.GetSectionPlane ();
00769 G4double a = sp.a();
00770 G4double b = sp.b();
00771 G4double c = sp.c();
00772 G4double d = sp.d();
00773 G4Transform3D transform = G4TranslateZ3D(-d);
00774 const G4Normal3D normal(a,b,c);
00775 if (normal != G4Normal3D(0,0,1)) {
00776 const G4double angle = std::acos(normal.dot(G4Normal3D(0,0,1)));
00777 const G4Vector3D axis = G4Normal3D(0,0,1).cross(normal);
00778 transform = G4Rotate3D(angle, axis) * transform;
00779 }
00780 sectioner = new G4DisplacedSolid
00781 ("_displaced_sectioning_box", sectionBox, transform);
00782 }
00783 return sectioner;
00784 }
00785
00786 G4VSolid* G4VSceneHandler::CreateCutawaySolid()
00787 {
00788 return 0;
00789 }
00790
00791 void G4VSceneHandler::LoadAtts(const G4Visible& visible, G4AttHolder* holder)
00792 {
00793
00794 const G4VisAttributes* va = visible.GetVisAttributes();
00795 if (va) {
00796 const std::map<G4String,G4AttDef>* vaDefs =
00797 va->GetAttDefs();
00798 if (vaDefs) {
00799 holder->AddAtts(visible.GetVisAttributes()->CreateAttValues(), vaDefs);
00800 }
00801 }
00802
00803 G4PhysicalVolumeModel* pPVModel =
00804 dynamic_cast<G4PhysicalVolumeModel*>(fpModel);
00805 if (pPVModel) {
00806
00807 const std::map<G4String,G4AttDef>* pvDefs = pPVModel->GetAttDefs();
00808 if (pvDefs) {
00809 holder->AddAtts(pPVModel->CreateCurrentAttValues(), pvDefs);
00810 }
00811 }
00812
00813 G4TrajectoriesModel* trajModel = dynamic_cast<G4TrajectoriesModel*>(fpModel);
00814 if (trajModel) {
00815
00816 const std::map<G4String,G4AttDef>* trajModelDefs = trajModel->GetAttDefs();
00817 if (trajModelDefs) {
00818 holder->AddAtts(trajModel->CreateCurrentAttValues(), trajModelDefs);
00819 }
00820
00821 const G4VTrajectory* traj = trajModel->GetCurrentTrajectory();
00822 const std::map<G4String,G4AttDef>* trajDefs = traj->GetAttDefs();
00823 if (trajDefs) {
00824 holder->AddAtts(traj->CreateAttValues(), trajDefs);
00825 }
00826 G4int nPoints = traj->GetPointEntries();
00827 for (G4int i = 0; i < nPoints; ++i) {
00828 G4VTrajectoryPoint* trajPoint = traj->GetPoint(i);
00829 const std::map<G4String,G4AttDef>* pointDefs = trajPoint->GetAttDefs();
00830 if (pointDefs) {
00831 holder->AddAtts(trajPoint->CreateAttValues(), pointDefs);
00832 }
00833 }
00834 }
00835
00836 G4HitsModel* hitsModel = dynamic_cast<G4HitsModel*>(fpModel);
00837 if (hitsModel) {
00838
00839 const G4VHit* hit = hitsModel->GetCurrentHit();
00840 const std::map<G4String,G4AttDef>* hitsDefs = hit->GetAttDefs();
00841 if (hitsDefs) {
00842 holder->AddAtts(hit->CreateAttValues(), hitsDefs);
00843 }
00844 }
00845 }
00846
00847 const G4Colour& G4VSceneHandler::GetColour (const G4Visible& visible) {
00848
00849 const G4Colour& colour = fpViewer ->
00850 GetApplicableVisAttributes (visible.GetVisAttributes ()) -> GetColour ();
00851 return colour;
00852 }
00853
00854 const G4Colour& G4VSceneHandler::GetTextColour (const G4Text& text) {
00855 const G4VisAttributes* pVA = text.GetVisAttributes ();
00856 if (!pVA) {
00857 pVA = fpViewer -> GetViewParameters (). GetDefaultTextVisAttributes ();
00858 }
00859 const G4Colour& colour = pVA -> GetColour ();
00860 return colour;
00861 }
00862
00863 G4double G4VSceneHandler::GetLineWidth(const G4VisAttributes* pVisAttribs)
00864 {
00865 G4double lineWidth = pVisAttribs->GetLineWidth();
00866 if (lineWidth < 1.) lineWidth = 1.;
00867 lineWidth *= fpViewer -> GetViewParameters().GetGlobalLineWidthScale();
00868 if (lineWidth < 1.) lineWidth = 1.;
00869 return lineWidth;
00870 }
00871
00872 G4ViewParameters::DrawingStyle G4VSceneHandler::GetDrawingStyle
00873 (const G4VisAttributes* pVisAttribs) {
00874
00875
00876
00877 G4ViewParameters::DrawingStyle style =
00878 fpViewer->GetViewParameters().GetDrawingStyle();
00879 if (pVisAttribs -> IsForceDrawingStyle ()) {
00880 G4VisAttributes::ForcedDrawingStyle forcedStyle =
00881 pVisAttribs -> GetForcedDrawingStyle ();
00882
00883
00884 switch (forcedStyle) {
00885 case (G4VisAttributes::solid):
00886 switch (style) {
00887 case (G4ViewParameters::hlr):
00888 style = G4ViewParameters::hlhsr;
00889 break;
00890 case (G4ViewParameters::wireframe):
00891 style = G4ViewParameters::hsr;
00892 break;
00893 case (G4ViewParameters::hlhsr):
00894 case (G4ViewParameters::hsr):
00895 default:
00896 break;
00897 }
00898 break;
00899 case (G4VisAttributes::wireframe):
00900 default:
00901
00902
00903
00904
00905 style = G4ViewParameters::wireframe;
00906 break;
00907 }
00908 }
00909 return style;
00910 }
00911
00912 G4bool G4VSceneHandler::GetAuxEdgeVisible (const G4VisAttributes* pVisAttribs) {
00913 G4bool isAuxEdgeVisible = fpViewer->GetViewParameters().IsAuxEdgeVisible ();
00914 if (pVisAttribs -> IsForceAuxEdgeVisible()) isAuxEdgeVisible = true;
00915 return isAuxEdgeVisible;
00916 }
00917
00918 G4double G4VSceneHandler::GetMarkerSize
00919 (const G4VMarker& marker,
00920 G4VSceneHandler::MarkerSizeType& markerSizeType)
00921 {
00922 G4bool userSpecified = marker.GetWorldSize() || marker.GetScreenSize();
00923 const G4VMarker& defaultMarker =
00924 fpViewer -> GetViewParameters().GetDefaultMarker();
00925 G4double size = userSpecified ?
00926 marker.GetWorldSize() : defaultMarker.GetWorldSize();
00927 if (size) {
00928
00929 markerSizeType = world;
00930 }
00931 else {
00932 size = userSpecified ?
00933 marker.GetScreenSize() : defaultMarker.GetScreenSize();
00934
00935 markerSizeType = screen;
00936 }
00937 size *= fpViewer -> GetViewParameters().GetGlobalMarkerScale();
00938 if (markerSizeType == screen && size < 1.) size = 1.;
00939 return size;
00940 }
00941
00942 G4int G4VSceneHandler::GetNoOfSides(const G4VisAttributes* pVisAttribs)
00943 {
00944
00945
00946
00947 G4int lineSegmentsPerCircle = fpViewer->GetViewParameters().GetNoOfSides();
00948 if (pVisAttribs) {
00949 if (pVisAttribs->IsForceLineSegmentsPerCircle())
00950 lineSegmentsPerCircle = pVisAttribs->GetForcedLineSegmentsPerCircle();
00951 const G4int nSegmentsMin = 12;
00952 if (lineSegmentsPerCircle < nSegmentsMin) {
00953 lineSegmentsPerCircle = nSegmentsMin;
00954 G4cout <<
00955 "G4VSceneHandler::GetNoOfSides: attempt to set the"
00956 "\nnumber of line segements per circle < " << nSegmentsMin
00957 << "; forced to " << lineSegmentsPerCircle << G4endl;
00958 }
00959 }
00960 return lineSegmentsPerCircle;
00961 }
00962
00963 std::ostream& operator << (std::ostream& os, const G4VSceneHandler& sh) {
00964
00965 os << "Scene handler " << sh.fName << " has "
00966 << sh.fViewerList.size () << " viewer(s):";
00967 for (size_t i = 0; i < sh.fViewerList.size (); i++) {
00968 os << "\n " << *(sh.fViewerList [i]);
00969 }
00970
00971 if (sh.fpScene) {
00972 os << "\n " << *sh.fpScene;
00973 }
00974 else {
00975 os << "\n This scene handler currently has no scene.";
00976 }
00977
00978 return os;
00979 }