Logo ROOT  
Reference Guide
RooBreitWigner.cxx
Go to the documentation of this file.
1/*****************************************************************************
2 * Project: RooFit *
3 * Package: RooFitModels *
4 * @(#)root/roofit:$Id$
5 * Authors: *
6 * AS, Abi Soffer, Colorado State University, abi@slac.stanford.edu *
7 * TS, Thomas Schietinger, SLAC, schieti@slac.stanford.edu *
8 * *
9 * Copyright (c) 2000-2005, Regents of the University of California *
10 * Colorado State University *
11 * and Stanford University. All rights reserved. *
12 * *
13 * Redistribution and use in source and binary forms, *
14 * with or without modification, are permitted according to the terms *
15 * listed in LICENSE (http://roofit.sourceforge.net/license.txt) *
16 *****************************************************************************/
17
18/** \class RooBreitWigner
19 \ingroup Roofit
20
21Class RooBreitWigner is a RooAbsPdf implementation
22that models a non-relativistic Breit-Wigner shape
23**/
24
25#include "RooFit.h"
26
27#include "Riostream.h"
28#include <math.h>
29
30#include "RooBreitWigner.h"
31#include "RooAbsReal.h"
32#include "RooRealVar.h"
33#include "RooBatchCompute.h"
34
35using namespace std;
36
38
39////////////////////////////////////////////////////////////////////////////////
40
41RooBreitWigner::RooBreitWigner(const char *name, const char *title,
42 RooAbsReal& _x, RooAbsReal& _mean,
43 RooAbsReal& _width) :
44 RooAbsPdf(name,title),
45 x("x","Dependent",this,_x),
46 mean("mean","Mean",this,_mean),
47 width("width","Width",this,_width)
48{
49}
50
51////////////////////////////////////////////////////////////////////////////////
52
54 RooAbsPdf(other,name), x("x",this,other.x), mean("mean",this,other.mean),
55 width("width",this,other.width)
56{
57}
58
59////////////////////////////////////////////////////////////////////////////////
60
62{
63 Double_t arg= x - mean;
64 return 1. / (arg*arg + 0.25*width*width);
65}
66
67////////////////////////////////////////////////////////////////////////////////
68/// Compute multiple values of BreitWigner distribution.
69void RooBreitWigner::computeBatch(cudaStream_t* stream, double* output, size_t nEvents, RooBatchCompute::DataMap& dataMap) const
70{
72 dispatch->compute(stream, RooBatchCompute::BreitWigner, output, nEvents, dataMap, {&*x,&*mean,&*width,&*_norm});
73}
74
75////////////////////////////////////////////////////////////////////////////////
76
77Int_t RooBreitWigner::getAnalyticalIntegral(RooArgSet& allVars, RooArgSet& analVars, const char* /*rangeName*/) const
78{
79 if (matchArgs(allVars,analVars,x)) return 1 ;
80 return 0 ;
81}
82
83////////////////////////////////////////////////////////////////////////////////
84
85Double_t RooBreitWigner::analyticalIntegral(Int_t code, const char* rangeName) const
86{
87 switch(code) {
88 case 1:
89 {
90 Double_t c = 2./width;
91 return c*(atan(c*(x.max(rangeName)-mean)) - atan(c*(x.min(rangeName)-mean)));
92 }
93 }
94
95 assert(0) ;
96 return 0 ;
97}
#define c(i)
Definition: RSha256.hxx:101
#define ClassImp(name)
Definition: Rtypes.h:364
include TDocParser_001 C image html pict1_TDocParser_001 png width
Definition: TDocParser.cxx:121
char name[80]
Definition: TGX11.cxx:110
RooAbsReal * _norm
Definition: RooAbsPdf.h:364
RooAbsReal is the common abstract base class for objects that represent a real value and implements f...
Definition: RooAbsReal.h:63
Bool_t matchArgs(const RooArgSet &allDeps, RooArgSet &numDeps, const RooArgProxy &a) const
Utility function for use in getAnalyticalIntegral().
RooArgSet is a container object that can hold multiple RooAbsArg objects.
Definition: RooArgSet.h:35
virtual void compute(cudaStream_t *, Computer, RestrictArr, size_t, const DataMap &, const VarVector &, const ArgVector &={})=0
Class RooBreitWigner is a RooAbsPdf implementation that models a non-relativistic Breit-Wigner shape.
Double_t evaluate() const
Evaluate this PDF / function / constant. Needs to be overridden by all derived classes.
void computeBatch(cudaStream_t *, double *output, size_t nEvents, RooBatchCompute::DataMap &) const
Compute multiple values of BreitWigner distribution.
Int_t getAnalyticalIntegral(RooArgSet &allVars, RooArgSet &analVars, const char *rangeName=0) const
Interface function getAnalyticalIntergral advertises the analytical integrals that are supported.
RooRealProxy width
RooRealProxy mean
Double_t analyticalIntegral(Int_t code, const char *rangeName=0) const
Implements the actual analytical integral(s) advertised by getAnalyticalIntegral.
RooRealProxy x
double min(const char *rname=0) const
Query lower limit of range. This requires the payload to be RooAbsRealLValue or derived.
double max(const char *rname=0) const
Query upper limit of range. This requires the payload to be RooAbsRealLValue or derived.
Double_t x[n]
Definition: legend1.C:17
R__EXTERN RooBatchComputeInterface * dispatchCUDA
std::map< DataKey, RooSpan< const double > > DataMap
R__EXTERN RooBatchComputeInterface * dispatchCPU
This dispatch pointer points to an implementation of the compute library, provided one has been loade...
static void output(int code)
Definition: gifencode.c:226