#include "postgres.h"
#include <math.h>
#include "utils/sampling.h"

Include dependency graph for sampling.c:

Functions
void	BlockSampler_Init (BlockSampler bs, BlockNumber nblocks, int samplesize, long randseed)

bool	BlockSampler_HasMore (BlockSampler bs)

BlockNumber	BlockSampler_Next (BlockSampler bs)

void	reservoir_init_selection_state (ReservoirState rs, int n)

double	reservoir_get_next_S (ReservoirState rs, double t, int n)

void	sampler_random_init_state (long seed, SamplerRandomState randstate)

double	sampler_random_fract (SamplerRandomState randstate)

double	anl_random_fract (void)

double	anl_init_selection_state (int n)

double	anl_get_next_S (double t, int n, double *stateptr)

Variables
static ReservoirStateData	oldrs

Function Documentation

double anl_get_next_S	(	double	t,
		int	n,
		double *	stateptr
	)

Definition at line 284 of file sampling.c.

References reservoir_get_next_S(), result, and ReservoirStateData::W.

 {
     double      result;
 
     oldrs.W = *stateptr;
     result = reservoir_get_next_S(&oldrs, t, n);
     *stateptr = oldrs.W;
     return result;
 }

double anl_init_selection_state ( int n )

Definition at line 273 of file sampling.c.

References random(), ReservoirStateData::randstate, sampler_random_fract(), and sampler_random_init_state().

 {
     /* initialize if first time through */
     if (oldrs.randstate[0] == 0)
         sampler_random_init_state(random(), oldrs.randstate);
 
     /* Initial value of W (for use when Algorithm Z is first applied) */
     return exp(-log(sampler_random_fract(oldrs.randstate)) / n);
 }

double anl_random_fract ( void )

Definition at line 262 of file sampling.c.

References random(), ReservoirStateData::randstate, sampler_random_fract(), and sampler_random_init_state().

 {
     /* initialize if first time through */
     if (oldrs.randstate[0] == 0)
         sampler_random_init_state(random(), oldrs.randstate);
 
     /* and compute a random fraction */
     return sampler_random_fract(oldrs.randstate);
 }

bool BlockSampler_HasMore ( BlockSampler bs )

Definition at line 54 of file sampling.c.

References BlockSamplerData::m, BlockSamplerData::N, BlockSamplerData::n, and BlockSamplerData::t.

Referenced by acquire_sample_rows(), and BlockSampler_Next().

 {
     return (bs->t < bs->N) && (bs->m < bs->n);
 }

void BlockSampler_Init	(	BlockSampler	bs,
		BlockNumber	nblocks,
		int	samplesize,
		long	randseed
	)

Definition at line 37 of file sampling.c.

References BlockSamplerData::m, BlockSamplerData::N, BlockSamplerData::n, BlockSamplerData::randstate, sampler_random_init_state(), and BlockSamplerData::t.

Referenced by acquire_sample_rows().

 {
     bs->N = nblocks;            /* measured table size */
 
     /*
      * If we decide to reduce samplesize for tables that have less or not much
      * more than samplesize blocks, here is the place to do it.
      */
     bs->n = samplesize;
     bs->t = 0;                  /* blocks scanned so far */
     bs->m = 0;                  /* blocks selected so far */
 
     sampler_random_init_state(randseed, bs->randstate);
 }

BlockNumber BlockSampler_Next ( BlockSampler bs )

Definition at line 60 of file sampling.c.

References Assert, BlockSampler_HasMore(), K, BlockSamplerData::m, BlockSamplerData::N, BlockSamplerData::n, BlockSamplerData::randstate, sampler_random_fract(), and BlockSamplerData::t.

Referenced by acquire_sample_rows().

 {
     BlockNumber K = bs->N - bs->t;  /* remaining blocks */
     int         k = bs->n - bs->m;  /* blocks still to sample */
     double      p;              /* probability to skip block */
     double      V;              /* random */
 
     Assert(BlockSampler_HasMore(bs));   /* hence K > 0 and k > 0 */
 
     if ((BlockNumber) k >= K)
     {
         /* need all the rest */
         bs->m++;
         return bs->t++;
     }
 
     /*----------
      * It is not obvious that this code matches Knuth's Algorithm S.
      * Knuth says to skip the current block with probability 1 - k/K.
      * If we are to skip, we should advance t (hence decrease K), and
      * repeat the same probabilistic test for the next block.  The naive
      * implementation thus requires a sampler_random_fract() call for each
      * block number.  But we can reduce this to one sampler_random_fract()
      * call per selected block, by noting that each time the while-test
      * succeeds, we can reinterpret V as a uniform random number in the range
      * 0 to p. Therefore, instead of choosing a new V, we just adjust p to be
      * the appropriate fraction of its former value, and our next loop
      * makes the appropriate probabilistic test.
      *
      * We have initially K > k > 0.  If the loop reduces K to equal k,
      * the next while-test must fail since p will become exactly zero
      * (we assume there will not be roundoff error in the division).
      * (Note: Knuth suggests a "<=" loop condition, but we use "<" just
      * to be doubly sure about roundoff error.)  Therefore K cannot become
      * less than k, which means that we cannot fail to select enough blocks.
      *----------
      */
     V = sampler_random_fract(bs->randstate);
     p = 1.0 - (double) k / (double) K;
     while (V < p)
     {
         /* skip */
         bs->t++;
         K--;                    /* keep K == N - t */
 
         /* adjust p to be new cutoff point in reduced range */
         p *= 1.0 - (double) k / (double) K;
     }
 
     /* select */
     bs->m++;
     return bs->t++;
 }

double reservoir_get_next_S	(	ReservoirState	rs,
		double	t,
		int	n
	)

Definition at line 142 of file sampling.c.

References ReservoirStateData::randstate, S, sampler_random_fract(), ReservoirStateData::W, and W.

Referenced by acquire_sample_rows(), analyze_row_processor(), anl_get_next_S(), and file_acquire_sample_rows().

 {
     double      S;
 
     /* The magic constant here is T from Vitter's paper */
     if (t <= (22.0 * n))
     {
         /* Process records using Algorithm X until t is large enough */
         double      V,
                     quot;
 
         V = sampler_random_fract(rs->randstate);    /* Generate V */
         S = 0;
         t += 1;
         /* Note: "num" in Vitter's code is always equal to t - n */
         quot = (t - (double) n) / t;
         /* Find min S satisfying (4.1) */
         while (quot > V)
         {
             S += 1;
             t += 1;
             quot *= (t - (double) n) / t;
         }
     }
     else
     {
         /* Now apply Algorithm Z */
         double      W = rs->W;
         double      term = t - (double) n + 1;
 
         for (;;)
         {
             double      numer,
                         numer_lim,
                         denom;
             double      U,
                         X,
                         lhs,
                         rhs,
                         y,
                         tmp;
 
             /* Generate U and X */
             U = sampler_random_fract(rs->randstate);
             X = t * (W - 1.0);
             S = floor(X);       /* S is tentatively set to floor(X) */
             /* Test if U <= h(S)/cg(X) in the manner of (6.3) */
             tmp = (t + 1) / term;
             lhs = exp(log(((U * tmp * tmp) * (term + S)) / (t + X)) / n);
             rhs = (((t + X) / (term + S)) * term) / t;
             if (lhs <= rhs)
             {
                 W = rhs / lhs;
                 break;
             }
             /* Test if U <= f(S)/cg(X) */
             y = (((U * (t + 1)) / term) * (t + S + 1)) / (t + X);
             if ((double) n < S)
             {
                 denom = t;
                 numer_lim = term + S;
             }
             else
             {
                 denom = t - (double) n + S;
                 numer_lim = t + 1;
             }
             for (numer = t + S; numer >= numer_lim; numer -= 1)
             {
                 y *= numer / denom;
                 denom -= 1;
             }
             W = exp(-log(sampler_random_fract(rs->randstate)) / n); /* Generate W in advance */
             if (exp(log(y) / n) <= (t + X) / t)
                 break;
         }
         rs->W = W;
     }
     return S;
 }

void reservoir_init_selection_state	(	ReservoirState	rs,
		int	n
	)

Definition at line 129 of file sampling.c.

References random(), ReservoirStateData::randstate, sampler_random_fract(), sampler_random_init_state(), and ReservoirStateData::W.

Referenced by acquire_sample_rows(), file_acquire_sample_rows(), and postgresAcquireSampleRowsFunc().

 {
     /*
      * Reservoir sampling is not used anywhere where it would need to return
      * repeatable results so we can initialize it randomly.
      */
     sampler_random_init_state(random(), rs->randstate);
 
     /* Initial value of W (for use when Algorithm Z is first applied) */
     rs->W = exp(-log(sampler_random_fract(rs->randstate)) / n);
 }

double sampler_random_fract ( SamplerRandomState randstate )

Definition at line 238 of file sampling.c.

References pg_erand48().

Referenced by acquire_sample_rows(), analyze_row_processor(), anl_init_selection_state(), anl_random_fract(), BlockSampler_Next(), file_acquire_sample_rows(), random_relative_prime(), reservoir_get_next_S(), reservoir_init_selection_state(), system_rows_nextsampleblock(), and system_time_nextsampleblock().

 {
     double      res;
 
     /* pg_erand48 returns a value in [0.0 - 1.0), so we must reject 0 */
     do
     {
         res = pg_erand48(randstate);
     } while (res == 0.0);
     return res;
 }

void sampler_random_init_state	(	long	seed,
		SamplerRandomState	randstate
	)

Definition at line 229 of file sampling.c.

Referenced by anl_init_selection_state(), anl_random_fract(), BlockSampler_Init(), reservoir_init_selection_state(), system_rows_nextsampleblock(), and system_time_nextsampleblock().

 {
     randstate[0] = 0x330e;      /* same as pg_erand48, but could be anything */
     randstate[1] = (unsigned short) seed;
     randstate[2] = (unsigned short) (seed >> 16);
 }

Variable Documentation

ReservoirStateData oldrs

static

Definition at line 259 of file sampling.c.

Functions

Variables

Function Documentation

Variable Documentation