path = '/Users/iljagolikov/Documents/ssh/study_Ilja/baka/period/fsum_area_'

pikkused1 = [];
aastad1 = 1993:1:2020;
algus1= [];
lopp1=[];
andmed1 = [];

for aasta = string(1993:1:2020)

ncf = append(path,aasta,'.nc');

data1 = ncread(ncf, 'siconc');
data1 = data1(:);

st1 = find(data1,1,'first');
ed1 = find(data1,1,'last');

if sum(data1)~=0
    
    algus1(end+1) = st1;
    lopp1(end+1) = ed1;
    pikkused1(end+1) = ed1 - st1 + 1 ;
    andmed1 = cat(2,andmed1,data1); 
    
else
    pikkused1(end+1) = NaN;     %kui polnud jaad siis pikkus = 0
    data1(1:243)=0;
    %data(1:48)= NaN;
    %data(end:-1:220)=NaN;
    algus1(end+1)=NaN;         % algus = NaN
    lopp1(end+1)=NaN;           % lopp = NaN
    andmed1 = cat(2,andmed1,data1);
end



k = plot(data1)
hold on
k.LineWidth = 0.8
title('Kogu jää pindala 1993/94-2020/21')
xlabel('Päevade arv alates 1-st oktoobrist') 
ylabel('Jää pindala [km^2]') 
disp(aasta);
end

pikkused1 = pikkused1(:);
aastad1 = aastad1(:);
algus1=algus1(:);
lopp1=lopp1(:);

tabel1 = table(aastad1,algus1,lopp1,pikkused1);
%%
path = '/Users/iljagolikov/Documents/ssh/study_Ilja/baka/period/fsum_icevol_'
pikkused = [];
aastad = 1993:1:2020;
algus= [];

lopp=[];

andmed = [];

for aasta = string(1993:1:2020)

ncf = append(path,aasta,'.nc');

data = ncread(ncf, 'sithick');
data = data(:);

st = find(data,1,'first');
ed = find(data,1,'last');

if sum(data)~=0
    %data(1:st)= NaN;
    %data(end:-1:ed)=NaN; 
    algus(end+1)=st;
    lopp(end+1)=ed;
    pikkused(end+1) = ed - st + 1;
    andmed = cat(2,andmed,data); 
else
    pikkused(end+1) = NaN;     %kui polnud jaad siis pikkus = 0
    data(1:243)=0;           % andmed  = NaN
    algus(end+1)=NaN;         % algus = NaN
    lopp(end+1)=NaN;           % lopp = NaN
    andmed = cat(2,andmed,data);
end



k = plot(data)
hold on
k.LineWidth = 0.8
title('Kogu jää ruumala 1993/94-2020/21')
xlabel('Päevade arv alates 1-st oktoobrist') 
ylabel('Jää ruumala [km^3]') 
disp(aasta);
end


pikkused = pikkused(:);
aastad = aastad(:);
algus=algus(:);
lopp=lopp(:);

tabel = table(aastad,algus,lopp,pikkused);

%%

%max S
data1 = andmed1;
data1(1:48,:) = NaN;
data1(end:-1:220,:) = NaN;

max_Ssum = max(data1);

%--------------------
%max V
data = andmed;
data(1:48,:) = NaN;
data(end:-1:220,:) = NaN;

max_Vsum = max(data);


%% Trendid
f = figure(2);
f.Position = [100 100 1200 500];

%pindala
lm1 = fitlm(aastad,max_Ssum)
subplot(2,1,1)
p1 = plot(lm1)

p1(1).LineWidth = 1.5;
p1(2).LineWidth = 2;
p1(3).LineWidth = 1.5;
p1(4).LineWidth = 1.5;

p1(1).MarkerSize = 8;
p1(3).LineStyle = '--';
p1(4).LineStyle = '--';

grid on
text(0.01, 0.9, '(a)', 'Units', 'normalized', 'FontSize', 24)
set(gca,'xtick',aastad(1:3:end));
set(gca,'ytick',0:3000:18000);
axis([1990 2022 -1000 19000])
xlabel('Aastad','fontweight','bold','FontSize', 20) 
ylabel('Jää pindala [km2]','fontweight','bold','FontSize', 20)

lgd = legend()
lgd.Location = 'northeast';
lgd.FontSize = 14;

ax = gca;
ax.XAxis.FontSize = 18;
ax.YAxis.FontSize = 18;
ax.LineWidth = 1;


%ruumala
sum_v = nansum(andmed)
lm2 = fitlm(aastad,max_Vsum)
subplot(2,1,2)
p2 = plot(lm2)

p2(1).LineWidth = 1.5;
p2(2).LineWidth = 2;
p2(3).LineWidth = 1.5;
p2(4).LineWidth = 1.5;

p2(1).MarkerSize = 8;
p2(3).LineStyle = '--';
p2(4).LineStyle = '--';


set(gca,'xtick',aastad(1:3:end));
set(gca,'ytick',0:1:7);
xlabel('Aastad','fontweight','bold','FontSize', 20) 
ylabel('Jää ruumala [km3]','fontweight','bold','FontSize', 20)
grid on
text(0.01, 0.9, '(b)', 'Units', 'normalized', 'FontSize', 24)
axis([1990 2022 -1 8])

lgd = legend()
lgd.Location = 'northeast';
lgd.FontSize = 14;

ax = gca;
ax.XAxis.FontSize = 18;
ax.YAxis.FontSize = 18;
ax.LineWidth = 1;



