residualR.m
上传用户:hfch80
上传日期:2007-10-25
资源大小:3637k
文件大小:1k
源码类别:

行业发展研究

开发平台:

Matlab

  1. function outIndex = residualR(inIndex,q);
  2. % PURPOSE : Performs the resampling stage of the SIR
  3. %           in order(number of samples) steps. It uses Liu's
  4. %           residual resampling algorithm and Niclas' magic line.
  5. % INPUTS  : - inIndex = Input particle indices.
  6. %           - q = Normalised importance ratios.
  7. % OUTPUTS : - outIndex = Resampled indices.
  8. % AUTHORS  : Arnaud Doucet and Nando de Freitas - Thanks for the acknowledgement.
  9. % DATE     : 08-09-98
  10. if nargin < 2, error('Not enough input arguments.'); end
  11. [S,arb] = size(q);  % S = Number of particles.
  12. % RESIDUAL RESAMPLING:
  13. % ====================
  14. N_babies= zeros(1,S);
  15. % first integer part
  16. q_res = S.*q'; %'
  17. N_babies = fix(q_res);
  18. % residual number of particles to sample
  19. N_res=S-sum(N_babies);
  20. if (N_res~=0)
  21.   q_res=(q_res-N_babies)/N_res;
  22.   cumDist= cumsum(q_res);   
  23.   % generate N_res ordered random variables uniformly distributed in [0,1]
  24.   u = fliplr(cumprod(rand(1,N_res).^(1./(N_res:-1:1))));
  25.   j=1;
  26.   for i=1:N_res
  27.     while (u(1,i)>cumDist(1,j))
  28.       j=j+1;
  29.     end
  30.     N_babies(1,j)=N_babies(1,j)+1;
  31.   end;
  32. end;
  33. % COPY RESAMPLED TRAJECTORIES:  
  34. % ============================
  35. index=1;
  36. for i=1:S
  37.   if (N_babies(1,i)>0)
  38.     for j=index:index+N_babies(1,i)-1
  39.       outIndex(j) = inIndex(i);
  40.     end;
  41.   end;   
  42.   index= index+N_babies(1,i);   
  43. end