%% JL20151014: TO BE REVISED FOR EFFICIENCY 
function PlygnMesh = PlygnMesh_Enrich2(PlygnMesh,BndryDescMat)
%% Enriching a polygonal mesh converted from PolyMesher 
% James Liu, Zhuoran Wang, ColoState; 2012/07--2017/08 

%%
NumEgs = PlygnMesh.NumEgs;
PlygnMesh.BndryEdge = zeros(NumEgs,1);
for ig=1:NumEgs 
  % Utilizing the info obtained the above 
  % PlygnMesh.edge2elem(ig,2)>0 if it is an interior edge 
  if PlygnMesh.edge2elem(ig,2)>0 
    continue;
  end
  % Assuming now RectMesh.edge2elem(2,jg)==0 
  x1 = PlygnMesh.node(PlygnMesh.edge(ig,1),1);
  x2 = PlygnMesh.node(PlygnMesh.edge(ig,2),1);
  y1 = PlygnMesh.node(PlygnMesh.edge(ig,1),2);
  y2 = PlygnMesh.node(PlygnMesh.edge(ig,2),2);
  % Screening 
  for k=1:size(BndryDescMat,1)  % Number of segments on the boundary
    X1 = BndryDescMat(k,1);  Y1 = BndryDescMat(k,2);  % Starting point
    X2 = BndryDescMat(k,3);  Y2 = BndryDescMat(k,4);  % Ending point
    if abs((x1-X1)*(Y2-Y1)-(X2-X1)*(y1-Y1))<1e-9 && abs((x2-X1)*(Y2-Y1)-(X2-X1)*(y2-Y1))<1e-9 
      PlygnMesh.BndryEdge(ig) = k;
      break;
    end
  end
end

return;