ITK  4.3.0
Insight Segmentation and Registration Toolkit
itkRankHistogram.h
Go to the documentation of this file.
1 /*=========================================================================
2  *
3  * Copyright Insight Software Consortium
4  *
5  * Licensed under the Apache License, Version 2.0 (the "License");
6  * you may not use this file except in compliance with the License.
7  * You may obtain a copy of the License at
8  *
9  * http://www.apache.org/licenses/LICENSE-2.0.txt
10  *
11  * Unless required by applicable law or agreed to in writing, software
12  * distributed under the License is distributed on an "AS IS" BASIS,
13  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14  * See the License for the specific language governing permissions and
15  * limitations under the License.
16  *
17  *=========================================================================*/
18 // histogram from the moving histogram operations
19 #ifndef __itkRankHistogram_h
20 #define __itkRankHistogram_h
21 
22 #include "itkIntTypes.h"
23 #include "itkNumericTraits.h"
24 
25 #include <map>
26 #include <vector>
27 
28 namespace itk
29 {
30 namespace Function
31 {
32 // a simple histogram class hierarchy. One subclass will be maps, the
33 // other vectors.
34 // This version is intended for keeping track of arbitrary ranks. It is
35 // based on the code from consolidatedMorphology.
36 //
37 // Support for different TCompare hasn't been tested, and shouldn't be
38 // necessary for the rank filters.
39 //
40 // This is a modified version for use with masks. Need to allow for
41 // the situation in which the map is empty
42 
43 /*
44  *
45  * This code was contributed in the Insight Journal paper:
46  * "Efficient implementation of kernel filtering"
47  * by Beare R., Lehmann G
48  * http://hdl.handle.net/1926/555
49  * http://www.insight-journal.org/browse/publication/160
50  *
51  */
52 
53 
54 template< class TInputPixel >
56 {
57 public:
58 
59  typedef std::less< TInputPixel > TCompare;
60 
62  {
63  m_Rank = 0.5;
64  m_Below = m_Entries = 0;
65  // can't set m_RankIt until something has been put in the histogram
66  m_Initialized = false;
69  {
71  }
72  else
73  {
75  }
77  m_RankIt = m_Map.begin(); // equivalent to setting to the initial value
78  }
79 
81  {}
82 
84  {
85  m_Map = hist.m_Map;
86  m_Rank = hist.m_Rank;
87  m_Below = hist.m_Below;
88  m_Entries = hist.m_Entries;
89  m_InitVal = hist.m_InitVal;
90  m_RankValue = hist.m_RankValue;
92  if ( m_Initialized )
93  {
94  m_RankIt = m_Map.find(m_RankValue);
95  }
96  return *this;
97  }
98 
99  void AddPixel(const TInputPixel & p)
100  {
101  m_Map[p]++;
102  if ( !m_Initialized )
103  {
104  m_Initialized = true;
105  m_RankIt = m_Map.begin();
106  m_Entries = m_Below = 0;
107  m_RankValue = p;
108  }
109  if ( m_Compare(p, m_RankValue) || p == m_RankValue )
110  {
111  ++m_Below;
112  }
113  ++m_Entries;
114  }
115 
116  void RemovePixel(const TInputPixel & p)
117  {
118  m_Map[p]--;
119  if ( m_Compare(p, m_RankValue) || p == m_RankValue )
120  {
121  --m_Below;
122  }
123  --m_Entries;
124  // this is the change that makes this version less efficient. The
125  // simplest approach I can think of with maps, though
126  if ( m_Entries <= 0 )
127  {
128  m_Initialized = false;
129  m_Below = 0;
130  m_Map.clear();
131  }
132  }
133 
134  bool IsValid()
135  {
136  return m_Initialized;
137  }
138 
139  TInputPixel GetValueBruteForce()
140  {
141  SizeValueType count = 0;
142  SizeValueType target = (int)( m_Rank * ( m_Entries - 1 ) ) + 1;
143  for( typename MapType::iterator it=m_Map.begin(); it != m_Map.end(); it++ )
144  {
145  count += it->second;
146  if( count >= target )
147  {
148  return it->first;
149  }
150  }
152  }
153 
154  TInputPixel GetValue(const TInputPixel &)
155  {
156  SizeValueType target = (SizeValueType)( m_Rank * ( m_Entries - 1 ) ) + 1;
157  SizeValueType total = m_Below;
158  SizeValueType ThisBin;
159  bool eraseFlag = false;
160 
161  if ( total < target )
162  {
163  typename MapType::iterator searchIt = m_RankIt;
164  typename MapType::iterator eraseIt;
165 
166  while ( searchIt != m_Map.end() )
167  {
168  // cleaning up the map - probably a better way of organising
169  // the loop. Currently makes sure that the search iterator is
170  // incremented before deleting
171  ++searchIt;
172  ThisBin = searchIt->second;
173  total += ThisBin;
174  if ( eraseFlag )
175  {
176  m_Map.erase(eraseIt);
177  eraseFlag = false;
178  }
179  if ( ThisBin <= 0 )
180  {
181  eraseFlag = true;
182  eraseIt = searchIt;
183  }
184  if ( total >= target )
185  {
186  break;
187  }
188  }
189  m_RankValue = searchIt->first;
190  m_RankIt = searchIt;
191  }
192  else
193  {
194  typename MapType::iterator searchIt = m_RankIt;
195  typename MapType::iterator eraseIt;
196 
197  while ( searchIt != m_Map.begin() )
198  {
199  ThisBin = searchIt->second;
200  unsigned int tbelow = total - ThisBin;
201  if ( tbelow < target ) // we've overshot
202  {
203  break;
204  }
205  if ( eraseFlag )
206  {
207  m_Map.erase(eraseIt);
208  eraseFlag = false;
209  }
210  if ( ThisBin <= 0 )
211  {
212  eraseIt = searchIt;
213  eraseFlag = true;
214  }
215  total = tbelow;
216 
217  --searchIt;
218  }
219  m_RankValue = searchIt->first;
220  m_RankIt = searchIt;
221  }
222 
223  m_Below = total;
224  itkAssertInDebugAndIgnoreInReleaseMacro( m_RankValue == GetValueBruteForce() );
225  return ( m_RankValue );
226  }
227 
228  void SetRank(float rank)
229  {
230  m_Rank = rank;
231  }
232 
233  void AddBoundary(){}
234 
235  void RemoveBoundary(){}
236 
238  {
239  return false;
240  }
241 
242 protected:
243  float m_Rank;
244 
245 private:
246  typedef typename std::map< TInputPixel, SizeValueType, TCompare > MapType;
247 
251  TInputPixel m_RankValue;
252  TInputPixel m_InitVal;
255  // This iterator will point at the desired rank value
256  typename MapType::iterator m_RankIt;
257 };
258 
259 
260 template< class TInputPixel >
262 {
263 public:
264  typedef std::less< TInputPixel > TCompare;
265 
267  {
269  m_Vec.resize(m_Size, 0);
272  {
274  }
275  else
276  {
278  }
279  m_Entries = m_Below = 0;
281  m_Rank = 0.5;
282  }
283 
285 
286  bool IsValid()
287  {
288  return m_Entries > 0;
289  }
290 
291  TInputPixel GetValueBruteForce()
292  {
293  SizeValueType count = 0;
294  SizeValueType target = (SizeValueType)( m_Rank * ( m_Entries - 1 ) ) + 1;
295  for( SizeValueType i=0; i<m_Size; i++ )
296  {
297  count += m_Vec[i];
298  if( count >= target )
299  {
301  }
302  }
304  }
305 
306  TInputPixel GetValue(const TInputPixel &)
307  {
308  return GetValueBruteForce();
309  }
310 
311  void AddPixel(const TInputPixel & p)
312  {
314 
315  m_Vec[q]++;
316  if ( m_Compare(p, m_RankValue) || p == m_RankValue )
317  {
318  ++m_Below;
319  }
320  ++m_Entries;
321  }
322 
323  void RemovePixel(const TInputPixel & p)
324  {
326 
327  itkAssertInDebugAndIgnoreInReleaseMacro( q >= 0 );
328  itkAssertInDebugAndIgnoreInReleaseMacro( q < (int)m_Vec.size() );
329  itkAssertInDebugAndIgnoreInReleaseMacro( m_Entries >= 1 );
330  itkAssertInDebugAndIgnoreInReleaseMacro( m_Vec[q] > 0 );
331 
332  m_Vec[q]--;
333  --m_Entries;
334 
335  if ( m_Compare(p, m_RankValue) || p == m_RankValue )
336  {
337  --m_Below;
338  }
339  }
340 
341  void SetRank(float rank)
342  {
343  m_Rank = rank;
344  }
345 
346  void AddBoundary(){}
347 
348  void RemoveBoundary(){}
349 
351  {
352  return true;
353  }
354 
355 protected:
356  float m_Rank;
357 
358 private:
359  typedef typename std::vector< SizeValueType > VecType;
360 
364  TInputPixel m_RankValue;
365  TInputPixel m_InitVal;
366  int m_Below;
368 };
369 
370 // now create MorphologicalGradientHistogram specilizations using the VectorMorphologicalGradientHistogram
371 // as base class
372 
375 template<>
376 class RankHistogram<unsigned char>:
377  public VectorRankHistogram<unsigned char>
378 {
379 };
380 
381 template<>
382 class RankHistogram<signed char>:
383  public VectorRankHistogram<signed char>
384 {
385 };
386 
387 template<>
388 class RankHistogram<bool>:
389  public VectorRankHistogram<bool>
390 {
391 };
392 
395 } // end namespace Function
396 } // end namespace itk
397 #endif
398