function [Q T p] = findQT(s,freq,p0)

% [Q T p] = findQT(s,freq)
% [Q T p] = findQT(s,freq,p0)
%
% Input Variables
%   s - ECG lead II segment approximately from the Q to the T
%   freq - sampling frequency

%
% Output Variables
%   Q - the start of the Q wave
%   T - the end of the T wave
%   p - parameters of the model
%
% Description
%   Uses a model based approach to compute the Q start and the T end.

% Created
%   4/20/2006, Richard J. Povinelli, Marquette University
%
% Modified
%   Mohamed A Mneimneh, Marquette University
% 
% This software is released under the terms of the GNU General
% Public License (http://www.gnu.org/copyleft/gpl.html).

s = s(:)'; %make sure s is a row vector
maxError = s(:)' * s(:); %sum of square error, straight line approximation

i = 1;
while i <= 30 %run optimization 10 times

 [x0 A b Aeq beq lb ub] = findQTlin1(s,freq); %get constraints

  if nargin == 3
    x0 = convertP(p0);
  end %if


  options = optimset('MaxIter',10000,'MaxFunEvals',100000,'TolFun',1e-3);
  [x{i} fval(i) exitflag] = fmincon(@(x) ecgError(x,s,freq),x0,A,b,Aeq,beq,lb,ub,@(x) findQTnonlin(x,freq), options);
  while (exitflag == 0)
    oldfval = fval(i); %save previous solution, because sometimes it gets worse
    oldx = x{i};
    [x{i} fval(i) exitflag] = fmincon(@(x) ecgError(x,s,freq),x0,A,b,Aeq,beq,lb,ub,@(x) findQTnonlin(x,freq));
    if fval(i) >= oldfval %bad run of fmincon
      fval(i) = oldfval;
      x{i} = oldx;
      exitflag = 100;
    end %if
  end %while

  if fval(i) < maxError %not a straight line solution
    i = i+1;
  end %if, otherwise try again

end %while

i1=1;
[temp index] = min(fval);
x1{i1} = x{index};
fval1(i1) = temp;
for i1 = 2:20
    options = optimset('MaxIter',10000,'MaxFunEvals',100000,'TolFun',1e-3);
    [x1{i1} fval1(i1) exitflag] = fmincon(@(x) ecgError(x,s,freq),x1{i1-1},A,b,Aeq,beq,lb,ub,@(x) findQTnonlin(x,freq), options);
%     i1 = i1 + 1;
end
[temp index] = min(fval1);

x = x1{index};
p = convertX(x);

