function [PID,gic] = c2d(PID,Ts,options)
%C2D  Conversion of continuous time PID to discrete time.

%   Author(s): Rong Chen
%   Copyright 2009-2017 The MathWorks, Inc.

gic = [];
if isa(PID,'ltipack.piddata2')
    % 2-DOF
    objFcn = @pidstd2;
else
    % 1-DOF
    objFcn = @pidstd;
end
switch lower(options.Method)
    case 'forwardeuler'
        % if method is ForwardEuler set Ts and Formulas directly.
        PID.Ts = Ts;
        PID.IFormula = 'F';
        PID.DFormula = 'F';
    case 'backwardeuler'
        % if method is BackwardEuler set Ts and Formulas directly.
        PID.Ts = Ts;
        PID.IFormula = 'B';
        PID.DFormula = 'B';
    case 'tustin'
        % if method is Tustin set Formulas to Trapezoidal.
        PID.Ts = Ts;
        PID.IFormula = 'T';
        PID.DFormula = 'T';
        % adjust Ti Td when prewarping is present
        w = options.PrewarpFrequency;
        if w ~= 0
            TsPre = tan(Ts*w/2)/(w/2);
            Ratio = TsPre/Ts;            
            PID.Ti = PID.Ti/Ratio;
            PID.Td = PID.Td/Ratio;
        end
    case 'zoh'
        % zoh(G(s),Ts) = (1-1/z)*Ztransform(G(s)/s) for Ts
        % s/(T*s+1) --> (1-1/z)*Z{1/(T*s+1)} = 1/T*(z-1)/(z-exp(-Ts/T))
        % default discretization method is Forward Euler
        % for Kp, Ti and N do not need to change
        if PID.Td>0
            if isinf(PID.N)
                ctrlMsgUtils.error('Control:transformation:NotSupportedImproperIZF')
            else
                EXP = exp(-Ts/PID.Td*PID.N);                
                PID.Td = Ts/(1-EXP)*PID.N;
            end
        end
        PID.Ts = Ts;
        PID.IFormula = 'F';
        PID.DFormula = 'F';
    case 'foh'
        % foh(G(s),Ts) = (z-1)^2/Ts/z*Ztransform(G(s)/s^2) for Ts
        % 1/s --> Ts*(z+1)/2/(z-1)
        % s/(T*s+1) --> (1-exp(-Ts/T))*(z-1)/(z-exp(-Ts/T))/Ts
        % default discretization method is Trapezoidal
        % for Kp, Ti and Td do not need to change
        if PID.Td>0
            if isinf(PID.N)
                ctrlMsgUtils.error('Control:transformation:NotSupportedImproperIZF')
            else
                EXP = exp(-Ts/PID.Td*PID.N);
                PID.N = PID.Td/(Ts/2*(1+EXP)/(1-EXP));
            end
        end
        PID.Ts = Ts;
        PID.IFormula = 'T';
        PID.DFormula = 'T';
    case 'impulse'
        if PID.Td>0 && isinf(PID.N)
            ctrlMsgUtils.error('Control:transformation:NotSupportedImproperIZF')
        end
        SSc = ss(PID);
        SSd = c2d(SSc,Ts,options);
        try
            PID = objFcn(SSd);
        catch
            ctrlMsgUtils.error('Control:transformation:pidConversionFailed','c2d');
        end
    case 'matched'
        ZPKc = zpk(PID);
        ZPKd = c2d(ZPKc,Ts,options);
        try
            PID = objFcn(ZPKd);
        catch
            ctrlMsgUtils.error('Control:transformation:pidConversionFailed','c2d');
        end
    case 'least-squares'
        if PID.Td>0 && isinf(PID.N)
            ctrlMsgUtils.error('Control:transformation:NotSupportedImproperL')
        end
        ZPKc = zpk(PID);
        ZPKd = c2d(ZPKc,Ts,options);
        try
            PID = objFcn(ZPKd);
        catch
            ctrlMsgUtils.error('Control:transformation:pidConversionFailed','c2d');
        end        
end
