forked from InsightSoftwareConsortium/ITK
-
Notifications
You must be signed in to change notification settings - Fork 0
/
ImageEntropy1.cxx
281 lines (230 loc) · 9.15 KB
/
ImageEntropy1.cxx
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
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
/*=========================================================================
*
* Copyright Insight Software Consortium
*
* Licensed under the Apache License, Version 2.0 (the "License");
* you may not use this file except in compliance with the License.
* You may obtain a copy of the License at
*
* http://www.apache.org/licenses/LICENSE-2.0.txt
*
* Unless required by applicable law or agreed to in writing, software
* distributed under the License is distributed on an "AS IS" BASIS,
* WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
* See the License for the specific language governing permissions and
* limitations under the License.
*
*=========================================================================*/
// Software Guide : BeginLatex
//
//
// This example shows how to compute the entropy of an image. More formally
// this should be said : The reduction in uncertainty gained when we measure
// the intensity of \emph{one} randomly selected pixel in this image, given
// that we already know the statistical distribution of the image intensity
// values.
//
// In practice it is almost never possible to know the real statistical
// distribution of intensities and we are forced to estimate it from the
// evaluation of the histogram from one or several images of similar nature.
// We can use the counts in histogram bins in order to compute frequencies and
// then consider those frequencies to be estimations of the probablility of a
// new value to belong to the intensity range of that bin.
//
// \index{Entropy!Images}
// \index{Image!Entropy}
// \index{Image!Amount of information}
// \index{Amount of information!Image}
//
// Software Guide : EndLatex
// Software Guide : BeginLatex
//
// Since the first stage in estimating the entropy of an image is to compute
// its histogram, we must start by including the headers of the classes that
// will perform such a computation. In this case, we are going to use a scalar
// image as input, therefore we need the
// \subdoxygen{Statistics}{ScalarImageToHistogramGenerator} class, as well as
// the image class.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
#include "itkScalarImageToHistogramGenerator.h"
#include "itkImage.h"
// Software Guide : EndCodeSnippet
#include "itkImageFileReader.h"
int main( int argc, char * argv [] )
{
if( argc < 3 )
{
std::cerr << "Missing command line arguments" << std::endl;
std::cerr << "Usage : ImageEntropy1 inputImageFileName ";
std::cerr << "numberOfHistogramBins" << std::endl;
return EXIT_FAILURE;
}
// Software Guide : BeginLatex
//
// The pixel type and dimension of the image are explicitly declared and then
// used for instantiating the image type.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
typedef unsigned char PixelType;
const unsigned int Dimension = 3;
typedef itk::Image< PixelType, Dimension > ImageType;
// Software Guide : EndCodeSnippet
typedef itk::ImageFileReader< ImageType > ReaderType;
ReaderType::Pointer reader = ReaderType::New();
reader->SetFileName( argv[1] );
try
{
reader->Update();
}
catch( itk::ExceptionObject & excp )
{
std::cerr << "Problem encoutered while reading image file : " << argv[1] << std::endl;
std::cerr << excp << std::endl;
return EXIT_FAILURE;
}
// Software Guide : BeginLatex
//
// The image type is used as template parameter for instantiating the histogram
// generator.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
typedef itk::Statistics::ScalarImageToHistogramGenerator<
ImageType > HistogramGeneratorType;
HistogramGeneratorType::Pointer histogramGenerator =
HistogramGeneratorType::New();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The parameters of the desired histogram are defined, including the
// number of bins and the marginal scale. For convenience in this example, we
// read the number of bins from the command line arguments. In this way we can
// easily experiment with different values for the number of bins and see how
// that choice affects the computation of the entropy.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
const unsigned int numberOfHistogramBins = atoi( argv[2] );
histogramGenerator->SetNumberOfBins( numberOfHistogramBins );
histogramGenerator->SetMarginalScale( 10.0 );
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We can then connect as input the output image from a reader and trigger the
// histogram computation by invoking the \code{Compute()} method in the
// generator.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
histogramGenerator->SetInput( reader->GetOutput() );
histogramGenerator->Compute();
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The resulting histogram can be recovered from the generator by using the
// \code{GetOutput()} method. A histogram class can be declared using the
// \code{HistogramType} trait from the generator.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
typedef HistogramGeneratorType::HistogramType HistogramType;
const HistogramType * histogram = histogramGenerator->GetOutput();
// Software Guide : EndCodeSnippet
const unsigned int histogramSize = histogram->Size();
std::cout << "Histogram size " << histogramSize << std::endl;
for (unsigned int bin=0; bin < histogramSize; ++bin)
{
std::cout << "bin = " << bin << " frequency = ";
std::cout << histogram->GetFrequency( bin, 0 ) << std::endl;
}
// Software Guide : BeginLatex
//
// We proceed now to compute the \emph{estimation} of entropy given the
// histogram. The first conceptual jump to be done here is to assume that
// the histogram, which is the simple count of frequency of occurrence for the
// gray scale values of the image pixels, can be normalized in order to estimate
// the probability density function \textbf{PDF} of the actual statistical
// distribution of pixel values.
//
// First we declare an iterator that will visit all the bins in the histogram.
// Then we obtain the total number of counts using the
// \code{GetTotalFrequency()} method, and we initialize the entropy variable
// to zero.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
HistogramType::ConstIterator itr = histogram->Begin();
HistogramType::ConstIterator end = histogram->End();
double Sum = histogram->GetTotalFrequency();
double Entropy = 0.0;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// We start now visiting every bin and estimating the probability of a pixel to
// have a value in the range of that bin. The base 2 logarithm of that
// probability is computed, and then weighted by the probability in order to
// compute the expected amount of information for any given pixel. Note that a
// minimum value is imposed for the probability in order to avoid computing
// logarithms of zeros.
//
// Note that the $\log{(2)}$ factor is used to convert the natural logarithm
// in to a logarithm of base 2, and makes it possible to report the entropy
// in its natural unit: the bit.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
while( itr != end )
{
const double probability = itr.GetFrequency() / Sum;
if( probability > 0.99 / Sum )
{
Entropy += - probability * std::log( probability ) / std::log( 2.0 );
}
++itr;
}
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// The result of this sum is considered to be our estimation of the image
// entropy. Note that the Entropy value will change depending on the number of
// histogram bins that we use for computing the histogram. This is particularly
// important when dealing with images whose pixel values have dynamic ranges so
// large that our number of bins will always underestimate the variability of
// the data.
//
// Software Guide : EndLatex
// Software Guide : BeginCodeSnippet
std::cout << "Image entropy = " << Entropy << " bits " << std::endl;
// Software Guide : EndCodeSnippet
// Software Guide : BeginLatex
//
// As an illustration, the application of this program to the image
//
// \begin{itemize}
// \item \code{Examples/Data/BrainProtonDensitySlice.png}
// \end{itemize}
//
// results in the following values of entropy for different values of number of
// histogram bins.
//
// \begin{center}
// \begin{tabular}{|l|r|r|r|r|r|}
// \hline
// Number of Histogram Bins & 16 & 32 & 64 & 128 & 255 \cr
// \hline
// Estimated Entropy (bits) & 3.02 & 3.98 & 4.92 & 5.89 & 6.88 \cr
// \hline
// \end{tabular}
// \end{center}
//
//
// This table highlights the importance of carefully considering the
// characteristics of the histograms used for estimating Information Theory
// measures such as the entropy.
//
// Software Guide : EndLatex
return EXIT_SUCCESS;
}