Index: /issm/trunk/src/c/Bamgx/objects/Triangles.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/objects/Triangles.cpp	(revision 2976)
+++ /issm/trunk/src/c/Bamgx/objects/Triangles.cpp	(revision 2977)
@@ -478,5 +478,5 @@
 	void Triangles::WriteMesh(BamgMesh* bamgmesh,BamgOpts* bamgopts){
 
-		int i,j;
+		int i,j,k,num;
 		int verbose;
 
@@ -507,4 +507,5 @@
 		bamgmesh->NumEdges=nbe;
 		xfree((void**)&bamgmesh->Edges);
+		int NumSegments=0;
 		if (nbe){
 			bamgmesh->Edges=(double*)xmalloc(3*nbe*sizeof(double));
@@ -513,11 +514,71 @@
 				bamgmesh->Edges[i*3+1]=Number(edges[i][1])+1; //back to M indexing
 				bamgmesh->Edges[i*3+2]=edges[i].ref;
-				/*if(edges[i].onGeometry){
-					printf("ref = %i [%i %i]\n",edges[i].ref,Number(edges[i][0])+1,Number(edges[i][1])+1);
-				}*/
+				if(edges[i].onGeometry){
+					NumSegments++;
+				}
 			}
 		}
 		else{
 			bamgmesh->Edges=NULL;
+		}
+
+		//Segments
+		if(verbose>3) printf("      writing Segments\n");
+		//chaining algorithm
+		int head_v[nbv];
+		int next_p[3*nbt];
+		for (i=0;i<nbv;i++) head_v[i]=-1;
+		k=0;
+		for (i=0;i<nbt;i++) {
+			//Do not take into account outside triangles (reft<0)
+			if (reft[i]>=0){
+				for (j=0;j<3;j++){
+					int v=Number(triangles[i][j]); //jth vertex of the ith triangle
+					next_p[k]=head_v[v];
+					head_v[v]=k++;
+				}
+			}
+		}
+		xfree((void**)&bamgmesh->Segments);
+		xfree((void**)&bamgmesh->SegmentsMarkers);
+		bamgmesh->NumSegments=NumSegments;
+		bamgmesh->Segments=(double*)xmalloc(3*NumSegments*sizeof(double));
+		bamgmesh->SegmentsMarkers=(double*)xmalloc(1*NumSegments*sizeof(double));
+		num=0;
+		for (i=0;i<nbe;i++){
+			if(edges[i].onGeometry){
+				//build segment
+				int i1=Number(edges[i][0]);
+				int i2=Number(edges[i][1]);
+				bool stop=false;
+				for(j=head_v[i1];j!=-1;j=next_p[j]){
+					for(k=0;k<3;k++){
+						if (Number(triangles[(int)j/3][k])==i1){
+							if (Number(triangles[(int)j/3][(int)((k+1)%3)])==i2){
+								bamgmesh->Segments[num*3+0]=Number(edges[i][0])+1; //back to M indexing
+								bamgmesh->Segments[num*3+1]=Number(edges[i][1])+1; //back to M indexing
+								bamgmesh->Segments[num*3+2]=(int)j/3+1;            //back to M indexing
+								bamgmesh->SegmentsMarkers[num]=edges[i].ref;
+								num+=1;
+								stop=true;
+								break;
+							}
+							if (Number(triangles[(int)j/3][(int)((k+2)%3)])==i2){
+								bamgmesh->Segments[num*3+0]=Number(edges[i][1])+1; //back to M indexing
+								bamgmesh->Segments[num*3+1]=Number(edges[i][0])+1; //back to M indexing
+								bamgmesh->Segments[num*3+2]=(int)j/3+1;            //back to M indexing
+								bamgmesh->SegmentsMarkers[num]=edges[i].ref;
+								num+=1;
+								stop=true;
+								break;
+							}
+						}
+					}
+					if(stop) break;
+				}
+				if (!stop){
+					throw ErrorException(__FUNCT__,exprintf("Element holding segment [%i %i] not found...",i1+1,i2+1));
+				}
+			}
 		}
 
@@ -539,6 +600,6 @@
 		//Triangles
 		if(verbose>3) printf("      writing Triangles\n");
-		Int4 k=nbInT-NbOfQuad*2;
-		Int4 num=0;
+		k=nbInT-NbOfQuad*2;
+		num=0;
 		bamgmesh->NumTriangles=k;
 		xfree((void**)&bamgmesh->Triangles);
Index: /issm/trunk/src/c/objects/BamgMesh.h
===================================================================
--- /issm/trunk/src/c/objects/BamgMesh.h	(revision 2976)
+++ /issm/trunk/src/c/objects/BamgMesh.h	(revision 2977)
@@ -16,4 +16,8 @@
 	int     NumEdges;
 	double* Edges;
+
+	int     NumSegments;
+	double* Segments;
+	double* SegmentsMarkers;
 
 	int     NumCrackedEdges;
Index: /issm/trunk/src/m/classes/public/bamg.m
===================================================================
--- /issm/trunk/src/m/classes/public/bamg.m	(revision 2976)
+++ /issm/trunk/src/m/classes/public/bamg.m	(revision 2977)
@@ -141,5 +141,5 @@
 
 %call Bamg
-[triangles vertices metric]=Bamg(bamg_mesh,bamg_geometry,bamg_options);
+[triangles vertices segments segmentsmarkers metric]=Bamg(bamg_mesh,bamg_geometry,bamg_options);
 
 % plug results onto model
@@ -147,4 +147,6 @@
 md.y=vertices(:,2);
 md.elements=triangles(:,1:3);
+md.segments=segments;
+md.segmentmarkers=segmentsmarkers;
 md.dummy=metric;
 
@@ -159,10 +161,4 @@
 md.elementonbed=ones(md.numberofelements,1);
 md.elementonsurface=ones(md.numberofelements,1);
-
-%Now, build the connectivity tables for this mesh.
-md.nodeconnectivity=NodeConnectivity(md.elements,md.numberofgrids);
-md.elementconnectivity=ElementConnectivity(md.elements,md.nodeconnectivity);
-
-%recreate segments
-md.segments=findsegments(md);
 md.gridonboundary=zeros(md.numberofgrids,1); md.gridonboundary(md.segments(:,1:2))=1;
+md.counter=1;
Index: /issm/trunk/src/mex/Bamg/Bamg.cpp
===================================================================
--- /issm/trunk/src/mex/Bamg/Bamg.cpp	(revision 2976)
+++ /issm/trunk/src/mex/Bamg/Bamg.cpp	(revision 2977)
@@ -107,4 +107,7 @@
 	bamgmesh.NumEdges=0;
 	bamgmesh.Edges=NULL;
+	bamgmesh.NumSegments=0;
+	bamgmesh.Segments=NULL;
+	bamgmesh.SegmentsMarkers=NULL;
 	FetchData(&NumCrackedEdgesMesh,mxGetField(BAMGMESH,0,"NumCrackedEdges"));
 	bamgmesh.NumCrackedEdges=NumCrackedEdgesMesh;
@@ -175,4 +178,6 @@
 	WriteData(TRIANGLESOUT,bamgmesh.Triangles,bamgmesh.NumTriangles,4);
 	WriteData(VERTICESOUT,bamgmesh.Vertices,bamgmesh.NumVertices,3);
+	WriteData(SEGMENTSOUT,bamgmesh.Segments,bamgmesh.NumSegments,3);
+	WriteData(SEGMENTSMARKERSOUT,bamgmesh.SegmentsMarkers,bamgmesh.NumSegments,1);
 	WriteData(METRICOUT,bamgopts.metric,nods,3);
 
Index: /issm/trunk/src/mex/Bamg/Bamg.h
===================================================================
--- /issm/trunk/src/mex/Bamg/Bamg.h	(revision 2976)
+++ /issm/trunk/src/mex/Bamg/Bamg.h	(revision 2977)
@@ -23,11 +23,13 @@
 
 /* serial output macros: */
-#define TRIANGLESOUT (mxArray**)&plhs[0]
-#define VERTICESOUT  (mxArray**)&plhs[1]
-#define METRICOUT    (mxArray**)&plhs[2]
+#define TRIANGLESOUT       (mxArray**)&plhs[0]
+#define VERTICESOUT        (mxArray**)&plhs[1]
+#define SEGMENTSOUT        (mxArray**)&plhs[2]
+#define SEGMENTSMARKERSOUT (mxArray**)&plhs[3]
+#define METRICOUT          (mxArray**)&plhs[4]
 
 /* serial arg counts: */
 #undef NLHS
-#define NLHS  3
+#define NLHS  5
 #undef NRHS
 #define NRHS  3