source: Revenant/geoworld/include/fluidyn.h@ 7ef8ec4

port/mars-tycoon
Last change on this file since 7ef8ec4 was 80a6a52, checked in by Jonathan Neufeld <support@…>, 3 years ago

Get to a compile state for terrain procedural generation
  • Property mode set to 100644
File size: 32.2 KB
Line 
1#pragma once
2
3#include <limits>
4
5#include "linalgadapt.h"
6
7#include <mars_calc.h>
8#include <mars_splines.h>
9#include <mars_geometry.h>
10#include <mars_util.h>
11#include <mars_ptr.h>
12#include <mars_streams.h>
13
14namespace fldyn
15{
16 template <class Derived, typename DP, typename HMPrec>
17 class IQueryBridge
18 {
19 public:
20 typedef DP DensityPrecision;
21
22 inline DensityPrecision density (const unsigned int x, const unsigned int y) const
23 { return static_cast <const Derived *> (this) -> density(x, y, 0); }
24
25 inline DensityPrecision density (const unsigned int x, const unsigned int y, const unsigned int z) const
26 { return static_cast <const Derived *> (this) -> density(x, y, z); }
27
28 inline HMPrec elevation (const unsigned int x, const unsigned int y) const
29 { return static_cast <const Derived *> (this) -> elevation(x, y); }
30 };
31
32 class Weights
33 {
34 public:
35 float
36 fwStraightness,
37 fwElevationInfluence,
38 fwDensityInfluence,
39 fwChaosInfluence,
40 fWeightSum;
41
42 inline Weights () : fwStraightness(0), fwElevationInfluence(0), fwDensityInfluence(0), fwChaosInfluence(0), fWeightSum(0) {}
43 inline Weights (
44 const float fwStraightness, const float fwElevationInfluence,
45 const float fwDensityInfluence, const float fwChaosInfluence
46 ) : fwStraightness(fwStraightness), fwElevationInfluence(fwElevationInfluence),
47 fwDensityInfluence(fwDensityInfluence), fwChaosInfluence(fwChaosInfluence),
48 fWeightSum(fwStraightness + fwElevationInfluence + fwDensityInfluence + fwChaosInfluence)
49 {}
50 };
51 static inline mars::ObjectStream & operator >> (mars::ObjectStream & ins, Weights & weights)
52 {
53 return ins >> weights.fwStraightness >> weights.fwElevationInfluence >> weights.fwDensityInfluence >> weights.fwChaosInfluence >> weights.fWeightSum;
54 }
55 static inline mars::ObjectStream & operator << (mars::ObjectStream & outs, const Weights & weights)
56 {
57 return outs << weights.fwStraightness << weights.fwElevationInfluence << weights.fwDensityInfluence << weights.fwChaosInfluence << weights.fWeightSum;
58 }
59
60 class FluiDynEx : public std::exception
61 {
62 };
63
64 template <class QB, typename T, typename H>
65 class SplineCPCandidate;
66
67 template <class QB, typename T, typename H>
68 class GuidedSplineCPCandidate;
69
70 template <class QB, typename T, typename H>
71 class FlowLine;
72
73 template <typename QB, typename T, typename H> // QB, Control point type, spline precision type
74 struct FluiDynTraits
75 {
76 typedef QB QueryBridgeType; // Bridge pattern for querying the world
77 typedef typename QB::DensityPrecision QBDensityPrec; // Precision used for sampling density
78 typedef SplineCPCandidate<QB, T, H> BasicCP; // Control-point type
79 typedef GuidedSplineCPCandidate<QB, T, H> GuidedCP;
80 typedef typename VectorTag< T >::V2 Vector2;
81 typedef typename VectorTag< T >::V3 Vector3;
82 typedef typename VectorTag< T >::PC PolarCoords;
83 typedef typename VectorTag< T >::CC CylindricalCoords;
84 typedef T Precision; // Precision of the n-D vector type previously mentioned
85 typedef H HMPrec;
86 typedef mars::Spline< mars::CubicBSplineBase <Vector2, Precision> > SplineType; // BBox computation depends on the Cubic B-Spline algorithm
87 typedef mars::TreeGuideNode< FlowLine <QB, T, H>, Vector3, T > GuideType; // Guide type used as the base
88 typedef FlowLine<QB, T, H> FlowLineType;
89
90 using SphericalCoords = typename mars::SphericalCoords< T >;
91 };
92
93 template< typename T >
94 struct StartPoint
95 {
96 typename VectorTag< T >::V3 position;
97 T direction;
98
99 StartPoint ()
100 : direction(0) {}
101 StartPoint (typename VectorTag< T >::V3 position, T direction)
102 : position(position), direction(direction) {}
103 };
104 template< typename T >
105 mars::ObjectStream & operator >> (mars::ObjectStream & ins, StartPoint< T > & sp)
106 {
107 return ins >> sp.position >> sp.direction;
108 }
109
110 template< typename T >
111 mars::ObjectStream & operator << (mars::ObjectStream & outs, const StartPoint< T > & sp)
112 {
113 return outs << sp.position << sp.direction;
114 }
115
116 template <class QB, typename T, typename H> // QB, Algebraic vector type
117 class SplineCPCandidate : public FluiDynTraits<QB, T, H>::Vector3, public FluiDynTraits<QB, T, H>
118 {
119 private:
120
121 template <typename J>
122 inline J nonzero (const J n) const
123 {
124 if (n == 0)
125 return std::numeric_limits< J >::epsilon();
126 else
127 return n;
128 }
129
130 protected:
131 using typename FluiDynTraits<QB, T, H>::QBDensityPrec;
132 using typename FluiDynTraits<QB, T, H>::HMPrec;
133 using typename FluiDynTraits<QB, T, H>::Vector2;
134 using typename FluiDynTraits<QB, T, H>::Vector3;
135 using typename FluiDynTraits<QB, T, H>::CylindricalCoords;
136 using typename FluiDynTraits<QB, T, H>::GuideType;
137 using typename FluiDynTraits<QB, T, H>::FlowLineType;
138 using typename FluiDynTraits<QB, T, H>::BasicCP;
139 using typename FluiDynTraits<QB, T, H>::GuidedCP;
140 using typename FluiDynTraits<QB, T, H>::PolarCoords;
141
142 QBDensityPrec _maxDens;
143 HMPrec _maxElev, _minElev;
144 float _fDot;
145
146 public:
147 inline explicit SplineCPCandidate (const Vector3 & p)
148 : Vector3(p), _maxDens(0), _minElev(0), _maxElev(0), _fDot(0) {}
149
150 inline SplineCPCandidate (
151 const Vector3 & p1, const Vector3 & p0, const Vector3 & dp0,
152 const mars::ptr< QB > & query
153 ) : Vector3(p1)
154 {
155 const Vector3 &
156 dp = p1 - p0,
157 pM = dp / 2 + p0;
158
159 // TODO: Check for negative/signed values of p1 and pM
160 _maxDens = std::max(
161 query->density(static_cast <unsigned int> (p1.x), static_cast <unsigned int> (p1.y)),
162 query->density(static_cast <unsigned int> (pM.x), static_cast <unsigned int> (pM.y))
163 );
164
165#ifdef _DEBUG
166 if (_maxDens == 0)
167 throw FluiDynEx();
168#endif
169
170 const unsigned short
171 nElevTerm = query->elevation(static_cast <unsigned int> (p1.x), static_cast <unsigned int> (p1.y)),
172 nElevMid = query->elevation(static_cast <unsigned int> (pM.x), static_cast <unsigned int> (pM.y));
173
174 if (nElevTerm > nElevMid)
175 {
176 _maxElev = nonzero(nElevTerm);
177 _minElev = nElevMid;
178 } else
179 {
180 _maxElev = nonzero(nElevMid);
181 _minElev = nElevTerm;
182 }
183
184 _fDot = dp0 * dp;
185#ifdef _DEBUG
186 if (_fDot == 0)
187 throw FluiDynEx();
188#endif
189 }
190
191 float compareTo (const SplineCPCandidate & other, const float fStuckFactor, const Weights & weights) const
192 {
193 return
194 static_cast <float> (_maxElev) / static_cast <float> (other._maxElev) * weights.fwElevationInfluence +
195 static_cast <float> (_maxDens / other._maxDens * weights.fwDensityInfluence) +
196 other._fDot / _fDot * weights.fwStraightness;
197 }
198
199 inline const QBDensityPrec & getDensity () const { return _maxDens; }
200 inline const HMPrec & getMaxElevation () const { return _maxElev; }
201 inline const HMPrec & getMinElevation () const { return _minElev; }
202 inline const float & getStraightness () const { return _fDot; }
203 };
204
205 template <class QB, typename T, typename H> // QB, Guide parent type, Algebraic vector type
206 class GuidedSplineCPCandidate : public SplineCPCandidate <QB, T, H>
207 {
208 protected:
209 using typename FluiDynTraits<QB, T, H>::QBDensityPrec;
210 using typename FluiDynTraits<QB, T, H>::HMPrec;
211 using typename FluiDynTraits<QB, T, H>::Vector2;
212 using typename FluiDynTraits<QB, T, H>::Vector3;
213 using typename FluiDynTraits<QB, T, H>::CylindricalCoords;
214 using typename FluiDynTraits<QB, T, H>::GuideType;
215 using typename FluiDynTraits<QB, T, H>::FlowLineType;
216 using typename FluiDynTraits<QB, T, H>::BasicCP;
217 using typename FluiDynTraits<QB, T, H>::GuidedCP;
218 using typename FluiDynTraits<QB, T, H>::PolarCoords;
219
220 typename GuideType::DistanceType _gdist;
221
222 public:
223 inline GuidedSplineCPCandidate (
224 const Vector3 & p1, const Vector3 & p0, const Vector3 & dp0, const typename GuideType::DistanceType & gdist,
225 const mars::ptr< QB > & query, const GuideType & guide
226 ) : SplineCPCandidate<QB, T, H> (p1, p0, dp0, query), _gdist(gdist) {}
227
228 float compareTo (const GuidedSplineCPCandidate & other, const float fStuckFactor, const Weights & weights) const
229 {
230 return
231 SplineCPCandidate<QB, T, H>::compareTo(other, fStuckFactor, weights) +
232 _gdist.distance / other._gdist.distance * weights.fwChaosInfluence * (fStuckFactor + 1); // XTODO: Boosts the anti-chaos influence by a straight factor of how many times the same guide segment was successively selected, this may be overly simplistic, see line below as well
233 }
234
235 const typename FlowLineType::SegmentType & getNearestGuideSegment () const { return _gdist.segment; }
236 const float & getGuideDistance () const { return _gdist; }
237 };
238
239 template <class QB, typename T, typename H> // QB, List of algebraic vector type, spline precision type
240 class FlowLine : public FluiDynTraits<QB, T, H>::GuideType, public FluiDynTraits<QB, T, H>
241 {
242 protected:
243 using typename FluiDynTraits<QB, T, H>::QBDensityPrec;
244 using typename FluiDynTraits<QB, T, H>::HMPrec;
245 using typename FluiDynTraits<QB, T, H>::Vector2;
246 using typename FluiDynTraits<QB, T, H>::Vector3;
247 using typename FluiDynTraits<QB, T, H>::CylindricalCoords;
248 using typename FluiDynTraits<QB, T, H>::GuideType;
249 using typename FluiDynTraits<QB, T, H>::FlowLineType;
250 using typename FluiDynTraits<QB, T, H>::BasicCP;
251 using typename FluiDynTraits<QB, T, H>::GuidedCP;
252 using typename FluiDynTraits<QB, T, H>::PolarCoords;
253 using typename FluiDynTraits<QB, T, H>::SplineType;
254 using typename FluiDynTraits<QB, T, H>::Precision;
255 using typename GuideType::DistanceType;
256
257 public:
258 template <typename R, typename E1, typename E2>
259 class Iterator : public mars::SplineZoneIterator<
260 SplineType,
261 E2, E1,
262 R, // The raster precision
263 typename VectorTag<R>::V2, // The raster version of the output vector type
264 typename VectorTag<Precision>::V3 // Used for "unormal"
265 >
266 {
267 private:
268 const FlowLineType & _parent;
269
270 inline Iterator (const FlowLineType & parent, const E1 & fnEaseIn, const E2 & fnEaseOut)
271 : mars::SplineZoneIterator<SplineType, E2, E1, R, typename VectorTag<R>::V2, typename VectorTag<Precision>::V3 > (parent.getSpline(), parent.radius, fnEaseOut, fnEaseIn, parent.unormal), _parent (parent) {}
272
273 friend FlowLineType;
274 public:
275
276 inline typename SplineType::VectorType stress () const
277 { return (this->laplace() - this->gradient()) / this->seglen(); }
278 };
279
280 public:
281 const unsigned int radius;
282 const typename VectorTag<Precision>::V3 unormal;
283
284 FlowLine (const typename GuideType::PointList & points, const unsigned int radius, const typename VectorTag<Precision>::V3 & unormal)
285 : GuideType (points), _spline (SplineType::template createFrom <typename GuideType::PointList> (points.begin(), points.end())), radius(radius), unormal(unormal) {}
286
287 FlowLine (const GuideType & parent, const typename GuideType::PointList & points, const unsigned int radius, const typename VectorTag<Precision>::V3 & unormal)
288 : GuideType (points, parent), _spline (SplineType::template createFrom <typename GuideType::PointList> (points.begin(), points.end())), radius(radius), unormal(unormal) {}
289
290 template <typename R, typename E1, typename E2>
291 inline Iterator< R, E1, E2 > iterate (const E1 & fnEaseIn = E1(), const E2 & fnEaseOut = E2()) const
292 {
293 return Iterator< R, E1, E2 > (*this, fnEaseIn, fnEaseOut);
294 }
295
296 // Computes the distance between this flow-line and an arbitrary point
297 inline float sdist (const Vector2 & p) const
298 {
299 DistanceType dist = *this - p;
300 return _spline.compute(dist.segment.index(), dist.time());
301 }
302
303 inline typename SplineType::real getLength () const { return _spline.len(); }
304 inline const SplineType & getSpline () const { return _spline; }
305
306 void offsetBy (const Vector3 & p)
307 {
308 GuideType::offsetBy(p);
309 _spline += p;
310 }
311
312 private:
313 SplineType _spline;
314 };
315
316 template <class QB, typename T, typename H>
317 class PathFinder : public FluiDynTraits< QB, T, H >
318 {
319 public:
320 using typename FluiDynTraits<QB, T, H>::FlowLineType;
321
322 private:
323 using typename FluiDynTraits<QB, T, H>::BasicCP;
324 using typename FluiDynTraits<QB, T, H>::GuidedCP;
325 using typename FluiDynTraits<QB, T, H>::Vector2;
326 using typename FluiDynTraits<QB, T, H>::Vector3;
327 using typename FluiDynTraits<QB, T, H>::GuideType;
328 using typename FluiDynTraits<QB, T, H>::CylindricalCoords;
329 using typename FluiDynTraits<QB, T, H>::Precision;
330 using typename FluiDynTraits<QB, T, H>::SphericalCoords;
331
332 typedef std::list <BasicCP> CandidatesList;
333 typedef std::list <GuidedCP> GuidedCandidatesList;
334
335 float
336 _frSnapToGuide,
337 _fMaxSegLenSQ,
338 _fGravitySlopeGrace;
339
340 mars::RangeX< float >
341 _mmfSplineSeg,
342 _mmfBendZenith;
343
344 unsigned int
345 _nCandidateSampleAmt,
346 _nBendZenithSteps;
347
348 mars::BBox< T >
349 _bbox;
350
351 Vector3 _vun3;
352 Weights _weights;
353
354 mars::ptr< QB > _query;
355
356 static Vector3 createUnitZ () { return Vector3(0, 0, 1); }
357
358 inline BasicCP createSplineCPCandidate(
359 const Vector3 & p1, const Vector3 & p0, const Vector3 & dp0
360 ) const { return BasicCP(p1, p0, dp0, _query); }
361
362 inline GuidedCP createGuidedSplineCPCandidate (
363 const Vector3 & p1, const Vector3 & p0, const Vector3 & dp0, const typename GuideType::DistanceType & gdist, const GuideType & guide
364 ) const { return GuidedCP(p1, p0, dp0, gdist, _query, guide); }
365
366 template <typename CL> typename CL::const_iterator findBest (const CL & candidates, const unsigned int nStuckCount) const
367 {
368 typename CL::const_iterator io = candidates.begin();
369 const float fStuckFactor = static_cast <float> (nStuckCount);
370
371 for (typename CL::const_iterator i = candidates.begin(); i != candidates.end(); ++i)
372 {
373 if (i->compareTo(*io, fStuckFactor, _weights) / (_weights.fWeightSum + fStuckFactor) < 1.0) // Factoring in the stuck factor
374 {
375 io = i;
376 }
377 }
378 return io;
379 }
380
381 template <typename CaDList>
382 inline void candidatesReset (CaDList & candidates, float & zenith) const
383 {
384 using namespace mars;
385
386 candidates.clear();
387 zenith = _mmfBendZenith.minimum;
388 }
389 inline void initCandidateSearch (
390 Vector3 & pCurr,
391 CylindricalCoords & ccDir, CylindricalCoords & ccDir0, CylindricalCoords & ccTest,
392 typename GuideType::PointList & points
393 ) const
394 {
395 using namespace mars;
396
397 ccDir.z =
398 ccDir0.z =
399 pCurr.z = _query->elevation(
400 static_cast< unsigned int > (mars::RNDi (pCurr.x)),
401 static_cast< unsigned int > (mars::RNDi (pCurr.y))
402 );
403
404 ccDir0.azimuth =
405 ccTest.azimuth = ccDir.azimuth;
406
407 points.push_back(BasicCP(pCurr));
408 }
409
410 inline bool computeCandidate (
411 const float zenith, const float pstep, const Vector3 & pCurr,
412 CylindricalCoords & ccTest, CylindricalCoords & ccDir0, Vector3 & pTest
413 ) const
414 {
415 using namespace mars;
416
417 // Calculate a candidate position at some random point within an angular threshold of the preset flow direction
418 ccTest.azimuth = ccDir0.azimuth + RANDf(zenith) - zenith / 2.0f;
419 ccTest.p = pstep + _mmfSplineSeg.minimum;
420
421 // Add the candidate point to the current point plus the density vector for a test
422 // Actually what does this do anyway again?
423 //+ _query->density(static_cast <HMPrec> (posTest.x), static_cast <HMPrec> (posTest.y));
424
425 pTest = pCurr + ccTest;
426
427 pTest.z =
428 ccTest.z = static_cast< float > (
429 _query->elevation(
430 mars::RNDi(pTest.x), // TODO: Verify use of RNDi
431 mars::RNDi(pTest.y) // TODO: Verify use of RNDi
432 )
433 );
434
435 return
436 static_cast< Precision > (pTest.z - pCurr.z)
437 / ccTest.p
438 < _fGravitySlopeGrace;
439 }
440
441 inline bool isBBoxApproved (const T x, const T y, const T fRadius) const
442 { return _bbox.inside(x - fRadius, y - fRadius) && _bbox.inside(x + fRadius, y + fRadius); }
443
444 void pathFind( const Vector2 & p, const float fDir, typename GuideType::PointList & points, const T fRadius ) const
445 {
446 using namespace mars;
447
448 CandidatesList candidates;
449 typename CandidatesList::const_iterator iCandidate;
450 Vector3 posTest, pCurr;
451 CylindricalCoords
452 ccDir (1, fDir, 0),
453 ccTest, // TODO: Make dirLast a scalar cuz it's only for angle
454 ccDir0 = ccDir; // Assign a non-zero magnitude initializer to avoid division by zero problems later
455 // Calculate maximum segment length
456 const float pstep = static_cast< float > (_mmfSplineSeg.delta);
457
458 pCurr.x = p.x;
459 pCurr.y = p.y;
460 initCandidateSearch(pCurr, ccDir, ccDir0, ccTest, points);
461
462 float zenith;
463 unsigned short iz;
464 do
465 {
466 candidatesReset(candidates, zenith);
467
468 for (iz = 0; iz < _nBendZenithSteps && candidates.size() < _nCandidateSampleAmt; ++iz, zenith += _mmfBendZenith.step)
469 {
470 if (computeCandidate(zenith, pstep, pCurr, ccTest, ccDir0, posTest) && // Returns false when elevation slopes upwards
471 isBBoxApproved(posTest.x, posTest.y, fRadius)/* && (posTest.x != pCurr.x || posTest.y != pCurr.y) */) // ccTest.p > 0 assumption guarantees the previous check is unnecessary
472 {
473 candidates.push_back(createSplineCPCandidate(posTest, pCurr, ccDir0));
474 }
475 }
476
477 iCandidate = findBest (candidates, 0);
478
479 if (iCandidate != candidates.end())
480 {
481 const BasicCP & best = *iCandidate;
482 points.push_back(best);
483 ccDir0 = static_cast< Vector3 > (best - pCurr);
484 pCurr = best;
485 }
486 } while (candidates.size() >= _nCandidateSampleAmt); // XTODO: Superficially tested: Tests that the distance between the current position and the end of the guide is still large enough
487 }
488
489 void pathFind( const FlowLineType & fl, typename GuideType::PointList & points, const T fRadius ) const
490 {
491 using namespace mars;
492 using PolarCoords = typename VectorTag< T >::PC;
493
494 GuidedCandidatesList candidates;
495 typename GuidedCandidatesList::const_iterator iCandidate;
496 Vector3 posTest, pCurr, pDir = fl.back() - fl.front();
497 CylindricalCoords
498 ccDir = pDir,
499 ccTest, // TODO: Make dirLast a scalar cuz it's only for angle
500 ccDirPrev = ccDir; // Assign a non-zero magnitude initializer to avoid division by zero problems later
501 typename FlowLineType::SegmentType igseg = fl.begin(); // Tracks the last accepted guide segment head
502 unsigned int nStuckCount = 0;
503 bool bExitCondDeviantChild = false;
504
505 pCurr = _bbox.template clampXY< typename Vector3::Precision > (fl.front() + (mars::U(pDir) & _vun3 * RANDf(std::min(_mmfSplineSeg.maximum, MAG(fl.grad(fl.begin()))))));
506 initCandidateSearch(pCurr, ccDir, ccDirPrev, ccTest, points);
507
508 unsigned short iz;
509 float zenith;
510 do
511 {
512 // Get the gradient for the current guide segment
513 const PolarCoords ccGuideSegDir = static_cast<Vector2> (fl.grad(igseg));
514 // Calculate maximum segment length relating to the current guide segment
515 const float pstep = std::min(static_cast< float > (_mmfSplineSeg.delta), ccGuideSegDir.p - _mmfSplineSeg.minimum);
516 const float pstepSQ = SQ(pstep);
517
518 // TODO: Optimize angle calculation
519 //ccGuideSegDir.azimuth = PolarCoords (guide.tail(igseg) - pCurr).azimuth;
520
521 candidatesReset <GuidedCandidatesList> (candidates, zenith);
522
523 for (iz = 0; iz < _nBendZenithSteps && candidates.size() < _nCandidateSampleAmt; ++iz, zenith += _mmfBendZenith.step)
524 {
525 if (computeCandidate(zenith, pstep, pCurr, ccTest, ccDirPrev, posTest) && // Returns false when elevation slopes upwards
526 isBBoxApproved(posTest.x, posTest.y, fRadius)
527 ) /*&& (posTest.x != pCurr.x || posTest.y != pCurr.y) (See above)*/
528 {
529 // Calculate median point between candidate position and source position
530 typename GuideType::DistanceType gdist = fl.distance (posTest, igseg);
531 const Vector3 pTail = fl.tail(gdist.segment);
532
533 if (MAGSQ(posTest - pTail) / pstepSQ < _frSnapToGuide)
534 posTest = pTail;
535
536 candidates.push_back(
537 createGuidedSplineCPCandidate(posTest, pCurr, ccDirPrev, gdist, fl)
538 );
539 }
540 }
541
542 iCandidate = findBest (candidates, nStuckCount);
543
544 if (iCandidate != candidates.end())
545 {
546 const GuidedCP & best = *iCandidate;
547 points.push_back(best);
548 ccDirPrev = static_cast< Vector3 > (best - pCurr);
549 pCurr = best;
550
551 // If the nearest guide segment was the same as last time, then increase the stuck count
552 // otherwise the guide segment was incremented and we're not stuck
553 if (igseg != best.getNearestGuideSegment())
554 {
555 igseg = best.getNearestGuideSegment();
556 nStuckCount = 0;
557 } else
558 {
559 ++nStuckCount;
560 }
561 }
562 bExitCondDeviantChild = MAGSQ(static_cast <Vector3> (pCurr - fl.back())) > _fMaxSegLenSQ;
563 } while (bExitCondDeviantChild && candidates.size() >= _nCandidateSampleAmt); // XTODO: Superficially tested: Tests that the distance between the current position and the end of the guide is still large enough
564
565 if (!bExitCondDeviantChild)
566 points.push_back(fl.back()); // HACK: Just add the guide's last point as this last point in the path find
567 }
568
569 public:
570 PathFinder ()
571 : _frSnapToGuide(0), _fMaxSegLenSQ(0), _fGravitySlopeGrace(0), _nCandidateSampleAmt(0), _nBendZenithSteps(0), _vun3(createUnitZ())
572 {}
573
574 PathFinder (
575 const mars::RangeX< float > & mmfBendZenith,
576 const unsigned int nBendZenithSteps,
577 const unsigned int nCandidateSampleAmt,
578 const mars::RangeX< float > & mmfSplineSeg,
579 const float fGravitySlopeGrace,
580 const float frSnapToGuide,
581 const float fwStraightness,
582 const float fwElevationInfluence,
583 const float fwDensityInfluence,
584 const float fwChaosInfluence
585 ) : _mmfBendZenith(mars::RangeX< float > (mmfBendZenith.minimum, mmfBendZenith.maximum, static_cast< float > (mmfBendZenith.delta / static_cast <float> (nBendZenithSteps)))),
586 _nBendZenithSteps(nBendZenithSteps), _nCandidateSampleAmt (nCandidateSampleAmt), _mmfSplineSeg(mmfSplineSeg), _fGravitySlopeGrace(fGravitySlopeGrace),
587 _frSnapToGuide(frSnapToGuide),
588 _weights(fwStraightness, fwElevationInfluence, fwDensityInfluence, fwChaosInfluence),
589 _fMaxSegLenSQ (mars::SQ(mmfSplineSeg.maximum)), _vun3(createUnitZ())
590 {}
591
592 inline const mars::RangeX< float > & getBendZenith () const { return _mmfBendZenith; }
593 inline unsigned int getBendZenithSteps () const { return _nBendZenithSteps; }
594 inline unsigned int getCandidateSampling () const { return _nCandidateSampleAmt; }
595 inline const mars::RangeX< float > & getSplineSegment () const { return _mmfSplineSeg; }
596 inline float getGravitySlopeGrace () const { return _fGravitySlopeGrace; }
597 inline float getSnapToGuide () const { return _frSnapToGuide; }
598 inline const Weights & getWeights () const { return _weights; }
599
600 inline FlowLineType * find (const Vector2 & p, const float fDir, const unsigned int nRadius) const
601 {
602 typename GuideType::PointList points;
603
604 pathFind(p, fDir, points, static_cast< T > (nRadius));
605 if (points.size() < 2)
606 return NULL;
607 else
608 return new FlowLineType(points, nRadius, _vun3);
609 }
610
611 inline FlowLineType * find (const FlowLineType & guide, const unsigned int nRadius) const
612 {
613 typename GuideType::PointList points;
614
615 pathFind(guide, points, static_cast< T > (nRadius));
616 if (points.size() < 2)
617 return NULL;
618 else
619 return new FlowLineType(guide, points, nRadius, _vun3);
620 }
621
622 inline void setBBoxLimit (const mars::BBox< T > & bbox)
623 {
624 _bbox = bbox;
625 }
626
627 inline void setQueryCallback (const mars::ptr< QB > & pqb)
628 {
629 _query = pqb;
630 }
631
632 inline mars::ObjectStream & operator >> (mars::ObjectStream & outs) const
633 {
634 return outs << _mmfBendZenith << _nBendZenithSteps << _nCandidateSampleAmt << _mmfSplineSeg << _fGravitySlopeGrace << _frSnapToGuide << _weights;
635 }
636 inline mars::ObjectStream & operator << (mars::ObjectStream & ins)
637 {
638 ins >> _mmfBendZenith >> _nBendZenithSteps >> _nCandidateSampleAmt >> _mmfSplineSeg >> _fGravitySlopeGrace >> _frSnapToGuide >> _weights;
639 _fMaxSegLenSQ = mars::SQ(_mmfSplineSeg.maximum);
640 return ins;
641 }
642 };
643
644 template <class QB, typename T, typename H, typename DS>
645 class ChannelSet : public FluiDynTraits< QB, T, H >
646 {
647 public:
648 using typename FluiDynTraits<QB, T, H>::FlowLineType;
649 using typename FluiDynTraits<QB, T, H>::Precision;
650 using typename FluiDynTraits<QB, T, H>::GuideType;
651 using typename FluiDynTraits<QB, T, H>::Vector2;
652 using typename FluiDynTraits<QB, T, H>::Vector3;
653
654 typedef std::list < mars::ptr< FlowLineType > > GuideList;
655 typedef PathFinder< QB, T, H > MyPathFinder;
656 typedef std::list < StartPoint< T > > StartingPoints;
657
658 protected:
659 virtual void applyFlowLine (const FlowLineType & fl, DS * pData) const = 0;
660 virtual void finishedTier (const GuideList & guides, DS * pData) const {}
661
662 private:
663 mars::RangeX< float >
664 _mmfFlowlineWidth;
665
666 mars::RangeX< unsigned int>
667 _mmnKidsPerLevel;
668
669 StartingPoints _starts;
670 MyPathFinder _finder;
671 GuideList _guides;
672 mars::BBox< Precision > _bbox;
673
674 void computeBBox ()
675 {
676 using namespace mars;
677
678 GuideList tier1, tier2,
679 * pTierSrc = &tier1,
680 * pTierKids = &tier2;
681
682 *pTierSrc = _guides;
683 _bbox = BBox< Precision >();
684 while (!pTierSrc->empty())
685 {
686 pTierKids->clear();
687 for (typename GuideList::iterator i = pTierSrc->begin(); i != pTierSrc->end(); ++i)
688 {
689 for (typename GuideType::SegmentType k = (*i)->begin(); k != (*i)->end(); ++k)
690 {
691 const Vector2
692 v1 = *k - static_cast< float > ((*i)->radius),
693 v2 = *k + static_cast< float > ((*i)->radius);
694
695 if (_bbox.empty())
696 _bbox = BBox< Precision > (v1, v2);
697 else
698 {
699 _bbox.add(v1);
700 _bbox.add(v2);
701 }
702 }
703 for (typename GuideType::KidList::iterator j = (*i)->kidsBegin(); j != (*i)->kidsEnd(); ++j)
704 pTierKids->push_back(*j);
705 }
706 std::swap(pTierKids, pTierSrc); // Flip the guide lists
707 }
708 }
709
710 // Traverses up the ancestry of flow-lines starting at the specified one and returns the ancestral flowline
711 // that is closest to the specified point based on an approximated spline calculation
712 const FlowLineType * closestFlowLine (const Vector3 & p, const FlowLineType & flBegin)
713 {
714 const Precision fMin = std::numeric_limits<Precision>::max;
715 const FlowLineType * pfl = &flBegin, *pflResult = NULL;
716
717 do
718 {
719 const Precision fDist = pfl->sdist(p);
720
721 if (fDist < fMin)
722 {
723 fMin = fDist;
724 pflResult = pfl;
725 }
726 } while ((pfl = pfl->parent()) != NULL);
727
728 return pflResult;
729 }
730
731 public:
732 ChannelSet() {}
733 ChannelSet(
734 const mars::RangeX< float > & mmfBendZenith,
735 const unsigned int nBendZenithSteps,
736 const unsigned int nCandidateSampleAmt,
737 const mars::RangeX< float > & mmfSplineSeg,
738 const float fGravitySlopeGrace,
739 const float frSnapToGuide,
740 const mars::RangeX< float > & mmfFlowlineWidth,
741 const mars::RangeX< unsigned int > & mmnKidsPerLevel,
742 const float fwStraightness,
743 const float fwElevationInfluence,
744 const float fwDensityInfluence,
745 const float fwChaosInfluence
746 )
747 : _finder(mmfBendZenith, nBendZenithSteps, nCandidateSampleAmt, mmfSplineSeg, fGravitySlopeGrace, frSnapToGuide, fwStraightness, fwElevationInfluence, fwDensityInfluence, fwChaosInfluence),
748 _mmfFlowlineWidth (mmfFlowlineWidth), _mmnKidsPerLevel(mmnKidsPerLevel) {}
749
750 inline const mars::RangeX< float > & getWidthRange() const { return _mmfFlowlineWidth; }
751 inline const mars::RangeX< unsigned int > & getKidsPerLevel () const { return _mmnKidsPerLevel; }
752 inline const MyPathFinder & getPathFinder () const { return _finder; }
753
754 inline void addStartPoint (const StartPoint< Precision > & spt)
755 { _starts.push_back(spt); }
756 inline void addStartPoint (const Vector3 & pos, const Precision fDir)
757 { addStartPoint(StartPoint< Precision >(pos, fDir)); }
758
759 // Query is not referenced beyond this call
760 void init (const mars::ptr< QB > & query, const mars::BBox< Precision > & bbox, const unsigned int nLevels)
761 {
762 using namespace mars;
763
764 _finder.setQueryCallback(query);
765 _guides.clear();
766
767 GuideList tier1, tier2,
768 * pTierSrc = &tier1,
769 * pTierKids = &tier2;
770
771 _finder.setBBoxLimit(bbox);
772
773 for (typename StartingPoints::const_iterator i = _starts.begin(); i != _starts.end(); ++i)
774 {
775 // HACK: No longer calculating radius based on configuration or level depth
776 const unsigned int nRandRadius = // HACK: Redundapendency in doLevel
777 static_cast <unsigned int> (_mmfFlowlineWidth.next());
778 //const unsigned int nRandRadius = mars::RAND(6) + 5;
779
780 mars::ptr< FlowLineType > pfl = _finder.find(i->position, i->direction, nRandRadius);
781
782 if (pfl != NULL)
783 {
784 pTierSrc->push_back (pfl);
785 _guides.push_back(pfl);
786 }
787 }
788
789 for (unsigned int l = nLevels; l > 0 && !pTierSrc->empty(); --l)
790 {
791 pTierKids->clear();
792 for (typename GuideList::iterator i = pTierSrc->begin(); i != pTierSrc->end(); ++i)
793 {
794 const unsigned int nIterations = _mmnKidsPerLevel.next();
795
796 // First step is to find paths for all the subsequent flow-lines at this order of magnitude
797 for (unsigned int k = 0; k < nIterations; ++k)
798 {
799 ptr< FlowLineType > flChild = _finder.find(*(*i), static_cast< unsigned int > (_mmfFlowlineWidth.next()));
800
801 if (flChild != NULL)
802 (*i)->addChild(
803 flChild
804 );
805 }
806
807 // Second step is to add sub-sequent flows to the next level-down list
808 pTierKids->insert(pTierKids->end(), (*i)->kidsBegin(), (*i)->kidsEnd());
809 }
810 std::swap(pTierKids, pTierSrc); // Flip the guide lists
811 }
812
813 computeBBox();
814
815 _finder.setQueryCallback(NULL);
816 }
817
818 void run (DS * pData) const
819 {
820 using namespace mars;
821
822 GuideList tier1, tier2,
823 * pTierSrc = &tier1,
824 * pTierKids = &tier2;
825
826 *pTierSrc = _guides;
827 while (!pTierSrc->empty())
828 {
829 pTierKids->clear();
830 for (typename GuideList::iterator i = pTierSrc->begin(); i != pTierSrc->end(); ++i)
831 {
832 applyFlowLine(*(*i), pData);
833 for (typename GuideType::KidList::iterator j = (*i)->kidsBegin(); j != (*i)->kidsEnd(); ++j)
834 pTierKids->push_back(*j);
835 }
836 finishedTier (*pTierSrc, pData);
837
838 std::swap(pTierKids, pTierSrc); // Flip the guide lists
839 }
840 }
841
842 mars::BBox< Precision > getBBox () const { return _bbox; }
843
844 inline typename GuideList::const_iterator begin() const { return _guides.begin(); }
845 inline typename GuideList::const_iterator end() const { return _guides.end(); }
846 inline typename GuideList::iterator begin() { return _guides.begin(); }
847 inline typename GuideList::iterator end() { return _guides.end(); }
848 inline size_t count () const { return _guides.size(); }
849
850 mars::ObjectStream & operator >> (mars::ObjectStream & outs) const
851 {
852 _finder >> outs;
853 outs << _mmfFlowlineWidth << _mmnKidsPerLevel << _starts << _bbox;
854
855 GuideList tier1(_guides), tier2,
856 * pTierSrc = &tier1,
857 * pTierKids = &tier2;
858
859 do
860 {
861 outs << static_cast< size_t > (pTierSrc->size());
862 pTierKids->clear();
863
864 for (typename GuideList::const_iterator i = pTierSrc->begin(); i != pTierSrc->end(); ++i)
865 {
866 outs << (*i)->radius << (*i)->unormal << (*i)->points() << (*i)->kidsCount();
867 pTierKids->insert(pTierKids->end(), (*i)->kidsBegin(), (*i)->kidsEnd());
868 }
869
870 std::swap(pTierKids, pTierSrc);
871 } while(!pTierSrc->empty());
872
873 return outs;
874 }
875
876 mars::ObjectStream & operator << (mars::ObjectStream & ins)
877 {
878 using namespace mars;
879
880 _finder << ins;
881 ins >> _mmfFlowlineWidth >> _mmnKidsPerLevel >> _starts >> _bbox;
882
883 struct GuideInfo
884 {
885 FlowLineType * parent, * line;
886 size_t childcount;
887
888 GuideInfo ()
889 : parent(NULL), childcount(0) {}
890 GuideInfo (FlowLineType * pParent)
891 : parent(pParent), childcount(0) {}
892
893 ObjectStream & operator << (ObjectStream & ins)
894 {
895 unsigned int nRadius;
896 Vector3 unorm;
897 typename GuideType::PointList points;
898
899 ins >> nRadius >> unorm >> points >> childcount;
900 if (parent == NULL)
901 {
902 line = new FlowLineType(points, nRadius, unorm);
903 }
904 else
905 {
906 line = new FlowLineType(*parent, points, nRadius, unorm);
907 parent->addChild(line);
908 }
909 return ins;
910 }
911 };
912
913 typedef std::vector< GuideInfo > GuideInfoList;
914
915 GuideInfoList
916 tier1, tier2,
917 * pSrcParentGuides = &tier1,
918 * pKidGuides = &tier2;
919
920 size_t nCount;
921 GuideInfo info;
922 size_t j;
923
924 ins >> nCount;
925 _guides.clear();
926 for (size_t i = 0; i < nCount; ++i)
927 {
928 info << ins;
929 pSrcParentGuides->push_back(info);
930 _guides.push_back(info.line);
931 }
932
933 do
934 {
935 pKidGuides->clear();
936 for (typename GuideInfoList::iterator i = pSrcParentGuides->begin(); i != pSrcParentGuides->end(); ++i)
937 {
938 GuideInfo info (i->line);
939
940 for (j = 0; j < i->childcount; ++j)
941 {
942 info << ins;
943 pKidGuides->push_back(info);
944 }
945 }
946
947 std::swap(pSrcParentGuides, pKidGuides);
948
949 } while (!pSrcParentGuides->empty());
950
951 return ins;
952 }
953 };
954}
Note: See TracBrowser for help on using the repository browser.