% The function dncm_filt.m creates a filter characterizing
% the direct neural coupling between Qlu and F.
%
% Function arguments:
%	th - current parameter values
%	S - 3x1 vector containing the desired granularity, duration
%	    of autonomic impulse responses, and averaging ratio of 
%           length of averaging window to Sgran
%	  - [Sgran; Slength; Sratio] 
%
% Function outputs:
%	h - vector representing dncm discrete-time filter
%

function h = dncm_filt(th,S)

% Assigning variables.
Sgran = S(1);
Slength = S(2);
Sratio = S(3);

% Creating unit-area beta-sympathetic impulse response.
sir = zeros(Slength,1);
st1 = ((2/Sgran)-((Sratio*Sgran/2)/Sgran))+1;
en1 = length(Sgran:Sgran:3.0)+st1-1;
st2 = en1+1;
en2 = length(Sgran:Sgran:25.0)+st2-1;
sir(st1:en1) = ((1.0/42.0)*(Sgran:Sgran:3.0))';
sir(st2:en2) = ((-1.0/350.0)*(Sgran:Sgran:25.0)+(1.0/14.0))';
sir = Sgran*sir;

% Creating unit-area parasympathetic impulse response.
vir = zeros(Slength,1);
numz = (0.5-(Sratio*Sgran/2))/Sgran;
st1 = 1+numz;
en1 = length(Sgran:Sgran:1)+st1-1;
st2 = en1+1;
en2 = length(Sgran:Sgran:1)+st2-1;
vir(st1:en1) = (1)*(Sgran:Sgran:1)';
vir(st2:en2) = (1)*(-(Sgran:Sgran:1)+1)';
vir = Sgran*vir;

% Scaling by supine static gain values.
%bres = th(38)*((7.2/9)/60000)*sir;
%pres = th(39)*((15/9)/60000)*vir;

% Scaling by standing static gain values.
bres = th(53)*th(38)*0.5*((18/9)/60000)*sir;
pres = th(52)*th(39)*0.3*((17.6/9)/60000)*vir;

% Forming impulse response between Qlu and F.  That is,
% parasympathetic withdrawal followed by beta-sympathetic
% withdrawal.
h = pres-bres;