clear;

data = xlsread('Hourly_electricity_data_2016.xlsx', 'H2:N8785');
windPower   = data(:, 3);
consumption = data(:, 7);

% Consumption minus wind profiles __________________________________________________________________________________________________________

windScaleFactor = [];                                               % TODO: add correct factors for different cases here

for j = 1:3 
    consumptionMinusWind(:, j) = ;                                  % TODO: add formula for consumption minus wind here
end

% Duration curves __________________________________________________________________________________________________________________________

sortedConsumptionMinusWind = zeros(length(consumption), 3);

for j = 1:3
    sortedConsumptionMinusWind(:,j) = ;                                    % TODO: sort the data to descending order
end

figure(1);
plot(sortedConsumptionMinusWind,'LineWidth',1.5);
colormap jet
xlabel('Hours (h)');
ylabel('Power (MW)');
title('Duration curves'); 
legend('Consumption - 1 x wind', 'Consumption - 3 x wind', 'Consumption - 10 x wind');
axis([0 9000 -4500 15500]);
grid on;
print(figure(1), '-dpng', '-r400', '2_duration');

% Variability histogram ____________________________________________________________________________________________________________________

ramps = zeros(length(consumption) - 1, 3);

for j = 1:3
    
    timestep = 1; % 1 hour
    for i = 1 : size(ramps, 1)
        ramps(i, j) = ;                                             % TODO: insert formula for calculating ramps
    end
    
    minBin      = 100 * floor(min(min(ramps)) / 100);
    maxBin      = 100 * ceil (max(max(ramps)) / 100);
    rampUpBins  = minBin : 20 : maxBin;
    frequencies = hist(ramps, rampUpBins);
    
end

figure(2);
plot(rampUpBins, frequencies, 'LineWidth', 1.5);
xlabel('Ramp-up/down (MW/h)');
ylabel('Frequency (number of 1 h timeslots)');
title('Variability histograms');
legend('Consumption - 1 x wind', 'Consumption - 3 x wind', 'Consumption - 10 x wind');
grid on;
print(figure(2), '-dpng', '-r400', '2_variability');

% Transition diagram _______________________________________________________________________________________________________________________

for j = 1:3
    
    initialPowerLevel  = ;                                          % TODO: insert variable according to problem description

    timestep = 1; % 1 hour
    ramps    = zeros(length(initialPowerLevel) - 1, 1);
    
    for i = 1 : length(ramps)
        ramps(i) = ;                                                % TODO: insert formula for calculating ramps
    end

    figure(2 + j);
    scatter(initialPowerLevel(1 : end-1), ramps, '.k');
    xlabel('Initial power level (MW)');
    ylabel('Ramp-up/down (MW/h)');
    title(['Transition diagram: Power consumption - ', num2str(windScaleFactor(j)), ' x Wind']);
    grid on;
    print(figure(2 + j), '-dpng', '-r400', ['2_transition_x', num2str(windScaleFactor(j))]);
    
end