function md=meshconvert(md,varargin)
%CONVERTMESH - convert mesh to bamg mesh
%
%   Usage:
%      md=meshconvert(md);
%      md=meshconvert(md,index,x,y);

if nargin~=1 & nargin~=4,
	help meshconvert
	error('meshconvert error message: bad usage');
end

if nargin==1,
	x=md.x;
	y=md.y;
	index=md.elements;
else
	x=varargin{1};
	y=varargin{2};
	index=varargin{3};
end

%call Bamg
[bamgmesh_out bamggeom_out]=BamgConvertMesh(index,x,y);

% plug results onto model
md.bamg=struct();
md.bamg.mesh=bamgmesh(bamgmesh_out);
md.bamg.geometry=bamggeom(bamggeom_out);
md.x=bamgmesh_out.Vertices(:,1);
md.y=bamgmesh_out.Vertices(:,2);
md.elements=bamgmesh_out.Triangles(:,1:3);
md.edges=bamgmesh_out.IssmEdges;
md.segments=bamgmesh_out.IssmSegments(:,1:3);
md.segmentmarkers=bamgmesh_out.IssmSegments(:,4);

%Fill in rest of fields:
md.dim=2;
md.numberofelements=size(md.elements,1);
md.numberofnodes=length(md.x);
md.z=zeros(md.numberofnodes,1);
md.nodeonbed=ones(md.numberofnodes,1);
md.nodeonwater=zeros(md.numberofnodes,1);
md.nodeonsurface=ones(md.numberofnodes,1);
md.elementonbed=ones(md.numberofelements,1);
md.elementonsurface=ones(md.numberofelements,1);
md.nodeonboundary=zeros(md.numberofnodes,1); md.nodeonboundary(md.segments(:,1:2))=1;