[IsectErr,Xclip,Yclip]=intrsect ...

function [IsectErr,Xclip,Yclip]=intrsect ...
  (Xclip,Yclip,Line,iStartPoint, ...
  iPtConsider,IdxNewPt,Epsilon)
%
% [IsectErr,Xclip,Yclip]=intrsect ...
%   (Xclip,Yclip,Line,iStartPoint, ...
%   iPtConsider,IdxNewPt,Epsilon)
% 
% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% This function determines the point of
% intersection between a specified line
% and reference polygon using linear
% interpolation.
%
% Xclip,Yclip  - coordinates of reference 
%                polygon on input
% Line         - vector defining the clipping
%                line in form:
%                Line(1)+Line(2)*X+Line(3)*Y=0
% iStartPoint  - starting node in reference
%                polygon segment
% iPtConsider  - ending node in reference
%                polygon segment
% IdxNewPt     - index into Xclip/Yclip to
%                store generated point
% Epsilon      - for determining zero values
%
% IsectErr     - error condition if set to one
% Xclip,Yclip  - updated coordinates of clipped
%                polygon on output
%
% User m functions called:  none
%----------------------------------------------

IsectErr=0;
C1=Line(2)*Xclip(iStartPoint)+ ...
   Line(3)*Yclip(iStartPoint);
C2=Line(2)*Xclip(iPtConsider)+ ...
   Line(3)*Yclip(iPtConsider)-C1;
C1=C1+Line(1);
if abs(C2) < Epsilon
  if abs(C1) < Epsilon
    %...Line is tangent to plane, return point
    %...halfway between
    xtmp=(Xclip(iStartPoint)+ ...
         Xclip(iPtConsider))/2;
    ytmp=(Yclip(iStartPoint)+ ...
         Yclip(iPtConsider))/2;
  else
    %...Line is parallel or in the plane
    IsectErr=1; return;
  end
else
  %...Ordinary situation
  C12=-C1/C2;
  xtmp=Xclip(iStartPoint)+ ...
       C12*(Xclip(iPtConsider)- ...
       Xclip(iStartPoint));
  ytmp=Yclip(iStartPoint)+ ...
       C12*(Yclip(iPtConsider)- ...
       Yclip(iStartPoint));
end
Xclip(IdxNewPt)=xtmp; Yclip(IdxNewPt)=ytmp;