im2bw_com(IM)

function BW = im2bw_com(IM)
% This function segments input image
% by detecting center of mass of input image
% input: IM  input image (color or grayscale)
% output: BW  segmented image
% example
% IM = imread('house.jpg');
% BW = im2bw_com(IM);
% imshow(BW) or imagesc(BW)

% size of input image
dim = size(IM);
% change to grayscale
if length(dim)==3
    IM = uint8((double(IM(:,:,1))*2+double(IM(:,:,2))*4+double(IM(:,:,3)))/7);
end
I = IM(:);

% histgram of input image
y = (histc(I,[0:255]))';

vec0 = [0:255];

for m = 1:length(y)-1
    % center of mass of histgram
    com = m-1;
    % distance from com
    vec = abs(com-vec0);
    % divide histgram to two classes 
    classa = y(1:m);
    veca = vec(1:m);
    classb = y(m+1:end);
    vecb = vec(m+1:end);
    % calculate moments of two classes
    vala = sum(classa.*veca);
    valb = sum(classb.*vecb);
    % calculate difference of two moments
    dif(m) = abs(vala-valb);
end

% find most balanced com
[minv,ind] = min(dif);
% threshold
thresh = ind-1
% segmentation based on threshold
BW = (IM>=thresh);