Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4TrajectoryDrawerUtils Namespace Reference

Enumerations

enum  TimesValidity { InvalidTimes, ValidTimes }
 

Functions

void GetPoints (const G4VTrajectory &traj, G4Polyline &trajectoryLine, G4Polymarker &auxiliaryPoints, G4Polymarker &stepPoints)
 
void DrawLineAndPoints (const G4VTrajectory &traj, const G4VisTrajContext &)
 
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)
 
static void SliceLine (G4double timeIncrement, G4Polyline &trajectoryLine, std::vector< G4double > &trajectoryLineTimes)
 
static void DrawWithoutTime (const G4VisTrajContext &myContext, G4Polyline &trajectoryLine, G4Polymarker &auxiliaryPoints, G4Polymarker &stepPoints)
 
static void DrawWithTime (const G4VisTrajContext &myContext, G4Polyline &trajectoryLine, G4Polymarker &auxiliaryPoints, G4Polymarker &stepPoints, std::vector< G4double > &trajectoryLineTimes, std::vector< G4double > &auxiliaryPointTimes, std::vector< G4double > &stepPointTimes)
 

Enumeration Type Documentation

Function Documentation

void G4TrajectoryDrawerUtils::DrawLineAndPoints ( const G4VTrajectory traj,
const G4VisTrajContext context 
)

Definition at line 314 of file G4TrajectoryDrawerUtils.cc.

315  {
316  // Return if don't need to do anything
317  if (!context.GetDrawLine() && !context.GetDrawAuxPts() && !context.GetDrawStepPts()) return;
318 
319  // Get points and times (times are returned only if time-slicing
320  // is requested).
321  G4Polyline trajectoryLine;
322  G4Polymarker stepPoints;
323  G4Polymarker auxiliaryPoints;
324  std::vector<G4double> trajectoryLineTimes;
325  std::vector<G4double> stepPointTimes;
326  std::vector<G4double> auxiliaryPointTimes;
327 
329  (traj, context,
330  trajectoryLine, auxiliaryPoints, stepPoints,
331  trajectoryLineTimes, auxiliaryPointTimes, stepPointTimes);
332 
333  if (validity == ValidTimes) {
334 
336  trajectoryLine, trajectoryLineTimes);
337 
338  DrawWithTime(context,
339  trajectoryLine, auxiliaryPoints, stepPoints,
340  trajectoryLineTimes, auxiliaryPointTimes, stepPointTimes);
341 
342  } else {
343 
344  DrawWithoutTime(context, trajectoryLine, auxiliaryPoints, stepPoints);
345 
346  }
347  }
G4bool GetDrawLine() const
G4double GetTimeSliceInterval() const
static void SliceLine(G4double timeIncrement, G4Polyline &trajectoryLine, std::vector< G4double > &trajectoryLineTimes)
static void DrawWithTime(const G4VisTrajContext &myContext, G4Polyline &trajectoryLine, G4Polymarker &auxiliaryPoints, G4Polymarker &stepPoints, std::vector< G4double > &trajectoryLineTimes, std::vector< G4double > &auxiliaryPointTimes, std::vector< G4double > &stepPointTimes)
G4bool GetDrawStepPts() const
static void DrawWithoutTime(const G4VisTrajContext &myContext, G4Polyline &trajectoryLine, G4Polymarker &auxiliaryPoints, G4Polymarker &stepPoints)
G4bool GetDrawAuxPts() 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)

Here is the call graph for this function:

Here is the caller graph for this function:

static void G4TrajectoryDrawerUtils::DrawWithoutTime ( const G4VisTrajContext myContext,
G4Polyline trajectoryLine,
G4Polymarker auxiliaryPoints,
G4Polymarker stepPoints 
)
static

Definition at line 208 of file G4TrajectoryDrawerUtils.cc.

212  {
213  // Draw without time slice information
214 
216  if (0 == pVVisManager) return;
217 
218  if (myContext.GetDrawLine()) {
219  G4VisAttributes trajectoryLineAttribs(myContext.GetLineColour());
220  trajectoryLineAttribs.SetVisibility(myContext.GetLineVisible());
221  trajectoryLine.SetVisAttributes(&trajectoryLineAttribs);
222 
223  pVVisManager->Draw(trajectoryLine);
224  }
225 
226  if (myContext.GetDrawAuxPts() && (auxiliaryPoints.size() > 0)) {
227  auxiliaryPoints.SetMarkerType(myContext.GetAuxPtsType());
228  auxiliaryPoints.SetSize(myContext.GetAuxPtsSizeType(), myContext.GetAuxPtsSize());
229  auxiliaryPoints.SetFillStyle(myContext.GetAuxPtsFillStyle());
230 
231  G4VisAttributes auxiliaryPointsAttribs(myContext.GetAuxPtsColour());
232  auxiliaryPointsAttribs.SetVisibility(myContext.GetAuxPtsVisible());
233  auxiliaryPoints.SetVisAttributes(&auxiliaryPointsAttribs);
234 
235  pVVisManager->Draw(auxiliaryPoints);
236  }
237 
238  if (myContext.GetDrawStepPts() && (stepPoints.size() > 0)) {
239  stepPoints.SetMarkerType(myContext.GetStepPtsType());
240  stepPoints.SetSize(myContext.GetStepPtsSizeType(), myContext.GetStepPtsSize());
241  stepPoints.SetFillStyle(myContext.GetStepPtsFillStyle());
242 
243  G4VisAttributes stepPointsAttribs(myContext.GetStepPtsColour());
244  stepPointsAttribs.SetVisibility(myContext.GetStepPtsVisible());
245  stepPoints.SetVisAttributes(&stepPointsAttribs);
246 
247  pVVisManager->Draw(stepPoints);
248  }
249  }
virtual void Draw(const G4Circle &, const G4Transform3D &objectTransformation=G4Transform3D())=0
G4VMarker::FillStyle GetAuxPtsFillStyle() const
G4double GetStepPtsSize() const
G4VMarker::FillStyle GetStepPtsFillStyle() const
G4bool GetDrawLine() const
static G4VVisManager * GetConcreteInstance()
G4double GetAuxPtsSize() const
void SetMarkerType(MarkerType)
void SetFillStyle(FillStyle)
G4Colour GetAuxPtsColour() const
G4Polymarker::MarkerType GetAuxPtsType() const
G4VMarker::SizeType GetAuxPtsSizeType() const
G4bool GetAuxPtsVisible() const
void SetVisibility(G4bool=true)
G4bool GetLineVisible() const
G4VMarker::SizeType GetStepPtsSizeType() const
void SetSize(SizeType, G4double)
Definition: G4VMarker.cc:119
G4bool GetDrawStepPts() const
void SetVisAttributes(const G4VisAttributes *)
Definition: G4Visible.cc:80
G4Colour GetLineColour() const
G4Colour GetStepPtsColour() const
G4Polymarker::MarkerType GetStepPtsType() const
G4bool GetDrawAuxPts() const
G4bool GetStepPtsVisible() const

Here is the call graph for this function:

Here is the caller graph for this function:

static void G4TrajectoryDrawerUtils::DrawWithTime ( const G4VisTrajContext myContext,
G4Polyline trajectoryLine,
G4Polymarker auxiliaryPoints,
G4Polymarker stepPoints,
std::vector< G4double > &  trajectoryLineTimes,
std::vector< G4double > &  auxiliaryPointTimes,
std::vector< G4double > &  stepPointTimes 
)
static

Definition at line 251 of file G4TrajectoryDrawerUtils.cc.

258  {
259  // Draw with time slice information
260 
262  if (0 == pVVisManager) return;
263 
264  if (myContext.GetDrawLine()) {
265  G4VisAttributes trajectoryLineAttribs(myContext.GetLineColour());
266  trajectoryLineAttribs.SetVisibility(myContext.GetLineVisible());
267 
268  for (size_t i = 1; i < trajectoryLine.size(); ++i ) {
269  G4Polyline slice;
270  slice.push_back(trajectoryLine[i -1]);
271  slice.push_back(trajectoryLine[i]);
272  trajectoryLineAttribs.SetStartTime(trajectoryLineTimes[i - 1]);
273  trajectoryLineAttribs.SetEndTime(trajectoryLineTimes[i]);
274  slice.SetVisAttributes(&trajectoryLineAttribs);
275  pVVisManager->Draw(slice);
276  }
277  }
278 
279  if (myContext.GetDrawAuxPts() && (auxiliaryPoints.size() > 0)) {
280  G4VisAttributes auxiliaryPointsAttribs(myContext.GetAuxPtsColour());
281  auxiliaryPointsAttribs.SetVisibility(myContext.GetAuxPtsVisible());
282 
283  for (size_t i = 0; i < auxiliaryPoints.size(); ++i ) {
284  G4Polymarker point;
285  point.push_back(auxiliaryPoints[i]);
286  point.SetMarkerType(myContext.GetAuxPtsType());
287  point.SetSize(myContext.GetAuxPtsSizeType(), myContext.GetAuxPtsSize());
288  point.SetFillStyle(myContext.GetAuxPtsFillStyle());
289  auxiliaryPointsAttribs.SetStartTime(auxiliaryPointTimes[i]);
290  auxiliaryPointsAttribs.SetEndTime(auxiliaryPointTimes[i]);
291  point.SetVisAttributes(&auxiliaryPointsAttribs);
292  pVVisManager->Draw(point);
293  }
294  }
295 
296  if (myContext.GetDrawStepPts() && (stepPoints.size() > 0)) {
297  G4VisAttributes stepPointsAttribs(myContext.GetStepPtsColour());
298  stepPointsAttribs.SetVisibility(myContext.GetStepPtsVisible());
299 
300  for (size_t i = 0; i < stepPoints.size(); ++i ) {
301  G4Polymarker point;
302  point.push_back(stepPoints[i]);
303  point.SetMarkerType(myContext.GetStepPtsType());
304  point.SetSize(myContext.GetStepPtsSizeType(), myContext.GetStepPtsSize());
305  point.SetFillStyle(myContext.GetStepPtsFillStyle());
306  stepPointsAttribs.SetStartTime(stepPointTimes[i]);
307  stepPointsAttribs.SetEndTime(stepPointTimes[i]);
308  point.SetVisAttributes(&stepPointsAttribs);
309  pVVisManager->Draw(point);
310  }
311  }
312  }
virtual void Draw(const G4Circle &, const G4Transform3D &objectTransformation=G4Transform3D())=0
G4VMarker::FillStyle GetAuxPtsFillStyle() const
G4double GetStepPtsSize() const
G4VMarker::FillStyle GetStepPtsFillStyle() const
G4bool GetDrawLine() const
static G4VVisManager * GetConcreteInstance()
G4double GetAuxPtsSize() const
void SetMarkerType(MarkerType)
void SetFillStyle(FillStyle)
G4Colour GetAuxPtsColour() const
G4Polymarker::MarkerType GetAuxPtsType() const
G4VMarker::SizeType GetAuxPtsSizeType() const
G4bool GetAuxPtsVisible() const
void SetVisibility(G4bool=true)
G4bool GetLineVisible() const
G4VMarker::SizeType GetStepPtsSizeType() const
void SetSize(SizeType, G4double)
Definition: G4VMarker.cc:119
G4bool GetDrawStepPts() const
void SetVisAttributes(const G4VisAttributes *)
Definition: G4Visible.cc:80
G4Colour GetLineColour() const
G4Colour GetStepPtsColour() const
G4Polymarker::MarkerType GetStepPtsType() const
G4bool GetDrawAuxPts() const
G4bool GetStepPtsVisible() const

Here is the call graph for this function:

Here is the caller graph for this function:

void G4TrajectoryDrawerUtils::GetPoints ( const G4VTrajectory traj,
G4Polyline trajectoryLine,
G4Polymarker auxiliaryPoints,
G4Polymarker stepPoints 
)
TimesValidity G4TrajectoryDrawerUtils::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 
)

Definition at line 47 of file G4TrajectoryDrawerUtils.cc.

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*************************************************************************"
120  << G4endl;
121  warnedTimesNotFound = true;
122  }
123  validity = InvalidTimes;
124  }
125  }
126  delete trajectoryPointAttValues; // (Must be deleted after use.)
127  }
128 
129  const std::vector<G4ThreeVector>* auxiliaries
130  = aTrajectoryPoint->GetAuxiliaryPoints();
131  if (0 != auxiliaries) {
132  for (size_t iAux=0; iAux<auxiliaries->size(); ++iAux) {
133  const G4ThreeVector& auxPointPosition = (*auxiliaries)[iAux];
134  if (positions.size() == 0 ||
135  auxPointPosition != positions[positions.size()-1]) {
136  // Only store if first or if different
137  positions.push_back(trajectoryPointPosition);
138  trajectoryLine.push_back(auxPointPosition);
139  auxiliaryPoints.push_back(auxPointPosition);
140  if (validity == ValidTimes) {
141  // Interpolate time for auxiliary points...
142  G4double s1 =
143  (auxPointPosition - lastTrajectoryPointPosition).mag();
144  G4double s2 =
145  (trajectoryPointPosition - auxPointPosition).mag();
146  G4double t = trajectoryPointPreTime +
147  (trajectoryPointPostTime - trajectoryPointPreTime) *
148  (s1 / (s1 + s2));
149  trajectoryLineTimes.push_back(t);
150  auxiliaryPointTimes.push_back(t);
151  }
152  }
153  }
154  }
155 
156  positions.push_back(trajectoryPointPosition);
157  trajectoryLine.push_back(trajectoryPointPosition);
158  stepPoints.push_back(trajectoryPointPosition);
159  if (validity == ValidTimes) {
160  trajectoryLineTimes.push_back(trajectoryPointPostTime);
161  stepPointTimes.push_back(trajectoryPointPostTime);
162  }
163  lastTrajectoryPointPosition = trajectoryPointPosition;
164  }
165  }
166  return validity;
167  }
virtual const std::vector< G4ThreeVector > * GetAuxiliaryPoints() const
virtual G4VTrajectoryPoint * GetPoint(G4int i) const =0
G4double GetTimeSliceInterval() const
int G4int
Definition: G4Types.hh:78
static G4double ConvertToDimensionedDouble(const char *st)
Definition: G4UIcommand.cc:463
virtual std::vector< G4AttValue > * CreateAttValues() const
virtual int GetPointEntries() const =0
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
virtual const G4ThreeVector GetPosition() const =0
T max(const T t1, const T t2)
brief Return the largest of the two arguments
#define G4endl
Definition: G4ios.hh:61
double G4double
Definition: G4Types.hh:76

Here is the call graph for this function:

Here is the caller graph for this function:

static void G4TrajectoryDrawerUtils::SliceLine ( G4double  timeIncrement,
G4Polyline trajectoryLine,
std::vector< G4double > &  trajectoryLineTimes 
)
static

Definition at line 169 of file G4TrajectoryDrawerUtils.cc.

172  {
173  // Assumes valid arguments from GetPoints and GetTimes.
174 
175  G4Polyline newTrajectoryLine;
176  std::vector<G4double> newTrajectoryLineTimes;
177 
178  newTrajectoryLine.push_back(trajectoryLine[0]);
179  newTrajectoryLineTimes.push_back(trajectoryLineTimes[0]);
180  size_t lineSize = trajectoryLine.size();
181  if (lineSize > 1) {
182  for (size_t i = 1; i < trajectoryLine.size(); ++i) {
183  G4double deltaT = trajectoryLineTimes[i] - trajectoryLineTimes[i - 1];
184  if (deltaT > 0.) {
185  G4double practicalTimeIncrement =
186  std::max(timeIncrement, deltaT / 100.);
187  for (G4double t =
188  (int(trajectoryLineTimes[i - 1]/practicalTimeIncrement) + 1) *
189  practicalTimeIncrement;
190  t <= trajectoryLineTimes[i];
191  t += practicalTimeIncrement) {
192  G4ThreeVector pos = trajectoryLine[i - 1] +
193  (trajectoryLine[i] - trajectoryLine[i - 1]) *
194  ((t - trajectoryLineTimes[i - 1]) / deltaT);
195  newTrajectoryLine.push_back(pos);
196  newTrajectoryLineTimes.push_back(t);
197  }
198  }
199  newTrajectoryLine.push_back(trajectoryLine[i]);
200  newTrajectoryLineTimes.push_back(trajectoryLineTimes[i]);
201  }
202  }
203 
204  trajectoryLine = newTrajectoryLine;
205  trajectoryLineTimes = newTrajectoryLineTimes;
206  }
T max(const T t1, const T t2)
brief Return the largest of the two arguments
double G4double
Definition: G4Types.hh:76
static const G4double pos

Here is the call graph for this function:

Here is the caller graph for this function: