% fMRI_SIMULATION

% APMA 0410
% Fall 2008

clear

% PARAMETERS
N = 100;            % linear size of array. (Number of voxels is N*N)
N_modes = 3;        % Number of modes in each condition
breadth = 6;        % Half-width of modes


centers_L = breadth + floor((N-2*breadth)*rand(N_modes,2));    % The modes in the LIE  condition
centers_T = breadth + floor((N-2*breadth)*rand(N_modes,2));    % The modes in the TRUE condition

[X,Y] = meshgrid(1:N,1:N);

N_trials = 100;     % Number of trials in each condition
a_peak = 1;         % peak voxel activation
a_noise = 10;       % noise amplitude


act_T = zeros(N,N); % activation in T condition
act_L = zeros(N,N); % activation in L condition

for mode = 1:N_modes
    
    center_X = centers_T(mode,1);
    center_Y = centers_T(mode,2);
    act_T = act_T + a_peak*exp(-((X-center_X).^2 + (Y-center_Y).^2) / breadth^2 );
    center_X = centers_L(mode,1);
    center_Y = centers_L(mode,2);
    act_L = act_L + a_peak*exp(-((X-center_X).^2 + (Y-center_Y).^2) / breadth^2 );
    
end

for trial = 1:N_trials
    data_T(:,:,trial) = act_T + a_noise*randn(N,N);
    data_L(:,:,trial) = act_L + a_noise*randn(N,N);
end

% [0,1]-NORMALIZATION:

min_data = min([data_L(:);data_T(:)]);
data_T = data_T - min_data;
data_L = data_L - min_data;

max_data = max([data_L(:);data_T(:)]);
data_T = data_T/max_data;
data_L = data_L/max_data;

save fMRI_dataset_1 data_T data_L


% Creating fake data (same modes for the two conditions):

centers = breadth + floor((N-2*breadth)*rand(N_modes,2));    % The modes in both conditions

act = zeros(N,N);   % activation in all conditions

for mode = 1:N_modes
    center_X = centers(mode,1);
    center_Y = centers(mode,2);
    act = act + a_peak*exp(-((X-center_X).^2 + (Y-center_Y).^2) / breadth^2 );
end

for trial = 1:N_trials
    data_T(:,:,trial) = act + a_noise*randn(N,N);
    data_L(:,:,trial) = act + a_noise*randn(N,N);
end

% [0,1]-NORMALIZATION:

min_data = min([data_L(:);data_T(:)]);
data_T = data_T - min_data;
data_L = data_L - min_data;

max_data = max([data_L(:);data_T(:)]);
data_T = data_T/max_data;
data_L = data_L/max_data;

save fMRI_dataset_2 data_T data_L
