function zDiscrete = discretiseTf(tfContinuous,rSampleTime,varargin)
% discretiseTf
%
% zDiscrete = discretiseTf(tfContinuous,rSampleTime,stMethod)
%
% discretiseTf returns a structure with the tranfer function whereby the dc
% gain of its discrete form is removed and then the delay is removed as
% well.  The sequence of the algorithm is:
%
% tfContinuous --(discretitation)--> tfInitialDiscreteTf
% tfInitialDiscreteTf minus rSeparateGain
% z(tfInitialDiscreteTf - rSeparateGain)---(continuos)--> tfDiscrete
%
% xDiscrete, a structure containing all the required information, is then
% created.
%
% The string stMethod selects the discretization method among the following:
%        'zoh'       Zero-order hold on the inputs
%        'foh'       Linear interpolation of inputs
%        'impulse'   Impulse-invariant discretization
%        'tustin'    Bilinear (Tustin) approximation.
%        'matched'   Matched pole-zero method (for SISO systems only).
% If no method is entered the tustin method is used.  Note that the third
% input must be named varargin rather than stMethod to achieve this.
%
% Example of use:
% G = tf([1 0.5],[1 0.3]);
% H = discretiseTf(G,0.02);
% H:
%     tfContinuousTf:   [1x1 tf]
%     arNumerator:      [0.001994017946162 0.199401794616152]
%     arDenominator:    [1 0.299999999999997]
%     rGain:            1.001994017946162

%   Project: IMPACT
%   Title of Code: complexToMagPhase
%   Relevant additional documentation:
%   Developer: Adam Stock (based on original by Carlos Gonzalez)
%   External Dependencies: NA
%
%   Disclaimer
%   CAUTION: It is strongly recommended that changes to this code
%   should be avoided without first consulting with authors.


%% Check the inputs
if(~isobject(tfContinuous))
    error('Input argument type not valid');
end
if(~isnumeric(rSampleTime))
    error('Input argument type not valid');
end
% If no method is specified use Tustin, otherwise use the method specified
if (nargin==2)
    stMethod = 'tustin';
elseif(nargin==3)
    stMethod=varargin{1};
else
    error('Wrong number of input parameters');
end

%% Discretise the system

%%
% Define the discrete operator z
z = tf('z',rSampleTime);

%%
% Ensure the transfer function has been simplified and then convert using
% the appropriate method
tfContinuousMinimal = minreal(tfContinuous);
tfInitialDiscreteTf=c2d(tfContinuousMinimal,rSampleTime,stMethod);

%%
% Find the separate gain
[arNumerator,~]=tfdata(tfInitialDiscreteTf,'v');
rSeparateGain=arNumerator(1);

%%
% Find the transfer function in continuous time that is paired with the
% separate gain
tfDiscreteDelayed = tfInitialDiscreteTf-rSeparateGain;
tfDiscrete = tfDiscreteDelayed*z;
tfContinuousOutput = d2c(tfDiscrete,stMethod);

%%
% Save the requisite values as a structure

[arContinuousNumerator,arContinuousDenominator] = tfdata(tfContinuousOutput,'v');
zDiscrete.tfContinuous= tfContinuousOutput;
zDiscrete.arNumerator = arContinuousNumerator;
zDiscrete.arDenominator = arContinuousDenominator;
zDiscrete.rGain = rSeparateGain;

















