Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4TrajectoryDrawerUtils.cc
이 파일의 문서화 페이지로 가기
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 // $Id: G4TrajectoryDrawerUtils.cc 107920 2017-12-13 08:06:00Z gcosmo $
27 //
28 // Jane Tinslay, John Allison, Joseph Perl November 2005
29 //
31 #include "G4Colour.hh"
32 #include "G4Polyline.hh"
33 #include "G4Polymarker.hh"
34 #include "G4VTrajectory.hh"
35 #include "G4VTrajectoryPoint.hh"
36 #include "G4VisAttributes.hh"
37 #include "G4VisTrajContext.hh"
38 #include "G4VVisManager.hh"
39 #include "G4UIcommand.hh"
40 #include "G4AttValue.hh"
41 
42 namespace G4TrajectoryDrawerUtils {
43 
45 
47  (const G4VTrajectory& traj,
49  G4Polyline& trajectoryLine,
50  G4Polymarker& auxiliaryPoints,
51  G4Polymarker& stepPoints,
52  std::vector<G4double>& trajectoryLineTimes,
53  std::vector<G4double>& auxiliaryPointTimes,
54  std::vector<G4double>& stepPointTimes)
55  {
56  TimesValidity validity = InvalidTimes;
57  if (context.GetTimeSliceInterval()) validity = ValidTimes;
58 
59  // Memory for last trajectory point position for auxiliary point
60  // time interpolation algorithm. There are no auxiliary points
61  // for the first trajectory point, so its initial value is
62  // immaterial.
63  G4ThreeVector lastTrajectoryPointPosition;
64 
65  // Keep positions. Don't store unless first or different.
66  std::vector<G4ThreeVector> positions;
67 
68  for (G4int iPoint=0; iPoint<traj.GetPointEntries(); iPoint++) {
69 
70  G4VTrajectoryPoint* aTrajectoryPoint = traj.GetPoint(iPoint);
71  const G4ThreeVector& trajectoryPointPosition =
72  aTrajectoryPoint->GetPosition();
73 
74  // Only store if first or if different
75  if (positions.size() == 0 ||
76  trajectoryPointPosition != positions[positions.size()-1]) {
77 
78  // Pre- and Post-Point times from the trajectory point...
79  G4double trajectoryPointPreTime = -std::numeric_limits<double>::max();
80  G4double trajectoryPointPostTime = std::numeric_limits<double>::max();
81 
82  if (context.GetTimeSliceInterval() && validity == ValidTimes) {
83 
84  std::vector<G4AttValue>* trajectoryPointAttValues =
85  aTrajectoryPoint->CreateAttValues();
86  if (!trajectoryPointAttValues) {
87  static G4bool warnedNoAttValues = false;
88  if (!warnedNoAttValues) {
89  G4cout <<
90  "*************************************************************************"
91  "\n* WARNING: G4TrajectoryDrawerUtils::GetPointsAndTimes: no att values."
92  "\n*************************************************************************"
93  << G4endl;
94  warnedNoAttValues = true;
95  }
96  validity = InvalidTimes;
97  } else {
98  G4bool foundPreTime = false, foundPostTime = false;
99  for (std::vector<G4AttValue>::iterator i =
100  trajectoryPointAttValues->begin();
101  i != trajectoryPointAttValues->end(); ++i) {
102  if (i->GetName() == "PreT") {
103  trajectoryPointPreTime =
105  foundPreTime = true;
106  }
107  if (i->GetName() == "PostT") {
108  trajectoryPointPostTime =
110  foundPostTime = true;
111  }
112  }
113  if (!foundPreTime || !foundPostTime) {
114  static G4bool warnedTimesNotFound = false;
115  if (!warnedTimesNotFound) {
116  G4cout <<
117  "*************************************************************************"
118  "\n* WARNING: G4TrajectoryDrawerUtils::GetPointsAndTimes: times not found."
119  "\n You need to specify \"/vis/scene/add/trajectories rich\""
120  "\n*************************************************************************"
121  << G4endl;
122  warnedTimesNotFound = true;
123  }
124  validity = InvalidTimes;
125  }
126  }
127  delete trajectoryPointAttValues; // (Must be deleted after use.)
128  }
129 
130  const std::vector<G4ThreeVector>* auxiliaries
131  = aTrajectoryPoint->GetAuxiliaryPoints();
132  if (0 != auxiliaries) {
133  for (size_t iAux=0; iAux<auxiliaries->size(); ++iAux) {
134  const G4ThreeVector& auxPointPosition = (*auxiliaries)[iAux];
135  if (positions.size() == 0 ||
136  auxPointPosition != positions[positions.size()-1]) {
137  // Only store if first or if different
138  positions.push_back(trajectoryPointPosition);
139  trajectoryLine.push_back(auxPointPosition);
140  auxiliaryPoints.push_back(auxPointPosition);
141  if (validity == ValidTimes) {
142  // Interpolate time for auxiliary points...
143  G4double s1 =
144  (auxPointPosition - lastTrajectoryPointPosition).mag();
145  G4double s2 =
146  (trajectoryPointPosition - auxPointPosition).mag();
147  G4double t = trajectoryPointPreTime +
148  (trajectoryPointPostTime - trajectoryPointPreTime) *
149  (s1 / (s1 + s2));
150  trajectoryLineTimes.push_back(t);
151  auxiliaryPointTimes.push_back(t);
152  }
153  }
154  }
155  }
156 
157  positions.push_back(trajectoryPointPosition);
158  trajectoryLine.push_back(trajectoryPointPosition);
159  stepPoints.push_back(trajectoryPointPosition);
160  if (validity == ValidTimes) {
161  trajectoryLineTimes.push_back(trajectoryPointPostTime);
162  stepPointTimes.push_back(trajectoryPointPostTime);
163  }
164  lastTrajectoryPointPosition = trajectoryPointPosition;
165  }
166  }
167  return validity;
168  }
169 
170  static void SliceLine(G4double timeIncrement,
171  G4Polyline& trajectoryLine,
172  std::vector<G4double>& trajectoryLineTimes)
173  {
174  // Assumes valid arguments from GetPoints and GetTimes.
175 
176  G4Polyline newTrajectoryLine;
177  std::vector<G4double> newTrajectoryLineTimes;
178 
179  newTrajectoryLine.push_back(trajectoryLine[0]);
180  newTrajectoryLineTimes.push_back(trajectoryLineTimes[0]);
181  size_t lineSize = trajectoryLine.size();
182  if (lineSize > 1) {
183  for (size_t i = 1; i < trajectoryLine.size(); ++i) {
184  G4double deltaT = trajectoryLineTimes[i] - trajectoryLineTimes[i - 1];
185  if (deltaT > 0.) {
186  G4double practicalTimeIncrement =
187  std::max(timeIncrement, deltaT / 100.);
188  for (G4double t =
189  (int(trajectoryLineTimes[i - 1]/practicalTimeIncrement) + 1) *
190  practicalTimeIncrement;
191  t <= trajectoryLineTimes[i];
192  t += practicalTimeIncrement) {
193  G4ThreeVector pos = trajectoryLine[i - 1] +
194  (trajectoryLine[i] - trajectoryLine[i - 1]) *
195  ((t - trajectoryLineTimes[i - 1]) / deltaT);
196  newTrajectoryLine.push_back(pos);
197  newTrajectoryLineTimes.push_back(t);
198  }
199  }
200  newTrajectoryLine.push_back(trajectoryLine[i]);
201  newTrajectoryLineTimes.push_back(trajectoryLineTimes[i]);
202  }
203  }
204 
205  trajectoryLine = newTrajectoryLine;
206  trajectoryLineTimes = newTrajectoryLineTimes;
207  }
208 
209  static void DrawWithoutTime(const G4VisTrajContext& myContext,
210  G4Polyline& trajectoryLine,
211  G4Polymarker& auxiliaryPoints,
212  G4Polymarker& stepPoints)
213  {
214  // Draw without time slice information
215 
217  if (0 == pVVisManager) return;
218 
219  if (myContext.GetDrawLine()) {
220  G4VisAttributes trajectoryLineAttribs(myContext.GetLineColour());
221  trajectoryLineAttribs.SetVisibility(myContext.GetLineVisible());
222  trajectoryLine.SetVisAttributes(&trajectoryLineAttribs);
223 
224  pVVisManager->Draw(trajectoryLine);
225  }
226 
227  if (myContext.GetDrawAuxPts() && (auxiliaryPoints.size() > 0)) {
228  auxiliaryPoints.SetMarkerType(myContext.GetAuxPtsType());
229  auxiliaryPoints.SetSize(myContext.GetAuxPtsSizeType(), myContext.GetAuxPtsSize());
230  auxiliaryPoints.SetFillStyle(myContext.GetAuxPtsFillStyle());
231 
232  G4VisAttributes auxiliaryPointsAttribs(myContext.GetAuxPtsColour());
233  auxiliaryPointsAttribs.SetVisibility(myContext.GetAuxPtsVisible());
234  auxiliaryPoints.SetVisAttributes(&auxiliaryPointsAttribs);
235 
236  pVVisManager->Draw(auxiliaryPoints);
237  }
238 
239  if (myContext.GetDrawStepPts() && (stepPoints.size() > 0)) {
240  stepPoints.SetMarkerType(myContext.GetStepPtsType());
241  stepPoints.SetSize(myContext.GetStepPtsSizeType(), myContext.GetStepPtsSize());
242  stepPoints.SetFillStyle(myContext.GetStepPtsFillStyle());
243 
244  G4VisAttributes stepPointsAttribs(myContext.GetStepPtsColour());
245  stepPointsAttribs.SetVisibility(myContext.GetStepPtsVisible());
246  stepPoints.SetVisAttributes(&stepPointsAttribs);
247 
248  pVVisManager->Draw(stepPoints);
249  }
250  }
251 
252  static void DrawWithTime(const G4VisTrajContext& myContext,
253  G4Polyline& trajectoryLine,
254  G4Polymarker& auxiliaryPoints,
255  G4Polymarker& stepPoints,
256  std::vector<G4double>& trajectoryLineTimes,
257  std::vector<G4double>& auxiliaryPointTimes,
258  std::vector<G4double>& stepPointTimes)
259  {
260  // Draw with time slice information
261 
263  if (0 == pVVisManager) return;
264 
265  if (myContext.GetDrawLine()) {
266  G4VisAttributes trajectoryLineAttribs(myContext.GetLineColour());
267  trajectoryLineAttribs.SetVisibility(myContext.GetLineVisible());
268 
269  for (size_t i = 1; i < trajectoryLine.size(); ++i ) {
270  G4Polyline slice;
271  slice.push_back(trajectoryLine[i -1]);
272  slice.push_back(trajectoryLine[i]);
273  trajectoryLineAttribs.SetStartTime(trajectoryLineTimes[i - 1]);
274  trajectoryLineAttribs.SetEndTime(trajectoryLineTimes[i]);
275  slice.SetVisAttributes(&trajectoryLineAttribs);
276  pVVisManager->Draw(slice);
277  }
278  }
279 
280  if (myContext.GetDrawAuxPts() && (auxiliaryPoints.size() > 0)) {
281  G4VisAttributes auxiliaryPointsAttribs(myContext.GetAuxPtsColour());
282  auxiliaryPointsAttribs.SetVisibility(myContext.GetAuxPtsVisible());
283 
284  for (size_t i = 0; i < auxiliaryPoints.size(); ++i ) {
285  G4Polymarker point;
286  point.push_back(auxiliaryPoints[i]);
287  point.SetMarkerType(myContext.GetAuxPtsType());
288  point.SetSize(myContext.GetAuxPtsSizeType(), myContext.GetAuxPtsSize());
289  point.SetFillStyle(myContext.GetAuxPtsFillStyle());
290  auxiliaryPointsAttribs.SetStartTime(auxiliaryPointTimes[i]);
291  auxiliaryPointsAttribs.SetEndTime(auxiliaryPointTimes[i]);
292  point.SetVisAttributes(&auxiliaryPointsAttribs);
293  pVVisManager->Draw(point);
294  }
295  }
296 
297  if (myContext.GetDrawStepPts() && (stepPoints.size() > 0)) {
298  G4VisAttributes stepPointsAttribs(myContext.GetStepPtsColour());
299  stepPointsAttribs.SetVisibility(myContext.GetStepPtsVisible());
300 
301  for (size_t i = 0; i < stepPoints.size(); ++i ) {
302  G4Polymarker point;
303  point.push_back(stepPoints[i]);
304  point.SetMarkerType(myContext.GetStepPtsType());
305  point.SetSize(myContext.GetStepPtsSizeType(), myContext.GetStepPtsSize());
306  point.SetFillStyle(myContext.GetStepPtsFillStyle());
307  stepPointsAttribs.SetStartTime(stepPointTimes[i]);
308  stepPointsAttribs.SetEndTime(stepPointTimes[i]);
309  point.SetVisAttributes(&stepPointsAttribs);
310  pVVisManager->Draw(point);
311  }
312  }
313  }
314 
315  void DrawLineAndPoints(const G4VTrajectory& traj, const G4VisTrajContext& context)
316  {
317  // Return if don't need to do anything
318  if (!context.GetDrawLine() && !context.GetDrawAuxPts() && !context.GetDrawStepPts()) return;
319 
320  // Get points and times (times are returned only if time-slicing
321  // is requested).
322  G4Polyline trajectoryLine;
323  G4Polymarker stepPoints;
324  G4Polymarker auxiliaryPoints;
325  std::vector<G4double> trajectoryLineTimes;
326  std::vector<G4double> stepPointTimes;
327  std::vector<G4double> auxiliaryPointTimes;
328 
330  (traj, context,
331  trajectoryLine, auxiliaryPoints, stepPoints,
332  trajectoryLineTimes, auxiliaryPointTimes, stepPointTimes);
333 
334  if (validity == ValidTimes) {
335 
337  trajectoryLine, trajectoryLineTimes);
338 
339  DrawWithTime(context,
340  trajectoryLine, auxiliaryPoints, stepPoints,
341  trajectoryLineTimes, auxiliaryPointTimes, stepPointTimes);
342 
343  } else {
344 
345  DrawWithoutTime(context, trajectoryLine, auxiliaryPoints, stepPoints);
346 
347  }
348  }
349 }
G4bool GetStepPtsVisible() const
void DrawLineAndPoints(const G4VTrajectory &traj, const G4VisTrajContext &)
static void DrawWithTime(const G4VisTrajContext &myContext, G4Polyline &trajectoryLine, G4Polymarker &auxiliaryPoints, G4Polymarker &stepPoints, std::vector< G4double > &trajectoryLineTimes, std::vector< G4double > &auxiliaryPointTimes, std::vector< G4double > &stepPointTimes)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
virtual const G4ThreeVector GetPosition() const =0
virtual void Draw(const G4Circle &, const G4Transform3D &objectTransformation=G4Transform3D())=0
static const G4double pos
G4bool GetDrawLine() const
G4VMarker::SizeType GetAuxPtsSizeType() const
G4bool GetLineVisible() const
G4double GetAuxPtsSize() const
#define G4endl
Definition: G4ios.hh:61
G4double GetStepPtsSize() const
G4Polymarker::MarkerType GetStepPtsType() const
TimesValidity GetPointsAndTimes(const G4VTrajectory &traj, const G4VisTrajContext &context, G4Polyline &trajectoryLine, G4Polymarker &auxiliaryPoints, G4Polymarker &stepPoints, std::vector< G4double > &trajectoryLineTimes, std::vector< G4double > &auxiliaryPointTimes, std::vector< G4double > &stepPointTimes)
G4Polymarker::MarkerType GetAuxPtsType() const
static G4VVisManager * GetConcreteInstance()
virtual std::vector< G4AttValue > * CreateAttValues() const
G4double GetTimeSliceInterval() const
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
G4VMarker::FillStyle GetStepPtsFillStyle() const
virtual int GetPointEntries() const =0
G4VMarker::SizeType GetStepPtsSizeType() const
static void DrawWithoutTime(const G4VisTrajContext &myContext, G4Polyline &trajectoryLine, G4Polymarker &auxiliaryPoints, G4Polymarker &stepPoints)
void SetSize(SizeType, G4double)
Definition: G4VMarker.cc:119
G4bool GetDrawAuxPts() const
G4bool GetAuxPtsVisible() const
G4Colour GetAuxPtsColour() const
static void SliceLine(G4double timeIncrement, G4Polyline &trajectoryLine, std::vector< G4double > &trajectoryLineTimes)
void SetVisibility(G4bool=true)
const XML_Char * context
Definition: expat.h:434
int G4int
Definition: G4Types.hh:78
G4Colour GetStepPtsColour() const
G4VMarker::FillStyle GetAuxPtsFillStyle() const
virtual const std::vector< G4ThreeVector > * GetAuxiliaryPoints() const
G4Colour GetLineColour() const
G4GLOB_DLL std::ostream G4cout
void SetMarkerType(MarkerType)
G4bool GetDrawStepPts() const
static G4double ConvertToDimensionedDouble(const char *st)
Definition: G4UIcommand.cc:465
virtual G4VTrajectoryPoint * GetPoint(G4int i) const =0
void SetVisAttributes(const G4VisAttributes *)
Definition: G4Visible.cc:80
void SetFillStyle(FillStyle)