xmsgridtrace  1.0
XmGridTrace.cpp
Go to the documentation of this file.
1 //------------------------------------------------------------------------------
6 //------------------------------------------------------------------------------
7 
8 //----- Included files ---------------------------------------------------------
9 
10 // 1. Precompiled header
11 
12 // 2. My own header
14 
15 // 3. Standard library headers
16 #include <sstream>
17 
18 // 4. External library headers
19 
20 // 5. Shared code headers
21 #include <xmscore/math/math.h>
22 #include <xmscore/misc/XmError.h>
23 #include <xmscore/misc/XmLog.h>
24 #include <xmscore/misc/xmstype.h> // XM_ZERO_TOL
28 
29 // 6. Non-shared code headers
30 
31 //----- Forward declarations ---------------------------------------------------
32 
33 //----- External globals -------------------------------------------------------
34 
35 //----- Namespace declaration --------------------------------------------------
36 
37 //----- Constants / Enumerations -----------------------------------------------
38 
39 //----- Classes / Structs ------------------------------------------------------
40 
41 //----- Internal functions -----------------------------------------------------
42 namespace xms
43 {
44 namespace
45 {
47 
48 //----- Class / Function definitions -------------------------------------------
49 
52 class XmGridTraceImpl : public XmGridTrace
53 {
54 public:
55  XmGridTraceImpl(BSHP<XmUGrid> a_ugrid);
56  ~XmGridTraceImpl(){};
57 
58  double GetVectorMultiplier() const final;
59  void SetVectorMultiplier(const double a_vectorMultiplier) final;
60 
61  double GetMaxTracingTime() const final;
62  void SetMaxTracingTime(const double a_maxTracingTime) final;
63 
64  double GetMaxTracingDistance() const final;
65  void SetMaxTracingDistance(const double a_maxTracingDistance) final;
66 
67  double GetMinDeltaTime() const final;
68  void SetMinDeltaTime(const double a_minDeltaTime) final;
69 
70  double GetMaxChangeDistance() const final;
71  void SetMaxChangeDistance(const double a_maxChangeDistance) final;
72 
73  double GetMaxChangeVelocity() const final;
74  void SetMaxChangeVelocity(const double a_maxChangeVelocity) final;
75 
76  double GetMaxChangeDirectionInRadians() const final;
77  void SetMaxChangeDirectionInRadians(const double a_maxChangeDirection) final;
78 
79  void AddGridScalarsAtTime(const VecPt3d& a_scalars,
80  DataLocationEnum a_scalarLoc,
81  xms::DynBitset& a_activity,
82  DataLocationEnum a_activityLoc,
83  double a_time) final;
84 
85  void TracePoint(const Pt3d& a_pt,
86  const double& a_ptTime,
87  VecPt3d& a_outTrace,
88  VecDbl& a_outTimes) final;
89 
90  std::string GetExitMessage() final;
91 
92 private:
93  bool GetVectorAtLocationAndTime(const xms::Pt3d& a_pt,
94  double a_currentTime,
95  xms::Pt3d& a_data) const;
96 
97  BSHP<XmUGrid> m_ugrid;
98  double m_vectorMultiplier=1;
99  double m_maxTracingTime=-1;
101  double m_minDeltaTime=1;
105 
110  double m_time1=-1;
115  double m_time2=-1;
116  double m_distTraveled=0;
117 
118  std::string m_exitMessage;
119 protected:
120 };
121 double iGetDirAsCosTheta(double a_vx0, double a_vy0, double a_vx1, double a_vy1)
122 {
123  // Should be the cosine of the angle
124  double mag0 = sqrt(a_vx0 * a_vx0 + a_vy0 * a_vy0);
125  double mag1 = sqrt(a_vx1 * a_vx1 + a_vy1 * a_vy1);
126  return (a_vx0 * a_vx1 + a_vy0 * a_vy1) / (mag0 * mag1);
127 }
128 //------------------------------------------------------------------------------
131 //------------------------------------------------------------------------------
132 XmGridTraceImpl::XmGridTraceImpl(BSHP<XmUGrid> a_ugrid)
133 : m_ugrid(a_ugrid)
134 {
135 }
136 
141 //------------------------------------------------------------------------------
143 //------------------------------------------------------------------------------
144 double XmGridTraceImpl::GetVectorMultiplier() const
145 {
146  return m_vectorMultiplier;
147 } // XmGridTraceImpl::GetVectorMultiplier
148 //------------------------------------------------------------------------------
151 //------------------------------------------------------------------------------
152 void XmGridTraceImpl::SetVectorMultiplier(const double a_vectorMultiplier)
153 {
154  m_vectorMultiplier = a_vectorMultiplier;
155 } // XmGridTraceImpl::GetVectorMultiplier
156 //------------------------------------------------------------------------------
158 //------------------------------------------------------------------------------
159 double XmGridTraceImpl::GetMaxTracingTime() const
160 {
161  return m_maxTracingTime;
162 } // XmGridTraceImpl::GetMaxTracingTime
163 //------------------------------------------------------------------------------
166 //------------------------------------------------------------------------------
167 void XmGridTraceImpl::SetMaxTracingTime(const double a_maxTracingTime)
168 {
169  m_maxTracingTime = a_maxTracingTime;
170 } // XmGridTraceImpl::SetMaxTracingTime
171 //------------------------------------------------------------------------------
173 //------------------------------------------------------------------------------
174 double XmGridTraceImpl::GetMaxTracingDistance() const
175 {
176  return m_maxTracingDistance;
177 } // XmGridTraceImpl::GetMaxTracingDistance
178 //------------------------------------------------------------------------------
181 //------------------------------------------------------------------------------
182 void XmGridTraceImpl::SetMaxTracingDistance(const double a_maxTracingDistance)
183 {
184  m_maxTracingDistance = a_maxTracingDistance;
185 } // XmGridTraceImpl::SetMaxTracingDistance
186 //------------------------------------------------------------------------------
188 //------------------------------------------------------------------------------
189 double XmGridTraceImpl::GetMinDeltaTime() const
190 {
191  return m_minDeltaTime;
192 } // XmGridTraceImpl::GetMinDeltaTime
193 //------------------------------------------------------------------------------
196 //------------------------------------------------------------------------------
197 void XmGridTraceImpl::SetMinDeltaTime(const double a_minDeltaTime)
198 {
199  if (a_minDeltaTime <= 0) // Must have an exit condition to avoid infinite loops
200  {
202  }
203  else
204  {
205  m_minDeltaTime = a_minDeltaTime;
206  }
207 } // XmGridTraceImpl::SetMinDeltaTime
208 //------------------------------------------------------------------------------
210 //------------------------------------------------------------------------------
211 double XmGridTraceImpl::GetMaxChangeDistance() const
212 {
213  return m_maxChangeDistance;
214 } // XmGridTraceImpl::GetMaxChangeDistance
215 //------------------------------------------------------------------------------
218 //------------------------------------------------------------------------------
219 void XmGridTraceImpl::SetMaxChangeDistance(const double a_maxChangeDistance)
220 {
221  m_maxChangeDistance = a_maxChangeDistance;
222 } // XmGridTraceImpl::SetMaxChangeDistance
223 //------------------------------------------------------------------------------
225 //------------------------------------------------------------------------------
226 double XmGridTraceImpl::GetMaxChangeVelocity() const
227 {
228  return m_maxChangeVelocity;
229 } // XmGridTraceImpl::GetMaxChangeVelocity
230 //------------------------------------------------------------------------------
233 //------------------------------------------------------------------------------
234 void XmGridTraceImpl::SetMaxChangeVelocity(const double a_maxChangeVelocity)
235 {
236  m_maxChangeVelocity = a_maxChangeVelocity;
237 } // XmGridTraceImpl::SetMaxChangeVelocity
238 //------------------------------------------------------------------------------
240 //------------------------------------------------------------------------------
241 double XmGridTraceImpl::GetMaxChangeDirectionInRadians() const
242 {
244 } // XmGridTraceImpl::GetMaxChangeDirectionInRadians
245 //------------------------------------------------------------------------------
248 //------------------------------------------------------------------------------
249 void XmGridTraceImpl::SetMaxChangeDirectionInRadians(const double a_maxChangeDirection)
250 {
251  if (a_maxChangeDirection > XM_PI || a_maxChangeDirection < 0)
252  {
253  }
254  m_maxChangeDirectionInRadians = a_maxChangeDirection;
255 } // XmGridTraceImpl::SetMaxChangeDirectionInRadians
256 //------------------------------------------------------------------------------
258 //------------------------------------------------------------------------------
259 std::string XmGridTraceImpl::GetExitMessage()
260 {
261  return m_exitMessage;
262 } // XmGridTraceImpl::GetExitMessage
263 //------------------------------------------------------------------------------
272 //------------------------------------------------------------------------------
273 void XmGridTraceImpl::AddGridScalarsAtTime(const VecPt3d& a_scalars,
274  DataLocationEnum a_scalarLoc,
275  xms::DynBitset& a_activity,
276  DataLocationEnum a_activityLoc,
277  double a_time)
278 {
280  {
283  m_time1 = m_time2;
284  }
285  m_extractor2x = XmUGrid2dDataExtractor::New(m_ugrid);
286  m_extractor2y = XmUGrid2dDataExtractor::New(m_ugrid);
287 
288  m_time2 = a_time;
289  std::vector<float> xx, yy;
290  for (auto& pt : a_scalars)
291  {
292  xx.push_back((float)pt.x);
293  yy.push_back((float)pt.y);
294  }
295  if (a_scalarLoc == DataLocationEnum::LOC_POINTS)
296  {
297  m_extractor2x->SetGridPointScalars(xx, a_activity, a_activityLoc);
298  m_extractor2y->SetGridPointScalars(yy, a_activity, a_activityLoc);
299  }
300  else
301  {
302  m_extractor2x->SetGridCellScalars(xx, a_activity, a_activityLoc);
303  m_extractor2y->SetGridCellScalars(yy, a_activity, a_activityLoc);
304  }
305 }
306 
307 //------------------------------------------------------------------------------
313 //------------------------------------------------------------------------------
314 void XmGridTraceImpl::TracePoint(const Pt3d& a_pt,
315  const double& a_ptTime,
316  VecPt3d& a_outTrace,
317  VecDbl& a_outTimes)
318 {
319  m_exitMessage.clear();
320  double deltaT = 1.00;
321  double mag0 = 0, mag1 = 0;
322  Pt3d pt0 = a_pt, pt1;
323  double vx0 = 0, vx1 = 0, vy0 = 0, vy1 = 0, elapsedTime = 0;
324  bool bContinue = true;
325  Pt3d vtkVec; // Rename this variable
326  Pt3d vtkPt;
327  Pt3d vector;
328 
329  m_distTraveled = 0;
330  a_outTrace.clear();
331  a_outTimes.clear();
332  if (a_ptTime > m_time2 || // Test if the time specified is after the time range
333  !GetVectorAtLocationAndTime(a_pt, a_ptTime, vector)) // Ensure nothing fails during extraction
334  {
335  m_exitMessage = "Error occurred while extracting point0.";
336  return;
337  }
338  if (EQ_TOL(vector.x, XM_NODATA, 1) || EQ_TOL(vector.y, XM_NODATA, 1))
339  {
340  m_exitMessage = "Point does not start inside an active cell.";
341  return;
342  }
343 
344  a_outTrace.push_back(a_pt);
345  a_outTimes.push_back(a_ptTime);
346 
347  vx0 = vector.x * m_vectorMultiplier;
348  vy0 = vector.y * m_vectorMultiplier;
349  mag0 = sqrt(vector.x * vector.x + vector.y * vector.y);
350 
351  double maxAngleChange = cos(m_maxChangeDirectionInRadians);
352 
353  while (bContinue)
354  {
355  if (m_maxChangeDistance > 0)
356  {
357  // make sure deltaT is small enough to not go past the max dist
359  double denom = (vx0 * vx0) + (vy0 * vy0) + (m_maxChangeDistance * XM_ZERO_TOL);
360  double dt = sqrt(d2 / denom);
361  if (deltaT > dt)
362  {
363  deltaT = dt;
364  m_exitMessage = "Change distance was greater than the max change distance.";
365  }
366  }
367  // If the change in DeltaT would push us beyond the time step, set it to hit the timestep
368  if (elapsedTime + deltaT + a_ptTime > m_time2)
369  {
370  deltaT = m_time2 - elapsedTime - a_ptTime;
371  bContinue = false; // This will be the last point traced
372  m_exitMessage = "The point has traveled beyond, or reached the second time step.";
373  }
374  // If the change in delta time would push beyond the max tracing time, set it to hit max tracing
375  // time
376  if (m_maxTracingTime > 0 && (elapsedTime + deltaT) > m_maxTracingTime)
377  {
378  deltaT = m_maxTracingTime - elapsedTime;
379  bContinue = false; // This will be the last point traced
380  m_exitMessage = "Exceeded or reached max tracing time.";
381  }
382 
383  // compute candidate point
384  pt1.x = pt0.x + deltaT * vx0;
385  pt1.y = pt0.y + deltaT * vy0;
386 
387  if (!GetVectorAtLocationAndTime(pt1, a_ptTime + elapsedTime + deltaT, vtkVec))
388  {
389  a_outTrace.clear();
390  a_outTimes.clear();
391  m_exitMessage = "Error occurred while extracting point1";
392  return;
393  }
394  // if pt1 outside of domain, compute new deltaT to get to boundary
395  if (EQ_TOL(vtkVec.x, XM_NODATA, 1) || EQ_TOL(vtkVec.y, XM_NODATA, 1))
396  {
397  m_exitMessage = "Point has traveled out of domain.";
398  VecPt3d points = {pt0, pt1};
399  // DataLocationEnum is irrelevant here.
400  BSHP<XmUGrid2dPolylineDataExtractor> polylineExtractor =
401  XmUGrid2dPolylineDataExtractor::New(m_ugrid, DataLocationEnum::LOC_POINTS);
402  polylineExtractor->SetPolyline(points);
403  points = polylineExtractor->GetExtractLocations();
404  if (points.size() < 3)
405  {
406  XM_LOG(xmlog::error, "Gridtracer failed to find an intersection when exiting grid.");
407  return;
408  }
409  double segDist = Mdist(pt0.x, pt0.y, pt1.x, pt1.y);
410  pt1 = points[points.size() - 2];
411  double newSegDist = Mdist(pt0.x, pt0.y, pt1.x, pt1.y);
412  deltaT *= (newSegDist / segDist);
413  bContinue = false;
414  if (!GetVectorAtLocationAndTime(pt1, a_ptTime + elapsedTime + deltaT, vtkVec) ||
415  vtkVec.x == XM_NODATA || vtkVec.y == XM_NODATA)
416  {
417  m_exitMessage = "Error occurred while extracting point1";
418  return;
419  }
420  }
421  vx1 = vtkVec.x;
422  vy1 = vtkVec.y;
423  vx1 *= m_vectorMultiplier;
424  vy1 *= m_vectorMultiplier;
425 
426  if (EQ_TOL(vx1, 0.0, .0001) && EQ_TOL(vy1, 0.0, .0001)) // No velocity
427  {
428  a_outTrace.push_back(pt1);
429  a_outTimes.push_back(a_ptTime + elapsedTime + deltaT);
430  m_exitMessage = "Velocity has gone to zero.";
431  return;
432  }
433  bool bSplit = false;
434 
435  mag1 = sqrt(vx1 * vx1 + vy1 * vy1);
436 
437  // do we subdivide?
438 
439  if (!bSplit && m_maxChangeVelocity > 0)
440  {
441  double changeVel = fabs(mag1 - mag0);
442  if (changeVel > m_maxChangeVelocity)
443  {
444  bSplit = true;
445  m_exitMessage = "Point has exceeded max change velocity.";
446  }
447  }
448  if (!bSplit && m_maxChangeDirectionInRadians > 0)
449  {
450  double dir = iGetDirAsCosTheta(vx0, vy0, vx1, vy1);
451  if (dir < maxAngleChange)
452  {
453  bSplit = true;
454  m_exitMessage = "Point has exceeded max change direction.";
455  }
456  }
457  if (bSplit)
458  {
459  bContinue = true;
460  deltaT /= 2;
461  if (m_minDeltaTime > 0 && deltaT < m_minDeltaTime)
462  {
463  // done, exit
464  bContinue = false;
465  m_exitMessage += " Delta time was less than min delta time.";
466  }
467  }
468  else
469  {
470  double segDist = Mdist(pt0.x, pt0.y, pt1.x, pt1.y);
471  m_distTraveled += segDist;
473  {
474  // because our last point exceeded the exitDistance
475  // find this point by linear calculations
476  double distancePast = m_distTraveled - m_maxTracingDistance;
477  double perc = distancePast / segDist;
478  Pt3d newPt;
479  newPt.x = (pt0.x * perc) + (pt1.x * (1 - perc));
480  newPt.y = (pt0.y * perc) + (pt1.y * (1 - perc));
481 
483  a_outTrace.push_back(newPt);
484  a_outTimes.push_back(a_ptTime + elapsedTime + deltaT * perc);
485  m_exitMessage = "Point has reached or exceeded the max tracing distance.";
486  return;
487  }
488 
489  // add new pt if not identical to last
490  int size = (int)a_outTrace.size();
491  if (size > 0)
492  {
493  if (!EQ_TOL(pt1.x, a_outTrace.at(size - 1).x, XM_ZERO_TOL) ||
494  !EQ_TOL(pt1.y, a_outTrace.at(size - 1).y, XM_ZERO_TOL))
495  {
496  a_outTrace.push_back(pt1);
497  }
498  }
499  else
500  {
501  a_outTrace.push_back(pt1);
502  }
503  pt0 = pt1;
504  elapsedTime += deltaT;
505  vx0 = vx1;
506  vy0 = vy1;
507  deltaT *= 1.2;
508  mag0 = mag1;
509  a_outTimes.push_back(a_ptTime + elapsedTime);
510  }
511  } // while ()
512 } // XmGridTraceImpl::TracePoint
513 //------------------------------------------------------------------------------
518 //------------------------------------------------------------------------------
519 bool XmGridTraceImpl::GetVectorAtLocationAndTime(const xms::Pt3d& a_pt,
520  double a_currentTime,
521  xms::Pt3d& a_data) const
522 {
523  xms::VecPt3d loc;
524  loc.push_back(a_pt);
525  m_extractor1x->SetExtractLocations(loc);
526  m_extractor1y->SetExtractLocations(loc);
527  xms::VecFlt dataOutx1;
528  xms::VecFlt dataOuty1;
529  m_extractor1x->ExtractData(dataOutx1);
530  m_extractor1y->ExtractData(dataOuty1);
531  if (dataOutx1.size() != 1 || dataOuty1.size() != 1)
532  {
533  XM_LOG(xmlog::error, "Gridtracer: An error occured when extracting data");
534  return false;
535  }
536 
537  m_extractor2x->SetExtractLocations(loc);
538  m_extractor2y->SetExtractLocations(loc);
539  xms::VecFlt dataOutx2;
540  xms::VecFlt dataOuty2;
541  m_extractor2x->ExtractData(dataOutx2);
542  m_extractor2y->ExtractData(dataOuty2);
543  if (dataOutx2.size() != 1 || dataOuty2.size() != 1)
544  {
545  XM_LOG(xmlog::error, "Gridtracer: An error occured when extracting data");
546  return false;
547  }
548 
549  if (a_currentTime < m_time1 - XM_ZERO_TOL)
550  {
551  XM_LOG(xmlog::warning, "Gridtracer: The given time is before the first time step.");
552  a_currentTime = m_time1;
553  }
554  double totalTime = fabs(m_time1 - m_time2);
555  double perc1 = fabs(a_currentTime - m_time1) / totalTime;
556  double perc2 = fabs(a_currentTime - m_time2) / totalTime;
557  a_data.x = dataOutx1[0] * perc1 + dataOutx2[0] * perc2;
558  a_data.y = dataOuty1[0] * perc1 + dataOuty2[0] * perc2;
559  return true;
560 } // XmGridTraceImpl::GetVectorAtLocationAndTime
561 } // namespace {}
566 //------------------------------------------------------------------------------
568 //------------------------------------------------------------------------------
569 XmGridTrace::XmGridTrace()
570 {
571 } // XmGridTrace::XmGridTrace
572 //------------------------------------------------------------------------------
574 //------------------------------------------------------------------------------
576 {
577 } // XmGridTrace::~XmGridTrace
578 //------------------------------------------------------------------------------
582 //------------------------------------------------------------------------------
583 BSHP<XmGridTrace> XmGridTrace::New(BSHP<XmUGrid> a_ugrid)
584 {
585  return BSHP<XmGridTraceImpl>(new XmGridTraceImpl(a_ugrid));
586 } // XmGridTrace::New
587 
588 } // namespace xms
589 #ifdef CXX_TEST
591 
593 #include <xmsgrid/ugrid/XmUGrid.h>
594 
595 using namespace xms;
596 namespace
597 {
598 //------------------------------------------------------------------------------
601 //------------------------------------------------------------------------------
602 void iCreateDefaultSingleCell(BSHP<XmGridTrace>& a_tracer)
603 {
604  // 3----2
605  // | 1 /|
606  // | / |
607  // | / |
608  // |/ 0 |
609  // 0----1
610  VecPt3d points = {{0, 0, 0}, {1, 0, 0}, {1, 1, 0}, {0, 1, 0}};
611  VecInt cells = {XMU_TRIANGLE, 3, 0, 1, 2, XMU_TRIANGLE, 3, 2, 3, 0};
612  BSHP<XmUGrid> ugrid = XmUGrid::New(points, cells);
613  a_tracer = XmGridTrace::New(ugrid);
614  const double vm = 1;
615  a_tracer->SetVectorMultiplier(vm);
616  TS_ASSERT_EQUALS(a_tracer->GetVectorMultiplier(), vm);
617 
618  const double tt = 100;
619  a_tracer->SetMaxTracingTime(tt);
620  TS_ASSERT_EQUALS(a_tracer->GetMaxTracingTime(), tt);
621 
622  const double td = 100;
623  a_tracer->SetMaxTracingDistance(td);
624  TS_ASSERT_EQUALS(a_tracer->GetMaxTracingDistance(), td);
625 
626  const double dt = .1;
627  a_tracer->SetMinDeltaTime(dt);
628  TS_ASSERT_EQUALS(a_tracer->GetMinDeltaTime(), dt);
629 
630  const double cd = 100;
631  a_tracer->SetMaxChangeDistance(cd);
632  TS_ASSERT_EQUALS(a_tracer->GetMaxChangeDistance(), cd);
633 
634  const double cv = 100;
635  a_tracer->SetMaxChangeVelocity(cv);
636  TS_ASSERT_EQUALS(a_tracer->GetMaxChangeVelocity(), cv);
637 
638  const double cdir = 1.5 * XM_PI;
639  a_tracer->SetMaxChangeDirectionInRadians(cdir);
640  TS_ASSERT_EQUALS(a_tracer->GetMaxChangeDirectionInRadians(), cdir);
641 
642  double time = 0;
643  VecPt3d scalars = {{1, 1, 0}, {1, 1, 0}, {1, 1, 0}, {1, 1, 0}};
644  DynBitset pointActivity;
645  for (int i = 0; i < 4; ++i)
646  {
647  pointActivity.push_back(true);
648  }
649  a_tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
650  DataLocationEnum::LOC_POINTS, time);
651 
652  time = 10;
653  // Uses exact same scalars/pointActivity
654  a_tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
655  DataLocationEnum::LOC_POINTS, time);
656 } // iCreateDefaultSingleCell
657 //------------------------------------------------------------------------------
660 //------------------------------------------------------------------------------
661 void iCreateDefaultTwoCell(BSHP<XmGridTrace>& a_tracer)
662 {
663  // clang-format off
664  // 3========2========5
665  // | | |
666  // | | |
667  // | | |
668  // 0--------1--------4
669  // clang-format on
670 
671  VecPt3d points = {{0, 0, 0}, {1, 0, 0}, {1, 1, 0}, {0, 1, 0}, {2, 0, 0}, {2, 1, 0}};
672  VecInt cells = {XMU_QUAD, 4, 0, 1, 2, 3, XMU_QUAD, 4, 1, 4, 5, 2};
673  BSHP<XmUGrid> ugrid = XmUGrid::New(points, cells);
674  a_tracer = XmGridTrace::New(ugrid);
675  const double vm = 1;
676  a_tracer->SetVectorMultiplier(vm);
677 
678  const double tt = 100;
679  a_tracer->SetMaxTracingTime(tt);
680 
681  const double td = 100;
682  a_tracer->SetMaxTracingDistance(td);
683 
684  const double dt = .1;
685  a_tracer->SetMinDeltaTime(dt);
686 
687  const double cd = 100;
688  a_tracer->SetMaxChangeDistance(cd);
689 
690  const double cv = 100;
691  a_tracer->SetMaxChangeVelocity(cv);
692 
693  const double cdir = 1.5 * XM_PI;
694  a_tracer->SetMaxChangeDirectionInRadians(cdir);
695 
696  double time = 0;
697  VecPt3d scalars = {{.1, 0, 0}, {.2, 0, 0}};
698  DynBitset pointActivity;
699  for (int i = 0; i < 2; ++i)
700  {
701  pointActivity.push_back(true);
702  }
703  a_tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_CELLS, pointActivity,
704  DataLocationEnum::LOC_CELLS, time);
705 
706  time = 10;
707  // Uses exact same scalars/pointActivity
708  a_tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_CELLS, pointActivity,
709  DataLocationEnum::LOC_CELLS, time);
710 } // iCreateDefaultTwoCell
711 }
716 //------------------------------------------------------------------------------
718 //------------------------------------------------------------------------------
720 {
721  BSHP<XmGridTrace> tracer;
722  iCreateDefaultSingleCell(tracer);
723 
724  VecPt3d outTrace;
725  VecDbl outTimes;
726  Pt3d startPoint = {.5, .5, 0};
727  double startTime = .5;
728 
729  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
730 
731  VecPt3d expectedOutTrace = {{.5, .5, 0}, {1, 1, 0}};
732  VecDbl expectedOutTimes = {.5, 1};
733  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
734  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
735 } // XmGridTraceUnitTests::testTracePoint
736 //------------------------------------------------------------------------------
738 //------------------------------------------------------------------------------
740 {
741  BSHP<XmGridTrace> tracer;
742  iCreateDefaultSingleCell(tracer);
743  tracer->SetMaxChangeDistance(.25);
744 
745  VecPt3d outTrace;
746  VecDbl outTimes;
747  Pt3d startPoint = {.5, .5, 0};
748  double startTime = .5;
749 
750  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
751 
752  VecPt3d expectedOutTrace = {{.5, .5, 0},
753  {0.67677668424809445, 0.67677668424809445, 0.00000000000000000},
754  {0.85355336849618890, 0.85355336849618890, 0.00000000000000000},
755  {1, 1, 0}};
756  VecDbl expectedOutTimes = {.5, 0.67677668424809445, 0.85355336849618890, 1};
757  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
758  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
759 } // XmGridTraceUnitTests::testMaxChangeDistance
760 //------------------------------------------------------------------------------
762 //------------------------------------------------------------------------------
764 {
765  BSHP<XmGridTrace> tracer;
766  iCreateDefaultSingleCell(tracer);
767 
768  // Push on different scalars
769  double time = 0;
770  VecPt3d scalars = {{.1, .1, 0}, {.1, .1, 0}, {.1, .1, 0}, {.1, .1, 0}};
771  DynBitset pointActivity;
772  for (int i = 0; i < 4; ++i)
773  {
774  pointActivity.push_back(true);
775  }
776  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
777  DataLocationEnum::LOC_POINTS, time);
778 
779  time = 10;
780  // Uses exact same scalars/pointActivity
781  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
782  DataLocationEnum::LOC_POINTS, time);
783 
784  VecPt3d outTrace;
785  VecDbl outTimes;
786  Pt3d startPoint = {.5, .5, 0};
787  double startTime = .5;
788  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
789 
790  VecPt3d expectedOutTrace = {{.5, .5, 0},
791  {0.60000000149011612, 0.60000000149011612, 0},
792  {0.72000000327825542, 0.72000000327825542, 0},
793  {0.86400000542402267, 0.86400000542402267, 0},
794  {1, 1, 0}};
795  VecDbl expectedOutTimes = {.5, 1.5, 2.7, 4.14, 5.5};
796 
797  if (expectedOutTrace.size() == outTrace.size())
798  {
799  for (int i = 0; i < expectedOutTrace.size(); ++i)
800  {
801  TS_ASSERT_DELTA(expectedOutTrace[i].x, outTrace[i].x, .001);
802  TS_ASSERT_DELTA(expectedOutTrace[i].y, outTrace[i].y, .001);
803  TS_ASSERT_DELTA(expectedOutTrace[i].z, outTrace[i].z, .001);
804  }
805  }
806  else
807  TS_FAIL("Expected trace size != actual trace size");
808  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
809  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, 0.001);
810 } // XmGridTraceUnitTests::testSmallScalarsTracePoint
811 //------------------------------------------------------------------------------
813 //------------------------------------------------------------------------------
815 {
816  BSHP<XmGridTrace> tracer;
817  iCreateDefaultSingleCell(tracer);
818 
819  tracer->SetMaxChangeDirectionInRadians(XM_PI * .2);
820  tracer->SetMinDeltaTime(-1);
821  // Push on different scalars
822  double time = 0;
823  VecPt3d scalars = {{0, 1, 0}, {-1, 0, 0}, {0, -1, 0}, {1, 0, 0}};
824  DynBitset pointActivity;
825  for (int i = 0; i < 4; ++i)
826  {
827  pointActivity.push_back(true);
828  }
829  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
830  DataLocationEnum::LOC_POINTS, time);
831 
832  time = 10;
833  // Uses exact same scalars/pointActivity
834  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
835  DataLocationEnum::LOC_POINTS, time);
836 
837  VecPt3d outTrace;
838  VecDbl outTimes;
839  Pt3d startPoint = {0, 0, 0};
840  double startTime = .5;
841  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
842 
843  VecPt3d expectedOutTrace = {{0, 0, 0},
844  {0.00000000000000000, 0.25000000000000000, 0.00000000000000000},
845  {0.074999999999999997, 0.47499999999999998, 0.00000000000000000},
846  {0.21900000214576720, 0.63699999570846555, 0.00000000000000000},
847  {0.30928799843788146, 0.66810399758815764, 0.00000000000000000},
848  {0.40229310507774352, 0.67396399235725402, 0.00000000000000000},
849  {0.48679361495018003, 0.65024498560905453, 0.00000000000000000},
850  {0.54780151323509219, 0.59909560095787040, 0.00000000000000000},
851  {0.55928876277122497, 0.56619817004051198, 0.00000000000000000},
852  {0.56114558691518779, 0.53247499044700608, 0.00000000000000000},
853  {0.55189971330840681, 0.50228363992173752, 0.00000000000000000},
854  {0.53269911067322617, 0.48131557500677169, 0.00000000000000000},
855  {0.52076836142536975, 0.47806150355091476, 0.00000000000000000},
856  {0.50886902895577013, 0.47838753608466128, 0.00000000000000000},
857  {0.49867742691962913, 0.48264835153512164, 0.00000000000000000},
858  {0.49224616907898289, 0.49014090685121131, 0.00000000000000000},
859  {0.49173935940609609, 0.49438094923206660, 0.00000000000000000},
860  {0.49250246625151450, 0.49839053740482164, 0.00000000000000000},
861  {0.49454361321306389, 0.50154755045413602, 0.00000000000000000},
862  {0.49745717820065949, 0.50317358562752867, 0.00000000000000000},
863  {0.49888395770889871, 0.50301615091938545, 0.00000000000000000},
864  {0.50012160117661586, 0.50244704462921241, 0.00000000000000000},
865  {0.50095740046883197, 0.50152383477622209, 0.00000000000000000},
866  {0.50107955145675120, 0.50098875888952354, 0.00000000000000000},
867  {0.50105605626599747, 0.50045352403892940, 0.00000000000000000},
868  {0.50086894918345870, 0.49998474718699493, 0.00000000000000000},
869  {0.50053945884260675, 0.49966662451478433, 0.00000000000000000},
870  {0.50034430627617277, 0.49962054739305783, 0.00000000000000000},
871  {0.50015012042108842, 0.49962997721910873, 0.00000000000000000},
872  {0.49998265395837810, 0.49970077747304897, 0.00000000000000000},
873  {0.49987374966305814, 0.49982308521808211, 0.00000000000000000},
874  {0.49986302487024292, 0.49988726006383088, 0.00000000000000000},
875  {0.49986815504448728, 0.49994012045071656, 0.00000000000000000}};
876  VecDbl expectedOutTimes = {.5,
877  0.75000000000000000,
878  1.0500000000000000,
879  1.4100000000000001,
880  1.6260000000000001,
881  1.8852000000000002,
882  2.1962400000000004,
883  2.5694880000000002,
884  2.7934368000000003,
885  3.0621753600000003,
886  3.3846616320000003,
887  3.7716451584000001,
888  4.0038352742400001,
889  4.2824634132480002,
890  4.6168171800576001,
891  5.0180417002291202,
892  5.2587764123320317,
893  5.5476580668555258,
894  5.8943160522837186,
895  6.3103056347975501,
896  6.5598993843058491,
897  6.8594118837158078,
898  7.2188268830077584,
899  7.4344758825829285,
900  7.6932546820731327,
901  8.0037892414613783,
902  8.3764307127272737,
903  8.6000155954868092,
904  8.8683174547982535,
905  9.1902796859719871,
906  9.5766343633804656,
907  9.7883171816902319,
908  10.000000000000000};
909  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
910  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
911 } // XmGridTraceUnitTests::testStrongDirectionChange
912 //------------------------------------------------------------------------------
914 //------------------------------------------------------------------------------
916 {
917  BSHP<XmGridTrace> tracer;
918  iCreateDefaultSingleCell(tracer);
919 
920  tracer->SetMaxChangeDirectionInRadians(XM_PI * .2);
921  tracer->SetMinDeltaTime(-1);
922  tracer->SetMaxTracingTime(5);
923  // Push on different scalars
924  double time = 0;
925  VecPt3d scalars = {{0, 1, 0}, {-1, 0, 0}, {0, -1, 0}, {1, 0, 0}};
926  DynBitset pointActivity;
927  for (int i = 0; i < 4; ++i)
928  {
929  pointActivity.push_back(true);
930  }
931  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
932  DataLocationEnum::LOC_POINTS, time);
933 
934  time = 10;
935  // Uses exact same scalars/pointActivity
936  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
937  DataLocationEnum::LOC_POINTS, time);
938 
939  VecPt3d outTrace;
940  VecDbl outTimes;
941  Pt3d startPoint = {0, 0, 0};
942  double startTime = .5;
943  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
944 
945  VecPt3d expectedOutTrace = {{0, 0, 0},
946  {0.00000000000000000, 0.25000000000000000, 0.00000000000000000},
947  {0.074999999999999997, 0.47499999999999998, 0.00000000000000000},
948  {0.21900000214576720, 0.63699999570846555, 0.00000000000000000},
949  {0.30928799843788146, 0.66810399758815764, 0.00000000000000000},
950  {0.40229310507774352, 0.67396399235725402, 0.00000000000000000},
951  {0.48679361495018003, 0.65024498560905453, 0.00000000000000000},
952  {0.54780151323509219, 0.59909560095787040, 0.00000000000000000},
953  {0.55928876277122497, 0.56619817004051198, 0.00000000000000000},
954  {0.56114558691518779, 0.53247499044700608, 0.00000000000000000},
955  {0.55189971330840681, 0.50228363992173752, 0.00000000000000000},
956  {0.53269911067322617, 0.48131557500677169, 0.00000000000000000},
957  {0.52076836142536975, 0.47806150355091476, 0.00000000000000000},
958  {0.50886902895577013, 0.47838753608466128, 0.00000000000000000},
959  {0.49867742691962913, 0.48264835153512164, 0.00000000000000000},
960  {0.49224616907898289, 0.49014090685121131, 0.00000000000000000},
961  {0.49173935940609609, 0.49438094923206660, 0.00000000000000000},
962  {0.49237657318600692, 0.49772905815126539, 0.00000000000000000}};
963  VecDbl expectedOutTimes = {.5,
964  0.75000000000000000,
965  1.0500000000000000,
966  1.4100000000000001,
967  1.6260000000000001,
968  1.8852000000000002,
969  2.1962400000000004,
970  2.5694880000000002,
971  2.7934368000000003,
972  3.0621753600000003,
973  3.3846616320000003,
974  3.7716451584000001,
975  4.0038352742400001,
976  4.2824634132480002,
977  4.6168171800576001,
978  5.0180417002291202,
979  5.2587764123320317,
980  5.5};
981  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
982  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
983 } // XmGridTraceUnitTests::testMaxTracingTime
984 //------------------------------------------------------------------------------
986 //------------------------------------------------------------------------------
988 {
989  BSHP<XmGridTrace> tracer;
990  iCreateDefaultSingleCell(tracer);
991 
992  tracer->SetMaxChangeDirectionInRadians(XM_PI * .2);
993  tracer->SetMinDeltaTime(-1);
994  tracer->SetMaxTracingDistance(1.0);
995  // Push on different scalars
996  double time = 0;
997  VecPt3d scalars = {{0, 1, 0}, {-1, 0, 0}, {0, -1, 0}, {1, 0, 0}};
998  DynBitset pointActivity;
999  for (int i = 0; i < 4; ++i)
1000  {
1001  pointActivity.push_back(true);
1002  }
1003  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
1004  DataLocationEnum::LOC_POINTS, time);
1005 
1006  time = 10;
1007  // Uses exact same scalars/pointActivity
1008  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
1009  DataLocationEnum::LOC_POINTS, time);
1010 
1011  VecPt3d outTrace;
1012  VecDbl outTimes;
1013  Pt3d startPoint = {0, 0, 0};
1014  double startTime = .5;
1015  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1016 
1017  VecPt3d expectedOutTrace = {{0, 0, 0},
1018  {0.00000000000000000, 0.25000000000000000, 0.00000000000000000},
1019  {0.074999999999999997, 0.47499999999999998, 0.00000000000000000},
1020  {0.21900000214576720, 0.63699999570846555, 0.00000000000000000},
1021  {0.30928799843788146, 0.66810399758815764, 0.00000000000000000},
1022  {0.40229310507774352, 0.67396399235725402, 0.00000000000000000},
1023  {0.48679361495018003, 0.65024498560905453, 0.00000000000000000},
1024  {0.50183556502673621, 0.63763372523131490, 0.00000000000000000}};
1025  VecDbl expectedOutTimes = {.5,
1026  0.75000000000000000,
1027  1.0500000000000000,
1028  1.4100000000000001,
1029  1.6260000000000001,
1030  1.8852000000000002,
1031  2.1962400000000004,
1032  2.4774609356360582};
1033  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1034  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1035 } // XmGridTraceUnitTests::testMaxTracingDistance
1036 //------------------------------------------------------------------------------
1038 //------------------------------------------------------------------------------
1040 {
1041  BSHP<XmGridTrace> tracer;
1042  iCreateDefaultSingleCell(tracer);
1043 
1044  VecPt3d outTrace;
1045  VecDbl outTimes;
1046  Pt3d startPoint = {-.1, 0, 0};
1047  double startTime = .5;
1048 
1049  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1050 
1051  VecPt3d expectedOutTrace = {};
1052  VecDbl expectedOutTimes = {};
1053  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1054  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1055 } // XmGridTraceUnitTests::testStartOutOfCell
1056 //------------------------------------------------------------------------------
1058 //------------------------------------------------------------------------------
1060 {
1061  TS_ASSERT_DELTA(iGetDirAsCosTheta(0, 1, 1, 0), 0, .1); // 90 degree angle
1062  TS_ASSERT_DELTA(iGetDirAsCosTheta(0, 1, 1, 1), .707, .1); // 45 degree angle
1063  TS_ASSERT_DELTA(iGetDirAsCosTheta(0, 1, 0, -1), -1, .1); // 180 degree angle
1064  TS_ASSERT_DELTA(iGetDirAsCosTheta(1, 1, -1, -1), -1, .1); // 180 degree angle
1065  TS_ASSERT_DELTA(iGetDirAsCosTheta(2, 5, 3, 6), .9965, .1); // almost 0 degree angle
1066  TS_ASSERT_DELTA(iGetDirAsCosTheta(0, 1, 0, 1), 1, .1); // 0 degree angle
1067 } // XmGridTraceUnitTests::testDotProduct
1068 //------------------------------------------------------------------------------
1070 //------------------------------------------------------------------------------
1072 {
1073  BSHP<XmGridTrace> tracer;
1074  iCreateDefaultSingleCell(tracer);
1075 
1076  VecPt3d outTrace;
1077  VecDbl outTimes;
1078  Pt3d startPoint = {.5, .5, 0};
1079  double startTime = 10.1; // Just beyond the 2nd time
1080 
1081  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1082 
1083  VecPt3d expectedOutTrace = {};
1084  VecDbl expectedOutTimes = {};
1085  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1086  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1087 } // XmGridTraceUnitTests::testBeyondTimestep
1088 //------------------------------------------------------------------------------
1090 //------------------------------------------------------------------------------
1092 {
1093  BSHP<XmGridTrace> tracer;
1094  iCreateDefaultSingleCell(tracer);
1095 
1096  VecPt3d outTrace;
1097  VecDbl outTimes;
1098  Pt3d startPoint = {.5, .5, 0};
1099  double startTime = -0.1; // Just beyond the 2nd time
1100 
1101  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1102 
1103  VecPt3d expectedOutTrace = {{.5, .5, 0}, {1, 1, 0}};
1104  VecDbl expectedOutTimes = {-.1, .4};
1105  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1106  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1107 } // XmGridTraceUnitTests::testBeforeTimestep
1108 //------------------------------------------------------------------------------
1110 //------------------------------------------------------------------------------
1112 {
1113  BSHP<XmGridTrace> tracer;
1114  iCreateDefaultSingleCell(tracer);
1115 
1116  tracer->SetMaxChangeDirectionInRadians(XM_PI * .2);
1117  tracer->SetMinDeltaTime(-1);
1118  tracer->SetVectorMultiplier(0.5);
1119  // Push on different scalars
1120  double time = 0;
1121  VecPt3d scalars = {{0, 1, 0}, {-1, 0, 0}, {0, -1, 0}, {1, 0, 0}};
1122  DynBitset pointActivity;
1123  for (int i = 0; i < 4; ++i)
1124  {
1125  pointActivity.push_back(true);
1126  }
1127  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
1128  DataLocationEnum::LOC_POINTS, time);
1129 
1130  time = 10;
1131  // Uses exact same scalars/pointActivity
1132  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_POINTS, pointActivity,
1133  DataLocationEnum::LOC_POINTS, time);
1134 
1135  VecPt3d outTrace;
1136  VecDbl outTimes;
1137  Pt3d startPoint = {0, 0, 0};
1138  double startTime = .5;
1139  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1140 
1141  VecPt3d expectedOutTrace = {{0, 0, 0},
1142  {0.00000000000000000, 0.25000000000000000, 0.00000000000000000},
1143  {0.074999999999999997, 0.47499999999999998, 0.00000000000000000},
1144  {0.21900000214576720, 0.63699999570846555, 0.00000000000000000},
1145  {0.30928799843788146, 0.66810399758815764, 0.00000000000000000},
1146  {0.40229310507774352, 0.67396399235725402, 0.00000000000000000},
1147  {0.48679361495018003, 0.65024498560905453, 0.00000000000000000},
1148  {0.54780151323509219, 0.59909560095787040, 0.00000000000000000},
1149  {0.55928876277122497, 0.56619817004051198, 0.00000000000000000},
1150  {0.56114558691518779, 0.53247499044700608, 0.00000000000000000},
1151  {0.55189971330840681, 0.50228363992173752, 0.00000000000000000},
1152  {0.53269911067322617, 0.48131557500677169, 0.00000000000000000},
1153  {0.52076836142536975, 0.47806150355091476, 0.00000000000000000},
1154  {0.50886902895577013, 0.47838753608466128, 0.00000000000000000},
1155  {0.49867742691962913, 0.48264835153512164, 0.00000000000000000},
1156  {0.49224616907898289, 0.49014090685121131, 0.00000000000000000},
1157  {0.49175783605462037, 0.49422637094165467, 0.00000000000000000}};
1158  VecDbl expectedOutTimes = {.5,
1159  1.0000000000000000,
1160  1.6000000000000001,
1161  2.3200000000000003,
1162  2.7520000000000002,
1163  3.2704000000000004,
1164  3.8924800000000004,
1165  4.6389760000000004,
1166  5.0868736000000006,
1167  5.6243507200000007,
1168  6.2693232640000005,
1169  7.0432903168000003,
1170  7.5076705484800001,
1171  8.0649268264960003,
1172  8.7336343601152002,
1173  9.5360834004582404,
1174  10.000000000000000};
1175  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1176  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1177 } // XmGridTraceUnitTests::testVectorMultiplier
1178 //------------------------------------------------------------------------------
1180 //------------------------------------------------------------------------------
1182 {
1183  BSHP<XmGridTrace> tracer;
1184  iCreateDefaultTwoCell(tracer);
1185 
1186  VecPt3d outTrace;
1187  VecDbl outTimes;
1188  Pt3d startPoint = {.5, .5, 0};
1189  double startTime = 0;
1190 
1191  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1192 
1193  VecPt3d expectedOutTrace = {{.5, .5, 0},
1194  {0.60000000149011612, 0.50000000000000000, 0.00000000000000000},
1195  {0.73200000077486038, 0.50000000000000000, 0.00000000000000000},
1196  {0.90940801054239273, 0.50000000000000000, 0.00000000000000000},
1197  {1.1529537134766579, 0.50000000000000000, 0.00000000000000000},
1198  {1.4957102079987525, 0.50000000000000000, 0.00000000000000000},
1199  {1.9923067892670629, 0.50000000000000000, 0.00000000000000000},
1200  {2, .5, 0}};
1201  VecDbl expectedOutTimes = {0,
1202  1.0000000000000000,
1203  2.2000000000000002,
1204  3.6400000000000001,
1205  5.3680000000000003,
1206  7.4416000000000002,
1207  9.9299199999999992,
1208  9.9683860530914945};
1209  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1210  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1211 } // XmGridTraceUnitTests::testMultiCell
1212 //------------------------------------------------------------------------------
1216 //------------------------------------------------------------------------------
1218 {
1219  BSHP<XmGridTrace> tracer;
1220  iCreateDefaultTwoCell(tracer);
1221  tracer->SetMaxChangeVelocity(.01);
1222  tracer->SetMinDeltaTime(0.001);
1223 
1224  VecPt3d outTrace;
1225  VecDbl outTimes;
1226  Pt3d startPoint = {.5, .5, 0};
1227  double startTime = 0;
1228 
1229  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1230 
1231  VecPt3d expectedOutTrace = {{.5, .5, 0},
1232  {0.60000000149011612, 0.50000000000000000, 0.00000000000000000},
1233  {0.66600000113248825, 0.50000000000000000, 0.00000000000000000},
1234  {0.74995200067758561, 0.50000000000000000, 0.00000000000000000},
1235  {0.80394992786645891, 0.50000000000000000, 0.00000000000000000},
1236  {0.87154669338464741, 0.50000000000000000, 0.00000000000000000},
1237  {0.95686786960840231, 0.50000000000000000, 0.00000000000000000},
1238  {1.0112451727318765, 0.50000000000000000, 0.00000000000000000},
1239  {1.0789334834771158, 0.50000000000000000, 0.00000000000000000},
1240  {1.1637975516948893, 0.50000000000000000, 0.00000000000000000},
1241  {1.2174527417415202, 0.50000000000000000, 0.00000000000000000},
1242  {1.2839153379250163, 0.50000000000000000, 0.00000000000000000},
1243  {1.3667568384715398, 0.50000000000000000, 0.00000000000000000},
1244  {1.4187699365302351, 0.50000000000000000, 0.00000000000000000},
1245  {1.4829247317645506, 0.50000000000000000, 0.00000000000000000},
1246  {1.5624845364724593, 0.50000000000000000, 0.00000000000000000},
1247  {1.6587784227485147, 0.50000000000000000, 0.00000000000000000},
1248  {1.7743310862797812, 0.50000000000000000, 0.00000000000000000},
1249  {1.9129942825173010, 0.50000000000000000, 0.00000000000000000},
1250  {2, .5, 0}};
1251  VecDbl expectedOutTimes = {0,
1252  1.0000000000000000,
1253  1.6000000000000001,
1254  2.3200000000000003,
1255  2.7520000000000002,
1256  3.2704000000000004,
1257  3.8924800000000004,
1258  4.2657280000000002,
1259  4.7136256000000003,
1260  5.2511027200000004,
1261  5.5735889920000004,
1262  5.9605725184000002,
1263  6.4249527500800001,
1264  6.7035808890880002,
1265  7.0379346558976001,
1266  7.4391591760691202,
1267  7.9206286002749442,
1268  8.4983919093219331,
1269  9.1917078801783187,
1270  9.6267364611093829};
1271  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1272  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1273 } // XmGridTraceUnitTests::testMaxChangeVelocity
1274 //------------------------------------------------------------------------------
1276 //------------------------------------------------------------------------------
1278 {
1279  BSHP<XmGridTrace> tracer;
1280  iCreateDefaultTwoCell(tracer);
1281 
1282  double time = 20;
1283  VecPt3d scalars = {{.2, 0, 0}, {.3, 0, 0}};
1284  DynBitset pointActivity;
1285  for (int i = 0; i < 2; ++i)
1286  {
1287  pointActivity.push_back(true);
1288  }
1289  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_CELLS, pointActivity,
1290  DataLocationEnum::LOC_CELLS, time);
1291 
1292  VecPt3d outTrace;
1293  VecDbl outTimes;
1294  Pt3d startPoint = {.5, .5, 0};
1295  double startTime = 10;
1296 
1297  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1298 
1299  VecPt3d expectedOutTrace = {{.5, .5, 0},
1300  {0.70000000298023224, 0.50000000000000000, 0.00000000000000000},
1301  {0.95200000226497650, 0.50000000000000000, 0.00000000000000000},
1302  {1.2734079944372176, 0.50000000000000000, 0.00000000000000000},
1303  {1.6897536998434066, 0.50000000000000000, 0.00000000000000000},
1304  {2, .5, 0}};
1305  VecDbl expectedOutTimes = {10,
1306  11.000000000000000,
1307  12.199999999999999,
1308  13.640000000000001,
1309  15.368000000000000,
1310  16.627525378316030};
1311  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1312  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1313 } // XmGridTraceUnitTests::testUniqueTimeSteps
1314 //------------------------------------------------------------------------------
1318 //------------------------------------------------------------------------------
1320 {
1321  BSHP<XmGridTrace> tracer;
1322  iCreateDefaultTwoCell(tracer);
1323 
1324  double time = 20;
1325  VecPt3d scalars = {{.2, 0, 0}, {99999, 0, 0}};
1326  DynBitset pointActivity;
1327  pointActivity.push_back(true);
1328  pointActivity.push_back(false);
1329  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_CELLS, pointActivity,
1330  DataLocationEnum::LOC_CELLS, time);
1331 
1332  VecPt3d outTrace;
1333  VecDbl outTimes;
1334  Pt3d startPoint = {.5, .5, 0};
1335  double startTime = 10;
1336 
1337  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1338 
1339  VecPt3d expectedOutTrace = {{.5, .5, 0},
1340  {0.70000000298023224, 0.50000000000000000, 0.00000000000000000},
1341  {0.93040000677108770, 0.50000000000000000, 0.00000000000000000},
1342  {0.99788877571821222, 0.50000000000000000, 0.00000000000000000}};
1343  VecDbl expectedOutTimes = {10, 11.000000000000000, 12.199999999999999, 12.560000000000000};
1344  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1345  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1346 } // XmGridTraceUnitTests::testInactiveCell
1347 //------------------------------------------------------------------------------
1349 //------------------------------------------------------------------------------
1351 {
1352  BSHP<XmGridTrace> tracer;
1353  iCreateDefaultTwoCell(tracer);
1354 
1355  double time = 20;
1356  VecPt3d scalars = {{1, 0, 0}, {99999, 0, 0}};
1357  DynBitset pointActivity;
1358  pointActivity.push_back(false);
1359  pointActivity.push_back(true);
1360  tracer->AddGridScalarsAtTime(scalars, DataLocationEnum::LOC_CELLS, pointActivity,
1361  DataLocationEnum::LOC_CELLS, time);
1362 
1363  VecPt3d outTrace;
1364  VecDbl outTimes;
1365  Pt3d startPoint = {.5, .5, 0};
1366  double startTime = 10;
1367 
1368  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1369 
1370  VecPt3d expectedOutTrace = {};
1371  VecDbl expectedOutTimes = {};
1372  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1373  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1374 } // XmGridTraceUnitTests::testStartInactiveCell
1375 //------------------------------------------------------------------------------
1379 //------------------------------------------------------------------------------
1382 {
1383  // Graph with scalar vectors indicated
1384  // -> ->
1385  // 6----7----8|
1386  // | | |v
1387  // | | |
1388  // | | |
1389  // ^| | |
1390  // |3----4----5|
1391  // | | |v
1392  // | | |
1393  // | | |
1394  // ^| | |
1395  // |0----1----2
1396  // <- <--
1397  // Step 1: Create the grid
1398  VecPt3d points = {{0, 0, 0}, {1, 0, 0}, {2, 0, 0}, {0, 1, 0}, {1, 1, 0},
1399  {2, 1, 0}, {0, 2, 0}, {1, 2, 0}, {2, 2, 0}};
1400  VecInt cells = {XMU_QUAD, 4, 0, 1, 4, 3, XMU_QUAD, 4, 1, 2, 5, 4,
1401  XMU_QUAD, 4, 3, 4, 7, 6, XMU_QUAD, 4, 4, 5, 8, 7};
1402  BSHP<XmUGrid> ugrid = XmUGrid::New(points, cells);
1403 
1404  // Step 2: Create the tracer from the grid
1405  BSHP<XmGridTrace> tracer = XmGridTrace::New(ugrid);
1406 
1407  // Step 3: Set up the constraints on the tracer
1408  tracer->SetVectorMultiplier(2);
1409  tracer->SetMaxTracingTime(-1);
1410  tracer->SetMaxTracingDistance(-1);
1411  tracer->SetMinDeltaTime(.01);
1412  tracer->SetMaxChangeDistance(-1);
1413  tracer->SetMaxChangeVelocity(-1);
1414  tracer->SetMaxChangeDirectionInRadians(XM_PI / 4);
1415 
1416  // Step 4: Set up the velocity vectors for both time steps. Insert timesteps sequentially
1417  double time = 0;
1418  // Scalars are set such that they circle around the edge of the graph in a clockwise direction
1419  // Z component is not used in scalars
1420  VecPt3d scalars1 = {{0, 1, 0}, {-.1, 0, 0}, {-1, 0, 0}, {0, .1, 0}, {0, 0, 0},
1421  {0, -.1, 0}, {1, 0, 0}, {.1, 0, 0}, {0, -1, 0}};
1422  DynBitset pointActivity;
1423  for (int i = 0; i < 9; ++i)
1424  {
1425  pointActivity.push_back(true);
1426  }
1427  tracer->AddGridScalarsAtTime(scalars1, DataLocationEnum::LOC_POINTS, pointActivity,
1428  DataLocationEnum::LOC_POINTS, time);
1429 
1430  // For the second timestep scalars are doubled to indicate an increase in magnitude
1431  VecPt3d scalars2 = {{0, 2, 0}, {-.2, 0, 0}, {-2, 0, 0}, {0, .2, 0}, {0, 0, 0},
1432  {0, -.2, 0}, {2, 0, 0}, {.2, 0, 0}, {0, -2, 0}};
1433  time = 20;
1434  // Uses exact same scalars/pointActivity
1435  tracer->AddGridScalarsAtTime(scalars2, DataLocationEnum::LOC_POINTS, pointActivity,
1436  DataLocationEnum::LOC_POINTS, time);
1437 
1438  VecPt3d outTrace;
1439  VecDbl outTimes;
1440  Pt3d startPoint = {.5, .5, 0};
1441  double startTime = 0;
1442 
1443  // Step 5: Trace the point
1444  tracer->TracePoint(startPoint, startTime, outTrace, outTimes);
1445  // show the cause for termination by calling GetExitMessage
1446  // std::cout << tracer->GetExitMessage();
1447 
1448  // Expected values for this simulation
1449  VecPt3d expectedOutTrace = {{0.50000000000000000, 0.50000000000000000, 0.00000000000000000},
1450  {0.50000000000000000, 1.2500000000000000, 0.00000000000000000},
1451  {0.54457812566426578, 1.3391562513285316, 0.00000000000000000},
1452  {0.61632493250262921, 1.4354984729093498, 0.00000000000000000},
1453  {0.72535406450374607, 1.5315533661126233, 0.00000000000000000},
1454  {0.88236797164001590, 1.6126801842666139, 0.00000000000000000},
1455  {0.98873181403598276, 1.6331015959080102, 0.00000000000000000},
1456  {1.0538503898747653, 1.6342606013582104, 0.00000000000000000},
1457  {1.1249433009705341, 1.5683006835455087, 0.00000000000000000},
1458  {1.1895097427498795, 1.3863448896225066, 0.00000000000000000},
1459  {1.2235242118635632, 1.0588590059131318, 0.00000000000000000},
1460  {1.2235242118635632, 0.90477286425654002, 0.00000000000000000},
1461  {1.2005336220528682, 0.85080764250970042, 0.00000000000000000},
1462  {1.1581790674742278, 0.79387770198395835, 0.00000000000000000},
1463  {1.0896874578697060, 0.74131697161132859, 0.00000000000000000},
1464  {0.98966250551038770, 0.70663752692174131, 0.00000000000000000},
1465  {0.95806149614159530, 0.71817980325332686, 0.00000000000000000},
1466  {0.92629620502521459, 0.77371504022050730, 0.00000000000000000},
1467  {0.90239412753251202, 0.88917318465162865, 0.00000000000000000},
1468  {0.89995172701803572, 1.0694875660697027, 0.00000000000000000},
1469  {0.91503139037776327, 1.0911992829869794, 0.00000000000000000},
1470  {0.93816744602651825, 1.1127546977629765, 0.00000000000000000},
1471  {0.97140028507849163, 1.1309789606067331, 0.00000000000000000},
1472  {0.99364912627842006, 1.1358370729524059, 0.00000000000000000},
1473  {1.0071524474802995, 1.1364684019706512, 0.00000000000000000},
1474  {1.0223447138862345, 1.1280655805979485, 0.00000000000000000},
1475  {1.0369737821057583, 1.0971462034407997, 0.00000000000000000},
1476  {1.0467397711865176, 1.0371377237101163, 0.00000000000000000},
1477  {1.0467397711865176, 0.96499504248441559, 0.00000000000000000},
1478  {1.0390576209755447, 0.95473758230148376, 0.00000000000000000},
1479  {1.0276444556154691, 0.94488898976070590, 0.00000000000000000},
1480  {1.0208791233912420, 0.94149540451099356, 0.00000000000000000}};
1481  VecDbl expectedOutTimes = {
1482  0.00000000000000000, 0.37500000000000000, 0.82499999999999996, 1.3649999999999998,
1483  2.0129999999999999, 2.7905999999999995, 3.2571599999999994, 3.5370959999999991,
1484  3.8730191999999990, 4.2761270399999987, 4.7598564479999981, 5.3403317375999979,
1485  6.0369020851199977, 6.8727865021439971, 7.8758478025727969, 9.0795213630873555,
1486  9.4406234312417237, 9.8739459130269651, 10.393932891169255, 11.017917264940003,
1487  11.766698513464901, 12.665236011694777, 13.743481009570628, 14.390428008296139,
1488  14.778596207531445, 15.244398046613812, 15.803360253512654, 16.474114901791264,
1489  17.279020479725595, 18.244907173246794, 19.403971205472232, 20.000000000000000};
1490  TS_ASSERT_DELTA_VECPT3D(expectedOutTrace, outTrace, .0001);
1491  TS_ASSERT_DELTA_VEC(expectedOutTimes, outTimes, .0001);
1492 } // XmGridTraceUnitTests::testTutorial
1494 
1495 #endif
BSHP< XmUGrid2dDataExtractor > m_extractor1x
data extractor for the x component for the first time step
double m_distTraveled
distance traveled in the last TracePoint operation
#define XM_LOG(A, B)
void testStrongDirectionChange()
test behavior when having large changes in direction
std::vector< int > VecInt
void testMultiCell()
test behavior of multiple cells
void testTutorial()
2nd cell is inactive in the 2nd time step. // Thus it does not pull as hard. Also once the point reac...
void testMaxTracingTime()
test setting max tracing time
double m_maxChangeDistance
maximum distance per trace step
std::vector< double > VecDbl
#define TS_ASSERT_DELTA_VECPT3D(a, b, delta)
void testUniqueTimeSteps()
Test behavior for unique timesteps.
double m_maxChangeVelocity
maximum change in velocity per trace step
void testInactiveCell()
2nd cell is inactive in the 2nd time step. // Thus it does not pull as hard. Also once the point reac...
BSHP< XmUGrid > m_ugrid
UGrid for the TracePoint operation.
Definition: XmGridTrace.cpp:97
void testVectorMultiplier()
test behavior of the vector multiplier
Contains the XmGridTrace Class and supporting data types.
static BSHP< XmUGrid > New()
double m_time1
#define XM_ZERO_TOL
double m_time2
time of the second time step
void testStartInactiveCell()
The point starts in an inactive cell, and doesnt move.
void testMaxTracingDistance()
test setting the max tracing distance
boost::dynamic_bitset< size_t > DynBitset
void testBeyondTimestep()
test behavior when starting beyond the second time step
BSHP< XmUGrid2dDataExtractor > m_extractor2x
data extractor for the x component for the second time step
void testBasicTracePoint()
test the basic functionality of trace point
DataLocationEnum
double m_maxTracingDistance
maximum distance for trace
#define XM_NODATA
bool EQ_TOL(const _T &A, const _U &B, const _V &tolerance)
void testDotProduct()
test the angle function
void testStartOutOfCell()
test behavior when starting point is out of the cell
virtual ~XmGridTrace()
Deconstruct XmGridTrace.
std::string m_exitMessage
exit message for the last TracePoint operation
void testMaxChangeVelocity()
Testing what happens when the maximum change in velocity is low / It reaches a point of high accelera...
XMS Namespace.
double m_minDeltaTime
minimum time per trace step
static BSHP< XmGridTrace > New(BSHP< XmUGrid > a_ugrid)
Construct XmGridTrace for a UGrid.
Traces points in an XmUGrid following a vector dataset.
Definition: XmGridTrace.h:39
#define XM_PI
double Mdist(_T x1, _U y1, _V x2, _W y2)
void testBeforeTimestep()
test the behavior when starting before the first timestep
void testSmallScalarsTracePoint()
test with small scalars to create more points
double m_maxChangeDirectionInRadians
maxmium change in direction per trace step
void testMaxChangeDistance()
Speed is limited to .25. It doesnt reach the edge because it goes below min delta time...
double m_vectorMultiplier
multiplier for all vectors in grid
Definition: XmGridTrace.cpp:98
BSHP< XmUGrid2dDataExtractor > m_extractor2y
data extractor for the y component for the second time step
BSHP< XmUGrid2dDataExtractor > m_extractor1y
data extractor for the y component for the first time step
double m_maxTracingTime
maximum time for trace
Definition: XmGridTrace.cpp:99
#define TS_ASSERT_DELTA_VEC(a, b, delta)
std::vector< Pt3d > VecPt3d