function [idx]=dynamickmeans(centknowed,var1,var2,var3,ik)

% dynamickmeans  routine to apply a dynamic k mean clustering
% [idx]=dynamickmeans(centknowed,var1,var2,var3,ik)
%
% input:
% centknowed  = parametre 0 or 1 to identify is the centroids are known or unknow 
%             
% var1 =      must be a 3D matrix containing a variable in time, lat and long
%             dimnension. To get the results of Oderiz et al.(2022) should be the wave
%             power.
% var2 =      must be a 3D matrix containing a variable in time, lat and long
%             dimnension. To get the results of Oderiz et al.(2022) should be
%             the cosine of the mean wave direction
% var3 =      must be a 3D matrix containing a variable in time, lat and long
%             dimnension. To get the results of Oderiz et al.(2022) should be
%             the sine of the mean wave direction.
% ik =        the optimal number of clusters, can be obtained by applying elbow method.
%           
% output:
% idx     = 3D dimenstion of the cluster assigned to each parameter 
%   --------------------------------------------------------------------
%% Dynamci K means 
% This method applies the k means algorithm to variables in 3D (x,y,t).
% For more detail of the method read
% Odériz, I., Silva, R., Mortlock, T. R., Mori, N., Shimura, T., et al. (2021). 
% Natural variability and warming signals in global ocean wave climates. Geophysical
% Research Letters, 48, e2021GL093622. https://doi.org/10.1029/2021GL093622

%   --------------------------------------------------------------------
%% Copyright notice
% Cite following references when you use this function.
% Odériz, I., Silva, R., Mortlock, T. R., Mori, N., Shimura, T., et al. (2021). 
% Natural variability and warming signals in global ocean wave climates. Geophysical
% Research Letters, 48, e2021GL093622. https://doi.org/10.1029/2021GL093622

% Odériz, I.,  Mori, N., Shimura, T., Webb, A., Silva, R., Mortlock, T. R., (2022). 
% Transitional Wave Climate Regions on Continental and Polar Coasts in a Warming World. Nature
% Climate Change

% Odériz, I., Mori, N., Shimura, T., Webb, A., Silva, R., Mortlock, T. R., Dataset for Transitional
% Wave Climate Regions on Continental and Polar Coasts in a Warming World (Version 1). Zenodo. DOI: 10.5281/zenodo.6482839 (2022).
% --------------------------------------------------------------------
%       Itxaso Odériz
%       itxaso.oderiz@gmail.com	
%   --------------------------------------------------------------------

%% size of the 3D matrix
[nr nc nt]=size(var1);

%% 3D to 1D
var1=var1(:);
var2=var2(:);
var3=var3(:);

%% normalized var1 (var-min)/(max-min)
        Tmin=min(min(var1));
        Tmax=max(max(var1));
        VAR1=(var1-Tmin)/(Tmax-Tmin);

%% normalized var2
        dcmin=min(min(var2));
        dcmax=max(max(var2));
        VAR2=(var2-dcmin)/(dcmax-dcmin);
        
%% normalized var3
        dsmin=min(min(var3));
        dsmax=max(max(var3));
        VAR3=(var3-dsmin)/(dsmax-dsmin);
        
%% matrix to 
        VAR=[VAR1,VAR2,VAR3];
        
%% elijo el comando si se los centroides 
if centknowed==1
        %% the centroids are known
        fileID = fopen('centroids.txt','rt');
        C1= textscan(fileID,'%14f %14f %14f \n');
        fclose(fileID);
        CC1=cell2mat(C1);  
       
        %%normalized centroids from 0 to 1
        CC1(:,1)=(C1(:,1)-Tmin)/(Tmax-Tmin);
        CC1(:,2)=(C1(:,2)-dcmin)/(dcmax-dcmin);
        CC1(:,3)=(C1(:,3)-dsmin)/(dsmax-dsmin);
        
        %%Calulo los clusters con los centroides que ya se tenían 
        [idx,C2] = kmeans(VAR,ik,'Start',CC1); 
        
elseif centknowed==0
        %% the centroids are computed
       [idx,C,sumd, D] = kmeans(VAR,ik); 

        %%des-normalized centroids 
        CC(:,1)=C(:,1)*(Tmax-Tmin)+Tmin;
        CC(:,2)=C(:,2)*(dcmax-dcmin)+dcmin;
        CC(:,3)=C(:,3)*(dsmax-dsmin)+dsmin;

        %%save centroids
        save('centroids.txt','CC','-ascii');
        save(centroidsNorm.txt','C','-ascii');
end

%% from 1D to 3D
        idx=reshape(idx,[nr nc nt]);
end