laura is hosted by Hepforge, IPPP Durham
Laura++  v1r0
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauCPFitModel.cc
Go to the documentation of this file.
1 
2 // Copyright University of Warwick 2004 - 2013.
3 // Distributed under the Boost Software License, Version 1.0.
4 // (See accompanying file LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
5 
6 // Authors:
7 // Thomas Latham
8 // John Back
9 // Paul Harrison
10 
15 #include <iostream>
16 #include <iomanip>
17 #include <fstream>
18 #include <vector>
19 using std::cout;
20 using std::cerr;
21 using std::endl;
22 using std::setprecision;
23 using std::vector;
24 
25 #include "TVirtualFitter.h"
26 #include "TSystem.h"
27 #include "TMinuit.h"
28 #include "TRandom.h"
29 #include "TFile.h"
30 #include "TMath.h"
31 #include "TH2.h"
32 
33 #include "LauAbsBkgndDPModel.hh"
34 #include "LauAbsCoeffSet.hh"
35 #include "LauAbsDPDynamics.hh"
36 #include "LauAbsPdf.hh"
37 #include "LauAsymmCalc.hh"
38 #include "LauComplex.hh"
39 #include "LauConstants.hh"
40 #include "LauCPFitModel.hh"
41 #include "LauDaughters.hh"
42 #include "LauEffModel.hh"
43 #include "LauEmbeddedData.hh"
44 #include "LauFitNtuple.hh"
45 #include "LauKinematics.hh"
46 #include "LauPrint.hh"
47 #include "LauRandom.hh"
48 #include "LauScfMap.hh"
49 
50 ClassImp(LauCPFitModel)
51 
52 
53 LauCPFitModel::LauCPFitModel(LauAbsDPDynamics* negModel, LauAbsDPDynamics* posModel, Bool_t tagged, const TString& tagVarName) : LauAbsFitModel(),
54  negSigModel_(negModel), posSigModel_(posModel),
55  negKinematics_(negModel ? negModel->getKinematics() : 0),
56  posKinematics_(posModel ? posModel->getKinematics() : 0),
57  usingBkgnd_(kFALSE),
58  nSigComp_(0),
59  nSigDPPar_(0),
60  nExtraPdfPar_(0),
61  nNormPar_(0),
62  negMeanEff_("negMeanEff",0.0,0.0,1.0), posMeanEff_("posMeanEff",0.0,0.0,1.0),
63  negDPRate_("negDPRate",0.0,0.0,100.0), posDPRate_("posDPRate",0.0,0.0,100.0),
64  signalEvents_("signalEvents",1.0,0.0,1.0),
65  signalAsym_("signalAsym",0.0,-1.0,1.0),
66  forceAsym_(kFALSE),
67  tagged_(tagged),
68  tagVarName_(tagVarName),
69  curEvtCharge_(0),
70  useSCF_(kFALSE),
71  useSCFHist_(kFALSE),
72  scfFrac_("scfFrac",0.0,0.0,1.0),
73  scfFracHist_(0),
74  scfMap_(0),
75  compareFitData_(kFALSE),
76  negParent_("B-"), posParent_("B+"),
77  negSignalTree_(0), posSignalTree_(0),
78  reuseSignal_(kFALSE),
79  useNegReweighting_(kFALSE), usePosReweighting_(kFALSE),
80  sigDPLike_(0.0),
81  scfDPLike_(0.0),
82  sigExtraLike_(0.0),
83  scfExtraLike_(0.0),
84  sigTotalLike_(0.0),
85  scfTotalLike_(0.0)
86 {
87  LauDaughters* negDaug = negSigModel_->getDaughters();
88  if (negDaug != 0) {negParent_ = negDaug->getNameParent();}
89  LauDaughters* posDaug = posSigModel_->getDaughters();
90  if (posDaug != 0) {posParent_ = posDaug->getNameParent();}
91 }
92 
93 LauCPFitModel::~LauCPFitModel()
94 {
95  delete negSignalTree_;
96  delete posSignalTree_;
97  for (LauBkgndEmbDataList::iterator iter = negBkgndTree_.begin(); iter != negBkgndTree_.end(); ++iter) {
98  delete (*iter);
99  }
100  for (LauBkgndEmbDataList::iterator iter = posBkgndTree_.begin(); iter != posBkgndTree_.end(); ++iter) {
101  delete (*iter);
102  }
103  delete scfFracHist_;
104 }
105 
106 void LauCPFitModel::setupBkgndVectors()
107 {
108  UInt_t nBkgnds = this->nBkgndClasses();
109  negBkgndDPModels_.resize( nBkgnds );
110  posBkgndDPModels_.resize( nBkgnds );
111  negBkgndPdfs_.resize( nBkgnds );
112  posBkgndPdfs_.resize( nBkgnds );
113  bkgndEvents_.resize( nBkgnds );
114  bkgndAsym_.resize( nBkgnds );
115  negBkgndTree_.resize( nBkgnds );
116  posBkgndTree_.resize( nBkgnds );
117  reuseBkgnd_.resize( nBkgnds );
118  bkgndDPLike_.resize( nBkgnds );
119  bkgndExtraLike_.resize( nBkgnds );
120  bkgndTotalLike_.resize( nBkgnds );
121 }
122 
123 void LauCPFitModel::setNSigEvents(Double_t nSigEvents, Bool_t fixSigEvents)
124 {
125  signalEvents_.range(-2.0*(TMath::Abs(nSigEvents)+1.0),2.0*(TMath::Abs(nSigEvents)+1.0));
126  signalEvents_.value(nSigEvents);
127  signalEvents_.initValue(nSigEvents);
128  signalEvents_.genValue(nSigEvents);
129  signalEvents_.fixed(fixSigEvents);
130 
131  signalAsym_.range(-1.0,1.0);
132  signalAsym_.value(0.0);
133  signalAsym_.initValue(0.0);
134  signalAsym_.genValue(0.0);
135  signalAsym_.fixed(kTRUE);
136 }
137 
138 void LauCPFitModel::setNSigEvents( Double_t nSigEvents, Bool_t fixSigEvents, Double_t sigAsym, Bool_t fixSigAsym, Bool_t forceAsym )
139 {
140  signalEvents_.range(-2.0*(TMath::Abs(nSigEvents)+1.0),2.0*(TMath::Abs(nSigEvents)+1.0));
141  signalEvents_.value(nSigEvents);
142  signalEvents_.initValue(nSigEvents);
143  signalEvents_.genValue(nSigEvents);
144  signalEvents_.fixed(fixSigEvents);
145 
146  signalAsym_.range(-1.0,1.0);
147  signalAsym_.value(sigAsym);
148  signalAsym_.initValue(sigAsym);
149  signalAsym_.genValue(sigAsym);
150  signalAsym_.fixed(fixSigAsym);
151 
152  forceAsym_ = forceAsym;
153 }
154 
155 void LauCPFitModel::setNBkgndEvents( const TString& bkgndClass, Double_t nBkgndEvents, Bool_t fixBkgndEvents )
156 {
157  if ( ! this->validBkgndClass( bkgndClass ) ) {
158  cerr << "ERROR in LauCPFitModel::setNBkgndEvents : Invalid background class \"" << bkgndClass << "\"." << std::endl;
159  cerr << " : Background class names must be provided in \"setBkgndClassNames\" before any other background-related actions can be performed." << endl;
160  return;
161  }
162 
163  UInt_t bkgndID = this->bkgndClassID( bkgndClass );
164 
165  bkgndEvents_[bkgndID].name( bkgndClass+"Events" );
166  bkgndEvents_[bkgndID].range(-2.0*(TMath::Abs(nBkgndEvents)+1.0), 2.0*(TMath::Abs(nBkgndEvents)+1.0));
167  bkgndEvents_[bkgndID].value(nBkgndEvents);
168  bkgndEvents_[bkgndID].initValue(nBkgndEvents);
169  bkgndEvents_[bkgndID].genValue(nBkgndEvents);
170  bkgndEvents_[bkgndID].fixed(fixBkgndEvents);
171 
172  bkgndAsym_[bkgndID].name( bkgndClass+"Asym" );
173  bkgndAsym_[bkgndID].range(-1.0,1.0);
174  bkgndAsym_[bkgndID].value(0.0);
175  bkgndAsym_[bkgndID].initValue(0.0);
176  bkgndAsym_[bkgndID].genValue(0.0);
177  bkgndAsym_[bkgndID].fixed(kTRUE);
178 }
179 
180 void LauCPFitModel::setNBkgndEvents(const TString& bkgndClass, Double_t nBkgndEvents, Bool_t fixBkgndEvents, Double_t bkgndAsym, Bool_t fixBkgndAsym)
181 {
182  if ( ! this->validBkgndClass( bkgndClass ) ) {
183  cerr << "ERROR in LauCPFitModel::setNBkgndEvents : Invalid background class \"" << bkgndClass << "\"." << std::endl;
184  cerr << " : Background class names must be provided in \"setBkgndClassNames\" before any other background-related actions can be performed." << endl;
185  return;
186  }
187 
188  UInt_t bkgndID = this->bkgndClassID( bkgndClass );
189 
190  bkgndEvents_[bkgndID].name( bkgndClass+"Events" );
191  bkgndEvents_[bkgndID].range(-2.0*(TMath::Abs(nBkgndEvents)+1.0), 2.0*(TMath::Abs(nBkgndEvents)+1.0));
192  bkgndEvents_[bkgndID].value(nBkgndEvents);
193  bkgndEvents_[bkgndID].initValue(nBkgndEvents);
194  bkgndEvents_[bkgndID].genValue(nBkgndEvents);
195  bkgndEvents_[bkgndID].fixed(fixBkgndEvents);
196 
197  bkgndAsym_[bkgndID].name( bkgndClass+"Asym" );
198  bkgndAsym_[bkgndID].range(-1.0,1.0);
199  bkgndAsym_[bkgndID].value(bkgndAsym);
200  bkgndAsym_[bkgndID].initValue(bkgndAsym);
201  bkgndAsym_[bkgndID].genValue(bkgndAsym);
202  bkgndAsym_[bkgndID].fixed(fixBkgndAsym);
203 }
204 
205 void LauCPFitModel::splitSignalComponent( const TH2* dpHisto, Bool_t upperHalf, LauScfMap* scfMap )
206 {
207  if ( useSCF_ == kTRUE ) {
208  cerr << "ERROR in LauCPFitModel::splitSignalComponent : Have already setup SCF." << endl;
209  return;
210  }
211 
212  if ( dpHisto == 0 ) {
213  cerr << "ERROR in LauCPFitModel::splitSignalComponent : The histogram pointer is null." << endl;
214  return;
215  }
216 
217  LauDaughters* daughters = negSigModel_->getDaughters();
218  scfFracHist_ = new LauEffModel( daughters, 0 );
219  scfFracHist_->setEffHisto( dpHisto, kTRUE, kFALSE, 0.0, 0.0, upperHalf, daughters->squareDP() );
220 
221  scfMap_ = scfMap;
222 
223  useSCF_ = kTRUE;
224  useSCFHist_ = kTRUE;
225 }
226 
227 void LauCPFitModel::splitSignalComponent( Double_t scfFrac, Bool_t fixed )
228 {
229  if ( useSCF_ == kTRUE ) {
230  cerr << "ERROR in LauCPFitModel::splitSignalComponent : Have already setup SCF." << endl;
231  return;
232  }
233 
234  scfFrac_.range( 0.0, 1.0 );
235  scfFrac_.value( scfFrac ); scfFrac_.initValue( scfFrac ); scfFrac_.genValue( scfFrac );
236  scfFrac_.fixed( fixed );
237 
238  useSCF_ = kTRUE;
239  useSCFHist_ = kFALSE;
240 }
241 
242 void LauCPFitModel::setBkgndDPModels(const TString& bkgndClass, LauAbsBkgndDPModel* negModel, LauAbsBkgndDPModel* posModel)
243 {
244  if ((negModel==0) || (posModel==0)) {
245  cerr << "ERROR in LauCPFitModel::setBkgndDPModels : One or both of the model pointers is null." << endl;
246  return;
247  }
248 
249  // check that this background name is valid
250  if ( ! this->validBkgndClass( bkgndClass) ) {
251  cerr << "ERROR in LauCPFitModel::setBkgndDPModel : Invalid background class \"" << bkgndClass << "\"." << std::endl;
252  cerr << " : Background class names must be provided in \"setBkgndClassNames\" before any other background-related actions can be performed." << endl;
253  return;
254  }
255 
256  UInt_t bkgndID = this->bkgndClassID( bkgndClass );
257  negBkgndDPModels_[bkgndID] = negModel;
258  posBkgndDPModels_[bkgndID] = posModel;
259 
260  usingBkgnd_ = kTRUE;
261 }
262 
263 void LauCPFitModel::setSignalPdfs(LauAbsPdf* negPdf, LauAbsPdf* posPdf)
264 {
265  if ( tagged_ ) {
266  if (negPdf==0 || posPdf==0) {
267  cerr << "ERROR in LauCPFitModel::setSignalPdfs : One or both of the PDF pointers is null." << endl;
268  return;
269  }
270  } else {
271  // if we're doing an untagged analysis we will only use the negative PDFs
272  if ( negPdf==0 ) {
273  cerr << "ERROR in LauCPFitModel::setSignalPdfs : The negative PDF pointer is null." << endl;
274  return;
275  }
276  if ( posPdf!=0 ) {
277  cerr << "WARNING in LauCPFitModel::setSignalPdfs : Doing an untagged fit so will not use the positive PDF." << endl;
278  }
279  }
280  negSignalPdfs_.push_back(negPdf);
281  posSignalPdfs_.push_back(posPdf);
282 }
283 
284 void LauCPFitModel::setSCFPdfs(LauAbsPdf* negPdf, LauAbsPdf* posPdf)
285 {
286  if ( tagged_ ) {
287  if (negPdf==0 || posPdf==0) {
288  cerr << "ERROR in LauCPFitModel::setSCFPdfs : One or both of the PDF pointers is null." << endl;
289  return;
290  }
291  } else {
292  // if we're doing an untagged analysis we will only use the negative PDFs
293  if ( negPdf==0 ) {
294  cerr << "ERROR in LauCPFitModel::setSCFPdfs : The negative PDF pointer is null." << endl;
295  return;
296  }
297  if ( posPdf!=0 ) {
298  cerr << "WARNING in LauCPFitModel::setSCFPdfs : Doing an untagged fit so will not use the positive PDF." << endl;
299  }
300  }
301  negScfPdfs_.push_back(negPdf);
302  posScfPdfs_.push_back(posPdf);
303 }
304 
305 void LauCPFitModel::setBkgndPdfs(const TString& bkgndClass, LauAbsPdf* negPdf, LauAbsPdf* posPdf)
306 {
307  if ( tagged_ ) {
308  if (negPdf==0 || posPdf==0) {
309  cerr << "ERROR in LauCPFitModel::setBkgndPdfs : One or both of the PDF pointers is null." << endl;
310  return;
311  }
312  } else {
313  // if we're doing an untagged analysis we will only use the negative PDFs
314  if ( negPdf==0 ) {
315  cerr << "ERROR in LauCPFitModel::setBkgndPdfs : The negative PDF pointer is null." << endl;
316  return;
317  }
318  if ( posPdf!=0 ) {
319  cerr << "WARNING in LauCPFitModel::setBkgndPdfs : Doing an untagged fit so will not use the positive PDF." << endl;
320  }
321  }
322 
323  // check that this background name is valid
324  if ( ! this->validBkgndClass( bkgndClass ) ) {
325  cerr << "ERROR in LauCPFitModel::setBkgndPdfs : Invalid background class \"" << bkgndClass << "\"." << std::endl;
326  cerr << " : Background class names must be provided in \"setBkgndClassNames\" before any other background-related actions can be performed." << endl;
327  return;
328  }
329 
330  UInt_t bkgndID = this->bkgndClassID( bkgndClass );
331  negBkgndPdfs_[bkgndID].push_back(negPdf);
332  posBkgndPdfs_[bkgndID].push_back(posPdf);
333 
334  usingBkgnd_ = kTRUE;
335 }
336 
337 void LauCPFitModel::setAmpCoeffSet(LauAbsCoeffSet* coeffSet)
338 {
339  // Is there a component called compName in the signal model?
340  TString compName(coeffSet->name());
341  Bool_t negOK = negSigModel_->hasResonance(compName);
342  TString conjName = negSigModel_->getConjResName(compName);
343  Bool_t posOK = posSigModel_->hasResonance(conjName);
344  if (!negOK) {
345  cout << "ERROR in LauCPFitModel::setMagPhase : " << negParent_ << " signal DP model doesn't contain component \"" << compName << "\"." << endl;
346  return;
347  }
348  if (!posOK) {
349  cout << "ERROR in LauCPFitModel::setMagPhase : " << posParent_ << " signal DP model doesn't contain component \"" << conjName << "\"." << endl;
350  return;
351  }
352 
353  // Do we already have it in our list of names?
354  for (std::vector<LauAbsCoeffSet*>::const_iterator iter=coeffPars_.begin(); iter!=coeffPars_.end(); ++iter) {
355  if ((*iter)->name() == compName) {
356  cerr << "ERROR in LauCPFitModel::setAmpCoeffSet : Have already set coefficients for \"" << compName << "\"." << endl;
357  return;
358  }
359  }
360 
361  coeffSet->index(nSigComp_);
362  coeffPars_.push_back(coeffSet);
363 
364  TString parName = coeffSet->baseName(); parName += "FitFracAsym";
365  fitFracAsymm_.push_back(LauParameter(parName, 0.0, -1.0, 1.0));
366 
367  acp_.push_back(coeffSet->acp());
368 
369  ++nSigComp_;
370 
371  cout << "Added coefficients for component \"" << compName << "\" to the fit model." << endl;
372 }
373 
374 void LauCPFitModel::initialise()
375 {
376  // First of all check that, if a given component is being used,
377  // we've got the PDFs for all the variables involved
378  if ( this->useDP() ) {
379  if ((negSigModel_ == 0) || (posSigModel_ == 0)) {
380  cerr << "ERROR in LauCPFitModel::initialise : the pointer to one (neg or pos) of the signal DP models is null.\n";
381  cerr << " : Removing the Dalitz Plot from the model." << endl;
382  this->useDP(kFALSE);
383  }
384  if ( usingBkgnd_ ) {
385  if ( negBkgndDPModels_.empty() || posBkgndDPModels_.empty() ) {
386  cerr << "ERROR in LauCPFitModel::initialise : No background DP models found.\n";
387  cerr << " : Removing the Dalitz plot from the model." << endl;
388  this->useDP(kFALSE);
389  }
390  for (LauBkgndDPModelList::const_iterator dpmodel_iter = negBkgndDPModels_.begin(); dpmodel_iter != negBkgndDPModels_.end(); ++dpmodel_iter ) {
391  if ( (*dpmodel_iter) == 0 ) {
392  cerr << "ERROR in LauCPFitModel::initialise : The pointer to one of the background DP models is null.\n";
393  cerr << " : Removing the Dalitz Plot from the model." << endl;
394  this->useDP(kFALSE);
395  break;
396  }
397  }
398  for (LauBkgndDPModelList::const_iterator dpmodel_iter = posBkgndDPModels_.begin(); dpmodel_iter != posBkgndDPModels_.end(); ++dpmodel_iter ) {
399  if ( (*dpmodel_iter) == 0 ) {
400  cerr << "ERROR in LauCPFitModel::initialise : The pointer to one of the background DP models is null.\n";
401  cerr << " : Removing the Dalitz Plot from the model." << endl;
402  this->useDP(kFALSE);
403  break;
404  }
405  }
406  }
407  }
408 
409  // Next check that, if a given component is being used we've got the
410  // right number of PDFs for all the variables involved
411  // TODO - should probably check variable names and so on as well
412 
413  UInt_t nsigpdfvars(0);
414  for ( LauPdfList::const_iterator pdf_iter = negSignalPdfs_.begin(); pdf_iter != negSignalPdfs_.end(); ++pdf_iter ) {
415  std::vector<TString> varNames = (*pdf_iter)->varNames();
416  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
417  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
418  ++nsigpdfvars;
419  }
420  }
421  }
422  if (useSCF_) {
423  UInt_t nscfpdfvars(0);
424  for ( LauPdfList::const_iterator pdf_iter = negScfPdfs_.begin(); pdf_iter != negScfPdfs_.end(); ++pdf_iter ) {
425  std::vector<TString> varNames = (*pdf_iter)->varNames();
426  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
427  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
428  ++nscfpdfvars;
429  }
430  }
431  }
432  if (nscfpdfvars != nsigpdfvars) {
433  cerr << "ERROR in LauCPFitModel::initialise : There are " << nsigpdfvars << " TM signal PDF variables but " << nscfpdfvars << " SCF signal PDF variables." << endl;
434  gSystem->Exit(EXIT_FAILURE);
435  }
436  }
437  if (usingBkgnd_) {
438  for (LauBkgndPdfsList::const_iterator bgclass_iter = negBkgndPdfs_.begin(); bgclass_iter != negBkgndPdfs_.end(); ++bgclass_iter) {
439  UInt_t nbkgndpdfvars(0);
440  const LauPdfList& pdfList = (*bgclass_iter);
441  for ( LauPdfList::const_iterator pdf_iter = pdfList.begin(); pdf_iter != pdfList.end(); ++pdf_iter ) {
442  std::vector<TString> varNames = (*pdf_iter)->varNames();
443  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
444  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
445  ++nbkgndpdfvars;
446  }
447  }
448  }
449  if (nbkgndpdfvars != nsigpdfvars) {
450  cerr << "ERROR in LauCPFitModel::initialise : There are " << nsigpdfvars << " signal PDF variables but " << nbkgndpdfvars << " bkgnd PDF variables." << endl;
451  gSystem->Exit(EXIT_FAILURE);
452  }
453  }
454  }
455 
456  // Clear the vectors of parameter information so we can start from scratch
457  this->clearFitParVectors();
458 
459  // Set the fit parameters for signal and background models
460  this->setSignalDPParameters();
461 
462  // Set the fit parameters for the various extra PDFs
463  this->setExtraPdfParameters();
464 
465  // Set the initial bg and signal events
466  this->setFitNEvents();
467 
468  // Check that we have the expected number of fit variables
469  const LauParameterPList& fitVars = this->fitPars();
470  if (fitVars.size() != (nSigDPPar_ + nExtraPdfPar_ + nNormPar_)) {
471  cerr << "ERROR in LauCPFitModel::initialise : Number of fit parameters not of expected size. Exiting" << endl;
472  gSystem->Exit(EXIT_FAILURE);
473  }
474 
475  // From the initial parameter values calculate the coefficients
476  // so they can be passed to the signal model
477  this->updateCoeffs();
478 
479  // Initialisation
480  if (this->useDP() == kTRUE) {
481  this->initialiseDPModels();
482  }
483 
484  if (!this->useDP() && negSignalPdfs_.empty()) {
485  cerr << "ERROR in LauCPFitModel::initialise : Signal model doesn't exist for any variable." << endl;
486  gSystem->Exit(EXIT_FAILURE);
487  }
488 
489  this->setExtraNtupleVars();
490 }
491 
492 void LauCPFitModel::initialiseDPModels()
493 {
494  cout << "INFO in LauCPFitModel::initialiseDPModels : Initialising signal DP model" << endl;
495  negSigModel_->initialise(negCoeffs_);
496  posSigModel_->initialise(posCoeffs_);
497 
498  if (usingBkgnd_ == kTRUE) {
499  for (LauBkgndDPModelList::iterator iter = negBkgndDPModels_.begin(); iter != negBkgndDPModels_.end(); ++iter) {
500  (*iter)->initialise();
501  }
502  for (LauBkgndDPModelList::iterator iter = posBkgndDPModels_.begin(); iter != posBkgndDPModels_.end(); ++iter) {
503  (*iter)->initialise();
504  }
505  }
506 }
507 
508 void LauCPFitModel::setSignalDPParameters()
509 {
510  // Set the fit parameters for the signal model.
511 
512  nSigDPPar_ = 0;
513  if ( ! this->useDP() ) {
514  return;
515  }
516 
517  cout << "INFO in LauCPFitModel::setSignalDPParameters : Setting the initial fit parameters for the signal DP model." << endl;
518 
519  // Need to check that the number of components we have and that the dynamics has matches up
520  UInt_t nNegAmp = negSigModel_->getnAmp();
521  UInt_t nPosAmp = posSigModel_->getnAmp();
522  if (nNegAmp != nSigComp_ || nPosAmp != nSigComp_) {
523  cerr << "ERROR in LauCPFitModel::setSignalDPParameters : Number of signal DP components with magnitude and phase set not right." << endl;
524  gSystem->Exit(EXIT_FAILURE);
525  }
526 
527 
528  // Place signal model parameters in vector of fit variables
529  LauParameterPList& fitVars = this->fitPars();
530  for (UInt_t i = 0; i < nSigComp_; i++) {
531  LauParameterPList pars = coeffPars_[i]->getParameters();
532  for (LauParameterPList::iterator iter = pars.begin(); iter != pars.end(); ++iter) {
533  if ( !(*iter)->clone() ) {
534  fitVars.push_back(*iter);
535  ++nSigDPPar_;
536  }
537  }
538  }
539 }
540 
541 void LauCPFitModel::setExtraPdfParameters()
542 {
543  // Include all the parameters of the PDF in the fit
544  // NB all of them are passed to the fit, even though some have been fixed through parameter.fixed(kTRUE)
545  // With the new "cloned parameter" scheme only "original" parameters are passed to the fit.
546  // Their clones are updated automatically when the originals are updated.
547 
548  nExtraPdfPar_ = 0;
549 
550  nExtraPdfPar_ += this->addFitParameters(negSignalPdfs_);
551  if ( tagged_ ) {
552  nExtraPdfPar_ += this->addFitParameters(posSignalPdfs_);
553  }
554 
555  if (useSCF_ == kTRUE) {
556  nExtraPdfPar_ += this->addFitParameters(negScfPdfs_);
557  if ( tagged_ ) {
558  nExtraPdfPar_ += this->addFitParameters(posScfPdfs_);
559  }
560  }
561 
562  if (usingBkgnd_ == kTRUE) {
563  for (LauBkgndPdfsList::iterator iter = negBkgndPdfs_.begin(); iter != negBkgndPdfs_.end(); ++iter) {
564  nExtraPdfPar_ += this->addFitParameters(*iter);
565  }
566  if ( tagged_ ) {
567  for (LauBkgndPdfsList::iterator iter = posBkgndPdfs_.begin(); iter != posBkgndPdfs_.end(); ++iter) {
568  nExtraPdfPar_ += this->addFitParameters(*iter);
569  }
570  }
571  }
572 }
573 
574 void LauCPFitModel::setFitNEvents()
575 {
576  nNormPar_ = 0;
577 
578  // initialise the total number of events to be the sum of all the hypotheses
579  Double_t nTotEvts = signalEvents_.value();
580  for (LauBkgndYieldList::const_iterator iter = bkgndEvents_.begin(); iter != bkgndEvents_.end(); ++iter) {
581  nTotEvts += *iter;
582  }
583  this->eventsPerExpt(TMath::FloorNint(nTotEvts));
584 
585  LauParameterPList& fitVars = this->fitPars();
586 
587  // if doing an extended ML fit add the number of signal events into the fit parameters
588  if (this->doEMLFit()) {
589  cout << "INFO in LauCPFitModel::setFitNEvents : Initialising number of events for signal and background components..." << endl;
590  // add the signal fraction to the list of fit parameters
591  fitVars.push_back(&signalEvents_);
592  ++nNormPar_;
593  } else {
594  cout << "INFO in LauCPFitModel::setFitNEvents : Initialising number of events for background components (and hence signal)..." << endl;
595  }
596 
597  // if not using the DP in the model we need an explicit signal asymmetry parameter
598  if (this->useDP() == kFALSE) {
599  fitVars.push_back(&signalAsym_);
600  ++nNormPar_;
601  }
602 
603  if (useSCF_ && !useSCFHist_) {
604  fitVars.push_back(&scfFrac_);
605  ++nNormPar_;
606  }
607 
608  if (usingBkgnd_ == kTRUE) {
609  for (LauBkgndYieldList::iterator iter = bkgndEvents_.begin(); iter != bkgndEvents_.end(); ++iter) {
610  LauParameter& parameter = (*iter);
611  fitVars.push_back(&parameter);
612  ++nNormPar_;
613  }
614  for (LauBkgndYieldList::iterator iter = bkgndAsym_.begin(); iter != bkgndAsym_.end(); ++iter) {
615  LauParameter& parameter = (*iter);
616  fitVars.push_back(&parameter);
617  ++nNormPar_;
618  }
619  }
620 }
621 
622 void LauCPFitModel::setExtraNtupleVars()
623 {
624  // Set-up other parameters derived from the fit results, e.g. fit fractions.
625 
626  if (this->useDP() != kTRUE) {
627  return;
628  }
629 
630  // First clear the vectors so we start from scratch
631  this->clearExtraVarVectors();
632 
633  // Add the positive and negative fit fractions for each signal component
634  negFitFrac_ = negSigModel_->getFitFractions();
635  if (negFitFrac_.size() != nSigComp_) {
636  cerr << "ERROR in LauCPFitModel::setExtraNtupleVars : Initial Fit Fraction array of unexpected dimension: " << negFitFrac_.size() << endl;
637  gSystem->Exit(EXIT_FAILURE);
638  }
639  for (UInt_t i(0); i<nSigComp_; ++i) {
640  if (negFitFrac_[i].size() != nSigComp_) {
641  cerr << "ERROR in LauCPFitModel::setExtraNtupleVars : Initial Fit Fraction array of unexpected dimension: " << negFitFrac_[i].size() << endl;
642  gSystem->Exit(EXIT_FAILURE);
643  }
644  }
645 
646  LauParameterList& extraVars = this->extraPars();
647 
648  for (UInt_t i(0); i<nSigComp_; ++i) {
649  for (UInt_t j(i); j<nSigComp_; ++j) {
650  TString name = negFitFrac_[i][j].name();
651  name.Insert( name.Index("FitFrac"), "Neg" );
652  negFitFrac_[i][j].name(name);
653  extraVars.push_back(negFitFrac_[i][j]);
654  }
655  }
656 
657  posFitFrac_ = posSigModel_->getFitFractions();
658  if (posFitFrac_.size() != nSigComp_) {
659  cerr << "ERROR in LauCPFitModel::setExtraNtupleVars : Initial Fit Fraction array of unexpected dimension: " << posFitFrac_.size() << endl;
660  gSystem->Exit(EXIT_FAILURE);
661  }
662  for (UInt_t i(0); i<nSigComp_; ++i) {
663  if (posFitFrac_[i].size() != nSigComp_) {
664  cerr << "ERROR in LauCPFitModel::setExtraNtupleVars : Initial Fit Fraction array of unexpected dimension: " << posFitFrac_[i].size() << endl;
665  gSystem->Exit(EXIT_FAILURE);
666  }
667  }
668  for (UInt_t i(0); i<nSigComp_; ++i) {
669  for (UInt_t j(i); j<nSigComp_; ++j) {
670  TString name = posFitFrac_[i][j].name();
671  name.Insert( name.Index("FitFrac"), "Pos" );
672  posFitFrac_[i][j].name(name);
673  extraVars.push_back(posFitFrac_[i][j]);
674  }
675  }
676 
677  // Store any extra parameters/quantities from the DP model (e.g. K-matrix total fit fractions)
678  std::vector<LauParameter> negExtraPars = negSigModel_->getExtraParameters();
679  std::vector<LauParameter>::iterator negExtraIter;
680  for (negExtraIter = negExtraPars.begin(); negExtraIter != negExtraPars.end(); ++negExtraIter) {
681  LauParameter negExtraParameter = (*negExtraIter);
682  extraVars.push_back(negExtraParameter);
683  }
684 
685  std::vector<LauParameter> posExtraPars = posSigModel_->getExtraParameters();
686  std::vector<LauParameter>::iterator posExtraIter;
687  for (posExtraIter = posExtraPars.begin(); posExtraIter != posExtraPars.end(); ++posExtraIter) {
688  LauParameter posExtraParameter = (*posExtraIter);
689  extraVars.push_back(posExtraParameter);
690  }
691 
692  // Now add in the DP efficiency value
693  Double_t initMeanEff = negSigModel_->getMeanEff().initValue();
694  negMeanEff_.value(initMeanEff);
695  negMeanEff_.genValue(initMeanEff);
696  negMeanEff_.initValue(initMeanEff);
697  extraVars.push_back(negMeanEff_);
698 
699  initMeanEff = posSigModel_->getMeanEff().initValue();
700  posMeanEff_.value(initMeanEff);
701  posMeanEff_.genValue(initMeanEff);
702  posMeanEff_.initValue(initMeanEff);
703  extraVars.push_back(posMeanEff_);
704 
705  // Also add in the DP rates
706  Double_t initDPRate = negSigModel_->getDPRate().initValue();
707  negDPRate_.value(initDPRate);
708  negDPRate_.genValue(initDPRate);
709  negDPRate_.initValue(initDPRate);
710  extraVars.push_back(negDPRate_);
711 
712  initDPRate = posSigModel_->getDPRate().initValue();
713  posDPRate_.value(initDPRate);
714  posDPRate_.genValue(initDPRate);
715  posDPRate_.initValue(initDPRate);
716  extraVars.push_back(posDPRate_);
717 
718  // Calculate the CPC and CPV Fit Fractions, ACPs and FitFrac asymmetries
719  this->calcExtraFractions(kTRUE);
720  this->calcAsymmetries(kTRUE);
721 
722  // Add the CP violating and CP conserving fit fractions for each signal component
723  for (UInt_t i = 0; i < nSigComp_; i++) {
724  for (UInt_t j = i; j < nSigComp_; j++) {
725  extraVars.push_back(CPVFitFrac_[i][j]);
726  }
727  }
728  for (UInt_t i = 0; i < nSigComp_; i++) {
729  for (UInt_t j = i; j < nSigComp_; j++) {
730  extraVars.push_back(CPCFitFrac_[i][j]);
731  }
732  }
733 
734  // Add the Fit Fraction asymmetry for each signal component
735  for (UInt_t i = 0; i < nSigComp_; i++) {
736  extraVars.push_back(fitFracAsymm_[i]);
737  }
738 
739  // Add the calculated CP asymmetry for each signal component
740  for (UInt_t i = 0; i < nSigComp_; i++) {
741  extraVars.push_back(acp_[i]);
742  }
743 }
744 
745 void LauCPFitModel::calcExtraFractions(Bool_t initValues)
746 {
747  // Calculate the CP-conserving and CP-violating fit fractions
748 
749  if (initValues) {
750  // create the structure
751  CPCFitFrac_.clear();
752  CPVFitFrac_.clear();
753  CPCFitFrac_.resize(nSigComp_);
754  CPVFitFrac_.resize(nSigComp_);
755  for (UInt_t i(0); i<nSigComp_; ++i) {
756  CPCFitFrac_[i].resize(nSigComp_);
757  CPVFitFrac_[i].resize(nSigComp_);
758 
759  for (UInt_t j(i); j<nSigComp_; ++j) {
760  TString name = negFitFrac_[i][j].name();
761  name.Replace( name.Index("Neg"), 3, "CPC" );
762  CPCFitFrac_[i][j].name( name );
763  CPCFitFrac_[i][j].valueAndRange( 0.0, -100.0, 100.0 );
764 
765  name = negFitFrac_[i][j].name();
766  name.Replace( name.Index("Neg"), 3, "CPV" );
767  CPVFitFrac_[i][j].name( name );
768  CPVFitFrac_[i][j].valueAndRange( 0.0, -100.0, 100.0 );
769  }
770  }
771  }
772 
773  Double_t denom = negDPRate_ + posDPRate_;
774 
775  for (UInt_t i(0); i<nSigComp_; ++i) {
776  for (UInt_t j(i); j<nSigComp_; ++j) {
777 
778  Double_t negTerm = negFitFrac_[i][j]*negDPRate_;
779  Double_t posTerm = posFitFrac_[i][j]*posDPRate_;
780 
781  Double_t cpcFitFrac = (negTerm + posTerm)/denom;
782  Double_t cpvFitFrac = (negTerm - posTerm)/denom;
783 
784  CPCFitFrac_[i][j] = cpcFitFrac;
785  CPVFitFrac_[i][j] = cpvFitFrac;
786 
787  if (initValues) {
788  CPCFitFrac_[i][j].genValue(cpcFitFrac);
789  CPCFitFrac_[i][j].initValue(cpcFitFrac);
790 
791  CPVFitFrac_[i][j].genValue(cpvFitFrac);
792  CPVFitFrac_[i][j].initValue(cpvFitFrac);
793  }
794  }
795  }
796 }
797 
798 void LauCPFitModel::calcAsymmetries(Bool_t initValues)
799 {
800  // Calculate the CP asymmetries
801  // Also calculate the fit fraction asymmetries
802 
803  for (UInt_t i = 0; i < nSigComp_; i++) {
804 
805  acp_[i] = coeffPars_[i]->acp();
806 
807  LauAsymmCalc asymmCalc(negFitFrac_[i][i].value(), posFitFrac_[i][i].value());
808  Double_t asym = asymmCalc.getAsymmetry();
809  fitFracAsymm_[i] = asym;
810  if (initValues) {
811  fitFracAsymm_[i].genValue(asym);
812  fitFracAsymm_[i].initValue(asym);
813  }
814  }
815 }
816 
817 void LauCPFitModel::finaliseFitResults(const TString& tablePrefixName)
818 {
819  // Retrieve parameters from the fit results for calculations and toy generation
820  // and eventually store these in output root ntuples/text files
821 
822  // Now take the fit parameters and update them as necessary
823  // i.e. to make mag > 0.0, phase in the right range.
824  // This function will also calculate any other values, such as the
825  // fit fractions, using any errors provided by fitParErrors as appropriate.
826  // Also obtain the pull values: (measured - generated)/(average error)
827 
828  if (this->useDP() == kTRUE) {
829  for (UInt_t i = 0; i < nSigComp_; ++i) {
830  // Check whether we have "a/b > 0.0", and phases in the right range
831  coeffPars_[i]->finaliseValues();
832  }
833  }
834 
835  // update the pulls on the event fractions and asymmetries
836  if (this->doEMLFit()) {
837  signalEvents_.updatePull();
838  }
839  if (this->useDP() == kFALSE) {
840  signalAsym_.updatePull();
841  }
842  if (useSCF_ && !useSCFHist_) {
843  scfFrac_.updatePull();
844  }
845  if (usingBkgnd_ == kTRUE) {
846  for (LauBkgndYieldList::iterator iter = bkgndEvents_.begin(); iter != bkgndEvents_.end(); ++iter) {
847  iter->updatePull();
848  }
849  for (LauBkgndYieldList::iterator iter = bkgndAsym_.begin(); iter != bkgndAsym_.end(); ++iter) {
850  iter->updatePull();
851  }
852  }
853 
854  // Update the pulls on all the extra PDFs' parameters
855  this->updateFitParameters(negSignalPdfs_);
856  this->updateFitParameters(posSignalPdfs_);
857  if (useSCF_ == kTRUE) {
858  this->updateFitParameters(negScfPdfs_);
859  this->updateFitParameters(posScfPdfs_);
860  }
861  if (usingBkgnd_ == kTRUE) {
862  for (LauBkgndPdfsList::iterator iter = negBkgndPdfs_.begin(); iter != negBkgndPdfs_.end(); ++iter) {
863  this->updateFitParameters(*iter);
864  }
865  for (LauBkgndPdfsList::iterator iter = posBkgndPdfs_.begin(); iter != posBkgndPdfs_.end(); ++iter) {
866  this->updateFitParameters(*iter);
867  }
868  }
869 
870  // Fill the fit results to the ntuple
871 
872  // update the coefficients and then calculate the fit fractions and ACP's
873  if (this->useDP() == kTRUE) {
874  this->updateCoeffs();
875  negSigModel_->updateCoeffs(negCoeffs_); negSigModel_->calcExtraInfo();
876  posSigModel_->updateCoeffs(posCoeffs_); posSigModel_->calcExtraInfo();
877 
878  LauParArray negFitFrac = negSigModel_->getFitFractions();
879  if (negFitFrac.size() != nSigComp_) {
880  cerr << "ERROR in LauCPFitModel::finaliseFitResults : Fit Fraction array of unexpected dimension: " << negFitFrac.size() << endl;
881  gSystem->Exit(EXIT_FAILURE);
882  }
883  for (UInt_t i(0); i<nSigComp_; ++i) {
884  if (negFitFrac[i].size() != nSigComp_) {
885  cerr << "ERROR in LauCPFitModel::finaliseFitResults : Fit Fraction array of unexpected dimension: " << negFitFrac[i].size() << endl;
886  gSystem->Exit(EXIT_FAILURE);
887  }
888  }
889 
890  LauParArray posFitFrac = posSigModel_->getFitFractions();
891  if (posFitFrac.size() != nSigComp_) {
892  cerr << "ERROR in LauCPFitModel::finaliseFitResults : Fit Fraction array of unexpected dimension: " << posFitFrac.size() << endl;
893  gSystem->Exit(EXIT_FAILURE);
894  }
895  for (UInt_t i(0); i<nSigComp_; ++i) {
896  if (posFitFrac[i].size() != nSigComp_) {
897  cerr << "ERROR in LauCPFitModel::finaliseFitResults : Fit Fraction array of unexpected dimension: " << posFitFrac[i].size() << endl;
898  gSystem->Exit(EXIT_FAILURE);
899  }
900  }
901 
902  for (UInt_t i(0); i<nSigComp_; ++i) {
903  for (UInt_t j(i); j<nSigComp_; ++j) {
904  negFitFrac_[i][j].value(negFitFrac[i][j].value());
905  posFitFrac_[i][j].value(posFitFrac[i][j].value());
906  }
907  }
908 
909  negMeanEff_.value(negSigModel_->getMeanEff().value());
910  posMeanEff_.value(posSigModel_->getMeanEff().value());
911  negDPRate_.value(negSigModel_->getDPRate().value());
912  posDPRate_.value(posSigModel_->getDPRate().value());
913 
914  this->calcExtraFractions();
915  this->calcAsymmetries();
916 
917  // Then store the final fit parameters, and any extra parameters for
918  // the signal model (e.g. fit fractions, FF asymmetries, ACPs, mean efficiency and DP rate)
919 
920  this->clearExtraVarVectors();
921  LauParameterList& extraVars = this->extraPars();
922 
923  // Add the positive and negative fit fractions for each signal component
924  for (UInt_t i(0); i<nSigComp_; ++i) {
925  for (UInt_t j(i); j<nSigComp_; ++j) {
926  extraVars.push_back(negFitFrac_[i][j]);
927  }
928  }
929  for (UInt_t i(0); i<nSigComp_; ++i) {
930  for (UInt_t j(i); j<nSigComp_; ++j) {
931  extraVars.push_back(posFitFrac_[i][j]);
932  }
933  }
934 
935  // Store any extra parameters/quantities from the DP model (e.g. K-matrix total fit fractions)
936  std::vector<LauParameter> negExtraPars = negSigModel_->getExtraParameters();
937  std::vector<LauParameter>::iterator negExtraIter;
938  for (negExtraIter = negExtraPars.begin(); negExtraIter != negExtraPars.end(); ++negExtraIter) {
939  LauParameter negExtraParameter = (*negExtraIter);
940  extraVars.push_back(negExtraParameter);
941  }
942 
943  std::vector<LauParameter> posExtraPars = posSigModel_->getExtraParameters();
944  std::vector<LauParameter>::iterator posExtraIter;
945  for (posExtraIter = posExtraPars.begin(); posExtraIter != posExtraPars.end(); ++posExtraIter) {
946  LauParameter posExtraParameter = (*posExtraIter);
947  extraVars.push_back(posExtraParameter);
948  }
949 
950  extraVars.push_back(negMeanEff_);
951  extraVars.push_back(posMeanEff_);
952  extraVars.push_back(negDPRate_);
953  extraVars.push_back(posDPRate_);
954 
955  for (UInt_t i = 0; i < nSigComp_; i++) {
956  for (UInt_t j(i); j<nSigComp_; ++j) {
957  extraVars.push_back(CPVFitFrac_[i][j]);
958  }
959  }
960  for (UInt_t i = 0; i < nSigComp_; i++) {
961  for (UInt_t j(i); j<nSigComp_; ++j) {
962  extraVars.push_back(CPCFitFrac_[i][j]);
963  }
964  }
965  for (UInt_t i(0); i<nSigComp_; ++i) {
966  extraVars.push_back(fitFracAsymm_[i]);
967  }
968  for (UInt_t i(0); i<nSigComp_; ++i) {
969  extraVars.push_back(acp_[i]);
970  }
971 
972  this->printFitFractions(cout);
973  this->printAsymmetries(cout);
974  }
975 
976  const LauParameterPList& fitVars = this->fitPars();
977  const LauParameterList& extraVars = this->extraPars();
978  LauFitNtuple* ntuple = this->fitNtuple();
979  ntuple->storeParsAndErrors(fitVars, extraVars);
980 
981  // find out the correlation matrix for the parameters
982  ntuple->storeCorrMatrix(this->iExpt(), this->nll(), this->fitStatus());
983 
984  // Fill the data into ntuple
985  ntuple->updateFitNtuple();
986 
987  // Print out the partial fit fractions, phases and the
988  // averaged efficiency, reweighted by the dynamics (and anything else)
989  if (this->writeLatexTable()) {
990  TString sigOutFileName(tablePrefixName);
991  sigOutFileName += "_"; sigOutFileName += this->iExpt(); sigOutFileName += "Expt.tex";
992  this->writeOutTable(sigOutFileName);
993  }
994 }
995 
996 void LauCPFitModel::printFitFractions(std::ostream& output)
997 {
998  // Print out Fit Fractions, total DP rate and mean efficiency
999  // First for the B- events
1000  for (UInt_t i = 0; i < nSigComp_; i++) {
1001  const TString compName(coeffPars_[i]->name());
1002  output << negParent_ << " FitFraction for component " << i << " (" << compName << ") = " << negFitFrac_[i][i] << endl;
1003  }
1004  output << negParent_ << " overall DP rate (integral of matrix element squared) = " << negDPRate_ << endl;
1005  output << negParent_ << " average efficiency weighted by whole DP dynamics = " << negMeanEff_ << endl;
1006 
1007  // Then for the positive sample
1008  for (UInt_t i = 0; i < nSigComp_; i++) {
1009  const TString compName(coeffPars_[i]->name());
1010  const TString conjName(negSigModel_->getConjResName(compName));
1011  output << posParent_ << " FitFraction for component " << i << " (" << conjName << ") = " << posFitFrac_[i][i] << endl;
1012  }
1013  output << posParent_ << " overall DP rate (integral of matrix element squared) = " << posDPRate_ << endl;
1014  output << posParent_ << " average efficiency weighted by whole DP dynamics = " << posMeanEff_ << endl;
1015 }
1016 
1017 void LauCPFitModel::printAsymmetries(std::ostream& output)
1018 {
1019  for (UInt_t i = 0; i < nSigComp_; i++) {
1020  const TString compName(coeffPars_[i]->name());
1021  output << "Fit Fraction asymmetry for component " << i << " (" << compName << ") = " << fitFracAsymm_[i] << endl;
1022  }
1023  for (UInt_t i = 0; i < nSigComp_; i++) {
1024  const TString compName(coeffPars_[i]->name());
1025  output << "ACP for component " << i << " (" << compName << ") = " << acp_[i] << endl;
1026  }
1027 }
1028 
1029 void LauCPFitModel::writeOutTable(const TString& outputFile)
1030 {
1031  // Write out the results of the fit to a tex-readable table
1032  // TODO - need to include the yields in this table
1033  std::ofstream fout(outputFile);
1034  LauPrint print;
1035 
1036  cout << "Writing out results of the fit to the tex file " << outputFile << endl;
1037 
1038  if (this->useDP() == kTRUE) {
1039  // print the fit coefficients in one table
1040  coeffPars_.front()->printTableHeading(fout);
1041  for (UInt_t i = 0; i < nSigComp_; i++) {
1042  coeffPars_[i]->printTableRow(fout);
1043  }
1044  fout << "\\hline" << endl;
1045  fout << "\\end{tabular}" << endl << endl;
1046 
1047  // print the fit fractions and asymmetries in another
1048  fout << "\\begin{tabular}{|l|c|c|c|c|}" << endl;
1049  fout << "\\hline" << endl;
1050  fout << "Component & " << negParent_ << " Fit Fraction & " << posParent_ << " Fit Fraction & Fit Fraction Asymmetry & ACP \\\\" << endl;
1051  fout << "\\hline" << endl;
1052 
1053  Double_t negFitFracSum(0.0);
1054  Double_t posFitFracSum(0.0);
1055 
1056  for (UInt_t i = 0; i < nSigComp_; i++) {
1057  TString resName = coeffPars_[i]->name();
1058  resName = resName.ReplaceAll("_", "\\_");
1059 
1060  Double_t negFitFrac = negFitFrac_[i][i].value();
1061  Double_t posFitFrac = posFitFrac_[i][i].value();
1062  negFitFracSum += negFitFrac;
1063  posFitFracSum += posFitFrac;
1064 
1065  Double_t fitFracAsymm = fitFracAsymm_[i].value();
1066 
1067  Double_t acp = acp_[i].value();
1068  Double_t acpErr = acp_[i].error();
1069 
1070  fout << resName << " & $";
1071 
1072  print.printFormat(fout, negFitFrac);
1073  fout << "$ & $";
1074  print.printFormat(fout, posFitFrac);
1075  fout << "$ & $";
1076 
1077  print.printFormat(fout, fitFracAsymm);
1078  fout << "$ & $";
1079 
1080  print.printFormat(fout, acp);
1081  fout << " \\pm ";
1082  print.printFormat(fout, acpErr);
1083  fout << "$ \\\\" << endl;
1084  }
1085  fout << "\\hline" << endl;
1086 
1087  // Also print out sum of fit fractions
1088  fout << "Fit Fraction Sum & $";
1089  print.printFormat(fout, negFitFracSum);
1090  fout << "$ & $";
1091  print.printFormat(fout, posFitFracSum);
1092  fout << "$ & & \\\\" << endl;
1093  fout << "\\hline" << endl;
1094 
1095  fout << "DP rate & $";
1096  print.printFormat(fout, negDPRate_.value());
1097  fout << "$ & $";
1098  print.printFormat(fout, posDPRate_.value());
1099  fout << "$ & & \\\\" << endl;
1100 
1101  fout << "$< \\varepsilon > $ & $";
1102  print.printFormat(fout, negMeanEff_.value());
1103  fout << "$ & $";
1104  print.printFormat(fout, posMeanEff_.value());
1105  fout << "$ & & \\\\" << endl;
1106  fout << "\\hline" << endl;
1107  fout << "\\end{tabular}" << endl << endl;
1108  }
1109 
1110  if (!negSignalPdfs_.empty()) {
1111  fout << "\\begin{tabular}{|l|c|}" << endl;
1112  fout << "\\hline" << endl;
1113  if (useSCF_ == kTRUE) {
1114  fout << "\\Extra TM Signal PDFs' Parameters: & \\\\" << endl;
1115  } else {
1116  fout << "\\Extra Signal PDFs' Parameters: & \\\\" << endl;
1117  }
1118  this->printFitParameters(negSignalPdfs_, fout);
1119  if ( tagged_ ) {
1120  this->printFitParameters(posSignalPdfs_, fout);
1121  }
1122  if (useSCF_ == kTRUE && !negScfPdfs_.empty()) {
1123  fout << "\\hline" << endl;
1124  fout << "\\Extra SCF Signal PDFs' Parameters: & \\\\" << endl;
1125  this->printFitParameters(negScfPdfs_, fout);
1126  if ( tagged_ ) {
1127  this->printFitParameters(posScfPdfs_, fout);
1128  }
1129  }
1130  if (usingBkgnd_ == kTRUE && !negBkgndPdfs_.empty()) {
1131  fout << "\\hline" << endl;
1132  fout << "\\Extra Background PDFs' Parameters: & \\\\" << endl;
1133  for (LauBkgndPdfsList::const_iterator iter = negBkgndPdfs_.begin(); iter != negBkgndPdfs_.end(); ++iter) {
1134  this->printFitParameters(*iter, fout);
1135  }
1136  if ( tagged_ ) {
1137  for (LauBkgndPdfsList::const_iterator iter = posBkgndPdfs_.begin(); iter != posBkgndPdfs_.end(); ++iter) {
1138  this->printFitParameters(*iter, fout);
1139  }
1140  }
1141  }
1142  fout << "\\hline \n\\end{tabular}" << endl << endl;
1143  }
1144 }
1145 
1146 void LauCPFitModel::checkInitFitParams()
1147 {
1148  // Update the number of signal events to be total-sum(background events)
1149  this->updateSigEvents();
1150 
1151  // Check whether we want to have randomised initial fit parameters for the signal model
1152  if (this->useRandomInitFitPars() == kTRUE) {
1153  cout << "Setting random parameters for the signal model" << endl;
1154  this->randomiseInitFitPars();
1155  }
1156 }
1157 
1158 void LauCPFitModel::randomiseInitFitPars()
1159 {
1160  // Only randomise those parameters that are not fixed!
1161  cout << "Randomising the initial fit magnitudes and phases of the components..." << endl;
1162 
1163  for (UInt_t i = 0; i < nSigComp_; i++) {
1164  coeffPars_[i]->randomiseInitValues();
1165  }
1166 }
1167 
1168 LauCPFitModel::LauGenInfo LauCPFitModel::eventsToGenerate()
1169 {
1170  // Determine the number of events to generate for each hypothesis
1171  // If we're smearing then smear each one individually
1172 
1173  LauGenInfo nEvtsGen;
1174 
1175  // Signal
1176  Double_t evtWeight(1.0);
1177  Double_t nEvts = signalEvents_.genValue();
1178  if ( nEvts < 0.0 ) {
1179  evtWeight = -1.0;
1180  nEvts = TMath::Abs( nEvts );
1181  }
1182  Double_t asym(0.0);
1183  Double_t sigAsym(0.0);
1184  // need to include this as an alternative in case the DP isn't in the model
1185  if ( !this->useDP() || forceAsym_ ) {
1186  sigAsym = signalAsym_.genValue();
1187  } else {
1188  Double_t negRate = negSigModel_->getDPNorm();
1189  Double_t posRate = posSigModel_->getDPNorm();
1190  if (negRate+posRate>1e-30) {
1191  sigAsym = (negRate-posRate)/(negRate+posRate);
1192  }
1193  }
1194  asym = sigAsym;
1195  Int_t nPosEvts = static_cast<Int_t>((nEvts/2.0 * (1.0 - asym)) + 0.5);
1196  Int_t nNegEvts = static_cast<Int_t>((nEvts/2.0 * (1.0 + asym)) + 0.5);
1197  if (this->doPoissonSmearing()) {
1198  nNegEvts = LauRandom::randomFun()->Poisson(nNegEvts);
1199  nPosEvts = LauRandom::randomFun()->Poisson(nPosEvts);
1200  }
1201  nEvtsGen[std::make_pair("signal",-1)] = std::make_pair(nNegEvts,evtWeight);
1202  nEvtsGen[std::make_pair("signal",+1)] = std::make_pair(nPosEvts,evtWeight);
1203 
1204  // backgrounds
1205  const UInt_t nBkgnds = this->nBkgndClasses();
1206  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
1207  const TString& bkgndClass = this->bkgndClassName(bkgndID);
1208  const LauParameter& evtsPar = bkgndEvents_[bkgndID];
1209  const LauParameter& asymPar = bkgndAsym_[bkgndID];
1210  evtWeight = 1.0;
1211  nEvts = TMath::FloorNint( evtsPar.genValue() );
1212  if ( nEvts < 0 ) {
1213  evtWeight = -1.0;
1214  nEvts = TMath::Abs( nEvts );
1215  }
1216  asym = asymPar.genValue();
1217  nPosEvts = static_cast<Int_t>((nEvts/2.0 * (1.0 - asym)) + 0.5);
1218  nNegEvts = static_cast<Int_t>((nEvts/2.0 * (1.0 + asym)) + 0.5);
1219  if (this->doPoissonSmearing()) {
1220  nNegEvts = LauRandom::randomFun()->Poisson(nNegEvts);
1221  nPosEvts = LauRandom::randomFun()->Poisson(nPosEvts);
1222  }
1223  nEvtsGen[std::make_pair(bkgndClass,-1)] = std::make_pair(nNegEvts,evtWeight);
1224  nEvtsGen[std::make_pair(bkgndClass,+1)] = std::make_pair(nPosEvts,evtWeight);
1225  }
1226 
1227  cout << "Generating toy MC with:" << endl;
1228  cout << "Signal asymmetry = " << sigAsym << " and number of signal events = " << signalEvents_.genValue() << endl;
1229  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
1230  const TString& bkgndClass = this->bkgndClassName(bkgndID);
1231  const LauParameter& evtsPar = bkgndEvents_[bkgndID];
1232  const LauParameter& asymPar = bkgndAsym_[bkgndID];
1233  cout << bkgndClass << " asymmetry = " << asymPar.genValue() << " and number of " << bkgndClass << " events = " << evtsPar.genValue() << endl;
1234  }
1235 
1236  return nEvtsGen;
1237 }
1238 
1239 Bool_t LauCPFitModel::genExpt()
1240 {
1241  // Routine to generate toy Monte Carlo events according to the various models we have defined.
1242 
1243  // Determine the number of events to generate for each hypothesis
1244  LauGenInfo nEvts = this->eventsToGenerate();
1245 
1246  Bool_t genOK(kTRUE);
1247  Int_t evtNum(0);
1248 
1249  const UInt_t nBkgnds = this->nBkgndClasses();
1250  std::vector<TString> bkgndClassNames(nBkgnds);
1251  std::vector<TString> bkgndClassNamesGen(nBkgnds);
1252  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
1253  TString name( this->bkgndClassName(iBkgnd) );
1254  bkgndClassNames[iBkgnd] = name;
1255  bkgndClassNamesGen[iBkgnd] = "gen"+name;
1256  }
1257 
1258  const Bool_t storeSCFTruthInfo = ( useSCF_ || ( this->enableEmbedding() &&
1259  negSignalTree_ != 0 && negSignalTree_->haveBranch("mcMatch") &&
1260  posSignalTree_ != 0 && posSignalTree_->haveBranch("mcMatch") ) );
1261 
1262  // Loop over the hypotheses and generate the requested number of events for each one
1263  for (LauGenInfo::const_iterator iter = nEvts.begin(); iter != nEvts.end(); ++iter) {
1264 
1265  const TString& type(iter->first.first);
1266  curEvtCharge_ = iter->first.second;
1267  Double_t evtWeight( iter->second.second );
1268  Int_t nEvtsGen( iter->second.first );
1269 
1270  for (Int_t iEvt(0); iEvt<nEvtsGen; ++iEvt) {
1271 
1272  this->setGenNtupleDoubleBranchValue( "evtWeight", evtWeight );
1273 
1274  if (type == "signal") {
1275  this->setGenNtupleIntegerBranchValue("genSig",1);
1276  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
1277  this->setGenNtupleIntegerBranchValue( bkgndClassNamesGen[iBkgnd], 0 );
1278  }
1279  genOK = this->generateSignalEvent();
1280  } else {
1281  this->setGenNtupleIntegerBranchValue("genSig",0);
1282  if ( storeSCFTruthInfo ) {
1283  this->setGenNtupleIntegerBranchValue("genTMSig",0);
1284  this->setGenNtupleIntegerBranchValue("genSCFSig",0);
1285  }
1286  UInt_t bkgndID(0);
1287  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
1288  Int_t gen(0);
1289  if ( bkgndClassNames[iBkgnd] == type ) {
1290  gen = 1;
1291  bkgndID = iBkgnd;
1292  }
1293  this->setGenNtupleIntegerBranchValue( bkgndClassNamesGen[iBkgnd], gen );
1294  }
1295  genOK = this->generateBkgndEvent(bkgndID);
1296  }
1297  if (!genOK) {
1298  // If there was a problem with the generation then break out and return.
1299  // The problem model will have adjusted itself so that all should be OK next time.
1300  break;
1301  }
1302  if (this->useDP() == kTRUE) {
1303  this->setDPBranchValues();
1304  }
1305 
1306  // Store the event charge
1307  this->setGenNtupleIntegerBranchValue(tagVarName_,curEvtCharge_);
1308 
1309  // Store the event number (within this experiment)
1310  // and then increment it
1311  this->setGenNtupleIntegerBranchValue("iEvtWithinExpt",evtNum);
1312  ++evtNum;
1313 
1314  this->fillGenNtupleBranches();
1315  if (iEvt%500 == 0) {cout << "Generated event number " << iEvt << " out of " << nEvtsGen << " " << type << " events." << endl;}
1316  }
1317 
1318  if (!genOK) {
1319  break;
1320  }
1321  }
1322 
1323  if (this->useDP() && genOK) {
1324 
1325  negSigModel_->checkToyMC(kTRUE,kTRUE);
1326  posSigModel_->checkToyMC(kTRUE,kTRUE);
1327 
1328  // Get the fit fractions if they're to be written into the latex table
1329  if (this->writeLatexTable()) {
1330  LauParArray negFitFrac = negSigModel_->getFitFractions();
1331  if (negFitFrac.size() != nSigComp_) {
1332  cerr << "ERROR in LauCPFitModel::genExpt : Fit Fraction array of unexpected dimension: " << negFitFrac.size() << endl;
1333  gSystem->Exit(EXIT_FAILURE);
1334  }
1335  for (UInt_t i(0); i<nSigComp_; ++i) {
1336  if (negFitFrac[i].size() != nSigComp_) {
1337  cerr << "ERROR in LauCPFitModel::genExpt : Fit Fraction array of unexpected dimension: " << negFitFrac[i].size() << endl;
1338  gSystem->Exit(EXIT_FAILURE);
1339  }
1340  }
1341  LauParArray posFitFrac = posSigModel_->getFitFractions();
1342  if (posFitFrac.size() != nSigComp_) {
1343  cerr << "ERROR in LauCPFitModel::genExpt : Fit Fraction array of unexpected dimension: " << posFitFrac.size() << endl;
1344  gSystem->Exit(EXIT_FAILURE);
1345  }
1346  for (UInt_t i(0); i<nSigComp_; ++i) {
1347  if (posFitFrac[i].size() != nSigComp_) {
1348  cerr << "ERROR in LauCPFitModel::genExpt : Fit Fraction array of unexpected dimension: " << posFitFrac[i].size() << endl;
1349  gSystem->Exit(EXIT_FAILURE);
1350  }
1351  }
1352 
1353  for (UInt_t i(0); i<nSigComp_; ++i) {
1354  for (UInt_t j(i); j<nSigComp_; ++j) {
1355  negFitFrac_[i][j].value(negFitFrac[i][j].value());
1356  posFitFrac_[i][j].value(posFitFrac[i][j].value());
1357  }
1358  }
1359  negMeanEff_.value(negSigModel_->getMeanEff().value());
1360  posMeanEff_.value(posSigModel_->getMeanEff().value());
1361  negDPRate_.value(negSigModel_->getDPRate().value());
1362  posDPRate_.value(posSigModel_->getDPRate().value());
1363  }
1364  }
1365 
1366  // If we're reusing embedded events or if the generation is being
1367  // reset then clear the lists of used events
1368  if (reuseSignal_ || !genOK) {
1369  if (negSignalTree_) {
1370  negSignalTree_->clearUsedList();
1371  }
1372  if (posSignalTree_) {
1373  posSignalTree_->clearUsedList();
1374  }
1375  }
1376  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
1377  LauEmbeddedData* data = negBkgndTree_[bkgndID];
1378  if (reuseBkgnd_[bkgndID] || !genOK) {
1379  if (data) {
1380  data->clearUsedList();
1381  }
1382  }
1383  }
1384  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
1385  LauEmbeddedData* data = posBkgndTree_[bkgndID];
1386  if (reuseBkgnd_[bkgndID] || !genOK) {
1387  if (data) {
1388  data->clearUsedList();
1389  }
1390  }
1391  }
1392 
1393  return genOK;
1394 }
1395 
1396 Bool_t LauCPFitModel::generateSignalEvent()
1397 {
1398  // Generate signal event
1399  Bool_t genOK(kTRUE);
1400  Bool_t genSCF(kFALSE);
1401 
1402  LauAbsDPDynamics* model(0);
1403  LauKinematics* kinematics(0);
1404  LauEmbeddedData* embeddedData(0);
1405  LauPdfList* sigPdfs(0);
1406  LauPdfList* scfPdfs(0);
1407 
1408  Bool_t doReweighting(kFALSE);
1409 
1410  if (curEvtCharge_<0) {
1411  model = negSigModel_;
1412  kinematics = negKinematics_;
1413  sigPdfs = &negSignalPdfs_;
1414  scfPdfs = &negScfPdfs_;
1415  if (this->enableEmbedding()) {
1416  embeddedData = negSignalTree_;
1417  doReweighting = useNegReweighting_;
1418  }
1419  } else {
1420  model = posSigModel_;
1421  kinematics = posKinematics_;
1422  if ( tagged_ ) {
1423  sigPdfs = &posSignalPdfs_;
1424  scfPdfs = &posScfPdfs_;
1425  } else {
1426  sigPdfs = &negSignalPdfs_;
1427  scfPdfs = &negScfPdfs_;
1428  }
1429  if (this->enableEmbedding()) {
1430  embeddedData = posSignalTree_;
1431  doReweighting = usePosReweighting_;
1432  }
1433  }
1434 
1435  if (this->useDP()) {
1436  if (embeddedData) {
1437  if (doReweighting) {
1438  // Select a (random) event from the generated data. Then store the
1439  // reconstructed DP co-ords, together with other pdf information,
1440  // as the embedded data.
1441  genOK = embeddedData->getReweightedEvent(model);
1442  } else {
1443  // Just get the information of a (randomly) selected event in the
1444  // embedded data
1445  embeddedData->getEmbeddedEvent(kinematics);
1446  }
1447  genSCF = this->storeSignalMCMatch( embeddedData );
1448  } else {
1449  genOK = model->generate();
1450  if ( genOK && useSCF_ ) {
1451  Double_t frac(0.0);
1452  if ( useSCFHist_ ) {
1453  frac = scfFracHist_->calcEfficiency( kinematics );
1454  } else {
1455  frac = scfFrac_.genValue();
1456  }
1457  if ( frac < LauRandom::randomFun()->Rndm() ) {
1458  this->setGenNtupleIntegerBranchValue("genTMSig",1);
1459  this->setGenNtupleIntegerBranchValue("genSCFSig",0);
1460  genSCF = kFALSE;
1461  } else {
1462  this->setGenNtupleIntegerBranchValue("genTMSig",0);
1463  this->setGenNtupleIntegerBranchValue("genSCFSig",1);
1464  genSCF = kTRUE;
1465 
1466  // Optionally smear the DP position
1467  // of the SCF event
1468  if ( scfMap_ != 0 ) {
1469  Double_t xCoord(0.0), yCoord(0.0);
1470  if ( kinematics->squareDP() ) {
1471  xCoord = kinematics->getmPrime();
1472  yCoord = kinematics->getThetaPrime();
1473  } else {
1474  xCoord = kinematics->getm13Sq();
1475  yCoord = kinematics->getm23Sq();
1476  }
1477 
1478  // Find the bin number where this event is generated
1479  Int_t binNo = scfMap_->binNumber( xCoord, yCoord );
1480 
1481  // Retrieve the migration histogram
1482  TH2* histo = scfMap_->trueHist( binNo );
1483 
1484  LauEffModel * effModel = model->getEffModel();
1485  do {
1486  // Get a random point from the histogram
1487  histo->GetRandom2( xCoord, yCoord );
1488 
1489  // Update the kinematics
1490  if ( kinematics->squareDP() ) {
1491  kinematics->updateSqDPKinematics( xCoord, yCoord );
1492  } else {
1493  kinematics->updateKinematics( xCoord, yCoord );
1494  }
1495  } while ( ! effModel->passVeto( kinematics ) );
1496  }
1497  }
1498  }
1499  }
1500  } else {
1501  if (embeddedData) {
1502  embeddedData->getEmbeddedEvent(0);
1503  genSCF = this->storeSignalMCMatch( embeddedData );
1504  } else if ( useSCF_ ) {
1505  Double_t frac = scfFrac_.genValue();
1506  if ( frac < LauRandom::randomFun()->Rndm() ) {
1507  this->setGenNtupleIntegerBranchValue("genTMSig",1);
1508  this->setGenNtupleIntegerBranchValue("genSCFSig",0);
1509  genSCF = kFALSE;
1510  } else {
1511  this->setGenNtupleIntegerBranchValue("genTMSig",0);
1512  this->setGenNtupleIntegerBranchValue("genSCFSig",1);
1513  genSCF = kTRUE;
1514  }
1515  }
1516  }
1517  if (genOK) {
1518  if ( useSCF_ ) {
1519  if ( genSCF ) {
1520  this->generateExtraPdfValues(scfPdfs, embeddedData);
1521  } else {
1522  this->generateExtraPdfValues(sigPdfs, embeddedData);
1523  }
1524  } else {
1525  this->generateExtraPdfValues(sigPdfs, embeddedData);
1526  }
1527  }
1528  // Check for problems with the embedding
1529  if (embeddedData && (embeddedData->nEvents() == embeddedData->nUsedEvents())) {
1530  cerr << "WARNING in LauCPFitModel::generateSignalEvent : Source of embedded signal events used up, clearing the list of used events." << endl;
1531  embeddedData->clearUsedList();
1532  }
1533 
1534  return genOK;
1535 }
1536 
1537 Bool_t LauCPFitModel::generateBkgndEvent(UInt_t bkgndID)
1538 {
1539  // Generate Bkgnd event
1540  Bool_t genOK(kTRUE);
1541 
1542  LauAbsBkgndDPModel* model(0);
1543  LauEmbeddedData* embeddedData(0);
1544  LauPdfList* extraPdfs(0);
1545  LauKinematics* kinematics(0);
1546 
1547  if (curEvtCharge_<0) {
1548  model = negBkgndDPModels_[bkgndID];
1549  if (this->enableEmbedding()) {
1550  embeddedData = negBkgndTree_[bkgndID];
1551  }
1552  extraPdfs = &negBkgndPdfs_[bkgndID];
1553  kinematics = negKinematics_;
1554  } else {
1555  model = posBkgndDPModels_[bkgndID];
1556  if (this->enableEmbedding()) {
1557  embeddedData = posBkgndTree_[bkgndID];
1558  }
1559  if ( tagged_ ) {
1560  extraPdfs = &posBkgndPdfs_[bkgndID];
1561  } else {
1562  extraPdfs = &negBkgndPdfs_[bkgndID];
1563  }
1564  kinematics = posKinematics_;
1565  }
1566 
1567  if (this->useDP()) {
1568  if (embeddedData) {
1569  embeddedData->getEmbeddedEvent(kinematics);
1570  } else {
1571  if (model == 0) {
1572  const TString& bkgndClass = this->bkgndClassName(bkgndID);
1573  cerr << "ERROR in LauCPFitModel::generateBkgndEvent : Can't find the DP model for background class \"" << bkgndClass << "\"." << endl;
1574  gSystem->Exit(EXIT_FAILURE);
1575  }
1576  genOK = model->generate();
1577  }
1578  } else {
1579  if (embeddedData) {
1580  embeddedData->getEmbeddedEvent(0);
1581  }
1582  }
1583  if (genOK) {
1584  this->generateExtraPdfValues(extraPdfs, embeddedData);
1585  }
1586 
1587  // Check for problems with the embedding
1588  if (embeddedData && (embeddedData->nEvents() == embeddedData->nUsedEvents())) {
1589  const TString& bkgndClass = this->bkgndClassName(bkgndID);
1590  cerr << "WARNING in LauCPFitModel::generateBkgndEvent : Source of embedded " << bkgndClass << " events used up, clearing the list of used events." << endl;
1591  embeddedData->clearUsedList();
1592  }
1593 
1594  return genOK;
1595 }
1596 
1597 void LauCPFitModel::setupGenNtupleBranches()
1598 {
1599  // Setup the required ntuple branches
1600  this->addGenNtupleDoubleBranch("evtWeight");
1601  this->addGenNtupleIntegerBranch("genSig");
1602  if ( useSCF_ || ( this->enableEmbedding() &&
1603  negSignalTree_ != 0 && negSignalTree_->haveBranch("mcMatch") &&
1604  posSignalTree_ != 0 && posSignalTree_->haveBranch("mcMatch") ) ) {
1605  this->addGenNtupleIntegerBranch("genTMSig");
1606  this->addGenNtupleIntegerBranch("genSCFSig");
1607  }
1608  const UInt_t nBkgnds = this->nBkgndClasses();
1609  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
1610  TString name( this->bkgndClassName(iBkgnd) );
1611  name.Prepend("gen");
1612  this->addGenNtupleIntegerBranch(name);
1613  }
1614  this->addGenNtupleIntegerBranch("charge");
1615  if (this->useDP() == kTRUE) {
1616  this->addGenNtupleDoubleBranch("m12");
1617  this->addGenNtupleDoubleBranch("m23");
1618  this->addGenNtupleDoubleBranch("m13");
1619  this->addGenNtupleDoubleBranch("m12Sq");
1620  this->addGenNtupleDoubleBranch("m23Sq");
1621  this->addGenNtupleDoubleBranch("m13Sq");
1622  this->addGenNtupleDoubleBranch("cosHel12");
1623  this->addGenNtupleDoubleBranch("cosHel23");
1624  this->addGenNtupleDoubleBranch("cosHel13");
1625  if (negKinematics_->squareDP() && posKinematics_->squareDP()) {
1626  this->addGenNtupleDoubleBranch("mPrime");
1627  this->addGenNtupleDoubleBranch("thPrime");
1628  }
1629  }
1630  for (LauPdfList::const_iterator pdf_iter = negSignalPdfs_.begin(); pdf_iter != negSignalPdfs_.end(); ++pdf_iter) {
1631  std::vector<TString> varNames = (*pdf_iter)->varNames();
1632  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
1633  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
1634  this->addGenNtupleDoubleBranch( (*var_iter) );
1635  }
1636  }
1637  }
1638 }
1639 
1640 void LauCPFitModel::setDPBranchValues()
1641 {
1642  LauKinematics* kinematics(0);
1643  if (curEvtCharge_<0) {
1644  kinematics = negKinematics_;
1645  } else {
1646  kinematics = posKinematics_;
1647  }
1648 
1649  // Store all the DP information
1650  this->setGenNtupleDoubleBranchValue("m12", kinematics->getm12());
1651  this->setGenNtupleDoubleBranchValue("m23", kinematics->getm23());
1652  this->setGenNtupleDoubleBranchValue("m13", kinematics->getm13());
1653  this->setGenNtupleDoubleBranchValue("m12Sq", kinematics->getm12Sq());
1654  this->setGenNtupleDoubleBranchValue("m23Sq", kinematics->getm23Sq());
1655  this->setGenNtupleDoubleBranchValue("m13Sq", kinematics->getm13Sq());
1656  this->setGenNtupleDoubleBranchValue("cosHel12", kinematics->getc12());
1657  this->setGenNtupleDoubleBranchValue("cosHel23", kinematics->getc23());
1658  this->setGenNtupleDoubleBranchValue("cosHel13", kinematics->getc13());
1659  if (kinematics->squareDP()) {
1660  this->setGenNtupleDoubleBranchValue("mPrime", kinematics->getmPrime());
1661  this->setGenNtupleDoubleBranchValue("thPrime", kinematics->getThetaPrime());
1662  }
1663 }
1664 
1665 void LauCPFitModel::generateExtraPdfValues(LauPdfList* extraPdfs, LauEmbeddedData* embeddedData)
1666 {
1667  LauKinematics* kinematics(0);
1668  if (curEvtCharge_<0) {
1669  kinematics = negKinematics_;
1670  } else {
1671  kinematics = posKinematics_;
1672  }
1673 
1674  if (!extraPdfs) {
1675  std::cerr << "ERROR in LauCPFitModel::generateExtraPdfValues : Null pointer to PDF list." << std::endl;
1676  gSystem->Exit(EXIT_FAILURE);
1677  }
1678 
1679  if (extraPdfs->empty()) {
1680  //std::cerr << "WARNING in LauCPFitModel::generateExtraPdfValues : PDF list is empty." << std::endl;
1681  return;
1682  }
1683 
1684  // Generate from the extra PDFs
1685  for (LauPdfList::iterator pdf_iter = extraPdfs->begin(); pdf_iter != extraPdfs->end(); ++pdf_iter) {
1686  LauFitData genValues;
1687  if (embeddedData) {
1688  genValues = embeddedData->getValues( (*pdf_iter)->varNames() );
1689  } else {
1690  genValues = (*pdf_iter)->generate(kinematics);
1691  }
1692  for ( LauFitData::const_iterator var_iter = genValues.begin(); var_iter != genValues.end(); ++var_iter ) {
1693  TString varName = var_iter->first;
1694  if ( varName != "m13Sq" && varName != "m23Sq" ) {
1695  Double_t value = var_iter->second;
1696  this->setGenNtupleDoubleBranchValue(varName,value);
1697  }
1698  }
1699  }
1700 }
1701 
1702 Bool_t LauCPFitModel::storeSignalMCMatch(LauEmbeddedData* embeddedData)
1703 {
1704  // Default to TM
1705  Bool_t genSCF(kFALSE);
1706  Int_t match(1);
1707 
1708  // Check that we have a valid pointer and that embedded data has
1709  // the mcMatch branch. If so then get the match value.
1710  if ( embeddedData && embeddedData->haveBranch("mcMatch") ) {
1711  match = TMath::Nint( embeddedData->getValue("mcMatch") );
1712  }
1713 
1714  // Set the variables accordingly.
1715  if (match) {
1716  this->setGenNtupleIntegerBranchValue("genTMSig",1);
1717  this->setGenNtupleIntegerBranchValue("genSCFSig",0);
1718  genSCF = kFALSE;
1719  } else {
1720  this->setGenNtupleIntegerBranchValue("genTMSig",0);
1721  this->setGenNtupleIntegerBranchValue("genSCFSig",1);
1722  genSCF = kTRUE;
1723  }
1724 
1725  return genSCF;
1726 }
1727 
1728 void LauCPFitModel::propagateParUpdates()
1729 {
1730  // Update the signal parameters and then the total normalisation for the signal likelihood
1731  if (this->useDP() == kTRUE) {
1732  this->updateCoeffs();
1733  negSigModel_->updateCoeffs(negCoeffs_);
1734  posSigModel_->updateCoeffs(posCoeffs_);
1735  }
1736 
1737  // Update the signal fraction from the background fractions if not doing an extended fit
1738  if ( !this->doEMLFit() && !signalEvents_.fixed() ) {
1739  this->updateSigEvents();
1740  }
1741 }
1742 
1743 void LauCPFitModel::updateSigEvents()
1744 {
1745  // The background parameters will have been set from Minuit.
1746  // We need to update the signal events using these.
1747  Double_t nTotEvts = this->eventsPerExpt();
1748 
1749  signalEvents_.range(-2.0*nTotEvts,2.0*nTotEvts);
1750  for (LauBkgndYieldList::iterator iter = bkgndEvents_.begin(); iter != bkgndEvents_.end(); ++iter) {
1751  (*iter).range(-2.0*nTotEvts,2.0*nTotEvts);
1752  }
1753 
1754  if (signalEvents_.fixed()) {
1755  return;
1756  }
1757 
1758  // Subtract background events (if any) from signal.
1759  Double_t signalEvents = nTotEvts;
1760  if (usingBkgnd_ == kTRUE) {
1761  for (LauBkgndYieldList::const_iterator iter = bkgndEvents_.begin(); iter != bkgndEvents_.end(); ++iter) {
1762  signalEvents -= (*iter);
1763  }
1764  }
1765  signalEvents_.value(signalEvents);
1766 }
1767 
1768 void LauCPFitModel::cacheInputFitVars()
1769 {
1770  // Fill the internal data trees of the signal and background models.
1771  // Note that we store the events of both charges in both the
1772  // negative and the positive models. It's only later, at the stage
1773  // when the likelihood is being calculated, that we separate them.
1774 
1775  LauFitDataTree* inputFitData = this->fitData();
1776 
1777  // First the Dalitz plot variables (m_ij^2)
1778  if (this->useDP() == kTRUE) {
1779 
1780  // need to append SCF smearing bins before caching DP amplitudes
1781  if ( scfMap_ != 0 ) {
1782  this->appendBinCentres( inputFitData );
1783  }
1784 
1785  negSigModel_->fillDataTree(*inputFitData);
1786  posSigModel_->fillDataTree(*inputFitData);
1787 
1788  if (usingBkgnd_ == kTRUE) {
1789  for (LauBkgndDPModelList::iterator iter = negBkgndDPModels_.begin(); iter != negBkgndDPModels_.end(); ++iter) {
1790  (*iter)->fillDataTree(*inputFitData);
1791  }
1792  for (LauBkgndDPModelList::iterator iter = posBkgndDPModels_.begin(); iter != posBkgndDPModels_.end(); ++iter) {
1793  (*iter)->fillDataTree(*inputFitData);
1794  }
1795  }
1796  }
1797 
1798  // ...and then the extra PDFs
1799  this->cacheInfo(negSignalPdfs_, *inputFitData);
1800  this->cacheInfo(negScfPdfs_, *inputFitData);
1801  for (LauBkgndPdfsList::iterator iter = negBkgndPdfs_.begin(); iter != negBkgndPdfs_.end(); ++iter) {
1802  this->cacheInfo((*iter), *inputFitData);
1803  }
1804  if ( tagged_ ) {
1805  this->cacheInfo(posSignalPdfs_, *inputFitData);
1806  this->cacheInfo(posScfPdfs_, *inputFitData);
1807  for (LauBkgndPdfsList::iterator iter = posBkgndPdfs_.begin(); iter != posBkgndPdfs_.end(); ++iter) {
1808  this->cacheInfo((*iter), *inputFitData);
1809  }
1810  }
1811 
1812  // the SCF fractions and jacobians
1813  if ( useSCF_ && useSCFHist_ ) {
1814  if ( !inputFitData->haveBranch( "m13Sq" ) || !inputFitData->haveBranch( "m23Sq" ) ) {
1815  cerr << "ERROR in LauCPFitModel::cacheInputFitVars : Input data does not contain DP branches and so can't cache the SCF fraction." << endl;
1816  gSystem->Exit(EXIT_FAILURE);
1817  }
1818 
1819  UInt_t nEvents = inputFitData->nEvents();
1820  recoSCFFracs_.clear();
1821  recoSCFFracs_.reserve( nEvents );
1822  if ( negKinematics_->squareDP() ) {
1823  recoJacobians_.clear();
1824  recoJacobians_.reserve( nEvents );
1825  }
1826  for (UInt_t iEvt = 0; iEvt < nEvents; iEvt++) {
1827  const LauFitData& dataValues = inputFitData->getData(iEvt);
1828  LauFitData::const_iterator m13_iter = dataValues.find("m13Sq");
1829  LauFitData::const_iterator m23_iter = dataValues.find("m23Sq");
1830  negKinematics_->updateKinematics( m13_iter->second, m23_iter->second );
1831  Double_t scfFrac = scfFracHist_->calcEfficiency( negKinematics_ );
1832  recoSCFFracs_.push_back( scfFrac );
1833  if ( negKinematics_->squareDP() ) {
1834  recoJacobians_.push_back( negKinematics_->calcSqDPJacobian() );
1835  }
1836  }
1837  }
1838 
1839  // finally cache the event charge
1840  evtCharges_.clear();
1841  if ( tagged_ ) {
1842  if ( !inputFitData->haveBranch( tagVarName_ ) ) {
1843  cerr << "ERROR in LauCPFitModel::cacheInputFitVars : Input data does not contain branch \"" << tagVarName_ << "\"." << endl;
1844  gSystem->Exit(EXIT_FAILURE);
1845  }
1846  UInt_t nEvents = inputFitData->nEvents();
1847  evtCharges_.reserve( nEvents );
1848  for (UInt_t iEvt = 0; iEvt < nEvents; iEvt++) {
1849  const LauFitData& dataValues = inputFitData->getData(iEvt);
1850  LauFitData::const_iterator iter = dataValues.find( tagVarName_ );
1851  curEvtCharge_ = static_cast<Int_t>( iter->second );
1852  evtCharges_.push_back( curEvtCharge_ );
1853  }
1854  }
1855 }
1856 
1857 void LauCPFitModel::appendBinCentres( LauFitDataTree* inputData )
1858 {
1859  // We'll be caching the DP amplitudes and efficiencies of the centres of the true bins.
1860  // To do so, we attach some fake points at the end of inputData, the number of the entry
1861  // minus the total number of events corresponding to the number of the histogram for that
1862  // given true bin in the LauScfMap object. (What this means is that when Laura is provided with
1863  // the LauScfMap object by the user, it's the latter who has to make sure that it contains the
1864  // right number of histograms and in exactly the right order!)
1865 
1866  // Get the x and y co-ordinates of the bin centres
1867  std::vector<Double_t> binCentresXCoords;
1868  std::vector<Double_t> binCentresYCoords;
1869  scfMap_->listBinCentres(binCentresXCoords, binCentresYCoords);
1870 
1871  // The SCF histograms could be in square Dalitz plot histograms.
1872  // The dynamics takes normal Dalitz plot coords, so we might have to convert them back.
1873  Bool_t sqDP = negKinematics_->squareDP();
1874  UInt_t nBins = binCentresXCoords.size();
1875  fakeSCFFracs_.clear();
1876  fakeSCFFracs_.reserve( nBins );
1877  if ( sqDP ) {
1878  fakeJacobians_.clear();
1879  fakeJacobians_.reserve( nBins );
1880  }
1881 
1882  for (UInt_t iBin = 0; iBin < nBins; ++iBin) {
1883 
1884  if ( sqDP ) {
1885  negKinematics_->updateSqDPKinematics(binCentresXCoords[iBin],binCentresYCoords[iBin]);
1886  binCentresXCoords[iBin] = negKinematics_->getm13Sq();
1887  binCentresYCoords[iBin] = negKinematics_->getm23Sq();
1888  fakeJacobians_.push_back( negKinematics_->calcSqDPJacobian() );
1889  } else {
1890  negKinematics_->updateKinematics(binCentresXCoords[iBin],binCentresYCoords[iBin]);
1891  }
1892 
1893  fakeSCFFracs_.push_back( scfFracHist_->calcEfficiency( negKinematics_ ) );
1894  }
1895 
1896  // Set up inputFitVars_ object to hold the fake events
1897  inputData->appendFakePoints(binCentresXCoords,binCentresYCoords);
1898 }
1899 
1900 Double_t LauCPFitModel::getTotEvtLikelihood(UInt_t iEvt)
1901 {
1902  // Find out whether we have B- or B+
1903  if ( tagged_ ) {
1904  curEvtCharge_ = evtCharges_[iEvt];
1905 
1906  // check that the charge is either +1 or -1
1907  if (TMath::Abs(curEvtCharge_)!=1) {
1908  cerr << "ERROR in LauCPFitModel::getTotEvtLikelihood : Charge/tag not accepted value: " << curEvtCharge_ << endl;
1909  if (curEvtCharge_>0) {
1910  curEvtCharge_ = +1;
1911  } else {
1912  curEvtCharge_ = -1;
1913  }
1914  cerr << " : Making it: " << curEvtCharge_ << "." << endl;
1915  }
1916  }
1917 
1918  // Get the DP likelihood for signal and backgrounds
1919  this->getEvtDPLikelihood(iEvt);
1920 
1921  // Get the combined extra PDFs likelihood for signal and backgrounds
1922  this->getEvtExtraLikelihoods(iEvt);
1923 
1924  // If appropriate, combine the TM and SCF likelihoods
1925  Double_t sigLike = sigDPLike_ * sigExtraLike_;
1926  if ( useSCF_ ) {
1927  Double_t scfFrac(0.0);
1928  if (useSCFHist_) {
1929  scfFrac = recoSCFFracs_[iEvt];
1930  } else {
1931  scfFrac = scfFrac_.value();
1932  }
1933  sigLike *= (1.0 - scfFrac);
1934  if ( (scfMap_ != 0) && (this->useDP() == kTRUE) ) {
1935  // if we're smearing the SCF DP PDF then the SCF frac
1936  // is already included in the SCF DP likelihood
1937  sigLike += (scfDPLike_ * scfExtraLike_);
1938  } else {
1939  sigLike += (scfFrac * scfDPLike_ * scfExtraLike_);
1940  }
1941  }
1942 
1943  // Get the correct event fractions depending on the charge
1944  // Signal asymmetry is built into the DP model... but when the DP
1945  // isn't in the fit we need an explicit parameter
1946  Double_t signalEvents = signalEvents_ * 0.5;
1947  if (this->useDP() == kFALSE) {
1948  signalEvents *= (1.0 - curEvtCharge_ * signalAsym_);
1949  }
1950 
1951  // Construct the total event likelihood
1952  Double_t likelihood(0.0);
1953  if (usingBkgnd_) {
1954  likelihood = sigLike*signalEvents;
1955  const UInt_t nBkgnds = this->nBkgndClasses();
1956  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
1957  Double_t bkgndEvents = bkgndEvents_[bkgndID] * 0.5 * (1.0 - curEvtCharge_ * bkgndAsym_[bkgndID]);
1958  likelihood += bkgndEvents*bkgndDPLike_[bkgndID]*bkgndExtraLike_[bkgndID];
1959  }
1960  } else {
1961  likelihood = sigLike*0.5;
1962  }
1963  return likelihood;
1964 }
1965 
1966 Double_t LauCPFitModel::getEventSum() const
1967 {
1968  Double_t eventSum(0.0);
1969  eventSum += signalEvents_;
1970  if (usingBkgnd_) {
1971  for (LauBkgndYieldList::const_iterator iter = bkgndEvents_.begin(); iter != bkgndEvents_.end(); ++iter) {
1972  eventSum += (*iter);
1973  }
1974  }
1975  return eventSum;
1976 }
1977 
1978 void LauCPFitModel::getEvtDPLikelihood(UInt_t iEvt)
1979 {
1980  // Function to return the signal and background likelihoods for the
1981  // Dalitz plot for the given event evtNo.
1982 
1983  if ( ! this->useDP() ) {
1984  // There's always going to be a term in the likelihood for the
1985  // signal, so we'd better not zero it.
1986  sigDPLike_ = 1.0;
1987  scfDPLike_ = 1.0;
1988 
1989  const UInt_t nBkgnds = this->nBkgndClasses();
1990  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
1991  if (usingBkgnd_ == kTRUE) {
1992  bkgndDPLike_[bkgndID] = 1.0;
1993  } else {
1994  bkgndDPLike_[bkgndID] = 0.0;
1995  }
1996  }
1997 
1998  return;
1999  }
2000 
2001  const UInt_t nBkgnds = this->nBkgndClasses();
2002  if ( tagged_ ) {
2003  if (curEvtCharge_==+1) {
2004  posSigModel_->calcLikelihoodInfo(iEvt);
2005  sigDPLike_ = posSigModel_->getEvtDPAmp().abs2();
2006  sigDPLike_ *= posSigModel_->getEvtEff();
2007 
2008  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
2009  if (usingBkgnd_ == kTRUE) {
2010  bkgndDPLike_[bkgndID] = posBkgndDPModels_[bkgndID]->getLikelihood(iEvt);
2011  } else {
2012  bkgndDPLike_[bkgndID] = 0.0;
2013  }
2014  }
2015  } else {
2016  negSigModel_->calcLikelihoodInfo(iEvt);
2017  sigDPLike_ = negSigModel_->getEvtDPAmp().abs2();
2018  sigDPLike_ *= negSigModel_->getEvtEff();
2019 
2020  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
2021  if (usingBkgnd_ == kTRUE) {
2022  bkgndDPLike_[bkgndID] = negBkgndDPModels_[bkgndID]->getLikelihood(iEvt);
2023  } else {
2024  bkgndDPLike_[bkgndID] = 0.0;
2025  }
2026  }
2027  }
2028  } else {
2029  posSigModel_->calcLikelihoodInfo(iEvt);
2030  negSigModel_->calcLikelihoodInfo(iEvt);
2031 
2032  sigDPLike_ = 0.5 * ( posSigModel_->getEvtDPAmp().abs2() * posSigModel_->getEvtEff() +
2033  negSigModel_->getEvtDPAmp().abs2() * negSigModel_->getEvtEff() );
2034 
2035  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
2036  if (usingBkgnd_ == kTRUE) {
2037  bkgndDPLike_[bkgndID] = 0.5 * ( posBkgndDPModels_[bkgndID]->getLikelihood(iEvt) +
2038  negBkgndDPModels_[bkgndID]->getLikelihood(iEvt) );
2039  } else {
2040  bkgndDPLike_[bkgndID] = 0.0;
2041  }
2042  }
2043  }
2044 
2045  if ( useSCF_ == kTRUE ) {
2046  if ( scfMap_ == 0 ) {
2047  // we're not smearing the SCF DP position
2048  // so the likelihood is the same as the TM
2049  scfDPLike_ = sigDPLike_;
2050  } else {
2051  // calculate the smeared SCF DP likelihood
2052  scfDPLike_ = this->getEvtSCFDPLikelihood(iEvt);
2053  }
2054  }
2055 
2056  // Calculate the signal normalisation
2057  // NB the 2.0 is there so that the 0.5 factor is applied to
2058  // signal and background in the same place otherwise you get
2059  // normalisation problems when you switch off the DP in the fit
2060  Double_t norm = negSigModel_->getDPNorm() + posSigModel_->getDPNorm();
2061  sigDPLike_ *= 2.0/norm;
2062  scfDPLike_ *= 2.0/norm;
2063 }
2064 
2065 Double_t LauCPFitModel::getEvtSCFDPLikelihood(UInt_t iEvt)
2066 {
2067  Double_t scfDPLike(0.0);
2068 
2069  Double_t recoJacobian(1.0);
2070  Double_t xCoord(0.0);
2071  Double_t yCoord(0.0);
2072 
2073  Bool_t squareDP = negKinematics_->squareDP();
2074  if ( squareDP ) {
2075  xCoord = negSigModel_->getEvtmPrime();
2076  yCoord = negSigModel_->getEvtthPrime();
2077  recoJacobian = recoJacobians_[iEvt];
2078  } else {
2079  xCoord = negSigModel_->getEvtm13Sq();
2080  yCoord = negSigModel_->getEvtm23Sq();
2081  }
2082 
2083  // Find the bin that our reco event falls in
2084  Int_t recoBin = scfMap_->binNumber( xCoord, yCoord );
2085 
2086  // Find out which true Bins contribute to the given reco bin
2087  const std::vector<Int_t>* trueBins = scfMap_->trueBins(recoBin);
2088 
2089  const Int_t nDataEvents = this->eventsPerExpt();
2090 
2091  // Loop over the true bins
2092  for (std::vector<Int_t>::const_iterator iter = trueBins->begin(); iter != trueBins->end(); ++iter)
2093  {
2094  Int_t trueBin = (*iter);
2095 
2096  // prob of a true event in the given true bin migrating to the reco bin
2097  Double_t pRecoGivenTrue = scfMap_->prob( recoBin, trueBin );
2098  Double_t pTrue(0.0);
2099 
2100  // We've cached the DP amplitudes and the efficiency for the
2101  // true bin centres, just after the data points
2102  if ( tagged_ ) {
2103  LauAbsDPDynamics* sigModel(0);
2104  if (curEvtCharge_<0) {
2105  sigModel = negSigModel_;
2106  } else {
2107  sigModel = posSigModel_;
2108  }
2109 
2110  sigModel->calcLikelihoodInfo( nDataEvents + trueBin );
2111 
2112  pTrue = sigModel->getEvtDPAmp().abs2() * sigModel->getEvtEff();
2113  } else {
2114  posSigModel_->calcLikelihoodInfo( nDataEvents + trueBin );
2115  negSigModel_->calcLikelihoodInfo( nDataEvents + trueBin );
2116 
2117  pTrue = 0.5 * ( posSigModel_->getEvtDPAmp().abs2() * posSigModel_->getEvtEff() +
2118  negSigModel_->getEvtDPAmp().abs2() * negSigModel_->getEvtEff() );
2119  }
2120 
2121  // Get the cached SCF fraction (and jacobian if we're using the square DP)
2122  Double_t scfFraction = fakeSCFFracs_[ trueBin ];
2123  Double_t jacobian(1.0);
2124  if ( squareDP ) {
2125  jacobian = fakeJacobians_[ trueBin ];
2126  }
2127 
2128  scfDPLike += pTrue * jacobian * scfFraction * pRecoGivenTrue;
2129  }
2130 
2131  // Divide by the reco jacobian
2132  scfDPLike /= recoJacobian;
2133 
2134  return scfDPLike;
2135 }
2136 
2137 void LauCPFitModel::getEvtExtraLikelihoods(UInt_t iEvt)
2138 {
2139  // Function to return the signal and background likelihoods for the
2140  // extra variables for the given event evtNo.
2141 
2142  sigExtraLike_ = 1.0;
2143 
2144  const UInt_t nBkgnds = this->nBkgndClasses();
2145 
2146  if ( ! tagged_ || curEvtCharge_<0 ) {
2147  sigExtraLike_ = this->prodPdfValue( negSignalPdfs_, iEvt );
2148 
2149  if (useSCF_) {
2150  scfExtraLike_ = this->prodPdfValue( negScfPdfs_, iEvt );
2151  }
2152 
2153  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
2154  if (usingBkgnd_) {
2155  bkgndExtraLike_[bkgndID] = this->prodPdfValue( negBkgndPdfs_[bkgndID], iEvt );
2156  } else {
2157  bkgndExtraLike_[bkgndID] = 0.0;
2158  }
2159  }
2160  } else {
2161  sigExtraLike_ = this->prodPdfValue( negSignalPdfs_, iEvt );
2162 
2163  if (useSCF_) {
2164  scfExtraLike_ = this->prodPdfValue( negScfPdfs_, iEvt );
2165  }
2166 
2167  for ( UInt_t bkgndID(0); bkgndID < nBkgnds; ++bkgndID ) {
2168  if (usingBkgnd_) {
2169  bkgndExtraLike_[bkgndID] = this->prodPdfValue( negBkgndPdfs_[bkgndID], iEvt );
2170  } else {
2171  bkgndExtraLike_[bkgndID] = 0.0;
2172  }
2173  }
2174  }
2175 }
2176 
2177 void LauCPFitModel::updateCoeffs()
2178 {
2179  negCoeffs_.clear(); posCoeffs_.clear();
2180  negCoeffs_.reserve(nSigComp_); posCoeffs_.reserve(nSigComp_);
2181  for (UInt_t i = 0; i < nSigComp_; i++) {
2182  negCoeffs_.push_back(coeffPars_[i]->antiparticleCoeff());
2183  posCoeffs_.push_back(coeffPars_[i]->particleCoeff());
2184  }
2185 }
2186 
2187 void LauCPFitModel::setupSPlotNtupleBranches()
2188 {
2189  // add branches for storing the experiment number and the number of
2190  // the event within the current experiment
2191  this->addSPlotNtupleIntegerBranch("iExpt");
2192  this->addSPlotNtupleIntegerBranch("iEvtWithinExpt");
2193 
2194  // Store the efficiency of the event (for inclusive BF calculations).
2195  if (this->storeDPEff()) {
2196  this->addSPlotNtupleDoubleBranch("efficiency");
2197  if ( negSigModel_->usingScfModel() && posSigModel_->usingScfModel() ) {
2198  this->addSPlotNtupleDoubleBranch("scffraction");
2199  }
2200  }
2201 
2202  // Store the total event likelihood for each species.
2203  if (useSCF_) {
2204  this->addSPlotNtupleDoubleBranch("sigTMTotalLike");
2205  this->addSPlotNtupleDoubleBranch("sigSCFTotalLike");
2206  this->addSPlotNtupleDoubleBranch("sigSCFFrac");
2207  } else {
2208  this->addSPlotNtupleDoubleBranch("sigTotalLike");
2209  }
2210  if (usingBkgnd_) {
2211  const UInt_t nBkgnds = this->nBkgndClasses();
2212  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2213  TString name( this->bkgndClassName(iBkgnd) );
2214  name += "TotalLike";
2215  this->addSPlotNtupleDoubleBranch(name);
2216  }
2217  }
2218 
2219  // Store the DP likelihoods
2220  if (this->useDP()) {
2221  if (useSCF_) {
2222  this->addSPlotNtupleDoubleBranch("sigTMDPLike");
2223  this->addSPlotNtupleDoubleBranch("sigSCFDPLike");
2224  } else {
2225  this->addSPlotNtupleDoubleBranch("sigDPLike");
2226  }
2227  if (usingBkgnd_) {
2228  const UInt_t nBkgnds = this->nBkgndClasses();
2229  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2230  TString name( this->bkgndClassName(iBkgnd) );
2231  name += "DPLike";
2232  this->addSPlotNtupleDoubleBranch(name);
2233  }
2234  }
2235  }
2236 
2237  // Store the likelihoods for each extra PDF
2238  if (useSCF_) {
2239  this->addSPlotNtupleBranches(&negSignalPdfs_, "sigTM");
2240  this->addSPlotNtupleBranches(&negScfPdfs_, "sigSCF");
2241  } else {
2242  this->addSPlotNtupleBranches(&negSignalPdfs_, "sig");
2243  }
2244  if (usingBkgnd_) {
2245  const UInt_t nBkgnds = this->nBkgndClasses();
2246  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2247  const TString& bkgndClass = this->bkgndClassName(iBkgnd);
2248  const LauPdfList* pdfList = &(negBkgndPdfs_[iBkgnd]);
2249  this->addSPlotNtupleBranches(pdfList, bkgndClass);
2250  }
2251  }
2252 }
2253 
2254 void LauCPFitModel::addSPlotNtupleBranches(const LauPdfList* extraPdfs, const TString& prefix)
2255 {
2256  if (extraPdfs) {
2257  // Loop through each of the PDFs
2258  for (LauPdfList::const_iterator pdf_iter = extraPdfs->begin(); pdf_iter != extraPdfs->end(); ++pdf_iter) {
2259 
2260  // Count the number of input variables that are not
2261  // DP variables (used in the case where there is DP
2262  // dependence for e.g. MVA)
2263  UInt_t nVars(0);
2264  std::vector<TString> varNames = (*pdf_iter)->varNames();
2265  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
2266  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
2267  ++nVars;
2268  }
2269  }
2270 
2271  if ( nVars == 1 ) {
2272  // If the PDF only has one variable then
2273  // simply add one branch for that variable
2274  TString varName = (*pdf_iter)->varName();
2275  TString name(prefix);
2276  name += varName;
2277  name += "Like";
2278  this->addSPlotNtupleDoubleBranch(name);
2279  } else if ( nVars == 2 ) {
2280  // If the PDF has two variables then we
2281  // need a branch for them both together and
2282  // branches for each
2283  TString allVars("");
2284  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
2285  allVars += (*var_iter);
2286  TString name(prefix);
2287  name += (*var_iter);
2288  name += "Like";
2289  this->addSPlotNtupleDoubleBranch(name);
2290  }
2291  TString name(prefix);
2292  name += allVars;
2293  name += "Like";
2294  this->addSPlotNtupleDoubleBranch(name);
2295  } else {
2296  cerr << "WARNING in LauCPFitModel::addSPlotNtupleBranches : Can't yet deal with 3D PDFs." << endl;
2297  }
2298  }
2299  }
2300 }
2301 
2302 Double_t LauCPFitModel::setSPlotNtupleBranchValues(LauPdfList* extraPdfs, const TString& prefix, UInt_t iEvt)
2303 {
2304  // Store the PDF value for each variable in the list
2305  Double_t totalLike(1.0);
2306  Double_t extraLike(0.0);
2307  if (extraPdfs) {
2308  for (LauPdfList::iterator pdf_iter = extraPdfs->begin(); pdf_iter != extraPdfs->end(); ++pdf_iter) {
2309 
2310  // calculate the likelihood for this event
2311  (*pdf_iter)->calcLikelihoodInfo(iEvt);
2312  extraLike = (*pdf_iter)->getLikelihood();
2313  totalLike *= extraLike;
2314 
2315  // Count the number of input variables that are not
2316  // DP variables (used in the case where there is DP
2317  // dependence for e.g. MVA)
2318  UInt_t nVars(0);
2319  std::vector<TString> varNames = (*pdf_iter)->varNames();
2320  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
2321  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
2322  ++nVars;
2323  }
2324  }
2325 
2326  if ( nVars == 1 ) {
2327  // If the PDF only has one variable then
2328  // simply store the value for that variable
2329  TString varName = (*pdf_iter)->varName();
2330  TString name(prefix);
2331  name += varName;
2332  name += "Like";
2333  this->setSPlotNtupleDoubleBranchValue(name, extraLike);
2334  } else if ( nVars == 2 ) {
2335  // If the PDF has two variables then we
2336  // store the value for them both together
2337  // and for each on their own
2338  TString allVars("");
2339  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
2340  allVars += (*var_iter);
2341  TString name(prefix);
2342  name += (*var_iter);
2343  name += "Like";
2344  Double_t indivLike = (*pdf_iter)->getLikelihood( (*var_iter) );
2345  this->setSPlotNtupleDoubleBranchValue(name, indivLike);
2346  }
2347  TString name(prefix);
2348  name += allVars;
2349  name += "Like";
2350  this->setSPlotNtupleDoubleBranchValue(name, extraLike);
2351  } else {
2352  cerr << "WARNING in LauCPFitModel::setSPlotNtupleBranchValues : Can't yet deal with 3D PDFs." << endl;
2353  }
2354  }
2355  }
2356  return totalLike;
2357 }
2358 
2359 LauSPlot::NameSet LauCPFitModel::variableNames() const
2360 {
2361  LauSPlot::NameSet nameSet;
2362  if (this->useDP()) {
2363  nameSet.insert("DP");
2364  }
2365  // Loop through all the signal PDFs
2366  for (LauPdfList::const_iterator pdf_iter = negSignalPdfs_.begin(); pdf_iter != negSignalPdfs_.end(); ++pdf_iter) {
2367  // Loop over the variables involved in each PDF
2368  std::vector<TString> varNames = (*pdf_iter)->varNames();
2369  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
2370  // If they are not DP coordinates then add them
2371  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
2372  nameSet.insert( (*var_iter) );
2373  }
2374  }
2375  }
2376  return nameSet;
2377 }
2378 
2379 LauSPlot::NumbMap LauCPFitModel::freeSpeciesNames() const
2380 {
2381  LauSPlot::NumbMap numbMap;
2382  if (!signalEvents_.fixed() && this->doEMLFit()) {
2383  numbMap["sig"] = signalEvents_.genValue();
2384  }
2385  if ( usingBkgnd_ ) {
2386  const UInt_t nBkgnds = this->nBkgndClasses();
2387  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2388  const TString& bkgndClass = this->bkgndClassName(iBkgnd);
2389  const LauParameter& par = bkgndEvents_[iBkgnd];
2390  if (!par.fixed()) {
2391  numbMap[bkgndClass] = par.genValue();
2392  }
2393  }
2394  }
2395  return numbMap;
2396 }
2397 
2398 LauSPlot::NumbMap LauCPFitModel::fixdSpeciesNames() const
2399 {
2400  LauSPlot::NumbMap numbMap;
2401  if (signalEvents_.fixed() && this->doEMLFit()) {
2402  numbMap["sig"] = signalEvents_.genValue();
2403  }
2404  if ( usingBkgnd_ ) {
2405  const UInt_t nBkgnds = this->nBkgndClasses();
2406  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2407  const TString& bkgndClass = this->bkgndClassName(iBkgnd);
2408  const LauParameter& par = bkgndEvents_[iBkgnd];
2409  if (par.fixed()) {
2410  numbMap[bkgndClass] = par.genValue();
2411  }
2412  }
2413  }
2414  return numbMap;
2415 }
2416 
2417 LauSPlot::TwoDMap LauCPFitModel::twodimPDFs() const
2418 {
2419  // This makes the assumption that the form of the positive and
2420  // negative PDFs are the same, which seems reasonable to me
2421 
2422  LauSPlot::TwoDMap twodimMap;
2423 
2424  for (LauPdfList::const_iterator pdf_iter = negSignalPdfs_.begin(); pdf_iter != negSignalPdfs_.end(); ++pdf_iter) {
2425  // Count the number of input variables that are not DP variables
2426  UInt_t nVars(0);
2427  std::vector<TString> varNames = (*pdf_iter)->varNames();
2428  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
2429  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
2430  ++nVars;
2431  }
2432  }
2433  if ( nVars == 2 ) {
2434  if (useSCF_) {
2435  twodimMap.insert( std::make_pair( "sigTM", std::make_pair( varNames[0], varNames[1] ) ) );
2436  } else {
2437  twodimMap.insert( std::make_pair( "sig", std::make_pair( varNames[0], varNames[1] ) ) );
2438  }
2439  }
2440  }
2441 
2442  if ( useSCF_ ) {
2443  for (LauPdfList::const_iterator pdf_iter = negScfPdfs_.begin(); pdf_iter != negScfPdfs_.end(); ++pdf_iter) {
2444  // Count the number of input variables that are not DP variables
2445  UInt_t nVars(0);
2446  std::vector<TString> varNames = (*pdf_iter)->varNames();
2447  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
2448  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
2449  ++nVars;
2450  }
2451  }
2452  if ( nVars == 2 ) {
2453  twodimMap.insert( std::make_pair( "sigSCF", std::make_pair( varNames[0], varNames[1] ) ) );
2454  }
2455  }
2456  }
2457 
2458  if (usingBkgnd_) {
2459  const UInt_t nBkgnds = this->nBkgndClasses();
2460  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2461  const TString& bkgndClass = this->bkgndClassName(iBkgnd);
2462  const LauPdfList& pdfList = negBkgndPdfs_[iBkgnd];
2463  for (LauPdfList::const_iterator pdf_iter = pdfList.begin(); pdf_iter != pdfList.end(); ++pdf_iter) {
2464  // Count the number of input variables that are not DP variables
2465  UInt_t nVars(0);
2466  std::vector<TString> varNames = (*pdf_iter)->varNames();
2467  for ( std::vector<TString>::const_iterator var_iter = varNames.begin(); var_iter != varNames.end(); ++var_iter ) {
2468  if ( (*var_iter) != "m13Sq" && (*var_iter) != "m23Sq" ) {
2469  ++nVars;
2470  }
2471  }
2472  if ( nVars == 2 ) {
2473  twodimMap.insert( std::make_pair( bkgndClass, std::make_pair( varNames[0], varNames[1] ) ) );
2474  }
2475  }
2476  }
2477  }
2478 
2479  return twodimMap;
2480 }
2481 
2482 void LauCPFitModel::storePerEvtLlhds()
2483 {
2484  cout << "Storing per-event likelihood values..." << endl;
2485 
2486  // if we've not been using the DP model then we need to cache all
2487  // the info here so that we can get the efficiency from it
2488 
2489  LauFitDataTree* inputFitData = this->fitData();
2490 
2491  if (!this->useDP() && this->storeDPEff()) {
2492  negSigModel_->initialise(negCoeffs_);
2493  posSigModel_->initialise(posCoeffs_);
2494  negSigModel_->fillDataTree(*inputFitData);
2495  posSigModel_->fillDataTree(*inputFitData);
2496  }
2497 
2498  UInt_t evtsPerExpt(this->eventsPerExpt());
2499 
2500  LauAbsDPDynamics* sigModel(0);
2501  LauPdfList* sigPdfs(0);
2502  LauPdfList* scfPdfs(0);
2503  LauBkgndPdfsList* bkgndPdfs(0);
2504 
2505  for (UInt_t iEvt = 0; iEvt < evtsPerExpt; ++iEvt) {
2506 
2507  this->setSPlotNtupleIntegerBranchValue("iExpt",this->iExpt());
2508  this->setSPlotNtupleIntegerBranchValue("iEvtWithinExpt",iEvt);
2509 
2510  // Find out whether we have B- or B+
2511  if ( tagged_ ) {
2512  const LauFitData& dataValues = inputFitData->getData(iEvt);
2513  LauFitData::const_iterator iter = dataValues.find("charge");
2514  curEvtCharge_ = static_cast<Int_t>(iter->second);
2515 
2516  if (curEvtCharge_==+1) {
2517  sigModel = posSigModel_;
2518  sigPdfs = &posSignalPdfs_;
2519  scfPdfs = &posScfPdfs_;
2520  bkgndPdfs = &posBkgndPdfs_;
2521  } else {
2522  sigModel = negSigModel_;
2523  sigPdfs = &negSignalPdfs_;
2524  scfPdfs = &negScfPdfs_;
2525  bkgndPdfs = &negBkgndPdfs_;
2526  }
2527  } else {
2528  sigPdfs = &negSignalPdfs_;
2529  scfPdfs = &negScfPdfs_;
2530  bkgndPdfs = &negBkgndPdfs_;
2531  }
2532 
2533  // the DP information
2534  this->getEvtDPLikelihood(iEvt);
2535  if (this->storeDPEff()) {
2536  if (!this->useDP()) {
2537  posSigModel_->calcLikelihoodInfo(iEvt);
2538  negSigModel_->calcLikelihoodInfo(iEvt);
2539  }
2540  if ( tagged_ ) {
2541  this->setSPlotNtupleDoubleBranchValue("efficiency",sigModel->getEvtEff());
2542  if ( negSigModel_->usingScfModel() && posSigModel_->usingScfModel() ) {
2543  this->setSPlotNtupleDoubleBranchValue("scffraction",sigModel->getEvtScfFraction());
2544  }
2545  } else {
2546  this->setSPlotNtupleDoubleBranchValue("efficiency",0.5*(posSigModel_->getEvtEff() + negSigModel_->getEvtEff()) );
2547  if ( negSigModel_->usingScfModel() && posSigModel_->usingScfModel() ) {
2548  this->setSPlotNtupleDoubleBranchValue("scffraction",0.5*(posSigModel_->getEvtScfFraction() + negSigModel_->getEvtScfFraction()));
2549  }
2550  }
2551  }
2552  if (this->useDP()) {
2553  sigTotalLike_ = sigDPLike_;
2554  if (useSCF_) {
2555  scfTotalLike_ = scfDPLike_;
2556  this->setSPlotNtupleDoubleBranchValue("sigTMDPLike",sigDPLike_);
2557  this->setSPlotNtupleDoubleBranchValue("sigSCFDPLike",scfDPLike_);
2558  } else {
2559  this->setSPlotNtupleDoubleBranchValue("sigDPLike",sigDPLike_);
2560  }
2561  if (usingBkgnd_) {
2562  const UInt_t nBkgnds = this->nBkgndClasses();
2563  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2564  TString name = this->bkgndClassName(iBkgnd);
2565  name += "DPLike";
2566  this->setSPlotNtupleDoubleBranchValue(name,bkgndDPLike_[iBkgnd]);
2567  }
2568  }
2569  } else {
2570  sigTotalLike_ = 1.0;
2571  if (useSCF_) {
2572  scfTotalLike_ = 1.0;
2573  }
2574  if (usingBkgnd_) {
2575  const UInt_t nBkgnds = this->nBkgndClasses();
2576  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2577  bkgndTotalLike_[iBkgnd] = 1.0;
2578  }
2579  }
2580  }
2581 
2582  // the signal PDF values
2583  if ( useSCF_ ) {
2584  sigTotalLike_ *= this->setSPlotNtupleBranchValues(sigPdfs, "sigTM", iEvt);
2585  scfTotalLike_ *= this->setSPlotNtupleBranchValues(scfPdfs, "sigSCF", iEvt);
2586  } else {
2587  sigTotalLike_ *= this->setSPlotNtupleBranchValues(sigPdfs, "sig", iEvt);
2588  }
2589 
2590  // the beckground PDF values
2591  if (usingBkgnd_) {
2592  const UInt_t nBkgnds = this->nBkgndClasses();
2593  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2594  const TString& bkgndClass = this->bkgndClassName(iBkgnd);
2595  LauPdfList& pdfs = (*bkgndPdfs)[iBkgnd];
2596  bkgndTotalLike_[iBkgnd] *= this->setSPlotNtupleBranchValues(&(pdfs), bkgndClass, iEvt);
2597  }
2598  }
2599 
2600  // the total likelihoods
2601  if (useSCF_) {
2602  Double_t scfFrac(0.0);
2603  if ( useSCFHist_ ) {
2604  scfFrac = recoSCFFracs_[iEvt];
2605  } else {
2606  scfFrac = scfFrac_.value();
2607  }
2608  this->setSPlotNtupleDoubleBranchValue("sigSCFFrac",scfFrac);
2609  sigTotalLike_ *= ( 1.0 - scfFrac );
2610  if ( scfMap_ == 0 ) {
2611  scfTotalLike_ *= scfFrac;
2612  }
2613  this->setSPlotNtupleDoubleBranchValue("sigTMTotalLike",sigTotalLike_);
2614  this->setSPlotNtupleDoubleBranchValue("sigSCFTotalLike",scfTotalLike_);
2615  } else {
2616  this->setSPlotNtupleDoubleBranchValue("sigTotalLike",sigTotalLike_);
2617  }
2618  if (usingBkgnd_) {
2619  const UInt_t nBkgnds = this->nBkgndClasses();
2620  for ( UInt_t iBkgnd(0); iBkgnd < nBkgnds; ++iBkgnd ) {
2621  TString name = this->bkgndClassName(iBkgnd);
2622  name += "TotalLike";
2623  this->setSPlotNtupleDoubleBranchValue(name,bkgndTotalLike_[iBkgnd]);
2624  }
2625  }
2626  // fill the tree
2627  this->fillSPlotNtupleBranches();
2628  }
2629  cout << "Finished storing per-event likelihood values." << endl;
2630 }
2631 
2632 void LauCPFitModel::embedNegSignal(const TString& fileName, const TString& treeName,
2633  Bool_t reuseEventsWithinEnsemble, Bool_t reuseEventsWithinExperiment,
2634  Bool_t useReweighting)
2635 {
2636  if (negSignalTree_) {
2637  cerr << "ERROR in LauCPFitModel::embedNegSignal : Already embedding signal from a file." << endl;
2638  return;
2639  }
2640 
2641  if (!reuseEventsWithinEnsemble && reuseEventsWithinExperiment) {
2642  cerr << "WARNING in LauCPFitModel::embedNegSignal : Conflicting options provided, will not reuse events at all." << endl;
2643  reuseEventsWithinExperiment = kFALSE;
2644  }
2645 
2646  negSignalTree_ = new LauEmbeddedData(fileName,treeName,reuseEventsWithinExperiment);
2647  Bool_t dataOK = negSignalTree_->findBranches();
2648  if (!dataOK) {
2649  delete negSignalTree_; negSignalTree_ = 0;
2650  cerr << "ERROR in LauCPFitModel::embedNegSignal : Problem creating data tree for embedding." << endl;
2651  return;
2652  }
2653  reuseSignal_ = reuseEventsWithinEnsemble;
2654  useNegReweighting_ = useReweighting;
2655  if (this->enableEmbedding() == kFALSE) {this->enableEmbedding(kTRUE);}
2656 }
2657 
2658 void LauCPFitModel::embedNegBkgnd(const TString& bkgndClass, const TString& fileName, const TString& treeName,
2659  Bool_t reuseEventsWithinEnsemble, Bool_t reuseEventsWithinExperiment)
2660 {
2661  if ( ! this->validBkgndClass( bkgndClass ) ) {
2662  cerr << "ERROR in LauCPFitModel::embedBkgnd : Invalid background class \"" << bkgndClass << "\"." << std::endl;
2663  cerr << " : Background class names must be provided in \"setBkgndClassNames\" before any other background-related actions can be performed." << endl;
2664  return;
2665  }
2666 
2667  UInt_t bkgndID = this->bkgndClassID( bkgndClass );
2668 
2669  if (negBkgndTree_[bkgndID]) {
2670  cerr << "ERROR in LauCPFitModel::embedNegBkgnd : Already embedding background from a file." << endl;
2671  return;
2672  }
2673 
2674  if (!reuseEventsWithinEnsemble && reuseEventsWithinExperiment) {
2675  cerr << "WARNING in LauCPFitModel::embedNegBkgnd : Conflicting options provided, will not reuse events at all." << endl;
2676  reuseEventsWithinExperiment = kFALSE;
2677  }
2678 
2679  negBkgndTree_[bkgndID] = new LauEmbeddedData(fileName,treeName,reuseEventsWithinExperiment);
2680  Bool_t dataOK = negBkgndTree_[bkgndID]->findBranches();
2681  if (!dataOK) {
2682  delete negBkgndTree_[bkgndID]; negBkgndTree_[bkgndID] = 0;
2683  cerr << "ERROR in LauCPFitModel::embedNegBkgnd : Problem creating data tree for embedding." << endl;
2684  return;
2685  }
2686  reuseBkgnd_[bkgndID] = reuseEventsWithinEnsemble;
2687  if (this->enableEmbedding() == kFALSE) {this->enableEmbedding(kTRUE);}
2688 }
2689 
2690 void LauCPFitModel::embedPosSignal(const TString& fileName, const TString& treeName,
2691  Bool_t reuseEventsWithinEnsemble, Bool_t reuseEventsWithinExperiment,
2692  Bool_t useReweighting)
2693 {
2694  if (posSignalTree_) {
2695  cerr << "ERROR in LauCPFitModel::embedPosSignal : Already embedding signal from a file." << endl;
2696  return;
2697  }
2698 
2699  if (!reuseEventsWithinEnsemble && reuseEventsWithinExperiment) {
2700  cerr << "WARNING in LauCPFitModel::embedPosSignal : Conflicting options provided, will not reuse events at all." << endl;
2701  reuseEventsWithinExperiment = kFALSE;
2702  }
2703 
2704  posSignalTree_ = new LauEmbeddedData(fileName,treeName,reuseEventsWithinExperiment);
2705  Bool_t dataOK = posSignalTree_->findBranches();
2706  if (!dataOK) {
2707  delete posSignalTree_; posSignalTree_ = 0;
2708  cerr << "ERROR in LauCPFitModel::embedPosSignal : Problem creating data tree for embedding." << endl;
2709  return;
2710  }
2711  reuseSignal_ = reuseEventsWithinEnsemble;
2712  usePosReweighting_ = useReweighting;
2713  if (this->enableEmbedding() == kFALSE) {this->enableEmbedding(kTRUE);}
2714 }
2715 
2716 void LauCPFitModel::embedPosBkgnd(const TString& bkgndClass, const TString& fileName, const TString& treeName,
2717  Bool_t reuseEventsWithinEnsemble, Bool_t reuseEventsWithinExperiment)
2718 {
2719  if ( ! this->validBkgndClass( bkgndClass ) ) {
2720  cerr << "ERROR in LauCPFitModel::embedBkgnd : Invalid background class \"" << bkgndClass << "\"." << std::endl;
2721  cerr << " : Background class names must be provided in \"setBkgndClassNames\" before any other background-related actions can be performed." << endl;
2722  return;
2723  }
2724 
2725  UInt_t bkgndID = this->bkgndClassID( bkgndClass );
2726 
2727  if (posBkgndTree_[bkgndID]) {
2728  cerr << "ERROR in LauCPFitModel::embedPosBkgnd : Already embedding background from a file." << endl;
2729  return;
2730  }
2731 
2732  if (!reuseEventsWithinEnsemble && reuseEventsWithinExperiment) {
2733  cerr << "WARNING in LauCPFitModel::embedPosBkgnd : Conflicting options provided, will not reuse events at all." << endl;
2734  reuseEventsWithinExperiment = kFALSE;
2735  }
2736 
2737  posBkgndTree_[bkgndID] = new LauEmbeddedData(fileName,treeName,reuseEventsWithinExperiment);
2738  Bool_t dataOK = posBkgndTree_[bkgndID]->findBranches();
2739  if (!dataOK) {
2740  delete posBkgndTree_[bkgndID]; posBkgndTree_[bkgndID] = 0;
2741  cerr << "ERROR in LauCPFitModel::embedPosBkgnd : Problem creating data tree for embedding." << endl;
2742  return;
2743  }
2744  reuseBkgnd_[bkgndID] = reuseEventsWithinEnsemble;
2745  if (this->enableEmbedding() == kFALSE) {this->enableEmbedding(kTRUE);}
2746 }
2747 
2748 void LauCPFitModel::weightEvents( const TString& /*dataFileName*/, const TString& /*dataTreeName*/ )
2749 {
2750  cerr << "ERROR in LauCPFitModel::weightEvents : Method not available for this fit model." << endl;
2751  return;
2752 }
2753 
LauParameter::range
Double_t range() const
The range allowed for the parameter.
Definition: LauParameter.hh:186
LauCPFitModel::generateExtraPdfValues
void generateExtraPdfValues(LauPdfList *extraPdfs, LauEmbeddedData *embeddedData)
Generate from the extra PDFs.
Definition: LauCPFitModel.cc:1665
LauCPFitModel::splitSignalComponent
void splitSignalComponent(const TH2 *dpHisto, Bool_t upperHalf=kFALSE, LauScfMap *scfMap=0)
Split the signal component into well-reconstructed and mis-reconstructed parts.
Definition: LauCPFitModel.cc:205
LauCPFitModel::generateSignalEvent
Bool_t generateSignalEvent()
Generate signal event.
Definition: LauCPFitModel.cc:1396
LauAbsFitModel::addSPlotNtupleIntegerBranch
virtual void addSPlotNtupleIntegerBranch(const TString &name)
Add a branch to the sPlot tree for storing an integer.
Definition: LauAbsFitModel.cc:523
LauAbsDPDynamics
Class for defining the abstract interface for signal Dalitz plot dynamics.
Definition: LauAbsDPDynamics.hh:39
LauAbsBkgndDPModel
The abstract interface for a background Dalitz plot model.
Definition: LauAbsBkgndDPModel.hh:31
LauFitNtuple::storeCorrMatrix
void storeCorrMatrix(UInt_t iExpt, Double_t NLL, Int_t fitStatus)
Store the correlation matrix and other fit information.
Definition: LauFitNtuple.cc:61
LauKinematics::squareDP
Bool_t squareDP() const
Are the square Dalitz plot co-ordinates being calculated?
Definition: LauKinematics.hh:51
LauCPFitModel::acp_
std::vector< LauParameter > acp_
A_CP parameter.
Definition: LauCPFitModel.hh:501
LauCPFitModel::variableNames
virtual LauSPlot::NameSet variableNames() const
Returns the names of all variables in the fit.
Definition: LauCPFitModel.cc:2359
LauKinematics::getc23
Double_t getc23() const
Get the cosine of the helicity angle theta23.
Definition: LauKinematics.hh:193
LauRandom::randomFun
TRandom * randomFun()
Access the singleton random number generator with a particular seed.
Definition: LauRandom.cc:20
LauAbsFitModel::LauParameterList
std::vector< LauParameter > LauParameterList
List of parameters.
Definition: LauAbsFitModel.hh:310
LauParameter::fixed
Bool_t fixed() const
Check whether the parameter is fixed or floated.
Definition: LauParameter.hh:192
LauCPFitModel::negKinematics_
LauKinematics * negKinematics_
The B- Dalitz plot kinematics object.
Definition: LauCPFitModel.hh:444
LauEmbeddedData::nEvents
UInt_t nEvents() const
Retrieve the number of events.
Definition: LauEmbeddedData.hh:63
LauAbsFitModel::eventsPerExpt
UInt_t eventsPerExpt() const
Obtain the total number of events in the current experiment.
Definition: LauAbsFitModel.hh:188
LauAbsDPDynamics::getEvtthPrime
virtual Double_t getEvtthPrime() const =0
Retrieve the square Dalitz plot coordinate, theta&#39;, for the current event.
LauCPFitModel::setNBkgndEvents
virtual void setNBkgndEvents(const TString &bkgndClass, Double_t nBkgndEvents, Bool_t fixBkgndEvents=kFALSE)
Set the background event yield(s)
Definition: LauCPFitModel.cc:155
LauCPFitModel::printAsymmetries
virtual void printAsymmetries(std::ostream &output)
Print the asymmetries.
Definition: LauCPFitModel.cc:1017
LauAbsFitModel::writeLatexTable
Bool_t writeLatexTable() const
Determine whether writing out of the latex table is enabled.
Definition: LauAbsFitModel.hh:148
LauCPFitModel::negSignalPdfs_
LauPdfList negSignalPdfs_
The B- signal PDFs.
Definition: LauCPFitModel.hh:450
LauCPFitModel::addSPlotNtupleBranches
void addSPlotNtupleBranches(const LauPdfList *extraPdfs, const TString &prefix)
Add sPlot branches for the extra PDFs.
Definition: LauCPFitModel.cc:2254
LauCPFitModel::forceAsym_
Bool_t forceAsym_
Option to force an asymmetry.
Definition: LauCPFitModel.hh:522
LauCPFitModel::scfTotalLike_
Double_t scfTotalLike_
Total SCF likelihood.
Definition: LauCPFitModel.hh:632
LauCPFitModel::setExtraPdfParameters
void setExtraPdfParameters()
Set the fit parameters for the extra PDFs.
Definition: LauCPFitModel.cc:541
LauKinematics::updateSqDPKinematics
void updateSqDPKinematics(Double_t mPrime, Double_t thetaPrime)
Update all kinematic quantities based on the square DP co-ordinates m&#39; and theta&#39;.
Definition: LauKinematics.cc:102
LauCPFitModel::posBkgndDPModels_
LauBkgndDPModelList posBkgndDPModels_
The B+ background Dalitz plot models.
Definition: LauCPFitModel.hh:441
LauEffModel::setEffHisto
void setEffHisto(const TH2 *effHisto, Bool_t useInterpolation=kTRUE, Bool_t fluctuateBins=kFALSE, Double_t avEff=-1.0, Double_t absError=-1.0, Bool_t useUpperHalfOnly=kFALSE, Bool_t squareDP=kFALSE)
Set the efficiency variation across the phase space using a predetermined 2D histogram.
Definition: LauEffModel.cc:60
LauAbsFitModel::cacheInfo
void cacheInfo(LauPdfList &pdfList, const LauFitDataTree &theData)
Have all PDFs in the list cache the data.
Definition: LauAbsFitModel.cc:1294
LauAbsFitModel::setGenNtupleIntegerBranchValue
virtual void setGenNtupleIntegerBranchValue(const TString &name, Int_t value)
Set the value of an integer branch in the gen tree.
Definition: LauAbsFitModel.cc:498
LauAsymmCalc.hh
File containing declaration of LauAsymmCalc class.
LauCPFitModel::nNormPar_
UInt_t nNormPar_
Number of parameters.
Definition: LauCPFitModel.hh:480
LauCPFitModel::posSignalTree_
LauEmbeddedData * posSignalTree_
The B+ signal event tree.
Definition: LauCPFitModel.hh:589
LauAbsCoeffSet::baseName
virtual TString baseName() const
Retrieve the base name of the coefficient set.
Definition: LauAbsCoeffSet.hh:117
LauParameter::LauParameter
LauParameter()
Default constructor.
Definition: LauParameter.cc:30
LauKinematics::getm23
Double_t getm23() const
Get the m23 invariant mass.
Definition: LauKinematics.hh:159
LauCPFitModel::curEvtCharge_
Int_t curEvtCharge_
Current event charge.
Definition: LauCPFitModel.hh:537
LauCPFitModel::calcExtraFractions
void calcExtraFractions(Bool_t initValues=kFALSE)
Calculate the CP-conserving and CP-violating fit fractions.
Definition: LauCPFitModel.cc:745
LauParameter::name
const TString & name() const
The parameter name.
Definition: LauParameter.hh:132
LauAbsCoeffSet::acp
virtual LauParameter acp()=0
Calculate the CP asymmetry.
LauKinematics::getmPrime
Double_t getmPrime() const
Get m&#39; value.
Definition: LauKinematics.hh:205
LauDaughters
Class that defines the particular 3-body decay under study.
Definition: LauDaughters.hh:33
LauAbsFitModel::addFitParameters
UInt_t addFitParameters(LauPdfList &pdfList)
Add parameters of the PDFs in the list to the list of all fit parameters.
Definition: LauAbsFitModel.cc:1246
LauAsymmCalc
Class for calculating the asymmetry between two variables.
Definition: LauAsymmCalc.hh:25
LauCPFitModel::negSignalTree_
LauEmbeddedData * negSignalTree_
The B- signal event tree.
Definition: LauCPFitModel.hh:586
LauAbsFitModel::useDP
Bool_t useDP() const
Is the Dalitz plot term in the likelihood.
Definition: LauAbsFitModel.hh:81
LauAbsDPDynamics::getDaughters
LauDaughters * getDaughters()
Retrieve the daughters.
Definition: LauAbsDPDynamics.hh:293
LauKinematics::getc13
Double_t getc13() const
Get the cosine of the helicity angle theta13.
Definition: LauKinematics.hh:199
LauAbsCoeffSet.hh
File containing declaration of LauAbsCoeffSet class.
LauEmbeddedData::clearUsedList
void clearUsedList()
Clear the list of used events.
Definition: LauEmbeddedData.hh:104
LauCPFitModel::useNegReweighting_
Bool_t useNegReweighting_
Boolean to use reweighting for B-.
Definition: LauCPFitModel.hh:601
LauScfMap::trueHist
TH2 * trueHist(Int_t trueBin)
Retrieve the migration histogram for trueBin.
Definition: LauScfMap.cc:176
LauCPFitModel::negCoeffs_
std::vector< LauComplex > negCoeffs_
The complex coefficients for B-.
Definition: LauCPFitModel.hh:579
LauCPFitModel::posBkgndPdfs_
LauBkgndPdfsList posBkgndPdfs_
The B+ background PDFs.
Definition: LauCPFitModel.hh:465
LauSPlot::TwoDMap
std::multimap< TString, std::pair< TString, TString > > TwoDMap
Type to associate the name of the species that have 2D PDFs with the names of the two variables invol...
Definition: LauSPlot.hh:68
LauCPFitModel::negFitFrac_
LauParArray negFitFrac_
The B- fit fractions.
Definition: LauCPFitModel.hh:486
LauScfMap::listBinCentres
void listBinCentres(std::vector< Double_t > &xCoords, std::vector< Double_t > &yCoords) const
Create lists of the co-ordinates of the centres of true bins.
Definition: LauScfMap.cc:114
LauCPFitModel::getEvtExtraLikelihoods
virtual void getEvtExtraLikelihoods(UInt_t iEvt)
Determine the signal and background likelihood for the extra variables for a given event...
Definition: LauCPFitModel.cc:2137
LauCPFitModel::sigDPLike_
Double_t sigDPLike_
Signal DP likelihood value.
Definition: LauCPFitModel.hh:611
LauDaughters.hh
File containing declaration of LauDaughters class.
LauScfMap::trueBins
const std::vector< Int_t > * trueBins(Int_t recoBin) const
Find which true bins contribute to the given reco bin.
Definition: LauScfMap.cc:155
LauCPFitModel::negScfPdfs_
LauPdfList negScfPdfs_
The B- SCF PDFs.
Definition: LauCPFitModel.hh:456
LauCPFitModel::useSCF_
Bool_t useSCF_
Is the signal split into TM and SCF.
Definition: LauCPFitModel.hh:543
LauCPFitModel::genExpt
virtual Bool_t genExpt()
Toy MC generation and fitting overloaded functions.
Definition: LauCPFitModel.cc:1239
LauAbsFitModel::setGenNtupleDoubleBranchValue
virtual void setGenNtupleDoubleBranchValue(const TString &name, Double_t value)
Set the value of a double branch in the gen tree.
Definition: LauAbsFitModel.cc:503
LauCPFitModel::negParent_
TString negParent_
Name of the parent particle.
Definition: LauCPFitModel.hh:573
LauCPFitModel::posParent_
TString posParent_
Name of the parent particle.
Definition: LauCPFitModel.hh:576
LauCPFitModel::embedNegSignal
void embedNegSignal(const TString &fileName, const TString &treeName, Bool_t reuseEventsWithinEnsemble, Bool_t reuseEventsWithinExperiment=kFALSE, Bool_t useReweighting=kFALSE)
Embed full simulation events for the B- signal, rather than generating toy from the PDFs...
Definition: LauCPFitModel.cc:2632
LauScfMap.hh
File containing declaration of LauScfMap class.
LauCPFitModel::fakeSCFFracs_
std::vector< Double_t > fakeSCFFracs_
The cached values of the SCF fraction for each bin centre.
Definition: LauCPFitModel.hh:561
LauParArray
std::vector< std::vector< LauParameter > > LauParArray
Type to define an array of parameters.
Definition: LauParameter.hh:510
LauAbsDPDynamics::getDPRate
LauParameter getDPRate() const
Retrieve the overall Dalitz plot rate.
Definition: LauAbsDPDynamics.hh:251
LauAbsFitModel::fitStatus
Int_t fitStatus() const
Access the fit status information.
Definition: LauAbsFitModel.hh:668
LauCPFitModel::storeSignalMCMatch
Bool_t storeSignalMCMatch(LauEmbeddedData *embeddedData)
Store the MC truth info on the TM/SCF nature of the embedded signal event.
Definition: LauCPFitModel.cc:1702
LauCPFitModel::usingBkgnd_
Bool_t usingBkgnd_
Background boolean.
Definition: LauCPFitModel.hh:468
LauPrint.hh
File containing declaration of LauPrint class.
LauAbsDPDynamics::calcLikelihoodInfo
virtual void calcLikelihoodInfo(Double_t m13Sq, Double_t m23Sq)=0
Calculate the likelihood (and all associated information) given values of the Dalitz plot coordinates...
LauCPFitModel::eventsToGenerate
LauGenInfo eventsToGenerate()
Determine the number of events to generate for each hypothesis.
Definition: LauCPFitModel.cc:1168
LauAbsFitModel::prodPdfValue
Double_t prodPdfValue(LauPdfList &pdfList, UInt_t iEvt)
Calculate the product of the per-event likelihoods of the PDFs in the list.
Definition: LauAbsFitModel.cc:1301
LauAbsFitModel::LauPdfList
std::vector< LauAbsPdf * > LauPdfList
List of Pdfs.
Definition: LauAbsFitModel.hh:306
LauFitData
std::map< TString, Double_t > LauFitData
Type for holding event data.
Definition: LauFitDataTree.hh:31
LauCPFitModel::bkgndAsym_
LauBkgndYieldList bkgndAsym_
Background asymmetries(s)
Definition: LauCPFitModel.hh:528
LauKinematics::getm23Sq
Double_t getm23Sq() const
Get the m23 invariant mass square.
Definition: LauKinematics.hh:175
LauDaughters::squareDP
Bool_t squareDP() const
Determine to use or not the square Dalitz plot.
Definition: LauDaughters.hh:72
LauAbsDPDynamics::hasResonance
virtual Bool_t hasResonance(const TString &resName) const
Check whether this model includes a named resonance.
Definition: LauAbsDPDynamics.cc:107
LauAbsDPDynamics::getExtraParameters
std::vector< LauParameter > getExtraParameters()
Retrieve any extra parameters/quantities (e.g. K-matrix total fit fractions)
Definition: LauAbsDPDynamics.hh:335
LauCPFitModel::negDPRate_
LauParameter negDPRate_
The average DP rate for B-.
Definition: LauCPFitModel.hh:510
LauPrint
Class to define various output print commands.
Definition: LauPrint.hh:29
LauEmbeddedData::getValue
Double_t getValue(const TString &name) const
Get the value of a specified branch.
Definition: LauEmbeddedData.cc:175
LauCPFitModel::bkgndTotalLike_
std::vector< Double_t > bkgndTotalLike_
Total background likelihood(s)
Definition: LauCPFitModel.hh:635
LauAbsFitModel::fillGenNtupleBranches
virtual void fillGenNtupleBranches()
Fill the gen tuple branches.
Definition: LauAbsFitModel.cc:518
LauCPFitModel::recoSCFFracs_
std::vector< Double_t > recoSCFFracs_
The cached values of the SCF fraction for each event.
Definition: LauCPFitModel.hh:558
LauCPFitModel::posDPRate_
LauParameter posDPRate_
The average DP rate for B+.
Definition: LauCPFitModel.hh:513
LauFitNtuple::updateFitNtuple
void updateFitNtuple()
Update the fit ntuple.
Definition: LauFitNtuple.cc:186
LauKinematics.hh
File containing declaration of LauKinematics class.
LauAbsFitModel::setSPlotNtupleDoubleBranchValue
virtual void setSPlotNtupleDoubleBranchValue(const TString &name, Double_t value)
Set the value of a double branch in the sPlot tree.
Definition: LauAbsFitModel.cc:538
LauAbsDPDynamics::getMeanEff
LauParameter getMeanEff() const
Retrieve the mean efficiency across the Dalitz plot.
Definition: LauAbsDPDynamics.hh:245
LauEmbeddedData.hh
File containing declaration of LauEmbeddedData class.
LauCPFitModel::bkgndExtraLike_
std::vector< Double_t > bkgndExtraLike_
Background likelihood value(s) from extra PDFs.
Definition: LauCPFitModel.hh:626
LauCPFitModel::setupGenNtupleBranches
void setupGenNtupleBranches()
Setup the required ntuple branches.
Definition: LauCPFitModel.cc:1597
LauAbsFitModel::printFitParameters
void printFitParameters(const LauPdfList &pdfList, std::ostream &fout) const
Print the fit parameters for all PDFs in the list.
Definition: LauAbsFitModel.cc:1273
LauCPFitModel::posSigModel_
LauAbsDPDynamics * posSigModel_
The B+ signal Dalitz plot model.
Definition: LauCPFitModel.hh:435
LauCPFitModel::scfExtraLike_
Double_t scfExtraLike_
SCF likelihood from extra PDFs.
Definition: LauCPFitModel.hh:623
LauEmbeddedData::haveBranch
Bool_t haveBranch(const TString &name) const
Boolean to determine whether branch exists.
Definition: LauEmbeddedData.hh:75
LauCPFitModel.hh
File containing declaration of LauCPFitModel class.
LauCPFitModel::reuseSignal_
Bool_t reuseSignal_
Boolean to reuse signal events.
Definition: LauCPFitModel.hh:598
LauCPFitModel::posFitFrac_
LauParArray posFitFrac_
The B+ fit fractions.
Definition: LauCPFitModel.hh:489
LauCPFitModel::cacheInputFitVars
virtual void cacheInputFitVars()
Read in the input fit data variables, e.g. m13Sq and m23Sq.
Definition: LauCPFitModel.cc:1768
LauCPFitModel::setSignalDPParameters
void setSignalDPParameters()
Set the fit parameters for the DP model.
Definition: LauCPFitModel.cc:508
LauCPFitModel::randomiseInitFitPars
void randomiseInitFitPars()
Randomise the initial fit parameters.
Definition: LauCPFitModel.cc:1158
LauCPFitModel::appendBinCentres
void appendBinCentres(LauFitDataTree *inputData)
Definition: LauCPFitModel.cc:1857
LauAbsCoeffSet::index
virtual UInt_t index() const
Retrieve the index number of the coefficient set.
Definition: LauAbsCoeffSet.hh:131
LauKinematics::updateKinematics
void updateKinematics(Double_t m13Sq, Double_t m23Sq)
Update all kinematic quantities based on the DP co-ordinates m13Sq and m23Sq.
Definition: LauKinematics.cc:85
LauCPFitModel::checkInitFitParams
virtual void checkInitFitParams()
Check the initial fit parameters.
Definition: LauCPFitModel.cc:1146
LauAbsFitModel::extraPars
const LauParameterList & extraPars() const
Access the extra variables.
Definition: LauAbsFitModel.hh:645
LauComplex::abs2
Double_t abs2() const
Obtain the square of the absolute value of the complex number.
Definition: LauComplex.hh:229
LauAbsFitModel::LauParameterPList
std::vector< LauParameter * > LauParameterPList
List of parameter pointers.
Definition: LauAbsFitModel.hh:308
LauAbsDPDynamics::fillDataTree
virtual void fillDataTree(const LauFitDataTree &inputFitTree)=0
Obtain data from a fit tree.
LauCPFitModel::tagged_
const Bool_t tagged_
IS the analysis tagged?
Definition: LauCPFitModel.hh:531
LauCPFitModel::nExtraPdfPar_
UInt_t nExtraPdfPar_
Number of extra PDF parameters.
Definition: LauCPFitModel.hh:477
LauCPFitModel::finaliseFitResults
virtual void finaliseFitResults(const TString &tablePrefixName)
Get the fit results and store them.
Definition: LauCPFitModel.cc:817
LauEffModel::passVeto
Bool_t passVeto(const LauKinematics *kinematics) const
Determine whether the given DP position is outside the vetoes.
Definition: LauEffModel.cc:134
LauCPFitModel::printFitFractions
virtual void printFitFractions(std::ostream &output)
Print the fit fractions, total DP rate and mean efficiency.
Definition: LauCPFitModel.cc:996
LauEffModel::calcEfficiency
Double_t calcEfficiency(const LauKinematics *kinematics) const
Determine the efficiency for a given point in the Dalitz plot.
Definition: LauEffModel.cc:95
LauAbsDPDynamics::calcExtraInfo
virtual void calcExtraInfo(Bool_t init=kFALSE)=0
Calculate the fit fractions, mean efficiency and total DP rate.
LauCPFitModel
Class for defining a CP fit model.
Definition: LauCPFitModel.hh:47
LauAbsFitModel
Abstract interface to the fitting and toy MC model.
Definition: LauAbsFitModel.hh:71
LauAbsFitModel::iExpt
UInt_t iExpt() const
Obtain the number of the current experiment.
Definition: LauAbsFitModel.hh:197
LauCPFitModel::tagVarName_
const TString tagVarName_
Event charge.
Definition: LauCPFitModel.hh:534
LauAbsDPDynamics::getEvtScfFraction
virtual Double_t getEvtScfFraction() const =0
Retrieve the fraction of events that are poorly reconstructed (the self cross feed fraction) in the D...
LauAbsFitModel::addGenNtupleIntegerBranch
virtual void addGenNtupleIntegerBranch(const TString &name)
Add a branch to the gen tree for storing an integer.
Definition: LauAbsFitModel.cc:488
LauCPFitModel::evtCharges_
std::vector< Int_t > evtCharges_
Vector to store event charges.
Definition: LauCPFitModel.hh:540
LauAbsFitModel::fitNtuple
const LauFitNtuple * fitNtuple() const
Access the fit ntuple.
Definition: LauAbsFitModel.hh:649
LauScfMap::binNumber
Int_t binNumber(Double_t xCoord, Double_t yCoord) const
Find the global bin number for the given co-ordinates.
Definition: LauScfMap.cc:144
LauAbsFitModel::doPoissonSmearing
Bool_t doPoissonSmearing() const
Determine whether Poisson smearing is enabled for the toy MC generation.
Definition: LauAbsFitModel.hh:109
LauKinematics::getc12
Double_t getc12() const
Get the cosine of the helicity angle theta12.
Definition: LauKinematics.hh:187
LauKinematics::getm13
Double_t getm13() const
Get the m13 invariant mass.
Definition: LauKinematics.hh:164
LauCPFitModel::LauBkgndPdfsList
std::vector< LauPdfList > LauBkgndPdfsList
Typedef for a vector of background PDFs.
Definition: LauCPFitModel.hh:226
LauCPFitModel::calcAsymmetries
void calcAsymmetries(Bool_t initValues=kFALSE)
Calculate the CP asymmetries.
Definition: LauCPFitModel.cc:798
LauSPlot::NumbMap
std::map< TString, Double_t > NumbMap
Type to associate a category name with a double precision number, e.g. a yield or PDF value for a giv...
Definition: LauSPlot.hh:62
LauAbsFitModel::doEMLFit
Bool_t doEMLFit() const
Determine whether an extended maximum likelihood fit it being performed.
Definition: LauAbsFitModel.hh:100
LauCPFitModel::scfMap_
LauScfMap * scfMap_
The smearing matrix for the SCF DP PDF.
Definition: LauCPFitModel.hh:555
LauComplex.hh
File containing declaration of LauComplex class.
LauCPFitModel::sigExtraLike_
Double_t sigExtraLike_
Signal likelihood from extra PDFs.
Definition: LauCPFitModel.hh:620
LauAbsDPDynamics::getConjResName
virtual TString getConjResName(const TString &resName) const
Retrieve the name of the charge conjugate of a named resonance.
Definition: LauAbsDPDynamics.cc:117
LauEmbeddedData::getValues
LauFitData getValues(const std::vector< TString > &names) const
Get values of specified branches.
Definition: LauEmbeddedData.cc:186
LauCPFitModel::initialise
virtual void initialise()
Initialise the fit.
Definition: LauCPFitModel.cc:374
LauCPFitModel::fixdSpeciesNames
virtual LauSPlot::NumbMap fixdSpeciesNames() const
Returns the names and yields of species that are fixed in the fit.
Definition: LauCPFitModel.cc:2398
LauCPFitModel::setNSigEvents
virtual void setNSigEvents(Double_t nSigEvents, Bool_t fixSigEvents=kFALSE)
Set the signal event yield.
Definition: LauCPFitModel.cc:123
LauAbsCoeffSet
Class for defining the abstract interface for complex coefficient classes.
Definition: LauAbsCoeffSet.hh:34
LauAbsFitModel::nll
Double_t nll() const
Access the current NLL value.
Definition: LauAbsFitModel.hh:665
LauCPFitModel::sigTotalLike_
Double_t sigTotalLike_
Total signal likelihood.
Definition: LauCPFitModel.hh:629
LauAbsFitModel::clearFitParVectors
void clearFitParVectors()
Clear the vectors containing fit parameters.
Definition: LauAbsFitModel.cc:378
LauCPFitModel::setAmpCoeffSet
void setAmpCoeffSet(LauAbsCoeffSet *coeffSet)
Set the DP amplitude coefficients.
Definition: LauCPFitModel.cc:337
LauParameter
Class for defining the fit parameter objects.
Definition: LauParameter.hh:31
LauEmbeddedData
Class to store the data for embedding in toy experiments.
Definition: LauEmbeddedData.hh:33
LauAbsDPDynamics::getEvtmPrime
virtual Double_t getEvtmPrime() const =0
Retrieve the square Dalitz plot coordinate, m&#39;, for the current event.
LauAbsFitModel::addSPlotNtupleDoubleBranch
virtual void addSPlotNtupleDoubleBranch(const TString &name)
Add a branch to the sPlot tree for storing a double.
Definition: LauAbsFitModel.cc:528
LauCPFitModel::nSigComp_
UInt_t nSigComp_
Number of signal components.
Definition: LauCPFitModel.hh:471
LauCPFitModel::setupSPlotNtupleBranches
virtual void setupSPlotNtupleBranches()
Add branches to store experiment number and the event number within the experiment.
Definition: LauCPFitModel.cc:2187
LauAbsFitModel::clearExtraVarVectors
void clearExtraVarVectors()
Clear the vectors containing extra ntuple variables.
Definition: LauAbsFitModel.cc:384
LauEmbeddedData::getEmbeddedEvent
void getEmbeddedEvent(LauKinematics *kinematics)
Retrieve an event from the data sample.
Definition: LauEmbeddedData.cc:140
LauCPFitModel::bkgndDPLike_
std::vector< Double_t > bkgndDPLike_
Background DP likelihood value(s)
Definition: LauCPFitModel.hh:617
LauCPFitModel::embedNegBkgnd
void embedNegBkgnd(const TString &bgClass, const TString &fileName, const TString &treeName, Bool_t reuseEventsWithinEnsemble, Bool_t reuseEventsWithinExperiment=kFALSE)
Embed full simulation events for the given background class, rather than generating toy from the PDFs...
Definition: LauCPFitModel.cc:2658
LauCPFitModel::negBkgndTree_
LauBkgndEmbDataList negBkgndTree_
The B- background event tree.
Definition: LauCPFitModel.hh:592
LauFitDataTree::appendFakePoints
void appendFakePoints(const std::vector< Double_t > &xCoords, const std::vector< Double_t > &yCoords)
Add fake events to the data.
Definition: LauFitDataTree.cc:182
LauAbsDPDynamics::checkToyMC
virtual ToyMCStatus checkToyMC(Bool_t printErrorMessages=kTRUE, Bool_t printInfoMessages=kFALSE)=0
Check the status of the toy MC generation.
LauAbsPdf.hh
File containing declaration of LauAbsPdf class.
LauAbsFitModel::storeDPEff
Bool_t storeDPEff() const
Determine whether the efficiency information should be stored in the sPlot ntuple.
Definition: LauAbsFitModel.hh:169
LauCPFitModel::reuseBkgnd_
LauBkgndReuseEventsList reuseBkgnd_
Vector of booleans to reuse background events.
Definition: LauCPFitModel.hh:607
LauCPFitModel::coeffPars_
std::vector< LauAbsCoeffSet * > coeffPars_
Magnitudes and Phases.
Definition: LauCPFitModel.hh:483
LauCPFitModel::posMeanEff_
LauParameter posMeanEff_
The mean efficiency for B+ model.
Definition: LauCPFitModel.hh:507
LauCPFitModel::negSigModel_
LauAbsDPDynamics * negSigModel_
The B- signal Dalitz plot model.
Definition: LauCPFitModel.hh:432
LauAbsFitModel::bkgndClassID
UInt_t bkgndClassID(const TString &className) const
The number assigned to a background class.
Definition: LauAbsFitModel.cc:348
LauFitDataTree::getData
const LauFitData & getData(UInt_t iEvt) const
Retrieve the data for a given event.
Definition: LauFitDataTree.cc:312
LauRandom.hh
File containing LauRandom namespace.
LauAbsDPDynamics::getnAmp
UInt_t getnAmp() const
Retrieve the number of amplitude components.
Definition: LauAbsDPDynamics.hh:269
LauAbsFitModel::useRandomInitFitPars
Bool_t useRandomInitFitPars() const
Determine whether the initial values of the fit parameters, in particular the isobar coefficient para...
Definition: LauAbsFitModel.hh:172
LauEffModel.hh
File containing declaration of LauEffModel class.
LauAbsFitModel::fitData
const LauFitDataTree * fitData() const
Access the data store.
Definition: LauAbsFitModel.hh:661
LauAbsDPDynamics.hh
File containing declaration of LauAbsDPDynamics class.
LauCPFitModel::CPCFitFrac_
LauParArray CPCFitFrac_
The CP conserving fit fraction.
Definition: LauCPFitModel.hh:495
LauCPFitModel::weightEvents
virtual void weightEvents(const TString &dataFileName, const TString &dataTreeName)
Weight events based on the DP model.
Definition: LauCPFitModel.cc:2748
LauAbsFitModel::validBkgndClass
Bool_t validBkgndClass(const TString &className) const
Check if the given background class is in the list.
Definition: LauAbsFitModel.cc:331
LauCPFitModel::usePosReweighting_
Bool_t usePosReweighting_
Boolean to use reweighting for B+.
Definition: LauCPFitModel.hh:604
LauCPFitModel::generateBkgndEvent
Bool_t generateBkgndEvent(UInt_t bgID)
Generate background event.
Definition: LauCPFitModel.cc:1537
LauAbsDPDynamics::getEvtm23Sq
virtual Double_t getEvtm23Sq() const =0
Retrieve the invariant mass squared of the second and third daughters in the current event...
LauCPFitModel::writeOutTable
virtual void writeOutTable(const TString &outputFile)
Write the fit results in latex table format.
Definition: LauCPFitModel.cc:1029
LauKinematics::getm12Sq
Double_t getm12Sq() const
Get the m12 invariant mass square.
Definition: LauKinematics.hh:170
LauAbsDPDynamics::getDPNorm
Double_t getDPNorm() const
Retrieve the normalisation factor for the log-likelihood function.
Definition: LauAbsDPDynamics.hh:275
LauEffModel
Class that implements the efficiency description across the signal Dalitz plot.
Definition: LauEffModel.hh:37
LauDaughters::getNameParent
TString getNameParent() const
Get name of the parent particle.
Definition: LauDaughters.cc:241
LauEmbeddedData::nUsedEvents
UInt_t nUsedEvents() const
Retrieve the number of events that have already been sampled.
Definition: LauEmbeddedData.hh:69
LauAbsDPDynamics::updateCoeffs
virtual void updateCoeffs(const std::vector< LauComplex > &coeffs)
Update the complex coefficients for the resonances.
Definition: LauAbsDPDynamics.cc:88
LauKinematics::getm13Sq
Double_t getm13Sq() const
Get the m13 invariant mass square.
Definition: LauKinematics.hh:180
LauFitNtuple::storeParsAndErrors
void storeParsAndErrors(const std::vector< LauParameter * > &fitVars, const std::vector< LauParameter > &extraVars)
Store parameters and their errors.
Definition: LauFitNtuple.cc:156
LauCPFitModel::updateCoeffs
virtual void updateCoeffs()
Update the coefficients.
Definition: LauCPFitModel.cc:2177
LauAbsFitModel::setSPlotNtupleIntegerBranchValue
virtual void setSPlotNtupleIntegerBranchValue(const TString &name, Int_t value)
Set the value of an integer branch in the sPlot tree.
Definition: LauAbsFitModel.cc:533
LauCPFitModel::setSPlotNtupleBranchValues
Double_t setSPlotNtupleBranchValues(LauPdfList *extraPdfs, const TString &prefix, UInt_t iEvt)
Set the branches for the sPlot ntuple with extra PDFs.
Definition: LauCPFitModel.cc:2302
LauParameter::initValue
Double_t initValue() const
The initial value of the parameter.
Definition: LauParameter.hh:168
LauCPFitModel::~LauCPFitModel
virtual ~LauCPFitModel()
Destructor.
Definition: LauCPFitModel.cc:93
LauCPFitModel::posCoeffs_
std::vector< LauComplex > posCoeffs_
The complex coefficients for B+.
Definition: LauCPFitModel.hh:582
LauCPFitModel::initialiseDPModels
virtual void initialiseDPModels()
Initialise the signal DP models.
Definition: LauCPFitModel.cc:492
LauCPFitModel::updateSigEvents
void updateSigEvents()
Update the signal events after Minuit sets background parameters.
Definition: LauCPFitModel.cc:1743
LauAbsDPDynamics::getEvtm13Sq
virtual Double_t getEvtm13Sq() const =0
Retrieve the invariant mass squared of the first and third daughters in the current event...
LauEmbeddedData::getReweightedEvent
Bool_t getReweightedEvent(LauAbsDPDynamics *dynamics)
Retrieve an event from the data sample, applying an accept/reject based on the given DP model...
Definition: LauEmbeddedData.cc:59
LauKinematics::getm12
Double_t getm12() const
Get the m12 invariant mass.
Definition: LauKinematics.hh:154
LauAbsFitModel::fillSPlotNtupleBranches
virtual void fillSPlotNtupleBranches()
Fill the sPlot tuple.
Definition: LauAbsFitModel.cc:543
LauCPFitModel::setupBkgndVectors
virtual void setupBkgndVectors()
Define the length of the background vectors.
Definition: LauCPFitModel.cc:106
LauConstants.hh
File containing LauConstants namespace.
LauCPFitModel::setSignalPdfs
void setSignalPdfs(LauAbsPdf *negPdf, LauAbsPdf *posPdf)
Set the signal PDFs.
Definition: LauCPFitModel.cc:263
LauCPFitModel::setFitNEvents
void setFitNEvents()
Set the initial yields.
Definition: LauCPFitModel.cc:574
LauCPFitModel::twodimPDFs
virtual LauSPlot::TwoDMap twodimPDFs() const
Returns the species and variables for all 2D PDFs in the fit.
Definition: LauCPFitModel.cc:2417
LauCPFitModel::freeSpeciesNames
virtual LauSPlot::NumbMap freeSpeciesNames() const
Returns the names and yields of species that are free in the fit.
Definition: LauCPFitModel.cc:2379
LauCPFitModel::scfFracHist_
LauEffModel * scfFracHist_
The histogram giving the DP-dependence of the SCF fraction.
Definition: LauCPFitModel.hh:552
LauPrint::printFormat
void printFormat(std::ostream &stream, Double_t value) const
Method to choose the printing format to a specified level of precision.
Definition: LauPrint.cc:32
LauFitNtuple
Class to store the results from the fit into an ntuple.
Definition: LauFitNtuple.hh:41
LauCPFitModel::signalAsym_
LauParameter signalAsym_
Signal asymmetry.
Definition: LauCPFitModel.hh:519
LauCPFitModel::nSigDPPar_
UInt_t nSigDPPar_
Number of signal DP parameters.
Definition: LauCPFitModel.hh:474
LauFitDataTree::haveBranch
Bool_t haveBranch(const TString &name) const
Check if the named branch is stored.
Definition: LauFitDataTree.cc:237
LauCPFitModel::propagateParUpdates
virtual void propagateParUpdates()
Definition: LauCPFitModel.cc:1728
LauScfMap
Class for representing the 4D smearing matrix for mis-reconstructed signal (self cross feed) ...
Definition: LauScfMap.hh:29
LauCPFitModel::storePerEvtLlhds
virtual void storePerEvtLlhds()
Store the per event likelihood values.
Definition: LauCPFitModel.cc:2482
LauCPFitModel::scfDPLike_
Double_t scfDPLike_
SCF DP likelihood value.
Definition: LauCPFitModel.hh:614
LauKinematics::getThetaPrime
Double_t getThetaPrime() const
Get theta&#39; value.
Definition: LauKinematics.hh:210
LauParameter::updatePull
void updatePull()
Call to update the bias and pull values.
Definition: LauParameter.cc:445
LauAbsFitModel::nBkgndClasses
UInt_t nBkgndClasses() const
Returns the number of background classes.
Definition: LauAbsFitModel.hh:206
LauAbsFitModel::bkgndClassName
const TString & bkgndClassName(UInt_t classID) const
Get the name of a background class from the number.
Definition: LauAbsFitModel.cc:366
LauCPFitModel::recoJacobians_
std::vector< Double_t > recoJacobians_
The cached values of the sqDP jacobians for each event.
Definition: LauCPFitModel.hh:564
LauAbsFitModel::enableEmbedding
Bool_t enableEmbedding() const
Determine whether embedding of events is enabled in the generation.
Definition: LauAbsFitModel.hh:118
LauCPFitModel::getEventSum
virtual Double_t getEventSum() const
Get the total number of events.
Definition: LauCPFitModel.cc:1966
LauAbsCoeffSet::name
virtual TString name() const
Retrieve the name of the coefficient set.
Definition: LauAbsCoeffSet.hh:102
LauAbsPdf
Class for defining the abstract interface for PDF classes.
Definition: LauAbsPdf.hh:40
LauAbsFitModel::fitPars
const LauParameterPList & fitPars() const
Access the fit variables.
Definition: LauAbsFitModel.hh:641
LauCPFitModel::posScfPdfs_
LauPdfList posScfPdfs_
The B+ SCF PDFs.
Definition: LauCPFitModel.hh:459
LauAbsDPDynamics::getEvtDPAmp
virtual const LauComplex & getEvtDPAmp() const =0
Retrieve the total amplitude of all amplitude components at the current point in the Dalitz plot...
LauCPFitModel::bkgndEvents_
LauBkgndYieldList bkgndEvents_
Background yield(s)
Definition: LauCPFitModel.hh:525
LauKinematics
Class for calculating 3-body kinematic quantities.
Definition: LauKinematics.hh:31
LauParameter::value
Double_t value() const
The value of the parameter.
Definition: LauParameter.hh:138
LauCPFitModel::negBkgndDPModels_
LauBkgndDPModelList negBkgndDPModels_
The B- background Dalitz plot models.
Definition: LauCPFitModel.hh:438
LauAbsDPDynamics::getFitFractions
const LauParArray & getFitFractions() const
Retrieve the fit fractions for the amplitude components.
Definition: LauAbsDPDynamics.hh:257
LauFitDataTree::nEvents
UInt_t nEvents() const
Retrieve the number of events.
Definition: LauFitDataTree.hh:102
LauCPFitModel::CPVFitFrac_
LauParArray CPVFitFrac_
The CP violating fit fraction.
Definition: LauCPFitModel.hh:492
LauAbsFitModel::addGenNtupleDoubleBranch
virtual void addGenNtupleDoubleBranch(const TString &name)
Add a branch to the gen tree for storing a double.
Definition: LauAbsFitModel.cc:493
LauCPFitModel::posSignalPdfs_
LauPdfList posSignalPdfs_
The B+ signal PDFs.
Definition: LauCPFitModel.hh:453
LauKinematics::calcSqDPJacobian
Double_t calcSqDPJacobian()
Calculate the Jacobian for the transformation m23^2, m13^2 -&gt; m&#39;, theta&#39; (square DP) ...
Definition: LauKinematics.cc:128
LauAbsFitModel::updateFitParameters
void updateFitParameters(LauPdfList &pdfList)
Update the fit parameters for all PDFs in the list.
Definition: LauAbsFitModel.cc:1266
LauCPFitModel::embedPosBkgnd
void embedPosBkgnd(const TString &bgClass, const TString &fileName, const TString &treeName, Bool_t reuseEventsWithinEnsemble, Bool_t reuseEventsWithinExperiment=kFALSE)
Embed full simulation events for the given background class, rather than generating toy from the PDFs...
Definition: LauCPFitModel.cc:2716
LauCPFitModel::fakeJacobians_
std::vector< Double_t > fakeJacobians_
The cached values of the sqDP jacobians for each true bin.
Definition: LauCPFitModel.hh:567
LauAbsDPDynamics::getEvtEff
virtual Double_t getEvtEff() const =0
Retrieve the efficiency for the current event.
LauAsymmCalc::getAsymmetry
Double_t getAsymmetry()
Retrieve the asymmetry.
Definition: LauAsymmCalc.hh:42
LauCPFitModel::getEvtSCFDPLikelihood
virtual Double_t getEvtSCFDPLikelihood(UInt_t iEvt)
Calculate the SCF likelihood for the DP for a given event.
Definition: LauCPFitModel.cc:2065
LauCPFitModel::signalEvents_
LauParameter signalEvents_
Signal yield.
Definition: LauCPFitModel.hh:516
LauAbsDPDynamics::usingScfModel
Bool_t usingScfModel()
Check whether a self cross feed fraction model is being used.
Definition: LauAbsDPDynamics.hh:329
LauEmbeddedData::findBranches
Bool_t findBranches()
Find and read the branches in data tree.
Definition: LauEmbeddedData.cc:41
LauAbsBkgndDPModel.hh
File containing declaration of LauAbsBkgndDPModel class.
LauCPFitModel::setSCFPdfs
void setSCFPdfs(LauAbsPdf *negPdf, LauAbsPdf *posPdf)
Set the SCF PDF for a given variable.
Definition: LauCPFitModel.cc:284
LauCPFitModel::LauGenInfo
std::map< std::pair< TString, Int_t >, std::pair< Int_t, Double_t > > LauGenInfo
Define a map to be used to store a category name and numbers.
Definition: LauCPFitModel.hh:220
LauCPFitModel::useSCFHist_
Bool_t useSCFHist_
Is the SCF fraction DP-dependent.
Definition: LauCPFitModel.hh:546
LauCPFitModel::negMeanEff_
LauParameter negMeanEff_
The mean efficiency for B- model.
Definition: LauCPFitModel.hh:504
LauCPFitModel::setBkgndPdfs
void setBkgndPdfs(const TString &bkgndClass, LauAbsPdf *negPdf, LauAbsPdf *posPdf)
Set the background PDFs.
Definition: LauCPFitModel.cc:305
LauScfMap::prob
Double_t prob(Int_t recoBin, Int_t trueBin) const
Probability of a true event in the given true bin migrating to the reco bin.
Definition: LauScfMap.cc:165
LauCPFitModel::getTotEvtLikelihood
virtual Double_t getTotEvtLikelihood(UInt_t iEvt)
Get the total likelihood for each event.
Definition: LauCPFitModel.cc:1900
LauCPFitModel::getEvtDPLikelihood
virtual void getEvtDPLikelihood(UInt_t iEvt)
Calculate the signal and background likelihoods for the DP for a given event.
Definition: LauCPFitModel.cc:1978
LauAbsDPDynamics::getEffModel
LauEffModel * getEffModel()
Retrieve the model for the efficiency across the Dalitz plot.
Definition: LauAbsDPDynamics.hh:311
LauParameter::genValue
Double_t genValue() const
The value generated for the parameter.
Definition: LauParameter.hh:162
LauCPFitModel::posBkgndTree_
LauBkgndEmbDataList posBkgndTree_
The B+ background event tree.
Definition: LauCPFitModel.hh:595
LauCPFitModel::setBkgndDPModels
void setBkgndDPModels(const TString &bkgndClass, LauAbsBkgndDPModel *negModel, LauAbsBkgndDPModel *posModel)
Set the background DP models.
Definition: LauCPFitModel.cc:242
LauCPFitModel::negBkgndPdfs_
LauBkgndPdfsList negBkgndPdfs_
The B- background PDFs.
Definition: LauCPFitModel.hh:462
LauAbsDPDynamics::initialise
virtual void initialise(const std::vector< LauComplex > &coeffs)=0
Initialise the Dalitz plot dynamics.
LauFitDataTree
Class to store the input fit variables.
Definition: LauFitDataTree.hh:37
LauCPFitModel::scfFrac_
LauParameter scfFrac_
The (global) SCF fraction parameter.
Definition: LauCPFitModel.hh:549
LauSPlot::NameSet
std::set< TString > NameSet
Type to store names, e.g. of the discriminating/control variables.
Definition: LauSPlot.hh:59
LauFitNtuple.hh
File containing declaration of LauFitNtuple class.
LauCPFitModel::setDPBranchValues
void setDPBranchValues()
Store all of the DP information.
Definition: LauCPFitModel.cc:1640
LauCPFitModel::posKinematics_
LauKinematics * posKinematics_
The B+ Dalitz plot kinematics object.
Definition: LauCPFitModel.hh:447
LauCPFitModel::embedPosSignal
void embedPosSignal(const TString &fileName, const TString &treeName, Bool_t reuseEventsWithinEnsemble, Bool_t reuseEventsWithinExperiment=kFALSE, Bool_t useReweighting=kFALSE)
Embed full simulation events for the B+ signal, rather than generating toy from the PDFs...
Definition: LauCPFitModel.cc:2690
LauCPFitModel::fitFracAsymm_
std::vector< LauParameter > fitFracAsymm_
The fit fraction asymmetries.
Definition: LauCPFitModel.hh:498
LauAbsBkgndDPModel::generate
virtual Bool_t generate()=0
Generate a toy MC event from the model.
LauCPFitModel::setExtraNtupleVars
void setExtraNtupleVars()
Set-up other parameters that are derived from the fit results, e.g. fit fractions.
Definition: LauCPFitModel.cc:622
LauAbsDPDynamics::generate
virtual Bool_t generate()=0
Generate a toy MC signal event.
Generated by   doxygen 1.8.5