function [self_action, disease_left_neighbor, disease_right_neighbor] = disease_outcome( disease )
	%% for a given disease (defined below) returns possible outcomes/actions
	% disease probabilities matrix
	% notice that probabilities of  self actions  prob_self
	% i.e. stay ill, die, and heal should add up to 1
	% probabilities distributed in a vector where each positions corresponds to
	% 1 die 
	% 2 heal
    % there is no 3rd element since probability of stay ill = 1 - p_die - p_heal
   
    % array prob_spread probability to spread disease to left(1) and right(2) per cycle 
    % for each disease
	% 1 disease left  neighbor (this one should be <= 1)
	% 2 disease right neighbor (this one should be <= 1)

	% normal dead member  (death is non contagious) 
	prob_self(1,:)=[0.0, 0.0];  % [prob_die, prob_heal]
    prob_spread(1,:)=[0.0, 0.0];

	% weakly contagious but high mortality, hard to heal
	prob_self(2,:)=[0.8, 0.0];
    prob_spread(2,:)=[0.1, 0.1];

	% highly contagious but low mortality, easy to heal
	prob_self(3,:)=[0.1, 0.1];
    prob_spread(3,:)=[0.4, 0.4];

	% healthy alive member  (life is a disease too) 
	prob_self(4,:)=[0.0, 0.0];
    prob_spread(4,:)=[0.0, 0.0];


    % some human readable index definitions
    left=1;
    right=2;
    die=1;
    heal=2;


	% ready to roll dices
	% 1st do we disease anyone
	% roll for dices for left neighbor
	p=rand();
	if ( p  <= prob_spread(disease, left) )
		disease_left_neighbor=true;
	else
		disease_left_neighbor=false;
	end

	p=rand(); % reroll dices for right neighbor
	if( p <= prob_spread(disease, right) )
		disease_right_neighbor=true;
	else
		disease_right_neighbor=false;
	end

	% 2nd what about yourself
	p=rand();

    if ( p <= prob_self(disease, die))
		self_action = 1; %death
    elseif (p <= (prob_self(disease, die) + prob_self(disease,heal)) )
        % notice the sum !
		self_action = 4; %healing,  recovery to normal
	else
		self_action =  disease; % keep what you have
	end
	
end
			




