File: cmtkFilterMask.h

package info (click to toggle)
cmtk 3.3.1p2%2Bdfsg-4
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid
  • size: 10,524 kB
  • sloc: cpp: 87,098; ansic: 23,347; sh: 3,896; xml: 1,551; perl: 707; makefile: 334
file content (213 lines) | stat: -rw-r--r-- 5,938 bytes parent folder | download | duplicates (5)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
/*
//
//  Copyright 1997-2009 Torsten Rohlfing
//
//  Copyright 2004-2010, 2013 SRI International
//
//  This file is part of the Computational Morphometry Toolkit.
//
//  http://www.nitrc.org/projects/cmtk/
//
//  The Computational Morphometry Toolkit is free software: you can
//  redistribute it and/or modify it under the terms of the GNU General Public
//  License as published by the Free Software Foundation, either version 3 of
//  the License, or (at your option) any later version.
//
//  The Computational Morphometry Toolkit is distributed in the hope that it
//  will be useful, but WITHOUT ANY WARRANTY; without even the implied
//  warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//  GNU General Public License for more details.
//
//  You should have received a copy of the GNU General Public License along
//  with the Computational Morphometry Toolkit.  If not, see
//  <http://www.gnu.org/licenses/>.
//
//  $Revision: 5436 $
//
//  $LastChangedDate: 2018-12-10 19:01:20 -0800 (Mon, 10 Dec 2018) $
//
//  $LastChangedBy: torstenrohlfing $
//
*/

#ifndef __cmtkFilterMask_h_included_
#define __cmtkFilterMask_h_included_

#include <cmtkconfig.h>

#include <Base/cmtkTypes.h>
#include <Base/cmtkUnits.h>
#include <Base/cmtkFixedVector.h>
#include <Base/cmtkMathUtil.h>

#include <System/cmtkSmartPtr.h>

#include <math.h>
#include <vector>

namespace
cmtk
{

/** \addtogroup Base */
//@{

/** Filter mask pixel entry.
 * This class handles a single entry in a pre-computed filter mask for
 * multidimensional images with relative coordinates and filter coeffiecient 
 * for a single pixel.
 */
template<int DIM>
class FilterMaskPixel
{
public:
  /// This class.
  typedef FilterMaskPixel Self;

  /// Smart pointer.
  typedef SmartPointer<Self> SmartPtr;

  /// Default constructor.
  FilterMaskPixel() {}

  /// Explicit constructor.
  FilterMaskPixel
  ( const FixedVector<DIM,int>& location, const int relativeIndex, const Types::DataItem coefficient ) 
    : Location( location ),
      RelativeIndex ( relativeIndex ),
      Coefficient( coefficient ),
      PixelIndex( 0 ),
      Valid( false )
  {}

  /// Relative location of this pixel.
  FixedVector<DIM,int> Location;

  /// Relative index of this pixel in source image from center of kernel.
  int RelativeIndex;

  /// Filter coefficient.
  Types::DataItem Coefficient;

  /// Cached source index of image pixel for this element.
  int PixelIndex;

  /** Active flag.
   * This flag can be used in client code to flag pixels in the filter mask as
   * either valid (i.e., inside domain) or invalid (outsside domain). This
   * eliminates repeated range checking in cases where there are multiple
   * iterations over the filter mask at one location.
   */
  bool Valid;
};

/** Filter mask.
 * This class handles pre-computed filter masks for multidimensional images
 * with relative pixel coordinates and filter coeffiecients.
 */
template<int DIM>
class FilterMask : 
  /// Inherit from STL container.
  public std::vector< FilterMaskPixel<DIM> >
{
public:
  /// This class type.
  typedef FilterMask<DIM> Self;

  /// Direct parent class.
  typedef std::vector< FilterMaskPixel<DIM> > Superclass;
  
  /// Smert pointer to this class.
  typedef SmartPointer<Self> SmartPtr;

  /// Constructor.
  template<class F>
  FilterMask( const FixedVector<DIM,int>& dims, const FixedVector<DIM,Types::Coordinate>& deltas, const Types::Coordinate radius, F filter ) 
  {
    FixedVector<DIM,int> pixel;
    FixedVector<DIM,int> width;
    FixedVector<DIM,Types::Coordinate> position;

    for ( int dim = 0; dim < DIM; ++dim ) 
      {
      pixel[dim] = - (width[dim] = 1+static_cast<int>( radius / deltas[dim] ));
      position[dim] = pixel[dim] * deltas[dim];
      }

    bool done = false;
    while ( ! done ) 
      {
      // increment the DIM-digit pixel index counter including overflow.
      for ( int dim = 0; dim < DIM; ++dim ) 
	{
	++pixel[dim];
	if ( pixel[dim] <= width[dim] ) 
	  {
	  // no overflow, leave for loop since we're done
	  dim = DIM;
	  } 
	else
	  {
	  if ( dim+1 == DIM ) 
	    // was this the last dimension? if so, leave while() loop
	    done = true;
	  else 
	    { 
	    // no, then reset this dimension and repeat loop to increment next
	    pixel[dim] = -width[dim];
	    }
	  }
	}
      // are we done with the kernel?
      if ( ! done ) 
	{
	// no, then compute Euclidean distance from center
	Types::Coordinate distance = 0.0;
	for ( int dim = 0; dim < DIM; ++dim ) 
	  {
	  position[dim] = pixel[dim] * deltas[dim];
	  distance += position[dim] * position[dim];
	  }
	distance = sqrt( distance );
	// if distance is within radius then add a pixel to the filter mask
	if ( distance < radius ) 
	  {
	  const int index = pixel[0] + dims[0] * (pixel[1] + dims[1] * pixel[2] );
	  this->push_back( FilterMaskPixel<DIM>( pixel, index, filter( position ) ) );
	  }
	}
      }
  }
  
  /// Gaussian filter as an example of a concrete filter implementation.
  class Gaussian 
  {
  public:
    /// Constructor.
    Gaussian( const Units::GaussianSigma& standardDeviation ) 
    {
      InvStandardDeviation = 1.0 / standardDeviation.Value();
      NormFactor = 1.0 / (sqrt(2.0 * M_PI) * standardDeviation.Value());
    }
    
    /// Get filter coefficient at relative location from filter center.
    Types::DataItem operator() ( const FixedVector<DIM,Types::Coordinate>& relativePosition ) 
    {
      Types::Coordinate distance = 0;
      for ( int i = 0; i < DIM; ++i ) 
	distance += relativePosition[i] * relativePosition[i];
      return static_cast<Types::DataItem>( NormFactor * exp( -distance * MathUtil::Square( InvStandardDeviation ) / 2 ) );
    }
    
  private:
    /// Standard deviation.
    Types::Coordinate InvStandardDeviation;

    /// Gaussian normalization factor.
    Types::Coordinate NormFactor;
  };
};

} // namespace  cmtk

#endif // #ifndef __cmtkFilterMask_h_included_