Index: /issm/trunk-jpl/src/c/modules/BamgConvertMeshx/BamgConvertMeshx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/BamgConvertMeshx/BamgConvertMeshx.cpp	(revision 12599)
+++ /issm/trunk-jpl/src/c/modules/BamgConvertMeshx/BamgConvertMeshx.cpp	(revision 12600)
@@ -12,30 +12,19 @@
 using namespace std;
 
-int BamgConvertMeshx(BamgMesh* bamgmesh,BamgGeom* bamggeom,double* index,double* x,double* y,int nods,int nels){
-
-	/*Intermediary*/
-	int i,j,k;
-	int verbose=0;
-	int noerr=1;
+int BamgConvertMeshx(BamgMesh* bamgmesh,BamgGeom* bamggeom,int* index,double* x,double* y,int nods,int nels){
 
 	/*Options*/
-	BamgOpts* bamgopts=NULL;
-	bamgopts=new BamgOpts();
+	BamgOpts* bamgopts=new BamgOpts();
 
-	// read mesh
-	if(verbose) _printLine_("Reading mesh");
+	/*read mesh*/
 	Mesh Th(index,x,y,nods,nels); 
 
-	//write mesh and geometry
-	if (verbose) _printLine_("Write Geometry");
+	/*write mesh and geometry*/
 	Th.Gh.WriteGeometry(bamggeom,bamgopts);
-	if (verbose) _printLine_("Write Mesh");
 	Th.WriteMesh(bamgmesh,bamgopts);
 
-	//clean up
+	/*clean up and return*/
 	delete bamgopts;
-
-	/*No error return*/
-	return noerr;
+	return 1;
 
 }
Index: /issm/trunk-jpl/src/c/modules/BamgConvertMeshx/BamgConvertMeshx.h
===================================================================
--- /issm/trunk-jpl/src/c/modules/BamgConvertMeshx/BamgConvertMeshx.h	(revision 12599)
+++ /issm/trunk-jpl/src/c/modules/BamgConvertMeshx/BamgConvertMeshx.h	(revision 12600)
@@ -9,5 +9,5 @@
 
 /* local prototypes: */
-int BamgConvertMeshx(BamgMesh* bamgmesh,BamgGeom* bamggeom,double* index,double* x,double* y,int nods,int nels);
+int BamgConvertMeshx(BamgMesh* bamgmesh,BamgGeom* bamggeom,int* index,double* x,double* y,int nods,int nels);
 
 #endif
Index: /issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.cpp	(revision 12599)
+++ /issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.cpp	(revision 12600)
@@ -13,6 +13,6 @@
 using namespace std;
 
-int InterpFromMeshToMesh2dx(double** pdata_interp,double* index_data,double* x_data,double* y_data,int nods_data,int nels_data,
-			double* data,int data_rows,int data_cols,double* x_interp,double* y_interp,int nods_interp,double* default_values,int num_default_values, DataSet* contours){
+int InterpFromMeshToMesh2dx(double** pdata_interp,int* index_data,double* x_data,double* y_data,int nods_data,int nels_data,
+			double* data,int M_data,int N_data,double* x_interp,double* y_interp,int N_interp,Options* options){
 	
 	/*Output*/
@@ -20,4 +20,6 @@
 
 	/*Intermediary*/
+	bool   isdefault;
+	double defaultvalue;
 	R2     r;
 	I2     I;
@@ -29,76 +31,56 @@
 	double data_value;
 	Icoor2 dete[3];
-	int verbose=0;
-
-	/*default values: */
-	Vector*    vec_incontour=NULL;
-	double*    incontour=NULL;
-	bool   skip_bamg=false;
 
 	/*Checks*/
-	if (data_cols<=0){
-		_error2_("data provided has a negative number of columns");
+	if (M_data!=nods_data && M_data!=nels_data){
+		_error2_("data provided should have either " << nods_data << " or " << nels_data << " lines (not " << M_data << ")");
 	}
-	if (data_rows!=nods_data && data_rows!=nels_data){
-		_error2_("data provided should have either " << nods_data << " or " << nels_data << " lines (not " << data_rows << ")");
+
+	/*Get default*/
+	if(options->GetOption("default")){
+		isdefault=true;
+		options->Get(&defaultvalue,"default");
 	}
-	if((num_default_values) && (data_cols>1)){
-		_error2_("data provided can only have 1 column if a default value is provided");
-	}
-	
-	/*If default values supplied, figure out which nodes are inside the contour, including the border of the contour: */
-	if(num_default_values){
-		ContourToNodesx( &vec_incontour,x_interp,y_interp,nods_interp,contours,1);
-		incontour=vec_incontour->ToMPISerial();
+	else{
+		isdefault=false;
 	}
 
 	/*Initialize output*/
-	if (verbose) _printLine_("Initializing output vector");
-	data_interp=xNew<double>(nods_interp*data_cols);
+	data_interp=xNew<double>(N_interp*N_data);
 
-	// read background mesh 
-	if (verbose) _printLine_("Reading mesh");
+	/*read background mesh*/
 	Mesh Th(index_data,x_data,y_data,nods_data,nels_data); 
 	Th.CreateSingleVertexToTriangleConnectivity();
 
-	//Loop over output nodes
-	if (verbose) _printLine_("Loop over the nodes");
-	for(i=0;i<nods_interp;i++){
-		
-		/*reset skip_bamg: */
-		skip_bamg=false;
+	/*Loop over output nodes*/
+	for(i=0;i<N_interp;i++){
 
-		/*figure out if we should skip bamg logic: */
-		if(num_default_values){
-			if(!incontour[i]){
-				/*This node is not inside the contour. Skip Bamg logic and apply default value.: */
-				skip_bamg=true;
+		//Get current point coordinates
+		r.x=x_interp[i]; r.y=y_interp[i];
+		I2 I=Th.R2ToI2(r);
+
+		//Find triangle holding r/I
+		Triangle &tb=*Th.TriangleFindFromCoord(I,dete);
+
+		// internal point 
+		if (tb.det>0){ 
+			//Area coordinate
+			areacoord[0]= (double) dete[0]/tb.det;
+			areacoord[1]= (double) dete[1]/tb.det;
+			areacoord[2]= (double) dete[2]/tb.det;
+			//3 vertices of the triangle
+			i0=Th.GetId(tb[0]);
+			i1=Th.GetId(tb[1]);
+			i2=Th.GetId(tb[2]);
+			//triangle number
+			it=Th.GetId(tb);
+		}
+		//external point
+		else{
+			if(isdefault){
+				for(j=0;j<N_data;j++) data_interp[i*N_data+j]=defaultvalue;
+				continue;
 			}
-		}
-
-		if(skip_bamg==false){
-
-			//Get current point coordinates
-			r.x=x_interp[i]; r.y=y_interp[i];
-			I2 I=Th.R2ToI2(r);
-
-			//Find triangle holding r/I
-			Triangle &tb=*Th.TriangleFindFromCoord(I,dete);
-
-			// internal point 
-			if (tb.det>0){ 
-				//Area coordinate
-				areacoord[0]= (double) dete[0]/tb.det;
-				areacoord[1]= (double) dete[1]/tb.det;
-				areacoord[2]= (double) dete[2]/tb.det;
-				//3 vertices of the triangle
-				i0=Th.GetId(tb[0]);
-				i1=Th.GetId(tb[1]);
-				i2=Th.GetId(tb[2]);
-				//triangle number
-				it=Th.GetId(tb);
-			}
-			//external point
-			else {
+			else{
 				//Get closest adjacent triangle (inside the mesh)
 				AdjacentTriangle ta=CloseBoundaryEdge(I,&tb,aa,bb).Adj();
@@ -116,27 +98,22 @@
 				it=Th.GetId(tc);
 			}
-			
-			if (data_rows==nods_data){
-				for (j=0;j<data_cols;j++){
-					data_interp[i*data_cols+j]=areacoord[0]*data[data_cols*i0+j]+areacoord[1]*data[data_cols*i1+j]+areacoord[2]*data[data_cols*i2+j];
-				}
-			}
-			else{
-				for (j=0;j<data_cols;j++){
-					if (it<0 || it>=nels_data){
-						_error2_("Triangle number " << it << " not in [0 " << nels_data << "], because not correctly implemented yet... interpolate on grid first");
-					}
-					data_interp[i*data_cols+j]=data[data_cols*it+j];
-				}
+		}
+
+		if (M_data==nods_data){
+			for (j=0;j<N_data;j++){
+				data_interp[i*N_data+j]=areacoord[0]*data[N_data*i0+j]+areacoord[1]*data[N_data*i1+j]+areacoord[2]*data[N_data*i2+j];
 			}
 		}
 		else{
-			if(num_default_values==1) data_interp[i]=default_values[0];
-			else data_interp[i]=default_values[i];
+			for (j=0;j<N_data;j++){
+				if (it<0 || it>=nels_data){
+					_error2_("Triangle number " << it << " not in [0 " << nels_data << "], report bug to developers");
+				}
+				data_interp[i*N_data+j]=data[N_data*it+j];
+			}
 		}
 	}
 
 	/*Assign output pointers:*/
-	if (verbose) _printLine_("Assigning output");
 	*pdata_interp=data_interp;
 
Index: /issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.h
===================================================================
--- /issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.h	(revision 12599)
+++ /issm/trunk-jpl/src/c/modules/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.h	(revision 12600)
@@ -6,9 +6,8 @@
 #define _INTERPFROMMESHTOMESH2DX_H
 
-#include "../../objects/objects.h"
+class Options;
 
-/* local prototypes: */
-int InterpFromMeshToMesh2dx(double** pdata_interp,double* index_data,double* x_data,double* y_data,int nods_data,int nels_data,
-			double* data,int data_rows,int data_cols,double* x_interp,double* y_interp,int nods_interp,double* default_values,int num_default_values,DataSet* contours);
+int InterpFromMeshToMesh2dx(double** pdata_interp,int* index_data,double* x_data,double* y_data,int nods_data,int nels_data,
+			double* data,int M_data,int N_data,double* x_interp,double* y_interp,int N_interp,Options* options);
 
 #endif
Index: /issm/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp	(revision 12599)
+++ /issm/trunk-jpl/src/c/modules/Krigingx/Krigingx.cpp	(revision 12600)
@@ -62,5 +62,4 @@
 		double *data  = NULL;
 		int    *index = NULL;
-		double *indexd= NULL;
 
 		observations->ObservationList(&x,&y,&data,&nobs);
@@ -68,11 +67,8 @@
 		_printLine_("Generation Delaunay Triangulation");
 		BamgTriangulatex(&index,&nel,x,y,nobs);
-		indexd =xNewZeroInit<double>(nel*3);
-		for(int i=0;i<nel*3;i++) indexd[i]=(double)index[i];
-		xDelete<int>(index);
 
 		_printLine_("Interpolating");
 		xDelete<double>(predictions);
-		InterpFromMeshToMesh2dx(&predictions,indexd,x,y,nobs,nel,data,nobs,1,x_interp,y_interp,n_interp,NULL,0,new DataSet());
+		InterpFromMeshToMesh2dx(&predictions,index,x,y,nobs,nel,data,nobs,1,x_interp,y_interp,n_interp,options);
 		xDelete<double>(x);
 		xDelete<double>(y);
Index: /issm/trunk-jpl/src/c/modules/TriaSearchx/TriaSearchx.cpp
===================================================================
--- /issm/trunk-jpl/src/c/modules/TriaSearchx/TriaSearchx.cpp	(revision 12599)
+++ /issm/trunk-jpl/src/c/modules/TriaSearchx/TriaSearchx.cpp	(revision 12600)
@@ -13,5 +13,5 @@
 using namespace std;
 
-void TriaSearchx(double** ptria,double* index,int nel, double* x, double* y, int nods,double* x0, double* y0,int numberofnodes){
+void TriaSearchx(double** ptria,int* index,int nel, double* x, double* y, int nods,double* x0, double* y0,int numberofnodes){
 
 	/*Output*/
Index: /issm/trunk-jpl/src/c/modules/TriaSearchx/TriaSearchx.h
===================================================================
--- /issm/trunk-jpl/src/c/modules/TriaSearchx/TriaSearchx.h	(revision 12599)
+++ /issm/trunk-jpl/src/c/modules/TriaSearchx/TriaSearchx.h	(revision 12600)
@@ -9,5 +9,5 @@
 
 /* local prototypes: */
-void TriaSearchx(double** ptria,double* index,int nel, double* x, double* y, int nods,double* x0, double* y0,int numberofnodes);
+void TriaSearchx(double** ptria,int* index,int nel, double* x, double* y, int nods,double* x0, double* y0,int numberofnodes);
 
 #endif
Index: /issm/trunk-jpl/src/c/objects/Bamg/Mesh.cpp
===================================================================
--- /issm/trunk-jpl/src/c/objects/Bamg/Mesh.cpp	(revision 12599)
+++ /issm/trunk-jpl/src/c/objects/Bamg/Mesh.cpp	(revision 12600)
@@ -41,6 +41,6 @@
 	}
 	/*}}}*/
-	/*FUNCTION Mesh::Mesh(double* index,double* x,double* y,int nods,int nels){{{*/
-	Mesh::Mesh(double* index,double* x,double* y,int nods,int nels):Gh(*(new Geometry())),BTh(*this){
+	/*FUNCTION Mesh::Mesh(int* index,double* x,double* y,int nods,int nels){{{*/
+	Mesh::Mesh(int* index,double* x,double* y,int nods,int nels):Gh(*(new Geometry())),BTh(*this){
 
 		Init(0);
@@ -265,6 +265,6 @@
 
 	/*IO*/
-	/*FUNCTION Mesh::ReadMesh(double* index,double* x,double* y,int nods,int nels){{{*/
-	void Mesh::ReadMesh(double* index,double* x,double* y,int nods,int nels){
+	/*FUNCTION Mesh::ReadMesh(int* index,double* x,double* y,int nods,int nels){{{*/
+	void Mesh::ReadMesh(int* index,double* x,double* y,int nods,int nels){
 
 		double Hmin = HUGE_VAL;// the infinie value 
Index: /issm/trunk-jpl/src/c/objects/Bamg/Mesh.h
===================================================================
--- /issm/trunk-jpl/src/c/objects/Bamg/Mesh.h	(revision 12599)
+++ /issm/trunk-jpl/src/c/objects/Bamg/Mesh.h	(revision 12600)
@@ -57,5 +57,5 @@
 			//Constructors/Destructors
 			Mesh(BamgGeom* bamggeom,BamgMesh* bamgmesh,BamgOpts* bamgopts);
-			Mesh(double* index,double* x,double* y,int nods,int nels);/*MeshConvert*/
+			Mesh(int* index,double* x,double* y,int nods,int nels);/*MeshConvert*/
 			Mesh(double* x,double* y,int nods); /*BamgTriangulate*/
 			Mesh(Mesh &,Geometry * pGh=0,Mesh* pBTh=0,long maxnbv_in=0 ); //copy operator
@@ -110,5 +110,5 @@
 			BamgVertex* NearestVertex(Icoor1 i,Icoor1 j) ;
 			Triangle* TriangleFindFromCoord(const I2 & ,Icoor2 [3],Triangle *tstart=0) const;
-			void ReadMesh(double* index,double* x,double* y,int nods,int nels);
+			void ReadMesh(int* index,double* x,double* y,int nods,int nels);
 			void ReadMesh(BamgMesh* bamgmesh, BamgOpts* bamgopts);
 			void WriteMesh(BamgMesh* bamgmesh,BamgOpts* bamgopts);