ITK  4.4.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  if(this != &hist)
86  {
87  m_Map = hist.m_Map;
88  m_Rank = hist.m_Rank;
89  m_Below = hist.m_Below;
90  m_Entries = hist.m_Entries;
91  m_InitVal = hist.m_InitVal;
92  m_RankValue = hist.m_RankValue;
94  if ( m_Initialized )
95  {
96  m_RankIt = m_Map.find(m_RankValue);
97  }
98  }
99  return *this;
100  }
101 
102  void AddPixel(const TInputPixel & p)
103  {
104  m_Map[p]++;
105  if ( !m_Initialized )
106  {
107  m_Initialized = true;
108  m_RankIt = m_Map.begin();
109  m_Entries = m_Below = 0;
110  m_RankValue = p;
111  }
112  if ( m_Compare(p, m_RankValue) || p == m_RankValue )
113  {
114  ++m_Below;
115  }
116  ++m_Entries;
117  }
118 
119  void RemovePixel(const TInputPixel & p)
120  {
121  m_Map[p]--;
122  if ( m_Compare(p, m_RankValue) || p == m_RankValue )
123  {
124  --m_Below;
125  }
126  --m_Entries;
127  // this is the change that makes this version less efficient. The
128  // simplest approach I can think of with maps, though
129  if ( m_Entries <= 0 )
130  {
131  m_Initialized = false;
132  m_Below = 0;
133  m_Map.clear();
134  }
135  }
136 
137  bool IsValid()
138  {
139  return m_Initialized;
140  }
141 
142  TInputPixel GetValueBruteForce()
143  {
144  SizeValueType count = 0;
145  SizeValueType target = (int)( m_Rank * ( m_Entries - 1 ) ) + 1;
146  for( typename MapType::iterator it=m_Map.begin(); it != m_Map.end(); it++ )
147  {
148  count += it->second;
149  if( count >= target )
150  {
151  return it->first;
152  }
153  }
155  }
156 
157  TInputPixel GetValue(const TInputPixel &)
158  {
159  SizeValueType target = (SizeValueType)( m_Rank * ( m_Entries - 1 ) ) + 1;
160  SizeValueType total = m_Below;
161  SizeValueType ThisBin;
162  bool eraseFlag = false;
163 
164  if ( total < target )
165  {
166  typename MapType::iterator searchIt = m_RankIt;
167  typename MapType::iterator eraseIt;
168 
169  while ( searchIt != m_Map.end() )
170  {
171  // cleaning up the map - probably a better way of organising
172  // the loop. Currently makes sure that the search iterator is
173  // incremented before deleting
174  ++searchIt;
175  ThisBin = searchIt->second;
176  total += ThisBin;
177  if ( eraseFlag )
178  {
179  m_Map.erase(eraseIt);
180  eraseFlag = false;
181  }
182  if ( ThisBin <= 0 )
183  {
184  eraseFlag = true;
185  eraseIt = searchIt;
186  }
187  if ( total >= target )
188  {
189  break;
190  }
191  }
192  m_RankValue = searchIt->first;
193  m_RankIt = searchIt;
194  }
195  else
196  {
197  typename MapType::iterator searchIt = m_RankIt;
198  typename MapType::iterator eraseIt;
199 
200  while ( searchIt != m_Map.begin() )
201  {
202  ThisBin = searchIt->second;
203  unsigned int tbelow = total - ThisBin;
204  if ( tbelow < target ) // we've overshot
205  {
206  break;
207  }
208  if ( eraseFlag )
209  {
210  m_Map.erase(eraseIt);
211  eraseFlag = false;
212  }
213  if ( ThisBin <= 0 )
214  {
215  eraseIt = searchIt;
216  eraseFlag = true;
217  }
218  total = tbelow;
219 
220  --searchIt;
221  }
222  m_RankValue = searchIt->first;
223  m_RankIt = searchIt;
224  }
225 
226  m_Below = total;
227  itkAssertInDebugAndIgnoreInReleaseMacro( m_RankValue == GetValueBruteForce() );
228  return ( m_RankValue );
229  }
230 
231  void SetRank(float rank)
232  {
233  m_Rank = rank;
234  }
235 
236  void AddBoundary(){}
237 
238  void RemoveBoundary(){}
239 
241  {
242  return false;
243  }
244 
245 protected:
246  float m_Rank;
247 
248 private:
249  typedef typename std::map< TInputPixel, SizeValueType, TCompare > MapType;
250 
254  TInputPixel m_RankValue;
255  TInputPixel m_InitVal;
258  // This iterator will point at the desired rank value
259  typename MapType::iterator m_RankIt;
260 };
261 
262 
263 template< class TInputPixel >
265 {
266 public:
267  typedef std::less< TInputPixel > TCompare;
268 
270  {
272  m_Vec.resize(m_Size, 0);
275  {
277  }
278  else
279  {
281  }
282  m_Entries = m_Below = 0;
284  m_Rank = 0.5;
285  }
286 
288 
289  bool IsValid()
290  {
291  return m_Entries > 0;
292  }
293 
294  TInputPixel GetValueBruteForce()
295  {
296  SizeValueType count = 0;
297  SizeValueType target = (SizeValueType)( m_Rank * ( m_Entries - 1 ) ) + 1;
298  for( SizeValueType i=0; i<m_Size; i++ )
299  {
300  count += m_Vec[i];
301  if( count >= target )
302  {
304  }
305  }
307  }
308 
309  TInputPixel GetValue(const TInputPixel &)
310  {
311  return GetValueBruteForce();
312  }
313 
314  void AddPixel(const TInputPixel & p)
315  {
317 
318  m_Vec[q]++;
319  if ( m_Compare(p, m_RankValue) || p == m_RankValue )
320  {
321  ++m_Below;
322  }
323  ++m_Entries;
324  }
325 
326  void RemovePixel(const TInputPixel & p)
327  {
329 
330  itkAssertInDebugAndIgnoreInReleaseMacro( q >= 0 );
331  itkAssertInDebugAndIgnoreInReleaseMacro( q < (int)m_Vec.size() );
332  itkAssertInDebugAndIgnoreInReleaseMacro( m_Entries >= 1 );
333  itkAssertInDebugAndIgnoreInReleaseMacro( m_Vec[q] > 0 );
334 
335  m_Vec[q]--;
336  --m_Entries;
337 
338  if ( m_Compare(p, m_RankValue) || p == m_RankValue )
339  {
340  --m_Below;
341  }
342  }
343 
344  void SetRank(float rank)
345  {
346  m_Rank = rank;
347  }
348 
349  void AddBoundary(){}
350 
351  void RemoveBoundary(){}
352 
354  {
355  return true;
356  }
357 
358 protected:
359  float m_Rank;
360 
361 private:
362  typedef typename std::vector< SizeValueType > VecType;
363 
367  TInputPixel m_RankValue;
368  TInputPixel m_InitVal;
369  int m_Below;
371 };
372 
373 // now create MorphologicalGradientHistogram specilizations using the VectorMorphologicalGradientHistogram
374 // as base class
375 
378 template<>
379 class RankHistogram<unsigned char>:
380  public VectorRankHistogram<unsigned char>
381 {
382 };
383 
384 template<>
385 class RankHistogram<signed char>:
386  public VectorRankHistogram<signed char>
387 {
388 };
389 
390 template<>
391 class RankHistogram<bool>:
392  public VectorRankHistogram<bool>
393 {
394 };
395 
398 } // end namespace Function
399 } // end namespace itk
400 #endif
401