%%
%
% Freek Scheel
% free.scheel@deltares.nl
% OpenEarthTools
%

clearvars; close all; clc;

%%
%

tic;

disp('Initializing...');

current_file = mfilename('fullpath');

addpath(genpath(current_file(1:max(strfind(current_file,filesep)))));

data_file = [current_file(1:max(strfind(current_file,filesep))) 'raw' filesep 'korrel.txt'];

plot_figures = 0; % only needed if not present yet

%%

[D,U]  = sedimentatlas_korrel_read(data_file);

hold off; close all;

no_of_points_to_comp = size(D.cumPHI,2);
% no_of_points_to_comp = 200;
disp('Computing mud contents...');
for ii=1:no_of_points_to_comp
    perc=10;
    kml_print_data.d10(ii,1) = max([D.D(end,1) (D.D(find(D.cumPHI(:,ii)>perc,1,'last')))-((D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last'),ii)-perc)/((D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last'),ii)-D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last')+1,ii))))*(D.D(find(D.cumPHI(:,ii)>perc,1,'last'))-D.D(find(D.cumPHI(:,ii)>perc,1,'last')+1))]);
    perc=50;
    kml_print_data.d50(ii,1) = min([D.D(1,1) (D.D(find(D.cumPHI(:,ii)>perc,1,'last')))-((D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last'),ii)-perc)/((D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last'),ii)-D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last')+1,ii))))*(D.D(find(D.cumPHI(:,ii)>perc,1,'last'))-D.D(find(D.cumPHI(:,ii)>perc,1,'last')+1))]);
    perc=90;
    kml_print_data.d90(ii,1) = min([D.D(1,1) (D.D(find(D.cumPHI(:,ii)>perc,1,'last')))-((D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last'),ii)-perc)/((D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last'),ii)-D.cumPHI(find(D.cumPHI(:,ii)>perc,1,'last')+1,ii))))*(D.D(find(D.cumPHI(:,ii)>perc,1,'last'))-D.D(find(D.cumPHI(:,ii)>perc,1,'last')+1))]);
    
    if plot_figures == 1;
        plot(D.D,D.cumPHI(:,ii),'k','linewidth',4);
        xlabel('Grain size'); ylabel('Cumulative percentage'); xlim([0.1 4000]); ylim([0 100]); title('Cumulative distribution of sediment diameter'); set(gca,'XScale','log'); axis on; box on; grid on; hold off;
        print(gcf,'-dpng','-painters',[current_file(1:max(strfind(current_file,filesep))) 'generate_kmz' filesep 'kml_temp_data' filesep 'Figure_' num2str(ii,'%05.f') '.png'],'-r70');
        disp(['Plotting figure ' num2str(ii) '/' num2str(no_of_points_to_comp)]);
    end
end
close all;

% j = 1;
% sedimentatlas_plot_diameter(D.diameter,D.phi(j,:),D.cumphi(j,:),'t.png')
%% Calculate derived parameters
% * Rijkswaterstaat, 2004. Zwarts, L. (2004). Bodemgesteldheid en mechanische kokkelvisserij in de
%   Waddenzee. Rijkswaterstaat Waterdienst RIZA report 2004.028. Available
%   from: http://www.rijkswaterstaat.nl (in Dutch).

   D.pmud16_from_16   = 0.2988.*sum(D.PHI(17:end,:),1)'; % mud if D < 16 um
   D.pmud16_from_63   = 0.1879.*sum(D.PHI(13:end,:),1)'; % mud if D < 63 um
   D.total_smaller_63 = sum(D.PHI(13:end,:),1)'; % total <63 um
   
%% Plot   
%    clf
%    caxis            (log10([.1 50]))
%    scatter          (D.lon,D.lat,20,log10(D.pmud16_from_16),'filled')
%    [ax,h]=colorbarwithtitle('percentage mud in bed');
%    set              (ax,'YTick',log10([.1 .2 .5 1 2 5 10 20 50]))
%    set(ax,'YTickLabel',{'.1','.2','.5','1','2','5','10','20','50','100'})   
%    axislat          (52)
%    tickmap          ('ll')
%    grid on
%    print2screensize ('sediment_atlas_waddenzee_EVAII_mud_content')
   
%% Plot in Google Earth

[D.lon,D.lat]=convertCoordinates(D.X_COORD,D.Y_COORD,'CS1.code',28992,'CS2.code',4326);

   OPT.fileName           = [];
   OPT.description        = '';
   OPT.kmlName            = 'Mud fraction in Waddensea';
   OPT.openInGE           = false;
   OPT.text               = '';
   OPT.name               = [];
   OPT.lon                = [];
   OPT.lat                = [];
   OPT.z                  = [];
   OPT.open               = 0;
   OPT.debug              = 0;
   
   OPT.overlayXY          = [0 0.00];
   OPT.screenXY           = [0 0.04];

   OPT.iconnormalState         = 'circle-white.png';
   OPT.iconhighlightState      = 'circle-white.png';
   OPT.iconnormalStateScale    = 0.6;
   OPT.iconhighlightStateScale = 1.5;
   
   OPT.varname            = []; % statistics
   OPT.preview            = 0;  % add 2 previews of - and hence requires - varname
   OPT.credit             = 'Plot created with: www.OpenEarth.eu';
   OPT.resolveName        = '';


   if OPT.debug
      eol           = char(10);
   else
      eol           = '';
   end
   
if isempty(OPT.lon);OPT.lon = mean([max(D.lon) min(D.lon)]); end
if isempty(OPT.lat);OPT.lat = mean([max(D.lat) min(D.lat)]); end
if isempty(OPT.z  );OPT.z   = 170e3; end

OPT_header = struct(...
      'kmlName',OPT.kmlName,...
  'description',OPT.description,...
         'open',OPT.open,...
    'cameralon',OPT.lon,...
    'cameralat',OPT.lat,...
      'cameraz',OPT.z);
   output = KML_header(OPT_header);

colormap_for_bar = [0.05 0.1 0.2 0.5 1 2 5 10 20 50]';
   disp('Generating color bar...');
KMLcolorbar('CBfileName',[pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'step_1_sediment_atlas_waddenzee_mud_content_colorbar.kml'],...
     'CBcolorMap'      ,@(m) jet(m),...
     'CBcolorSteps'    ,length(jet),...
     'CBcLim'          ,log10(colormap_for_bar),...
     'CBcolorTick'     ,log10(colormap_for_bar),...
     'CBcolorTickLabel',{['  0'; num2str(colormap_for_bar(2:end,1))]},...
     'CBcolorTitle'    ,'Mud fraction [%] (based on < 63 \mum)',...
     'CBtemplateHor'   ,[],'CBtemplateVer',[],...
     'CBthickness'     ,.1,'CBalignmargin',0.02,...
     'CBalpha'         ,1,...
     'CBfontrgb'       ,[1 1 1],'CBbgcolor',[0 0 0],'CBtitlergb',[0 0 0],'CBframergb',[1 1 1],...
     'CBcolorMap',jet(size(colormap_for_bar,1))); % white letters with black halo
   close all;
   colormap_for_bar(1,1) = 0;

for ii=1:size(colormap_for_bar,1);
    colormap_jet_num = jet(size(colormap_for_bar,1));
    temp = dec2hex(round([ones(1,1) colormap_jet_num(ii,:)].*255));
    dot_colors_hex(ii,:) = [temp(1,:) temp(4,:) temp(3,:) temp(2,:)];
end

% ORIGINAL:

% output = [output ...
%     '<Style id="normalState"   >'...
%     '  <IconStyle><scale>' num2str(OPT.iconnormalStateScale) '</scale><color>' dot_colors_hex(1,:) '</color><Icon><href>'...
%     '  ' OPT.iconnormalState ''...
%     '  </href></Icon></IconStyle>'...
%     '  <LabelStyle><scale>0</scale></LabelStyle>'...
%     '  </Style>'...
%     '<Style id="highlightState">'...
%     '  <IconStyle><scale>' num2str(OPT.iconhighlightStateScale) '</scale><Icon><href>'...
%     '  ' OPT.iconhighlightState ''...
%     '  </href></Icon></IconStyle>'...
%     '  <BalloonStyle>'...
%     '  <text><![CDATA[<h3>$[name]</h3>'...
%     '  $[description]<br><br><br><center>Powered by:<br></center>'...
%     '  <center><img src="OpenEarth-logo-blurred-white-background.png" width="80"/><br></center>]]>'...
%     '  <center><font size="1">And a little bit of Freek</font></center></text>'...
%     '  </BalloonStyle></Style>'...
%     '<StyleMap id="MarkerBalloonStyle">'...
%     '  <Pair><key>normal</key><styleUrl>#normalState</styleUrl></Pair> '...
%     '  <Pair><key>highlight</key><styleUrl>#highlightState</styleUrl></Pair> '...
%     '  </StyleMap>'];
disp('Generating kml styles...');
for ii=1:size(dot_colors_hex,1)

    output = [output ...
        ['<Style id="normalState' num2str(ii,'%02.0f') '"   >']...
        '  <IconStyle><scale>' num2str(OPT.iconnormalStateScale) '</scale><color>' dot_colors_hex(ii,:) '</color><Icon><href>'...
        '  ' OPT.iconnormalState ''...
        '  </href></Icon></IconStyle>'...
        '  <LabelStyle><scale>0</scale></LabelStyle>'...
        '  </Style>'...
        '<Style id="highlightState' num2str(ii,'%02.0f') '">'...
        '  <IconStyle><scale>' num2str(OPT.iconhighlightStateScale) '</scale><color>' dot_colors_hex(ii,:) '</color><Icon><href>'...
        '  ' OPT.iconhighlightState ''...
        '  </href></Icon></IconStyle>'...
        '  <BalloonStyle>'...
        '  <text><![CDATA[<h3>$[name]</h3>'...
        '  $[description]<br><br><br><center>Powered by:<br></center>'...
        '  <center><img src="OpenEarth-logo-blurred-white-background.png" width="80"/><br></center>]]>'...
        '  <center><font size="1">And a little bit of Freek</font></center></text>'...
        '  </BalloonStyle></Style>'...
        ['<StyleMap id="MarkerBalloonStyle' num2str(ii,'%02.0f') '">']...
        ['  <Pair><key>normal</key><styleUrl>#normalState' num2str(ii,'%02.0f') '</styleUrl></Pair> ']...
        ['  <Pair><key>highlight</key><styleUrl>#highlightState' num2str(ii,'%02.0f') '</styleUrl></Pair> ']...
        '  </StyleMap>'];

end

%% labels

   output = [output, sprintf('\n%s','<Folder>')];
   output = [output, sprintf('\n%s','<open>',num2str(OPT.open),'</open>')];
   output = [output, sprintf('\n%s','<name>Platforms</name>')];

   color_indices = nearestpoint(D.pmud16_from_63,colormap_for_bar);
   
disp(['Adding ' num2str(no_of_points_to_comp) ' locations to kml (might take some time)...']);
   for ll=1:no_of_points_to_comp
   
locstring = [];
locstring = [locstring, sprintf([...
     '<tr><td colspan="2" bgcolor="#999999"><div style="color:#FFFFFF;">Where:</div></td></tr>' eol ...
     '<tr><td bgcolor="#FFFFFF">Lon.:   </td><td bgcolor="#FFFFFF">%s</td></tr>' eol ...
     '<tr><td bgcolor="#FFFFFF">Lat.:   </td><td bgcolor="#FFFFFF">%s</td></tr>' eol ...
     ],...
     [num2str(D.lon(1)) ' deg.'],...
     [num2str(D.lat(1)) ' deg.'])];

whatstring = [];
whatstring = [whatstring sprintf([...
     '<tr><td colspan="2" bgcolor="#999999"><div style="color:#FFFFFF;">Characteristics:</div></td></tr>' eol ...
     '<tr><td bgcolor="#FFFFFF">d50:   </td><td bgcolor="#FFFFFF">%s</td></tr>' eol ...
     '<tr><td bgcolor="#FFFFFF">d90:   </td><td bgcolor="#FFFFFF">%s</td></tr>' eol ...
     '<tr><td bgcolor="#FFFFFF">Pm:  </td><td bgcolor="#FFFFFF">%s</td></tr>' eol ...
     '<tr><td bgcolor="#FFFFFF">Pm:  </td><td bgcolor="#FFFFFF">%s</td></tr>' eol ...
     ],...
     [num2str(round(kml_print_data.d50(ll,1))) ' &#181;m'],...
     [num2str(round(kml_print_data.d90(ll,1))) ' &#181;m'],...
     [num2str(D.pmud16_from_16(ll,1),'%9.1f') '&#1642; (mud content based on < 16 &#181;m)'],...
     [num2str(D.pmud16_from_63(ll,1),'%9.1f') '&#1642; (mud content based on < 63 &#181;m)'])];  

   if ~isempty(OPT.text)
       tableheader = ['<tr><td colspan="2" bgcolor="#666666"><div style="color:#FFFFFF;">',char(OPT.text),'</div></td></tr>'];
    else
       tableheader = '';        
   end
   
   png_file = ['Figure_' num2str(ll,'%05.0f') '.png'];
    
   
   previewhtml  = ['<img src="',png_file,'" alt="Preview">'];
   
   
%    previewhtml  = ['<img src="',png_file,'" alt="Preview">'...
%                           '<img src="',png_file,'" alt="Preview"><hr>'];
   
   
   table = [tableheader...
        '<center><table bgcolor="#333333" cellpadding="3" cellspacing="1"><tbody>'...
        sprintf([locstring,...
        whatstring])...
        '</table></center><hr>' previewhtml];
   
   
   
   
   output = [output, sprintf([...
        '\n<Placemark>'...
        '<snippet></snippet>'...
        '<description><![CDATA[<hr>'...
        '%s'...                                                        % 3 table with links, doubel white spaces disappear
        ']]></description>'...
        '<styleUrl>#MarkerBalloonStyle</styleUrl>'...
        '<Point><coordinates>%3.8f,%3.8f,0</coordinates></Point>'... % 4 lon lat
        '<styleUrl>#MarkerBalloonStyle%02.0f</styleUrl></Placemark>'],...
         table,...
         [D.lon(ll) D.lat(ll) color_indices(ll,1)])];
     
     
     
     
   end
   
   disp('Writing kml file...');
   
   fid     = fopen([pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'step_1_sed_distr_and_mud_content_in_waddensea.kml'],'w');
   fprintf(fid,'%c',output); % handle any % in output, cannot handle \n any more
   fprintf(fid,'\n%s\n','</Folder>'); % handle any % in output, cannot handle \n any more
   output  = [];
   
   fprintf(fid,KML_footer,'%s');
   
   fclose(fid);
  
   disp('Done! KML generated');
   
   
disp('Combining kml, figures & colorbar into one kmz file (might take some time)...');
figure_names   = dir([pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'Figure_*.png']);
kml_names      = dir([pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'step_1_*.kml']);
KMLmerge_files('sourceFiles',{kml_names.name}','fileName',[pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'sed_distr_and_mud_content_in_waddensea.kml'],'cameralon',OPT.lon,'cameralat',OPT.lat,'cameraz',OPT.z);
KML2kmz([pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'sed_distr_and_mud_content_in_waddensea.kml'],[pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'OpenEarth-logo-blurred-white-background.png'],[pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'step_1_sediment_atlas_waddenzee_mud_content_colorbar_ver_lft.png'],OPT.iconnormalState,OPT.iconhighlightState,{figure_names.name}'); % oet logo, OPT.iconnormalState,OPT.iconhighlightState, are cached because they're googles
movefile([pwd filesep 'generate_kmz' filesep 'kml_temp_data' filesep 'sed_distr_and_mud_content_in_waddensea.kmz'],[pwd filesep 'generate_kmz' filesep 'sed_distr_and_mud_content_in_waddensea_incl_png.kmz'])
copyfile([pwd filesep 'generate_kmz' filesep 'sed_distr_and_mud_content_in_waddensea_incl_png.kmz'],[pwd filesep 'sed_distr_and_mud_content_in_waddensea_incl_png.kmz'])

disp(' ');
disp(['Script was finished in ' num2str(toc,'%5.0f') ' seconds']);
   

   
   
   %% End Freek
   


description =   ['<hr> Sediment Atlas WaddenZee: size distribution curves from 7502 sediment surface grabs.',...
                 'Mud content is derived from sediment distribution curves usingh EVA II formulations.',...
                 'This is a proof-of-concept demo of how scatter point could'...
                 'be presented in Google Earth for easy navigation in time space. '...
                 'The data in this proof-of-concept demo is a cache that is updated every now and then.'...
                 'For data and meta-data please visit the original source: <a href="http://www.waddenzee.nl/Sedimentatlas.729.0.html">www.waddenzee.nl</a>.'....
                 'or ask the publisher Rojkswaterstaat <a href="http://www.helpdeskwater.nl">www.helpdeskwater.nl</a>.'...
                 '<hr><table bgcolor="#333333" cellpadding="3" cellspacing="1"><tbody><tr><td colspan="2" bgcolor="#666666"><div style="color:#FFFFFF;">Credits:</div></td></tr>',...
                 '<tr><td    bgcolor="#FFFFFF">data source     </td><td bgcolor="#FFFFFF">Rijkswaterstaat</td></tr>',...
                 '<tr><td    bgcolor="#FFFFFF">data source url </td><td bgcolor="#FFFFFF">http://www.waddenzee.nl</td></tr>',...
                 '<tr><td    bgcolor="#FFFFFF">data provider   </td><td bgcolor="#FFFFFF">OpenEarth</td>',...
                 '<tr><td    bgcolor="#FFFFFF">data distributor</td><td bgcolor="#FFFFFF">http://www.OpenEarth.eu</td>',...
                 '</tr></tbody></table><hr>',...
                 '<p><font size="1" face="courier" color="gray">$HeadURL$ $Id$</font></p>'];


% 
% 
%    KMLscatter(D.lat,D.lon,log10(D.pmud16_from_16),...
%           'name',cellstr((num2str(D.pmud16_from_16'))),...
%       'fileName','sediment_atlas_waddenzee_EVAII_mud_content.kml',...
%    'description',description,...
%        'kmlName','Sediment Atlas WaddenZee: mud content',...
%   'CBcolorTitle','mud content [%]',...
%    'markerAlpha',.8,...
%      'CBcolorMap'      ,@(m) jet(m),...
%      'CBcolorSteps'    ,length(jet),...
%      'CBcLim'          ,log10([.1 50]),...
%      'CBcolorTick'     ,log10([.1 .2 .5 1 2 5 10 20 50]),...
%      'CBcolorTickLabel',{'.1','.2','.5','1','2','5','10','20','50'},...
%      'CBcolorTitle'    ,'mud fraction [%]'); % white letters with black halo)
%     
% %% test for log-scales on colorbar
% 
%   KMLcolorbar('CBfileName','sediment_atlas_waddenzee_EVAII_mud_content_colorbar.kml',...
%      'CBcolorMap'      ,@(m) jet(m),...
%      'CBcolorSteps'    ,length(jet),...
%      'CBcLim'          ,log10([.1 50]),...
%      'CBcolorTick'     ,log10([.1 .2 .5 1 2 5 10 20 50]),...
%      'CBcolorTickLabel',{'.1','.2','.5','1','2','5','10','20','50'},...
%      'CBcolorTitle'    ,'mud fraction [%]',...
%      'CBtemplateHor'   ,[],'CBtemplateVer',[],...
%      'CBthickness'     ,.05,'CBalignmargin',0.02,...
%      'CBfontrgb'       ,[1 1 1],'CBbgcolor',[0 0 0],'CBtitlergb',[1 1 1],'CBframergb',[1 1 1]); % white letters with black halo