Developer documentation
Version 3.0.3-105-gd3941f44
fixel_map.h
Go to the documentation of this file.
1/* Copyright (c) 2008-2022 the MRtrix3 contributors.
2 *
3 * This Source Code Form is subject to the terms of the Mozilla Public
4 * License, v. 2.0. If a copy of the MPL was not distributed with this
5 * file, You can obtain one at http://mozilla.org/MPL/2.0/.
6 *
7 * Covered Software is provided under this License on an "as is"
8 * basis, without warranty of any kind, either expressed, implied, or
9 * statutory, including, without limitation, warranties that the
10 * Covered Software is free of defects, merchantable, fit for a
11 * particular purpose or non-infringing.
12 * See the Mozilla Public License v. 2.0 for more details.
13 *
14 * For more details, see http://www.mrtrix.org/.
15 */
16
17#ifndef __dwi_fixel_map_h__
18#define __dwi_fixel_map_h__
19
20
21#include "header.h"
22#include "image.h"
23#include "algo/loop.h"
24#include "dwi/fmls.h"
25#include "dwi/directions/mask.h"
26
27
28namespace MR
29{
30 namespace DWI
31 {
32
33
34
35 template <class Fixel>
36 class Fixel_map { MEMALIGN(Fixel_map<Fixel>)
37
38 public:
39 Fixel_map (const Header& H) :
40 _header (H),
41 _accessor (Image<MapVoxel*>::scratch (_header, "fixel map voxels"))
42 {
43 for (auto l = Loop() (_accessor); l; ++l)
44 _accessor.value() = nullptr;
45 // fixels[0] is an invalid fixel, as provided by the relevant empty constructor
46 // This allows index 0 to be used as an error code, simplifying the implementation of MapVoxel and Iterator
47 fixels.push_back (Fixel());
48 }
49 Fixel_map (const Fixel_map&) = delete;
50
51 class MapVoxel;
53
54 virtual ~Fixel_map()
55 {
56 for (auto l = Loop() (_accessor); l; ++l) {
57 if (_accessor.value()) {
58 delete _accessor.value();
59 _accessor.value() = nullptr;
60 }
61 }
62 }
63
64
65 class Iterator;
66 class ConstIterator;
67
68 Iterator begin (VoxelAccessor& v) { return Iterator (v.value(), *this); }
69 ConstIterator begin (VoxelAccessor& v) const { return ConstIterator (v.value(), *this); }
70
71 Fixel& operator[] (const size_t i) { return fixels[i]; }
72 const Fixel& operator[] (const size_t i) const { return fixels[i]; }
73
74 virtual bool operator() (const FMLS::FOD_lobes& in);
75
76 // Functions can copy-construct their own voxel accessor from this and retain const-ness:
77 const VoxelAccessor& accessor() const { return _accessor; }
78
79 const ::MR::Header& header() const { return _header; }
80
81 protected:
83
84 private:
85 const class HeaderHelper : public ::MR::Header { MEMALIGN(HeaderHelper)
86 public:
87 HeaderHelper (const ::MR::Header& H) :
89 {
90 ndim() = 3;
91 }
92 HeaderHelper() = default;
93 } _header;
94
95 VoxelAccessor _accessor;
96
97 };
98
99
100
101
102 template <class Fixel>
103 class Fixel_map<Fixel>::MapVoxel
104 { MEMALIGN(Fixel_map<Fixel>)
105 public:
106 MapVoxel (const FMLS::FOD_lobes& in, const size_t first) :
107 first_fixel_index (first),
108 count (in.size()),
109 lookup_table (new uint8_t[in.lut.size()])
110 {
111 memcpy (lookup_table, &in.lut[0], in.lut.size() * sizeof (uint8_t));
112 }
113
114 MapVoxel (const size_t first, const size_t size) :
115 first_fixel_index (first),
116 count (size),
117 lookup_table (nullptr) { }
118
119 ~MapVoxel()
120 {
121 if (lookup_table) {
122 delete[] lookup_table;
123 lookup_table = nullptr;
124 }
125 }
126
127 size_t first_index() const { return first_fixel_index; }
128 size_t num_fixels() const { return count; }
129 bool empty() const { return !count; }
130
131 // Direction must have been assigned to a histogram bin first
132 size_t dir2fixel (const size_t dir) const
133 {
134 assert (lookup_table);
135 const size_t offset = lookup_table[dir];
136 return ((offset == count) ? 0 : (first_fixel_index + offset));
137 }
138
139
140 private:
141 size_t first_fixel_index, count;
142 uint8_t* lookup_table;
143
144 };
145
146
147
148 template <class Fixel>
150 friend class Fixel_map<Fixel>::ConstIterator;
151 public:
152 Iterator (const MapVoxel* const voxel, Fixel_map<Fixel>& parent) :
153 index (voxel ? voxel->first_index() : 0),
154 last (voxel ? (index + voxel->num_fixels()) : 0),
155 fixel_map (parent) { }
156 Iterator (const Iterator&) = default;
157 Iterator& operator++ () { ++index; return *this; }
158 Fixel& operator() () const { return fixel_map.fixels[index]; }
159 operator bool () const { return (index != last); }
160 operator size_t () const { return index; }
161 private:
162 size_t index, last;
163 Fixel_map<Fixel>& fixel_map;
164 };
165
166 template <class Fixel>
168 public:
169 ConstIterator (const MapVoxel* const voxel, const Fixel_map& parent) :
170 index (voxel ? voxel->first_index() : 0),
171 last (voxel ? (index + voxel->num_fixels()) : 0),
172 fixel_map (parent) { }
173 ConstIterator (const ConstIterator&) = default;
174 ConstIterator& operator++ () { ++index; return *this; }
175 const Fixel& operator() () const { return fixel_map.fixels[index]; }
176 operator bool () const { return (index != last); }
177 operator size_t () const { return index; }
178 private:
179 size_t index, last;
180 const Fixel_map<Fixel>& fixel_map;
181 };
182
183
184
185
186 template <class Fixel>
188 {
189 if (in.empty())
190 return true;
191 auto v = accessor();
192 assign_pos_of (in.vox).to (v);
193 if (is_out_of_bounds (v))
194 return false;
195 if (v.value())
196 throw Exception ("FIXME: FOD_map has received multiple segmentations for the same voxel!");
197 v.value() = new MapVoxel (in, fixels.size());
198 for (const auto& i : in)
199 fixels.push_back (Fixel (i));
200 return true;
201 }
202
203
204
205
206 }
207}
208
209
210
211#endif
ConstIterator(const MapVoxel *const voxel, const Fixel_map &parent)
Definition: fixel_map.h:169
ConstIterator(const ConstIterator &)=default
Iterator(const Iterator &)=default
vector< Fixel > fixels
Definition: fixel_map.h:82
functions and classes related to image data input/output
Definition: image.h:41
index_type offset
Definition: loop.h:33
FORCE_INLINE LoopAlongAxes Loop()
Definition: loop.h:419
#define NOMEMALIGN
Definition: memory.h:22
Definition: base.h:24
Eigen::MatrixXd H
size_t index
#define MEMALIGN(...)
Definition: types.h:185