laura is hosted by Hepforge, IPPP Durham
Laura++  v3r3
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauRealImagGammaCPCoeffSet.cc
Go to the documentation of this file.
1 
2 // Copyright University of Warwick 2013 - 2014.
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 <fstream>
17 #include <vector>
18 
19 #include "TMath.h"
20 #include "TRandom.h"
21 
23 #include "LauComplex.hh"
24 #include "LauConstants.hh"
25 #include "LauParameter.hh"
26 #include "LauPrint.hh"
27 
29 
30 
31 LauRealImagGammaCPCoeffSet::LauRealImagGammaCPCoeffSet(const TString& compName, const Double_t x, const Double_t y, const Double_t xCP, const Double_t yCP, const Double_t xbarCP, const Double_t ybarCP,
32  const Bool_t xFixed, const Bool_t yFixed, const Bool_t xCPFixed, const Bool_t yCPFixed, const Bool_t xbarCPFixed, const Bool_t ybarCPFixed) :
33  LauAbsCoeffSet(compName),
34  x_(new LauParameter("X", x, minRealImagPart_, maxRealImagPart_, xFixed)),
35  y_(new LauParameter("Y", y, minRealImagPart_, maxRealImagPart_, yFixed)),
36  xCP_(new LauParameter("XCP", xCP, minRealImagPart_, maxRealImagPart_, xCPFixed)),
37  yCP_(new LauParameter("YCP", yCP, minRealImagPart_, maxRealImagPart_, yCPFixed)),
38  xbarCP_(new LauParameter("XbarCP", xbarCP, minRealImagPart_, maxRealImagPart_, xbarCPFixed)),
39  ybarCP_(new LauParameter("YbarCP", ybarCP, minRealImagPart_, maxRealImagPart_, ybarCPFixed)),
40  nonCPPart_( x, y),
41  cpPart_( 1+xCP, yCP ),
42  cpAntiPart_( 1+xbarCP, ybarCP ),
43  particleCoeff_( nonCPPart_ * cpPart_ ),
44  antiparticleCoeff_( nonCPPart_ * cpAntiPart_ ),
45  acp_("ACP", (antiparticleCoeff_.abs2()-particleCoeff_.abs2())/(antiparticleCoeff_.abs2()+particleCoeff_.abs2()), -1.0, 1.0, xCPFixed&&yCPFixed&&xbarCPFixed&&ybarCPFixed)
46 {
47 }
48 
50  x_(0),
51  y_(0),
52  xCP_(0),
53  yCP_(0),
54  xbarCP_(0),
55  ybarCP_(0),
56  nonCPPart_( rhs.nonCPPart_ ),
57  cpPart_( rhs.cpPart_ ),
58  cpAntiPart_( rhs.cpAntiPart_ ),
59  particleCoeff_( rhs.particleCoeff_ ),
60  antiparticleCoeff_( rhs.antiparticleCoeff_ ),
61  acp_( rhs.acp_ )
62 {
63  if ( cloneOption == All || cloneOption == TieRealPart ) {
64  x_ = rhs.x_->createClone(constFactor);
65  } else {
66  x_ = new LauParameter("X", rhs.x_->value(), minRealImagPart_, maxRealImagPart_, rhs.x_->fixed());
67  if ( rhs.x_->blind() ) {
68  const LauBlind* blinder = rhs.x_->blinder();
69  x_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
70  }
71  }
72 
73  if ( cloneOption == All || cloneOption == TieImagPart ) {
74  y_ = rhs.y_->createClone(constFactor);
75  } else {
76  y_ = new LauParameter("Y", rhs.y_->value(), minRealImagPart_, maxRealImagPart_, rhs.y_->fixed());
77  if ( rhs.y_->blind() ) {
78  const LauBlind* blinder = rhs.y_->blinder();
79  y_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
80  }
81  }
82 
83  if ( cloneOption == All || cloneOption == TieCPPars ) {
84  xCP_ = rhs.xCP_->createClone(constFactor);
85  yCP_ = rhs.yCP_->createClone(constFactor);
86  xbarCP_ = rhs.xbarCP_->createClone(constFactor);
87  ybarCP_ = rhs.ybarCP_->createClone(constFactor);
88  } else {
89  xCP_ = new LauParameter("XCP", rhs.xCP_->value(), minRealImagPart_, maxRealImagPart_, rhs.xCP_->fixed());
90  if ( rhs.xCP_->blind() ) {
91  const LauBlind* blinder = rhs.xCP_->blinder();
92  xCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
93  }
94  yCP_ = new LauParameter("YCP", rhs.yCP_->value(), minRealImagPart_, maxRealImagPart_, rhs.yCP_->fixed());
95  if ( rhs.yCP_->blind() ) {
96  const LauBlind* blinder = rhs.yCP_->blinder();
97  yCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
98  }
100  if ( rhs.xbarCP_->blind() ) {
101  const LauBlind* blinder = rhs.xbarCP_->blinder();
102  xbarCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
103  }
105  if ( rhs.ybarCP_->blind() ) {
106  const LauBlind* blinder = rhs.ybarCP_->blinder();
107  ybarCP_->blindParameter( blinder->blindingString(), blinder->blindingWidth() );
108  }
109  }
110 }
111 
112 std::vector<LauParameter*> LauRealImagGammaCPCoeffSet::getParameters()
113 {
114  std::vector<LauParameter*> pars;
115  pars.push_back(x_);
116  pars.push_back(y_);
117  if(!xCP_->fixed()) pars.push_back(xCP_);
118  if(!yCP_->fixed()) pars.push_back(yCP_);
119  if(!xbarCP_->fixed()) pars.push_back(xbarCP_);
120  if(!ybarCP_->fixed()) pars.push_back(ybarCP_);
121  return pars;
122 }
123 
125 {
126  std::cout << "INFO in LauRealImagGammaCPCoeffSet::printParValues : Component \"" << this->name() << "\" has ";
127  std::cout << "x = " << x_->value() << ",\t";
128  std::cout << "y = " << y_->value() << ",\t";
129  std::cout << "xCP = " << xCP_->value() << ",\t";
130  std::cout << "yCP = " << yCP_->value() << ",\t";
131  std::cout << "xbarCP = " << xbarCP_->value() << ",\t";
132  std::cout << "ybarCP = " << ybarCP_->value() << "." << std::endl;
133 }
134 
135 void LauRealImagGammaCPCoeffSet::printTableHeading(std::ostream& stream) const
136 {
137  stream<<"\\begin{tabular}{|l|c|c|c|c|c|c|}"<<std::endl;
138  stream<<"\\hline"<<std::endl;
139  stream<<"Component & Real Part & Imaginary Part & Particle CP Real Part & Particle CP Imaginary Part & Antiparticle CP Real Part & Antiparticle CP Imaginary Part \\\\"<<std::endl;
140  stream<<"\\hline"<<std::endl;
141 }
142 
143 void LauRealImagGammaCPCoeffSet::printTableRow(std::ostream& stream) const
144 {
145  LauPrint print;
146  TString resName = this->name();
147  resName = resName.ReplaceAll("_", "\\_");
148  stream<<resName<<" & $";
149  print.printFormat(stream, x_->value());
150  stream<<" \\pm ";
151  print.printFormat(stream, x_->error());
152  stream<<"$ & $";
153  print.printFormat(stream, y_->value());
154  stream<<" \\pm ";
155  print.printFormat(stream, y_->error());
156  stream<<"$ & $";
157  print.printFormat(stream, xCP_->value());
158  stream<<" \\pm ";
159  print.printFormat(stream, xCP_->error());
160  stream<<"$ & $";
161  print.printFormat(stream, yCP_->value());
162  stream<<" \\pm ";
163  print.printFormat(stream, yCP_->error());
164  stream<<"$ & $";
165  print.printFormat(stream, xbarCP_->value());
166  stream<<" \\pm ";
167  print.printFormat(stream, xbarCP_->error());
168  stream<<"$ & $";
169  print.printFormat(stream, ybarCP_->value());
170  stream<<" \\pm ";
171  print.printFormat(stream, ybarCP_->error());
172  stream<<"$ \\\\"<<std::endl;
173 }
174 
176 {
177  if (x_->fixed() == kFALSE) {
178  // Choose a value for "X" between -3.0 and 3.0
179  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
180  x_->initValue(value); x_->value(value);
181  }
182  if (y_->fixed() == kFALSE) {
183  // Choose a value for "Y" between -3.0 and 3.0
184  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
185  y_->initValue(value); y_->value(value);
186  }
187  if (xCP_->fixed() == kFALSE) {
188  // Choose a value for "XCP" between -3.0 and 3.0
189  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
190  xCP_->initValue(value); xCP_->value(value);
191  }
192  if (yCP_->fixed() == kFALSE) {
193  // Choose a value for "YCP" between -3.0 and 3.0
194  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
195  yCP_->initValue(value); yCP_->value(value);
196  }
197  if (xbarCP_->fixed() == kFALSE) {
198  // Choose a value for "XbarCP" between -3.0 and 3.0
199  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
200  xbarCP_->initValue(value); xbarCP_->value(value);
201  }
202  if (ybarCP_->fixed() == kFALSE) {
203  // Choose a value for "YbarCP" between -3.0 and 3.0
204  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm()*6.0 - 3.0;
205  ybarCP_->initValue(value); ybarCP_->value(value);
206  }
207 }
208 
210 {
211  // update the pulls
212  x_->updatePull();
213  y_->updatePull();
214  xCP_->updatePull();
215  yCP_->updatePull();
216  xbarCP_->updatePull();
217  ybarCP_->updatePull();
218 }
219 
221 {
225  return particleCoeff_;
226 }
227 
229 {
233  return antiparticleCoeff_;
234 }
235 
237 {
238  std::cerr << "ERROR in LauCartesianGammaCPCoeffSet::setCoeffValues : Method not supported by this class - too many parameters" << std::endl;
239 }
240 
242 {
243  // set the name
244  TString parName(this->baseName()); parName += "_ACP";
245  acp_.name(parName);
246 
247  // work out the ACP value
248  const LauComplex nonCPPart( x_->value(), y_->value() );
249  const LauComplex cpPart( 1.0+xCP_->value(), yCP_->value() );
250  const LauComplex cpAntiPart( 1.0+xbarCP_->value(), ybarCP_->value() );
251  const LauComplex partCoeff = nonCPPart * cpPart;
252  const LauComplex antiCoeff = nonCPPart * cpAntiPart;
253 
254  const Double_t numer = antiCoeff.abs2() - partCoeff.abs2();
255  const Double_t denom = antiCoeff.abs2() + partCoeff.abs2();
256  const Double_t value = numer/denom;
257 
258  // is it fixed?
259  const Bool_t fixed = xCP_->fixed() && yCP_->fixed() && xbarCP_->fixed() && ybarCP_->fixed();
260  acp_.fixed(fixed);
261 
262  // we can't work out the error without the covariance matrix
263  const Double_t error(0.0);
264 
265  // set the value and error
266  acp_.valueAndErrors(value,error);
267 
268  return acp_;
269 }
270 
271 LauAbsCoeffSet* LauRealImagGammaCPCoeffSet::createClone(const TString& newName, CloneOption cloneOption, Double_t constFactor)
272 {
273  LauAbsCoeffSet* clone(0);
274  if ( cloneOption == All || cloneOption == TieRealPart || cloneOption == TieImagPart || cloneOption == TieCPPars ) {
275  clone = new LauRealImagGammaCPCoeffSet( *this, cloneOption, constFactor );
276  clone->name( newName );
277  } else {
278  std::cerr << "ERROR in LauRealImagGammaCPCoeffSet::createClone : Invalid clone option" << std::endl;
279  }
280  return clone;
281 }
282 
LauRealImagGammaCPCoeffSet(const TString &compName, const Double_t x, const Double_t y, const Double_t xCP, const Double_t yCP, const Double_t xbarCP, const Double_t ybarCP, const Bool_t xFixed, const Bool_t yFixed, const Bool_t xCPFixed, const Bool_t yCPFixed, const Bool_t xbarCPFixed, const Bool_t ybarCPFixed)
Constructor.
Bool_t fixed() const
Check whether the parameter is fixed or floated.
virtual LauAbsCoeffSet * createClone(const TString &newName, CloneOption cloneOption=All, Double_t constFactor=1.0)
Create a clone of the coefficient set.
LauParameter * xbarCP_
The real CP part for the antiparticle.
LauComplex cpAntiPart_
The CP part of the complex coefficient for the antiparticle.
ClassImp(LauAbsCoeffSet)
virtual void printParValues() const
Print the current values of the parameters.
LauParameter()
Default constructor.
Definition: LauParameter.cc:30
const TString & name() const
The parameter name.
Class for defining a complex coefficient using a Cartesian nonCP part multiplied by a simple Cartesia...
static Double_t maxRealImagPart_
Maximum allowed value of real/imaginary part parameters.
virtual void printTableRow(std::ostream &stream) const
Print the parameters of the complex coefficient as a row in the results table.
File containing declaration of LauPrint class.
Class to define various output print commands.
Definition: LauPrint.hh:29
LauParameter * ybarCP_
The imaginary CP part for the antiparticle.
virtual const LauComplex & particleCoeff()
Retrieve the complex coefficient for a particle.
LauComplex particleCoeff_
The particle complex coefficient.
LauComplex antiparticleCoeff_
The antiparticle complex coefficient.
CloneOption
Options for cloning operation.
Double_t abs2() const
Obtain the square of the absolute value of the complex number.
Definition: LauComplex.hh:232
Bool_t clone() const
Check whether is a clone or not.
LauComplex cpPart_
The CP part of the complex coefficient for the particle.
Bool_t blind() const
The blinding state.
LauParameter * yCP_
The imaginary CP part for the particle.
LauParameter * xCP_
The real CP part for the particle.
File containing declaration of LauParameter class.
File containing declaration of LauComplex class.
const TString & blindingString() const
Obtain the blinding string.
Definition: LauBlind.hh:62
Double_t error() const
The error on the parameter.
Class for defining the abstract interface for complex coefficient classes.
virtual void setCoeffValues(const LauComplex &coeff, const LauComplex &coeffBar, Bool_t init)
Set the parameters based on the complex coefficients for particles and antiparticles.
File containing declaration of LauRealImagGammaCPCoeffSet class.
Class for defining the fit parameter objects.
Definition: LauParameter.hh:35
void valueAndErrors(Double_t newValue, Double_t newError, Double_t newNegError=0.0, Double_t newPosError=0.0)
Set the value and errors on the parameter.
virtual const LauComplex & antiparticleCoeff()
Retrieve the complex coefficient for an antiparticle.
virtual std::vector< LauParameter * > getParameters()
Retrieve the parameters of the coefficient, e.g. so that they can be loaded into a fit...
static Double_t minRealImagPart_
Minimum allowed value of real/imaginary part parameters.
LauParameter * y_
The imaginary nonCP part.
const LauBlind * blinder() const
Access the blinder object.
LauParameter acp_
The CP asymmetry.
void setRealImagPart(Double_t realpart, Double_t imagpart)
Set both real and imaginary part.
Definition: LauComplex.hh:314
virtual void printTableHeading(std::ostream &stream) const
Print the column headings for a results table.
Double_t initValue() const
The initial value of the parameter.
void blindParameter(const TString &blindingString, const Double_t width)
Blind the parameter.
File containing LauConstants namespace.
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
Double_t unblindValue() const
The unblinded value of the parameter.
Class for defining a complex number.
Definition: LauComplex.hh:47
void updatePull()
Call to update the bias and pull values.
LauParameter * createClone(Double_t constFactor=1.0)
Method to create a clone from the parent parameter using the copy constructor.
virtual TString name() const
Retrieve the name of the coefficient set.
LauParameter * x_
The real nonCP part.
Double_t value() const
The value of the parameter.
Double_t blindingWidth() const
Obtain the Gaussian width.
Definition: LauBlind.hh:68
virtual const TString & baseName() const
Retrieve the base name of the coefficient set.
virtual void randomiseInitValues()
Randomise the starting values of the parameters for a fit.
Class for blinding and unblinding a number based on a blinding string.
Definition: LauBlind.hh:28
static TRandom * getRandomiser()
Access the randomiser.
LauComplex nonCPPart_
The nonCP part of the complex coefficient.
virtual void finaliseValues()
Make sure values are in &quot;standard&quot; ranges, e.g. phases should be between -pi and pi.
virtual LauParameter acp()
Calculate the CP asymmetry.