function imfeat = compute_imfeat(img, feat_id)
% extact color histogram, edge histogram and other things from image(s)
%
% 2007-09-16, xlx@ee.columbia.edu

if ischar(img)
    rgb = imread(img); % input can be an RGB matrix or a filename
else
    rgb = img;
end
hsv = rgb2hsv(rgb);
imgry = rgb2gray(rgb);
imsiz = size(hsv);

switch feat_id
    case 'color_hist',
        % color histogram
        nbins = [8 4 4]; % 8 bins for hue, 4 bins each for saturation and lightness
        num_shift = [16 4 1 128]; % multiplier for the linear numbering of bins
        
        qimg = zeros(imsiz(1:2));
        for i = 1 : 3
            tmp = floor( (hsv(:,:,i)-eps)*nbins(i)) ;
            qimg = qimg + tmp*num_shift(i);
        end
        clr_hist = hist(qimg(:), 0:num_shift(end)-1 );
        clr_hist = clr_hist/sum(clr_hist);
        feat = clr_hist;
    case 'edge_hist'
        % edge histogram
        % first compute image edges and gradient
        bw = edge(imgry, 'canny');
        [gx, gy] = gradient(double(imgry));

        ebins = pi/8*[-4:4]; % edge direction centers, -pi/2 to pi/2 in pi/8 increments
        edge_hist = hist(atan(gy(bw==1)./ (gx(bw==1) + gx(bw==1)==0)), ebins);
        edge_hist = edge_hist/sum(edge_hist);
        
        feat = edge_hist;
    case 'glcm'
        % co-occurrence texture
        os_inc = 1:2:5 ;
        offsets = [ zeros(1,3), os_inc, os_inc, -(os_inc); ...
            os_inc, zeros(1,3), os_inc, os_inc ]';

        imrsz = imresize(imgry, [NaN, 256]);
        % compute gray-level co-occurrence matrix
        glcm = graycomatrix(imrsz, 'Offset', offsets, 'symmetric', true);
        % stats of the GLCM
        stats = graycoprops(glcm);
        feat = stats;
        
    case 'clms' 
        % color moments
        clms = zeros(3, 3);
        for i = 1 : 3
            chnl = double(rgb(:,:,i)); % values in the current channel
            %h = hist(chnl(:), 256);
            clms(:, i) = [mean(chnl(:)), std(chnl(:)), skewness(chnl(:))];
        end
        feat = clms(:); % stretch the 3x3 colormoments to a 9x1 vector
    otherwise
        disp('do not know what feature to comptue');
        feat = [];
end

% done, return results
imfeat = struct('name', feat_id, 'feat', feat);