View Javadoc

1   /**
2    * Copyright (c) 2011, The University of Southampton and the individual contributors.
3    * All rights reserved.
4    *
5    * Redistribution and use in source and binary forms, with or without modification,
6    * are permitted provided that the following conditions are met:
7    *
8    *   * 	Redistributions of source code must retain the above copyright notice,
9    * 	this list of conditions and the following disclaimer.
10   *
11   *   *	Redistributions in binary form must reproduce the above copyright notice,
12   * 	this list of conditions and the following disclaimer in the documentation
13   * 	and/or other materials provided with the distribution.
14   *
15   *   *	Neither the name of the University of Southampton nor the names of its
16   * 	contributors may be used to endorse or promote products derived from this
17   * 	software without specific prior written permission.
18   *
19   * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
20   * ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
21   * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
22   * DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
23   * ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
24   * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
25   * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON
26   * ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
27   * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
28   * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29   */
30  package org.openimaj.image.processing.threshold;
31  
32  import org.openimaj.citation.annotation.Reference;
33  import org.openimaj.citation.annotation.ReferenceType;
34  import org.openimaj.image.FImage;
35  import org.openimaj.image.processor.ImageProcessor;
36  import org.openimaj.util.array.ArrayUtils;
37  import org.openimaj.util.pair.FloatFloatPair;
38  
39  /**
40   * Otsu's adaptive thresholding algorithm.
41   * 
42   * @see <a
43   *      href="http://en.wikipedia.org/wiki/Otsu's_method">http://en.wikipedia.org/wiki/Otsu&apos;s_method</a>
44   * 
45   * @author Jonathon Hare (jsh2@ecs.soton.ac.uk)
46   */
47  @Reference(
48  		type = ReferenceType.Article,
49  		author = { "Nobuyuki Otsu" },
50  		title = "A Threshold Selection Method from Gray-Level Histograms",
51  		year = "1979",
52  		journal = "Systems, Man and Cybernetics, IEEE Transactions on",
53  		pages = { "62", "66" },
54  		number = "1",
55  		volume = "9",
56  		customData = {
57  				"keywords",
58  				"Displays;Gaussian distribution;Histograms;Least squares approximation;Marine vehicles;Q measurement;Radar tracking;Sea measurements;Surveillance;Target tracking",
59  				"doi", "10.1109/TSMC.1979.4310076",
60  				"ISSN", "0018-9472"
61  		})
62  public class OtsuThreshold implements ImageProcessor<FImage> {
63  	private static final int DEFAULT_NUM_BINS = 256;
64  	int numBins = DEFAULT_NUM_BINS;
65  
66  	/**
67  	 * Default constructor
68  	 */
69  	public OtsuThreshold() {
70  
71  	}
72  
73  	/**
74  	 * Construct with the given number of histogram bins
75  	 * 
76  	 * @param numBins
77  	 *            the number of histogram bins
78  	 */
79  	public OtsuThreshold(int numBins) {
80  		this.numBins = numBins;
81  	}
82  
83  	protected static int[] makeHistogram(FImage fimg, int numBins) {
84  		final int[] histData = new int[numBins];
85  
86  		// Calculate histogram
87  		for (int r = 0; r < fimg.height; r++) {
88  			for (int c = 0; c < fimg.width; c++) {
89  				final int h = (int) (fimg.pixels[r][c] * (numBins - 1));
90  				histData[h]++;
91  			}
92  		}
93  
94  		return histData;
95  	}
96  
97  	protected static int[] makeHistogram(float[] data, int numBins, float min, float max) {
98  		final int[] histData = new int[numBins];
99  
100 		// Calculate histogram
101 		for (int c = 0; c < data.length; c++) {
102 			final float d = (data[c] - min) / (max - min);
103 			final int h = (int) (d * (numBins - 1));
104 			histData[h]++;
105 		}
106 
107 		return histData;
108 	}
109 
110 	/**
111 	 * Estimate the threshold for the given image.
112 	 * 
113 	 * @param img
114 	 *            the image
115 	 * @param numBins
116 	 *            the number of histogram bins
117 	 * @return the estimated threshold
118 	 */
119 	public static float calculateThreshold(FImage img, int numBins) {
120 		final int[] histData = makeHistogram(img, numBins);
121 
122 		// Total number of pixels
123 		final int total = img.getWidth() * img.getHeight();
124 
125 		return computeThresholdFromHistogram(histData, total);
126 	}
127 
128 	/**
129 	 * Estimate the threshold for the given data.
130 	 * <p>
131 	 * Internally, the data will be min-max normalised before the histogram is
132 	 * built, and the specified number of bins will cover the entire
133 	 * <code>max-min</code> range. The returned threshold will have
134 	 * <code>min</code> added to it to return it to the original range.
135 	 * 
136 	 * @param data
137 	 *            the data
138 	 * @param numBins
139 	 *            the number of histogram bins
140 	 * @return the estimated threshold
141 	 */
142 	public static float calculateThreshold(float[] data, int numBins) {
143 		final float min = ArrayUtils.minValue(data);
144 		final float max = ArrayUtils.maxValue(data);
145 		final int[] histData = makeHistogram(data, numBins, min, max);
146 
147 		return computeThresholdFromHistogram(histData, data.length) + min;
148 	}
149 
150 	/**
151 	 * Estimate the threshold and inter-class variance for the given data.
152 	 * <p>
153 	 * Internally, the data will be min-max normalised before the histogram is
154 	 * built, and the specified number of bins will cover the entire
155 	 * <code>max-min</code> range. The returned threshold will have
156 	 * <code>min</code> added to it to return it to the original range.
157 	 * 
158 	 * @param data
159 	 *            the data
160 	 * @param numBins
161 	 *            the number of histogram bins
162 	 * @return the estimated threshold and variance
163 	 */
164 	public static FloatFloatPair calculateThresholdAndVariance(float[] data, int numBins) {
165 		final float min = ArrayUtils.minValue(data);
166 		final float max = ArrayUtils.maxValue(data);
167 		final int[] histData = makeHistogram(data, numBins, min, max);
168 
169 		final FloatFloatPair result = computeThresholdAndVarianceFromHistogram(histData, data.length);
170 		result.first += min;
171 		return result;
172 	}
173 
174 	/**
175 	 * Estimate the threshold for the given histogram.
176 	 * 
177 	 * @param histData
178 	 *            the histogram
179 	 * @param total
180 	 *            the total number of items in the histogram
181 	 * @return the estimated threshold
182 	 */
183 	public static float computeThresholdFromHistogram(int[] histData, int total) {
184 		return computeThresholdAndVarianceFromHistogram(histData, total).first;
185 	}
186 
187 	/**
188 	 * Estimate the threshold and inter-class variance for the given histogram.
189 	 * 
190 	 * @param histData
191 	 *            the histogram
192 	 * @param total
193 	 *            the total number of items in the histogram
194 	 * @return the estimated threshold and variance
195 	 */
196 	public static FloatFloatPair computeThresholdAndVarianceFromHistogram(int[] histData, int total) {
197 		final int numBins = histData.length;
198 		float sum = 0;
199 		for (int t = 0; t < numBins; t++)
200 			sum += t * histData[t];
201 
202 		float sumB = 0;
203 		int wB = 0;
204 		int wF = 0;
205 
206 		float varMax = 0;
207 		float threshold = 0;
208 
209 		for (int t = 0; t < numBins; t++) {
210 			wB += histData[t]; // Weight Background
211 			if (wB == 0)
212 				continue;
213 
214 			wF = total - wB; // Weight Foreground
215 			if (wF == 0)
216 				break;
217 
218 			sumB += (t * histData[t]);
219 
220 			final float mB = sumB / wB; // Mean Background
221 			final float mF = (sum - sumB) / wF; // Mean Foreground
222 
223 			// Calculate Between Class Variance
224 			final float varBetween = (float) wB * (float) wF * (mB - mF) * (mB - mF);
225 
226 			// Check if new maximum found
227 			if (varBetween > varMax) {
228 				varMax = varBetween;
229 				threshold = t;
230 			}
231 		}
232 
233 		return new FloatFloatPair(threshold / (numBins - 1), varMax / total / total);
234 	}
235 
236 	@Override
237 	public void processImage(FImage image) {
238 		final float threshold = calculateThreshold(image, numBins);
239 
240 		image.threshold(threshold);
241 	}
242 }