% post_sampling_f
% have sampled from the prior and have observed number of failures. Use the
% observations to give weight to each of the sampled values
for k=1:KK      %for each sampled value...
        
    l_weight=0;
    s_r = samples(k,2) + samples(k,1)*exp(-samples(k,3).*(0:I));    %compute the failure rate for the sampled values 
    expected_fails = expected_function(samples(k,1),samples(k,2),samples(k,3));
    variance_fails = variances_func(samples(k,1),samples(k,2),samples(k,3));
    for i=2:I   %first obs is zero since components are working when they are produced.
        l_weight = l_weight + log_pnorm(obs(i),expected_fails(i),variance_fails(i));
    end
    
    weight(k) =exp(l_weight);
end
clear l_weight s_r expected_fails variance_fails k
display('computing weights')

maxt=max(weight);
temp = sum(exp(log(weight(1:KK)-log(maxt))));
temp2(1:KK) = log(weight(1:KK)) - log(maxt);
temp3=exp(temp2);
clear temp2;
sum = exp(log(maxt) + log(sum(temp3)));


        s=sum(weight)
weight = weight/s;
s_cdf=cumsum(weight);

for i=1:MM
    u=rand(1);  
    sample_x = 1;
        while (u > s_cdf(sample_x))
            sample_x=sample_x+1;
        end
    posterior_sim(i,:)=samples(sample_x,:);
end

posterior_sim_stored(count,:,:) = posterior_sim(:,:);
clear sample_x u i