% define the directory where the raw files are stored
datadir=fullfile('..','raw');
% check which files are in the data directory
files=dir(fullfile(datadir, 'topo_all_2009_pre.txt'));
dataset=load(fullfile(datadir, files(1).name));

% the file contains 6 columns
id=dataset(:,1);
x92=dataset(:,2);
y92=dataset(:,3);
z=dataset(:,4);
x84=dataset(:,5);
y84=dataset(:,6);

% plot3(z,x92,y92,'ro');
% split the data into the different transects
% compute distances between sequential points
d = (x92(1:(end-1)) - x92(2:end)).^2 + (y92(1:(end-1)) - y92(2:end)).^2;
[b, ix] = sort(d);

%plot(x92,y92);
% check the number of transects visualy (plot them x,y) (or use some
% cutoff)
ntransects = 20; 
% total number of points
npointstotal = size(x92,1)
splitpoints=ix(end-(ntransects-1):end);
splitpoints=[sort(splitpoints)' npointstotal];
 
% Compute the maximum number of points per transect  
npoints = max(diff([0 splitpoints]))

% Create a structure for each transect. 
transects = struct('x92', [], 'y92', [], 'lat', [], 'lon', [], 'z',[]);
start = 1
for i=1:ntransects
    transects(i).x92 = x92(start:splitpoints(i));
    transects(i).y92 = y92(start:splitpoints(i));
    transects(i).z = z(start:splitpoints(i));
    transects(i).lat = y84(start:splitpoints(i));
    transects(i).lon = x84(start:splitpoints(i));
    start = splitpoints(i)+1
end

% this is the file we're saving results to 
filename=('topo_all_2009_pre.nc');
nc_create_empty(filename)
nc_add_dimension(filename,'points',npoints)
nc_add_dimension(filename,'transects',ntransects)

varstruct.Name='x92';
varstruct.Nctype=nc_double;
varstruct.Dimension={'transects','points'};
nc_addvar(filename,varstruct)
nc_attput(filename,'x92','units','meters')
nc_attput(filename, 'x92', 'standard_name', 'projection_x_coordinate')

varstruct.Name='y92';
varstruct.Nctype=nc_double;
varstruct.Dimension={'transects','points'};
nc_addvar(filename,varstruct)
nc_attput(filename,'y92','units','meters')
nc_attput(filename, 'y92', 'standard_name', 'projection_y_coordinate')

varstruct.Name='z';
varstruct.Nctype=nc_double;
varstruct.Dimension={'transects','points'};
nc_addvar(filename,varstruct)
nc_attput(filename,'z','units','meters')
nc_attput(filename, 'z', 'standard_name', 'mean_high')

varstruct.Name='lat';
varstruct.Nctype=nc_double;
varstruct.Dimension={'transects','points'};
nc_addvar(filename,varstruct)
nc_attput(filename,'lat','units','degrees_east')
nc_attput(filename, 'lat', 'standard_name', 'latitude')

varstruct.Name='lon';
varstruct.Nctype=nc_double;
varstruct.Dimension={'transects','points'};
nc_addvar(filename,varstruct)
nc_attput(filename,'lon','units','degrees_north')
nc_attput(filename, 'lon', 'standard_name', 'longitude')

for i=1:ntransects
    nc_varput(filename, 'x92', transects(i).x92, [i-1,0], [1,size(transects(i).x92,1)])
    nc_varput(filename, 'y92', transects(i).y92, [i-1,0], [1,size(transects(i).x92,1)])
    nc_varput(filename, 'z', transects(i).z, [i-1,0], [1,size(transects(i).x92,1)])
    nc_varput(filename, 'lon', transects(i).lon, [i-1,0], [1,size(transects(i).x92,1)])
    nc_varput(filename, 'lat', transects(i).lat, [i-1,0], [1,size(transects(i).x92,1)])
end

nc_dump(filename)


