function I = rgb2gray(X) %#codegen

% Copyright 2013-2017 The MathWorks, Inc.

narginchk(1,1);

% Determine if input includes a 2-D or 3-D array
is2D = coder.const(numel(size(X)) == 2);
is3D = coder.const(numel(size(X)) == 3);
coder.internal.assert(is2D || is3D, 'MATLAB:images:rgb2gray:invalidInputSize');

if is2D
    coder.internal.errorIf((size(X,2) ~= 3 || size(X,1) < 1),'MATLAB:images:rgb2gray:invalidSizeForColormap');
    coder.internal.errorIf(~isa(X,'double'),'MATLAB:images:rgb2gray:notAValidColormap');
else % is3D
    coder.internal.errorIf((size(X,3) ~= 3),'MATLAB:images:rgb2gray:invalidInputSizeRGB');
end

%allow only uint8, uint16, single and double.
validateattributes(X, {'uint8', 'uint16', 'double', 'single'}, {'real'}, mfilename, 'RGB');

origSize = size(X);

% Calculate transformation matrix
T    = inv([1.0 0.956 0.621; 1.0 -0.272 -0.647; 1.0 -1.106 1.703]);
coef = T(1,:);

%coef = [0.2989 0.5870 0.1140];
coder.const(coef);

coder.extrinsic('images.internal.coder.useSharedLibrary');

if is3D
    %RGB
    sizeOutput = [origSize(1), origSize(2)];

    singleThread = images.internal.coder.useSingleThread();

    useSharedLibrary = coder.const(coder.internal.isTargetMATLABHost()) && ...
        coder.const(images.internal.coder.useSharedLibrary()) && ...
        coder.const(~singleThread);

    I = coder.nullcopy(ones(sizeOutput, 'like', X));

    if(useSharedLibrary)

        fcnName = ['rgb2gray_tbb_', images.internal.coder.getCtype(I)];
        I = images.internal.coder.buildable.Rgb2grayBuildable.rgb2graycore(...
            fcnName, ...
            X, ...
            size(X,1)* size(X,2), ...
            I);
    else
        if coder.isColumnMajor
            % Do transformation
            if isa(X, 'double') || isa(X, 'single')
                sizeI = numel(X)/3;
                coder.gpu.internal.kernelImpl(false);
                parfor i = 1:sizeI
                    I(i)    = [X(i) X(i+sizeI) X(i+2*sizeI)]*coef'; %#ok<PFBNS>
                    if I(i) > 1
                        I(i) = 1;
                    elseif I(i) < 0
                        I(i) = 0;
                    end
                end
            else
                %uint8 or uint16
                sizeI = numel(X)/3;
                if isa(X,'uint8')
                    coder.gpu.internal.kernelImpl(false);
                    parfor i = 1:sizeI
                        I(i)   = uint8(double([X(i) X(i+sizeI) X(i+2*sizeI)])*coef'); %#ok<PFBNS>
                    end
                elseif isa(X,'uint16')
                    coder.gpu.internal.kernelImpl(false);
                    parfor i = 1:sizeI
                        I(i)   = uint16(double([X(i) X(i+sizeI) X(i+2*sizeI)])*coef'); %#ok<PFBNS>
                    end
                end
            end

        else  % coder.isRowMajor
            if isa(X, 'double') || isa(X, 'single')
                coder.gpu.internal.kernelImpl(false);
                parfor m = 1:size(X,1)
                    for n = 1:size(X,2)
                        I(m,n)    = [X(m,n,1) X(m,n,2) X(m,n,3)]*coef';
                        if I(m,n) > 1
                            I(m,n) = 1;
                        elseif I(m,n) < 0
                            I(m,n) = 0;
                        end
                    end
                end
            else
                %uint8 or uint16
                if isa(X,'uint8')
                    coder.gpu.internal.kernelImpl(false);
                    parfor m = 1:size(X,1)
                        for n = 1:size(X,2)
                            I(m,n)   = uint8(double([X(m,n,1) X(m,n,2) X(m,n,3)])*coef');
                        end
                    end
                elseif isa(X,'uint16')
                    coder.gpu.internal.kernelImpl(false);
                    parfor m = 1:size(X,1)
                        for n = 1:size(X,2)
                            I(m,n)   = uint16(double([X(m,n,1) X(m,n,2) X(m,n,3)])*coef');
                        end
                    end
                end
            end
        end
    end
else % input is 2D
    I1 = X * coef';
    I2 = min(max(I1,0),1);
    I = [I2,I2,I2];
end
