%% load the raw data from synapses output file
raw = load('synapse.out');


%% get mean and variance from the 1000 neurons
ns = 999;
l = floor(length(raw)/ns);
raw2 = zeros(l,1);
raw2_var = zeros(l,1);
for i=1:l,
    % hint: adapted to read only the first 100 out of 1000 neurons
    raw2(i,1) =    mean(raw((i*ns-ns+1):(i*ns-900), 4), 1);
    raw2_var(i,1) = var(raw((i*ns-ns+1):(i*ns-900), 4), 1);
end

%% erase duplicate lines (the silence period of simulation)
l2 = floor(l/2) - 1;
raw3 = zeros(l2,1);
raw3_var = zeros(l2,1);
for i=0:l2,
    raw3(i+1,1) = raw2((2*i+1),1);
    raw3_var(i+1,1) = raw2_var((2*i+1),1);
end

%% display graphs
nx = 4;
nt = 100;
res = zeros(nx,nx);
k =1;
for i=0:(nx-1),
    for j=0:(nx-1),
        cur = raw3( ((i*nx + j)*nt + 1):((i*nx + j)*nt + nt - 1), 1);
        cur_var = raw3_var( ((i*nx + j)*nt + 1):((i*nx + j)*nt + nt - 1), 1);
        tr = (cur(nt - 1) > 0) + 2 * (sum(cur < 0) > 0);
        res(i+1,j+1) = tr;
        
         %if ((mod(i,5) == 0) && (mod(j,5) == 0)),
         %if (tr == 2),
         %if (i<=10 && j>10),
            subplot(nx, nx, k); 
            k = k + 1;
            %if (tr == 3),
                plot(  [ 1:0.5:(nt/2) ]', [ cur, zeros(nt-1, 1), cur-cur_var, cur+cur_var ]);
                title( sprintf('i=%d, j=%d', i, j));
            %end;
            %axis([-1 (nt+1) -1 1])
            set(gca,'xtick',0:10:(nt/2))
            %set(gca,'ytick',[])
         %end;
    end;
end;