32 #ifndef VSMC_RESAMPLE_RESIDUAL_STRATIFIED_HPP
33 #define VSMC_RESAMPLE_RESIDUAL_STRATIFIED_HPP
41 typedef cxx11::integral_constant<ResampleScheme, ResidualStratified>
53 template <
typename IntType,
typename RngType>
54 void operator() (std::size_t M, std::size_t N, RngType &rng,
55 const double *weight, IntType *replication)
61 double *
const rptr = &residual_[0];
62 double *
const iptr = &integral_[0];
63 for (std::size_t i = 0; i != M; ++i)
64 rptr[i] = modf(N * weight[i], iptr + i);
66 for (std::size_t i = 0; i != M; ++i)
68 for (std::size_t i = 0; i != M; ++i)
72 for (std::size_t i = 0; i != M; ++i)
73 R += static_cast<IntType>(iptr[i]);
74 std::size_t NN = N -
static_cast<std::size_t
>(R);
76 double *
const uptr = &u01_[0];
77 cxx11::uniform_real_distribution<double> runif(0, 1);
78 const double delta = 1.0 / N;
79 for (std::size_t i = 0; i != NN; ++i)
80 uptr[i] = runif(rng) * delta + i * delta;
81 inversion_(M, NN, rptr, uptr, replication,
true);
83 for (std::size_t i = 0; i != N; ++i)
84 replication[i] += static_cast<IntType>(iptr[i]);
90 std::vector<double, AlignedAllocator<double> > residual_;
91 std::vector<double, AlignedAllocator<double> > integral_;
92 std::vector<double, AlignedAllocator<double> > u01_;
97 #endif // VSMC_RESAMPLE_RESIDUAL_STRATIFIED_HPP
cxx11::integral_constant< ResampleScheme, ResidualStratified > ResampleResidualStratified
Resample forward decleration.