THESE/codes/graphe/graphcuts.m

   1 function [Ncut] = graphcuts(I,pad,MAXVAL)\r
   2 % function [Ncut] = graphcuts(I)\r
   3 % Input: I image\r
   4 % pad: spatial connectivity; eg. 3\r
   5 % MAXVAL: maximum image value\r
   6 % Output: Ncut: Binary map 0 or 1 corresponding to image segmentation\r
   7 I = double(I); [H,W] = size(I);\r
   8 % Find weights between nodes I1 and I2, w = exp(a*abs(I1-I2));\r
   9 % Set a to have a weight of 0.01 for diff = MAXVAL\r
  10 a = log(0.01)/MAXVAL; x = [0:MAXVAL/100:MAXVAL]'; y = exp(a*x);\r
  11 figure;plot(x,y);xlabel('intensity diff');ylabel('weights'); title('weights')\r
  12 ws = 2*pad + 1;\r
  13 if(ws <= 3)\r
  14     ws = 3;\r
  15 end\r
  16 %Build the weight matrix\r
  17 disp('Building Weight Matrix'); close all; tic\r
  18 WM = zeros(H*W,H*W); countWM = 0;\r
  19 for kk = 1:W\r
  20 for jj = 1:H\r
  21 mask = logical(zeros(H,W));\r
  22 cs = kk-pad; ce = kk+pad; rs = jj-pad; re = jj+pad;\r
  23 if(cs<1) \r
  24     cs = 1; \r
  25 end;\r
  26 if(ce>W)\r
  27     ce = W; \r
  28 end;\r
  29 if(rs<1) \r
  30     rs = 1;\r
  31 end;\r
  32 if(re>H) \r
  33     re = H; \r
  34 end;\r
  35 mask(rs:re,cs:ce) = 1;\r
  36 idx = find(mask==1);\r
  37 p = abs(I(idx) - I(jj,kk)); p = exp(a*p);\r
  38 countWM = countWM + 1; WM(countWM,idx) = p(:)';\r
  39 end\r
  40 end\r
  41 ttime = toc; disp(sprintf('Time for generating weight matrix = %f',ttime)); clear countWM\r
  42 % Weight between a node and iteself is 0\r
  43 for jj = 1:H*W \r
  44     WM(jj,jj) = 0; \r
  45 end; \r
  46 WM = sparse(WM);\r
  47 % Shi and Malik Algorithm: second smallest eigen vector\r
  48 disp('Finding Eigen Vector');\r
  49 d = sum(WM,2); D = diag(d); tic\r
  50 B = (D-WM); B = (B+B')/2; OPTS.disp = 0;\r
  51 [v,d,flag] = eigs(B,D,2,'SA',OPTS); ttime = toc;\r
  52 disp(sprintf('Time for finding eigen vector = %f',ttime)); clear OPTS\r
  53 y = v(:,2);\r
  54 Ncut = reshape(y,H,W);\r
  55 Ncut = Ncut > 0;