001/*******************************************************************************
002 * This software is provided as a supplement to the authors' textbooks on digital
003 * image processing published by Springer-Verlag in various languages and editions.
004 * Permission to use and distribute this software is granted under the BSD 2-Clause
005 * "Simplified" License (see http://opensource.org/licenses/BSD-2-Clause).
006 * Copyright (c) 2006-2023 Wilhelm Burger, Mark J. Burge. All rights reserved.
007 * Visit https://imagingbook.com for additional details.
008 ******************************************************************************/
009
010package imagingbook.common.threshold.global;
011
012/**
013 * <p>
014 * This is an implementation of the global thresholder proposed by Ridler and Calvard [1]. See Sec. 9.1.3 (Alg. 9.3) of
015 * [2] for a detailed description.
016 * </p>
017 * <p>
018 * [1] T.W. Ridler, S. Calvard, Picture thresholding using an iterative selection method, IEEE Transactions on Systems,
019 * Man and Cybernetics, SMC-8 (August 1978) 630-632. <br> [2] W. Burger, M.J. Burge, <em>Digital Image Processing
020 * &ndash; An Algorithmic Introduction</em>, 3rd ed, Springer (2022).
021 * </p>
022 *
023 * @author WB
024 * @version 2022/04/02
025 */
026public class IsodataThresholder implements GlobalThresholder {
027        
028        private static int MAX_ITERATIONS = 100;
029
030        private double[] M0;            // table of background means
031        private double[] M1;            // table of foreground means
032        
033        public IsodataThresholder() {
034                super();
035        }
036
037        @Override
038        public float getThreshold(int[] h) {
039                makeMeanTables(h);
040                int K = h.length;
041                int q = (int) M0[K-1];  // start with total mean
042                int qq;
043                
044                int i = 0;
045                do {
046                        i++; 
047                        if (M0[q] < 0 || M1[q] < 0)  // background or foreground is empty
048                                return NoThreshold;
049                        qq = q;
050                        q = (int) ((M0[q] + M1[q]) / 2);
051                } while (q != qq && i < MAX_ITERATIONS);
052                
053                return q;
054        }
055        
056        private void makeMeanTables(int[] h) {
057                int K = h.length;
058                M0 = new double[K];
059                M1 = new double[K];
060                
061                int n0 = 0;
062                long s0 = 0;
063                for (int q = 0; q < K; q++) {
064                        n0 = n0 + h[q];
065                        s0 = s0 + q * h[q];
066                        M0[q] = (n0 > 0) ? ((double) s0)/n0 : -1;
067                }
068                
069                int n1 = 0;
070                long s1 = 0;
071                M1[K-1] = 0;
072                for (int q = h.length-2; q >= 0; q--) {
073                        n1 = n1 + h[q+1];
074                        s1 = s1 + (q+1) * h[q+1];
075                        M1[q] = (n1 > 0) ? ((double) s1)/n1 : -1;
076                }
077        }
078        
079}