A.I.M.S algorithms


convol.h
Go to the documentation of this file.
1 /* This software and supporting documentation are distributed by
2  * Institut Federatif de Recherche 49
3  * CEA/NeuroSpin, Batiment 145,
4  * 91191 Gif-sur-Yvette cedex
5  * France
6  *
7  * This software is governed by the CeCILL-B license under
8  * French law and abiding by the rules of distribution of free software.
9  * You can use, modify and/or redistribute the software under the
10  * terms of the CeCILL-B license as circulated by CEA, CNRS
11  * and INRIA at the following URL "http://www.cecill.info".
12  *
13  * As a counterpart to the access to the source code and rights to copy,
14  * modify and redistribute granted by the license, users are provided only
15  * with a limited warranty and the software's author, the holder of the
16  * economic rights, and the successive licensors have only limited
17  * liability.
18  *
19  * In this respect, the user's attention is drawn to the risks associated
20  * with loading, using, modifying and/or developing or reproducing the
21  * software by the user in light of its specific status of free software,
22  * that may mean that it is complicated to manipulate, and that also
23  * therefore means that it is reserved for developers and experienced
24  * professionals having in-depth computer knowledge. Users are therefore
25  * encouraged to load and test the software's suitability as regards their
26  * requirements in conditions enabling the security of their systems and/or
27  * data to be ensured and, more generally, to use and operate it in the
28  * same conditions as regards security.
29  *
30  * The fact that you are presently reading this means that you have had
31  * knowledge of the CeCILL-B license and that you accept its terms.
32  */
33 
34 
35 #ifndef AIMS_SIGNALFILTER_CONVOL_H
36 #define AIMS_SIGNALFILTER_CONVOL_H
37 
39 #include <aims/data/data.h>
40 #include <cartobase/type/limits.h>
41 #include <float.h>
42 
43 
51 template <class T>
53 {
54 public:
55 
56  AimsConvolution(bool safe=true) : _safe(safe) {
57  set_safe_status(safe);
58  }
59  virtual ~AimsConvolution() { }
60 
62 
63  void set_safe_status(bool safe=true) {
64  if (safe)
66  else
68  _safe = safe;
69  }
70 
71 protected:
73  T doit_voxel(AimsData< T >& img, AimsData< T >& kernel,
74  int x, int y, int z, int t);
76  virtual T doit_voxel_safe(AimsData< T >& img, AimsData< T >& kernel,
77  int x, int y, int z, int t);
79  virtual T doit_voxel_unsafe(AimsData< T >& img, AimsData< T >& kernel,
80  int x, int y, int z, int t);
81  T min() const;
82  T max() const;
83 
85  AimsData< T >& kernel, int x, int y, int z, int t);
86 
87  bool _safe;
88 };
89 
97 template <class T>
99 {
100 public:
101 
103  bool safe=true) : AimsConvolution< T >(safe),
104  _sources(sources), _mask(mask) { }
107  virtual T doit_voxel_safe(AimsData< T >& img, AimsData< T >& kernel,
108  int x, int y, int z, int t);
110  virtual T doit_voxel_unsafe(AimsData< T >& img, AimsData< T >& kernel,
111  int x, int y, int z, int t);
112 
113 protected:
114  short _sources;
116 };
117 
118 
119 template< class T > inline
121  AimsData< T >& kernel )
122 {
123  ASSERT( kernel.dimT() == 1 );
124  if (_safe == false) ASSERT(img.borderWidth() >= 1);
125  else if (img.borderWidth() >= 1) set_safe_status(false);
126 
127  int x, y, z, t;
128 
129  int dx = img.dimX();
130  int dy = img.dimY();
131  int dz = img.dimZ();
132  int dt = img.dimT();
133 
134  AimsData< T > res( dx, dy, dz, dt );
135  res.setSizeXYZT( img.sizeX(), img.sizeY(), img.sizeZ(), img.sizeT() );
136 
137  for ( t=0; t<dt; t++ )
138  for ( z=0; z<dz; z++ )
139  for ( y=0; y<dy; y++ )
140  for ( x=0; x<dx; x++ )
141  {
142  T val = doit_voxel(img, kernel, x, y, z, t);
143  if ( val > max() ) val = max();
144  else if ( val < min() ) val = min();
145  res( x, y, z, t ) = val;
146  }
147 
148  return res;
149 }
150 
151 template< class T > inline
153  int x, int y, int z, int t)
154 {
155  return (this->*_doit_voxel_method)(img, kernel, x, y, z, t);
156 };
157 
158 template< class T > inline
160  AimsData< T >& kernel, int x, int y, int z, int t)
161 {
162  int i, j, k, idx, idy, idz;
163 
164  int mx = kernel.dimX();
165  int my = kernel.dimY();
166  int mz = kernel.dimZ();
167 
168  int mx2 = mx / 2;
169  int my2 = my / 2;
170  int mz2 = mz / 2;
171 
172  int dx = img.dimX();
173  int dy = img.dimY();
174  int dz = img.dimZ();
175 
176  T val = (T)0;
177 
178  for ( k=0; k<mz; k++ )
179  for ( j=0; j<my; j++ )
180  for ( i=0; i<mx; i++ )
181  {
182  idx = x - i + mx2;
183  idy = y - j + my2;
184  idz = z - k + mz2;
185 
186  if ( idx >= 0 && idy >= 0 && idz >= 0 &&
187  idx < dx && idy < dy && idz < dz )
188  val += kernel( i, j, k ) * img( idx, idy, idz, t );
189  }
190 
191  return val;
192 }
193 
194 template< class T > inline
196  AimsData< T >& kernel, int x, int y, int z, int t)
197 {
198  int i, j, k, idx, idy, idz;
199 
200  int mx = kernel.dimX();
201  int my = kernel.dimY();
202  int mz = kernel.dimZ();
203 
204  int mx2 = mx / 2;
205  int my2 = my / 2;
206  int mz2 = mz / 2;
207 
208  T val = (T)0;
209 
210  for ( k=0; k<mz; k++ )
211  for ( j=0; j<my; j++ )
212  for ( i=0; i<mx; i++ )
213  {
214  idx = x - i + mx2;
215  idy = y - j + my2;
216  idz = z - k + mz2;
217 
218  val += kernel( i, j, k ) * img( idx, idy, idz, t );
219  }
220 
221  return val;
222 }
223 
224 template< class T > inline
226  AimsData< T >& kernel, int x, int y, int z, int t)
227 {
228  if (_mask(x, y, z, t) == _sources) return img(x, y, z, t);
229  return AimsConvolution<T>::doit_voxel_safe(img, kernel, x, y, z, t);
230 }
231 
232 template< class T > inline
234  AimsData< T >& kernel, int x, int y, int z, int t)
235 {
236  if (_mask(x, y, z, t) == _sources) return img(x, y, z, t);
237  return AimsConvolution<T>::doit_voxel_unsafe(img, kernel, x, y, z, t);
238 }
239 
240 template <>
241 inline
243 {
245 }
246 
247 
248 template <>
249 inline
251 {
253 }
254 
255 
256 template <>
257 inline
259 {
260  return 0;
261 }
262 
263 
264 template <>
265 inline
267 {
269 }
270 
271 
272 template <>
273 inline
275 {
277 }
278 
279 
280 template <>
281 inline
283 {
285 }
286 
287 
288 template <>
289 inline
291 {
292  return 0;
293 }
294 
295 
296 template <>
297 inline
299 {
301 }
302 
303 
304 template <>
305 inline
307 {
309 }
310 
311 
312 template <>
313 inline
315 {
317 }
318 
319 
320 template <>
321 inline
323 {
324  return 0;
325 }
326 
327 
328 template <>
329 inline
331 {
333 }
334 
335 
336 template <>
337 inline
339 {
340  return (-FLT_MAX);
341 }
342 
343 
344 template <>
345 inline
347 {
348  return FLT_MAX;
349 }
350 
351 
352 template <>
353 inline
355 {
356  return DBL_MIN;
357 }
358 
359 
360 template <>
361 inline
363 {
364  return DBL_MAX;
365 }
366 
367 
368 // ### remove after everything has been moved to intN_t/uintN_t
369 #include <limits.h>
370 #if !defined(__sun__) || !defined(_CHAR_IS_SIGNED)
371 template<>
372 inline
374 {
376 }
377 
378 
379 template <>
380 inline
382 {
384 }
385 #endif
386 
387 
388 template <>
389 inline
391 {
392  return LONG_MIN;
393 }
394 
395 
396 template <>
397 inline
399 {
400  return LONG_MAX;
401 }
402 
403 
404 template <>
405 inline
407 {
408  return 0;
409 }
410 
411 
412 template <>
413 inline
415 {
416  return ULONG_MAX;
417 }
418 
419 
420 #endif
float min(float x, float y)
Definition: thickness.h:105
float max(float x, float y)
Definition: thickness.h:97
virtual T doit_voxel_safe(AimsData< T > &img, AimsData< T > &kernel, int x, int y, int z, int t)
called for each voxel (safe version)
Definition: convol.h:159
Make convolution only on a specified mask.
Definition: convol.h:98
int dimX() const
unsigned long ulong
float sizeY() const
int dimT() const
virtual ~AimsMaskedConvolution()
Definition: convol.h:105
T(AimsConvolution::* _doit_voxel_method)(AimsData< T > &img, AimsData< T > &kernel, int x, int y, int z, int t)
Definition: convol.h:84
virtual T doit_voxel_unsafe(AimsData< T > &img, AimsData< T > &kernel, int x, int y, int z, int t)
called for each voxel (unsafe version)
Definition: convol.h:233
The template class to make convolutions.
Definition: convol.h:52
int borderWidth() const
virtual T doit_voxel_safe(AimsData< T > &img, AimsData< T > &kernel, int x, int y, int z, int t)
called for each voxel (safe version)
Definition: convol.h:225
T doit_voxel(AimsData< T > &img, AimsData< T > &kernel, int x, int y, int z, int t)
called for each voxel (user selected function)
Definition: convol.h:152
static _Tp min()
AimsData< short > & _mask
Definition: convol.h:115
virtual T doit_voxel_unsafe(AimsData< T > &img, AimsData< T > &kernel, int x, int y, int z, int t)
called for each voxel (unsafe version)
Definition: convol.h:195
AimsConvolution(bool safe=true)
Definition: convol.h:56
float sizeT() const
void set_safe_status(bool safe=true)
Definition: convol.h:63
AimsMaskedConvolution(AimsData< short > &mask, short sources=0, bool safe=true)
Definition: convol.h:102
float sizeZ() const
virtual ~AimsConvolution()
Definition: convol.h:59
T max() const
T min() const
float sizeX() const
int dimZ() const
bool _safe
Definition: convol.h:87
void setSizeXYZT(float sizex=1.0f, float sizey=1.0f, float sizez=1.0f, float sizet=1.0f)
int dimY() const
#define ASSERT(EX)
AimsData< T > doit(AimsData< T > &, AimsData< T > &)
Definition: convol.h:120
static _Tp max()