Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
20 changes: 9 additions & 11 deletions roofit/roofit/inc/RooVoigtian.h
Original file line number Diff line number Diff line change
Expand Up @@ -16,27 +16,25 @@
#ifndef ROO_VOIGTIAN
#define ROO_VOIGTIAN

#include "RooAbsPdf.h"
#include "RooRealProxy.h"

class RooRealVar;
#include <RooAbsPdf.h>
#include <RooRealProxy.h>

class RooVoigtian : public RooAbsPdf {
public:
RooVoigtian() {} ;
RooVoigtian() {}
RooVoigtian(const char *name, const char *title,
RooAbsReal& _x, RooAbsReal& _mean,
RooAbsReal& _width, RooAbsReal& _sigma,
bool doFast = false);
RooVoigtian(const RooVoigtian& other, const char* name=nullptr) ;
TObject* clone(const char* newname) const override { return new RooVoigtian(*this,newname); }
inline ~RooVoigtian() override { }

// These methods allow the user to select the fast evaluation
// of the complex error function using look-up tables
// (default is the "slow" CERNlib algorithm)

/// Enable the fast evaluation of the complex error function using look-up
/// tables (default is the "slow" CERNlib algorithm).
inline void selectFastAlgorithm() { _doFast = true; }

/// Disable the fast evaluation of the complex error function using look-up
/// tables (default is the "slow" CERNlib algorithm).
inline void selectDefaultAlgorithm() { _doFast = false; }

protected:
Expand All @@ -52,7 +50,7 @@ class RooVoigtian : public RooAbsPdf {

private:

bool _doFast;
bool _doFast = false;
ClassDefOverride(RooVoigtian,2) // Voigtian PDF (Gauss (x) BreitWigner)
};

Expand Down
26 changes: 20 additions & 6 deletions roofit/roofit/src/RooVoigtian.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -23,20 +23,34 @@ function (the default CERNlib C335 algorithm, and a faster, look-up-table
based method). By default, RooVoigtian employs the default (CERNlib)
algorithm. Select the faster algorithm either in the constructor, or with
the selectFastAlgorithm() method.

\Note The "width" parameter that determines the Breit-Wigner shape
represents the **full width at half maximum (FWHM)** of the
Breit-Wigner (often referred to as \f$\Gamma\f$ or \f$2\gamma\f$).
**/

#include "RooVoigtian.h"
#include "RooRealVar.h"
#include "RooMath.h"
#include "RooBatchCompute.h"
#include <RooVoigtian.h>

#include <RooMath.h>
#include <RooBatchCompute.h>

#include <cmath>
#include <complex>
using namespace std;

ClassImp(RooVoigtian);

////////////////////////////////////////////////////////////////////////////////
/// Construct a RooVoigtian PDF, which represents the convolution of a
/// Breit-Wigner with a Gaussian.
/// \param name Name that identifies the PDF in computations.
/// \param title Title for plotting.
/// \param _x The observable for the PDF.
/// \param _mean The mean of the distribution.
/// \param _width The **full width at half maximum (FWHM)** of the Breit-Wigner
/// (often referred to as \f$\Gamma\f$ or \f$2\gamma\f$).
/// \param _sigma The width of the Gaussian distribution.
/// \param doFast Use the faster look-up-table-based method for the evaluation
/// of the complex error function.

RooVoigtian::RooVoigtian(const char *name, const char *title,
RooAbsReal& _x, RooAbsReal& _mean,
Expand Down Expand Up @@ -79,7 +93,7 @@ double RooVoigtian::evaluate() const
if (s==0.) return (1./(arg*arg+0.25*w*w));

// Gauss for zero width
if (w==0.) return exp(coef*arg*arg);
if (w==0.) return std::exp(coef*arg*arg);

// actual Voigtian for non-trivial width and sigma
double c = 1./(sqrt(2.)*s);
Expand Down