--- a/inst/graythresh.m
+++ b/inst/graythresh.m
@@ -9,76 +9,76 @@
 ## but WITHOUT ANY WARRANTY; without even the implied warranty of
 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 ## GNU General Public License for more details.
-
-## -*- texinfo -*-
-## @deftypefn {Function File} {@var{level}=} graythresh (@var{I})
-## Compute global image threshold using Otsu's method.
-##
-## The output @var{level} is a global threshold (level) that can be used to convert
-## an intensity image to a binary image with @code{im2bw}.
-## @var{level} is a normalized intensity value that lies in the range [0, 1]. 
-## 
-## The function uses Otsu's method, which chooses the threshold to
-## minimize the intraclass variance of the black and white pixels. 
-## 
-## Color images are converted grayscale before @var{level} is computed.
-## @seealso{im2bw}
-## @end deftypefn
-
-## Note:
-## This function is taken from
-## http://www.irit.fr/recherches/SAMOVA/MEMBERS/JOLY/Homepage_files/IRR05/Barre-Piquot/graythresh.m
-## I added texinfo documentation, error checking and sanitised the code.
-##    -- Søren Hauberg
-
-function level = graythresh (I)
-    ## Input checking
-    if (nargin != 1)
-      print_usage();
-    endif
-    if (!isgray(I) && !isrgb(I))
-      error("graythresh: input must be an image");
-    endif
-    
-    ## If the image is RGB convert it to grayscale
-    if (isrgb(I))
-      I = rgb2gray(I);
-    endif
-
-    ## Calculation of the normalized histogram
-    n = 256;
-    h = hist(I(:), n);        
-    h = h/(length(I(:))+1);
-    
-    ## Calculation of the cumulated histogram and the mean values
-    w = cumsum(h);
-    mu = zeros(n, 1); mu(1) = h(1);
-    for i=2:n
-        mu(i) = mu(i-1) + i*h(i);
-    end    
-         
-    ## Initialisation of the values used for the threshold calculation
-    w0 = w(1);
-    w1 = 1-w0;
-    mu0 = mu(1)/w0;
-    mu1 = (mu(end)-mu(1))/w1;
-    max = w0*w1*(mu1-mu0)*(mu1-mu0);
-    level = 1;
-    
-    ## For each step of the histogram, calculation of the threshold and storing of the maximum
-    for i = 2:n
-        w0 = w(i);
-        w1 = 1-w0;
-        mu0 = mu(i)/w0;
-        mu1 = (mu(end)-mu(i))/w1;
-        s = w0*w1*(mu1-mu0)*(mu1-mu0);
-        if (s > max)
-            max = s;
-            level = i;
-        endif
-    endfor
-    
-    ## Normalisation of the threshold        
-    level /= n;
-endfunction
-
+
+## -*- texinfo -*-
+## @deftypefn {Function File} {@var{level}=} graythresh (@var{I})
+## Compute global image threshold using Otsu's method.
+##
+## The output @var{level} is a global threshold (level) that can be used to convert
+## an intensity image to a binary image with @code{im2bw}.
+## @var{level} is a normalized intensity value that lies in the range [0, 1]. 
+## 
+## The function uses Otsu's method, which chooses the threshold to
+## minimize the intraclass variance of the black and white pixels. 
+## 
+## Color images are converted grayscale before @var{level} is computed.
+## @seealso{im2bw}
+## @end deftypefn
+
+## Note:
+## This function is taken from
+## http://www.irit.fr/recherches/SAMOVA/MEMBERS/JOLY/Homepage_files/IRR05/Barre-Piquot/graythresh.m
+## I added texinfo documentation, error checking and sanitised the code.
+##    -- Søren Hauberg
+
+function level = graythresh (I)
+    ## Input checking
+    if (nargin != 1)
+      print_usage();
+    endif
+    if (!isgray(I) && !isrgb(I))
+      error("graythresh: input must be an image");
+    endif
+    
+    ## If the image is RGB convert it to grayscale
+    if (isrgb(I))
+      I = rgb2gray(I);
+    endif
+
+    ## Calculation of the normalized histogram
+    n = 256;
+    h = hist(I(:), 1:n);        
+    h = h/(length(I(:))+1);
+    
+    ## Calculation of the cumulated histogram and the mean values
+    w = cumsum(h);
+    mu = zeros(n, 1); mu(1) = h(1);
+    for i=2:n
+        mu(i) = mu(i-1) + i*h(i);
+    end    
+         
+    ## Initialisation of the values used for the threshold calculation
+    level = find (h > 0, 1);
+    w0 = w(level);
+    w1 = 1-w0;
+    mu0 = mu(level)/w0;
+    mu1 = (mu(end)-mu(level))/w1;
+    max = w0*w1*(mu1-mu0)*(mu1-mu0);
+    
+    ## For each step of the histogram, calculation of the threshold and storing of the maximum
+    for i = find (h > 0)
+        w0 = w(i);
+        w1 = 1-w0;
+        mu0 = mu(i)/w0;
+        mu1 = (mu(end)-mu(i))/w1;
+        s = w0*w1*(mu1-mu0)*(mu1-mu0);
+        if (s > max)
+            max = s;
+            level = i;
+        endif
+    endfor
+    
+    ## Normalisation of the threshold        
+    level /= n;
+endfunction
+