laura is hosted by Hepforge, IPPP Durham
Laura++  3.6.0
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauCartesianGammaCPCoeffSet.cc
Go to the documentation of this file.
1 
2 /*
3 Copyright 2014 University of Warwick
4 
5 Licensed under the Apache License, Version 2.0 (the "License");
6 you may not use this file except in compliance with the License.
7 You may obtain a copy of the License at
8 
9  http://www.apache.org/licenses/LICENSE-2.0
10 
11 Unless required by applicable law or agreed to in writing, software
12 distributed under the License is distributed on an "AS IS" BASIS,
13 WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 See the License for the specific language governing permissions and
15 limitations under the License.
16 */
17 
18 /*
19 Laura++ package authors:
20 John Back
21 Paul Harrison
22 Thomas Latham
23 */
24 
30 
31 #include "LauComplex.hh"
32 #include "LauConstants.hh"
33 #include "LauParameter.hh"
34 #include "LauPrint.hh"
35 
36 #include "TMath.h"
37 #include "TRandom.h"
38 
39 #include <fstream>
40 #include <iostream>
41 #include <vector>
42 
44  const Double_t x,
45  const Double_t y,
46  const Double_t xCP,
47  const Double_t yCP,
48  const Double_t deltaXCP,
49  const Double_t deltaYCP,
50  const Bool_t xFixed,
51  const Bool_t yFixed,
52  const Bool_t xCPFixed,
53  const Bool_t yCPFixed,
54  const Bool_t deltaXCPFixed,
55  const Bool_t deltaYCPFixed,
56  const Bool_t deltaXCPSecondStage,
57  const Bool_t deltaYCPSecondStage ) :
58  LauAbsCoeffSet( compName ),
59  x_( new LauParameter( "X", x, minRealImagPart_, maxRealImagPart_, xFixed ) ),
60  y_( new LauParameter( "Y", y, minRealImagPart_, maxRealImagPart_, yFixed ) ),
61  xCP_( new LauParameter( "XCP", xCP, minRealImagPart_, maxRealImagPart_, xCPFixed ) ),
62  yCP_( new LauParameter( "YCP", yCP, minRealImagPart_, maxRealImagPart_, yCPFixed ) ),
63  deltaXCP_( new LauParameter( "DeltaXCP", deltaXCP, minDelta_, maxDelta_, deltaXCPFixed ) ),
64  deltaYCP_( new LauParameter( "DeltaYCP", deltaYCP, minDelta_, maxDelta_, deltaYCPFixed ) ),
65  nonCPPart_( x, y ),
66  cpPart_( 1 + xCP + deltaXCP, yCP + deltaYCP ),
67  cpAntiPart_( 1 + xCP - deltaXCP, yCP - deltaYCP ),
68  particleCoeff_( nonCPPart_ * cpPart_ ),
69  antiparticleCoeff_( nonCPPart_ * cpAntiPart_ ),
70  acp_( "ACP",
71  ( antiparticleCoeff_.abs2() - particleCoeff_.abs2() ) /
72  ( antiparticleCoeff_.abs2() + particleCoeff_.abs2() ),
73  -1.0,
74  1.0,
75  deltaXCPFixed && deltaYCPFixed )
76 {
77  if ( deltaXCPSecondStage && ! deltaXCPFixed ) {
78  deltaXCP_->secondStage( kTRUE );
79  deltaXCP_->initValue( 0.0 );
80  }
81  if ( deltaYCPSecondStage && ! deltaYCPFixed ) {
82  deltaYCP_->secondStage( kTRUE );
83  deltaYCP_->initValue( 0.0 );
84  }
85 }
86 
88  CloneOption cloneOption,
89  Double_t constFactor ) :
90  LauAbsCoeffSet( rhs.name() ),
91  x_( 0 ),
92  y_( 0 ),
93  xCP_( 0 ),
94  yCP_( 0 ),
95  deltaXCP_( 0 ),
96  deltaYCP_( 0 ),
97  nonCPPart_( rhs.nonCPPart_ ),
98  cpPart_( rhs.cpPart_ ),
99  cpAntiPart_( rhs.cpAntiPart_ ),
100  particleCoeff_( rhs.particleCoeff_ ),
101  antiparticleCoeff_( rhs.antiparticleCoeff_ ),
102  acp_( rhs.acp_ )
103 {
104  if ( cloneOption == All || cloneOption == TieRealPart ) {
105  x_ = rhs.x_->createClone( constFactor );
106  } else {
107  x_ = new LauParameter( "X",
108  rhs.x_->value(),
111  rhs.x_->fixed() );
112  if ( rhs.x_->blind() ) {
113  const LauBlind* blinder = rhs.x_->blinder();
115  }
116  }
117 
118  if ( cloneOption == All || cloneOption == TieImagPart ) {
119  y_ = rhs.y_->createClone( constFactor );
120  } else {
121  y_ = new LauParameter( "Y",
122  rhs.y_->value(),
125  rhs.y_->fixed() );
126  if ( rhs.y_->blind() ) {
127  const LauBlind* blinder = rhs.y_->blinder();
129  }
130  }
131 
132  if ( cloneOption == All || cloneOption == TieCPPars ) {
133  xCP_ = rhs.xCP_->createClone( constFactor );
134  yCP_ = rhs.yCP_->createClone( constFactor );
135  deltaXCP_ = rhs.deltaXCP_->createClone( constFactor );
136  deltaYCP_ = rhs.deltaYCP_->createClone( constFactor );
137  } else {
138  xCP_ = new LauParameter( "XCP",
139  rhs.xCP_->value(),
142  rhs.xCP_->fixed() );
143  if ( rhs.xCP_->blind() ) {
144  const LauBlind* blinder = rhs.xCP_->blinder();
146  }
147  yCP_ = new LauParameter( "YCP",
148  rhs.yCP_->value(),
151  rhs.yCP_->fixed() );
152  if ( rhs.yCP_->blind() ) {
153  const LauBlind* blinder = rhs.yCP_->blinder();
155  }
156  deltaXCP_ = new LauParameter( "DeltaXCP",
157  rhs.deltaXCP_->value(),
158  minDelta_,
159  maxDelta_,
160  rhs.deltaXCP_->fixed() );
161  deltaYCP_ = new LauParameter( "DeltaYCP",
162  rhs.deltaYCP_->value(),
163  minDelta_,
164  maxDelta_,
165  rhs.deltaYCP_->fixed() );
166  if ( rhs.deltaXCP_->secondStage() && ! rhs.deltaXCP_->fixed() ) {
167  deltaXCP_->secondStage( kTRUE );
168  deltaXCP_->initValue( 0.0 );
169  }
170  if ( rhs.deltaYCP_->secondStage() && ! rhs.deltaYCP_->fixed() ) {
171  deltaYCP_->secondStage( kTRUE );
172  deltaYCP_->initValue( 0.0 );
173  }
174  if ( rhs.deltaXCP_->blind() ) {
175  const LauBlind* blinder = rhs.deltaXCP_->blinder();
177  }
178  if ( rhs.deltaYCP_->blind() ) {
179  const LauBlind* blinder = rhs.deltaYCP_->blinder();
181  }
182  }
183 }
184 
185 std::vector<LauParameter*> LauCartesianGammaCPCoeffSet::getParameters()
186 {
187  std::vector<LauParameter*> pars;
188  pars.push_back( x_ );
189  pars.push_back( y_ );
190  if ( ! xCP_->fixed() )
191  pars.push_back( xCP_ );
192  if ( ! yCP_->fixed() )
193  pars.push_back( yCP_ );
194  if ( ! deltaXCP_->fixed() )
195  pars.push_back( deltaXCP_ );
196  if ( ! deltaYCP_->fixed() )
197  pars.push_back( deltaYCP_ );
198  return pars;
199 }
200 
202 {
203  std::cout << "INFO in LauCartesianGammaCPCoeffSet::printParValues : Component \""
204  << this->name() << "\" has ";
205  std::cout << "x = " << x_->value() << ",\t";
206  std::cout << "y = " << y_->value() << ",\t";
207  std::cout << "xCP = " << xCP_->value() << ",\t";
208  std::cout << "yCP = " << yCP_->value() << ",\t";
209  std::cout << "Delta xCP = " << deltaXCP_->value() << ",\t";
210  std::cout << "Delta yCP = " << deltaYCP_->value() << "." << std::endl;
211 }
212 
213 void LauCartesianGammaCPCoeffSet::printTableHeading( std::ostream& stream ) const
214 {
215  stream << "\\begin{tabular}{|l|c|c|c|c|c|c|}" << std::endl;
216  stream << "\\hline" << std::endl;
217  stream << "Component & Real Part & Imaginary Part & Real CP Part & Imaginary CP Part & $\\Delta$ Real CP Part & $\\Delta$ Imaginary CP Part \\\\"
218  << std::endl;
219  stream << "\\hline" << std::endl;
220 }
221 
222 void LauCartesianGammaCPCoeffSet::printTableRow( std::ostream& stream ) const
223 {
224  LauPrint print;
225  TString resName = this->name();
226  resName = resName.ReplaceAll( "_", "\\_" );
227  stream << resName << " & $";
228  print.printFormat( stream, x_->value() );
229  stream << " \\pm ";
230  print.printFormat( stream, x_->error() );
231  stream << "$ & $";
232  print.printFormat( stream, y_->value() );
233  stream << " \\pm ";
234  print.printFormat( stream, y_->error() );
235  stream << "$ & $";
236  print.printFormat( stream, xCP_->value() );
237  stream << " \\pm ";
238  print.printFormat( stream, xCP_->error() );
239  stream << "$ & $";
240  print.printFormat( stream, yCP_->value() );
241  stream << " \\pm ";
242  print.printFormat( stream, yCP_->error() );
243  stream << "$ & $";
244  print.printFormat( stream, deltaXCP_->value() );
245  stream << " \\pm ";
246  print.printFormat( stream, deltaXCP_->error() );
247  stream << "$ & $";
248  print.printFormat( stream, deltaYCP_->value() );
249  stream << " \\pm ";
250  print.printFormat( stream, deltaYCP_->error() );
251  stream << "$ \\\\" << std::endl;
252 }
253 
255 {
256  if ( x_->fixed() == kFALSE ) {
257  // Choose a value for "X" between -3.0 and 3.0
258  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm() * 6.0 - 3.0;
259  x_->initValue( value );
260  x_->value( value );
261  }
262  if ( y_->fixed() == kFALSE ) {
263  // Choose a value for "Y" between -3.0 and 3.0
264  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm() * 6.0 - 3.0;
265  y_->initValue( value );
266  y_->value( value );
267  }
268  if ( xCP_->fixed() == kFALSE ) {
269  // Choose a value for "XCP" between -3.0 and 3.0
270  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm() * 6.0 - 3.0;
271  xCP_->initValue( value );
272  xCP_->value( value );
273  }
274  if ( yCP_->fixed() == kFALSE ) {
275  // Choose a value for "YCP" between -3.0 and 3.0
276  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm() * 6.0 - 3.0;
277  yCP_->initValue( value );
278  yCP_->value( value );
279  }
280  if ( deltaXCP_->fixed() == kFALSE && deltaXCP_->secondStage() == kFALSE ) {
281  // Choose a value for "Delta XCP" between -0.5 and 0.5
282  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm() * 1.0 - 0.5;
284  deltaXCP_->value( value );
285  }
286  if ( deltaYCP_->fixed() == kFALSE && deltaYCP_->secondStage() == kFALSE ) {
287  // Choose a value for "Delta YCP" between -0.5 and 0.5
288  Double_t value = LauAbsCoeffSet::getRandomiser()->Rndm() * 1.0 - 0.5;
290  deltaYCP_->value( value );
291  }
292 }
293 
295 {
296  // update the pulls
297  x_->updatePull();
298  y_->updatePull();
299  xCP_->updatePull();
300  yCP_->updatePull();
303 }
304 
306 {
311  return particleCoeff_;
312 }
313 
315 {
320  return antiparticleCoeff_;
321 }
322 
324 {
325  std::cerr << "ERROR in LauCartesianGammaCPCoeffSet::setCoeffValues : Method not supported by this class - too many parameters"
326  << std::endl;
327 }
328 
330 {
331  // set the name
332  TString parName( this->baseName() );
333  parName += "_ACP";
334  acp_.name( parName );
335 
336  // work out the ACP value
337  const LauComplex nonCPPart( x_->value(), y_->value() );
338  const LauComplex cpPart( 1.0 + xCP_->value() + deltaXCP_->value(),
339  yCP_->value() + deltaYCP_->value() );
340  const LauComplex cpAntiPart( 1.0 + xCP_->value() - deltaXCP_->value(),
341  yCP_->value() - deltaYCP_->value() );
342  const LauComplex partCoeff = nonCPPart * cpPart;
343  const LauComplex antiCoeff = nonCPPart * cpAntiPart;
344 
345  const Double_t numer = antiCoeff.abs2() - partCoeff.abs2();
346  const Double_t denom = antiCoeff.abs2() + partCoeff.abs2();
347  const Double_t value = numer / denom;
348 
349  // is it fixed?
350  const Bool_t fixed = deltaXCP_->fixed() && deltaYCP_->fixed();
351  acp_.fixed( fixed );
352 
353  // we can't work out the error without the covariance matrix
354  const Double_t error( 0.0 );
355 
356  // set the value and error
358 
359  return acp_;
360 }
361 
363  CloneOption cloneOption,
364  Double_t constFactor )
365 {
366  LauAbsCoeffSet* clone( 0 );
367  if ( cloneOption == All || cloneOption == TieRealPart || cloneOption == TieImagPart ||
368  cloneOption == TieCPPars ) {
369  clone = new LauCartesianGammaCPCoeffSet( *this, cloneOption, constFactor );
370  clone->name( newName );
371  } else {
372  std::cerr << "ERROR in LauCartesianGammaCPCoeffSet::createClone : Invalid clone option"
373  << std::endl;
374  }
375  return clone;
376 }
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.
Class for defining the abstract interface for complex coefficient classes.
Class for defining the fit parameter objects.
Definition: LauParameter.hh:49
virtual void printTableRow(std::ostream &stream) const
Print the parameters of the complex coefficient as a row in the results table.
Double_t unblindValue() const
The unblinded value of the parameter.
virtual void printTableHeading(std::ostream &stream) const
Print the column headings for a results table.
LauComplex nonCPPart_
The nonCP part of the complex coefficient.
Double_t value() const
The value of the parameter.
LauComplex cpAntiPart_
The CP part of the complex coefficient for the antiparticle.
const LauBlind * blinder() const
Access the blinder object.
LauParameter * yCP_
The average CP imaginary part.
virtual TString name() const
Retrieve the name of the coefficient set.
File containing declaration of LauCartesianGammaCPCoeffSet class.
LauParameter * deltaXCP_
The asymmetric CP real part.
File containing declaration of LauParameter class.
LauParameter * x_
The nonCP real part.
virtual void finaliseValues()
Make sure values are in "standard" ranges, e.g. phases should be between -pi and pi.
virtual const LauComplex & particleCoeff()
Retrieve the complex coefficient for a particle.
LauComplex cpPart_
The CP part of the complex coefficient for the particle.
static Double_t minDelta_
Minimum allowed value of CP-violating real/imaginary part parameters.
virtual const TString & baseName() const
Retrieve the base name of the coefficient set.
Bool_t secondStage() const
Check whether the parameter should be floated only in the second stage of a two stage fit.
Class for defining a complex number.
Definition: LauComplex.hh:61
LauComplex particleCoeff_
The particle complex coefficient.
Bool_t blind() const
The blinding state.
Class for blinding and unblinding a number based on a blinding string.
Definition: LauBlind.hh:42
void setRealImagPart(Double_t realpart, Double_t imagpart)
Set both real and imaginary part.
Definition: LauComplex.hh:324
Bool_t clone() const
Check whether is a clone or not.
virtual void randomiseInitValues()
Randomise the starting values of the parameters for a fit.
File containing declaration of LauComplex class.
LauParameter * createClone(Double_t constFactor=1.0)
Method to create a clone from the parent parameter using the copy constructor.
void blindParameter(const TString &blindingString, const Double_t width)
Blind the parameter.
LauParameter acp_
The CP asymmetry.
virtual LauAbsCoeffSet * createClone(const TString &newName, CloneOption cloneOption=All, Double_t constFactor=1.0)
Create a clone of the coefficient set.
Double_t blindingWidth() const
Obtain the Gaussian width.
Definition: LauBlind.hh:82
Double_t abs2() const
Obtain the square of the absolute value of the complex number.
Definition: LauComplex.hh:260
virtual std::vector< LauParameter * > getParameters()
Retrieve the parameters of the coefficient, e.g. so that they can be loaded into a fit.
Bool_t fixed() const
Check whether the parameter is fixed or floated.
LauParameter()
Default constructor.
Definition: LauParameter.cc:40
virtual void setCoeffValues(const LauComplex &coeff, const LauComplex &coeffBar, Bool_t init)
Set the parameters based on the complex coefficients for particles and antiparticles.
virtual void printParValues() const
Print the current values of the parameters.
LauParameter * xCP_
The average CP real part.
File containing LauConstants namespace.
const TString & name() const
The parameter name.
const TString & blindingString() const
Obtain the blinding string.
Definition: LauBlind.hh:76
Double_t error() const
The error on the parameter.
static TRandom * getRandomiser()
Access the randomiser.
LauParameter * y_
The nonCP imaginary part.
void updatePull()
Call to update the bias and pull values.
LauParameter * deltaYCP_
The asymmetric CP imaginary part.
File containing declaration of LauPrint class and LauOutputLevel enum.
CloneOption
Options for cloning operation.
static Double_t maxRealImagPart_
Maximum allowed value of real/imaginary part parameters.
virtual const LauComplex & antiparticleCoeff()
Retrieve the complex coefficient for an antiparticle.
LauCartesianGammaCPCoeffSet(const TString &compName, const Double_t x, const Double_t y, const Double_t xCP, const Double_t yCP, const Double_t deltaXCP, const Double_t deltaYCP, const Bool_t xFixed, const Bool_t yFixed, const Bool_t xCPFixed, const Bool_t yCPFixed, const Bool_t deltaXCPFixed, const Bool_t deltaYCPFixed, const Bool_t deltaXCPSecondStage=kFALSE, const Bool_t deltaYCPSecondStage=kFALSE)
Constructor.
Class to define various output print commands.
Definition: LauPrint.hh:54
LauComplex antiparticleCoeff_
The antiparticle complex coefficient.
Class for defining a complex coefficient using the Cartesian gamma CP convention.
static Double_t minRealImagPart_
Minimum allowed value of real/imaginary part parameters.
void printFormat(std::ostream &stream, Double_t value) const
Method to choose the printing format to a specified level of precision.
Definition: LauPrint.cc:43
virtual LauParameter acp()
Calculate the CP asymmetry.
static Double_t maxDelta_
Maximum allowed value of CP-violating real/imaginary part parameters.
Double_t initValue() const
The initial value of the parameter.