laura is hosted by Hepforge, IPPP Durham
Laura++  v3r4
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauNovosibirskPdf.cc
Go to the documentation of this file.
1 
2 /*
3 Copyright 2008 University of Warwick
4 
5 Licensed under the Apache License, Version 2.0 (the "License");
6 you may not use this file except in compliance with the License.
7 You may obtain a copy of the License at
8 
9  http://www.apache.org/licenses/LICENSE-2.0
10 
11 Unless required by applicable law or agreed to in writing, software
12 distributed under the License is distributed on an "AS IS" BASIS,
13 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 See the License for the specific language governing permissions and
15 limitations under the License.
16 */
17 
18 /*
19 Laura++ package authors:
20 John Back
21 Paul Harrison
22 Thomas Latham
23 */
24 
29 #include <iostream>
30 
31 using std::cout;
32 using std::cerr;
33 using std::endl;
34 
35 #include "TMath.h"
36 #include "TSystem.h"
37 
38 #include "LauNovosibirskPdf.hh"
39 #include "LauConstants.hh"
40 
42 
43 LauNovosibirskPdf::LauNovosibirskPdf(const TString& theVarName, const vector<LauAbsRValue*>& params, Double_t minAbscissa, Double_t maxAbscissa) :
44  LauAbsPdf(theVarName, params, minAbscissa, maxAbscissa),
45  mean_(0),
46  sigma_(0),
47  tail_(0)
48 {
49  // Constructor for the Novosibirsk PDF.
50  //
51  // The parameters in params are the mean, sigma and tail.
52  // The last two arguments specify the range in which the PDF is defined, and the PDF
53  // will be normalised w.r.t. these limits.
54 
55  mean_ = this->findParameter("mean");
56  sigma_ = this->findParameter("sigma");
57  tail_ = this->findParameter("tail");
58 
59  if ((this->nParameters() != 3) || (mean_ == 0) || (sigma_ == 0) || (tail_ == 0)) {
60  cerr<<"ERROR in LauNovosibirskPdf constructor: LauNovosibirskPdf requires 3 parameters: \"mean\", \"sigma\"and \"tail\" "<<endl;
61  gSystem->Exit(EXIT_FAILURE);
62  }
63 
64  // Cache the normalisation factor
65  this->calcNorm();
66 }
67 
69 {
70  // Destructor
71 }
72 
74 {
75  // Check that the given abscissa is within the allowed range
76  if (!this->checkRange(abscissas)) {
77  gSystem->Exit(EXIT_FAILURE);
78  }
79 
80  // Get our abscissa
81  Double_t abscissa = abscissas[0];
82 
83  // Get the up to date parameter values
84  Double_t mean = mean_->unblindValue();
85  Double_t sigma = sigma_->unblindValue();
86  Double_t tail = tail_->unblindValue();
87 
88  // Evaluate the Novosibirsk PDF value
89 
90 
91  Double_t qa(0.0), qb(0.0), qx(0.0), qy(0.0);
92  Double_t arg(0.0);
93  Double_t value(0.0);
94 
95 
96  if(TMath::Abs(tail) < 1.e-7)
97  arg = 0.5*((abscissa - mean)/sigma)*((abscissa - mean)/sigma);
98  else {
99  qa = tail*TMath::Sqrt(LauConstants::log4);
100  qb = TMath::SinH(qa)/qa;
101  qx = (abscissa - mean)/sigma*qb;
102  qy = 1.0+ tail*qx;
103 
104  //---- Cutting curve from right side
105 
106  if( qy > 1.E-7) {
107  arg = 0.5*( (log(qy)/tail)*(log(qy)/tail) + tail*tail);
108  }else{
109  arg = 15.0;
110  }
111  }
112 
113  value = TMath::Exp(-arg);
114 
115  // if the parameters are floating then we
116  // need to recalculate the normalisation
117  if (!this->cachePDF() && !this->withinNormCalc() && !this->withinGeneration()) {
118  this->calcNorm();
119  }
120  this->setUnNormPDFVal(value);
121 }
122 
123 void LauNovosibirskPdf::calcPDFHeight( const LauKinematics* /*kinematics*/ )
124 {
125  if (this->heightUpToDate()) {
126  return;
127  }
128 
129  // Get the up to date parameter values
130  Double_t mean = mean_->unblindValue();
131 
132  LauAbscissas maxPoint(1);
133  maxPoint[0] = mean;
134 
135 
136  // Calculate the PDF height for the Bifurcated Gaussian function.
137 
138  if (mean < this->getMinAbscissa()) {
139  maxPoint[0] = this->getMinAbscissa();
140  } else if (mean > this->getMaxAbscissa()) {
141  maxPoint[0] = this->getMaxAbscissa();
142  }
143 
144  this->calcLikelihoodInfo(maxPoint);
145  Double_t height = this->getUnNormLikelihood();
146 
147  // Multiply by a small factor to avoid problems from rounding errors
148  height *= (1.0 + 1e-1);
149 
150  this->setMaxHeight(height);
151 }
152 
virtual void setUnNormPDFVal(Double_t unNormPDFVal)
Set the unnormalised likelihood.
Definition: LauAbsPdf.hh:383
virtual Double_t getMinAbscissa() const
Retrieve the minimum value of the (primary) abscissa.
Definition: LauAbsPdf.hh:131
Class for defining a Novosibirsk function PDF.
virtual Bool_t heightUpToDate() const
Check if the maximum height of the PDF is up to date.
Definition: LauAbsPdf.hh:278
ClassImp(LauAbsCoeffSet)
virtual void calcPDFHeight(const LauKinematics *kinematics)
Calculate the PDF height.
virtual Double_t getUnNormLikelihood() const
Retrieve the unnormalised likelihood value.
Definition: LauAbsPdf.hh:210
File containing declaration of LauNovosibirskPdf class.
virtual Bool_t checkRange(const LauAbscissas &abscissas) const
Check that all abscissas are within their allowed ranges.
Definition: LauAbsPdf.cc:227
virtual ~LauNovosibirskPdf()
Destructor.
virtual Bool_t withinNormCalc() const
Check whether the calcNorm method is running.
Definition: LauAbsPdf.hh:437
virtual Double_t unblindValue() const =0
The unblinded value of the parameter.
LauAbsRValue * sigma_
Gaussian sigma.
virtual Double_t getMaxAbscissa() const
Retrieve the maximum value of the (primary) abscissa.
Definition: LauAbsPdf.hh:137
virtual void setMaxHeight(Double_t maxHeight)
Set the maximum height.
Definition: LauAbsPdf.hh:345
LauAbsRValue * tail_
Gaussian tail.
const Double_t log4
Logarithm of four.
virtual Bool_t withinGeneration() const
Check whether the generate method is running.
Definition: LauAbsPdf.hh:449
virtual Bool_t cachePDF() const
Check if the PDF is to be cached.
Definition: LauAbsPdf.hh:290
virtual void calcLikelihoodInfo(const LauAbscissas &abscissas)
Calculate the likelihood (and intermediate info) for a given abscissa.
File containing LauConstants namespace.
virtual void calcNorm()
Calculate the normalisation factor of the PDF.
Definition: LauAbsPdf.cc:444
Class for defining the abstract interface for PDF classes.
Definition: LauAbsPdf.hh:55
LauAbsRValue * mean_
Gaussian mean.
Class for calculating 3-body kinematic quantities.
Double_t value() const
The value of the parameter.
Pure abstract base class for defining a parameter containing an R value.
Definition: LauAbsRValue.hh:43
std::vector< Double_t > LauAbscissas
The type used for containing multiple abscissa values.
Definition: LauAbsPdf.hh:59