function Ct = ath(p,t)
% ATH model
% Fusco et al. Journal of Biomedical Graphics and computing, 2(2);2012

Fp = p(1); % mmol/L
Tc = p(2); % min
E  = p(3); % fraction
ve = p(4); % fraction

for k = 1:length(t)
	[Ct(k)] = quad("convolution",0,t(k));
end
Ct(t<0)=0; % before injection Ct is 0

	function rc = convolution(tc)
		if (tc<=Tc)
		rc = weinmann(tc) .* Fp;
		end
		if (tc>Tc)
		rc = weinmann(tc) .* E*Fp*exp(-(E*Fp/ve)*(t(k)-tc-Tc));
		end
	endfunction

endfunction

