function varargout = ws_shapefile2kml(varargin)
%WS_SHAPEFILE2KML  from shape file to Google Earth file
% for the "ingrepen" of the Western Scheldt
% 
%   ws_shapefile2kml(<keyword,value>)   
%
% For the <keyword,value> pairs and their defaults call
%
%    OPT = ws_shapefile2kml()
%
% where the following keywords are required: 
% 'filename', 'filename_dir', 'epsg', 'color' 
%
% Example_1: GE plot 2D of baggeren, 'Amersfoort / RD New' in blue
%  
%    ws_shapefile2kml('filename','baggeren\',...
%                 'filename_dir','bagger',...
%                         'epsg', 28992,...
%                        'color',[0 0 1]); 
%
% Example_2: GE plot 3D of storten, 'Amersfoort / RD New' in red
%  
%    ws_shapefile2kml('filename','storten\',...
%                 'filename_dir','stort',...
%                         'epsg', 28992,...
%                        'color',[1 0 0]); 
%
% Example_3: GE plot of macrocells "rvaksas3_vaknr" (2D), 'Amersfoort / RD New' in white
%  
%    ws_shapefile2kml('filename','rvaksas3\',...
%                 'filename_dir','rvaksas3_vaknr',...
%                         'epsg', 28992,...
%                        'color',[1 1 1]); 
%
%See also: CONVERTCOORDINATES, ARC_SHAPE2KML, ARC_SHAPE_READ 


%% settings

OPT.directory_out     = 'P:\1202022-ltv-zandhuishouding\data\shapefiles';
OPT.filename          = ''; % 'storten\', 'baggeren\', 'zandwinning_rijk\', 'zandwinning_handel\', 'rvaksas3\'
OPT.filename_dir      = ''; % 'stort', 'bagger', 'win_derden', 'winc', 'rvaksas3', 'rvaksas3_vaknr'
OPT.shpfile           = '*.shp';
OPT.epsg              = ''; % 28992
OPT.color             = []; % [0 0.5 0]: win_derden; [0 1 0]: winc; [0 0 1]: bagger; [1 0 0]: stort 
OPT.cleanup           = 0;
OPT.pause             = 0;
OPT.plot              = 1;
OPT.fldname           = ''; % *.dbf fieldname entry to be used for name (default sequence number)
OPT.deleteSourceFiles = 1;
OPT.extrude           = 0; % 2D / 3D polygons
OPT.zScaleFun         = @(z) z.*100;

if nargin==0
  varargout = {OPT};
  return
end

OPT = setproperty(OPT,varargin{:});

%% for rvaksas3 and rvaksas3_vaknr only
if strcmp(OPT.filename_dir,'rvaksas3')
    arc_shape2kml('rvaksas3\rvaksas3.shp',...
                  'epsg',OPT.epsg,...
                 'color',OPT.color) % default: 'color',[1 1 0.5]
    return
elseif strcmp(OPT.filename_dir,'rvaksas3_vaknr')
    arc_shape2kml('rvaksas3\rvaksas3_vaknr.shp',...
                  'epsg',OPT.epsg,...
                 'color',OPT.color) % default: 'color',[1 1 1]
    return
end

%% list of the .shp files
fss = dir(fullfile(OPT.directory_out, OPT.filename, OPT.shpfile));

%% define volume names
volNames  = {'HVH__M3_'; 'HVH_M3'; 'KUUBS'; 'BEUNKUUBS'; 'hoev_m3'; };
areaNames = {'AREA'; 'AREAB'; 'AREABB'; 'HECTARES';};

%% cycle

% pre-allocation
date      = NaN(length(fss),1);
%shapename = NaN(length(fss),1);

%% define date
for ii = 1:length(fss)
shpname     = fss(ii).name;
date(ii)    = str2double(cellstr(shpname(end-7:end-4)));
shapename{ii} = shpname;
%% for storten only
if strcmp(fss(ii).name,'stort_2006_jan_jun.shp')
    shpname     = fss(ii).name;
    date(ii)    = str2double(cellstr(shpname(end-15:end-12)));
elseif strcmp(fss(ii).name,'stort_2006_jul_dec.shp')
    shpname     = fss(ii).name;
    date(ii)    = str2double(cellstr(shpname(end-15:end-12))) + 0.5;
    jj          = ii;
end
end
shapename = shapename(~isnan(date));
date      = datenum(num2str(date(~isnan(date))), 'yyyy');

%% for storten only
if strcmp(OPT.filename_dir,'stort') && date(17) == date(18)
    date(18) = datenum(2006, 7, 1);
end

%% extend last time interval for visibility period in Google Earth
date4GE        = date;
date4GE(end+1,1) = date(end)+365; % the duration of visualization is set to the whole reference year

%% maketiles
% cycle for different years
for ii = 1:length(date)
    disp(['date: ', num2str(datestr(date(ii)))])
    % define name
    OPT.fname         = fullfile(OPT.filename, char(shapename(ii)));

%% read shape file
    [R,M]  = arc_shape_read(OPT.fname);
        n  = length(R);
        
    % cycle for the polygons
    for i=1:n
        disp([num2str(i,'%0.3d'),' / ',num2str(n,'%0.3d')])
        x = R(i).X;
        y = R(i).Y;

        %% for zandwinning_rijk\ only
        if strcmp(OPT.filename_dir,'win_derden')
        R(i).AREA = 97167;
        end
        
        %% find z for each polygon
        volName   = isfield(R,volNames);
        areaName  = isfield(R,areaNames);
        if any(volName) && any(areaName)
            if ismember('AREA',fieldnames(R(i))) && R(i).AREA ~= 0
                z = eval(['R(i).',char(intersect(volNames,fieldnames(R)))])...
                    /(R(i).AREA);
            elseif ismember('AREAB',fieldnames(R(i))) &&...
                    R(i).AREA == 0 && R(i).AREAB ~= 0
                z = eval(['R(i).',char(intersect(volNames,fieldnames(R)))])...
                    /(R(i).AREAB);
            elseif ismember('AREABB',fieldnames(R(i))) &&...
                    R(i).AREA == 0 && R(i).AREAB == 0
                z = eval(['R(i).',char(intersect(volNames,fieldnames(R)))])...
                    /(R(i).AREABB);
            elseif ~ismember('AREA',fieldnames(R(i))) &&...
                    ismember('HECTARES',fieldnames(R(i)))
                z = eval(['R(i).',char(intersect(volNames,fieldnames(R)))])...
                    /(R(i).HECTARES*1e4);
            end
        end        
        
        if isempty(x) % 'zuid holland\waterkeringen\Waterkeringen.shp' no 23
            warning(['shape with length 0: ',OPT.fname,' no ',num2str(i)])
        else

        [lon,lat]           = convertCoordinates(R(i).X,R(i).Y,'CS1.code',OPT.epsg,'CS2.code',4326);

        % for patch-holes remove the NaNs
        mask = isnan(R(i).X);
        lonp = lon;
        latp = lat;
        lonp(mask) = lonp(1);
        latp(mask) = latp(1);
   
   %% plot in matlab
   
        if OPT.plot
            h(i) = patch(lonp,latp,[0 1 0]);
            hold on
            p(i) = plot (lon,lat,'k');
         
            axislat(52)
            %axis equal
            hold on
         
            set(h(i),'edgecolor','none')
        end
      
   %% plot as google earth files
      
        if ~exist([OPT.directory_out,filesep,fileparts(OPT.fname)])
        mkpath([OPT.directory_out,filesep,fileparts(OPT.fname)])
        end
   
        fillnames{i} = [OPT.directory_out,filesep,fileparts(OPT.fname),filesep,filename(OPT.fname),num2str(i,'%0.3d'),'_fill.kml'];
        edgenames{i} = [OPT.directory_out,filesep,fileparts(OPT.fname),filesep,filename(OPT.fname),num2str(i,'%0.3d'),'_edge.kml'];
        textnames{i} = [OPT.directory_out,filesep,fileparts(OPT.fname),filesep,filename(OPT.fname),num2str(i,'%0.3d'),'_text.kml'];
      
      %% name
      if ~isempty(OPT.fldname)
         try
          namestr = R(i).(OPT.fldname);
         catch
             disp(char(fieldnames(R)))
             error(['unknown field:',OPT.fldname])
             break
         end
      else
         namestr = num2str(i);
      end

      %% thin back line
      if OPT.extrude
      KMLline(lat,lon,z*ones(1,length(lat)),...
                       'fileName',edgenames{i},...
                        'kmlName',namestr,...
                         'timeIn',date4GE(ii:end-1),...
                        'timeOut',date4GE((ii+1):end),...
                           'fill',false,...
                      'zScaleFun',OPT.zScaleFun);      
      else
      KMLline(lat,lon,'fileName',edgenames{i},...
                        'kmlName',namestr,...
                         'timeIn',date4GE(ii:end-1),...
                        'timeOut',date4GE((ii+1):end));
      end
                    
%                      KMLline(lat_0', lon_0', z_KML', 'timeIn',timeIn,'timeOut',timeOut,...
%     'fileName', fname, 'lineColor',jet(size(z_KML,1)),'lineWidth',2,...
%     'fillColor',jet(size(z_KML,1)),'zScaleFun', @(z_KML)(z_KML+60)*10);

      %% colored fill
      if     strcmpi(R(i).Geometry,'Point')
      elseif strcmpi(R(i).Geometry,'Line') 
      KMLline(lat,lon,'fileName',fillnames{i},...
             'kmlName',namestr,...
              'timeIn',date4GE(ii:end-1),...
             'timeOut',date4GE((ii+1):end),...
           'lineAlpha',0.6,...
           'lineColor',OPT.color,...
           'lineWidth',5);
      elseif strcmpi(R(i).Geometry,'Line')|| ...
             strcmpi(R(i).Geometry,'Polygon')
         if OPT.extrude
      KMLpatch3(latp',lonp',z*ones(length(latp),1),...
            'fileName',fillnames{i},...
             'kmlName',namestr,...
              'timeIn',date4GE(ii:end-1),...
             'timeOut',date4GE((ii+1):end),...
           'lineAlpha',0,...
           'fillColor',OPT.color,...
           'fillAlpha',0.6,...
           'zScaleFun',OPT.zScaleFun);         
         else
      KMLpatch(latp',lonp','fileName',fillnames{i},...
             'kmlName',namestr,...
              'timeIn',date4GE(ii:end-1),...
             'timeOut',date4GE((ii+1):end),...
           'lineAlpha',0,...
           'fillColor',OPT.color,...
           'fillAlpha',0.6);
         end
      end

      fldnames = setxor(fieldnames(R),{'BoundingBox','X','Y','Geometry'});
         
      html = ['<table border="1"><tr>'];
      for j=1:length(fldnames)
         fldname = fldnames{j};
         value   = R(i).(fldname);
         if isnumeric(value)
         value = num2str(value);
         end
         html = [html '<td>' fldname '=</td><td>' value '</td><tr>'];
      end
      html = [html '</tr>'];
      
      KMLmarker(lat(1),lon(1),'fileName',textnames{i},...
                               'kmlName',namestr,...
                                'timeIn',date4GE(ii:end-1),...
                               'timeOut',date4GE((ii+1):end),...               
                                  'name',namestr,...
                                  'html',html,...
                      'colornormalState',OPT.color,...
                   'colorhighlightState',OPT.color);                    
               
        end        
   
    if OPT.pause
    pausedisp
    end
   
    end
    
  %%  mergefiles
  
     edgename = [OPT.directory_out,filesep,fileparts(OPT.fname),filesep,filename(OPT.fname),'_edge.kml'];
   KMLmerge_files('fileName',edgename,'sourceFiles',edgenames,...
                   'kmlName','edge',...
         'deleteSourceFiles',OPT.deleteSourceFiles);

   fillname = [OPT.directory_out,filesep,fileparts(OPT.fname),filesep,filename(OPT.fname),'_fill.kml'];
   KMLmerge_files('fileName',fillname,'sourceFiles',fillnames,...
                   'kmlName','fill',...
         'deleteSourceFiles',OPT.deleteSourceFiles);
   
   textname = [OPT.directory_out,filesep,fileparts(OPT.fname),filesep,filename(OPT.fname),'_text.kml'];
   KMLmerge_files('fileName',textname,'sourceFiles',textnames,...
                   'kmlName','text',...
         'deleteSourceFiles',OPT.deleteSourceFiles);
       
   try % the XML meta-data files are generally mess: they change per arcGIS version and per user.
   if ~isfield(M.dataset_description.identification,'dataset_title')
   M.dataset_description.identification.dataset_title = M.dataset_description.identification.alternative_title;
   end

   if ~isfield(M.dataset_description,'overview')
   M.dataset_description.overview.summary = '';
   else
      if ~isfield(M.dataset_description.overview,'summary')
      M.dataset_description.overview.summary = '';
      end
   end
   catch
      M.dataset_description                              = struct();
      M.dataset_description.identification.dataset_title = '';
      M.dataset_description.overview.summary             = '';
   end

   %if ~isempty(OPT.html)
   %  %M.dataset_description.overview.summary = ['<![CDATA[',str2line(textread(OPT.html,'%s')),']]>'];
   %   M.dataset_description.overview.summary = ['<A href="',filenameext(OPT.html),'" >metadata</A>'];
   %end

   fullname = [OPT.directory_out,filesep,fileparts(OPT.fname),filesep,filename(OPT.fname),'.kml'];
   KMLmerge_files('fileName',fullname,...
                  'sourceFiles',{fillname edgename textname},...
                  'foldernames',{'areas','edges','texts'},...
       'kmlName',M.dataset_description.identification.dataset_title,...
   'description',M.dataset_description.overview.summary,...
   'deleteSourceFiles',OPT.deleteSourceFiles);
   
%% clean-up   
   
   if OPT.cleanup
    for i=1:n
     delete(edgenames{i})   
     delete(fillnames{i})   
     delete(textnames{i})   
    end
    delete(edgename)   
    delete(fillname)   
    delete(textname)   
   end

%% clear variables

    clear edgenames fillnames textnames
    clear edgename  fillname  textname
   
end

for f = 1:length(shapename) % polyname to be defined first
polyname{f} = [OPT.directory_out,filesep,fileparts(OPT.fname),filesep,filename(char(shapename(f))),'.kml'];
end

KMLmerge_files('fileName',[OPT.filename_dir, '.kml'],...
            'sourceFiles', polyname,...% {fillname edgename textname},...
                'kmlName',OPT.filename_dir,...
            'description',fileparts(OPT.filename),...
      'deleteSourceFiles',0);

%% EOF
