function [x_f y_f MatData]=TracksToMatrix(xmin,ymin,nx,ny,posting,track_coord,track_values,theta,varargin),
%TRACKSTOMATRIX - Create a map from tracks
%
%   This routine creates a map from tracks using a Kriging algorithm given by the package DACE.
%   - xmin and ymin are two scalars that are the extreme values of the data matrix one wants to create
%   - posting is a scalar giving the posting of the matrix to be created in meters
%   - trac_coord is an Argus file containing the coordinates of the tracks (x and y)
%   - trav_values is a vector with the values along the track coordinates
%   - theta is a parameter of the correlation function stricly positive
%   varagin can contain a matrix of data to constrain the matrix to be created
%   it shoud contain indices x and y and the matrix
%
%   Usage:
%      [x_ y_f MatData]=TracksToMatrix(xmin,ymin,nx,ny,posting,track_coord,track_values,theta),
%
%   Example:
%      [x_f y_f Thickness]=TracksToMatrix(-10^6,10^6,200,180,1000,'trackcoord.exp',thickness_values,theta),
%      [x_f y_f Thickness]=TracksToMatrix(-10^6,10^6,200,180,1000,'trackcoord.exp',thickness_values,theta,x_m,y_m,globalthickness),
%
%   See also: CREATEDATABOUNDARIES, CREATEDATAMATRIX

%some checks
if nargin~=8 & nargin~=11,
	error(' TracksToMatrix error message: wrong number of arguments')
end

%Add path to dace
addpath([issmdir() '/externalpackages/dace/install'])

%First create the x_m and y_m fot the matrix
x_f=[xmin:posting:xmin+posting*nx]';
y_f=[ymin:posting:ymin+posting*ny]';

%Now create a bigger map we will then truncate
x_m=[xmin-posting*nx/2:posting:xmin+posting*nx*3/2]';
y_m=[ymin-posting*ny/2:posting:ymin+posting*ny*3/2]';

%Create DataMatrix with local averaging of tracks
[Mvalue Mx My]=CreateDataMatrix(x_m,y_m,track_coord,track_values);

%Create vector for these coordinates and values
Svalue=sparse(Mvalue);
Sx=sparse(Mx);
Sy=sparse(My);
[i,j,Values]=find(Svalue);
[i,j,X]=find(Sx);
[i,j,Y]=find(Sy);
Coord=[X,Y];

%Create boundaries for the track if specified
if nargin==11,
	disp('      constain the border of matrix with global matrix');
	%Find out the elements to counstrain the border
	x_matrix=varargin{1};
	y_matrix=varargin{2};
	datamatrix=varargin{3};

	%Create values on the border of the smaller map
	[Xedge,Yedge,EdgeValues]=CreateDataBoundaries(x_m,y_m,x_matrix,y_matrix,datamatrix);

	%Create values on the border of the larger map
	[Xedge2,Yedge2,EdgeValues2]=CreateDataBoundaries(x_f,y_f,x_matrix,y_matrix,datamatrix);

	%Add these values to the track values
	Values=[Values;EdgeValues;EdgeValues2];
	Coord=[X,Y;Xedge,Yedge;Xedge2,Yedge2];
end

%Create model for data
disp('      create the statiscal model of the data');
[dmodel,perf]=dacefit(Coord,Values,@regpoly1,@corrgauss,theta);

%Create design site(points where we are looking for the data)
Points=nodesamp([x_f(1)+posting/2,y_f(1)+posting/2;x_f(end)-posting/2,y_f(end)-posting/2],[length(x_f)-1;length(y_f)-1]);

%Compute data on these points
disp('      predict the values');
VecData=predictor(Points,dmodel);

%Reshape to get a matrix
MatData=reshape(VecData,ny,nx);

%remove DACE path
rmpath([issmdir() '/externalpackages/dace/install']);
