Index: /issm/trunk/src/c/Bamgx/Mesh2.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/Mesh2.cpp	(revision 2782)
+++ /issm/trunk/src/c/Bamgx/Mesh2.cpp	(revision 2783)
@@ -4131,4 +4131,5 @@
 			Triangles * OldCurrentTh =CurrentTh;
 			CurrentTh=this;
+			int verbosity=1;
 			//  Int4 NbTold = nbt;
 			// generation of the integer coor
@@ -4440,6 +4441,7 @@
 						 // 
 
-						 delete edge4;
-						 delete [] st;
+						 //delete edge4; //TESTESTESTSTSTSTSTST
+						 //delete [] st; //TESTESTESTSTSTSTSTST
+						 printf("ok1\n"); 
 						 for (i=0;i<nbv;i++)
 						  quadtree->Add(vertices[i]);
@@ -4463,4 +4465,5 @@
 							  }
 			  }
+			printf("ok2\n"); 
 			CurrentTh=OldCurrentTh;
 		}
Index: /issm/trunk/src/c/Bamgx/MeshRead.cpp
===================================================================
--- /issm/trunk/src/c/Bamgx/MeshRead.cpp	(revision 2782)
+++ /issm/trunk/src/c/Bamgx/MeshRead.cpp	(revision 2783)
@@ -41,46 +41,229 @@
 void Triangles::ReadFromMatlabMesh(BamgMesh* bamgmesh, BamgArgs* bamgargs){
 
-	int verbosity;
-	Real8 hmin = HUGE_VAL;// the infinie value 
+	int verbose;
+	Real8 Hmin = HUGE_VAL;// the infinie value 
 	Int4 i1,i2,i3,iref;
-	Int4 i;
+	Int4 i,j;
 	Int4 hvertices =0;
 	Int4 ifgeom=0;
 	Metric M1(1);
-	Int4 dim=2;
+	int Version=1,dim=2;
+
+	verbose=bamgargs->verbose;
 
 	nbv=bamgmesh->numberofnodes;
-	nbvx=bamgargs->maxnbv;
+	nbvx=nbv;
 	nbt=bamgmesh->numberofelements;
-	nbtx=2*nbvx-2; // for filling The Holes and quadrilaterals 
 	nbiv=nbvx;
 
-	vertices=(Vertex*)xmalloc(nbv*sizeof(Vertex));
-	ordre=(Vertex**)xmalloc(nbv*sizeof(Vertex*));
-
-	/*Create Vertices*/
-	for (i=0;i<bamgmesh->numberofnodes;i++){
-		vertices[i].r.x=bamgmesh->x[i];
-		vertices[i].r.y=bamgmesh->y[i];
-		vertices[i].ReferenceNumber=1;
-		vertices[i].DirOfSearch =NoDirOfSearch;
-		vertices[i].m=M1;
-		vertices[i].color =0;
-	}
-
-	/*Create triangles*/
-	triangles=(Triangle*)xmalloc(nbt*sizeof(Triangle));
-	for (i=0;i<bamgmesh->numberofelements;i++){
-		Triangle & t = triangles[i];
-		i1=(Int4)bamgmesh->index[i*3+0]-1; //for C indexing
-		i2=(Int4)bamgmesh->index[i*3+1]-1; //for C indexing
-		i3=(Int4)bamgmesh->index[i*3+2]-1; //for C indexing
-		t=Triangle(this,i1,i2,i3);
-		t.color=1;
-	}
-
-	/*Recreate geometry: */
-	ConsGeometry(-1);	//MaximalAngleOfCorner is actually in BamgGeom...
-	Gh.AfterRead();
+	//Vertices
+	if(bamgmesh->x){
+		if(verbose>3) printf("Reading Vertices\n");
+
+		vertices=(Vertex*)xmalloc(nbv*sizeof(Vertex));
+		ordre=(Vertex**)xmalloc(nbv*sizeof(Vertex*));
+
+		for (i=0;i<nbv;i++){
+			vertices[i].r.x=bamgmesh->x[i];
+			vertices[i].r.y=bamgmesh->y[i];
+			vertices[i].ReferenceNumber=1;
+			vertices[i].DirOfSearch =NoDirOfSearch;
+			vertices[i].m=M1;
+			vertices[i].color =0;
+		}
+		nbtx=2*nbvx-2; // for filling The Holes and quadrilaterals 
+	}
+	else{
+		if(verbose>3) throw ErrorException(__FUNCT__,exprintf("no Vertices found in the initial mesh"));
+	}
+
+	//Triangles
+	if(bamgmesh->index){
+		if(verbose>3) printf("Reading index\n");
+		triangles =new Triangle[nbt];
+		for (i=0;i<bamgmesh->numberofelements;i++){
+			Triangle & t = triangles[i];
+			i1=(Int4)bamgmesh->index[i*3+0]-1; //for C indexing
+			i2=(Int4)bamgmesh->index[i*3+1]-1; //for C indexing
+			i3=(Int4)bamgmesh->index[i*3+2]-1; //for C indexing
+			t=Triangle(this,i1,i2,i3);
+			t.color=1; //reference = 1 for all triangles since it has not been specified
+		}
+	}
+	else{
+		if(verbose>3) throw ErrorException(__FUNCT__,exprintf("no Vertices found in the initial mesh"));
+	}
+
+	//Quadrilaterals
+	if(bamgmesh->Quadrilaterals){
+		if(verbose>3) printf("Reading Quadrilaterals\n");
+		Int4 i1,i2,i3,i4,iref;
+		triangles =new Triangle[nbt];
+		for (i=0;i<bamgmesh->NumQuadrilaterals;i++){
+			Triangle & t1 = triangles[2*i];
+			Triangle & t2 = triangles[2*i+1];
+			i1=(Int4)bamgmesh->Quadrilaterals[i*4+0]-1; //for C indexing
+			i2=(Int4)bamgmesh->Quadrilaterals[i*4+1]-1; //for C indexing
+			i3=(Int4)bamgmesh->Quadrilaterals[i*4+2]-1; //for C indexing
+			i4=(Int4)bamgmesh->Quadrilaterals[i*4+3]-1; //for C indexing
+			t1=Triangle(this,i1,i2,i3);
+			t2=Triangle(this,i3,i4,i1);
+			t1.color=1; //reference = 1 for all triangles since it has not been specified
+			t2.color=1; //reference = 1 for all triangles since it has not been specified
+			t1.SetHidden(OppositeEdge[1]); // two times  because the adj was not created 
+			t2.SetHidden(OppositeEdge[1]); 
+		}
+	}
+	else{
+		if(verbose>3) printf("No Quadrilaterals found\n");
+	}
+
+	//hVertices
+	if(bamgmesh->hVertices){
+		if(verbose>3) printf("Reading hVertices\n");
+		hvertices=1;
+		for (i=0;i< nbv;i++){
+			vertices[i].m=Metric((Real4)bamgmesh->hVertices[i]);
+		}
+	}
+	else{
+		if(verbose>3) printf("No hVertices found\n");
+	}
+
+	//VertexOnGeometricEdge
+	if(bamgmesh->VertexOnGeometricEdge){
+		if(verbose>3) printf("Reading VertexOnGeometricEdge\n");
+		NbVerticesOnGeomEdge=bamgmesh->NumVertexOnGeometricEdge;
+		VerticesOnGeomEdge= new  VertexOnGeom[NbVerticesOnGeomEdge] ;
+		for (i=0;i<NbVerticesOnGeomEdge;i++){
+			Int4  i1,i2;
+			Real8 s;
+			//VertexOnGeom & v =VerticesOnGeomVertex[i0];
+			i1=(Int4)bamgmesh->VertexOnGeometricEdge[i*3+0]-1; //for C indexing
+			i2=(Int4)bamgmesh->VertexOnGeometricEdge[i*3+1]-1; //for C indexing
+			s =(Int4)bamgmesh->VertexOnGeometricEdge[i*3+2];
+			VerticesOnGeomEdge[i]=VertexOnGeom(vertices[i1],Gh.edges[i2],s);
+		}
+	}
+	else{
+		if(verbose>3) printf("No VertexOnGeometricEdge found\n");
+	}
+
+	//Edges
+	if (bamgmesh->Edges){
+		int i1,i2;
+		R2 zero2(0,0);
+		Real4 *len =0;
+
+		if(verbose>3) printf("Reading Edges\n");
+		nbe=bamgmesh->NumEdges;
+		edges = new Edge[nbe];
+
+		if (!hvertices) {
+			len = new Real4[nbv];
+			for(i=0;i<nbv;i++)
+			 len[i]=0;
+		}
+
+		for (i=0;i<nbe;i++){
+			i1=(int)bamgmesh->Edges[i*3+0]-1; //-1 for C indexing
+			i2=(int)bamgmesh->Edges[i*3+1]-1; //-1 for C indexing
+			edges[i].v[0]= vertices +i1;
+			edges[i].v[1]= vertices +i2;
+			edges[i].adj[0]=0;
+			edges[i].adj[1]=0;
+			R2 x12 = vertices[i2].r-vertices[i1].r;
+			Real8 l12=sqrt( (x12,x12));        
+
+			if (!hvertices) {
+				vertices[i1].color++;
+				vertices[i2].color++;
+				len[i1]+= l12;
+				len[i2] += l12;
+			}
+			Hmin = Min(Hmin,l12);
+		}
+
+		// definition  the default of the given mesh size 
+		if (!hvertices){
+			for (i=0;i<nbv;i++) 
+			 if (vertices[i].color > 0) 
+			  vertices[i].m=  Metric(len[i] /(Real4) vertices[i].color);
+			 else 
+			  vertices[i].m=  Metric(Hmin);
+			delete [] len;
+		}
+
+		// construction of edges[].adj 
+		for (i=0;i<nbv;i++){ 
+		 vertices[i].color = (vertices[i].color ==2) ? -1 : -2;
+		}
+		for (i=0;i<nbe;i++){
+			for (j=0;j<2;j++) { 
+				Vertex *v=edges[i].v[j];
+				Int4 i0=v->color,j0;
+				if(i0==-1){
+					v->color=i*2+j;
+				}
+				else if (i0>=0) {// i and i0 edge are adjacent by the vertex v
+					j0 =  i0%2;
+					i0 =  i0/2;
+					assert( v ==  edges[i0 ].v[j0]);
+					edges[i ].adj[ j ] =edges +i0;
+					edges[i0].adj[ j0] =edges +i ;
+					assert(edges[i0].v[j0] == v);
+					v->color = -3;
+				}
+			}
+		}
+	}
+	else{
+		if(verbose>3) printf("No Edges found\n");
+	}
+
+	//EdgeOnGeometricEdge
+	if(bamgmesh->EdgeOnGeometricEdge){
+		if(verbose>3) printf("Reading EdgeOnGeometricEdge\n");
+		int i1,i2,i,j;
+		i2=bamgmesh->NumEdgeOnGeometricEdge;
+		for (i1=0;i1<i2;i1++) {
+			i=(int)bamgmesh->EdgeOnGeometricEdge[i*2+0];
+			j=(int)bamgmesh->EdgeOnGeometricEdge[i*2+1];
+			if(!(i>0 && j >0 && i <= nbe && j <= Gh.nbe)) {
+				throw ErrorException(__FUNCT__,exprintf("EdgeOnGeometricEdge error: We must have : (i>0 && j >0 && i <= nbe && j <= Gh.nbe)"));
+			}
+			edges[i-1].on = Gh.edges + j-1;
+		}
+	}
+	else{
+		if(verbose>3) printf("No EdgeOnGeometricEdge found\n");
+	}
+
+	//SubDomain
+	if(bamgmesh->SubDomain){
+		Int4 i3,head,sens;
+		if(verbose>3) printf("Reading SubDomain\n");
+		NbSubDomains=bamgmesh->NumSubDomain;
+		subdomains = new SubDomain [ NbSubDomains ];
+		for (i=0;i<NbSubDomains;i++) {
+			i3  =(int)bamgmesh->SubDomain[i*3+0];
+			head=(int)bamgmesh->SubDomain[i*3+1]-1;//C indexing
+			sens=(int)bamgmesh->SubDomain[i*3+2];
+			if (i3!=23) throw ErrorException(__FUNCT__,exprintf("Bad Subdomain definition: first number should be 3"));
+			if (head<0 || head>=nbt) throw ErrorException(__FUNCT__,exprintf("Bad Subdomain definition: head should in [1 %i] (triangle number)",nbt));
+			//subdomains[i].head=i3;   //TESTTESTETETETETETETETEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEEE
+		}
+	}
+	else{
+		if(verbose>3) printf("No SubDomain found\n");
+	}
+
+	/*Recreate geometry if needed*/
+	if(1!=0) {
+		printf("Warning: mesh present but no geometry found. Reconstructing...\n");
+		/*Recreate geometry: */
+		ConsGeometry(-1);	//MaximalAngleOfCorner is actually in BamgGeom...
+		Gh.AfterRead();
+	}
 }
 
@@ -419,6 +602,4 @@
       }	  
 }
-
-
 
 
@@ -965,16 +1146,10 @@
 
   ReadFromMatlabMesh(bamgmesh,bamgargs);
-  printf("ok2\n");
-
   SetIntCoor();
-
   FillHoleInMesh();
-  throw ErrorException(__FUNCT__,exprintf("test test")); 
-
-
+  throw ErrorException(__FUNCT__,exprintf("TEST"));
 }
 
-void Geometry::ReadGeometry(const char * filename)
-{
+void Geometry::ReadGeometry(const char * filename) {
 	long int verbosity=0;
 
@@ -1238,6 +1413,4 @@
 	}
 }
-
-
 
 void Geometry::ReadGeometry(MeshIstream & f_in,const char * filename)