Index: /issm/trunk/src/c/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.cpp
===================================================================
--- /issm/trunk/src/c/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.cpp	(revision 2359)
+++ /issm/trunk/src/c/InterpFromMeshToMesh2dx/InterpFromMeshToMesh2dx.cpp	(revision 2360)
@@ -21,6 +21,8 @@
 	double area_1,area_2,area_3;
 	double data_value;
-	double xmin,xmax;
-	double ymin,ymax;
+	double x_prime_min,x_prime_max;
+	double y_prime_min,y_prime_max;
+	double x_tria_min,y_tria_min;
+	double x_tria_max,y_tria_max;
 
 	/*some checks*/
@@ -44,10 +46,10 @@
 
 	/*Get prime mesh extrema coordinates*/
-	xmin=x_prime[0]; xmax=x_prime[0];ymin=y_prime[0]; ymax=y_prime[0];
+	x_prime_min=x_prime[0]; x_prime_max=x_prime[0];y_prime_min=y_prime[0]; y_prime_max=y_prime[0];
 	for (i=1;i<nods_prime;i++){
-		if (x_prime[i]<xmin) xmin=x_prime[i];
-		if (x_prime[i]>xmax) xmax=x_prime[i];
-		if (y_prime[i]<ymin) ymin=y_prime[i];
-		if (y_prime[i]>ymax) ymax=y_prime[i];
+		if (x_prime[i]<x_prime_min) x_prime_min=x_prime[i];
+		if (x_prime[i]>x_prime_max) x_prime_max=x_prime[i];
+		if (y_prime[i]<y_prime_min) y_prime_min=y_prime[i];
+		if (y_prime[i]>y_prime_max) y_prime_max=y_prime[i];
 	}
 
@@ -64,8 +66,18 @@
 
 		/*if there is no point inside the domain, go to next iteration*/
-		if ( (x_data[(int)index_data[3*i+0]-1]<xmin) && (x_data[(int)index_data[3*i+1]-1]<xmin) && (x_data[(int)index_data[3*i+2]-1]<xmin)) continue;
-		if ( (x_data[(int)index_data[3*i+0]-1]>xmax) && (x_data[(int)index_data[3*i+1]-1]>xmax) && (x_data[(int)index_data[3*i+2]-1]>xmax)) continue;
-		if ( (y_data[(int)index_data[3*i+0]-1]<ymin) && (y_data[(int)index_data[3*i+1]-1]<ymin) && (y_data[(int)index_data[3*i+2]-1]<ymin)) continue;
-		if ( (y_data[(int)index_data[3*i+0]-1]>ymax) && (y_data[(int)index_data[3*i+1]-1]>ymax) && (y_data[(int)index_data[3*i+2]-1]>ymax)) continue;
+		if ( (x_data[(int)index_data[3*i+0]-1]<x_prime_min) && (x_data[(int)index_data[3*i+1]-1]<x_prime_min) && (x_data[(int)index_data[3*i+2]-1]<x_prime_min)) continue;
+		if ( (x_data[(int)index_data[3*i+0]-1]>x_prime_max) && (x_data[(int)index_data[3*i+1]-1]>x_prime_max) && (x_data[(int)index_data[3*i+2]-1]>x_prime_max)) continue;
+		if ( (y_data[(int)index_data[3*i+0]-1]<y_prime_min) && (y_data[(int)index_data[3*i+1]-1]<y_prime_min) && (y_data[(int)index_data[3*i+2]-1]<y_prime_min)) continue;
+		if ( (y_data[(int)index_data[3*i+0]-1]>y_prime_max) && (y_data[(int)index_data[3*i+1]-1]>y_prime_max) && (y_data[(int)index_data[3*i+2]-1]>y_prime_max)) continue;
+
+		/*Get extrema coordinates of current elements*/
+		x_tria_min=x_data[(int)index_data[3*i+0]-1]; x_tria_max=x_tria_min;
+		y_tria_min=y_data[(int)index_data[3*i+0]-1]; y_tria_max=y_tria_min;
+		for (j=1;j<3;j++){
+			if(x_data[(int)index_data[3*i+j]-1]<x_tria_min) x_tria_min=x_data[(int)index_data[3*i+j]-1];
+			if(x_data[(int)index_data[3*i+j]-1]>x_tria_max) x_tria_max=x_data[(int)index_data[3*i+j]-1];
+			if(y_data[(int)index_data[3*i+j]-1]<y_tria_min) y_tria_min=y_data[(int)index_data[3*i+j]-1];
+			if(y_data[(int)index_data[3*i+j]-1]>y_tria_max) y_tria_max=y_data[(int)index_data[3*i+j]-1];
+		}
 
 		/*get area of the current element (Jacobian = 2 * area)*/
@@ -77,4 +89,7 @@
 		/*loop over the prime nodes*/
 		for (j=0;j<nods_prime;j++){
+
+			/*if the current point is not in the triangle, continue*/
+			if ( (x_tria_min>x_prime[j]) || (x_tria_max<x_prime[j]) || (y_tria_min>y_prime[j]) || (y_tria_max<y_prime[j]) ) continue;
 
 			/*Get first area coordinate = det(x-x3  x2-x3 ; y-y3   y2-y3)/area*/