00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030 #include "G4TrajectoryDrawerUtils.hh"
00031 #include "G4Colour.hh"
00032 #include "G4Polyline.hh"
00033 #include "G4Polymarker.hh"
00034 #include "G4VTrajectory.hh"
00035 #include "G4VTrajectoryPoint.hh"
00036 #include "G4VisAttributes.hh"
00037 #include "G4VisTrajContext.hh"
00038 #include "G4VVisManager.hh"
00039 #include "G4UIcommand.hh"
00040 #include "G4AttValue.hh"
00041
00042 namespace G4TrajectoryDrawerUtils {
00043
00044 enum TimesValidity {InvalidTimes, ValidTimes};
00045
00046 TimesValidity GetPointsAndTimes
00047 (const G4VTrajectory& traj,
00048 const G4VisTrajContext& context,
00049 G4Polyline& trajectoryLine,
00050 G4Polymarker& auxiliaryPoints,
00051 G4Polymarker& stepPoints,
00052 std::vector<G4double>& trajectoryLineTimes,
00053 std::vector<G4double>& auxiliaryPointTimes,
00054 std::vector<G4double>& stepPointTimes)
00055 {
00056 TimesValidity validity = InvalidTimes;
00057 if (context.GetTimeSliceInterval()) validity = ValidTimes;
00058
00059
00060
00061
00062
00063 G4ThreeVector lastTrajectoryPointPosition;
00064
00065
00066 std::vector<G4ThreeVector> positions;
00067
00068 for (G4int iPoint=0; iPoint<traj.GetPointEntries(); iPoint++) {
00069
00070 G4VTrajectoryPoint* aTrajectoryPoint = traj.GetPoint(iPoint);
00071 const G4ThreeVector& trajectoryPointPosition =
00072 aTrajectoryPoint->GetPosition();
00073
00074
00075 if (positions.size() == 0 ||
00076 trajectoryPointPosition != positions[positions.size()-1]) {
00077
00078
00079 G4double trajectoryPointPreTime = -std::numeric_limits<double>::max();
00080 G4double trajectoryPointPostTime = std::numeric_limits<double>::max();
00081
00082 if (context.GetTimeSliceInterval() && validity == ValidTimes) {
00083
00084 std::vector<G4AttValue>* trajectoryPointAttValues =
00085 aTrajectoryPoint->CreateAttValues();
00086 if (!trajectoryPointAttValues) {
00087 static G4bool warnedNoAttValues = false;
00088 if (!warnedNoAttValues) {
00089 G4cout <<
00090 "*************************************************************************"
00091 "\n* WARNING: G4TrajectoryDrawerUtils::GetPointsAndTimes: no att values."
00092 "\n*************************************************************************"
00093 << G4endl;
00094 warnedNoAttValues = true;
00095 }
00096 validity = InvalidTimes;
00097 } else {
00098 G4bool foundPreTime = false, foundPostTime = false;
00099 for (std::vector<G4AttValue>::iterator i =
00100 trajectoryPointAttValues->begin();
00101 i != trajectoryPointAttValues->end(); ++i) {
00102 if (i->GetName() == "PreT") {
00103 trajectoryPointPreTime =
00104 G4UIcommand::ConvertToDimensionedDouble(i->GetValue());
00105 foundPreTime = true;
00106 }
00107 if (i->GetName() == "PostT") {
00108 trajectoryPointPostTime =
00109 G4UIcommand::ConvertToDimensionedDouble(i->GetValue());
00110 foundPostTime = true;
00111 }
00112 }
00113 if (!foundPreTime || !foundPostTime) {
00114 static G4bool warnedTimesNotFound = false;
00115 if (!warnedTimesNotFound) {
00116 G4cout <<
00117 "*************************************************************************"
00118 "\n* WARNING: G4TrajectoryDrawerUtils::GetPointsAndTimes: times not found."
00119 "\n*************************************************************************"
00120 << G4endl;
00121 warnedTimesNotFound = true;
00122 }
00123 validity = InvalidTimes;
00124 }
00125 }
00126 }
00127
00128 const std::vector<G4ThreeVector>* auxiliaries
00129 = aTrajectoryPoint->GetAuxiliaryPoints();
00130 if (0 != auxiliaries) {
00131 for (size_t iAux=0; iAux<auxiliaries->size(); ++iAux) {
00132 const G4ThreeVector& auxPointPosition = (*auxiliaries)[iAux];
00133 if (positions.size() == 0 ||
00134 auxPointPosition != positions[positions.size()-1]) {
00135
00136 positions.push_back(trajectoryPointPosition);
00137 trajectoryLine.push_back(auxPointPosition);
00138 auxiliaryPoints.push_back(auxPointPosition);
00139 if (validity == ValidTimes) {
00140
00141 G4double s1 =
00142 (auxPointPosition - lastTrajectoryPointPosition).mag();
00143 G4double s2 =
00144 (trajectoryPointPosition - auxPointPosition).mag();
00145 G4double t = trajectoryPointPreTime +
00146 (trajectoryPointPostTime - trajectoryPointPreTime) *
00147 (s1 / (s1 + s2));
00148 trajectoryLineTimes.push_back(t);
00149 auxiliaryPointTimes.push_back(t);
00150 }
00151 }
00152 }
00153 }
00154
00155 positions.push_back(trajectoryPointPosition);
00156 trajectoryLine.push_back(trajectoryPointPosition);
00157 stepPoints.push_back(trajectoryPointPosition);
00158 if (validity == ValidTimes) {
00159 trajectoryLineTimes.push_back(trajectoryPointPostTime);
00160 stepPointTimes.push_back(trajectoryPointPostTime);
00161 }
00162 lastTrajectoryPointPosition = trajectoryPointPosition;
00163 }
00164 }
00165 return validity;
00166 }
00167
00168 static void SliceLine(G4double timeIncrement,
00169 G4Polyline& trajectoryLine,
00170 std::vector<G4double>& trajectoryLineTimes)
00171 {
00172
00173
00174 G4Polyline newTrajectoryLine;
00175 std::vector<G4double> newTrajectoryLineTimes;
00176
00177 newTrajectoryLine.push_back(trajectoryLine[0]);
00178 newTrajectoryLineTimes.push_back(trajectoryLineTimes[0]);
00179 size_t lineSize = trajectoryLine.size();
00180 if (lineSize > 1) {
00181 for (size_t i = 1; i < trajectoryLine.size(); ++i) {
00182 G4double deltaT = trajectoryLineTimes[i] - trajectoryLineTimes[i - 1];
00183 if (deltaT > 0.) {
00184 G4double practicalTimeIncrement =
00185 std::max(timeIncrement, deltaT / 100.);
00186 for (G4double t =
00187 (int(trajectoryLineTimes[i - 1]/practicalTimeIncrement) + 1) *
00188 practicalTimeIncrement;
00189 t <= trajectoryLineTimes[i];
00190 t += practicalTimeIncrement) {
00191 G4ThreeVector pos = trajectoryLine[i - 1] +
00192 (trajectoryLine[i] - trajectoryLine[i - 1]) *
00193 ((t - trajectoryLineTimes[i - 1]) / deltaT);
00194 newTrajectoryLine.push_back(pos);
00195 newTrajectoryLineTimes.push_back(t);
00196 }
00197 }
00198 newTrajectoryLine.push_back(trajectoryLine[i]);
00199 newTrajectoryLineTimes.push_back(trajectoryLineTimes[i]);
00200 }
00201 }
00202
00203 trajectoryLine = newTrajectoryLine;
00204 trajectoryLineTimes = newTrajectoryLineTimes;
00205 }
00206
00207 static void DrawWithoutTime(const G4VisTrajContext& myContext,
00208 G4Polyline& trajectoryLine,
00209 G4Polymarker& auxiliaryPoints,
00210 G4Polymarker& stepPoints)
00211 {
00212
00213
00214 G4VVisManager* pVVisManager = G4VVisManager::GetConcreteInstance();
00215 if (0 == pVVisManager) return;
00216
00217 if (myContext.GetDrawLine()) {
00218 G4VisAttributes trajectoryLineAttribs(myContext.GetLineColour());
00219 trajectoryLineAttribs.SetVisibility(myContext.GetLineVisible());
00220 trajectoryLine.SetVisAttributes(&trajectoryLineAttribs);
00221
00222 pVVisManager->Draw(trajectoryLine);
00223 }
00224
00225 if (myContext.GetDrawAuxPts() && (auxiliaryPoints.size() > 0)) {
00226 auxiliaryPoints.SetMarkerType(myContext.GetAuxPtsType());
00227 auxiliaryPoints.SetSize(myContext.GetAuxPtsSizeType(), myContext.GetAuxPtsSize());
00228 auxiliaryPoints.SetFillStyle(myContext.GetAuxPtsFillStyle());
00229
00230 G4VisAttributes auxiliaryPointsAttribs(myContext.GetAuxPtsColour());
00231 auxiliaryPointsAttribs.SetVisibility(myContext.GetAuxPtsVisible());
00232 auxiliaryPoints.SetVisAttributes(&auxiliaryPointsAttribs);
00233
00234 pVVisManager->Draw(auxiliaryPoints);
00235 }
00236
00237 if (myContext.GetDrawStepPts() && (stepPoints.size() > 0)) {
00238 stepPoints.SetMarkerType(myContext.GetStepPtsType());
00239 stepPoints.SetSize(myContext.GetStepPtsSizeType(), myContext.GetStepPtsSize());
00240 stepPoints.SetFillStyle(myContext.GetStepPtsFillStyle());
00241
00242 G4VisAttributes stepPointsAttribs(myContext.GetStepPtsColour());
00243 stepPointsAttribs.SetVisibility(myContext.GetStepPtsVisible());
00244 stepPoints.SetVisAttributes(&stepPointsAttribs);
00245
00246 pVVisManager->Draw(stepPoints);
00247 }
00248 }
00249
00250 static void DrawWithTime(const G4VisTrajContext& myContext,
00251 G4Polyline& trajectoryLine,
00252 G4Polymarker& auxiliaryPoints,
00253 G4Polymarker& stepPoints,
00254 std::vector<G4double>& trajectoryLineTimes,
00255 std::vector<G4double>& auxiliaryPointTimes,
00256 std::vector<G4double>& stepPointTimes)
00257 {
00258
00259
00260 G4VVisManager* pVVisManager = G4VVisManager::GetConcreteInstance();
00261 if (0 == pVVisManager) return;
00262
00263 if (myContext.GetDrawLine()) {
00264 G4VisAttributes trajectoryLineAttribs(myContext.GetLineColour());
00265 trajectoryLineAttribs.SetVisibility(myContext.GetLineVisible());
00266
00267 for (size_t i = 1; i < trajectoryLine.size(); ++i ) {
00268 G4Polyline slice;
00269 slice.push_back(trajectoryLine[i -1]);
00270 slice.push_back(trajectoryLine[i]);
00271 trajectoryLineAttribs.SetStartTime(trajectoryLineTimes[i - 1]);
00272 trajectoryLineAttribs.SetEndTime(trajectoryLineTimes[i]);
00273 slice.SetVisAttributes(&trajectoryLineAttribs);
00274 pVVisManager->Draw(slice);
00275 }
00276 }
00277
00278 if (myContext.GetDrawAuxPts() && (auxiliaryPoints.size() > 0)) {
00279 G4VisAttributes auxiliaryPointsAttribs(myContext.GetAuxPtsColour());
00280 auxiliaryPointsAttribs.SetVisibility(myContext.GetAuxPtsVisible());
00281
00282 for (size_t i = 0; i < auxiliaryPoints.size(); ++i ) {
00283 G4Polymarker point;
00284 point.push_back(auxiliaryPoints[i]);
00285 point.SetMarkerType(myContext.GetAuxPtsType());
00286 point.SetSize(myContext.GetAuxPtsSizeType(), myContext.GetAuxPtsSize());
00287 point.SetFillStyle(myContext.GetAuxPtsFillStyle());
00288 auxiliaryPointsAttribs.SetStartTime(auxiliaryPointTimes[i]);
00289 auxiliaryPointsAttribs.SetEndTime(auxiliaryPointTimes[i]);
00290 point.SetVisAttributes(&auxiliaryPointsAttribs);
00291 pVVisManager->Draw(point);
00292 }
00293 }
00294
00295 if (myContext.GetDrawStepPts() && (stepPoints.size() > 0)) {
00296 G4VisAttributes stepPointsAttribs(myContext.GetStepPtsColour());
00297 stepPointsAttribs.SetVisibility(myContext.GetStepPtsVisible());
00298
00299 for (size_t i = 0; i < stepPoints.size(); ++i ) {
00300 G4Polymarker point;
00301 point.push_back(stepPoints[i]);
00302 point.SetMarkerType(myContext.GetStepPtsType());
00303 point.SetSize(myContext.GetStepPtsSizeType(), myContext.GetStepPtsSize());
00304 point.SetFillStyle(myContext.GetStepPtsFillStyle());
00305 stepPointsAttribs.SetStartTime(stepPointTimes[i]);
00306 stepPointsAttribs.SetEndTime(stepPointTimes[i]);
00307 point.SetVisAttributes(&stepPointsAttribs);
00308 pVVisManager->Draw(point);
00309 }
00310 }
00311 }
00312
00313 void DrawLineAndPoints(const G4VTrajectory& traj, const G4VisTrajContext& context, const G4int& i_mode)
00314 {
00315 static G4bool warnedAboutIMode = false;
00316 G4ExceptionDescription ed;
00317 ed << "WARNING: DEPRECATED use of i_mode (i_mode: " << i_mode
00318 << "). Feature will be removed at a future major release.";
00319 if (!warnedAboutIMode) {
00320 G4Exception
00321 ("G4TrajectoryDrawerUtils::DrawLineAndPoints(traj, context, i_mode)",
00322 "modeling0125", JustWarning, ed);
00323 warnedAboutIMode = true;
00324 }
00325
00326
00327 G4VisTrajContext myContext(context);
00328
00329 if (i_mode != 0) {
00330 const G4double markerSize = std::abs(i_mode)/1000;
00331 G4bool lineRequired (i_mode >= 0);
00332 G4bool markersRequired (markerSize > 0.);
00333
00334 myContext.SetDrawLine(lineRequired);
00335 myContext.SetDrawAuxPts(markersRequired);
00336 myContext.SetDrawStepPts(markersRequired);
00337
00338 myContext.SetAuxPtsSize(markerSize);
00339 myContext.SetStepPtsSize(markerSize);
00340 }
00341
00342
00343 if (!myContext.GetDrawLine() && !myContext.GetDrawAuxPts() && !myContext.GetDrawStepPts()) return;
00344
00345
00346
00347 G4Polyline trajectoryLine;
00348 G4Polymarker stepPoints;
00349 G4Polymarker auxiliaryPoints;
00350 std::vector<G4double> trajectoryLineTimes;
00351 std::vector<G4double> stepPointTimes;
00352 std::vector<G4double> auxiliaryPointTimes;
00353
00354 TimesValidity validity = GetPointsAndTimes
00355 (traj, context,
00356 trajectoryLine, auxiliaryPoints, stepPoints,
00357 trajectoryLineTimes, auxiliaryPointTimes, stepPointTimes);
00358
00359 if (validity == ValidTimes) {
00360
00361 SliceLine(context.GetTimeSliceInterval(),
00362 trajectoryLine, trajectoryLineTimes);
00363
00364 DrawWithTime(context,
00365 trajectoryLine, auxiliaryPoints, stepPoints,
00366 trajectoryLineTimes, auxiliaryPointTimes, stepPointTimes);
00367
00368 } else {
00369
00370 DrawWithoutTime(context, trajectoryLine, auxiliaryPoints, stepPoints);
00371
00372 }
00373 }
00374
00375 void DrawLineAndPoints(const G4VTrajectory& traj, const G4VisTrajContext& context)
00376 {
00377
00378 if (!context.GetDrawLine() && !context.GetDrawAuxPts() && !context.GetDrawStepPts()) return;
00379
00380
00381
00382 G4Polyline trajectoryLine;
00383 G4Polymarker stepPoints;
00384 G4Polymarker auxiliaryPoints;
00385 std::vector<G4double> trajectoryLineTimes;
00386 std::vector<G4double> stepPointTimes;
00387 std::vector<G4double> auxiliaryPointTimes;
00388
00389 TimesValidity validity = GetPointsAndTimes
00390 (traj, context,
00391 trajectoryLine, auxiliaryPoints, stepPoints,
00392 trajectoryLineTimes, auxiliaryPointTimes, stepPointTimes);
00393
00394 if (validity == ValidTimes) {
00395
00396 SliceLine(context.GetTimeSliceInterval(),
00397 trajectoryLine, trajectoryLineTimes);
00398
00399 DrawWithTime(context,
00400 trajectoryLine, auxiliaryPoints, stepPoints,
00401 trajectoryLineTimes, auxiliaryPointTimes, stepPointTimes);
00402
00403 } else {
00404
00405 DrawWithoutTime(context, trajectoryLine, auxiliaryPoints, stepPoints);
00406
00407 }
00408 }
00409 }