Index: /issm/trunk/src/m/utils/Interp/InterpFromFile.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/InterpFromFile.m	(revision 1182)
+++ /issm/trunk/src/m/utils/Interp/InterpFromFile.m	(revision 1183)
@@ -48,16 +48,14 @@
 
 	%x -> name begins with "x" and is a vector
-	if (strcmpi(A(i).name(1),'x') & min(A(i).size)==1),
+	if (strncmpi(A(i).name,'x',1) & min(A(i).size)==1),
 		xenum=i;
 
 	%y -> name begins with "y" and is a vector
-	elseif (strcmpi(A(i).name(1),'y') & min(A(i).size)==1),
+	elseif (strncmpi(A(i).name,'y',1) & min(A(i).size)==1),
 		yenum=i;
 
 	%index -> name begins with "index" or "elements" and 3 columns
 	elseif length(A)==4
-		if (( (length(A(i).name)>=5 & strcmpi(A(i).name(1:5),'index')) ...
-				| (length(A(i).name)>=8 & strcmpi(A(i).name(1:8),'elements')) )...
-				& A(i).size(2)==3),
+		if ((strncmpi(A(i).name,'index',5) | strncmpi(A(i).name,'elements',8)) & A(i).size(2)==3),
 			indexenum=i;
 		end
@@ -200,5 +198,5 @@
 
 	%interpolate
-	data_out=InterpFromMesh2d(index_data,x_data,y_data,data,x,y,default_value);
+	data_out=InterpFromMesh2d(index_data,x_data(:),y_data(:),data(:),x(:),y(:),default_value);
 
 else
@@ -218,4 +216,4 @@
 
 	%interpolate
-	data_out=InterpFromGrid(x_data,y_data,data,x,y,default_value);
-end
+	data_out=InterpFromGrid(x_data(:),y_data(:),data,x(:),y(:),default_value);
+end
Index: sm/trunk/src/m/utils/Interp/griddata_grid_to_mesh.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_grid_to_mesh.m	(revision 1182)
+++ 	(revision )
@@ -1,95 +1,0 @@
-function zp=griddata_grid_to_mesh(x_m,y_m,u,index,x,y);
-%GRIDDATA_GRID_TO_MESH - interpolates data from a grid to a mesh   
-%
-%   Usage:
-%      zp=griddata_grid_to_mesh(x_m,y_m,u,index,x,y)
-
-pos=find(isnan(u));
-for n=1:10,
-   pos=find(isnan(u));
-   u(pos)=u(pos+1);
-end
-
-zp=zeros(length(x),1);
-for n=1:length(x),
-   if mod(n,100)==0,
-      disp(['Count is ' num2str(n/length(x)*100)]);
-   end
-   
-   pos1=find((x_m-x(n))<0);
-   pos2=find((x_m-x(n))>0);
-   pos3=find((y_m-y(n))<0);
-   pos4=find((y_m-y(n))>0);
-   
-   if (isempty(pos1) | isempty(pos2) | isempty(pos3) | isempty(pos4)),
-      out_of_bounds=1;
-   else
-      out_of_bounds=0;
-   end
-   
-   if out_of_bounds==0,   
-   %first find where (x(n),y(n)) lies.
-   pos1=find((x_m-x(n))<0);
-   pos2=find((x_m-x(n))>0);
-   x_bound=[pos1(length(pos1)) pos2(1)];   
-   
-   pos1=find((y_m-y(n))<0);
-   pos2=find((y_m-y(n))>0);
-   y_bound=[pos1(length(pos1)) pos2(1)];
-  
-   
-   
-   %we look for a linear interpolation, therefore have to decide 
-   % in which triangle our point lies.
-   A=[x_m(x_bound(1)) y_m(y_bound(1))]';
-   u_A=u(x_bound(1),y_bound(1));
-   B=[x_m(x_bound(2)) y_m(y_bound(1))]';
-   u_B=u(x_bound(2),y_bound(1));
-   C=[x_m(x_bound(2)) y_m(y_bound(2))]';
-   u_C=u(x_bound(2),y_bound(2));
-   D=[x_m(x_bound(1)) y_m(y_bound(2))]';
-   u_D=u(x_bound(1),y_bound(2));
-   
-   
-   
-   res=(D(2)-B(2))*x(n)+(B(1)-D(1))*y(n)+D(1)*B(2)-D(2)*B(1);
-   if res>0,
-      xx=[B(1) C(1) D(1)];
-      yy=[B(2) C(2) D(2)];
-      uu=[u_B u_C u_D];
-      ind=[1 2 3];
-   else
-      xx=[A(1) B(1) D(1)];
-      yy=[A(2) B(2) D(2)];
-      uu=[u_A u_B u_D];
-      ind=[1 2 3];
-   end
-   
-   
-	alpha=zeros(1,3);
-	beta=zeros(1,3);
-	gamma=zeros(1,3);
-	X=inv([xx(ind(1,:))' yy(ind(1,:))' ones(3,1)]);
-	alpha(1,:)=X(1,:);
-	beta(1,:)=X(2,:);
-	gamma(1,:)=X(3,:);
-   
-   if length(find(isnan(uu)))==0,
-      zp(n)=sum(uu.*(alpha.*x(n)+beta.*y(n)+gamma));
-   else
-      if length(find(~isnan(uu)))==0,
-         zp(n)=NaN;
-         else
-         pos=find(~isnan(uu));
-         zp(n)=uu(pos(1));
-      end
-   end
-   
-   else
-   zp(n)=NaN;
-end
-
-   
-end
-%disp('Filling final gaps');
-%zp=void_fill(zp,x,y);
Index: sm/trunk/src/m/utils/Interp/griddata_grid_to_mesh_completer.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_grid_to_mesh_completer.m	(revision 1182)
+++ 	(revision )
@@ -1,98 +1,0 @@
-function zp=griddata_grid_to_mesh_completer(x_m,y_m,u,index,x,y,zp);
-%GRIDDATA_GRID_TO_MESH_COMPLETER - interpolates data from a grid to a mesh   
-%
-%   Usage:
-%      zp=griddata_grid_to_mesh_completer(x_m,y_m,u,index,x,y,zp)
-
-pos=find(isnan(u));
-for n=1:100,
-   pos=find(isnan(u));
-   u(pos)=u(pos+1);
-end
-
-imagesc(u)
-break
-
-zp=zeros(length(x),1);
-for n=1:length(x),
-   if mod(n,100)==0,
-      disp(['Count is ' num2str(n/length(x)*100)]);
-   end
-   
-   pos1=find((x_m-x(n))<0);
-   pos2=find((x_m-x(n))>0);
-   pos3=find((y_m-y(n))<0);
-   pos4=find((y_m-y(n))>0);
-   
-   if (isempty(pos1) | isempty(pos2) | isempty(pos3) | isempty(pos4)),
-      out_of_bounds=1;
-   else
-      out_of_bounds=0;
-   end
-   
-   if out_of_bounds==0,   
-   %first find where (x(n),y(n)) lies.
-   pos1=find((x_m-x(n))<0);
-   pos2=find((x_m-x(n))>0);
-   x_bound=[pos1(length(pos1)) pos2(1)];   
-   
-   pos1=find((y_m-y(n))<0);
-   pos2=find((y_m-y(n))>0);
-   y_bound=[pos1(length(pos1)) pos2(1)];   
-   
-   
-   %we look for a linear interpolation, therefore have to decide 
-   % in which triangle our point lies.
-   
-   A=[x_m(x_bound(1)) y_m(y_bound(1))]';
-   u_A=u(x_bound(1),y_bound(1));
-   B=[x_m(x_bound(2)) y_m(y_bound(1))]';
-   u_B=u(x_bound(2),y_bound(1));
-   C=[x_m(x_bound(2)) y_m(y_bound(2))]';
-   u_C=u(x_bound(2),y_bound(2));
-   D=[x_m(x_bound(1)) y_m(y_bound(2))]';
-   u_D=u(x_bound(1),y_bound(2));
-   
-   
-   
-   res=(D(2)-B(2))*x(n)+(B(1)-D(1))*y(n)+D(1)*B(2)-D(2)*B(1);
-   if res>0,
-      xx=[B(1) C(1) D(1)];
-      yy=[B(2) C(2) D(2)];
-      uu=[u_B u_C u_D];
-      ind=[1 2 3];
-   else
-      xx=[A(1) B(1) D(1)];
-      yy=[A(2) B(2) D(2)];
-      uu=[u_A u_B u_D];
-      ind=[1 2 3];
-   end
-   
-   
-	alpha=zeros(1,3);
-	beta=zeros(1,3);
-	gamma=zeros(1,3);
-	X=inv([xx(ind(1,:))' yy(ind(1,:))' ones(3,1)]);
-	alpha(1,:)=X(1,:);
-	beta(1,:)=X(2,:);
-	gamma(1,:)=X(3,:);
-   
-   if length(find(isnan(uu)))==0,
-      zp(n)=sum(uu.*(alpha.*x(n)+beta.*y(n)+gamma));
-   else
-      if length(find(~isnan(uu)))==0,
-         zp(n)=NaN;
-         else
-         pos=find(~isnan(uu));
-         zp(n)=uu(pos(1));
-      end
-   end
-   
-   else
-   zp(n)=NaN;
-end
-
-   
-end
-%disp('Filling final gaps');
-%zp=void_fill(zp,x,y);
Index: sm/trunk/src/m/utils/Interp/griddata_grid_to_mesh_nan.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_grid_to_mesh_nan.m	(revision 1182)
+++ 	(revision )
@@ -1,90 +1,0 @@
-function zp=griddata_grid_to_mesh_nan(x_m,y_m,u,index,x,y);
-%GRIDDATA_GRID_TO_MESH_NAN - interpolates data from a grid to a mesh   
-%
-%   Usage:
-%      zp=griddata_grid_to_mesh_nan(x_m,y_m,u,index,x,y)
-
-zp=zeros(length(x),1);
-for n=1:length(x),
-   if mod(n,100)==0,
-      disp(['Count is ' num2str(n/length(x)*100)]);
-   end
-   
-   pos1=find((x_m-x(n))<0);
-   pos2=find((x_m-x(n))>0);
-   pos3=find((y_m-y(n))<0);
-   pos4=find((y_m-y(n))>0);
-   
-   if (isempty(pos1) | isempty(pos2) | isempty(pos3) | isempty(pos4)),
-      out_of_bounds=1;
-   else
-      out_of_bounds=0;
-   end
-   
-   if out_of_bounds==0,   
-   %first find where (x(n),y(n)) lies.
-   pos1=find((x_m-x(n))<0);
-   pos2=find((x_m-x(n))>0);
-   x_bound=[pos1(length(pos1)) pos2(1)];   
-   
-   pos1=find((y_m-y(n))<0);
-   pos2=find((y_m-y(n))>0);
-   y_bound=[pos1(length(pos1)) pos2(1)];   
-   
-   
-   %we look for a linear interpolation, therefore have to decide 
-   % in which triangle our point lies.
-   
-   A=[x_m(x_bound(1)) y_m(y_bound(1))]';
-   u_A=u(x_bound(1),y_bound(1));
-   B=[x_m(x_bound(2)) y_m(y_bound(1))]';
-   u_B=u(x_bound(2),y_bound(1));
-   C=[x_m(x_bound(2)) y_m(y_bound(2))]';
-   u_C=u(x_bound(2),y_bound(2));
-   D=[x_m(x_bound(1)) y_m(y_bound(2))]';
-   u_D=u(x_bound(1),y_bound(2));
-   
-   
-   
-   res=(D(2)-B(2))*x(n)+(B(1)-D(1))*y(n)+D(1)*B(2)-D(2)*B(1);
-   if res>0,
-      xx=[B(1) C(1) D(1)];
-      yy=[B(2) C(2) D(2)];
-      uu=[u_B u_C u_D];
-      ind=[1 2 3];
-   else
-      xx=[A(1) B(1) D(1)];
-      yy=[A(2) B(2) D(2)];
-      uu=[u_A u_B u_D];
-      ind=[1 2 3];
-   end
-   
-   
-	alpha=zeros(1,3);
-	beta=zeros(1,3);
-	gamma=zeros(1,3);
-	X=inv([xx(ind(1,:))' yy(ind(1,:))' ones(3,1)]);
-	alpha(1,:)=X(1,:);
-	beta(1,:)=X(2,:);
-	gamma(1,:)=X(3,:);
-   
-   if length(find(isnan(uu)))==0,
-      zp(n)=sum(uu.*(alpha.*x(n)+beta.*y(n)+gamma));
-   else
-      if length(find(~isnan(uu)))==0,
-         zp(n)=NaN;
-         else
-         pos=find(~isnan(uu));
-         zp(n)=uu(pos(1));
-      end
-   end
-   
-   else
-   zp(n)=NaN;
-
-end
-
-   
-end
-%disp('Filling final gaps');
-%zp=void_fill(zp,x,y);
Index: sm/trunk/src/m/utils/Interp/griddata_mesh_to_grid.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_mesh_to_grid.m	(revision 1182)
+++ 	(revision )
@@ -1,49 +1,0 @@
-function grid_data=griddata_mesh_to_grid(x_m,y_m,grid_data,x,y,mesh_data,varargin)
-%GRIDDATA_MESH_TO_GRID - interpolates data from a mesh to a grid
-%
-%   This routine takes a mesh (x,y) and values for every grid of the mesh (mesh_data),
-%   and interpolates it (using nearest interpolation) onto the gridded data (grid_data). 
-%   x_m and y_m (vectors of coordinates) are used to find interpolation locations onto grids.
-%
-%   Usage:
-%      grid_data=griddata_mesh_to_grid(x_m,y_m,grid_data,x,y,mesh_data,varargin)
-
-if nargin==7,
-	roundup=varargin{1};
-else
-	roundup=0;
-end
-
-s=size(grid_data);
-if length(x_m) ~= (s(2)+1),
-	error('x_m should be wider than grid_data by 1');
-end
-
-if length(y_m) ~= (s(1)+1),
-	error('y_m should be larger than grid_data by 1');
-end
-
-x_m0=x_m(1:end-1);
-x_m1=x_m(2:end);
-y_m0=y_m(1:end-1);
-y_m1=y_m(2:end);
-
-for i=1:length(x),
-	if mod(i,100000)==0,
-		disp(['   progress ' num2str(i/length(x)*100) '%']);
-	end
-	posx=find( (x(i)-x_m1)<=0 & (x(i)-x_m0)>0 );
-	posy=find( (y(i)-y_m1)<=0 & (y(i)-y_m0)>0 );
-	
-	y0=posy-roundup;
-	y1=posy+roundup;
-	if y0<=0, y0=1; end
-	if y1>s(1), y1=s(1);end
-
-	x0=posx-roundup;
-	x1=posx+roundup;
-	if x0<=0, x0=1; end
-	if x1>s(2), x1=s(2);end
-
-	grid_data(y0:y1,x0:x1)=mesh_data(i);
-end
Index: sm/trunk/src/m/utils/Interp/griddata_mesh_to_mesh.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_mesh_to_mesh.m	(revision 1182)
+++ 	(revision )
@@ -1,99 +1,0 @@
-function data_prime=griddata_mesh_to_mesh(index,x,y,data,x_prime,y_prime,varargin);
-%GRIDDATA_MESH_TO_MESH - interpolates data from a  2d mesh to a 2d mesh
-%
-%   griddata_mesh_to_mesh: interpolate data (from mesh with triangulation index,x,y)
-%   into mesh with triangulation (x_prime,y_prime).
-%
-%   Usage:
-%      data_prime=griddata_mesh_to_mesh(index,x,y,data,x_prime,y_prime,varargin)
-
-	global client_server_mode IMDataCounter uset
-
-	%restrict mesh_prime grids to the ones which are inside the mesh hull 
-	in=inhull([x_prime y_prime],[x y]);
-	pos_hull=find(in);
-
-	%keep only mesh_prime grids inside the hull
-	x_prime_hull=x_prime(pos_hull);
-	y_prime_hull=y_prime(pos_hull);
-
-	%call griddata_mesh_to_mesh_hulled on this new sub mesh from mesh 2
-	if strcmpi(client_server_mode,'yes'),
-		r_data_prime_hull=GriddataMeshToMeshHulled_client(index,x,y,data,x_prime_hull,y_prime_hull);
-		data_prime_hull=IMdb('select matrix from r_data_prime_hull');
-		IMdb('drop r_data_prime_hull');
-	else
-		data_prime_hull=griddata_mesh_to_mesh_hulled(index,x,y,data,x_prime_hull,y_prime_hull,varargin);
-	end
-
-	%return output
-	data_prime=zeros(length(x_prime),1); data_prime(:)=NaN;
-	data_prime(pos_hull)=data_prime_hull;
-end
-
-function data_prime=griddata_mesh_to_mesh_hulled(index,x,y,data,x_prime,y_prime,varargin);
-
-	%Analyse data: is it given per node or per element?
-	if isempty(varargin{1})
-		element_data=(length(data)==size(index,1));%data is given for a given element
-		grid_data=(length(data)==length(x));       %data is given for a given grid
-		if element_data==0 & grid_data==0, error('griddata_mesh_to_mesh_hulled error message: unknown size of data, specify element or node'), end
-	else
-		if strcmpi(varargin{1},'node')
-			element_data=0;
-			grid_data=1;
-		elseif strcmpi(varargin{1},'element')
-			element_data=1;
-			grid_data=0;
-		else
-			error('griddata_mesh_to_mesh_hulled error message: wrong input argument. Must be ''node'' or ''element''')
-		end
-	end
-
-	%Some parameters and initialization of data_prime
-	nel=length(index);
-	data_prime=zeros(length(x_prime),1); data_prime(:)=NaN;
-
-	%compute some quantities that will speed up the process
-	x3x1=x(index(:,1))-x(index(:,3));
-	y3y1=y(index(:,1))-y(index(:,3));
-	x3x2=x(index(:,2))-x(index(:,3));
-	y3y2=y(index(:,2))-y(index(:,3));
-	x3=x(index(:,3));
-	y3=y(index(:,3));
-
-	%triangles jacobian determinant (2*area)
-	delta=x(index(:,2)).*y(index(:,3))-y(index(:,2)).*x(index(:,3))-x(index(:,1)).*y(index(:,3))+y(index(:,1)).*x(index(:,3))+x(index(:,1)).*y(index(:,2))-y(index(:,1)).*x(index(:,2));
-
-	if nel>100,
-		fprintf('%s','      interpolation progress:   0.00 %');
-	end
-
-	%Get area coordinates:
-	for n=1:nel,
-		if mod(n,100)==0,
-			fprintf('\b\b\b\b\b\b\b')
-			fprintf('%5.2f%s',n/nel*100,' %');
-		end
-		%first area coordinate
-		area_1=(y3y2(n)*(x_prime-x3(n))-x3x2(n)*(y_prime-y3(n)))/delta(n);
-		%second area coordinate
-		area_2=(x3x1(n)*(y_prime-y3(n))-y3y1(n)*(x_prime-x3(n)))/delta(n);
-		%third area coordinate
-		area_3=1-area_1-area_2;
-		%get grids in mesh 2 that have all area coordinates positive (meaning they belong to mesh 1 triangle n:
-		pos=find((area_1>=0) & (area_2>=0) & (area_3>=0));
-		if ~isempty(pos),
-			if grid_data
-				%mesh 2 grids in pos belong to mesh 1 element n. interpolate value linearly.
-				data_prime(pos)=area_1(pos)*data(index(n,1))+area_2(pos)*data(index(n,2))+area_3(pos)*data(index(n,3));
-			elseif element_data
-				data_prime(pos)=data(n);
-			end
-		end
-	end
-	if nel>100,
-		fprintf('\b\b\b\b\b\b\b\b')
-		fprintf('%4.2f%s\n',100,' %');
-	end
-end
Index: sm/trunk/src/m/utils/Interp/griddata_mesh_to_mesh_3d.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_mesh_to_mesh_3d.m	(revision 1182)
+++ 	(revision )
@@ -1,113 +1,0 @@
-function data_prime=griddata_mesh_to_mesh_3d(index,x,y,z,data,x_prime,y_prime,z_prime,varargin);
-%GRIDDATA_MESH_TO_MESH_3D - interpolates data from a 3d mesh to a 3d mesh
-%
-%   griddata_mesh_to_mesh: interpolate data (from mesh with triangulation index,x,y,z)
-%   into mesh with triangulation (x_prime,y_prime,z_prime).
-%
-%   Usage:
-%      data_prime=griddata_mesh_to_mesh_3d(index,x,y,z,data,x_prime,y_prime,z_prime,varargin)
-
-	global client_server_mode IMDataCounter uset
-
-	%restrict mesh_prime grids to the ones which are inside the mesh hull 
-	in=inhull([x_prime y_prime z_prime],[x y z]);
-	pos_hull=find(in);
-
-	%keep only mesh_prime grids inside the hull
-	x_prime_hull=x_prime(pos_hull);
-	y_prime_hull=y_prime(pos_hull);
-	z_prime_hull=z_prime(pos_hull);
-
-	%call griddata_mesh_to_mesh_hulled on this new sub mesh from mesh 2
-	if strcmpi(client_server_mode,'yes'),
-		r_data_prime_hull=GriddataMeshToMeshHulled_client(index,x,y,data,x_prime_hull,y_prime_hull);
-		data_prime_hull=IMdb('select matrix from r_data_prime_hull');
-		IMdb('drop r_data_prime_hull');
-	else
-		data_prime_hull=griddata_mesh_to_mesh_hulled(index,x,y,z,data,x_prime_hull,y_prime_hull,z_prime_hull,varargin);
-	end
-
-	%return output
-	data_prime=zeros(length(x_prime),1); data_prime(:)=NaN;
-	data_prime(pos_hull)=data_prime_hull;
-end
-
-function data_prime=griddata_mesh_to_mesh_hulled(index,x,y,z,data,x_prime,y_prime,z_prime,varargin);
-
-	%Analyse data: is it given per node or per element?
-	if isempty(varargin{1})
-		element_data=(length(data)==size(index,1));%data is given for a given element
-		grid_data=(length(data)==length(x));       %data is given for a given grid
-		if element_data==0 & grid_data==0, error('griddata_mesh_to_mesh_hulled error message: unknown size of data, specify element or node'), end
-	else
-		if strcmpi(varargin{1},'node')
-			element_data=0;
-			grid_data=1;
-		elseif strcmpi(varargin{1},'element')
-			element_data=1;
-			grid_data=0;
-		else
-			error('griddata_mesh_to_mesh_hulled error message: wrong input argument. Must be ''node'' or ''element''')
-		end
-	end
-
-	%Some parameters and initialization of data_prime
-	nel=length(index);
-	data_prime=zeros(length(x_prime),1); data_prime(:)=NaN;
-
-	%compute some quantities that will speed up the process
-	x3x1=x(index(:,1))-x(index(:,3));
-	y3y1=y(index(:,1))-y(index(:,3));
-	x3x2=x(index(:,2))-x(index(:,3));
-	y3y2=y(index(:,2))-y(index(:,3));
-	x3=x(index(:,3));
-	y3=y(index(:,3));
-
-	%triangles jacobian determinant (2*area)
-	delta=x(index(:,2)).*y(index(:,3))-y(index(:,2)).*x(index(:,3))-x(index(:,1)).*y(index(:,3))+y(index(:,1)).*x(index(:,3))+x(index(:,1)).*y(index(:,2))-y(index(:,1)).*x(index(:,2));
-
-	if nel>100,
-		fprintf('%s','      interpolation progress:   0.00 %');
-	end
-
-	%Get area coordinates:
-	for n=1:nel,
-		if mod(n,100)==0,
-			fprintf('\b\b\b\b\b\b\b')
-			fprintf('%5.2f%s',n/nel*100,' %');
-		end
-		%first area coordinate
-		area_1=(y3y2(n)*(x_prime-x3(n))-x3x2(n)*(y_prime-y3(n)))/delta(n);
-		%second area coordinate
-		area_2=(x3x1(n)*(y_prime-y3(n))-y3y1(n)*(x_prime-x3(n)))/delta(n);
-		%third area coordinate
-		area_3=1-area_1-area_2;
-
-		%get grids in mesh prime that have all area coordinates positive (meaning they are on the same extruded column):
-		pos_horiz=find((area_1>=0) & (area_2>=0) & (area_3>=0));
-		if ~isempty(pos_horiz),
-
-			%Figure out in which extruded penta it is located
-			lower=area_1(pos_horiz)*z(index(n,1))+area_2(pos_horiz)*z(index(n,2))+area_3(pos_horiz)*z(index(n,3));
-			upper=area_1(pos_horiz)*z(index(n,4))+area_2(pos_horiz)*z(index(n,5))+area_3(pos_horiz)*z(index(n,6));
-			xi=2*(z_prime(pos_horiz)-lower)./(upper-lower)-1;
-			pos=find(z_prime(pos_horiz)>=lower & z_prime(pos_horiz)<=upper);
-			pos_vert=pos_horiz(pos);
-			xi=xi(pos);
-
-			if ~isempty(pos_vert)
-				if grid_data
-					%mesh 2 grids in pos belong to mesh 1 element n. interpolate value linearly.
-					data_prime(pos_vert)=area_1(pos_vert).*(1-xi)/2*data(index(n,1))+area_2(pos_vert).*(1-xi)/2*data(index(n,2))+area_3(pos_vert).*(1-xi)/2*data(index(n,3))...
-							+area_1(pos_vert).*(1+xi)/2*data(index(n,4))+area_2(pos_vert).*(1+xi)/2*data(index(n,5))+area_3(pos_vert).*(1+xi)/2*data(index(n,6));
-				elseif element_data
-					data_prime(pos_vert)=data(n);
-				end
-			end
-		end
-	end
-	if nel>100,
-		fprintf('\b\b\b\b\b\b\b\b')
-		fprintf('%4.2f%s\n',100,' %');
-	end
-end
Index: sm/trunk/src/m/utils/Interp/griddata_nearest.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_nearest.m	(revision 1182)
+++ 	(revision )
@@ -1,24 +1,0 @@
-function [u]=griddata_nearest(x_m,y_m,u_m,x,y)
-%GRIDDATA_NEAREST - returns interpolated u_m function of x_m and y_m  on a grid
-%
-%   This function returns interpolated u_m function of x_m and y_m,
-%   on the grid defined by x and y.
-%   See griddata_nearest_matrix for the case where x_m,y_m and u_m are matrices 
-%   and x,y are vectores.
-%   See griddata_perso for a user defined version of griddata ( with fuzz coeff).
-%
-%   Usage:
-%      [u]=griddata_nearest(x_m,y_m,u_m,x,y)
-
-u=zeros(length(x),1);
-count=0;
-for n=1:length(x),
-   
-   if n>length(x)/10*count,
-      disp(n/length(x)*100);
-      count=count+1;
-   end
-   
-   [posx,posy]=find_coord(x_m,y_m,x(n),y(n));
-   u(n)=u_m(posx,posy);
-end
Index: sm/trunk/src/m/utils/Interp/griddata_nearest_matrix.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_nearest_matrix.m	(revision 1182)
+++ 	(revision )
@@ -1,25 +1,0 @@
-function [u]=griddata_nearest_matrix(x_m,y_m,u_m,x,y)
-%GRIDDATA_NEAREST_MATRIX - interpolate a function on a grid
-%
-%   This function returns interpolated u_m function of x_m and y_m,
-%   on the grid defined by x and y.
-%   INPUT x_m,y_m,u_M, x,y where x_m,y_m and u_m are 2d matrices, 
-%   and x,y are one dimensional vectors,
-%   See griddata_nearest_matrix for the case where x_m,y_m and u_m are matrices 
-%   and x,y are vectores.
-%   See griddata_perso for a user defined version of griddata ( with fuzz coeff).
-%
-%   Usage:
-%      [u]=griddata_nearest_matrix(x_m,y_m,u_m,x,y)
-
-u=zeros(length(x),1);
-count=0;
-for n=1:length(x),
-    if n>length(x)/10*count,
-      disp(n/length(x)*100);
-      count=count+1;
-   end
-
-   [posx,posy]=find_coord_matrix(x_m,y_m,x(n),y(n));
-   u(n)=u_m(posx,posy);
-end
Index: sm/trunk/src/m/utils/Interp/griddata_nearest_nan.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_nearest_nan.m	(revision 1182)
+++ 	(revision )
@@ -1,21 +1,0 @@
-function [u]=griddata_nearest_nan(x_m,y_m,u_m,x,y)
-%GRIDDATA_NEAREST_NAN - interpolates a function on a grid
-%
-%   This function returns interpolated u_m function of x_m and y_m,
-%   on the grid defined by x and y.
-%   INPUT x_m,y_m,u_M, x,y where x_m,y_m, u_m, x,y are all one dimensional vectors,
-%   See griddata_nearest_matrix for the case where x_m,y_m and u_m are matrices 
-%   and x,y are vectores.
-%   See griddata_perso for a user defined version of griddata ( with fuzz coeff).
-%
-%   Usage:
-%      [u]=griddata_nearest_nan(x_m,y_m,u_m,x,y)
-
-u=zeros(length(x),1);
-
-pos=find(~isnan(u_m));
-
-for n=1:length(x),
-   [posx,posy]=find_coord(x_m,y_m,x(n),y(n));
-   u(n)=u_m(posx,posy);
-end
Index: sm/trunk/src/m/utils/Interp/griddata_nearest_vector.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_nearest_vector.m	(revision 1182)
+++ 	(revision )
@@ -1,28 +1,0 @@
-function [u]=griddata_nearest_vector(x_m,y_m,u_m,x,y)
-%GRIDDATA_NEAREST_VECTOR - interpolates a function on a grid
-%
-%   This function returns interpolated u_m function of x_m and y_m,
-%   on the grid defined by x and y.
-%   INPUT x_m,y_m,u_M, x,y where x_m,y_m, u_m, x,y are all one dimensional vectors,
-%   See griddata_nearest_matrix for the case where x_m,y_m and u_m are matrices 
-%   and x,y are vectores.
-%   See griddata_perso for a user defined version of griddata ( with fuzz coeff).
-%
-%   Usage:
-%      [u]=griddata_nearest_vector(x_m,y_m,u_m,x,y)
-
-u=zeros(length(x),1);
-for n=1:length(x),
-   if mod(n,length(x)/10)==0,
-      disp(n/length(x)*100);
-   end
-   
-   dist=sqrt((x(n)-x_m).^2+(y(n)-y_m).^2);
-   pos=find(dist==min(dist));
-   if isnan(u_m(pos)),
-      pos
-   end
-   
-   u(n)=u_m(pos);
-   
-end
Index: sm/trunk/src/m/utils/Interp/griddata_perso.m
===================================================================
--- /issm/trunk/src/m/utils/Interp/griddata_perso.m	(revision 1182)
+++ 	(revision )
@@ -1,26 +1,0 @@
-function zp=griddata_perso(index,x,y,u,xp,yp,alpha,beta,gamma);
-%GRIDDATA_PERSO - interpolates linearly a profile on a mesh
-%
-%   Usage:
-%      zp=griddata_perso(index,x,y,u,xp,yp,alpha,beta,gamma)
-
-zp=zeros(length(xp),1);
-for n=1:length(xp),
-   P=[xp(n) yp(n)];
-   dist=sqrt((x-P(1)).^2+(y-P(2)).^2);
-   posQ=find(dist==min(dist));
-   if min(dist)==0,
-      zp(n)=u(n);
-   else
-   	posQ=posQ(1);
-      elems=find(index(:,1)==posQ | index(:,2)==posQ | index(:,3)==posQ);
-   	for i=1:length(elems),
-         el=elems(i);
-         if isinpoly(P(1),P(2),x(index(el,:)),y(index(el,:)))==1,
-            zp(n)=sum((alpha(el,:)).*u(index(el,:))')*P(1)+...
-            sum((beta(el,:)).*u(index(el,:))')*P(2)+...
-            sum((gamma(el,:)).*u(index(el,:))');
-			end
-      end
-   end
-end