Geant4  10.01
G4VisCommandsSceneAdd.cc
Go to the documentation of this file.
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: G4VisCommandsSceneAdd.cc 87360 2014-12-01 16:07:16Z gcosmo $
28 // /vis/scene/add commands - John Allison 9th August 1998
29 
30 #include "G4VisCommandsSceneAdd.hh"
31 
33 #include "G4LogicalVolumeStore.hh"
34 #include "G4PhysicalVolumeModel.hh"
35 #include "G4LogicalVolumeModel.hh"
36 #include "G4ModelingParameters.hh"
37 #include "G4HitsModel.hh"
38 #include "G4DigiModel.hh"
39 #include "G4MagneticFieldModel.hh"
40 #include "G4PSHitsModel.hh"
41 #include "G4TrajectoriesModel.hh"
42 #include "G4ScaleModel.hh"
43 #include "G4TextModel.hh"
44 #include "G4ArrowModel.hh"
45 #include "G4AxesModel.hh"
47 #include "G4ParticleTable.hh"
48 #include "G4ParticleDefinition.hh"
49 #include "G4ApplicationState.hh"
50 #include "G4VUserVisAction.hh"
51 #include "G4CallbackModel.hh"
52 #include "G4UnionSolid.hh"
53 #include "G4SubtractionSolid.hh"
54 #include "G4Polyhedron.hh"
55 #include "G4UImanager.hh"
56 #include "G4UIcommand.hh"
57 #include "G4UIcmdWithAString.hh"
59 #include "G4UIcmdWithAnInteger.hh"
60 #include "G4Tokenizer.hh"
61 #include "G4RunManager.hh"
62 #ifdef G4MULTITHREADED
63 #include "G4MTRunManager.hh"
64 #endif
65 #include "G4StateManager.hh"
66 #include "G4Run.hh"
67 #include "G4Event.hh"
68 #include "G4Trajectory.hh"
69 #include "G4TrajectoryPoint.hh"
70 #include "G4RichTrajectory.hh"
71 #include "G4RichTrajectoryPoint.hh"
72 #include "G4SmoothTrajectory.hh"
74 #include "G4AttDef.hh"
75 #include "G4AttCheck.hh"
76 #include "G4Polyline.hh"
77 #include "G4UnitsTable.hh"
78 #include "G4PhysicalConstants.hh"
79 #include "G4SystemOfUnits.hh"
80 
81 #include <sstream>
82 
83 // Local function with some frequently used error printing...
86  if (verbosity >= G4VisManager::warnings) {
87  G4cout <<
88  "WARNING: For some reason, possibly mentioned above, it has not been"
89  "\n possible to add to the scene."
90  << G4endl;
91  }
92 }
93 
95 
97  fpCommand = new G4UIcommand("/vis/scene/add/arrow", this);
98  fpCommand -> SetGuidance ("Adds arrow to current scene.");
99  G4bool omitable;
100  G4UIparameter* parameter;
101  parameter = new G4UIparameter ("x1", 'd', omitable = false);
102  fpCommand -> SetParameter (parameter);
103  parameter = new G4UIparameter ("y1", 'd', omitable = false);
104  fpCommand -> SetParameter (parameter);
105  parameter = new G4UIparameter ("z1", 'd', omitable = false);
106  fpCommand -> SetParameter (parameter);
107  parameter = new G4UIparameter ("x2", 'd', omitable = false);
108  fpCommand -> SetParameter (parameter);
109  parameter = new G4UIparameter ("y2", 'd', omitable = false);
110  fpCommand -> SetParameter (parameter);
111  parameter = new G4UIparameter ("z2", 'd', omitable = false);
112  fpCommand -> SetParameter (parameter);
113  parameter = new G4UIparameter ("unit", 's', omitable = true);
114  parameter->SetDefaultValue ("m");
115  fpCommand->SetParameter (parameter);
116 }
117 
119  delete fpCommand;
120 }
121 
123  return "";
124 }
125 
127 {
129  G4bool warn(verbosity >= G4VisManager::warnings);
130 
131  G4Scene* pScene = fpVisManager->GetCurrentScene();
132  if (!pScene) {
133  if (verbosity >= G4VisManager::errors) {
134  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
135  }
136  return;
137  }
138 
139  G4String unitString;
140  G4double x1, y1, z1, x2, y2, z2;
141  std::istringstream is(newValue);
142  is >> x1 >> y1 >> z1 >> x2 >> y2 >> z2 >> unitString;
143  G4double unit = G4UIcommand::ValueOf(unitString);
144  x1 *= unit; y1 *= unit; z1 *= unit;
145  x2 *= unit; y2 *= unit; z2 *= unit;
146 
147  // Consult scene for arrow width.
148  const G4VisExtent& sceneExtent = pScene->GetExtent();
149  G4double arrowWidth =
150  0.005 * fCurrentLineWidth * sceneExtent.GetExtentRadius();
151 
152  G4VModel* model = new G4ArrowModel
153  (x1, y1, z1, x2, y2, z2,
154  arrowWidth, fCurrentColour, newValue);
155 
156  const G4String& currentSceneName = pScene -> GetName ();
157  G4bool successful = pScene -> AddRunDurationModel (model, warn);
158  if (successful) {
159  if (verbosity >= G4VisManager::confirmations) {
160  G4cout << "Arrow has been added to scene \""
161  << currentSceneName << "\"."
162  << G4endl;
163  }
164  }
165  else G4VisCommandsSceneAddUnsuccessful(verbosity);
166  UpdateVisManagerScene (currentSceneName);
167 }
168 
170 
172  fpCommand = new G4UIcommand("/vis/scene/add/arrow2D", this);
173  fpCommand -> SetGuidance ("Adds 2D arrow to current scene.");
174  G4bool omitable;
175  G4UIparameter* parameter;
176  parameter = new G4UIparameter ("x1", 'd', omitable = false);
177  fpCommand -> SetParameter (parameter);
178  parameter = new G4UIparameter ("y1", 'd', omitable = false);
179  fpCommand -> SetParameter (parameter);
180  parameter = new G4UIparameter ("x2", 'd', omitable = false);
181  fpCommand -> SetParameter (parameter);
182  parameter = new G4UIparameter ("y2", 'd', omitable = false);
183  fpCommand -> SetParameter (parameter);
184 }
185 
187  delete fpCommand;
188 }
189 
191  return "";
192 }
193 
195 {
197  G4bool warn(verbosity >= G4VisManager::warnings);
198 
199  G4Scene* pScene = fpVisManager->GetCurrentScene();
200  if (!pScene) {
201  if (verbosity >= G4VisManager::errors) {
202  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
203  }
204  return;
205  }
206 
207  G4double x1, y1, x2, y2;
208  std::istringstream is(newValue);
209  is >> x1 >> y1 >> x2 >> y2;
210 
211  Arrow2D* arrow2D = new Arrow2D
212  (x1, y1, x2, y2, fCurrentLineWidth, fCurrentColour);
213  G4VModel* model =
215  model->SetType("Arrow2D");
216  model->SetGlobalTag("Arrow2D");
217  model->SetGlobalDescription("Arrow2D: " + newValue);
218  const G4String& currentSceneName = pScene -> GetName ();
219  G4bool successful = pScene -> AddRunDurationModel (model, warn);
220  if (successful) {
221  if (verbosity >= G4VisManager::confirmations) {
222  G4cout << "A 2D arrow has been added to scene \""
223  << currentSceneName << "\"."
224  << G4endl;
225  }
226  }
227  else G4VisCommandsSceneAddUnsuccessful(verbosity);
228  UpdateVisManagerScene (currentSceneName);
229 }
230 
233  G4double x2, G4double y2,
234  G4double width, const G4Colour& colour):
235  fWidth(width), fColour(colour)
236 {
237  fShaftPolyline.push_back(G4Point3D(x1,y1,0));
238  fShaftPolyline.push_back(G4Point3D(x2,y2,0));
239  G4Vector3D arrowDirection = G4Vector3D(x2-x1,y2-y1,0).unit();
240  G4Vector3D arrowPointLeftDirection(arrowDirection);
241  arrowPointLeftDirection.rotateZ(150.*deg);
242  G4Vector3D arrowPointRightDirection(arrowDirection);
243  arrowPointRightDirection.rotateZ(-150.*deg);
244  fHeadPolyline.push_back(G4Point3D(x2,y2,0)+0.04*arrowPointLeftDirection);
245  fHeadPolyline.push_back(G4Point3D(x2,y2,0));
246  fHeadPolyline.push_back(G4Point3D(x2,y2,0)+0.04*arrowPointRightDirection);
247  G4VisAttributes va;
248  va.SetLineWidth(fWidth);
249  va.SetColour(fColour);
250  fShaftPolyline.SetVisAttributes(va);
251  fHeadPolyline.SetVisAttributes(va);
252 }
253 
254 void G4VisCommandSceneAddArrow2D::Arrow2D::operator()
255  (G4VGraphicsScene& sceneHandler, const G4Transform3D&)
256 {
257  sceneHandler.BeginPrimitives2D();
258  sceneHandler.AddPrimitive(fShaftPolyline);
259  sceneHandler.AddPrimitive(fHeadPolyline);
260  sceneHandler.EndPrimitives2D();
261 }
262 
264 
266  G4bool omitable;
267  fpCommand = new G4UIcommand ("/vis/scene/add/axes", this);
268  fpCommand -> SetGuidance ("Add axes.");
269  fpCommand -> SetGuidance
270  ("Draws axes at (x0, y0, z0) of given length and colour.");
271  fpCommand -> SetGuidance
272  ("If \"unitcolour\" is \"auto\", x, y and z will be red, green and blue"
273  "\n respectively. Otherwise choose from the pre-defined text-specified"
274  "\n colours - see information printed by the vis manager at start-up or"
275  "\n use \"/vis/list\".");
276  fpCommand -> SetGuidance
277  ("If \"length\" is negative, it is set to about 25% of scene extent.");
278  fpCommand -> SetGuidance
279  ("If \"showtext\" is false, annotations are suppressed.");
280  G4UIparameter* parameter;
281  parameter = new G4UIparameter ("x0", 'd', omitable = true);
282  parameter->SetDefaultValue (0.);
283  fpCommand->SetParameter (parameter);
284  parameter = new G4UIparameter ("y0", 'd', omitable = true);
285  parameter->SetDefaultValue (0.);
286  fpCommand->SetParameter (parameter);
287  parameter = new G4UIparameter ("z0", 'd', omitable = true);
288  parameter->SetDefaultValue (0.);
289  fpCommand->SetParameter (parameter);
290  parameter = new G4UIparameter ("length", 'd', omitable = true);
291  parameter->SetDefaultValue (-1.);
292  fpCommand->SetParameter (parameter);
293  parameter = new G4UIparameter ("unit", 's', omitable = true);
294  parameter->SetDefaultValue ("m");
295  fpCommand->SetParameter (parameter);
296  parameter = new G4UIparameter ("unitcolour", 's', omitable = true);
297  parameter->SetDefaultValue ("auto");
298  fpCommand->SetParameter (parameter);
299  parameter = new G4UIparameter ("showtext", 'b', omitable = true);
300  parameter->SetDefaultValue ("true");
301  fpCommand->SetParameter (parameter);
302 }
303 
305  delete fpCommand;
306 }
307 
309  return "";
310 }
311 
313 
315  G4bool warn(verbosity >= G4VisManager::warnings);
316 
317  G4Scene* pScene = fpVisManager->GetCurrentScene();
318  if (!pScene) {
319  if (verbosity >= G4VisManager::errors) {
320  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
321  }
322  return;
323  } else {
324  if (pScene->GetExtent().GetExtentRadius() <= 0.) {
325  if (verbosity >= G4VisManager::errors) {
326  G4cerr
327  << "ERROR: Scene has no extent. Add volumes or use \"/vis/scene/add/extent\"."
328  << G4endl;
329  }
330  return;
331  }
332  }
333 
334  G4String unitString, colourString, showTextString;
335  G4double x0, y0, z0, length;
336  std::istringstream is (newValue);
337  is >> x0 >> y0 >> z0 >> length >> unitString
338  >> colourString >> showTextString;
339  G4bool showText = G4UIcommand::ConvertToBool(showTextString);
340 
341 
342  G4double unit = G4UIcommand::ValueOf(unitString);
343  x0 *= unit; y0 *= unit; z0 *= unit;
344  const G4VisExtent& sceneExtent = pScene->GetExtent(); // Existing extent.
345  if (length < 0.) {
346  const G4double lengthMax = 0.5 * sceneExtent.GetExtentRadius();
347  const G4double intLog10Length = std::floor(std::log10(lengthMax));
348  length = std::pow(10,intLog10Length);
349  if (5.*length < lengthMax) length *= 5.;
350  else if (2.*length < lengthMax) length *= 2.;
351  } else {
352  length *= unit;
353  }
354 
355  // Consult scene for arrow width...
356  G4double arrowWidth =
357  0.005 * fCurrentLineWidth * sceneExtent.GetExtentRadius();
358  // ...but limit it to length/50.
359  if (arrowWidth > length/50.) arrowWidth = length/50.;
360 
361  G4VModel* model = new G4AxesModel
362  (x0, y0, z0, length, arrowWidth, colourString, newValue,
363  showText, fCurrentTextSize);
364 
365  G4bool successful = pScene -> AddRunDurationModel (model, warn);
366  const G4String& currentSceneName = pScene -> GetName ();
367  if (successful) {
368  if (verbosity >= G4VisManager::confirmations) {
369  G4cout << "Axes of length " << G4BestUnit(length,"Length")
370  << "have been added to scene \"" << currentSceneName << "\"."
371  << G4endl;
372  }
373  }
374  else G4VisCommandsSceneAddUnsuccessful(verbosity);
375  UpdateVisManagerScene (currentSceneName);
376 }
377 
379 
381  G4bool omitable;
382  fpCommand = new G4UIcommand ("/vis/scene/add/date", this);
383  fpCommand -> SetGuidance ("Adds date to current scene.");
384  fpCommand -> SetGuidance
385  ("If \"date\"is omitted, the current date and time is drawn."
386  "\nOtherwise, the string, including the rest of the line, is drawn.");
387  G4UIparameter* parameter;
388  parameter = new G4UIparameter ("size", 'i', omitable = true);
389  parameter -> SetGuidance ("Screen size of text in pixels.");
390  parameter -> SetDefaultValue (18);
391  fpCommand -> SetParameter (parameter);
392  parameter = new G4UIparameter ("x-position", 'd', omitable = true);
393  parameter -> SetGuidance ("x screen position in range -1 < x < 1.");
394  parameter -> SetDefaultValue (0.95);
395  fpCommand -> SetParameter (parameter);
396  parameter = new G4UIparameter ("y-position", 'd', omitable = true);
397  parameter -> SetGuidance ("y screen position in range -1 < y < 1.");
398  parameter -> SetDefaultValue (0.9);
399  fpCommand -> SetParameter (parameter);
400  parameter = new G4UIparameter ("layout", 's', omitable = true);
401  parameter -> SetGuidance ("Layout, i.e., adjustment: left|centre|right.");
402  parameter -> SetDefaultValue ("right");
403  fpCommand -> SetParameter (parameter);
404  parameter = new G4UIparameter ("date", 's', omitable = true);
405  parameter -> SetDefaultValue ("-");
406  fpCommand -> SetParameter (parameter);
407 }
408 
410  delete fpCommand;
411 }
412 
414  return "";
415 }
416 
418 {
420  G4bool warn(verbosity >= G4VisManager::warnings);
421 
422  G4Scene* pScene = fpVisManager->GetCurrentScene();
423  if (!pScene) {
424  if (verbosity >= G4VisManager::errors) {
425  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
426  }
427  return;
428  }
429 
430  G4int size;
431  G4double x, y;
432  G4String layoutString, dateString;
433  std::istringstream is(newValue);
434  is >> size >> x >> y >> layoutString >> dateString;
435  // Read rest of line, if any.
436  const size_t NREMAINDER = 100;
437  char remainder[NREMAINDER];
438  remainder[0]='\0'; // In case there is nothing remaining.
439  is.getline(remainder, NREMAINDER);
440  dateString += remainder;
441  G4Text::Layout layout = G4Text::right;
442  if (layoutString(0) == 'l') layout = G4Text::left;
443  else if (layoutString(0) == 'c') layout = G4Text::centre;
444  else if (layoutString(0) == 'r') layout = G4Text::right;
445 
446  Date* date = new Date(fpVisManager, size, x, y, layout, dateString);
447  G4VModel* model =
449  model->SetType("Date");
450  model->SetGlobalTag("Date");
451  model->SetGlobalDescription("Date");
452  const G4String& currentSceneName = pScene -> GetName ();
453  G4bool successful = pScene -> AddRunDurationModel (model, warn);
454  if (successful) {
455  if (verbosity >= G4VisManager::confirmations) {
456  G4cout << "Date has been added to scene \""
457  << currentSceneName << "\"."
458  << G4endl;
459  }
460  }
461  else G4VisCommandsSceneAddUnsuccessful(verbosity);
462  UpdateVisManagerScene (currentSceneName);
463 }
464 
465 void G4VisCommandSceneAddDate::Date::operator()
466  (G4VGraphicsScene& sceneHandler, const G4Transform3D&)
467 {
468  G4String time;
469  if (fDate == "-") {
470  time = fTimer.GetClockTime();
471  } else {
472  time = fDate;
473  }
474  // Check for \n, starting from back, and erase.
475  std::string::size_type i = time.rfind('\n');
476  if (i != std::string::npos) time.erase(i);
477  G4Text text(time, G4Point3D(fX, fY, 0.));
478  text.SetScreenSize(fSize);
479  text.SetLayout(fLayout);
480  G4VisAttributes textAtts(G4Colour(0.,1.,1));
481  text.SetVisAttributes(textAtts);
482  sceneHandler.BeginPrimitives2D();
483  sceneHandler.AddPrimitive(text);
484  sceneHandler.EndPrimitives2D();
485 }
486 
488 
490  fpCommand = new G4UIcmdWithoutParameter ("/vis/scene/add/digis", this);
491  fpCommand -> SetGuidance ("Adds digis to current scene.");
492  fpCommand -> SetGuidance
493  ("Digis are drawn at end of event when the scene in which"
494  "\nthey are added is current.");
495 }
496 
498  delete fpCommand;
499 }
500 
502  return "";
503 }
504 
506 
508  G4bool warn(verbosity >= G4VisManager::warnings);
509 
510  G4Scene* pScene = fpVisManager->GetCurrentScene();
511  if (!pScene) {
512  if (verbosity >= G4VisManager::errors) {
513  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
514  }
515  return;
516  }
517 
518  G4DigiModel* model = new G4DigiModel;
519  const G4String& currentSceneName = pScene -> GetName ();
520  G4bool successful = pScene -> AddEndOfEventModel (model, warn);
521  if (successful) {
522  if (verbosity >= G4VisManager::confirmations) {
523  G4cout << "Digis, if any, will be drawn at end of run in scene \""
524  << currentSceneName << "\"."
525  << G4endl;
526  }
527  }
528  else G4VisCommandsSceneAddUnsuccessful(verbosity);
529  UpdateVisManagerScene (currentSceneName);
530 }
531 
533 
535  G4bool omitable;
536  fpCommand = new G4UIcommand ("/vis/scene/add/eventID", this);
537  fpCommand -> SetGuidance ("Adds eventID to current scene.");
538  fpCommand -> SetGuidance
539  ("Run and event numbers are drawn at end of event or run when"
540  "\n the scene in which they are added is current.");
541  G4UIparameter* parameter;
542  parameter = new G4UIparameter ("size", 'i', omitable = true);
543  parameter -> SetGuidance ("Screen size of text in pixels.");
544  parameter -> SetDefaultValue (18);
545  fpCommand -> SetParameter (parameter);
546  parameter = new G4UIparameter ("x-position", 'd', omitable = true);
547  parameter -> SetGuidance ("x screen position in range -1 < x < 1.");
548  parameter -> SetDefaultValue (-0.95);
549  fpCommand -> SetParameter (parameter);
550  parameter = new G4UIparameter ("y-position", 'd', omitable = true);
551  parameter -> SetGuidance ("y screen position in range -1 < y < 1.");
552  parameter -> SetDefaultValue (0.9);
553  fpCommand -> SetParameter (parameter);
554  parameter = new G4UIparameter ("layout", 's', omitable = true);
555  parameter -> SetGuidance ("Layout, i.e., adjustment: left|centre|right.");
556  parameter -> SetDefaultValue ("left");
557  fpCommand -> SetParameter (parameter);
558 }
559 
561  delete fpCommand;
562 }
563 
565  return "";
566 }
567 
569 {
571  G4bool warn(verbosity >= G4VisManager::warnings);
572 
573  G4Scene* pScene = fpVisManager->GetCurrentScene();
574  if (!pScene) {
575  if (verbosity >= G4VisManager::errors) {
576  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
577  }
578  return;
579  }
580 
581  G4int size;
582  G4double x, y;
583  G4String layoutString;
584  std::istringstream is(newValue);
585  is >> size >> x >> y >> layoutString;
586 
587  G4Text::Layout layout = G4Text::right;
588  if (layoutString(0) == 'l') layout = G4Text::left;
589  else if (layoutString(0) == 'c') layout = G4Text::centre;
590  else if (layoutString(0) == 'r') layout = G4Text::right;
591 
592  EventID* eventID = new EventID(fpVisManager, size, x, y, layout);
593  G4VModel* model =
595  model->SetType("EventID");
596  model->SetGlobalTag("EventID");
597  model->SetGlobalDescription("EventID");
598  const G4String& currentSceneName = pScene -> GetName ();
599  G4bool successful = pScene -> AddEndOfEventModel (model, warn);
600  if (successful) {
601  if (verbosity >= G4VisManager::confirmations) {
602  G4cout << "EventID has been added to scene \""
603  << currentSceneName << "\"."
604  << G4endl;
605  }
606  }
607  else G4VisCommandsSceneAddUnsuccessful(verbosity);
608  UpdateVisManagerScene (currentSceneName);
609 }
610 
611 void G4VisCommandSceneAddEventID::EventID::operator()
612  (G4VGraphicsScene& sceneHandler, const G4Transform3D&)
613 {
614  const G4Run* currentRun = 0;
616 #ifdef G4MULTITHREADED
619  }
620 #endif
621  if (runManager) currentRun = runManager->GetCurrentRun();
622 
623  G4VModel* model = fpVisManager->GetCurrentSceneHandler()->GetModel();
624  const G4ModelingParameters* mp = 0;
625  const G4Event* currentEvent = 0;
626  if (model) {
627  mp = model->GetModelingParameters();
628  currentEvent = mp->GetEvent();
629  } else {
630  G4VisManager::Verbosity verbosity = fpVisManager->GetVerbosity();
631  if (verbosity >= G4VisManager::errors) {
632  G4cerr << "ERROR: No model defined for this SceneHandler : "
633  << fpVisManager->GetCurrentSceneHandler()->GetName()
634  << G4endl;
635  }
636  }
637  if (currentRun && currentEvent) {
638  G4int runID = currentRun->GetRunID();
639  G4int eventID = currentEvent->GetEventID();
640  std::ostringstream oss;
641  if (fpVisManager->GetCurrentScene()->GetRefreshAtEndOfEvent()) {
642  oss << "Run " << runID << " Event " << eventID;
643  } else {
644  G4int nEvents = 0;
646  G4ApplicationState state = stateManager->GetCurrentState();
647  if (state == G4State_EventProc) {
648  nEvents = currentRun->GetNumberOfEventToBeProcessed();
649  } else {
650  const std::vector<const G4Event*>* events =
651  currentRun->GetEventVector();
652  if (events) nEvents = events->size();
653  }
654  if (eventID < nEvents - 1) return; // Not last event.
655  else {
656  oss << "Run " << runID << " (" << nEvents << " event";
657  if (nEvents != 1) oss << 's';
658  oss << ')';
659  }
660  }
661  G4Text text(oss.str(), G4Point3D(fX, fY, 0.));
662  text.SetScreenSize(fSize);
663  text.SetLayout(fLayout);
664  G4VisAttributes textAtts(G4Colour(0.,1.,1));
665  text.SetVisAttributes(textAtts);
666  sceneHandler.BeginPrimitives2D();
667  sceneHandler.AddPrimitive(text);
668  sceneHandler.EndPrimitives2D();
669  }
670 }
671 
673 
675  fpCommand = new G4UIcommand("/vis/scene/add/extent", this);
676  fpCommand -> SetGuidance
677  ("Adds a dummy model with given extent to the current scene."
678  "\nThis can be used to provide an extent to the scene even if"
679  "\nno other models with extent are available. For example,"
680  "\neven if there is no geometry.");
681  G4bool omitable;
682  G4UIparameter* parameter;
683  parameter = new G4UIparameter ("xmin", 'd', omitable = true);
684  parameter -> SetDefaultValue (0.);
685  fpCommand -> SetParameter (parameter);
686  parameter = new G4UIparameter ("xmax", 'd', omitable = true);
687  parameter -> SetDefaultValue (0.);
688  fpCommand -> SetParameter (parameter);
689  parameter = new G4UIparameter ("ymin", 'd', omitable = true);
690  parameter -> SetDefaultValue (0.);
691  fpCommand -> SetParameter (parameter);
692  parameter = new G4UIparameter ("ymax", 'd', omitable = true);
693  parameter -> SetDefaultValue (0.);
694  fpCommand -> SetParameter (parameter);
695  parameter = new G4UIparameter ("zmin", 'd', omitable = true);
696  parameter -> SetDefaultValue (0.);
697  fpCommand -> SetParameter (parameter);
698  parameter = new G4UIparameter ("zmax", 'd', omitable = true);
699  parameter -> SetDefaultValue (0.);
700  fpCommand -> SetParameter (parameter);
701  parameter = new G4UIparameter ("unit", 's', omitable = true);
702  parameter -> SetDefaultValue ("m");
703  fpCommand -> SetParameter (parameter);
704 }
705 
707  delete fpCommand;
708 }
709 
711  return "";
712 }
713 
715 {
717  G4bool warn(verbosity >= G4VisManager::warnings);
718 
719  G4Scene* pScene = fpVisManager->GetCurrentScene();
720  if (!pScene) {
721  if (verbosity >= G4VisManager::errors) {
722  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
723  }
724  return;
725  }
726 
727  G4double xmin, xmax, ymin, ymax, zmin, zmax;
728  G4String unitString;
729  std::istringstream is(newValue);
730  is >> xmin >> xmax >> ymin >> ymax >> zmin >> zmax >> unitString;
731  G4double unit = G4UIcommand::ValueOf(unitString);
732  xmin *= unit; xmax *= unit;
733  ymin *= unit; ymax *= unit;
734  zmin *= unit; zmax *= unit;
735 
736  G4VisExtent visExtent(xmin, xmax, ymin, ymax, zmin, zmax);
737  Extent* extent = new Extent(xmin, xmax, ymin, ymax, zmin, zmax);
738  G4VModel* model =
740  model->SetType("Extent");
741  model->SetGlobalTag("Extent");
742  model->SetGlobalDescription("Extent: " + newValue);
743  model->SetExtent(visExtent);
744  const G4String& currentSceneName = pScene -> GetName ();
745  G4bool successful = pScene -> AddRunDurationModel (model, warn);
746  if (successful) {
747  if (verbosity >= G4VisManager::confirmations) {
748  G4cout << "A benign model with extent "
749  << visExtent
750  << " has been added to scene \""
751  << currentSceneName << "\"."
752  << G4endl;
753  }
754  }
755  else G4VisCommandsSceneAddUnsuccessful(verbosity);
756  UpdateVisManagerScene (currentSceneName);
757 }
758 
760 (G4double xmin, G4double xmax,
761  G4double ymin, G4double ymax,
762  G4double zmin, G4double zmax):
763 fExtent(xmin,xmax,ymin,ymax,zmin,zmax)
764 {}
765 
766 void G4VisCommandSceneAddExtent::Extent::operator()
768 {}
769 
771 
773  fpCommand = new G4UIcommand("/vis/scene/add/frame", this);
774  fpCommand -> SetGuidance ("Adds frame to current scene.");
775  G4bool omitable;
776  G4UIparameter* parameter;
777  parameter = new G4UIparameter ("size", 'd', omitable = true);
778  parameter -> SetGuidance ("Size of frame. 1 = full window.");
779  parameter -> SetParameterRange ("size > 0 && size <=1");
780  parameter -> SetDefaultValue (0.97);
781  fpCommand -> SetParameter (parameter);
782 }
783 
785  delete fpCommand;
786 }
787 
789  return "";
790 }
791 
793 {
795  G4bool warn(verbosity >= G4VisManager::warnings);
796 
797  G4Scene* pScene = fpVisManager->GetCurrentScene();
798  if (!pScene) {
799  if (verbosity >= G4VisManager::errors) {
800  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
801  }
802  return;
803  }
804 
805  G4double size;
806  std::istringstream is(newValue);
807  is >> size;
808 
809  Frame* frame = new Frame(size, fCurrentLineWidth, fCurrentColour);
810  G4VModel* model =
812  model->SetType("Frame");
813  model->SetGlobalTag("Frame");
814  model->SetGlobalDescription("Frame: " + newValue);
815  const G4String& currentSceneName = pScene -> GetName ();
816  G4bool successful = pScene -> AddRunDurationModel (model, warn);
817  if (successful) {
818  if (verbosity >= G4VisManager::confirmations) {
819  G4cout << "Frame has been added to scene \""
820  << currentSceneName << "\"."
821  << G4endl;
822  }
823  }
824  else G4VisCommandsSceneAddUnsuccessful(verbosity);
825  UpdateVisManagerScene (currentSceneName);
826 }
827 
828 void G4VisCommandSceneAddFrame::Frame::operator()
829  (G4VGraphicsScene& sceneHandler, const G4Transform3D&)
830 {
831  G4Polyline frame;
832  frame.push_back(G4Point3D( fSize, fSize, 0.));
833  frame.push_back(G4Point3D(-fSize, fSize, 0.));
834  frame.push_back(G4Point3D(-fSize, -fSize, 0.));
835  frame.push_back(G4Point3D( fSize, -fSize, 0.));
836  frame.push_back(G4Point3D( fSize, fSize, 0.));
837  G4VisAttributes va;
838  va.SetLineWidth(fWidth);
839  va.SetColour(fColour);
840  frame.SetVisAttributes(va);
841  sceneHandler.BeginPrimitives2D();
842  sceneHandler.AddPrimitive(frame);
843  sceneHandler.EndPrimitives2D();
844 }
845 
847 
849  fpCommand = new G4UIcmdWithoutParameter ("/vis/scene/add/hits", this);
850  fpCommand -> SetGuidance ("Adds hits to current scene.");
851  fpCommand -> SetGuidance
852  ("Hits are drawn at end of event when the scene in which"
853  "\nthey are added is current.");
854 }
855 
857  delete fpCommand;
858 }
859 
861  return "";
862 }
863 
865 
867  G4bool warn(verbosity >= G4VisManager::warnings);
868 
869  G4Scene* pScene = fpVisManager->GetCurrentScene();
870  if (!pScene) {
871  if (verbosity >= G4VisManager::errors) {
872  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
873  }
874  return;
875  }
876 
877  G4HitsModel* model = new G4HitsModel;
878  const G4String& currentSceneName = pScene -> GetName ();
879  G4bool successful = pScene -> AddEndOfEventModel (model, warn);
880  if (successful) {
881  if (verbosity >= G4VisManager::confirmations) {
882  G4cout << "Hits, if any, will be drawn at end of run in scene \""
883  << currentSceneName << "\"."
884  << G4endl;
885  }
886  }
887  else G4VisCommandsSceneAddUnsuccessful(verbosity);
888  UpdateVisManagerScene (currentSceneName);
889 }
890 
892 
894  fpCommand = new G4UIcommand("/vis/scene/add/line", this);
895  fpCommand -> SetGuidance ("Adds line to current scene.");
896  G4bool omitable;
897  G4UIparameter* parameter;
898  parameter = new G4UIparameter ("x1", 'd', omitable = false);
899  fpCommand -> SetParameter (parameter);
900  parameter = new G4UIparameter ("y1", 'd', omitable = false);
901  fpCommand -> SetParameter (parameter);
902  parameter = new G4UIparameter ("z1", 'd', omitable = false);
903  fpCommand -> SetParameter (parameter);
904  parameter = new G4UIparameter ("x2", 'd', omitable = false);
905  fpCommand -> SetParameter (parameter);
906  parameter = new G4UIparameter ("y2", 'd', omitable = false);
907  fpCommand -> SetParameter (parameter);
908  parameter = new G4UIparameter ("z2", 'd', omitable = false);
909  fpCommand -> SetParameter (parameter);
910  parameter = new G4UIparameter ("unit", 's', omitable = true);
911  parameter->SetDefaultValue ("m");
912  fpCommand->SetParameter (parameter);
913 }
914 
916  delete fpCommand;
917 }
918 
920  return "";
921 }
922 
924 {
926  G4bool warn(verbosity >= G4VisManager::warnings);
927 
928  G4Scene* pScene = fpVisManager->GetCurrentScene();
929  if (!pScene) {
930  if (verbosity >= G4VisManager::errors) {
931  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
932  }
933  return;
934  }
935 
936  G4String unitString;
937  G4double x1, y1, z1, x2, y2, z2;
938  std::istringstream is(newValue);
939  is >> x1 >> y1 >> z1 >> x2 >> y2 >> z2 >> unitString;
940  G4double unit = G4UIcommand::ValueOf(unitString);
941  x1 *= unit; y1 *= unit; z1 *= unit;
942  x2 *= unit; y2 *= unit; z2 *= unit;
943 
944  Line* line = new Line(x1, y1, z1, x2, y2, z2,
946  G4VModel* model =
948  model->SetType("Line");
949  model->SetGlobalTag("Line");
950  model->SetGlobalDescription("Line: " + newValue);
951  const G4String& currentSceneName = pScene -> GetName ();
952  G4bool successful = pScene -> AddRunDurationModel (model, warn);
953  if (successful) {
954  if (verbosity >= G4VisManager::confirmations) {
955  G4cout << "Line has been added to scene \""
956  << currentSceneName << "\"."
957  << G4endl;
958  }
959  }
960  else G4VisCommandsSceneAddUnsuccessful(verbosity);
961  UpdateVisManagerScene (currentSceneName);
962 }
963 
966  G4double x2, G4double y2, G4double z2,
967  G4double width, const G4Colour& colour):
968  fWidth(width), fColour(colour)
969 {
970  fPolyline.push_back(G4Point3D(x1,y1,z1));
971  fPolyline.push_back(G4Point3D(x2,y2,z2));
972  G4VisAttributes va;
973  va.SetLineWidth(fWidth);
974  va.SetColour(fColour);
975  fPolyline.SetVisAttributes(va);
976 }
977 
978 void G4VisCommandSceneAddLine::Line::operator()
979  (G4VGraphicsScene& sceneHandler, const G4Transform3D&)
980 {
981  sceneHandler.BeginPrimitives();
982  sceneHandler.AddPrimitive(fPolyline);
983  sceneHandler.EndPrimitives();
984 }
985 
987 
989  fpCommand = new G4UIcommand("/vis/scene/add/line2D", this);
990  fpCommand -> SetGuidance ("Adds 2D line to current scene.");
991  G4bool omitable;
992  G4UIparameter* parameter;
993  parameter = new G4UIparameter ("x1", 'd', omitable = false);
994  fpCommand -> SetParameter (parameter);
995  parameter = new G4UIparameter ("y1", 'd', omitable = false);
996  fpCommand -> SetParameter (parameter);
997  parameter = new G4UIparameter ("x2", 'd', omitable = false);
998  fpCommand -> SetParameter (parameter);
999  parameter = new G4UIparameter ("y2", 'd', omitable = false);
1000  fpCommand -> SetParameter (parameter);
1001 }
1002 
1004  delete fpCommand;
1005 }
1006 
1008  return "";
1009 }
1010 
1012 {
1014  G4bool warn(verbosity >= G4VisManager::warnings);
1015 
1016  G4Scene* pScene = fpVisManager->GetCurrentScene();
1017  if (!pScene) {
1018  if (verbosity >= G4VisManager::errors) {
1019  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
1020  }
1021  return;
1022  }
1023 
1024  G4double x1, y1, x2, y2;
1025  std::istringstream is(newValue);
1026  is >> x1 >> y1 >> x2 >> y2;
1027 
1028  Line2D* line2D = new Line2D
1029  (x1, y1, x2, y2, fCurrentLineWidth, fCurrentColour);
1030  G4VModel* model =
1032  model->SetType("Line2D");
1033  model->SetGlobalTag("Line2D");
1034  model->SetGlobalDescription("Line2D: " + newValue);
1035  const G4String& currentSceneName = pScene -> GetName ();
1036  G4bool successful = pScene -> AddRunDurationModel (model, warn);
1037  if (successful) {
1038  if (verbosity >= G4VisManager::confirmations) {
1039  G4cout << "A 2D line has been added to scene \""
1040  << currentSceneName << "\"."
1041  << G4endl;
1042  }
1043  }
1044  else G4VisCommandsSceneAddUnsuccessful(verbosity);
1045  UpdateVisManagerScene (currentSceneName);
1046 }
1047 
1050  G4double x2, G4double y2,
1051  G4double width, const G4Colour& colour):
1052  fWidth(width), fColour(colour)
1053 {
1054  fPolyline.push_back(G4Point3D(x1,y1,0));
1055  fPolyline.push_back(G4Point3D(x2,y2,0));
1056  G4VisAttributes va;
1057  va.SetLineWidth(fWidth);
1058  va.SetColour(fColour);
1059  fPolyline.SetVisAttributes(va);
1060 }
1061 
1062 void G4VisCommandSceneAddLine2D::Line2D::operator()
1063  (G4VGraphicsScene& sceneHandler, const G4Transform3D&)
1064 {
1065  sceneHandler.BeginPrimitives2D();
1066  sceneHandler.AddPrimitive(fPolyline);
1067  sceneHandler.EndPrimitives2D();
1068 }
1069 
1071 
1073  G4bool omitable;
1074  fpCommand = new G4UIcommand ("/vis/scene/add/logicalVolume", this);
1075  fpCommand -> SetGuidance ("Adds a logical volume to the current scene,");
1076  fpCommand -> SetGuidance
1077  ("Shows boolean components (if any), voxels (if any) and readout geometry"
1078  "\n(if any). Note: voxels are not constructed until start of run -"
1079  "\n \"/run/beamOn\". (For voxels without a run, \"/run/beamOn 0\".)");
1080  G4UIparameter* parameter;
1081  parameter = new G4UIparameter ("logical-volume-name", 's', omitable = false);
1082  fpCommand -> SetParameter (parameter);
1083  parameter = new G4UIparameter ("depth-of-descent", 'i', omitable = true);
1084  parameter -> SetGuidance ("Depth of descent of geometry hierarchy.");
1085  parameter -> SetDefaultValue (1);
1086  fpCommand -> SetParameter (parameter);
1087  parameter = new G4UIparameter ("booleans-flag", 'b', omitable = true);
1088  parameter -> SetDefaultValue (true);
1089  fpCommand -> SetParameter (parameter);
1090  parameter = new G4UIparameter ("voxels-flag", 'b', omitable = true);
1091  parameter -> SetDefaultValue (true);
1092  fpCommand -> SetParameter (parameter);
1093  parameter = new G4UIparameter ("readout-flag", 'b', omitable = true);
1094  parameter -> SetDefaultValue (true);
1095  fpCommand -> SetParameter (parameter);
1096  parameter = new G4UIparameter ("axes-flag", 'b', omitable = true);
1097  parameter -> SetDefaultValue (true);
1098  parameter -> SetGuidance ("Set \"false\" to suppress axes.");
1099  fpCommand -> SetParameter (parameter);
1100 }
1101 
1103  delete fpCommand;
1104 }
1105 
1107  return "";
1108 }
1109 
1111  G4String newValue) {
1112 
1114  G4bool warn(verbosity >= G4VisManager::warnings);
1115 
1116  G4Scene* pScene = fpVisManager->GetCurrentScene();
1117  if (!pScene) {
1118  if (verbosity >= G4VisManager::errors) {
1119  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
1120  }
1121  return;
1122  }
1123 
1124  G4String name;
1125  G4int requestedDepthOfDescent;
1126  G4String booleansString, voxelsString, readoutString, axesString;
1127  std::istringstream is (newValue);
1128  is >> name >> requestedDepthOfDescent
1129  >> booleansString >> voxelsString >> readoutString >> axesString;
1130  G4bool booleans = G4UIcommand::ConvertToBool(booleansString);
1131  G4bool voxels = G4UIcommand::ConvertToBool(voxelsString);
1132  G4bool readout = G4UIcommand::ConvertToBool(readoutString);
1133  G4bool axes = G4UIcommand::ConvertToBool(axesString);
1134 
1136  int nLV = pLVStore -> size ();
1137  int iLV;
1138  G4LogicalVolume* pLV = 0;
1139  for (iLV = 0; iLV < nLV; iLV++ ) {
1140  pLV = (*pLVStore) [iLV];
1141  if (pLV -> GetName () == name) break;
1142  }
1143  if (iLV == nLV) {
1144  if (verbosity >= G4VisManager::errors) {
1145  G4cerr << "ERROR: Logical volume " << name
1146  << " not found in logical volume store." << G4endl;
1147  }
1148  return;
1149  }
1150 
1151  const std::vector<G4Scene::Model>& rdModelList =
1152  pScene -> GetRunDurationModelList();
1153  std::vector<G4Scene::Model>::const_iterator i;
1154  for (i = rdModelList.begin(); i != rdModelList.end(); ++i) {
1155  if (i->fpModel->GetGlobalDescription().find("Volume") != std::string::npos) break;
1156  }
1157  if (i != rdModelList.end()) {
1158  if (verbosity >= G4VisManager::errors) {
1159  G4cout << "There is already a volume, \""
1160  << i->fpModel->GetGlobalDescription()
1161  << "\",\n in the run-duration model list of scene \""
1162  << pScene -> GetName()
1163  << "\".\n Your logical volume must be the only volume in the scene."
1164  << "\n Create a new scene and try again:"
1165  << "\n /vis/specify " << name
1166  << "\n or"
1167  << "\n /vis/scene/create"
1168  << "\n /vis/scene/add/logicalVolume " << name
1169  << "\n /vis/sceneHandler/attach"
1170  << "\n (and also, if necessary, /vis/viewer/flush)"
1171  << G4endl;
1172  }
1173  return;
1174  }
1175 
1177  (pLV, requestedDepthOfDescent, booleans, voxels, readout);
1178  const G4String& currentSceneName = pScene -> GetName ();
1179  G4bool successful = pScene -> AddRunDurationModel (model, warn);
1180 
1181  if (successful) {
1182 
1183  G4bool axesSuccessful = false;
1184  if (axes) {
1185  const G4double radius = model->GetExtent().GetExtentRadius();
1186  const G4double axisLengthMax = radius / 2.;
1187  const G4double intLog10Length = std::floor(std::log10(axisLengthMax));
1188  G4double axisLength = std::pow(10,intLog10Length);
1189  if (5.*axisLength < axisLengthMax) axisLength *= 5.;
1190  else if (2.*axisLength < axisLengthMax) axisLength *= 2.;
1191  const G4double axisWidth = axisLength / 20.;
1192  G4VModel* axesModel = new G4AxesModel(0.,0.,0.,axisLength,axisWidth);
1193  axesSuccessful = pScene -> AddRunDurationModel (axesModel, warn);
1194  }
1195 
1196 // if (verbosity >= G4VisManager::warnings) {
1197 // const std::map<G4String,G4AttDef>* attDefs = model->GetAttDefs();
1198 // std::vector<G4AttValue>* attValues = model->CreateCurrentAttValues();
1199 // G4cout << G4AttCheck(attValues, attDefs);
1200 // delete attValues;
1201 // }
1202 
1203  if (verbosity >= G4VisManager::confirmations) {
1204  G4cout << "Logical volume \"" << pLV -> GetName ()
1205  << " with requested depth of descent "
1206  << requestedDepthOfDescent
1207  << ",\n with";
1208  if (!booleans) G4cout << "out";
1209  G4cout << " boolean components, with";
1210  if (!voxels) G4cout << "out";
1211  G4cout << " voxels and with";
1212  if (!readout) G4cout << "out";
1213  G4cout << " readout geometry,"
1214  << "\n has been added to scene \"" << currentSceneName << "\".";
1215  if (axes) {
1216  if (axesSuccessful) {
1217  G4cout <<
1218  "\n Axes have also been added at the origin of local cooordinates.";
1219  } else {
1220  G4cout <<
1221  "\n Axes have not been added for some reason possibly stated above.";
1222  }
1223  }
1224  G4cout << G4endl;
1225  }
1226  }
1227  else {
1229  return;
1230  }
1231 
1232  UpdateVisManagerScene (currentSceneName);
1233 }
1234 
1235 
1237 
1239  G4bool omitable;
1240  fpCommand = new G4UIcommand ("/vis/scene/add/logo", this);
1241  fpCommand -> SetGuidance ("Adds a G4 logo to the current scene.");
1242  fpCommand -> SetGuidance
1243  ("If \"unit\" is \"auto\", height is roughly one tenth of scene extent.");
1244  fpCommand -> SetGuidance
1245  ("\"direction\" is that of outward-facing normal to front face of logo."
1246  "\nIf \"direction\" is \"auto\", logo faces the user in the current viewer.");
1247  fpCommand -> SetGuidance
1248  ("\nIf \"placement\" is \"auto\", logo is placed at bottom right of screen"
1249  "\n when viewed from logo direction.");
1250  G4UIparameter* parameter;
1251  parameter = new G4UIparameter ("height", 'd', omitable = true);
1252  parameter->SetDefaultValue (1.);
1253  fpCommand->SetParameter (parameter);
1254  parameter = new G4UIparameter ("unit", 's', omitable = true);
1255  parameter->SetDefaultValue ("auto");
1256  fpCommand->SetParameter (parameter);
1257  parameter = new G4UIparameter ("direction", 's', omitable = true);
1258  parameter->SetGuidance ("auto|[-]x|[-]y|[-]z");
1259  parameter->SetDefaultValue ("auto");
1260  fpCommand->SetParameter (parameter);
1261  parameter = new G4UIparameter ("red", 'd', omitable = true);
1262  parameter->SetDefaultValue (0.);
1263  fpCommand->SetParameter (parameter);
1264  parameter = new G4UIparameter ("green", 'd', omitable = true);
1265  parameter->SetDefaultValue (1.);
1266  fpCommand->SetParameter (parameter);
1267  parameter = new G4UIparameter ("blue", 'd', omitable = true);
1268  parameter->SetDefaultValue (0.);
1269  fpCommand->SetParameter (parameter);
1270  parameter = new G4UIparameter ("placement", 's', omitable = true);
1271  parameter -> SetParameterCandidates("auto manual");
1272  parameter->SetDefaultValue ("auto");
1273  fpCommand->SetParameter (parameter);
1274  parameter = new G4UIparameter ("xmid", 'd', omitable = true);
1275  parameter->SetDefaultValue (0.);
1276  fpCommand->SetParameter (parameter);
1277  parameter = new G4UIparameter ("ymid", 'd', omitable = true);
1278  parameter->SetDefaultValue (0.);
1279  fpCommand->SetParameter (parameter);
1280  parameter = new G4UIparameter ("zmid", 'd', omitable = true);
1281  parameter->SetDefaultValue (0.);
1282  fpCommand->SetParameter (parameter);
1283  parameter = new G4UIparameter ("unit", 's', omitable = true);
1284  parameter->SetDefaultValue ("m");
1285  fpCommand->SetParameter (parameter);
1286 }
1287 
1289  delete fpCommand;
1290 }
1291 
1293  return "";
1294 }
1295 
1297 
1299  G4bool warn = verbosity >= G4VisManager::warnings;
1300 
1301  G4Scene* pScene = fpVisManager->GetCurrentScene();
1302  if (!pScene) {
1303  if (verbosity >= G4VisManager::errors) {
1304  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
1305  }
1306  return;
1307  } else {
1308  if (pScene->GetExtent().GetExtentRadius() <= 0.) {
1309  if (verbosity >= G4VisManager::errors) {
1310  G4cerr
1311  << "ERROR: Scene has no extent. Add volumes or use \"/vis/scene/add/extent\"."
1312  << G4endl;
1313  }
1314  return;
1315  }
1316  }
1317 
1318  G4VViewer* pViewer = fpVisManager->GetCurrentViewer();
1319  if (!pViewer) {
1320  if (verbosity >= G4VisManager::errors) {
1321  G4cerr <<
1322  "ERROR: G4VisCommandSceneAddLogo::SetNewValue: no viewer."
1323  "\n Auto direction needs a viewer."
1324  << G4endl;
1325  }
1326  return;
1327  }
1328 
1329  G4double userHeight, red, green, blue, xmid, ymid, zmid;
1330  G4String userHeightUnit, direction, auto_manual, positionUnit;
1331  std::istringstream is (newValue);
1332  is >> userHeight >> userHeightUnit >> direction
1333  >> red >> green >> blue
1334  >> auto_manual
1335  >> xmid >> ymid >> zmid >> positionUnit;
1336 
1337  G4double height = userHeight;
1338  const G4VisExtent& sceneExtent = pScene->GetExtent(); // Existing extent.
1339  if (userHeightUnit == "auto") {
1340  height *= 0.2 * sceneExtent.GetExtentRadius();
1341  } else {
1342  height *= G4UIcommand::ValueOf(userHeightUnit);
1343  }
1344 
1345  G4double unit = G4UIcommand::ValueOf(positionUnit);
1346  xmid *= unit; ymid *= unit; zmid *= unit;
1347 
1348  Direction logoDirection = X; // Initialise to keep some compilers happy.
1349  if (direction == "auto") {
1350  // Take cue from viewer
1351  const G4Vector3D& vp =
1353  if (vp.x() > vp.y() && vp.x() > vp.z()) logoDirection = X;
1354  else if (vp.x() < vp.y() && vp.x() < vp.z()) logoDirection = minusX;
1355  else if (vp.y() > vp.x() && vp.y() > vp.z()) logoDirection = Y;
1356  else if (vp.y() < vp.x() && vp.y() < vp.z()) logoDirection = minusY;
1357  else if (vp.z() > vp.x() && vp.z() > vp.y()) logoDirection = Z;
1358  else if (vp.z() < vp.x() && vp.z() < vp.y()) logoDirection = minusZ;
1359  }
1360  else if (direction(0) == 'x') logoDirection = X;
1361  else if (direction(0) == 'y') logoDirection = Y;
1362  else if (direction(0) == 'z') logoDirection = Z;
1363  else if (direction(0) == '-') {
1364  if (direction(1) == 'x') logoDirection = minusX;
1365  else if (direction(1) == 'y') logoDirection = minusY;
1366  else if (direction(1) == 'z') logoDirection = minusZ;
1367  } else {
1368  if (verbosity >= G4VisManager::errors) {
1369  G4cerr << "ERROR: Unrecogniseed direction: \""
1370  << direction << "\"." << G4endl;
1371  return;
1372  }
1373  }
1374 
1375  G4bool autoPlacing = false; if (auto_manual == "auto") autoPlacing = true;
1376  // Parameters read and interpreted.
1377 
1378  // Current scene extent
1379  const G4double xmin = sceneExtent.GetXmin();
1380  const G4double xmax = sceneExtent.GetXmax();
1381  const G4double ymin = sceneExtent.GetYmin();
1382  const G4double ymax = sceneExtent.GetYmax();
1383  const G4double zmin = sceneExtent.GetZmin();
1384  const G4double zmax = sceneExtent.GetZmax();
1385 
1386  // Test existing extent and issue warnings...
1387  G4bool worried = false;
1388  if (sceneExtent.GetExtentRadius() == 0) {
1389  worried = true;
1390  if (verbosity >= G4VisManager::warnings) {
1391  G4cout <<
1392  "WARNING: Existing scene does not yet have any extent."
1393  "\n Maybe you have not yet added any geometrical object."
1394  << G4endl;
1395  }
1396  }
1397 
1398  // Useful constants, etc...
1399  const G4double halfHeight(height / 2.);
1400  const G4double comfort(0.01); // 0.15 seems too big. 0.05 might be better.
1401  const G4double freeHeightFraction (1. + 2. * comfort);
1402 
1403  // Test existing scene for room...
1404  G4bool room = true;
1405  switch (logoDirection) {
1406  case X:
1407  case minusX:
1408  if (freeHeightFraction * (xmax - xmin) < height) room = false; break;
1409  case Y:
1410  case minusY:
1411  if (freeHeightFraction * (ymax - ymin) < height) room = false; break;
1412  case Z:
1413  case minusZ:
1414  if (freeHeightFraction * (zmax - zmin) < height) room = false; break;
1415  }
1416  if (!room) {
1417  worried = true;
1418  if (verbosity >= G4VisManager::warnings) {
1419  G4cout <<
1420  "WARNING: Not enough room in existing scene. Maybe logo is too large."
1421  << G4endl;
1422  }
1423  }
1424  if (worried) {
1425  if (verbosity >= G4VisManager::warnings) {
1426  G4cout <<
1427  "WARNING: The logo you have asked for is bigger than the existing"
1428  "\n scene. Maybe you have added it too soon. It is recommended that"
1429  "\n you add the logo last so that it can be correctly auto-positioned"
1430  "\n so as not to be obscured by any existing object and so that the"
1431  "\n view parameters can be correctly recalculated."
1432  << G4endl;
1433  }
1434  }
1435 
1436  G4double sxmid(xmid), symid(ymid), szmid(zmid);
1437  if (autoPlacing) {
1438  // Aim to place at bottom right of screen when viewed from logoDirection.
1439  // Give some comfort zone.
1440  const G4double xComfort = comfort * (xmax - xmin);
1441  const G4double yComfort = comfort * (ymax - ymin);
1442  const G4double zComfort = comfort * (zmax - zmin);
1443  switch (logoDirection) {
1444  case X: // y-axis up, z-axis to left?
1445  sxmid = xmax + halfHeight + xComfort;
1446  symid = ymin - yComfort;
1447  szmid = zmin - zComfort;
1448  break;
1449  case minusX: // y-axis up, z-axis to right?
1450  sxmid = xmin - halfHeight - xComfort;
1451  symid = ymin - yComfort;
1452  szmid = zmax + zComfort;
1453  break;
1454  case Y: // z-axis up, x-axis to left?
1455  sxmid = xmin - xComfort;
1456  symid = ymax + halfHeight + yComfort;
1457  szmid = zmin - zComfort;
1458  break;
1459  case minusY: // z-axis up, x-axis to right?
1460  sxmid = xmax + xComfort;
1461  symid = ymin - halfHeight - yComfort;
1462  szmid = zmin - zComfort;
1463  break;
1464  case Z: // y-axis up, x-axis to right?
1465  sxmid = xmax + xComfort;
1466  symid = ymin - yComfort;
1467  szmid = zmax + halfHeight + zComfort;
1468  break;
1469  case minusZ: // y-axis up, x-axis to left?
1470  sxmid = xmin - xComfort;
1471  symid = ymin - yComfort;
1472  szmid = zmin - halfHeight - zComfort;
1473  break;
1474  }
1475  }
1476 
1477  G4Transform3D transform;
1478  switch (logoDirection) {
1479  case X: // y-axis up, z-axis to left?
1480  transform = G4RotateY3D(halfpi);
1481  break;
1482  case minusX: // y-axis up, z-axis to right?
1483  transform = G4RotateY3D(-halfpi);
1484  break;
1485  case Y: // z-axis up, x-axis to left?
1486  transform = G4RotateX3D(-halfpi) * G4RotateZ3D(pi);
1487  break;
1488  case minusY: // z-axis up, x-axis to right?
1489  transform = G4RotateX3D(halfpi);
1490  break;
1491  case Z: // y-axis up, x-axis to right?
1492  // No transformation required.
1493  break;
1494  case minusZ: // y-axis up, x-axis to left?
1495  transform = G4RotateY3D(pi);
1496  break;
1497  }
1498  transform = G4Translate3D(sxmid,symid,szmid) * transform;
1499 
1500  G4VisAttributes visAtts(G4Colour(red, green, blue));
1501  visAtts.SetForceSolid(true); // Always solid.
1502 
1503  G4Logo* logo = new G4Logo(height,visAtts);
1504  G4VModel* model =
1506  model->SetType("G4Logo");
1507  model->SetGlobalTag("G4Logo");
1508  model->SetGlobalDescription("G4Logo: " + newValue);
1509  model->SetTransformation(transform);
1510  // Note: it is the responsibility of the model to act upon this, but
1511  // the extent is in local coordinates...
1512  G4double& h = height;
1513  G4double h2 = h/2.;
1514  G4VisExtent extent(-h,h,-h2,h2,-h2,h2);
1515  model->SetExtent(extent);
1516  // This extent gets "added" to existing scene extent in
1517  // AddRunDurationModel below.
1518  const G4String& currentSceneName = pScene -> GetName ();
1519  G4bool successful = pScene -> AddRunDurationModel (model, warn);
1520  if (successful) {
1521  if (verbosity >= G4VisManager::confirmations) {
1522  G4cout << "G4 Logo of height " << userHeight << ' ' << userHeightUnit
1523  << ", " << direction << "-direction, added to scene \""
1524  << currentSceneName << "\"";
1525  if (verbosity >= G4VisManager::parameters) {
1526  G4cout << "\n with extent " << extent
1527  << "\n at " << transform.getRotation()
1528  << transform.getTranslation();
1529  }
1530  G4cout << G4endl;
1531  }
1532  }
1533  else G4VisCommandsSceneAddUnsuccessful(verbosity);
1534  UpdateVisManagerScene (currentSceneName);
1535 }
1536 
1538 (G4double height, const G4VisAttributes& visAtts):
1539  fVisAtts(visAtts)
1540  {
1541  const G4double& h = height;
1542  const G4double h2 = 0.5 * h; // Half height.
1543  const G4double ri = 0.25 * h; // Inner radius.
1544  const G4double ro = 0.5 * h; // Outer radius.
1545  const G4double ro2 = 0.5 * ro; // Half outer radius.
1546  const G4double w = ro - ri; // Width.
1547  const G4double w2 = 0.5 * w; // Half width.
1548  const G4double d2 = 0.2 * h; // Half depth.
1549  const G4double f1 = 0.05 * h; // left edge of stem of "4".
1550  const G4double f2 = -0.3 * h; // bottom edge of cross of "4".
1551  const G4double e = 1.e-4 * h; // epsilon.
1552  const G4double xt = f1, yt = h2; // Top of slope.
1553  const G4double xb = -h2, yb = f2 + w; // Bottom of slope.
1554  const G4double dx = xt - xb, dy = yt - yb;
1555  const G4double angle = std::atan2(dy,dx);
1556  G4RotationMatrix rm;
1557  rm.rotateZ(angle*rad);
1558  const G4double d = std::sqrt(dx * dx + dy * dy);
1559  const G4double ss = h; // Half height of square subtractor
1560  const G4double y8 = ss; // Choose y of subtractor for outer slope.
1561  const G4double x8 = ((-ss * d - dx * (yt - y8)) / dy) + xt;
1562  G4double y9 = ss; // Choose y of subtractor for inner slope.
1563  G4double x9 = ((-(ss - w) * d - dx * (yt - y8)) / dy) + xt;
1564  // But to get inner, we make a triangle translated by...
1565  const G4double xtr = ss - f1, ytr = -ss - f2 -w;
1566  x9 += xtr; y9 += ytr;
1567 
1568  // G...
1569  G4Tubs tG("tG",ri,ro,d2,0.15*pi,1.85*pi);
1570  G4Box bG("bG",w2,ro2,d2);
1571  G4UnionSolid logoG("logoG",&tG,&bG,G4Translate3D(ri+w2,-ro2,0.));
1572  fpG = logoG.CreatePolyhedron();
1573  fpG->SetVisAttributes(&fVisAtts);
1574  fpG->Transform(G4Translate3D(-0.55*h,0.,0.));
1575 
1576  // 4...
1577  G4Box b1("b1",h2,h2,d2);
1578  G4Box bS("bS",ss,ss,d2+e); // Subtractor.
1579  G4Box bS2("bS2",ss,ss,d2+2.*e); // 2nd Subtractor.
1580  G4SubtractionSolid s1("s1",&b1,&bS,G4Translate3D(f1-ss,f2-ss,0.));
1581  G4SubtractionSolid s2("s2",&s1,&bS,G4Translate3D(f1+ss+w,f2-ss,0.));
1582  G4SubtractionSolid s3("s3",&s2,&bS,G4Translate3D(f1+ss+w,f2+ss+w,0.));
1584  ("s4",&s3,&bS,G4Transform3D(rm,G4ThreeVector(x8,y8,0.)));
1585  G4SubtractionSolid s5 // Triangular hole.
1586  ("s5",&bS,&bS2,G4Transform3D(rm,G4ThreeVector(x9,y9,0.)));
1587  G4SubtractionSolid logo4("logo4",&s4,&s5,G4Translate3D(-xtr,-ytr,0.));
1588  fp4 = logo4.CreatePolyhedron();
1589  /* Experiment with creating own polyhedron...
1590  int nNodes = 4;
1591  int nFaces = 4;
1592  double xyz[][3] = {{0,0,0},{1*m,0,0},{0,1*m,0},{0,0,1*m}};
1593  int faces[][4] = {{1,3,2,0},{1,2,4,0},{1,4,3,0},{2,3,4,0}};
1594  fp4 = new G4Polyhedron();
1595  fp4->createPolyhedron(nNodes,nFaces,xyz,faces);
1596  */
1597  fp4->SetVisAttributes(&fVisAtts);
1598  fp4->Transform(G4Translate3D(0.55*h,0.,0.));
1599 }
1600 
1602  delete fpG;
1603  delete fp4;
1604 }
1605 
1606 void G4VisCommandSceneAddLogo::G4Logo::operator()
1607  (G4VGraphicsScene& sceneHandler, const G4Transform3D& transform) {
1608  sceneHandler.BeginPrimitives(transform);
1609  sceneHandler.AddPrimitive(*fpG);
1610  sceneHandler.AddPrimitive(*fp4);
1611  sceneHandler.EndPrimitives();
1612 }
1613 
1615 
1617  G4bool omitable;
1618  fpCommand = new G4UIcommand ("/vis/scene/add/logo2D", this);
1619  fpCommand -> SetGuidance ("Adds 2D logo to current scene.");
1620  G4UIparameter* parameter;
1621  parameter = new G4UIparameter ("size", 'i', omitable = true);
1622  parameter -> SetGuidance ("Screen size of text in pixels.");
1623  parameter -> SetDefaultValue (48);
1624  fpCommand -> SetParameter (parameter);
1625  parameter = new G4UIparameter ("x-position", 'd', omitable = true);
1626  parameter -> SetGuidance ("x screen position in range -1 < x < 1.");
1627  parameter -> SetDefaultValue (-0.9);
1628  fpCommand -> SetParameter (parameter);
1629  parameter = new G4UIparameter ("y-position", 'd', omitable = true);
1630  parameter -> SetGuidance ("y screen position in range -1 < y < 1.");
1631  parameter -> SetDefaultValue (-0.9);
1632  fpCommand -> SetParameter (parameter);
1633  parameter = new G4UIparameter ("layout", 's', omitable = true);
1634  parameter -> SetGuidance ("Layout, i.e., adjustment: left|centre|right.");
1635  parameter -> SetDefaultValue ("left");
1636  fpCommand -> SetParameter (parameter);
1637 }
1638 
1640  delete fpCommand;
1641 }
1642 
1644  return "";
1645 }
1646 
1648 {
1650  G4bool warn(verbosity >= G4VisManager::warnings);
1651 
1652  G4Scene* pScene = fpVisManager->GetCurrentScene();
1653  if (!pScene) {
1654  if (verbosity >= G4VisManager::errors) {
1655  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
1656  }
1657  return;
1658  }
1659 
1660  G4int size;
1661  G4double x, y;
1662  G4String layoutString;
1663  std::istringstream is(newValue);
1664  is >> size >> x >> y >> layoutString;
1665  G4Text::Layout layout = G4Text::right;
1666  if (layoutString(0) == 'l') layout = G4Text::left;
1667  else if (layoutString(0) == 'c') layout = G4Text::centre;
1668  else if (layoutString(0) == 'r') layout = G4Text::right;
1669 
1670  Logo2D* logo2D = new Logo2D(fpVisManager, size, x, y, layout);
1671  G4VModel* model =
1673  model->SetType("G4Logo2D");
1674  model->SetGlobalTag("G4Logo2D");
1675  model->SetGlobalDescription("G4Logo2D: " + newValue);
1676  const G4String& currentSceneName = pScene -> GetName ();
1677  G4bool successful = pScene -> AddRunDurationModel (model, warn);
1678  if (successful) {
1679  if (verbosity >= G4VisManager::confirmations) {
1680  G4cout << "2D logo has been added to scene \""
1681  << currentSceneName << "\"."
1682  << G4endl;
1683  }
1684  }
1685  else G4VisCommandsSceneAddUnsuccessful(verbosity);
1686  UpdateVisManagerScene (currentSceneName);
1687 }
1688 
1689 void G4VisCommandSceneAddLogo2D::Logo2D::operator()
1690  (G4VGraphicsScene& sceneHandler, const G4Transform3D&)
1691 {
1692  G4Text text("Geant4", G4Point3D(fX, fY, 0.));
1693  text.SetScreenSize(fSize);
1694  text.SetLayout(fLayout);
1695  G4VisAttributes textAtts(G4Colour::Brown());
1696  text.SetVisAttributes(textAtts);
1697  sceneHandler.BeginPrimitives2D();
1698  sceneHandler.AddPrimitive(text);
1699  sceneHandler.EndPrimitives2D();
1700 }
1701 
1703 
1705  G4bool ommitable;
1706  fpCommand = new G4UIcmdWithAnInteger ("/vis/scene/add/magneticField", this);
1707  fpCommand->SetParameterName("nDataPointsPerHalfScene",ommitable=true);
1708  fpCommand->SetDefaultValue(10);
1709  fpCommand -> SetGuidance
1710  ("Adds magnetic field representation to current scene.");
1711  fpCommand -> SetGuidance
1712  ("The parameter is no. of data points per half scene. So, possibly, at"
1713  "\nmaximum, the number of data points sampled is (2*n+1)^3, which can grow"
1714  "\nlarge--be warned!"
1715  "\nYou might find that your scene is cluttered by thousands of arrows for"
1716  "\nthe default number of data points, so try reducing to 2 or 3, e.g:"
1717  "\n /vis/scene/add/magneticField 3"
1718  "\nor, if only a small part of the scene has a field:"
1719  "\n /vis/scene/add/magneticField 50 # or more");
1720  fpCommand -> SetGuidance
1721  ("In the arrow representation, the length of the arrow is proportional"
1722  "\nto the magnitude of the field and the colour is mapped onto the range"
1723  "\nas a fraction of the maximum magnitude: 0->0.5->1 is blue->green->red.");
1724 }
1725 
1727  delete fpCommand;
1728 }
1729 
1731  return "";
1732 }
1733 
1735 (G4UIcommand*, G4String newValue) {
1736 
1738  G4bool warn(verbosity >= G4VisManager::warnings);
1739 
1740  G4Scene* pScene = fpVisManager->GetCurrentScene();
1741  if (!pScene) {
1742  if (verbosity >= G4VisManager::errors) {
1743  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
1744  }
1745  return;
1746  }
1747 
1748  G4int nDataPointsPerHalfScene = fpCommand->GetNewIntValue(newValue);
1749 
1750  G4MagneticFieldModel* model =
1751  new G4MagneticFieldModel(nDataPointsPerHalfScene);
1752  const G4String& currentSceneName = pScene -> GetName ();
1753  G4bool successful = pScene -> AddRunDurationModel (model, warn);
1754  if (successful) {
1755  if (verbosity >= G4VisManager::confirmations) {
1756  G4cout << "Magnetic field, if any, will be drawn in scene \""
1757  << currentSceneName
1758  << "\"\n with "
1759  << nDataPointsPerHalfScene << " data points per half scene."
1760  << G4endl;
1761  }
1762  }
1763  else G4VisCommandsSceneAddUnsuccessful(verbosity);
1764  UpdateVisManagerScene (currentSceneName);
1765 }
1766 
1768 
1770  G4bool omitable;
1771  fpCommand = new G4UIcmdWithAString ("/vis/scene/add/psHits", this);
1772  fpCommand -> SetGuidance
1773  ("Adds Primitive Scorer Hits (PSHits) to current scene.");
1774  fpCommand -> SetGuidance
1775  ("PSHits are drawn at end of run when the scene in which"
1776  "\nthey are added is current.");
1777  fpCommand -> SetGuidance
1778  ("Optional parameter specifies name of scoring map. By default all"
1779  "\nscoring maps registered with the G4ScoringManager are drawn.");
1780  fpCommand -> SetParameterName ("mapname", omitable = true);
1781  fpCommand -> SetDefaultValue ("all");
1782 }
1783 
1785  delete fpCommand;
1786 }
1787 
1789  return "";
1790 }
1791 
1793 (G4UIcommand*, G4String newValue)
1794 {
1796  G4bool warn(verbosity >= G4VisManager::warnings);
1797 
1798  G4Scene* pScene = fpVisManager->GetCurrentScene();
1799  if (!pScene) {
1800  if (verbosity >= G4VisManager::errors) {
1801  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
1802  }
1803  return;
1804  }
1805 
1806  G4PSHitsModel* model = new G4PSHitsModel(newValue);
1807  const G4String& currentSceneName = pScene -> GetName ();
1808  G4bool successful = pScene -> AddEndOfRunModel (model, warn);
1809  if (successful) {
1810  if (verbosity >= G4VisManager::confirmations) {
1811  if (newValue == "all") {
1812  G4cout << "All Primitive Scorer hits";
1813  } else {
1814  G4cout << "Hits of Primitive Scorer \"" << newValue << '"';
1815  }
1816  G4cout << " will be drawn at end of run in scene \""
1817  << currentSceneName << "\"."
1818  << G4endl;
1819  }
1820  }
1821  else G4VisCommandsSceneAddUnsuccessful(verbosity);
1822  UpdateVisManagerScene (currentSceneName);
1823 }
1824 
1826 
1828  G4bool omitable;
1829  fpCommand = new G4UIcommand ("/vis/scene/add/scale", this);
1830  fpCommand -> SetGuidance
1831  ("Adds an annotated scale line to the current scene.");
1832  fpCommand -> SetGuidance
1833  ("If \"unit\" is \"auto\", length is roughly one tenth of the scene extent.");
1834  fpCommand -> SetGuidance
1835  ("If \"direction\" is \"auto\", scale is roughly in the plane of the current view.");
1836  fpCommand -> SetGuidance
1837  ("If \"placement\" is \"auto\", scale is placed at bottom left of current view."
1838  "\n Otherwise placed at (xmid,ymid,zmid).");
1839  fpCommand -> SetGuidance (G4Scale::GetGuidanceString());
1840  G4UIparameter* parameter;
1841  parameter = new G4UIparameter ("length", 'd', omitable = true);
1842  parameter->SetDefaultValue (1.);
1843  fpCommand->SetParameter (parameter);
1844  parameter = new G4UIparameter ("unit", 's', omitable = true);
1845  parameter->SetDefaultValue ("auto");
1846  fpCommand->SetParameter (parameter);
1847  parameter = new G4UIparameter ("direction", 's', omitable = true);
1848  parameter->SetGuidance ("auto|x|y|z");
1849  parameter->SetDefaultValue ("auto");
1850  fpCommand->SetParameter (parameter);
1851  parameter = new G4UIparameter ("red", 'd', omitable = true);
1852  parameter->SetDefaultValue (1.);
1853  fpCommand->SetParameter (parameter);
1854  parameter = new G4UIparameter ("green", 'd', omitable = true);
1855  parameter->SetDefaultValue (0.);
1856  fpCommand->SetParameter (parameter);
1857  parameter = new G4UIparameter ("blue", 'd', omitable = true);
1858  parameter->SetDefaultValue (0.);
1859  fpCommand->SetParameter (parameter);
1860  parameter = new G4UIparameter ("placement", 's', omitable = true);
1861  parameter -> SetParameterCandidates("auto manual");
1862  parameter->SetDefaultValue ("auto");
1863  fpCommand->SetParameter (parameter);
1864  parameter = new G4UIparameter ("xmid", 'd', omitable = true);
1865  parameter->SetDefaultValue (0.);
1866  fpCommand->SetParameter (parameter);
1867  parameter = new G4UIparameter ("ymid", 'd', omitable = true);
1868  parameter->SetDefaultValue (0.);
1869  fpCommand->SetParameter (parameter);
1870  parameter = new G4UIparameter ("zmid", 'd', omitable = true);
1871  parameter->SetDefaultValue (0.);
1872  fpCommand->SetParameter (parameter);
1873  parameter = new G4UIparameter ("unit", 's', omitable = true);
1874  parameter->SetDefaultValue ("m");
1875  fpCommand->SetParameter (parameter);
1876 }
1877 
1879  delete fpCommand;
1880 }
1881 
1883  return "";
1884 }
1885 
1887 
1889  G4bool warn = verbosity >= G4VisManager::warnings;
1890 
1891  G4Scene* pScene = fpVisManager->GetCurrentScene();
1892  if (!pScene) {
1893  if (verbosity >= G4VisManager::errors) {
1894  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
1895  }
1896  return;
1897  } else {
1898  if (pScene->GetExtent().GetExtentRadius() <= 0.) {
1899  if (verbosity >= G4VisManager::errors) {
1900  G4cerr
1901  << "ERROR: Scene has no extent. Add volumes or use \"/vis/scene/add/extent\"."
1902  << G4endl;
1903  }
1904  return;
1905  }
1906  }
1907 
1908  G4double userLength, red, green, blue, xmid, ymid, zmid;
1909  G4String userLengthUnit, direction, auto_manual, positionUnit;
1910  std::istringstream is (newValue);
1911  is >> userLength >> userLengthUnit >> direction
1912  >> red >> green >> blue
1913  >> auto_manual
1914  >> xmid >> ymid >> zmid >> positionUnit;
1915 
1916  G4double length = userLength;
1917  const G4VisExtent& sceneExtent = pScene->GetExtent(); // Existing extent.
1918  if (userLengthUnit == "auto") {
1919  const G4double lengthMax = 0.5 * sceneExtent.GetExtentRadius();
1920  const G4double intLog10Length = std::floor(std::log10(lengthMax));
1921  length = std::pow(10,intLog10Length);
1922  if (5.*length < lengthMax) length *= 5.;
1923  else if (2.*length < lengthMax) length *= 2.;
1924  } else {
1925  length *= G4UIcommand::ValueOf(userLengthUnit);
1926  }
1927  G4String annotation = G4BestUnit(length,"Length");
1928 
1929  G4double unit = G4UIcommand::ValueOf(positionUnit);
1930  xmid *= unit; ymid *= unit; zmid *= unit;
1931 
1932  G4Scale::Direction scaleDirection (G4Scale::x);
1933  if (direction(0) == 'y') scaleDirection = G4Scale::y;
1934  if (direction(0) == 'z') scaleDirection = G4Scale::z;
1935 
1936  G4VViewer* pViewer = fpVisManager->GetCurrentViewer();
1937  if (!pViewer) {
1938  if (verbosity >= G4VisManager::errors) {
1939  G4cerr <<
1940  "ERROR: G4VisCommandSceneAddScale::SetNewValue: no viewer."
1941  "\n Auto direction needs a viewer."
1942  << G4endl;
1943  }
1944  return;
1945  }
1946 
1947  const G4Vector3D& vp =
1949  const G4Vector3D& up =
1950  pViewer->GetViewParameters().GetUpVector();
1951 
1952  if (direction == "auto") { // Takes cue from viewer.
1953  if (std::abs(vp.x()) > std::abs(vp.y()) &&
1954  std::abs(vp.x()) > std::abs(vp.z())) { // x viewpoint
1955  if (std::abs(up.y()) > std::abs(up.z())) scaleDirection = G4Scale::z;
1956  else scaleDirection = G4Scale::y;
1957  }
1958  else if (std::abs(vp.y()) > std::abs(vp.x()) &&
1959  std::abs(vp.y()) > std::abs(vp.z())) { // y viewpoint
1960  if (std::abs(up.x()) > std::abs(up.z())) scaleDirection = G4Scale::z;
1961  else scaleDirection = G4Scale::x;
1962  }
1963  else if (std::abs(vp.z()) > std::abs(vp.x()) &&
1964  std::abs(vp.z()) > std::abs(vp.y())) { // z viewpoint
1965  if (std::abs(up.y()) > std::abs(up.x())) scaleDirection = G4Scale::x;
1966  else scaleDirection = G4Scale::y;
1967  }
1968  }
1969 
1970  G4bool autoPlacing = false; if (auto_manual == "auto") autoPlacing = true;
1971  // Parameters read and interpreted.
1972 
1973  // Useful constants, etc...
1974  const G4double halfLength(length / 2.);
1975  const G4double comfort(0.01); // 0.15 seems too big. 0.05 might be better.
1976  const G4double freeLengthFraction (1. + 2. * comfort);
1977 
1978  const G4double xmin = sceneExtent.GetXmin();
1979  const G4double xmax = sceneExtent.GetXmax();
1980  const G4double ymin = sceneExtent.GetYmin();
1981  const G4double ymax = sceneExtent.GetYmax();
1982  const G4double zmin = sceneExtent.GetZmin();
1983  const G4double zmax = sceneExtent.GetZmax();
1984 
1985  // Test existing extent and issue warnings...
1986  G4bool worried = false;
1987  if (sceneExtent.GetExtentRadius() == 0) {
1988  worried = true;
1989  if (verbosity >= G4VisManager::warnings) {
1990  G4cout <<
1991  "WARNING: Existing scene does not yet have any extent."
1992  "\n Maybe you have not yet added any geometrical object."
1993  << G4endl;
1994  }
1995  }
1996  // Test existing scene for room...
1997  G4bool room = true;
1998  switch (scaleDirection) {
1999  case G4Scale::x:
2000  if (freeLengthFraction * (xmax - xmin) < length) room = false; break;
2001  case G4Scale::y:
2002  if (freeLengthFraction * (ymax - ymin) < length) room = false; break;
2003  case G4Scale::z:
2004  if (freeLengthFraction * (zmax - zmin) < length) room = false; break;
2005  }
2006  if (!room) {
2007  worried = true;
2008  if (verbosity >= G4VisManager::warnings) {
2009  G4cout <<
2010  "WARNING: Not enough room in existing scene. Maybe scale is too long."
2011  << G4endl;
2012  }
2013  }
2014  if (worried) {
2015  if (verbosity >= G4VisManager::warnings) {
2016  G4cout <<
2017  "WARNING: The scale you have asked for is bigger than the existing"
2018  "\n scene. Maybe you have added it too soon. It is recommended that"
2019  "\n you add the scale last so that it can be correctly auto-positioned"
2020  "\n so as not to be obscured by any existing object and so that the"
2021  "\n view parameters can be correctly recalculated."
2022  << G4endl;
2023  }
2024  }
2025 
2026  // Let's go ahead a construct a scale and a scale model. Since the
2027  // placing is done here, this G4Scale is *not* auto-placed...
2028  G4Scale scale(length, annotation, scaleDirection,
2029  false, xmid, ymid, zmid,
2031  G4VisAttributes visAttr(G4Colour(red, green, blue));
2032  scale.SetVisAttributes(visAttr);
2033  G4VModel* model = new G4ScaleModel(scale);
2034  G4String globalDescription = model->GetGlobalDescription();
2035  globalDescription += " (" + newValue + ")";
2036  model->SetGlobalDescription(globalDescription);
2037 
2038  // Now figure out the extent...
2039  //
2040  // From the G4Scale.hh:
2041  //
2042  // This creates a representation of annotated line in the specified
2043  // direction with tick marks at the end. If autoPlacing is true it
2044  // is required to be centred at the front, right, bottom corner of
2045  // the world space, comfortably outside the existing bounding
2046  // box/sphere so that existing objects do not obscure it. Otherwise
2047  // it is required to be drawn with mid-point at (xmid, ymid, zmid).
2048  //
2049  // The auto placing algorithm might be:
2050  // x = xmin + (1 + comfort) * (xmax - xmin)
2051  // y = ymin - comfort * (ymax - ymin)
2052  // z = zmin + (1 + comfort) * (zmax - zmin)
2053  // if direction == x then (x - length,y,z) to (x,y,z)
2054  // if direction == y then (x,y,z) to (x,y + length,z)
2055  // if direction == z then (x,y,z - length) to (x,y,z)
2056  //
2057  // End of clip from G4Scale.hh:
2058  //
2059  // Implement this in two parts. Here, use the scale's extent to
2060  // "expand" the scene's extent. Then rendering - in
2061  // G4VSceneHandler::AddPrimitive(const G4Scale&) - simply has to
2062  // ensure it's within the new extent.
2063  //
2064 
2065  G4double sxmid(xmid), symid(ymid), szmid(zmid);
2066  if (autoPlacing) {
2067  // Aim to place at bottom right of screen in current view.
2068  // Give some comfort zone.
2069  const G4double xComfort = comfort * (xmax - xmin);
2070  const G4double yComfort = comfort * (ymax - ymin);
2071  const G4double zComfort = comfort * (zmax - zmin);
2072  switch (scaleDirection) {
2073  case G4Scale::x:
2074  if (vp.z() > 0.) {
2075  sxmid = xmax + xComfort;
2076  symid = ymin - yComfort;
2077  szmid = zmin - zComfort;
2078  } else {
2079  sxmid = xmin - xComfort;
2080  symid = ymin - yComfort;
2081  szmid = zmax + zComfort;
2082  }
2083  break;
2084  case G4Scale::y:
2085  if (vp.x() > 0.) {
2086  sxmid = xmin - xComfort;
2087  symid = ymax + yComfort;
2088  szmid = zmin - zComfort;
2089  } else {
2090  sxmid = xmax + xComfort;
2091  symid = ymin - yComfort;
2092  szmid = zmin - zComfort;
2093  }
2094  break;
2095  case G4Scale::z:
2096  if (vp.x() > 0.) {
2097  sxmid = xmax + xComfort;
2098  symid = ymin - yComfort;
2099  szmid = zmax + zComfort;
2100  } else {
2101  sxmid = xmin - xComfort;
2102  symid = ymin - yComfort;
2103  szmid = zmax + zComfort;
2104  }
2105  break;
2106  }
2107  }
2108 
2109  /* Old code - kept for future reference.
2110  G4double sxmid(xmid), symid(ymid), szmid(zmid);
2111  if (autoPlacing) {
2112  sxmid = xmin + onePlusComfort * (xmax - xmin);
2113  symid = ymin - comfort * (ymax - ymin);
2114  szmid = zmin + onePlusComfort * (zmax - zmin);
2115  switch (scaleDirection) {
2116  case G4Scale::x:
2117  sxmid -= halfLength;
2118  break;
2119  case G4Scale::y:
2120  symid += halfLength;
2121  break;
2122  case G4Scale::z:
2123  szmid -= halfLength;
2124  break;
2125  }
2126  }
2127  */
2128 
2129  /* sxmin, etc., not actually used. Comment out to prevent compiler
2130  warnings but keep in case need in future. Extract transform and
2131  scaleExtent into reduced code below.
2132  G4double sxmin(sxmid), sxmax(sxmid);
2133  G4double symin(symid), symax(symid);
2134  G4double szmin(szmid), szmax(szmid);
2135  G4Transform3D transform;
2136  G4VisExtent scaleExtent;
2137  switch (scaleDirection) {
2138  case G4Scale::x:
2139  sxmin = sxmid - halfLength;
2140  sxmax = sxmid + halfLength;
2141  scaleExtent = G4VisExtent(-halfLength,halfLength,0,0,0,0);
2142  break;
2143  case G4Scale::y:
2144  symin = symid - halfLength;
2145  symax = symid + halfLength;
2146  transform = G4RotateZ3D(halfpi);
2147  scaleExtent = G4VisExtent(0,0,-halfLength,halfLength,0,0);
2148  break;
2149  case G4Scale::z:
2150  szmin = szmid - halfLength;
2151  szmax = szmid + halfLength;
2152  transform = G4RotateY3D(halfpi);
2153  scaleExtent = G4VisExtent(0,0,0,0,-halfLength,halfLength);
2154  break;
2155  }
2156  */
2157  G4Transform3D transform;
2158  G4VisExtent scaleExtent;
2159  switch (scaleDirection) {
2160  case G4Scale::x:
2161  scaleExtent = G4VisExtent(-halfLength,halfLength,0,0,0,0);
2162  break;
2163  case G4Scale::y:
2164  transform = G4RotateZ3D(halfpi);
2165  scaleExtent = G4VisExtent(0,0,-halfLength,halfLength,0,0);
2166  break;
2167  case G4Scale::z:
2168  transform = G4RotateY3D(halfpi);
2169  scaleExtent = G4VisExtent(0,0,0,0,-halfLength,halfLength);
2170  break;
2171  }
2172  transform = G4Translate3D(sxmid,symid,szmid) * transform;
2174 
2175 
2176  model->SetTransformation(transform);
2177  // Note: it is the responsibility of the model to act upon this, but
2178  // the extent is in local coordinates...
2179  model->SetExtent(scaleExtent);
2180  // This extent gets "added" to existing scene extent in
2181  // AddRunDurationModel below.
2182 
2183  const G4String& currentSceneName = pScene -> GetName ();
2184  G4bool successful = pScene -> AddRunDurationModel (model, warn);
2185  if (successful) {
2186  if (verbosity >= G4VisManager::confirmations) {
2187  G4cout << "Scale of " << annotation
2188  << " added to scene \"" << currentSceneName << "\".";
2189  if (verbosity >= G4VisManager::parameters) {
2190  G4cout << "\n with extent " << scaleExtent
2191  << "\n at " << transform.getRotation()
2192  << transform.getTranslation();
2193  }
2194  G4cout << G4endl;
2195  }
2196  }
2197  else G4VisCommandsSceneAddUnsuccessful(verbosity);
2198  UpdateVisManagerScene (currentSceneName);
2199 }
2200 
2201 
2203 
2205  G4bool omitable;
2206  fpCommand = new G4UIcommand ("/vis/scene/add/text", this);
2207  fpCommand -> SetGuidance ("Adds text to current scene.");
2208  fpCommand -> SetGuidance
2209  ("Use \"/vis/set/textColour\" to set colour.");
2210  fpCommand -> SetGuidance
2211  ("Use \"/vis/set/textLayout\" to set layout:");
2212  G4UIparameter* parameter;
2213  parameter = new G4UIparameter ("x", 'd', omitable = true);
2214  parameter->SetDefaultValue (0);
2215  fpCommand->SetParameter (parameter);
2216  parameter = new G4UIparameter ("y", 'd', omitable = true);
2217  parameter->SetDefaultValue (0);
2218  fpCommand->SetParameter (parameter);
2219  parameter = new G4UIparameter ("z", 'd', omitable = true);
2220  parameter->SetDefaultValue (0);
2221  fpCommand->SetParameter (parameter);
2222  parameter = new G4UIparameter ("unit", 's', omitable = true);
2223  parameter->SetDefaultValue ("m");
2224  fpCommand->SetParameter (parameter);
2225  parameter = new G4UIparameter ("font_size", 'd', omitable = true);
2226  parameter->SetDefaultValue (12);
2227  parameter->SetGuidance ("pixels");
2228  fpCommand->SetParameter (parameter);
2229  parameter = new G4UIparameter ("x_offset", 'd', omitable = true);
2230  parameter->SetDefaultValue (0);
2231  parameter->SetGuidance ("pixels");
2232  fpCommand->SetParameter (parameter);
2233  parameter = new G4UIparameter ("y_offset", 'd', omitable = true);
2234  parameter->SetDefaultValue (0);
2235  parameter->SetGuidance ("pixels");
2236  fpCommand->SetParameter (parameter);
2237  parameter = new G4UIparameter ("text", 's', omitable = true);
2238  parameter->SetGuidance ("The rest of the line is text.");
2239  parameter->SetDefaultValue ("Hello G4");
2240  fpCommand->SetParameter (parameter);
2241 }
2242 
2244  delete fpCommand;
2245 }
2246 
2248  return "";
2249 }
2250 
2252 
2254  G4bool warn = verbosity >= G4VisManager::warnings;
2255 
2256  G4Scene* pScene = fpVisManager->GetCurrentScene();
2257  if (!pScene) {
2258  if (verbosity >= G4VisManager::errors) {
2259  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
2260  }
2261  return;
2262  }
2263 
2264  G4Tokenizer next(newValue);
2265  G4double x = StoD(next());
2266  G4double y = StoD(next());
2267  G4double z = StoD(next());
2268  G4String unitString = next();
2269  G4double font_size = StoD(next());
2270  G4double x_offset = StoD(next());
2271  G4double y_offset = StoD(next());
2272  G4String text = next("\n");
2273 
2274  G4double unit = G4UIcommand::ValueOf(unitString);
2275  x *= unit; y *= unit; z *= unit;
2276 
2277  G4Text g4text(text, G4Point3D(x,y,z));
2279  g4text.SetVisAttributes(visAtts);
2280  g4text.SetLayout(fCurrentTextLayout);
2281  g4text.SetScreenSize(font_size);
2282  g4text.SetOffset(x_offset,y_offset);
2283  G4VModel* model = new G4TextModel(g4text);
2284  const G4String& currentSceneName = pScene -> GetName ();
2285  G4bool successful = pScene -> AddRunDurationModel (model, warn);
2286  if (successful) {
2287  if (verbosity >= G4VisManager::confirmations) {
2288  G4cout << "Text \"" << text
2289  << "\" has been added to scene \"" << currentSceneName << "\"."
2290  << G4endl;
2291  }
2292  }
2293  else G4VisCommandsSceneAddUnsuccessful(verbosity);
2294  UpdateVisManagerScene (currentSceneName);
2295 }
2296 
2297 
2299 
2301  G4bool omitable;
2302  fpCommand = new G4UIcommand ("/vis/scene/add/text2D", this);
2303  fpCommand -> SetGuidance ("Adds 2D text to current scene.");
2304  fpCommand -> SetGuidance
2305  ("Use \"/vis/set/textColour\" to set colour.");
2306  fpCommand -> SetGuidance
2307  ("Use \"/vis/set/textLayout\" to set layout:");
2308  G4UIparameter* parameter;
2309  parameter = new G4UIparameter ("x", 'd', omitable = true);
2310  parameter->SetDefaultValue (0);
2311  fpCommand->SetParameter (parameter);
2312  parameter = new G4UIparameter ("y", 'd', omitable = true);
2313  parameter->SetDefaultValue (0);
2314  fpCommand->SetParameter (parameter);
2315  parameter = new G4UIparameter ("font_size", 'd', omitable = true);
2316  parameter->SetDefaultValue (12);
2317  parameter->SetGuidance ("pixels");
2318  fpCommand->SetParameter (parameter);
2319  parameter = new G4UIparameter ("x_offset", 'd', omitable = true);
2320  parameter->SetDefaultValue (0);
2321  parameter->SetGuidance ("pixels");
2322  fpCommand->SetParameter (parameter);
2323  parameter = new G4UIparameter ("y_offset", 'd', omitable = true);
2324  parameter->SetDefaultValue (0);
2325  parameter->SetGuidance ("pixels");
2326  fpCommand->SetParameter (parameter);
2327  parameter = new G4UIparameter ("text", 's', omitable = true);
2328  parameter->SetGuidance ("The rest of the line is text.");
2329  parameter->SetDefaultValue ("Hello G4");
2330  fpCommand->SetParameter (parameter);
2331 }
2332 
2334  delete fpCommand;
2335 }
2336 
2338  return "";
2339 }
2340 
2342 
2344  G4bool warn = verbosity >= G4VisManager::warnings;
2345 
2346  G4Scene* pScene = fpVisManager->GetCurrentScene();
2347  if (!pScene) {
2348  if (verbosity >= G4VisManager::errors) {
2349  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
2350  }
2351  return;
2352  }
2353 
2354  G4Tokenizer next(newValue);
2355  G4double x = StoD(next());
2356  G4double y = StoD(next());
2357  G4double font_size = StoD(next());
2358  G4double x_offset = StoD(next());
2359  G4double y_offset = StoD(next());
2360  G4String text = next("\n");
2361 
2362  G4Text g4text(text, G4Point3D(x,y,0.));
2364  g4text.SetVisAttributes(visAtts);
2365  g4text.SetLayout(fCurrentTextLayout);
2366  g4text.SetScreenSize(font_size);
2367  g4text.SetOffset(x_offset,y_offset);
2368  G4Text2D* g4text2D = new G4Text2D(g4text);
2369  G4VModel* model =
2371  model->SetType("Text2D");
2372  model->SetGlobalTag("Text2D");
2373  model->SetGlobalDescription("Text2D: " + newValue);
2374  const G4String& currentSceneName = pScene -> GetName ();
2375  G4bool successful = pScene -> AddRunDurationModel (model, warn);
2376  if (successful) {
2377  if (verbosity >= G4VisManager::confirmations) {
2378  G4cout << "2D text \"" << text
2379  << "\" has been added to scene \"" << currentSceneName << "\"."
2380  << G4endl;
2381  }
2382  }
2383  else G4VisCommandsSceneAddUnsuccessful(verbosity);
2384  UpdateVisManagerScene (currentSceneName);
2385 }
2386 
2388  fText(text)
2389 {}
2390 
2391 void G4VisCommandSceneAddText2D::G4Text2D::operator()
2392  (G4VGraphicsScene& sceneHandler, const G4Transform3D& transform) {
2393  sceneHandler.BeginPrimitives2D(transform);
2394  sceneHandler.AddPrimitive(fText);
2395  sceneHandler.EndPrimitives2D();
2396 }
2397 
2398 
2400 
2402  G4bool omitable;
2404  ("/vis/scene/add/trajectories", this);
2405  fpCommand -> SetGuidance
2406  ("Adds trajectories to current scene.");
2407  fpCommand -> SetGuidance
2408  ("Causes trajectories, if any, to be drawn at the end of processing an"
2409  "\nevent. Switches on trajectory storing and sets the"
2410  "\ndefault trajectory type.");
2411  fpCommand -> SetGuidance
2412  ("The command line parameter list determines the default trajectory type."
2413  "\nIf it contains the string \"smooth\", auxiliary inter-step points will"
2414  "\nbe inserted to improve the smoothness of the drawing of a curved"
2415  "\ntrajectory."
2416  "\nIf it contains the string \"rich\", significant extra information will"
2417  "\nbe stored in the trajectory (G4RichTrajectory) amenable to modeling"
2418  "\nand filtering with \"/vis/modeling/trajectories/create/drawByAttribute\""
2419  "\nand \"/vis/filtering/trajectories/create/attributeFilter\" commands."
2420  "\nIt may contain both strings in any order.");
2421  fpCommand -> SetGuidance
2422  ("\nTo switch off trajectory storing: \"/tracking/storeTrajectory 0\"."
2423  "\nSee also \"/vis/scene/endOfEventAction\".");
2424  fpCommand -> SetGuidance
2425  ("Note: This only sets the default. Independently of the result of this"
2426  "\ncommand, a user may instantiate a trajectory that overrides this default"
2427  "\nin PreUserTrackingAction.");
2428  fpCommand -> SetParameterName ("default-trajectory-type", omitable = true);
2429  fpCommand -> SetDefaultValue ("");
2430 }
2431 
2433  delete fpCommand;
2434 }
2435 
2437  return "";
2438 }
2439 
2441  G4String newValue) {
2442 
2444  G4bool warn = verbosity >= G4VisManager::warnings;
2445 
2446  G4Scene* pScene = fpVisManager->GetCurrentScene();
2447  if (!pScene) {
2448  if (verbosity >= G4VisManager::errors) {
2449  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
2450  }
2451  return;
2452  }
2453 
2454  G4bool smooth = false;
2455  G4bool rich = false;
2456  if (newValue.find("smooth") != std::string::npos) smooth = true;
2457  if (newValue.find("rich") != std::string::npos) rich = true;
2458  if (newValue.size() && !(rich || smooth)) {
2459  if (verbosity >= G4VisManager::errors) {
2460  G4cerr << "ERROR: Unrecognised parameter \"" << newValue << "\""
2461  "\n No action taken."
2462  << G4endl;
2463  }
2464  return;
2465  }
2466 
2467  G4UImanager* UImanager = G4UImanager::GetUIpointer();
2468  G4int keepVerbose = UImanager->GetVerboseLevel();
2469  G4int newVerbose = 2;
2470  UImanager->SetVerboseLevel(newVerbose);
2471  G4String defaultTrajectoryType;
2472  if (smooth && rich) {
2473  UImanager->ApplyCommand("/tracking/storeTrajectory 4");
2474  defaultTrajectoryType = "G4RichTrajectory configured for smooth steps";
2475  } else if (smooth) {
2476  UImanager->ApplyCommand("/tracking/storeTrajectory 2");
2477  defaultTrajectoryType = "G4SmoothTrajectory";
2478  } else if (rich) {
2479  UImanager->ApplyCommand("/tracking/storeTrajectory 3");
2480  defaultTrajectoryType = "G4RichTrajectory";
2481  } else {
2482  UImanager->ApplyCommand("/tracking/storeTrajectory 1");
2483  defaultTrajectoryType = "G4Trajectory";
2484  }
2485  UImanager->SetVerboseLevel(keepVerbose);
2486 
2487  if (verbosity >= G4VisManager::errors) {
2488  G4cout <<
2489  "Attributes available for modeling and filtering with"
2490  "\n \"/vis/modeling/trajectories/create/drawByAttribute\" and"
2491  "\n \"/vis/filtering/trajectories/create/attributeFilter\" commands:"
2492  << G4endl;
2494  if (rich) {
2497  } else if (smooth) {
2500  } else {
2502  << *G4TrajectoryPoint().GetAttDefs();
2503  }
2504  }
2505 
2507  const G4String& currentSceneName = pScene -> GetName ();
2508  pScene -> AddEndOfEventModel (model, warn);
2509 
2510  if (verbosity >= G4VisManager::confirmations) {
2511  G4cout << "Default trajectory type " << defaultTrajectoryType
2512  << "\n will be used to store trajectories for scene \""
2513  << currentSceneName << "\"."
2514  << G4endl;
2515  }
2516 
2517  if (verbosity >= G4VisManager::warnings) {
2518  G4cout <<
2519  "WARNING: Trajectory storing has been requested. This action may be"
2520  "\n reversed with \"/tracking/storeTrajectory 0\"."
2521  << G4endl;
2522  }
2523  UpdateVisManagerScene (currentSceneName);
2524 }
2525 
2527 
2529  G4bool omitable;
2530  fpCommand = new G4UIcmdWithAString("/vis/scene/add/userAction",this);
2531  fpCommand -> SetGuidance
2532  ("Add named Vis User Action to current scene.");
2533  fpCommand -> SetGuidance
2534  ("Attempts to match search string to name of action - use unique sub-string.");
2535  fpCommand -> SetGuidance
2536  ("(Use /vis/list to see names of registered actions.)");
2537  fpCommand -> SetGuidance
2538  ("If name == \"all\" (default), all actions are added.");
2539  fpCommand -> SetParameterName("action-name", omitable = true);
2540  fpCommand -> SetDefaultValue("all");
2541 }
2542 
2544  delete fpCommand;
2545 }
2546 
2548  return "";
2549 }
2550 
2552 (G4UIcommand*, G4String newValue) {
2553 
2555 
2556  G4Scene* pScene = fpVisManager->GetCurrentScene();
2557  if (!pScene) {
2558  if (verbosity >= G4VisManager::errors) {
2559  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
2560  }
2561  return;
2562  }
2563 
2564  G4bool any = false;
2565 
2566  const std::vector<G4VisManager::UserVisAction>& runDurationUserVisActions =
2568  for (size_t i = 0; i < runDurationUserVisActions.size(); i++) {
2569  const G4String& name = runDurationUserVisActions[i].fName;
2570  G4VUserVisAction* visAction = runDurationUserVisActions[i].fpUserVisAction;
2571  if (newValue == "all" || name.find(newValue) != std::string::npos) {
2572  any = true;
2573  AddVisAction(name,visAction,pScene,runDuration,verbosity);
2574  }
2575  }
2576 
2577  const std::vector<G4VisManager::UserVisAction>& endOfEventUserVisActions =
2579  for (size_t i = 0; i < endOfEventUserVisActions.size(); i++) {
2580  const G4String& name = endOfEventUserVisActions[i].fName;
2581  G4VUserVisAction* visAction = endOfEventUserVisActions[i].fpUserVisAction;
2582  if (newValue == "all" || name.find(newValue) != std::string::npos) {
2583  any = true;
2584  AddVisAction(name,visAction,pScene,endOfEvent,verbosity);
2585  }
2586  }
2587 
2588  const std::vector<G4VisManager::UserVisAction>& endOfRunUserVisActions =
2590  for (size_t i = 0; i < endOfRunUserVisActions.size(); i++) {
2591  const G4String& name = endOfRunUserVisActions[i].fName;
2592  G4VUserVisAction* visAction = endOfRunUserVisActions[i].fpUserVisAction;
2593  if (newValue == "all" || name.find(newValue) != std::string::npos) {
2594  any = true;
2595  AddVisAction(name,visAction,pScene,endOfRun,verbosity);
2596  }
2597  }
2598 
2599  if (!any) {
2600  if (verbosity >= G4VisManager::warnings) {
2601  G4cout << "WARNING: No User Vis Action registered." << G4endl;
2602  }
2603  return;
2604  }
2605 
2606  const G4String& currentSceneName = pScene -> GetName ();
2607  UpdateVisManagerScene (currentSceneName);
2608 }
2609 
2611 (const G4String& name,
2612  G4VUserVisAction* visAction,
2613  G4Scene* pScene,
2615  G4VisManager::Verbosity verbosity)
2616 {
2617  G4bool warn = verbosity >= G4VisManager::warnings;
2618 
2619  const std::map<G4VUserVisAction*,G4VisExtent>& visExtentMap =
2621  G4VisExtent extent;
2622  std::map<G4VUserVisAction*,G4VisExtent>::const_iterator i =
2623  visExtentMap.find(visAction);
2624  if (i != visExtentMap.end()) extent = i->second;
2625  if (warn) {
2626  if (extent.GetExtentRadius() <= 0.) {
2627  G4cout << "WARNING: User Vis Action extent is null." << G4endl;
2628  }
2629  }
2630 
2631  G4VModel* model = new G4CallbackModel<G4VUserVisAction>(visAction);
2632  model->SetType("User Vis Action");
2633  model->SetGlobalTag(name);
2634  model->SetGlobalDescription(name);
2635  model->SetExtent(extent);
2636  G4bool successful = false;;
2637  switch (type) {
2638  case runDuration:
2639  successful = pScene -> AddRunDurationModel (model, warn);
2640  break;
2641  case endOfEvent:
2642  successful = pScene -> AddEndOfEventModel (model, warn);
2643  break;
2644  case endOfRun:
2645  successful = pScene -> AddEndOfRunModel (model, warn);
2646  break;
2647  }
2648  if (successful && verbosity >= G4VisManager::confirmations) {
2649  const G4String& currentSceneName = pScene -> GetName ();
2650  G4cout << "User Vis Action added to scene \""
2651  << currentSceneName << "\"";
2652  if (verbosity >= G4VisManager::parameters) {
2653  G4cout << "\n with extent " << extent;
2654  }
2655  G4cout << G4endl;
2656  }
2657  else G4VisCommandsSceneAddUnsuccessful(verbosity);
2658 }
2659 
2661 
2663  G4bool omitable;
2664  fpCommand = new G4UIcommand ("/vis/scene/add/volume", this);
2665  fpCommand -> SetGuidance
2666  ("Adds a physical volume to current scene, with optional clipping volume.");
2667  fpCommand -> SetGuidance
2668  ("If physical-volume-name is \"world\" (the default), the top of the"
2669  "\nmain geometry tree (material world) is added. If \"worlds\", the"
2670  "\ntop of all worlds - material world and parallel worlds, if any - are"
2671  "\nadded. Otherwise a search of all worlds is made, taking the first"
2672  "\nmatching occurence only. To see a representation of the geometry"
2673  "\nhierarchy of the worlds, try \"/vis/drawTree [worlds]\" or one of the"
2674  "\ndriver/browser combinations that have the required functionality, e.g., HepRep.");
2675  fpCommand -> SetGuidance
2676  ("If clip-volume-type is specified, the subsequent parameters are used to"
2677  "\nto define a clipping volume. For example,"
2678  "\n\"vis/scene/add/volume ! ! ! -box km 0 1 0 1 0 1\" will draw the world"
2679  "\nwith the positive octant cut away.");
2680  fpCommand -> SetGuidance
2681  ("If clip-volume-type is prepended with '-', the clip-volume is subtracted"
2682  "\n(cutaway). (This is the default if there is no prepended character.)"
2683  "\nIf '*' is prepended, the intersection of the physical-volume and the"
2684  "\nclip-volume is made. (You can make a section/DCUT with a thin box, for"
2685  "\nexample).");
2686  fpCommand -> SetGuidance
2687  ("For \"box\", the parameters are xmin,xmax,ymin,ymax,zmin,zmax."
2688  "\nOnly \"box\" is programmed at present.");
2689  G4UIparameter* parameter;
2690  parameter = new G4UIparameter ("physical-volume-name", 's', omitable = true);
2691  parameter -> SetDefaultValue ("world");
2692  fpCommand -> SetParameter (parameter);
2693  parameter = new G4UIparameter ("copy-no", 'i', omitable = true);
2694  parameter -> SetGuidance
2695  ("If negative, matches any copy no. First name match is taken.");
2696  parameter -> SetDefaultValue (-1);
2697  fpCommand -> SetParameter (parameter);
2698  parameter = new G4UIparameter ("depth-of-descent", 'i', omitable = true);
2699  parameter -> SetGuidance
2700  ("Depth of descent of geometry hierarchy. Default = unlimited depth.");
2701  parameter -> SetDefaultValue (G4Scene::UNLIMITED);
2702  fpCommand -> SetParameter (parameter);
2703  parameter = new G4UIparameter ("clip-volume-type", 's', omitable = true);
2704  parameter -> SetParameterCandidates("none box -box *box");
2705  parameter -> SetDefaultValue ("none");
2706  parameter -> SetGuidance("[-|*]type. See general guidance.");
2707  fpCommand -> SetParameter (parameter);
2708  parameter = new G4UIparameter ("parameter-unit", 's', omitable = true);
2709  parameter -> SetDefaultValue ("m");
2710  fpCommand -> SetParameter (parameter);
2711  parameter = new G4UIparameter ("parameter-1", 'd', omitable = true);
2712  parameter -> SetDefaultValue (0.);
2713  fpCommand -> SetParameter (parameter);
2714  parameter = new G4UIparameter ("parameter-2", 'd', omitable = true);
2715  parameter -> SetDefaultValue (0.);
2716  fpCommand -> SetParameter (parameter);
2717  parameter = new G4UIparameter ("parameter-3", 'd', omitable = true);
2718  parameter -> SetDefaultValue (0.);
2719  fpCommand -> SetParameter (parameter);
2720  parameter = new G4UIparameter ("parameter-4", 'd', omitable = true);
2721  parameter -> SetDefaultValue (0.);
2722  fpCommand -> SetParameter (parameter);
2723  parameter = new G4UIparameter ("parameter-5", 'd', omitable = true);
2724  parameter -> SetDefaultValue (0.);
2725  fpCommand -> SetParameter (parameter);
2726  parameter = new G4UIparameter ("parameter-6", 'd', omitable = true);
2727  parameter -> SetDefaultValue (0.);
2728  fpCommand -> SetParameter (parameter);
2729 }
2730 
2732  delete fpCommand;
2733 }
2734 
2736  return "world 0 -1";
2737 }
2738 
2740  G4String newValue) {
2741 
2743  G4bool warn = verbosity >= G4VisManager::warnings;
2744 
2745  G4Scene* pScene = fpVisManager->GetCurrentScene();
2746  if (!pScene) {
2747  if (verbosity >= G4VisManager::errors) {
2748  G4cerr << "ERROR: No current scene. Please create one." << G4endl;
2749  }
2750  return;
2751  }
2752 
2753  G4String name, clipVolumeType, parameterUnit;
2754  G4int copyNo, requestedDepthOfDescent;
2755  G4double param1, param2, param3, param4, param5, param6;
2756  std::istringstream is (newValue);
2757  is >> name >> copyNo >> requestedDepthOfDescent
2758  >> clipVolumeType >> parameterUnit
2759  >> param1 >> param2 >> param3 >> param4 >> param5 >> param6;
2761  G4PhysicalVolumeModel::subtraction; // Default subtraction mode.
2762  if (clipVolumeType[size_t(0)] == '-') {
2763  clipVolumeType = clipVolumeType.substr(1); // Remove first character.
2764  } else if (clipVolumeType[size_t(0)] == '*') {
2765  clippingMode = G4PhysicalVolumeModel::intersection;
2766  clipVolumeType = clipVolumeType.substr(1);
2767  }
2768  G4double unit = G4UIcommand::ValueOf(parameterUnit);
2769  param1 *= unit; param2 *= unit; param3 *= unit;
2770  param4 *= unit; param5 *= unit; param6 *= unit;
2771 
2772  G4TransportationManager* transportationManager =
2774 
2775  size_t nWorlds = transportationManager->GetNoWorlds();
2776  if (nWorlds > 1) { // Parallel worlds in operation...
2777  if (verbosity >= G4VisManager::warnings) {
2778  static G4bool warned = false;
2779  if (!warned && name != "worlds") {
2780  G4cout <<
2781  "WARNING: Parallel worlds in operation. To visualise, specify"
2782  "\n \"worlds\" or the parallel world volume or sub-volume name"
2783  "\n and control visibility with /vis/geometry."
2784  << G4endl;
2785  std::vector<G4VPhysicalVolume*>::iterator iterWorld =
2786  transportationManager->GetWorldsIterator();
2787  for (size_t i = 0; i < nWorlds; ++i, ++iterWorld) {
2788  G4cout << " World " << i << ": " << (*iterWorld)->GetName()
2789  << G4endl;
2790  warned = true;
2791  }
2792  }
2793  }
2794  }
2795 
2796  G4VPhysicalVolume* world = *(transportationManager->GetWorldsIterator());
2797 
2798  if (!world) {
2799  if (verbosity >= G4VisManager::errors) {
2800  G4cerr <<
2801  "ERROR: G4VisCommandSceneAddVolume::SetNewValue:"
2802  "\n No world. Maybe the geometry has not yet been defined."
2803  "\n Try \"/run/initialize\""
2804  << G4endl;
2805  }
2806  return;
2807  }
2808 
2809  const std::vector<G4Scene::Model>& rdModelList = pScene -> GetRunDurationModelList();
2810  std::vector<G4Scene::Model>::const_iterator it;
2811  for (it = rdModelList.begin(); it != rdModelList.end(); ++it) {
2812  if (it->fpModel->GetGlobalDescription().find("G4PhysicalVolumeModel")
2813  != std::string::npos) {
2814  if (((G4PhysicalVolumeModel*)(it->fpModel))->GetTopPhysicalVolume () == world) break;
2815  }
2816  }
2817  if (it != rdModelList.end()) {
2818  if (verbosity >= G4VisManager::warnings) {
2819  G4cout << "WARNING: There is already a volume, \""
2820  << it -> fpModel -> GetGlobalDescription()
2821  << "\",\n in the run-duration model list of scene \""
2822  << pScene -> GetName()
2823  << "\".\n To get a clean scene:"
2824  << "\n /vis/drawVolume " << name
2825  << "\n or"
2826  << "\n /vis/scene/create"
2827  << "\n /vis/scene/add/volume " << name
2828  << "\n /vis/sceneHandler/attach"
2829  << "\n (and also, if necessary, /vis/viewer/flush)"
2830  << G4endl;
2831  }
2832  return;
2833  }
2834 
2835  std::vector<G4PhysicalVolumeModel*> models;
2836  std::vector<G4VPhysicalVolume*> foundVolumes;
2837  G4VPhysicalVolume* foundWorld = 0;
2839  typedef std::vector<PVNodeID> PVPath;
2840  PVPath foundFullPVPath;
2841  std::vector<G4int> foundDepths;
2842  std::vector<G4Transform3D> transformations;
2843 
2844  if (name == "world") {
2845 
2846  models.push_back
2847  (new G4PhysicalVolumeModel (world, requestedDepthOfDescent));
2848  foundVolumes.push_back(world);
2849  foundDepths.push_back(0);
2850  transformations.push_back(G4Transform3D());
2851 
2852  } else if (name == "worlds") {
2853 
2854  if (nWorlds == 0) {
2855  if (verbosity >= G4VisManager::warnings) {
2856  G4cout <<
2857  "WARNING: G4VisCommandSceneAddVolume::SetNewValue:"
2858  "\n Parallel worlds requested but none exist."
2859  "\n Just adding material world."
2860  << G4endl;
2861  }
2862  }
2863  std::vector<G4VPhysicalVolume*>::iterator iterWorld =
2864  transportationManager->GetWorldsIterator();
2865  for (size_t i = 0; i < nWorlds; ++i, ++iterWorld) {
2866  models.push_back
2867  (new G4PhysicalVolumeModel (*iterWorld, requestedDepthOfDescent));
2868  foundVolumes.push_back(*iterWorld);
2869  foundDepths.push_back(0);
2870  transformations.push_back(G4Transform3D());
2871  }
2872 
2873  } else { // Search all worlds...
2874 
2875  std::vector<G4VPhysicalVolume*>::iterator iterWorld =
2876  transportationManager->GetWorldsIterator();
2877  for (size_t i = 0; i < nWorlds; ++i, ++iterWorld) {
2878  G4PhysicalVolumeModel searchModel (*iterWorld); // Unlimited depth.
2879  G4ModelingParameters mp; // Default - no culling.
2880  searchModel.SetModelingParameters (&mp);
2881  G4PhysicalVolumeSearchScene searchScene (&searchModel, name, copyNo);
2882  searchModel.DescribeYourselfTo (searchScene); // Initiate search.
2883  G4VPhysicalVolume* foundVolume = searchScene.GetFoundVolume ();
2884  if (foundVolume) {
2885  foundWorld = *iterWorld;
2886  foundVolumes.push_back(foundVolume);
2887  foundFullPVPath = searchScene.GetFoundFullPVPath();
2888  foundDepths.push_back(searchScene.GetFoundDepth());
2889  transformations.push_back(searchScene.GetFoundTransformation());
2890  break;
2891  }
2892  }
2893 
2894  if (foundVolumes.size()) {
2895  for (size_t i = 0; i < foundVolumes.size(); ++i) {
2896  G4PhysicalVolumeModel* foundPVModel = new G4PhysicalVolumeModel
2897  (foundVolumes[i], requestedDepthOfDescent, transformations[i]);
2898  foundFullPVPath.pop_back(); // "Base" is "Found - 1".
2899  foundPVModel->SetBaseFullPVPath(foundFullPVPath);
2900  models.push_back(foundPVModel);
2901  }
2902  } else {
2903  if (verbosity >= G4VisManager::errors) {
2904  G4cerr << "ERROR: Volume \"" << name << "\"";
2905  if (copyNo >= 0) {
2906  G4cerr << ", copy no. " << copyNo << ",";
2907  }
2908  G4cerr << " not found." << G4endl;
2909  }
2910  return;
2911  }
2912  }
2913 
2914  if (clipVolumeType == "box") {
2915  const G4double dX = (param2 - param1) / 2.;
2916  const G4double dY = (param4 - param3) / 2.;
2917  const G4double dZ = (param6 - param5) / 2.;
2918  const G4double x0 = (param2 + param1) / 2.;
2919  const G4double y0 = (param4 + param3) / 2.;
2920  const G4double z0 = (param6 + param5) / 2.;
2921  G4VSolid* clippingSolid =
2922  new G4DisplacedSolid
2923  ("_displaced_clipping_box",
2924  new G4Box("_clipping_box",dX,dY,dZ),
2925  G4Translate3D(x0,y0,z0));
2926  for (size_t i = 0; i < foundVolumes.size(); ++i) {
2927  models[i]->SetClippingSolid(clippingSolid);
2928  models[i]->SetClippingMode(clippingMode);
2929  }
2930  } // If any other shape consider NumberOfRotationSides!!!!!!!!!!!
2931 
2932  const G4String& currentSceneName = pScene -> GetName ();
2933  G4bool failure = true;
2934  for (size_t i = 0; i < foundVolumes.size(); ++i) {
2935  G4bool successful = pScene -> AddRunDurationModel (models[i], warn);
2936  if (successful) {
2937  failure = false;
2938  if (verbosity >= G4VisManager::confirmations) {
2939  G4cout << "First occurrence of \""
2940  << foundVolumes[i] -> GetName ()
2941  << "\"";
2942  if (copyNo >= 0) {
2943  G4cout << ", copy no. " << copyNo << ",";
2944  }
2945  G4cout << "\n found ";
2946  if (foundWorld)
2947  G4cout << "in world \"" << foundWorld->GetName() << "\" ";
2948  G4cout << "at depth " << foundDepths[i]
2949  << ",\n with a requested depth of further descent of ";
2950  if (requestedDepthOfDescent < 0) {
2951  G4cout << "<0 (unlimited)";
2952  }
2953  else {
2954  G4cout << requestedDepthOfDescent;
2955  }
2956  G4cout << ",\n has been added to scene \"" << currentSceneName << "\"."
2957  << G4endl;
2958  }
2959  }
2960  }
2961 
2962  if (failure) {
2964  return;
2965  }
2966 
2967  UpdateVisManagerScene (currentSceneName);
2968 }
void SetNewValue(G4UIcommand *command, G4String newValue)
void SetParameter(G4UIparameter *const newParameter)
Definition: G4UIcommand.hh:152
G4String GetCurrentValue(G4UIcommand *command)
void SetColour(const G4Colour &)
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
G4Logo(G4double height, const G4VisAttributes &)
void SetGlobalTag(const G4String &)
void SetNewValue(G4UIcommand *command, G4String newValue)
HepGeom::RotateX3D G4RotateX3D
Definition: G4Text.hh:73
G4String GetCurrentValue(G4UIcommand *command)
HepGeom::RotateY3D G4RotateY3D
Definition: test07.cc:36
G4double GetXmin() const
Definition: G4VisExtent.hh:89
static const G4double f2
CLHEP::Hep3Vector G4ThreeVector
CLHEP::HepRotation G4RotationMatrix
static G4Colour fCurrentColour
void SetNewValue(G4UIcommand *command, G4String newValue)
G4String GetCurrentValue(G4UIcommand *command)
std::vector< G4VPhysicalVolume * >::iterator GetWorldsIterator()
void SetNewValue(G4UIcommand *command, G4String newValue)
void SetLayout(Layout)
G4String GetCurrentValue(G4UIcommand *command)
G4double z
Definition: TRTMaterials.hh:39
Definition: G4Box.hh:64
G4String GetCurrentValue(G4UIcommand *command)
void SetNewValue(G4UIcommand *command, G4String newValue)
void UpdateVisManagerScene(const G4String &sceneName="")
G4String name
Definition: TRTMaterials.hh:40
G4String GetCurrentValue(G4UIcommand *command)
const G4double pi
void SetLineWidth(G4double)
Definition: G4Tubs.hh:85
static G4Text::Layout fCurrentTextLayout
Line(G4double x1, G4double y1, G4double z1, G4double x2, G4double y2, G4double z2, G4double width, const G4Colour &colour)
const G4ViewParameters & GetViewParameters() const
G4double GetXmax() const
Definition: G4VisExtent.hh:90
const G4ModelingParameters * GetModelingParameters() const
void SetDefaultValue(const char *theDefaultValue)
virtual void BeginPrimitives(const G4Transform3D &objectTransformation=G4Transform3D())=0
#define width
const std::vector< const G4Event * > * GetEventVector() const
Definition: G4Run.hh:115
void SetNewValue(G4UIcommand *command, G4String newValue)
void SetNewValue(G4UIcommand *command, G4String newValue)
void SetNewValue(G4UIcommand *command, G4String newValue)
HepGeom::Point3D< G4double > G4Point3D
Definition: G4Point3D.hh:35
HepGeom::Vector3D< G4double > G4Vector3D
Definition: G4Vector3D.hh:35
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
static G4Colour Brown()
Definition: G4Colour.hh:147
void SetForceSolid(G4bool)
Direction
Definition: G4Scale.hh:42
void SetExtent(const G4VisExtent &)
const G4VisExtent & GetExtent() const
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
G4int GetVerboseLevel() const
Definition: G4UImanager.hh:227
G4String GetCurrentValue(G4UIcommand *command)
void SetNewValue(G4UIcommand *command, G4String newValue)
static G4double fCurrentTextSize
Definition: test07.cc:36
int G4int
Definition: G4Types.hh:78
void SetModelingParameters(const G4ModelingParameters *)
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
HepGeom::RotateZ3D G4RotateZ3D
const G4Run * GetCurrentRun() const
void SetVerboseLevel(G4int val)
Definition: G4UImanager.hh:225
void SetTransformation(const G4Transform3D &)
static void G4VisCommandsSceneAddUnsuccessful(G4VisManager::Verbosity verbosity)
virtual void AddPrimitive(const G4Polyline &)=0
static G4UImanager * GetUIpointer()
Definition: G4UImanager.cc:58
G4int GetEventID() const
Definition: G4Event.hh:151
void SetNewValue(G4UIcommand *command, G4String newValue)
G4String GetCurrentValue(G4UIcommand *command)
G4String GetCurrentValue(G4UIcommand *command)
G4String GetCurrentValue(G4UIcommand *command)
Line2D(G4double x1, G4double y1, G4double x2, G4double y2, G4double width, const G4Colour &colour)
G4String GetCurrentValue(G4UIcommand *command)
static G4StateManager * GetStateManager()
void AddVisAction(const G4String &name, G4VUserVisAction *, G4Scene *, ActionType, G4VisManager::Verbosity)
G4double GetExtentRadius() const
Definition: G4VisExtent.cc:73
const G4Vector3D & GetViewpointDirection() const
void SetNewValue(G4UIcommand *command, G4String newValue)
G4GLOB_DLL std::ostream G4cout
const G4Event * GetEvent() const
void SetNewValue(G4UIcommand *command, G4String newValue)
const G4VisExtent & GetExtent() const
const G4String & GetGlobalDescription() const
const G4String & GetName() const
void SetNewValue(G4UIcommand *command, G4String newValue)
static G4bool ConvertToBool(const char *st)
Definition: G4UIcommand.cc:425
void SetNewValue(G4UIcommand *command, G4String newValue)
static const double deg
Definition: G4SIunits.hh:133
G4double GetYmax() const
Definition: G4VisExtent.hh:92
bool G4bool
Definition: G4Types.hh:79
G4String GetCurrentValue(G4UIcommand *command)
G4PhysicalVolumeModel::G4PhysicalVolumeNodeID PVNodeID
G4Polyhedron * CreatePolyhedron() const
Extent(G4double xmin, G4double xmax, G4double ymin, G4double ymax, G4double zmin, G4double zmax)
G4String GetCurrentValue(G4UIcommand *command)
static G4MTRunManager * GetMasterRunManager()
void SetType(const G4String &)
G4int GetRunID() const
Definition: G4Run.hh:76
Definition: G4Run.hh:46
static const G4String & GetGuidanceString()
Definition: G4Scale.cc:66
static G4LogicalVolumeStore * GetInstance()
const std::vector< UserVisAction > & GetEndOfEventUserVisActions() const
HepGeom::Transform3D G4Transform3D
std::vector< PVNodeID > PVPath
G4double GetZmax() const
Definition: G4VisExtent.hh:94
G4ApplicationState GetCurrentState() const
const std::map< G4VUserVisAction *, G4VisExtent > & GetUserVisActionExtents() const
void SetBaseFullPVPath(const std::vector< G4PhysicalVolumeNodeID > &baseFullPVPath)
void SetNewValue(G4UIcommand *command, G4String newValue)
G4bool IsMultithreadedApplication()
Definition: G4Threading.cc:134
static const G4double f1
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
static const double rad
Definition: G4SIunits.hh:130
void SetVisAttributes(const G4VisAttributes *)
Definition: G4Visible.cc:80
const std::vector< UserVisAction > & GetRunDurationUserVisActions() const
static G4TransportationManager * GetTransportationManager()
G4String GetCurrentValue(G4UIcommand *command)
void SetNewValue(G4UIcommand *command, G4String newValue)
static G4RunManager * GetRunManager()
Definition: G4RunManager.cc:79
const std::vector< G4PhysicalVolumeModel::G4PhysicalVolumeNodeID > & GetFoundFullPVPath() const
virtual void EndPrimitives2D()=0
G4String GetCurrentValue(G4UIcommand *command)
void SetOffset(double dx, double dy)
void SetNewValue(G4UIcommand *command, G4String newValue)
virtual void BeginPrimitives2D(const G4Transform3D &objectTransformation=G4Transform3D())=0
G4double StoD(G4String s)
static G4double ValueOf(const char *unitName)
Definition: G4UIcommand.cc:308
Layout
Definition: G4Text.hh:77
G4VPhysicalVolume * GetFoundVolume() const
const std::map< G4String, G4AttDef > * GetAttDefs() const
G4String GetCurrentValue(G4UIcommand *command)
static G4double fCurrentLineWidth
G4double GetZmin() const
Definition: G4VisExtent.hh:93
static G4Colour fCurrentTextColour
static Verbosity GetVerbosity()
static const G4double b1
HepGeom::Translate3D G4Translate3D
G4double GetYmin() const
Definition: G4VisExtent.hh:91
G4UIcmdWithoutParameter * fpCommand
G4int GetNumberOfEventToBeProcessed() const
Definition: G4Run.hh:83
size_t GetNoWorlds() const
G4String GetCurrentValue(G4UIcommand *command)
#define G4endl
Definition: G4ios.hh:61
void SetNewValue(G4UIcommand *command, G4String newValue)
G4UIcmdWithoutParameter * fpCommand
void SetNewValue(G4UIcommand *command, G4String newValue)
G4VViewer * GetCurrentViewer() const
virtual void EndPrimitives()=0
const G4Transform3D & GetFoundTransformation() const
const std::vector< UserVisAction > & GetEndOfRunUserVisActions() const
G4String GetCurrentValue(G4UIcommand *command)
G4String GetCurrentValue(G4UIcommand *command)
double G4double
Definition: G4Types.hh:76
G4String GetCurrentValue(G4UIcommand *command)
void SetNewValue(G4UIcommand *command, G4String newValue)
void SetNewValue(G4UIcommand *command, G4String newValue)
void SetGuidance(const char *theGuidance)
const G4Vector3D & GetUpVector() const
static const G4double d2
G4String GetCurrentValue(G4UIcommand *command)
Arrow2D(G4double x1, G4double y1, G4double x2, G4double y2, G4double width, const G4Colour &colour)
G4ApplicationState
void SetGlobalDescription(const G4String &)
G4String GetCurrentValue(G4UIcommand *command)
G4int ApplyCommand(const char *aCommand)
Definition: G4UImanager.cc:432
G4Scene * GetCurrentScene() const
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const
void SetNewValue(G4UIcommand *command, G4String newValue)
G4GLOB_DLL std::ostream G4cerr
void SetScreenSize(G4double)
void DescribeYourselfTo(G4VGraphicsScene &)
static G4VisManager * fpVisManager
virtual const std::map< G4String, G4AttDef > * GetAttDefs() const