// \begin{verbatim}
function w=distorsion(p,ggamma)
w= p^ggamma ./ ( ( p^ggamma + (1-p)^ggamma )^{1/ggamma} )
endfunction

ggamma_plus= 0.61;
ggamma_minus= 0.69;

function pi=weight(proba,ggamma)
cumul_proba=cumsum(proba);
distort_cumul_proba = distorsion( cumul_proba , ggamma ) ;
ip = distort_cumul_proba - [0 distort_cumul_proba(1:($-1))];
pi=ip($:(-1):1);
endfunction

alpha=0.88;
bbeta=0.88;
lambda=2.25;

function [H] = Heavyside(x)
     // Heavyside function 
     H = bool2s(x>=0)
     // notice that Heavyside(0)=1
 endfunction

function u=utility(x)
  u=Heavyside(x).*x^alpha - lambda* Heavyside(-x).* (-x)^bbeta
endfunction

function v=CPT_eval(lotery,anchor)
values=lotery(1,:)-anchor;
probabilities=lotery(2,:);
//
ind_plus=find(values > 0) ;
values_plus=values(values > 0) ;
proba_plus = probabilities(ind_plus);
aborp_plus=proba_plus($:(-1):1);
utility_plus=utility(values_plus) ;
weighting_plus= weight ( aborp_plus , ggamma_plus) ;
//
ind_minus=find(values <= 0) ;
values_minus=values(values <= 0) ;
proba_minus = probabilities(ind_minus);
utility_minus=utility(values_minus) ;
weighting_minus= weight ( proba_minus , ggamma_minus) ;
// 
v= sum(utility_plus .* weighting_plus) + ...
   sum(utility_minus .* weighting_minus)
endfunction

//\end{verbatim}