OpenWalnut  1.3.1
WThreadedTrackingFunction.h
1 //---------------------------------------------------------------------------
2 //
3 // Project: OpenWalnut ( http://www.openwalnut.org )
4 //
5 // Copyright 2009 OpenWalnut Community, BSV@Uni-Leipzig and CNCF@MPI-CBS
6 // For more information see http://www.openwalnut.org/copying
7 //
8 // This file is part of OpenWalnut.
9 //
10 // OpenWalnut is free software: you can redistribute it and/or modify
11 // it under the terms of the GNU Lesser General Public License as published by
12 // the Free Software Foundation, either version 3 of the License, or
13 // (at your option) any later version.
14 //
15 // OpenWalnut is distributed in the hope that it will be useful,
16 // but WITHOUT ANY WARRANTY; without even the implied warranty of
17 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
18 // GNU Lesser General Public License for more details.
19 //
20 // You should have received a copy of the GNU Lesser General Public License
21 // along with OpenWalnut. If not, see <http://www.gnu.org/licenses/>.
22 //
23 //---------------------------------------------------------------------------
24 
25 #ifndef WTHREADEDTRACKINGFUNCTION_H
26 #define WTHREADEDTRACKINGFUNCTION_H
27 
28 #include <vector>
29 #include <utility>
30 
31 #include <boost/array.hpp>
32 
33 #include "../common/math/linearAlgebra/WLinearAlgebra.h"
34 #include "../common/WSharedObject.h"
35 #include "../common/WThreadedJobs.h"
36 
37 #include "WDataSetSingle.h"
38 
39 class WThreadedTrackingFunctionTest; //! forward declaration
40 
41 
42 /**
43  * For tracking related functionality.
44  */
45 namespace wtracking // note that this is not final
46 {
47  // an epsilon value for various floating point comparisons
48 #define TRACKING_EPS 0.0000001
49 
50  /**
51  * \class WTrackingUtility
52  *
53  * A class that provides untility functions and typedefs for tracking algorithms.
54  */
56  {
57  public:
58  //! define a job type for tracking algorithms
59  typedef std::pair< WVector3d, WVector3d > JobType;
60 
61  //! the dataset type
63 
64  //! a pointer to a dataset
65  typedef boost::shared_ptr< DataSetType const > DataSetPtr;
66 
67  //! a function that calculates a direction to continue tracking
68  typedef boost::function< WVector3d ( DataSetPtr, JobType const& ) > DirFunc;
69 
70  //! a pointer to a regular 3d grid
71  // other grid types are not supported at the moment
72  typedef boost::shared_ptr< WGridRegular3D > Grid3DPtr;
73 
74  /**
75  * A function that follows a direction until leaving the current voxel.
76  *
77  * \param dataset A pointer to the input dataset.
78  * \param job A pair of vectors, the position and the direction of the last integration.
79  * \param dirFunc A function that computes the next direction.
80  *
81  * \return true, iff the calculated point is a valid position inside the grid
82  */
83  static bool followToNextVoxel( DataSetPtr dataset, JobType& job, DirFunc const& dirFunc );
84 
85  // one could add a runge-kutta-integrator too
86 
87  /**
88  * Check if a point is on the boundary of the given grid, where boundary
89  * means a distance less then TRACKING_EPS from any plane between
90  * voxels. This does not check if the position is actually inside the grid.
91  *
92  * \param grid The grid.
93  * \param pos The position to test.
94  *
95  * \return true, iff the position is on any voxel boundary
96  */
97  static bool onBoundary( Grid3DPtr grid, WVector3d const& pos );
98 
99  /**
100  * Calculate the distance from a given position to the nearest voxel boundary
101  * on the ray from the position along the given direction.
102  *
103  * \param grid The grid.
104  * \param pos The starting position of the ray.
105  * \param dir The normalized direction of the ray.
106  *
107  * \return The distance to the next voxel boundary.
108  *
109  * \note pos + getDistanceToBoundary( grid, pos, dir ) * dir will be a position on a voxel boundary
110  */
111  static double getDistanceToBoundary( Grid3DPtr grid, WVector3d const& pos, WVector3d const& dir );
112  };
113 
114  //////////////////////////////////////////////////////////////////////////////////////////
115 
116  /**
117  * \class WThreadedTrackingFunction
118  *
119  * Implements a generalized multithreaded tracking algorithm. A function that calculates the direction
120  * and a function that calculates a new position have to be provided.
121  *
122  * Output values can be retrieved via two visitor functions that get called per fiber tracked and
123  * per point calculated respectively.
124  *
125  * There are a certain number n of seeds per direction, this meens n*n*n seeds per voxel. For every
126  * seed, m fibers get integrated. These two parameters are the seedPositions and seedsPerVoxel parameters
127  * of the constructor, respectively.
128  *
129  * A 'cubic' region of the grid can be chosen for seeding. The v0 and v1 parameters of the constructor
130  * are the starting/target voxel coords. Example:
131  *
132  * v0: 1, 1, 1
133  * v1: 4, 5, 3
134  *
135  * In this case, only voxels between coords 1 to 3 in the x-direction, the voxels 1 to 4 in y- and the voxels 1 to 2
136  * in z-direction are used for seeding.
137  *
138  * Note that voxels at the first (0) and last (grid->getNbCoords*()) position in any direction are
139  * invalid seeding voxels as they are partially outside of the grid.
140  */
141  class WThreadedTrackingFunction : public WThreadedJobs< WTrackingUtility::DataSetType, WTrackingUtility::JobType >
142  {
143  //! make the test a friend
144  friend class ::WThreadedTrackingFunctionTest;
145 
146  //! the job type
148 
149  //! the dataset type
151 
152  //! a pointer to a dataset
154 
155  //! the grid type
157 
158  //! a pointer to the grid
159  typedef boost::shared_ptr< GridType > GridPtr;
160 
161  //! the direction calculation function
163 
164  //! the path integration function
165  typedef boost::function< bool ( DataSetPtr, JobType&, DirFunc const& ) > NextPositionFunc;
166 
167  //! a visitor function for fibers
168  typedef boost::function< void ( std::vector< WVector3d > const& ) > FiberVisitorFunc;
169 
170  //! a visitor function type for points
171  typedef boost::function< void ( WVector3d const& ) > PointVisitorFunc;
172 
173  //! the base class, a threaded job function
175 
176  //! this type
178 
179  public:
180  /**
181  * Constructor.
182  *
183  * \param dataset A pointer to a dataset.
184  * \param dirFunc A direction calculation function.
185  * \param nextFunc A position integration function.
186  * \param fiberVst A visitor for fibers.
187  * \param pointVst A visitor for points.
188  * \param seedPositions The number of seed positions in every direction per voxel.
189  * \param seedsPerPos The number of fibers startet from every seed position.
190  * \param v0 A vector of starting voxel indices for every direction.
191  * \param v1 A vector of target voxel indices for every direction.
192  */
193  WThreadedTrackingFunction( DataSetPtr dataset, DirFunc dirFunc, NextPositionFunc nextFunc,
194  FiberVisitorFunc fiberVst, PointVisitorFunc pointVst,
195  std::size_t seedPositions = 1, std::size_t seedsPerPos = 1,
196  std::vector< int > v0 = std::vector< int >(),
197  std::vector< int > v1 = std::vector< int >() );
198 
199  /**
200  * Destructor.
201  */
202  virtual ~WThreadedTrackingFunction();
203 
204  /**
205  * The job generator.
206  *
207  * \param job The next job (output).
208  *
209  * \return false, iff there are no more jobs.
210  */
211  virtual bool getJob( JobType& job ); // NOLINT
212 
213  /**
214  * The calculation per job.
215  *
216  * \param input The input dataset.
217  * \param job The job.
218  */
219  virtual void compute( DataSetPtr input, JobType const& job );
220 
221  private:
222  /**
223  * \class IndexType
224  *
225  * An index for seed positions.
226  */
227  class IndexType
228  {
229  friend class ::WThreadedTrackingFunctionTest;
230  public:
231  /**
232  * Construct an invalid index.
233  */
234  IndexType();
235 
236  /**
237  * Construct an index.
238  *
239  * \param grid The grid.
240  * \param v0 A vector of starting voxel indices for every direction.
241  * \param v1 A vector of target voxel indices for every direction.
242  * \param seedPositions The number of seed positions in every direction per voxel.
243  * \param seedsPerPosition The number of fibers startet from every seed position.
244  */
245  IndexType( GridPtr grid, std::vector< int > const& v0,
246  std::vector< int > const& v1, std::size_t seedPositions,
247  std::size_t seedsPerPosition );
248 
249  /**
250  * Increase the index by one, effectively generating the next seed position.
251  *
252  * \return *this
253  */
255 
256  /**
257  * Check if there aren't any more seed positions.
258  *
259  * \return true, iff there aren't any more seed positions.
260  */
261  bool done();
262 
263  /**
264  * Create a job from this index.
265  *
266  * \return The job that is the current position.
267  */
268  JobType job();
269 
270  private:
271  //! a pointer to the grid
273 
274  //! true, iff there are no more seeds
275  bool m_done;
276 
277  //! the position in the seed space
278  boost::array< std::size_t, 4 > m_pos;
279 
280  //! the minimum position in the seed space
281  boost::array< std::size_t, 4 > m_min;
282 
283  //! the maximum position in the seed space
284  boost::array< std::size_t, 4 > m_max;
285 
286  //! the relative (to the size of a voxel) distance between seeds
287  double m_offset;
288  };
289 
290  //! a pointer to the grid
292 
293  //! a function that returns the next direction
295 
296  //! a function that calculates the next position
298 
299  //! the fiber visitor
301 
302  //! the point visitor
304 
305  //! the maximum number of points per forward/backward integration of a fiber
306  std::size_t m_maxPoints;
307 
308  //! the current index/seed position
310  };
311 
312 } /* namespace wtracking */
313 
314 #endif // WTHREADEDTRACKINGFUNCTION_H