function [phi2,DD] = sweep3(phi2,n,h,N,phi,iix,iiy) 
%% second order correction 
%% update based on quadrant splitting with 1rst order direction 
%% gradient 

DD = zeros(size(phi)) ;   %% save Determinant 

 for j = iiy
  for i = iix
%%%% shutof source data (exact)
 if (phi(i,j) > 0.1)  

%% find the second order stencil using first order solution  
%% upwind quadrant does not work for grid lines 
%% 
if ( phi(i-1, j)  < phi(i,j) ) 
ix = -1 ; 
a1  = phi(i,j)-phi(i-1,j) ; 
elseif  ( phi(i+1, j) < phi(i,j) ) 
ix  = 1 ; 
a1  = phi(i,j)-phi(i+1,j) ; 
else 
ix = 0 ; 
a1 =  0 ; 
end 

if ( phi(i, j-1) < phi(i,j) ) 
b1 = phi(i,j) - phi(i,j-1)  ; 
jy = -1 ; 
elseif  ( phi(i,j) > phi(i,j+1)  ) 
b1 = phi(i,j) - phi(i,j+1)  ; 
jy  = 1 ; 
else 
b1 = 0  ; 
jy =  0 ; 
end 

alpha = max(eps,b1/max(eps,a1)) ;  %% this is the tangent with machine precision cutoff 

%%%%%%%%  stencil region descisison paramter  3 zones center and close to grid lines 
cp =  pi/12 ;   %% angle of cutoff 
lcut = tan(cp) ; ucut = tan(pi/2-cp) ; 

%%% central stencil 
fx = phi2(i+ix,j) ; 
fy = phi2(i,j+jy) ; 
fxy = phi2(i+ix,j+jy) ; 

if (alpha  < lcut) 
%% right  
d2 = phi2(i+ix,j+1)+phi2(i+ix,j-1)-2*fx ;     
B = 1/2*(phi2(i+ix,j+1)-phi2(i+ix,j-1)) + jy*alpha*d2  ; 
A = -fx + 1/2*(alpha)^2*d2 ; 
c1 = 1 ; 
c2 = 2*A ; 
c3 = B^2 + A^2 - h^2*n(i,j)^2 ; 
D = c2^2-4*c1*c3 ; 

elseif (alpha > ucut) 

d2 = phi2(i+1,j+jy)+phi2(i-1,j+jy)-2*fy ;     
A = 1/2*(phi2(i+1,j+jy)-phi2(i-1,j+jy)) + ix/alpha*d2  ; 
B = -fy + 1/2*(1/alpha)^2*d2 ; 
c1 = 1 ; 
c2 = 2*B ; 
c3 = B^2 + A^2 - h^2*n(i,j)^2 ; 
D = c2^2-4*c1*c3 ; 

else 

A = (alpha/2*(fy-fxy) - (1-alpha/2)*fx)  ; 
B = (1/(2*alpha)*(fx-fxy) - (1-1/(2*alpha))*fy)  ; 
c1 = (1-alpha/2)^2 + (1-1/(2*alpha))^2 ; 
c2 = 2*A*(1-alpha/2) + 2*B*(1-1/(2*alpha)) ; 
c3 = B^2 + A^2 - h^2*n(i,j)^2 ; 
D = c2^2-4*c1*c3 ; 

end 

DD(i,j) = D ; %% store determinant
if ( D >=0)  %% must be positive 
phi2(i,j) =(-c2+sqrt(D))/(2*c1) ;  
end 

%%%%%%%%%%%%
%end %% if 
end %% for 
end %% for  
end %% if 
