gsl_randist.h

Go to the documentation of this file.

/*
 * $Id: gsl_randist.h 17393 2014-08-23 20:19:14Z predoehl $
 */

#ifndef GSL_RANDIST_H_INCLUDED_LIBKJB_UOFARIZONAVISION
#define GSL_RANDIST_H_INCLUDED_LIBKJB_UOFARIZONAVISION

#include <l_cpp/l_exception.h>
#include <m_cpp/m_vector.h>
#include <gsl_cpp/gsl_rng.h>

#ifdef KJB_HAVE_GSL
// these are GSL headers, not our wrapper
#include "gsl/gsl_randist.h"
#else
#warning "Compiling GNU GSL wrapper without GNU GSL; it will not run properly"
#endif


namespace kjb {

class Gsl_ran_discrete {

#ifdef KJB_HAVE_GSL
    gsl_ran_discrete_t* m_opaque;
#endif

public:

    Gsl_ran_discrete( size_t event_count, const double* event_probs )
#ifdef KJB_HAVE_GSL
    :   m_opaque( gsl_ran_discrete_preproc( event_count, event_probs ) )
    {
        ETX_2( 00 == m_opaque, "Gsl_ran_discrete ctor:  bad alloc" );
    }
#else
    {
        KJB_THROW_2( Missing_dependency, "GNU GSL" );
    }
#endif

    Gsl_ran_discrete( const Vector& event_probs )
#ifdef KJB_HAVE_GSL
    :   m_opaque( 00 )
    {
        Gsl_ran_discrete grd(
                            event_probs.size(),
                            event_probs.get_c_vector() -> elements
                        );
        swap( grd );
    }
#else
    {
        KJB_THROW_2( Missing_dependency, "GNU GSL" );
    }
#endif


    Gsl_ran_discrete( const std::vector< double >& event_probs )
#ifdef KJB_HAVE_GSL
    :   m_opaque( 00 )
    {
        // Herb Sutter tells me the following code is legit:
        // http://herbsutter.com/2008/04/07/cringe-not-vectors-are-guaranteed-to-be-contiguous/
        Gsl_ran_discrete grd( event_probs.size(), & event_probs[ 0 ] );
        swap( grd );
    }
#else
    {
        KJB_THROW_2( Missing_dependency, "GNU GSL" );
    }
#endif

    double sample( const gsl_rng* rng )
    {
#ifdef KJB_HAVE_GSL
        return gsl_ran_discrete( rng, m_opaque );
#endif
    }

    double pdf( size_t k )
    {
#ifdef KJB_HAVE_GSL
        return gsl_ran_discrete_pdf( k, m_opaque );
#endif
    }


    void swap( Gsl_ran_discrete& that )
    {
#ifdef KJB_HAVE_GSL
        using std::swap;

        swap( m_opaque, that.m_opaque );
#endif
    }


#ifdef KJB_HAVE_GSL
    ~Gsl_ran_discrete()
    {
        gsl_ran_discrete_free( m_opaque );
    }
#endif

};


} // end ns kjb

#endif /* GSL_RANDIST_H_INCLUDED_LIBKJB_UOFARIZONAVISION */