clear all
close all

%%
load rustdata.mat

A = double(rustdata);
%A = log(A+1);
% for i = 2 : 23
%     A(:,:,i) = A(:,:,i) - A(:,:,1);
% end


B = reshape(A,308*308,23);
B = B(:,1:23);

%% Transform B to standard normal
for i = 1 : size(B,1)
    %B(:,i) = ( B(:,i) ) ./ std( B(:,i) );
    B(i,:) = ( B(i,:) - mean(B(i,1:5)) ) / std(B(i,1:5));
end
%%
Signal = B(:,1:end);
T = size(Signal, 2); 
d = size(Signal, 1);
w = 200;
Pre_mean = zeros(d,1); % the mean vector before the change-point
Gamma = 10000;
%% Choose the threshold
% b_CUSUM = Threshold_normal_mean('CUSUM', T, d, Pre_mean, Gamma, w);
% b_GLR = Threshold_normal_mean('GLR', T, d, Pre_mean, Gamma, w);
% % b_SRRS = Threshold_normal_mean('SRRS', T, d, Pre_mean, Gamma, w);
% % b_ACM = Threshold_normal_mean('ACM', T, d, Pre_mean, Gamma, w); 
% % b_WLCUSUM = Threshold_normal_mean('WLCUSUM', T, d, Pre_mean, Gamma, w);
% b_WLGLR = Threshold_normal_mean('WLGLR', T, d, Pre_mean, Gamma, w);
% b_WLSRRS = Threshold_normal_mean('WLSRRS', T, d, Pre_mean, Gamma, w);
% b_WLACM = Threshold_normal_mean('WLACM', T, d, Pre_mean, Gamma, w); 
%b_WLSRRS_L1 = Threshold_normal_mean('WLSRRS_L1', 100, d, Pre_mean, Gamma, w);
%b_WLACM_L1 = Threshold_normal_mean('WLACM_L1', 100, d, Pre_mean, Gamma, w); 

%%
b_CUSUM = 0; 
b_WLACM = 0;
b_WLGLR = 48365;
b_WLSRRS = 0;
b_WLSRRS_L1 = 0;
b_WLACM_L1 = 0;



%% Compute the Detection delay

EDD_CUSUM = EDD_normal_mean_rust('CUSUM', T, d, Pre_mean, Signal, b_CUSUM, w);
disp(mean(EDD_CUSUM));
disp(std(EDD_CUSUM));
% EDD_GLR = EDD_normal_mean('GLR', Len, d, Post_mean, b_GLR, w);
% EDD_SRRS = EDD_normal_mean('SRRS', Len, d, Post_mean, b_SRRS, w);
% EDD_ACM = EDD_normal_mean('ACM', Len, d, Post_mean, b_ACM, w);
% EDD_WLCUSUM = EDD_normal_mean('WLCUSUM', Len, d, Post_mean, b_WLCUSUM, w);
EDD_WLGLR = EDD_normal_mean_rust('WLGLR', T, d, Pre_mean, Signal, b_WLGLR, w);
disp(mean(EDD_WLGLR));
disp(std(EDD_WLGLR));
EDD_WLSRRS = EDD_normal_mean_rust('WLSRRS', T, d, Pre_mean, Signal, b_WLSRRS, w);
disp(mean(EDD_WLSRRS));
disp(std(EDD_WLSRRS));
EDD_WLACM = EDD_normal_mean_rust('WLACM', T, d, Pre_mean, Signal, b_WLACM, w);
disp(mean(EDD_WLACM));
disp(std(EDD_WLACM));
EDD_WLSRRS_L1 = EDD_normal_mean_rust('WLSRRS_L1', T, d, Pre_mean, Signal, b_WLSRRS_L1, w);
disp(mean(EDD_WLSRRS_L1))
disp(std(EDD_WLSRRS_L1));
EDD_WLACM_L1 = EDD_normal_mean_rust('WLACM_L1', T, d, Pre_mean, Signal, b_WLACM_L1, w);
disp(mean(EDD_WLACM_L1));
disp(std(EDD_WLACM_L1));
%%
%figure;
%imagesc(A(:,:,1));

