clear;

data = xlsread('Hourly_electricity_data_2016.xlsx','H2:N8785');
production = data(:,1:6); % Columns 1-6, all rows

% Duration curves
sorted_production = ; % TODO: sort the data to descending order

% Ramp-up/down

timestep = 1; % 1 hour
ramps = zeros(size(production,1)-1,6);
for i = 1:length(ramps) % Calculating delta P / delta t
    ramps(i,:) = ; % TODO: insert formula for calculating ramps
end   

rbins = 100*floor(min(min(ramps))/100):20:100*ceil(max(max(ramps))/100); % Ramp-up bins
freq = hist(ramps,rbins); % Frequencies

% Plots duration curve

figure(1);
plot(sorted_production,'LineWidth',1.5);
colormap jet(6)
xlabel('Hours (h)');
ylabel('Power (MW)');
title('Power duration curve'); 
legend('Nuclear power','Hydropower','Wind power',...
    'CHP','Separate thermal production','Imports');
axis([0 9000 0 5000]);

% Plots ramp histogram

figure(2);
plot(rbins,freq,'LineWidth',1.5);
xlabel('Ramp-up/down (MW/h)');
ylabel('Frequency (number of 1 h timeslots)');
title('Variability analysis');
legend('Nuclear power','Hydropower','Wind power',...
    'CHP','Separate thermal production','Imports');

% Saves figures as .png files

print(figure(1), '-dpng', '-r400', 'a_duration');
print(figure(2), '-dpng', '-r400', 'a_variability');