clear all
close all

praw=fullfile('..','Raw');
pconv=fullfile('..','Processed');

load(fullfile(praw,'sorted_profiles'))
load(fullfile(praw,'origins'))

origin.x0=origin.x;
origin.y0=origin.y;

[x2,y2] = convertCoordinates(origin.x,origin.y,'CS1.code',27493,'CS2.code',4326);

origin.x=x2;
origin.y=y2;

fangle=deg2rad(90-(origin.azimuth-180))

outputfile=fullfile(pconv,'micore_ualg_topo.nc');

% find maximum length
mxl=max([profs.xl]);
  
dx=0.5;

% create crosshore and longshore variables
xx=0:dx:mxl;

yy=[0 500 1000 1500 2000 700];

nx1=length(xx);
ny1=length(yy);

% create time variables
tmp=[];

for i=1:length(profs)
    
    tmp=[tmp ;profs(i).sdates(:)];
    
end

sdates=unique(tmp);

nt=length(sdates);

% make the final elevation cube
inds='ABCDEP';

mat=nan(nx1,ny1,nt);
lonn=nan(nx1,ny1);
latt=nan(nx1,ny1);

error_log=[];

for ii=1:length(sdates)
    
    sd=sdates(ii);
    
    for i=1:length(profs)
    
        [~,ia,ib]=intersect(sd,profs(i).sdates);
        
        if ~isempty(ib)
            
            x1=profs(i).data(ib).xcext;
            y1=profs(i).data(ib).zcext;

            [~,ia2,ib2]=intersect(xx,x1);

    %         double check which profile it is
            ip=find(inds==profs(i).data(ib).profile);
            mat(ia2,ip,ii)=y1(ib2);
            
            tmpx1=origin.x0(ip)+x1.*sin(fangle);
            tmpy1=origin.y0(ip)+x1.*cos(fangle);

            [x3,y3] = convertCoordinates(tmpx1,tmpy1,'CS1.code',27493,'CS2.code',4326);

            lonn(ia2,ip,1)=x3(ib2);
            latt(ia2,ip,1)=y3(ib2);


    %         control
            if length(ib2)~=length(y1)

                error_log(end+1,:)=[ii,i];

            end
            
        end
        
    end
    
end

% quick visual check
for i=1:nt
    
    figure
    
    imagesc(mat(:,:,i))
    
end

nc_create_empty(outputfile);

nc_add_dimension(outputfile, 'longshore_distance', ny1);
longVariable = struct(...
    'Name', 'long_distance', ...
    'Nctype', 'float', ...
    'Dimension', {{'longshore_distance'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'metre') ...
    );
nc_addvar(outputfile, longVariable);

nc_add_dimension(outputfile, 'crosshore_distance', nx1);
crossVariable = struct(...
    'Name', 'cross_distance', ...
    'Nctype', 'float', ...
    'Dimension', {{'crosshore_distance'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'metre') ...
    );
nc_addvar(outputfile, crossVariable);

lonVariable = struct(...
    'Name', 'lon', ...
    'Nctype', 'float', ...
    'Dimension', {{'crosshore_distance','longshore_distance'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'degrees_west') ...
    );
nc_addvar(outputfile, lonVariable);

latVariable = struct(...
    'Name', 'lat', ...
    'Nctype', 'float', ...
    'Dimension', {{'crosshore_distance','longshore_distance'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'degrees_north') ...
    );
nc_addvar(outputfile, latVariable);

nc_add_dimension(outputfile, 'time', nt);
timeVariable = struct(...
    'Name', 'time', ...
    'Nctype', 'float', ...
    'Dimension', {{'time'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'serial_date') ...
    );
nc_addvar(outputfile, timeVariable);

heightVariable = struct(...
    'Name', 'elevation', ...
    'Nctype', 'float', ...
    'Dimension', {{'crosshore_distance','longshore_distance', 'time'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'metre') ...
    );
nc_addvar(outputfile, heightVariable);

azimuthVariable = struct(...
    'Name', 'azimuth', ...
    'Nctype', 'float', ...
    'Dimension', {{'longshore_distance'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'degrees') ...
    );
nc_addvar(outputfile, azimuthVariable);

x_originvariable = struct(...
    'Name', 'origin_lon', ...
    'Nctype', 'float', ...
    'Dimension', {{'longshore_distance'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'degrees_west') ...
    );
nc_addvar(outputfile, x_originvariable);

y_originvariable = struct(...
    'Name', 'origin_lat', ...
    'Nctype', 'float', ...
    'Dimension', {{'longshore_distance'}}, ...
    'Attribute', struct('Name', 'unit', 'Value', 'degrees_north') ...
    );
nc_addvar(outputfile, y_originvariable);

nc_varput(outputfile, 'elevation', mat)
nc_varput(outputfile, 'time', sdates)
nc_varput(outputfile, 'long_distance', yy)
nc_varput(outputfile, 'cross_distance', xx)
nc_varput(outputfile, 'lon', lonn)
nc_varput(outputfile, 'lat', latt)
nc_varput(outputfile, 'origin_lon', origin.x)
nc_varput(outputfile, 'origin_lat', origin.y)
nc_varput(outputfile, 'azimuth', origin.azimuth)


% surface(nc_varget(outputfile, 'height')')

nc_dump(outputfile)
