  We present an unbiased numerical integration algorithm that handles both
low-frequency regions and high frequency details of multidimensional integrals.
It combines quadrature and Monte Carlo integration, by using a quadrature-base
approximation as a control variate of the signal. We adaptively build the
control variate constructed as a piecewise polynomial, which can be
analytically integrated, and accurately reconstructs the low frequency regions
of the integrand. We then recover the high-frequency details missed by the
control variate by using Monte Carlo integration of the residual. Our work
leverages importance sampling techniques by working in primary space, allowing
the combination of multiple mappings; this enables multiple importance sampling
in quadrature-based integration. Our algorithm is generic, and can be applied
to any complex multidimensional integral. We demonstrate its effectiveness with
four applications with low dimensionality: transmittance estimation in
heterogeneous participating media, low-order scattering in homogeneous media,
direct illumination computation, and rendering of distributed effects. Finally,
we show how our technique is extensible to integrands of higher dimensionality,
by computing the control variate on Monte Carlo estimates of the
high-dimensional signal, and accounting for such additional dimensionality on
the residual as well. In all cases, we show accurate results and faster
convergence compared to previous approaches.
