clear all
theDir = 'd:\checkouts\kustlijnkaartenboek\raw\';
files = dir([theDir 'tkl*.xls']);


%% read XLS files
KLdat = [];
for ii = 1:length(files);
    disp(['processing ' files(ii).name])
    [data, txt, rawData] = xlsread([theDir files(ii).name]);
    text{ii,1}           = files(ii).name(4:7);
    text(ii,2:7)         = txt(1,1:6);
    year                 = str2double(files(ii).name(4:7));
    KL(:,2) = data(:,1); % 2 'section'
    KL(:,3) = data(:,2); % 3 'raai'
    KL(:,4) = data(:,3); % 4 'BKL'
    KL(:,5) = data(:,4); % 5 'TKLyear'
    KL(:,6) = data(:,5); % 6 'Verschil'
    KL(:,7) = data(:,6); % 7 'Trend'   
    KL(:,1) = year;      % 1 'year'
    KL(isnan(KL(:,2)),:)=[];
    KLdat = [KLdat;KL];
    KL    = [];
end

%% some data preprocessing

year       = KLdat(:,1);
section    = KLdat(:,2);
transectNr = KLdat(:,3);
trID       = KLdat(:,2)*1e6+KLdat(:,3);
BKL        = KLdat(:,4);
TKL        = KLdat(:,5);
TND        = KLdat(:,7);

transectNr(section==13&transectNr==1528)=1525;
transectNr(section==17&transectNr==1242)=1241;

disp('converting coordinates');
    [xBKL,yBKL] =             xRSP2xyRD(BKL,section,transectNr); % jarkus_rsp2xy
[LonBKL,LatBKL] = convertCoordinates(xBKL,yBKL,STD,'CS1.name','Amersfoort / RD New','CS2.name','WGS 84','CS2.type','geo');
    [xTKL,yTKL] =             xRSP2xyRD(TKL,section,transectNr);
[LonTKL,LatTKL] = convertCoordinates(xTKL,yTKL,STD,'CS1.name','Amersfoort / RD New','CS2.name','WGS 84','CS2.type','geo');
    [xTND,yTND] =             xRSP2xyRD(TND,section,transectNr);
[LonTND,LatTND] = convertCoordinates(xTND,yTND,STD,'CS1.name','Amersfoort / RD New','CS2.name','WGS 84','CS2.type','geo');

%% put the data in structure D
disp('pre allocating  data in structure');
%pre allocating
D.year = unique(year);
D.trID = unique(trID);
D.BKL = nan(length(D.trID),length(D.year));

   D.BKL = nan(length(D.trID),length(D.year));
  D.xBKL = nan(length(D.trID),length(D.year));
  D.yBKL = nan(length(D.trID),length(D.year));
D.LatBKL = nan(length(D.trID),length(D.year));
D.LonBKL = nan(length(D.trID),length(D.year));

   D.TKL = nan(length(D.trID),length(D.year));
  D.xTKL = nan(length(D.trID),length(D.year));
  D.yTKL = nan(length(D.trID),length(D.year));
D.LatTKL = nan(length(D.trID),length(D.year));
D.LonTKL = nan(length(D.trID),length(D.year));

   D.TND = nan(length(D.trID),length(D.year));
  D.xTND = nan(length(D.trID),length(D.year));
  D.yTND = nan(length(D.trID),length(D.year));
D.LatTND = nan(length(D.trID),length(D.year));
D.LonTND = nan(length(D.trID),length(D.year));

% loop through the data
for jj = 1:length(D.year);
    disp(['reshaping ' num2str(D.year(jj))])
    for ii = 1:length(D.trID) 
        ind = find(trID==D.trID(ii)&year==D.year(jj));
        if ~isempty(ind)
               D.BKL(ii,jj) =    BKL(ind);
              D.xBKL(ii,jj) =   xBKL(ind);
              D.yBKL(ii,jj) =   yBKL(ind);
            D.LatBKL(ii,jj) = LatBKL(ind);
            D.LonBKL(ii,jj) = LonBKL(ind);
               
               D.TKL(ii,jj) =    TKL(ind);
              D.xTKL(ii,jj) =   xTKL(ind);
              D.yTKL(ii,jj) =   yTKL(ind);
            D.LatTKL(ii,jj) = LatTKL(ind);
            D.LonTKL(ii,jj) = LonTKL(ind);
                           
               D.TND(ii,jj) =    TND(ind);
              D.xTND(ii,jj) =   xTND(ind);
              D.yTND(ii,jj) =   yTND(ind);
            D.LatTND(ii,jj) = LatTND(ind);
            D.LonTND(ii,jj) = LonTND(ind);
        end
    end
end

save('KustLijnData.mat','D');