Index: /issm/trunk/src/c/objects/BamgGeom.cpp
===================================================================
--- /issm/trunk/src/c/objects/BamgGeom.cpp	(revision 3272)
+++ /issm/trunk/src/c/objects/BamgGeom.cpp	(revision 3273)
@@ -25,2 +25,101 @@
 
 }
+
+void BamgGeomWrite(mxArray** pbamggeom_mat, BamgGeom* bamggeom){
+
+	/*Output*/
+	mxArray*    bamggeom_mat=NULL;
+
+	/*Intermediary*/
+	int         i,j;
+	mxArray*    pfield=NULL;
+	mxArray*    pfield2=NULL;
+	int         numfields=14;
+	const char* fnames[numfields];
+	mwSize      ndim=2;
+	mwSize      dimensions[2]={1,1};
+
+	fnames[0] = "NumVertices";
+	fnames[1] = "Vertices";
+	fnames[2] = "NumEdges";
+	fnames[3] = "Edges";
+	fnames[4] = "NumTangentAtEdges";
+	fnames[5] = "TangentAtEdges";
+	fnames[6] = "NumRequiredVertices";
+	fnames[7] = "RequiredVertices";
+	fnames[8] = "NumRequiredEdges";
+	fnames[9] = "RequiredEdges";
+	fnames[10]= "NumCrackedEdges";
+	fnames[11]= "CrackedEdges";
+	fnames[12]= "NumSubDomains";
+	fnames[13]= "SubDomains";
+
+	bamggeom_mat=mxCreateStructArray(ndim,dimensions,numfields,fnames);
+
+	mxSetField(bamggeom_mat,0,"NumVertices",mxCreateDoubleScalar(bamggeom->NumVertices)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,3);
+	mxSetN(pfield,bamggeom->NumVertices);
+	mxSetPr(pfield,bamggeom->Vertices);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamggeom_mat,0,"Vertices",pfield2);
+
+	mxSetField(bamggeom_mat,0,"NumEdges",mxCreateDoubleScalar(bamggeom->NumEdges)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,3);
+	mxSetN(pfield,bamggeom->NumEdges);
+	mxSetPr(pfield,bamggeom->Edges);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamggeom_mat,0,"Edges",pfield2);
+
+	mxSetField(bamggeom_mat,0,"NumTangentAtEdges",mxCreateDoubleScalar(bamggeom->NumTangentAtEdges)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,4);
+	mxSetN(pfield,bamggeom->NumTangentAtEdges);
+	mxSetPr(pfield,bamggeom->TangentAtEdges);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamggeom_mat,0,"TangentAtEdges",pfield2);
+
+	mxSetField(bamggeom_mat,0,"NumRequiredVertices",mxCreateDoubleScalar(bamggeom->NumRequiredVertices)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,1);
+	mxSetN(pfield,bamggeom->NumRequiredVertices);
+	mxSetPr(pfield,bamggeom->RequiredVertices);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamggeom_mat,0,"RequiredVertices",pfield2);
+
+	mxSetField(bamggeom_mat,0,"NumRequiredEdges",mxCreateDoubleScalar(bamggeom->NumRequiredEdges)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,1);
+	mxSetN(pfield,bamggeom->NumRequiredEdges);
+	mxSetPr(pfield,bamggeom->RequiredEdges);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamggeom_mat,0,"RequiredEdges",pfield2);
+
+	mxSetField(bamggeom_mat,0,"NumCrackedEdges",mxCreateDoubleScalar(bamggeom->NumCrackedEdges)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,2);
+	mxSetN(pfield,bamggeom->NumCrackedEdges);
+	mxSetPr(pfield,bamggeom->CrackedEdges);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamggeom_mat,0,"CrackedEdges",pfield2);
+
+	mxSetField(bamggeom_mat,0,"NumSubDomains",mxCreateDoubleScalar(bamggeom->NumSubDomains)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,4);
+	mxSetN(pfield,bamggeom->NumSubDomains);
+	mxSetPr(pfield,bamggeom->SubDomains);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamggeom_mat,0,"SubDomains",pfield2);
+
+	/*Assign output pointer*/
+	*pbamggeom_mat=bamggeom_mat;
+
+}
Index: /issm/trunk/src/c/objects/BamgGeom.h
===================================================================
--- /issm/trunk/src/c/objects/BamgGeom.h	(revision 3272)
+++ /issm/trunk/src/c/objects/BamgGeom.h	(revision 3273)
@@ -4,4 +4,6 @@
 #ifndef BAMGGEOM_H_
 #define BAMGGEOM_H_
+
+#include <mex.h>
 
 struct BamgGeom{
@@ -37,4 +39,5 @@
 
 void BamgGeomInit(BamgGeom* bamggeom);
+void BamgGeomWrite(mxArray** bamggeom_mat,BamgGeom* bamggeom);
 
 #endif
Index: /issm/trunk/src/c/objects/BamgMesh.cpp
===================================================================
--- /issm/trunk/src/c/objects/BamgMesh.cpp	(revision 3272)
+++ /issm/trunk/src/c/objects/BamgMesh.cpp	(revision 3273)
@@ -30,2 +30,152 @@
 
 }
+
+void BamgMeshWrite(mxArray** pbamgmesh_mat, BamgMesh* bamgmesh){
+
+	/*Output*/
+	mxArray*    bamgmesh_mat=NULL;
+
+	/*Intermediary*/
+	int         i,j;
+	mxArray*    pfield=NULL;
+	mxArray*    pfield2=NULL;
+	int         numfields=23;
+	const char* fnames[numfields];
+	mwSize      ndim=2;
+	mwSize      dimensions[2]={1,1};
+
+	fnames[0] = "NumTriangles";
+	fnames[1] = "Triangles";
+	fnames[2] = "NumVertices";
+	fnames[3] = "Vertices";
+	fnames[4] = "NumEdges";
+	fnames[5] = "Edges";
+	fnames[6] = "NumSegments";
+	fnames[7] = "Segments";
+	fnames[8] = "SegmentsMarkers";
+	fnames[9] = "NumCrackedEdges";
+	fnames[10] = "CrackedEdges";
+	fnames[11] = "NumQuadrilaterals";
+	fnames[12] = "Quadrilaterals";
+	fnames[13] = "NumVerticesOnGeometricVertex";
+	fnames[14] = "VerticesOnGeometricVertex";
+	fnames[15] = "NumVerticesOnGeometricEdge";
+	fnames[16] = "VerticesOnGeometricEdge";
+	fnames[17] = "NumEdgesOnGeometricEdge";
+	fnames[18] = "EdgesOnGeometricEdge";
+	fnames[19] = "NumSubDomains";
+	fnames[20] = "SubDomains";
+	fnames[21] = "NumSubDomainsFromGeom";
+	fnames[22] = "SubDomainsFromGeom";
+
+	bamgmesh_mat=mxCreateStructArray(ndim,dimensions,numfields,fnames);
+
+	mxSetField(bamgmesh_mat,0,"NumTriangles",mxCreateDoubleScalar(bamgmesh->NumTriangles)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,4);
+	mxSetN(pfield,bamgmesh->NumTriangles);
+	mxSetPr(pfield,bamgmesh->Triangles);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"Triangles",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumVertices",mxCreateDoubleScalar(bamgmesh->NumVertices)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,3);
+	mxSetN(pfield,bamgmesh->NumVertices);
+	mxSetPr(pfield,bamgmesh->Vertices);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"Vertices",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumEdges",mxCreateDoubleScalar(bamgmesh->NumEdges)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,3);
+	mxSetN(pfield,bamgmesh->NumEdges);
+	mxSetPr(pfield,bamgmesh->Edges);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"Edges",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumSegments",mxCreateDoubleScalar(bamgmesh->NumSegments)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,3);
+	mxSetN(pfield,bamgmesh->NumSegments);
+	mxSetPr(pfield,bamgmesh->Segments);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"Segments",pfield2);
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,1);
+	mxSetN(pfield,bamgmesh->NumSegments);
+	mxSetPr(pfield,bamgmesh->SegmentsMarkers);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"SegmentsMarkers",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumCrackedEdges",mxCreateDoubleScalar(bamgmesh->NumCrackedEdges)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,2);
+	mxSetN(pfield,bamgmesh->NumCrackedEdges);
+	mxSetPr(pfield,bamgmesh->CrackedEdges);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"CrackedEdges",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumQuadrilaterals",mxCreateDoubleScalar(bamgmesh->NumQuadrilaterals)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,5);
+	mxSetN(pfield,bamgmesh->NumQuadrilaterals);
+	mxSetPr(pfield,bamgmesh->Quadrilaterals);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"Quadrilaterals",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumVerticesOnGeometricVertex",mxCreateDoubleScalar(bamgmesh->NumVerticesOnGeometricVertex)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,2);
+	mxSetN(pfield,bamgmesh->NumVerticesOnGeometricVertex);
+	mxSetPr(pfield,bamgmesh->VerticesOnGeometricVertex);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"VerticesOnGeometricVertex",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumVerticesOnGeometricEdge",mxCreateDoubleScalar(bamgmesh->NumVerticesOnGeometricEdge)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,2);
+	mxSetN(pfield,bamgmesh->NumVerticesOnGeometricEdge);
+	mxSetPr(pfield,bamgmesh->VerticesOnGeometricEdge);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"VerticesOnGeometricEdge",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumEdgesOnGeometricEdge",mxCreateDoubleScalar(bamgmesh->NumEdgesOnGeometricEdge)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,2);
+	mxSetN(pfield,bamgmesh->NumEdgesOnGeometricEdge);
+	mxSetPr(pfield,bamgmesh->EdgesOnGeometricEdge);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"EdgesOnGeometricEdge",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumSubDomains",mxCreateDoubleScalar(bamgmesh->NumSubDomains)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,4);
+	mxSetN(pfield,bamgmesh->NumSubDomains);
+	mxSetPr(pfield,bamgmesh->SubDomains);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"SubDomains",pfield2);
+
+	mxSetField(bamgmesh_mat,0,"NumSubDomainsFromGeom",mxCreateDoubleScalar(bamgmesh->NumSubDomainsFromGeom)); 
+
+	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
+	mxSetM(pfield,4);
+	mxSetN(pfield,bamgmesh->NumSubDomainsFromGeom);
+	mxSetPr(pfield,bamgmesh->SubDomainsFromGeom);
+	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
+	mxSetField(bamgmesh_mat,0,"SubDomainsFromGeom",pfield2);
+
+	/*Assign output pointer*/
+	*pbamgmesh_mat=bamgmesh_mat;
+}
Index: /issm/trunk/src/c/objects/BamgMesh.h
===================================================================
--- /issm/trunk/src/c/objects/BamgMesh.h	(revision 3272)
+++ /issm/trunk/src/c/objects/BamgMesh.h	(revision 3273)
@@ -4,4 +4,6 @@
 #ifndef BAMGMESH_H_
 #define BAMGMESH_H_
+
+#include <mex.h>
 
 struct BamgMesh{
@@ -45,4 +47,5 @@
 
 void BamgMeshInit(BamgMesh* bamgmesh);
+void BamgMeshWrite(mxArray** bamgmesh_mat,BamgMesh* bamgmesh);
 
 #endif
Index: /issm/trunk/src/m/classes/@model/model.m
===================================================================
--- /issm/trunk/src/m/classes/@model/model.m	(revision 3272)
+++ /issm/trunk/src/m/classes/@model/model.m	(revision 3273)
@@ -22,4 +22,5 @@
 	
 	%Mesh
+	md.bamg=NaN;
 	md.type='';
 	md.numberofelements=0;
Index: /issm/trunk/src/m/classes/public/bamg.m
===================================================================
--- /issm/trunk/src/m/classes/public/bamg.m	(revision 3272)
+++ /issm/trunk/src/m/classes/public/bamg.m	(revision 3273)
@@ -276,4 +276,7 @@
 
 % plug results onto model
+md.bamg=struct();
+md.bamg.mesh=bamgmesh_out;
+md.bamg.geometry=bamggeom_out;
 md.x=bamgmesh_out.Vertices(:,1);
 md.y=bamgmesh_out.Vertices(:,2);
Index: /issm/trunk/src/m/classes/public/display/displaymesh.m
===================================================================
--- /issm/trunk/src/m/classes/public/display/displaymesh.m	(revision 3272)
+++ /issm/trunk/src/m/classes/public/display/displaymesh.m	(revision 3273)
@@ -39,4 +39,5 @@
 disp(sprintf('\n      Properties:'));
 fielddisplay(md,'type','mesh type');
+fielddisplay(md,'bamg','Geometry and 2d mesh properties (if generated by Bamg)');
 fielddisplay(md,'penalties','penalties list');
 fielddisplay(md,'gridonbed','lower nodes flags list');
Index: /issm/trunk/src/mex/Bamg/Bamg.cpp
===================================================================
--- /issm/trunk/src/mex/Bamg/Bamg.cpp	(revision 3272)
+++ /issm/trunk/src/mex/Bamg/Bamg.cpp	(revision 3273)
@@ -5,4 +5,8 @@
 
 void mexFunction( int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]){
+
+	/*Outputs*/
+	mxArray* bamgmesh_mat=NULL;
+	mxArray* bamggeom_mat=NULL;
 
 	/*diverse: */
@@ -83,234 +87,10 @@
 	Bamgx(&bamgmesh_out,&bamggeom_out,&bamgmesh_in,&bamggeom_in,&bamgopts);
 
-	/*Variables*/
-	mxArray*    bamgmesh_mat=NULL;
-	mxArray*    pfield=NULL;
-	mxArray*    pfield2=NULL;
-	int         numfields=23;
-	const char* fnames[numfields];
-	mwSize      ndim=2;
-	mwSize      dimensions[2]={1,1};
-
-	fnames[0] = "NumTriangles";
-	fnames[1] = "Triangles";
-	fnames[2] = "NumVertices";
-	fnames[3] = "Vertices";
-	fnames[4] = "NumEdges";
-	fnames[5] = "Edges";
-	fnames[6] = "NumSegments";
-	fnames[7] = "Segments";
-	fnames[8] = "SegmentsMarkers";
-	fnames[9] = "NumCrackedEdges";
-	fnames[10] = "CrackedEdges";
-	fnames[11] = "NumQuadrilaterals";
-	fnames[12] = "Quadrilaterals";
-	fnames[13] = "NumVerticesOnGeometricVertex";
-	fnames[14] = "VerticesOnGeometricVertex";
-	fnames[15] = "NumVerticesOnGeometricEdge";
-	fnames[16] = "VerticesOnGeometricEdge";
-	fnames[17] = "NumEdgesOnGeometricEdge";
-	fnames[18] = "EdgesOnGeometricEdge";
-	fnames[19] = "NumSubDomains";
-	fnames[20] = "SubDomains";
-	fnames[21] = "NumSubDomainsFromGeom";
-	fnames[22] = "SubDomainsFromGeom";
-
-	bamgmesh_mat=mxCreateStructArray(ndim,dimensions,numfields,fnames);
-
-	mxSetField(bamgmesh_mat,0,"NumTriangles",mxCreateDoubleScalar(bamgmesh_out.NumTriangles)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,4);
-	mxSetN(pfield,bamgmesh_out.NumTriangles);
-	mxSetPr(pfield,bamgmesh_out.Triangles);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"Triangles",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumVertices",mxCreateDoubleScalar(bamgmesh_out.NumVertices)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,3);
-	mxSetN(pfield,bamgmesh_out.NumVertices);
-	mxSetPr(pfield,bamgmesh_out.Vertices);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"Vertices",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumEdges",mxCreateDoubleScalar(bamgmesh_out.NumEdges)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,3);
-	mxSetN(pfield,bamgmesh_out.NumEdges);
-	mxSetPr(pfield,bamgmesh_out.Edges);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"Edges",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumSegments",mxCreateDoubleScalar(bamgmesh_out.NumSegments)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,3);
-	mxSetN(pfield,bamgmesh_out.NumSegments);
-	mxSetPr(pfield,bamgmesh_out.Segments);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"Segments",pfield2);
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,1);
-	mxSetN(pfield,bamgmesh_out.NumSegments);
-	mxSetPr(pfield,bamgmesh_out.SegmentsMarkers);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"SegmentsMarkers",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumCrackedEdges",mxCreateDoubleScalar(bamgmesh_out.NumCrackedEdges)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,2);
-	mxSetN(pfield,bamgmesh_out.NumCrackedEdges);
-	mxSetPr(pfield,bamgmesh_out.CrackedEdges);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"CrackedEdges",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumQuadrilaterals",mxCreateDoubleScalar(bamgmesh_out.NumQuadrilaterals)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,5);
-	mxSetN(pfield,bamgmesh_out.NumQuadrilaterals);
-	mxSetPr(pfield,bamgmesh_out.Quadrilaterals);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"Quadrilaterals",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumVerticesOnGeometricVertex",mxCreateDoubleScalar(bamgmesh_out.NumVerticesOnGeometricVertex)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,2);
-	mxSetN(pfield,bamgmesh_out.NumVerticesOnGeometricVertex);
-	mxSetPr(pfield,bamgmesh_out.VerticesOnGeometricVertex);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"VerticesOnGeometricVertex",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumVerticesOnGeometricEdge",mxCreateDoubleScalar(bamgmesh_out.NumVerticesOnGeometricEdge)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,2);
-	mxSetN(pfield,bamgmesh_out.NumVerticesOnGeometricEdge);
-	mxSetPr(pfield,bamgmesh_out.VerticesOnGeometricEdge);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"VerticesOnGeometricEdge",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumEdgesOnGeometricEdge",mxCreateDoubleScalar(bamgmesh_out.NumEdgesOnGeometricEdge)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,2);
-	mxSetN(pfield,bamgmesh_out.NumEdgesOnGeometricEdge);
-	mxSetPr(pfield,bamgmesh_out.EdgesOnGeometricEdge);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"EdgesOnGeometricEdge",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumSubDomains",mxCreateDoubleScalar(bamgmesh_out.NumSubDomains)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,4);
-	mxSetN(pfield,bamgmesh_out.NumSubDomains);
-	mxSetPr(pfield,bamgmesh_out.SubDomains);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"SubDomains",pfield2);
-
-	mxSetField(bamgmesh_mat,0,"NumSubDomainsFromGeom",mxCreateDoubleScalar(bamgmesh_out.NumSubDomainsFromGeom)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,4);
-	mxSetN(pfield,bamgmesh_out.NumSubDomainsFromGeom);
-	mxSetPr(pfield,bamgmesh_out.SubDomainsFromGeom);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamgmesh_mat,0,"SubDomainsFromGeom",pfield2);
+	/*Generate output Matlab Structures*/
+	BamgMeshWrite(&bamgmesh_mat,&bamgmesh_out);
+	BamgGeomWrite(&bamggeom_mat,&bamggeom_out);
 
 	/*assign output datasets: */
 	*BAMGMESHOUT=bamgmesh_mat;
-
-	/*Variables*/
-	mxArray*    bamggeom_mat=NULL;
-	int         numfields2=14;
-	const char* fnames2[numfields];
-
-	fnames[0] = "NumVertices";
-	fnames[1] = "Vertices";
-	fnames[2] = "NumEdges";
-	fnames[3] = "Edges";
-	fnames[4] = "NumTangentAtEdges";
-	fnames[5] = "TangentAtEdges";
-	fnames[6] = "NumRequiredVertices";
-	fnames[7] = "RequiredVertices";
-	fnames[8] = "NumRequiredEdges";
-	fnames[9] = "RequiredEdges";
-	fnames[10]= "NumCrackedEdges";
-	fnames[11]= "CrackedEdges";
-	fnames[12]= "NumEdges";
-	fnames[13]= "SubDomains";
-
-	bamggeom_mat=mxCreateStructArray(ndim,dimensions,numfields,fnames);
-
-	mxSetField(bamggeom_mat,0,"NumVertices",mxCreateDoubleScalar(bamggeom_out.NumVertices)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,3);
-	mxSetN(pfield,bamggeom_out.NumVertices);
-	mxSetPr(pfield,bamggeom_out.Vertices);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamggeom_mat,0,"Vertices",pfield2);
-
-	mxSetField(bamggeom_mat,0,"NumEdges",mxCreateDoubleScalar(bamggeom_out.NumEdges)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,3);
-	mxSetN(pfield,bamggeom_out.NumEdges);
-	mxSetPr(pfield,bamggeom_out.Edges);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamggeom_mat,0,"Edges",pfield2);
-
-	mxSetField(bamggeom_mat,0,"NumTangentAtEdges",mxCreateDoubleScalar(bamggeom_out.NumTangentAtEdges)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,4);
-	mxSetN(pfield,bamggeom_out.NumTangentAtEdges);
-	mxSetPr(pfield,bamggeom_out.TangentAtEdges);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamggeom_mat,0,"TangentAtEdges",pfield2);
-
-	mxSetField(bamggeom_mat,0,"NumRequiredVertices",mxCreateDoubleScalar(bamggeom_out.NumRequiredVertices)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,1);
-	mxSetN(pfield,bamggeom_out.NumRequiredVertices);
-	mxSetPr(pfield,bamggeom_out.RequiredVertices);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamggeom_mat,0,"RequiredVertices",pfield2);
-
-	mxSetField(bamggeom_mat,0,"NumRequiredEdges",mxCreateDoubleScalar(bamggeom_out.NumRequiredEdges)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,1);
-	mxSetN(pfield,bamggeom_out.NumRequiredEdges);
-	mxSetPr(pfield,bamggeom_out.RequiredEdges);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamggeom_mat,0,"RequiredEdges",pfield2);
-
-	mxSetField(bamggeom_mat,0,"NumCrackedEdges",mxCreateDoubleScalar(bamggeom_out.NumCrackedEdges)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,2);
-	mxSetN(pfield,bamggeom_out.NumCrackedEdges);
-	mxSetPr(pfield,bamggeom_out.CrackedEdges);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamggeom_mat,0,"CrackedEdges",pfield2);
-
-	mxSetField(bamggeom_mat,0,"NumSubDomains",mxCreateDoubleScalar(bamggeom_out.NumSubDomains)); 
-
-	pfield=mxCreateDoubleMatrix(0,0,mxREAL);
-	mxSetM(pfield,4);
-	mxSetN(pfield,bamggeom_out.NumSubDomains);
-	mxSetPr(pfield,bamggeom_out.SubDomains);
-	mexCallMATLAB(1,&pfield2,1,&pfield,"'");//transpose
-	mxSetField(bamggeom_mat,0,"SubDomains",pfield2);
-
-	/*assign output datasets: */
 	*BAMGGEOMOUT=bamggeom_mat;
 
@@ -322,8 +102,7 @@
 }
 
-void BamgUsage(void)
-{
+void BamgUsage(void){
 	_printf_("\n");
-	_printf_("   usage: [elements,vertices,segments,segmentmarkers,metric]=%s(bamgmesh,bamggeom,bamgoptions);\n",__FUNCT__);
+	_printf_("   usage: [bamgmesh,bamggeom]=%s(bamgmesh,bamggeom,bamgoptions);\n",__FUNCT__);
 	_printf_("\n");
 }
