laura is hosted by Hepforge, IPPP Durham
Laura++  v3r3
A maximum likelihood fitting package for performing Dalitz-plot analysis.
LauResonanceMaker.cc
Go to the documentation of this file.
1 
2 // Copyright University of Warwick 2004 - 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 
17 #include "LauAbsResonance.hh"
18 #include "LauBelleNR.hh"
19 #include "LauBelleSymNR.hh"
20 #include "LauBreitWignerRes.hh"
21 #include "LauDabbaRes.hh"
22 #include "LauDaughters.hh"
23 #include "LauEFKLLMRes.hh"
24 #include "LauFlatteRes.hh"
25 #include "LauFlatNR.hh"
26 #include "LauGaussIncohRes.hh"
27 #include "LauGounarisSakuraiRes.hh"
28 #include "LauKappaRes.hh"
29 #include "LauLASSRes.hh"
30 #include "LauLASSBWRes.hh"
31 #include "LauLASSNRRes.hh"
34 #include "LauNRAmplitude.hh"
35 #include "LauPolNR.hh"
36 #include "LauRelBreitWignerRes.hh"
37 #include "LauResonanceInfo.hh"
38 #include "LauResonanceMaker.hh"
39 #include "LauRhoOmegaMix.hh"
40 #include "LauSigmaRes.hh"
41 
43 
44 
46 
47 
49  nResDefMax_(0),
50  bwBarrierType_(LauBlattWeisskopfFactor::BWPrimeBarrier),
51  bwRestFrame_(LauBlattWeisskopfFactor::ResonanceFrame),
52  spinFormalism_(LauAbsResonance::Zemach_P),
53  summaryPrinted_(kFALSE)
54 {
55  this->createResonanceVector();
57 }
58 
60 {
61  for ( std::vector<LauBlattWeisskopfFactor*>::iterator iter = bwIndepFactors_.begin(); iter != bwIndepFactors_.end(); ++iter ) {
62  delete *iter;
63  }
64  bwIndepFactors_.clear();
65  for ( BWFactorCategoryMap::iterator iter = bwFactors_.begin(); iter != bwFactors_.end(); ++iter ) {
66  delete iter->second.bwFactor_;
67  }
68  bwFactors_.clear();
69 }
70 
72 {
73  if ( resonanceMaker_ == 0 ) {
75  }
76 
77  return *resonanceMaker_;
78 }
79 
81 {
82  // Function to create all possible resonances that this class supports.
83  // Also add in the sigma and kappa - but a special paramterisation is used
84  // instead of the PDG "pole mass and width" values.
85 
86  std::cout << "INFO in LauResonanceMaker::createResonanceVector : Setting up possible resonance states..." << std::endl;
87 
88  LauResonanceInfo* neutral(0);
89  LauResonanceInfo* positve(0);
90  LauResonanceInfo* negatve(0);
91 
92  // Define the resonance names and store them in the array
93  resInfo_.clear();
94  resInfo_.reserve(100);
95  // rho resonances name, mass, width, spin, charge, default BW category, BW radius parameter (defaults to 4.0)
96  // rho(770)
97  neutral = new LauResonanceInfo("rho0(770)", 0.77526, 0.1478, 1, 0, LauBlattWeisskopfFactor::Light, 5.3);
98  positve = new LauResonanceInfo("rho+(770)", 0.77511, 0.1491, 1, 1, LauBlattWeisskopfFactor::Light, 5.3);
99  negatve = positve->createChargeConjugate();
100  resInfo_.push_back( neutral );
101  resInfo_.push_back( positve );
102  resInfo_.push_back( negatve );
103  // The following two lines of code are placed here in order to allow the following, rather niche, scenario:
104  // The LauRhoOmegaMix code permits (through the use of the optional independentPar argument of LauResonanceInfo::addExtraParameter) the magnitude and phase of the rho/omega mixing to potentially differ between the decay of the parent particle to rho0 X and the parent antiparticle to rho0 Xbar.
105  // This can be acheived by using the rho0(770) record in one case and the rho0(770)_COPY record in the other.
106  neutral = neutral->createSharedParameterRecord("rho0(770)_COPY");
107  resInfo_.push_back( neutral );
108  // rho(1450)
109  neutral = new LauResonanceInfo("rho0(1450)", 1.465, 0.400, 1, 0, LauBlattWeisskopfFactor::Light );
110  positve = new LauResonanceInfo("rho+(1450)", 1.465, 0.400, 1, 1, LauBlattWeisskopfFactor::Light );
111  negatve = positve->createChargeConjugate();
112  resInfo_.push_back( neutral );
113  resInfo_.push_back( positve );
114  resInfo_.push_back( negatve );
115  // rho_3(1690)
116  neutral = new LauResonanceInfo("rho0_3(1690)", 1.686, 0.186, 3, 0, LauBlattWeisskopfFactor::Light );
117  positve = new LauResonanceInfo("rho+_3(1690)", 1.686, 0.186, 3, 1, LauBlattWeisskopfFactor::Light );
118  negatve = positve->createChargeConjugate();
119  resInfo_.push_back( neutral );
120  resInfo_.push_back( positve );
121  resInfo_.push_back( negatve );
122  // rho(1700)
123  neutral = new LauResonanceInfo("rho0(1700)", 1.720, 0.250, 1, 0, LauBlattWeisskopfFactor::Light );
124  positve = new LauResonanceInfo("rho+(1700)", 1.720, 0.250, 1, 1, LauBlattWeisskopfFactor::Light );
125  negatve = positve->createChargeConjugate();
126  resInfo_.push_back( neutral );
127  resInfo_.push_back( positve );
128  resInfo_.push_back( negatve );
129  // rho(1900)
130  neutral = new LauResonanceInfo("rho0(1900)", 1.909, 0.130, 1, 0, LauBlattWeisskopfFactor::Light );
131  positve = new LauResonanceInfo("rho+(1900)", 1.909, 0.130, 1, 1, LauBlattWeisskopfFactor::Light );
132  negatve = positve->createChargeConjugate();
133  resInfo_.push_back( neutral );
134  resInfo_.push_back( positve );
135  resInfo_.push_back( negatve );
136  // rho_3(1990)
137  neutral = new LauResonanceInfo("rho0_3(1990)", 1.982, 0.188, 3, 0, LauBlattWeisskopfFactor::Light );
138  positve = new LauResonanceInfo("rho+_3(1990)", 1.982, 0.188, 3, 1, LauBlattWeisskopfFactor::Light );
139  negatve = positve->createChargeConjugate();
140  resInfo_.push_back( neutral );
141  resInfo_.push_back( positve );
142  resInfo_.push_back( negatve );
143 
144  // K* resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
145  // K*(892)
146  neutral = new LauResonanceInfo("K*0(892)", 0.89581, 0.0474, 1, 0, LauBlattWeisskopfFactor::Kstar, 3.0);
147  positve = new LauResonanceInfo("K*+(892)", 0.89166, 0.0508, 1, 1, LauBlattWeisskopfFactor::Kstar, 3.0);
148  negatve = positve->createChargeConjugate();
149  resInfo_.push_back( neutral );
150  resInfo_.push_back( positve );
151  resInfo_.push_back( negatve );
152  // K*(1410)
153  neutral = new LauResonanceInfo("K*0(1410)", 1.414, 0.232, 1, 0, LauBlattWeisskopfFactor::Kstar );
154  positve = new LauResonanceInfo("K*+(1410)", 1.414, 0.232, 1, 1, LauBlattWeisskopfFactor::Kstar );
155  negatve = positve->createChargeConjugate();
156  resInfo_.push_back( neutral );
157  resInfo_.push_back( positve );
158  resInfo_.push_back( negatve );
159  // K*_0(1430)
160  neutral = new LauResonanceInfo("K*0_0(1430)", 1.425, 0.270, 0, 0, LauBlattWeisskopfFactor::Kstar );
161  positve = new LauResonanceInfo("K*+_0(1430)", 1.425, 0.270, 0, 1, LauBlattWeisskopfFactor::Kstar );
162  negatve = positve->createChargeConjugate();
163  resInfo_.push_back( neutral );
164  resInfo_.push_back( positve );
165  resInfo_.push_back( negatve );
166  // LASS nonresonant model
167  neutral = neutral->createSharedParameterRecord("LASSNR0");
168  positve = positve->createSharedParameterRecord("LASSNR+");
169  negatve = negatve->createSharedParameterRecord("LASSNR-");
170  resInfo_.push_back( neutral );
171  resInfo_.push_back( positve );
172  resInfo_.push_back( negatve );
173  // K*_2(1430)
174  neutral = new LauResonanceInfo("K*0_2(1430)", 1.4324, 0.109, 2, 0, LauBlattWeisskopfFactor::Kstar );
175  positve = new LauResonanceInfo("K*+_2(1430)", 1.4256, 0.0985, 2, 1, LauBlattWeisskopfFactor::Kstar );
176  negatve = positve->createChargeConjugate();
177  resInfo_.push_back( neutral );
178  resInfo_.push_back( positve );
179  resInfo_.push_back( negatve );
180  // K*(1680)
181  neutral = new LauResonanceInfo("K*0(1680)", 1.717, 0.322, 1, 0, LauBlattWeisskopfFactor::Kstar );
182  positve = new LauResonanceInfo("K*+(1680)", 1.717, 0.322, 1, 1, LauBlattWeisskopfFactor::Kstar );
183  negatve = positve->createChargeConjugate();
184  resInfo_.push_back( neutral );
185  resInfo_.push_back( positve );
186  resInfo_.push_back( negatve );
187  // K*(1950)
188  neutral = new LauResonanceInfo("K*0_0(1950)", 1.945, 0.201, 0, 0, LauBlattWeisskopfFactor::Kstar );
189  positve = new LauResonanceInfo("K*+_0(1950)", 1.945, 0.201, 0, 1, LauBlattWeisskopfFactor::Kstar );
190  negatve = positve->createChargeConjugate();
191  resInfo_.push_back( neutral );
192  resInfo_.push_back( positve );
193  resInfo_.push_back( negatve );
194 
195  // phi resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
196  // phi(1020)
197  neutral = new LauResonanceInfo("phi(1020)", 1.019461, 0.004266, 1, 0, LauBlattWeisskopfFactor::Light );
198  resInfo_.push_back( neutral );
199  // phi(1680)
200  neutral = new LauResonanceInfo("phi(1680)", 1.680, 0.150, 1, 0, LauBlattWeisskopfFactor::Light );
201  resInfo_.push_back( neutral );
202 
203  // f resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
204  // f_0(980)
205  neutral = new LauResonanceInfo("f_0(980)", 0.990, 0.070, 0, 0, LauBlattWeisskopfFactor::Light );
206  resInfo_.push_back( neutral );
207  // f_2(1270)
208  neutral = new LauResonanceInfo("f_2(1270)", 1.2751, 0.1851, 2, 0, LauBlattWeisskopfFactor::Light );
209  resInfo_.push_back( neutral );
210  // f_0(1370)
211  neutral = new LauResonanceInfo("f_0(1370)", 1.370, 0.350, 0, 0, LauBlattWeisskopfFactor::Light );
212  resInfo_.push_back( neutral );
213  // f'_0(1300), from Belle's Kspipi paper
214  neutral = new LauResonanceInfo("f'_0(1300)", 1.449, 0.126, 0, 0, LauBlattWeisskopfFactor::Light );
215  resInfo_.push_back( neutral );
216  // f_2(1430)
217  neutral = new LauResonanceInfo("f_2(1430)", 1.430, 0.150, 2, 0, LauBlattWeisskopfFactor::Light ); // PDG width in the range 13 - 150
218  resInfo_.push_back( neutral );
219  // f_0(1500)
220  neutral = new LauResonanceInfo("f_0(1500)", 1.505, 0.109, 0, 0, LauBlattWeisskopfFactor::Light );
221  resInfo_.push_back( neutral );
222  // f'_2(1525)
223  neutral = new LauResonanceInfo("f'_2(1525)", 1.525, 0.073, 2, 0, LauBlattWeisskopfFactor::Light );
224  resInfo_.push_back( neutral );
225  // f_2(1565)
226  neutral = new LauResonanceInfo("f_2(1565)", 1.562, 0.134, 2, 0, LauBlattWeisskopfFactor::Light );
227  resInfo_.push_back( neutral );
228  // f_2(1640)
229  neutral = new LauResonanceInfo("f_2(1640)", 1.639, 0.099, 2, 0, LauBlattWeisskopfFactor::Light );
230  resInfo_.push_back( neutral );
231  // f_0(1710)
232  neutral = new LauResonanceInfo("f_0(1710)", 1.722, 0.135, 0, 0, LauBlattWeisskopfFactor::Light );
233  resInfo_.push_back( neutral );
234  // f_2(1810)
235  neutral = new LauResonanceInfo("f_2(1810)", 1.816, 0.197, 2, 0, LauBlattWeisskopfFactor::Light );
236  resInfo_.push_back( neutral );
237  // f_2(1910)
238  neutral = new LauResonanceInfo("f_2(1910)", 1.903, 0.196, 2, 0, LauBlattWeisskopfFactor::Light );
239  resInfo_.push_back( neutral );
240  // f_2(1950)
241  neutral = new LauResonanceInfo("f_2(1950)", 1.944, 0.472, 2, 0, LauBlattWeisskopfFactor::Light );
242  resInfo_.push_back( neutral );
243  // f_2(2010)
244  neutral = new LauResonanceInfo("f_2(2010)", 2.011, 0.202, 2, 0, LauBlattWeisskopfFactor::Light );
245  resInfo_.push_back( neutral );
246  // f_0(2020)
247  neutral = new LauResonanceInfo("f_0(2020)", 1.992, 0.442, 0, 0, LauBlattWeisskopfFactor::Light );
248  resInfo_.push_back( neutral );
249  // f_4(2050)
250  neutral = new LauResonanceInfo("f_4(2050)", 2.018, 0.237, 4, 0, LauBlattWeisskopfFactor::Light );
251  resInfo_.push_back( neutral );
252  // f_0(2100)
253  neutral = new LauResonanceInfo("f_0(2100)", 2.101, 0.224, 0, 0, LauBlattWeisskopfFactor::Light );
254  resInfo_.push_back( neutral );
255 
256  // omega resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
257  // omega(782)
258  neutral = new LauResonanceInfo("omega(782)", 0.78265, 0.00849, 1, 0, LauBlattWeisskopfFactor::Light );
259  resInfo_.push_back( neutral );
260 
261  // a resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
262  // a_0(980)
263  neutral = new LauResonanceInfo("a0_0(980)", 0.980, 0.092, 0, 0, LauBlattWeisskopfFactor::Light );
264  positve = new LauResonanceInfo("a+_0(980)", 0.980, 0.092, 0, 1, LauBlattWeisskopfFactor::Light );
265  negatve = positve->createChargeConjugate();
266  resInfo_.push_back( neutral );
267  resInfo_.push_back( positve );
268  resInfo_.push_back( negatve );
269  // a_0(1450)
270  neutral = new LauResonanceInfo("a0_0(1450)", 1.474, 0.265, 0, 0, LauBlattWeisskopfFactor::Light );
271  positve = new LauResonanceInfo("a+_0(1450)", 1.474, 0.265, 0, 1, LauBlattWeisskopfFactor::Light );
272  negatve = positve->createChargeConjugate();
273  resInfo_.push_back( neutral );
274  resInfo_.push_back( positve );
275  resInfo_.push_back( negatve );
276  // a_2(1320)
277  neutral = new LauResonanceInfo("a0_2(1320)", 1.3190, 0.1050, 2, 0, LauBlattWeisskopfFactor::Light );
278  positve = new LauResonanceInfo("a+_2(1320)", 1.3190, 0.1050, 2, 1, LauBlattWeisskopfFactor::Light );
279  negatve = positve->createChargeConjugate();
280  resInfo_.push_back( neutral );
281  resInfo_.push_back( positve );
282  resInfo_.push_back( negatve );
283 
284  // charmonium resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
285  // chi_c0
286  neutral = new LauResonanceInfo("chi_c0", 3.41475, 0.0105, 0, 0, LauBlattWeisskopfFactor::Charmonium );
287  resInfo_.push_back( neutral );
288  // chi_c1
289  neutral = new LauResonanceInfo("chi_c1", 3.51066, 0.00084, 0, 0, LauBlattWeisskopfFactor::Charmonium );
290  resInfo_.push_back( neutral );
291  // chi_c2
292  neutral = new LauResonanceInfo("chi_c2", 3.55620, 0.00193, 2, 0, LauBlattWeisskopfFactor::Charmonium );
293  resInfo_.push_back( neutral );
294  // X(3872)
295  neutral = new LauResonanceInfo("X(3872)", 3.87169, 0.0012, 1, 0, LauBlattWeisskopfFactor::Charmonium );
296  resInfo_.push_back( neutral );
297 
298  // unknown scalars name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
299  // sigma
300  neutral = new LauResonanceInfo("sigma0", 0.475, 0.550, 0, 0, LauBlattWeisskopfFactor::Light );
301  positve = new LauResonanceInfo("sigma+", 0.475, 0.550, 0, 1, LauBlattWeisskopfFactor::Light );
302  negatve = positve->createChargeConjugate();
303  resInfo_.push_back( neutral );
304  resInfo_.push_back( positve );
305  resInfo_.push_back( negatve );
306  // kappa
307  neutral = new LauResonanceInfo("kappa0", 0.682, 0.547, 0, 0, LauBlattWeisskopfFactor::Kstar );
308  positve = new LauResonanceInfo("kappa+", 0.682, 0.547, 0, 1, LauBlattWeisskopfFactor::Kstar );
309  negatve = positve->createChargeConjugate();
310  resInfo_.push_back( neutral );
311  resInfo_.push_back( positve );
312  resInfo_.push_back( negatve );
313  // dabba
314  neutral = new LauResonanceInfo("dabba0", 2.098, 0.520, 0, 0, LauBlattWeisskopfFactor::Charm );
315  positve = new LauResonanceInfo("dabba+", 2.098, 0.520, 0, 1, LauBlattWeisskopfFactor::Charm );
316  negatve = positve->createChargeConjugate();
317  resInfo_.push_back( neutral );
318  resInfo_.push_back( positve );
319  resInfo_.push_back( negatve );
320 
321  // excited charm states name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
322  // D*
323  neutral = new LauResonanceInfo("D*0", 2.00696, 0.0021, 1, 0, LauBlattWeisskopfFactor::Charm );
324  positve = new LauResonanceInfo("D*+", 2.01026, 83.4e-6, 1, 1, LauBlattWeisskopfFactor::Charm );
325  negatve = positve->createChargeConjugate();
326  resInfo_.push_back( neutral );
327  resInfo_.push_back( positve );
328  resInfo_.push_back( negatve );
329  // D*_0
330  neutral = new LauResonanceInfo("D*0_0", 2.318, 0.267, 0, 0, LauBlattWeisskopfFactor::Charm );
331  positve = new LauResonanceInfo("D*+_0", 2.403, 0.283, 0, 1, LauBlattWeisskopfFactor::Charm );
332  negatve = positve->createChargeConjugate();
333  resInfo_.push_back( neutral );
334  resInfo_.push_back( positve );
335  resInfo_.push_back( negatve );
336  // D*_2
337  //AVERAGE--neutral = new LauResonanceInfo("D*0_2", 2.4618, 0.049, 2, 0 );
338  neutral = new LauResonanceInfo("D*0_2", 2.4626, 0.049, 2, 0, LauBlattWeisskopfFactor::Charm );
339  positve = new LauResonanceInfo("D*+_2", 2.4643, 0.037, 2, 1, LauBlattWeisskopfFactor::Charm );
340  negatve = positve->createChargeConjugate();
341  resInfo_.push_back( neutral );
342  resInfo_.push_back( positve );
343  resInfo_.push_back( negatve );
344  // D1(2420)
345  neutral = new LauResonanceInfo("D0_1(2420)", 2.4214, 0.0274, 1, 0, LauBlattWeisskopfFactor::Charm );
346  positve = new LauResonanceInfo("D+_1(2420)", 2.4232, 0.025, 1, 1, LauBlattWeisskopfFactor::Charm );
347  negatve = positve->createChargeConjugate();
348  resInfo_.push_back( neutral );
349  resInfo_.push_back( positve );
350  resInfo_.push_back( negatve );
351  // D(2600)
352  //OLD--neutral = new LauResonanceInfo("D0(2600)", 2.6087, 0.093, 0, 0 );
353  //OLD--positve = new LauResonanceInfo("D+(2600)", 2.6213, 0.093, 0, 1 );
354  neutral = new LauResonanceInfo("D0(2600)", 2.612, 0.093, 0, 0, LauBlattWeisskopfFactor::Charm );
355  positve = new LauResonanceInfo("D+(2600)", 2.612, 0.093, 0, 1, LauBlattWeisskopfFactor::Charm );
356  negatve = positve->createChargeConjugate();
357  resInfo_.push_back( neutral );
358  resInfo_.push_back( positve );
359  resInfo_.push_back( negatve );
360  // D(2760)
361  //OLD-- neutral = new LauResonanceInfo("D0(2760)", 2.7633, 0.061, 1, 0 );
362  //OLD-- positve = new LauResonanceInfo("D+(2760)", 2.7697, 0.061, 1, 1 );
363  neutral = new LauResonanceInfo("D0(2760)", 2.761, 0.063, 1, 0, LauBlattWeisskopfFactor::Charm );
364  positve = new LauResonanceInfo("D+(2760)", 2.761, 0.063, 1, 1, LauBlattWeisskopfFactor::Charm );
365  negatve = positve->createChargeConjugate();
366  resInfo_.push_back( neutral );
367  resInfo_.push_back( positve );
368  resInfo_.push_back( negatve );
369  // D(2900)
370  neutral = new LauResonanceInfo("D0(3000)", 3.0, 0.15, 0, 0, LauBlattWeisskopfFactor::Charm );
371  resInfo_.push_back( neutral );
372  // D(3400)
373  neutral = new LauResonanceInfo("D0(3400)", 3.4, 0.15, 0, 0, LauBlattWeisskopfFactor::Charm );
374  resInfo_.push_back( neutral );
375 
376  // excited strange charm name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
377  // Ds*
378  positve = new LauResonanceInfo("Ds*+", 2.1121, 0.0019, 1, 1, LauBlattWeisskopfFactor::StrangeCharm );
379  negatve = positve->createChargeConjugate();
380  resInfo_.push_back( positve );
381  resInfo_.push_back( negatve );
382  // Ds0*(2317)
383  positve = new LauResonanceInfo("Ds*+_0(2317)", 2.3177, 0.0038, 0, 1, LauBlattWeisskopfFactor::StrangeCharm );
384  negatve = positve->createChargeConjugate();
385  resInfo_.push_back( positve );
386  resInfo_.push_back( negatve );
387  // Ds2*(2573)
388  positve = new LauResonanceInfo("Ds*+_2(2573)", 2.5719, 0.017, 2, 1, LauBlattWeisskopfFactor::StrangeCharm );
389  negatve = positve->createChargeConjugate();
390  resInfo_.push_back( positve );
391  resInfo_.push_back( negatve );
392  // Ds1*(2700)
393  positve = new LauResonanceInfo("Ds*+_1(2700)", 2.709, 0.117, 1, 1, LauBlattWeisskopfFactor::StrangeCharm );
394  negatve = positve->createChargeConjugate();
395  resInfo_.push_back( positve );
396  resInfo_.push_back( negatve );
397  // Ds1*(2860)
398  positve = new LauResonanceInfo("Ds*+_1(2860)", 2.862, 0.180, 1, 1, LauBlattWeisskopfFactor::StrangeCharm );
399  negatve = positve->createChargeConjugate();
400  resInfo_.push_back( positve );
401  resInfo_.push_back( negatve );
402  // Ds3*(2860)
403  positve = new LauResonanceInfo("Ds*+_3(2860)", 2.862, 0.058, 3, 1, LauBlattWeisskopfFactor::StrangeCharm );
404  negatve = positve->createChargeConjugate();
405  resInfo_.push_back( positve );
406  resInfo_.push_back( negatve );
407 
408  // excited bottom states name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
409  // B*
410  neutral = new LauResonanceInfo("B*0", 5.3252, 0.00, 1, 0, LauBlattWeisskopfFactor::Beauty, 6.0);
411  positve = new LauResonanceInfo("B*+", 5.3252, 0.00, 1, 1, LauBlattWeisskopfFactor::Beauty, 6.0);
412  negatve = positve->createChargeConjugate();
413  resInfo_.push_back( neutral );
414  resInfo_.push_back( positve );
415  resInfo_.push_back( negatve );
416 
417  // excited strange bottom name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
418  // Bs*
419  neutral = new LauResonanceInfo("Bs*0", 5.4154, 0.00, 1, 0, LauBlattWeisskopfFactor::StrangeBeauty, 6.0);
420  resInfo_.push_back( neutral );
421 
422  // nonresonant models name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
423  // Phase-space nonresonant model
424  neutral = new LauResonanceInfo("NonReson", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
425  resInfo_.push_back( neutral );
426  // Theory-based nonresonant model
427  neutral = new LauResonanceInfo("NRModel", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
428  resInfo_.push_back( neutral );
429  // Belle nonresonant models
430  neutral = new LauResonanceInfo("BelleSymNR", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
431  resInfo_.push_back( neutral );
432  neutral = new LauResonanceInfo("BelleNR", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
433  resInfo_.push_back( neutral );
434  positve = new LauResonanceInfo("BelleNR+", 0.0, 0.0, 0, 1, LauBlattWeisskopfFactor::Light );
435  negatve = positve->createChargeConjugate();
436  resInfo_.push_back( positve );
437  resInfo_.push_back( negatve );
438  neutral = new LauResonanceInfo("BelleNR_Swave", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
439  resInfo_.push_back( neutral );
440  positve = new LauResonanceInfo("BelleNR_Swave+",0.0, 0.0, 0, 1, LauBlattWeisskopfFactor::Light );
441  negatve = positve->createChargeConjugate();
442  resInfo_.push_back( positve );
443  resInfo_.push_back( negatve );
444  neutral = new LauResonanceInfo("BelleNR_Pwave", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light );
445  resInfo_.push_back( neutral );
446  positve = new LauResonanceInfo("BelleNR_Pwave+",0.0, 0.0, 1, 1, LauBlattWeisskopfFactor::Light );
447  negatve = positve->createChargeConjugate();
448  resInfo_.push_back( positve );
449  resInfo_.push_back( negatve );
450  neutral = new LauResonanceInfo("BelleNR_Dwave", 0.0, 0.0, 2, 0, LauBlattWeisskopfFactor::Light );
451  resInfo_.push_back( neutral );
452  positve = new LauResonanceInfo("BelleNR_Dwave+",0.0, 0.0, 2, 1, LauBlattWeisskopfFactor::Light );
453  negatve = positve->createChargeConjugate();
454  resInfo_.push_back( positve );
455  resInfo_.push_back( negatve );
456  neutral = new LauResonanceInfo("BelleNR_Fwave", 0.0, 0.0, 3, 0, LauBlattWeisskopfFactor::Light );
457  resInfo_.push_back( neutral );
458  positve = new LauResonanceInfo("BelleNR_Fwave+",0.0, 0.0, 3, 1, LauBlattWeisskopfFactor::Light );
459  negatve = positve->createChargeConjugate();
460  resInfo_.push_back( positve );
461  resInfo_.push_back( negatve );
462  // Taylor expansion nonresonant model
463  neutral = new LauResonanceInfo("NRTaylor", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
464  resInfo_.push_back( neutral );
465  // Polynomial nonresonant models
466  neutral = new LauResonanceInfo("PolNR_S0", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
467  resInfo_.push_back( neutral );
468  neutral = new LauResonanceInfo("PolNR_S1", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
469  resInfo_.push_back( neutral );
470  neutral = new LauResonanceInfo("PolNR_S2", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
471  resInfo_.push_back( neutral );
472  neutral = new LauResonanceInfo("PolNR_P0", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light );
473  resInfo_.push_back( neutral );
474  neutral = new LauResonanceInfo("PolNR_P1", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light );
475  resInfo_.push_back( neutral );
476  neutral = new LauResonanceInfo("PolNR_P2", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light );
477  resInfo_.push_back( neutral );
478 
479  // Fake resonances for S-Wave splines
480  neutral = new LauResonanceInfo("Spline_S0", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
481  resInfo_.push_back( neutral );
482  neutral = new LauResonanceInfo("Spline_S0_Bar", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
483  resInfo_.push_back( neutral );
484 
485  nResDefMax_ = resInfo_.size();
486 }
487 
489 {
490  // Check whether any BW factors have been created and bail out if so
491  if ( ! bwIndepFactors_.empty() ) {
492  std::cerr << "ERROR in LauResonanceMaker::setBWType : some barrier factors have already been created - cannot change the barrier type now!" << std::endl;
493  return;
494  }
495  for ( BWFactorCategoryMap::const_iterator iter = bwFactors_.begin(); iter != bwFactors_.end(); ++iter ) {
496  if ( iter->second.bwFactor_ != 0 ) {
497  std::cerr << "ERROR in LauResonanceMaker::setBWType : some barrier factors have already been created - cannot change the barrier type now!" << std::endl;
498  return;
499  }
500  }
501 
502  bwBarrierType_ = bwType;
503 }
504 
506 {
507  // Check whether any BW factors have been created and bail out if so
508  if ( ! bwIndepFactors_.empty() ) {
509  std::cerr << "ERROR in LauResonanceMaker::setBWBachelorRestFrame : some barrier factors have already been created - cannot change the rest frame now!" << std::endl;
510  return;
511  }
512  for ( BWFactorCategoryMap::const_iterator iter = bwFactors_.begin(); iter != bwFactors_.end(); ++iter ) {
513  if ( iter->second.bwFactor_ != 0 ) {
514  std::cerr << "ERROR in LauResonanceMaker::setBWBachelorRestFrame : some barrier factors have already been created - cannot change the rest frame now!" << std::endl;
515  return;
516  }
517  }
518 
519  bwRestFrame_ = restFrame;
520 }
521 
523 {
524  if ( summaryPrinted_ ) {
525  std::cerr << "ERROR in LauResonanceMaker::setSpinFormalism : cannot redefine the spin formalism after creating one or more resonances" << std::endl;
526  return;
527  }
528  spinFormalism_ = spinType;
529 }
530 
532 {
533  if ( bwCategory == LauBlattWeisskopfFactor::Default || bwCategory == LauBlattWeisskopfFactor::Indep ) {
534  std::cerr << "WARNING in LauResonanceMaker::setDefaultBWRadius : cannot set radius values for Default or Indep categories" << std::endl;
535  return;
536  }
537 
538  // Check if we have an information object for this category
539  BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory );
540  if ( factor_iter != bwFactors_.end() ) {
541  // If so, we can set the value in the information object
542  BlattWeisskopfCategoryInfo& categoryInfo = factor_iter->second;
543  categoryInfo.defaultRadius_ = bwRadius;
544 
545  // Then we can check if a LauBlattWeisskopfFactor object has been created for this category
546  LauBlattWeisskopfFactor* bwFactor = categoryInfo.bwFactor_;
547  if ( bwFactor != 0 ) {
548  // If it has then we can also set its radius value directly
549  LauParameter* radius = bwFactor->getRadiusParameter();
550  radius->value(bwRadius);
551  radius->initValue(bwRadius);
552  radius->genValue(bwRadius);
553  }
554  } else {
555  // If not then we just store the value to be used later
556  BlattWeisskopfCategoryInfo& categoryInfo = bwFactors_[bwCategory];
557  categoryInfo.bwFactor_ = 0;
558  categoryInfo.defaultRadius_ = bwRadius;
559  categoryInfo.radiusFixed_ = kTRUE;
560  }
561 }
562 
564 {
565  if ( bwCategory == LauBlattWeisskopfFactor::Default || bwCategory == LauBlattWeisskopfFactor::Indep ) {
566  std::cerr << "WARNING in LauResonanceMaker::fixBWRadius : cannot fix/float radius values for Default or Indep categories" << std::endl;
567  return;
568  }
569 
570  // Check if we have an information object for this category
571  BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory );
572  if ( factor_iter != bwFactors_.end() ) {
573  // If so, we can set the value in the information object
574  BlattWeisskopfCategoryInfo& categoryInfo = factor_iter->second;
575  categoryInfo.radiusFixed_ = fixRadius;
576 
577  // Then we can check if a LauBlattWeisskopfFactor object has been created for this category
578  LauBlattWeisskopfFactor* bwFactor = categoryInfo.bwFactor_;
579  if ( bwFactor != 0 ) {
580  // If it has then we can also fix/float its radius value directly
581  LauParameter* radius = bwFactor->getRadiusParameter();
582  radius->fixed(fixRadius);
583  }
584  } else {
585  // If not then we just store the value to be used later
586  BlattWeisskopfCategoryInfo& categoryInfo = bwFactors_[bwCategory];
587  categoryInfo.bwFactor_ = 0;
588  categoryInfo.defaultRadius_ = -1.0;
589  categoryInfo.radiusFixed_ = fixRadius;
590  }
591 }
592 
594 {
595  LauBlattWeisskopfFactor* bwFactor(0);
596 
597  // If this is an independent factor, create it and add it to the list of independent factors, then return it
598  if ( bwCategory == LauBlattWeisskopfFactor::Indep ) {
599  bwFactor = new LauBlattWeisskopfFactor( *resInfo, bwBarrierType_, bwRestFrame_, bwCategory );
600  bwIndepFactors_.push_back(bwFactor);
601  return bwFactor;
602  }
603 
604  // Otherwise, look up the category in the category information map
605  BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory );
606 
607  if ( factor_iter == bwFactors_.end() ) {
608  // If the category is currently undefined we need to create it
609  bwFactor = new LauBlattWeisskopfFactor( *resInfo, bwBarrierType_, bwRestFrame_, bwCategory );
610 
611  BlattWeisskopfCategoryInfo& categoryInfo = bwFactors_[bwCategory];
612  categoryInfo.bwFactor_ = bwFactor;
613  categoryInfo.defaultRadius_ = bwFactor->getRadiusParameter()->value();
614  categoryInfo.radiusFixed_ = kTRUE;
615  } else {
616  // If it exists, we can check if the factor object has been created
617  BlattWeisskopfCategoryInfo& categoryInfo = factor_iter->second;
618 
619  if ( categoryInfo.bwFactor_ != 0 ) {
620  // If so, simply clone it
621  const UInt_t resSpin = resInfo->getSpin();
622  bwFactor = categoryInfo.bwFactor_->createClone( resSpin );
623  } else {
624  // Otherwise we need to create it, using the default value if it has been set
625  if ( categoryInfo.defaultRadius_ >= 0.0 ) {
626  bwFactor = new LauBlattWeisskopfFactor( *resInfo, categoryInfo.defaultRadius_, bwBarrierType_, bwRestFrame_, bwCategory );
627  } else {
628  bwFactor = new LauBlattWeisskopfFactor( *resInfo, bwBarrierType_, bwRestFrame_, bwCategory );
629  }
630  categoryInfo.bwFactor_ = bwFactor;
631 
632  // Set whether the radius should be fixed/floated
633  LauParameter* radius = bwFactor->getRadiusParameter();
634  radius->fixed( categoryInfo.radiusFixed_ );
635  }
636  }
637 
638  return bwFactor;
639 }
640 
641 LauAbsResonance* LauResonanceMaker::getResonance(const LauDaughters* daughters, const TString& resName, const Int_t resPairAmpInt, const LauAbsResonance::LauResonanceModel resType, const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory)
642 {
643  // Routine to return the appropriate LauAbsResonance object given the resonance
644  // name (resName), which daughter is the bachelor track (resPairAmpInt = 1,2 or 3),
645  // and the resonance type ("BW" = Breit-Wigner, "Flatte" = Flatte distribution).
646 
647  // If this is the first resonance we are making, first print a summary of the formalism
648  if ( ! summaryPrinted_ ) {
649  std::cout << "INFO in LauResonanceMaker::getResonance : Freezing amplitude formalism:" << std::endl;
650  switch ( spinFormalism_ ) {
652  std::cout << " : Spin factors use Zemach spin tensors, with bachelor momentum in resonance rest frame" << std::endl;
653  break;
655  std::cout << " : Spin factors use Zemach spin tensors, with bachelor momentum in parent rest frame" << std::endl;
656  break;
658  std::cout << " : Spin factors use Covariant spin tensors" << std::endl;
659  break;
661  std::cout << " : Spin factors are just Legendre polynomials" << std::endl;
662  break;
663  }
664  switch ( bwBarrierType_ ) {
666  std::cout << " : Blatt-Weisskopf barrier factors are the 'non-primed' form" << std::endl;
667  break;
669  std::cout << " : Blatt-Weisskopf barrier factors are the 'primed' form" << std::endl;
670  break;
672  std::cout << " : Blatt-Weisskopf barrier factors are the exponential form" << std::endl;
673  break;
674  }
675  switch ( bwRestFrame_ ) {
677  std::cout << " : Blatt-Weisskopf barrier factors use bachelor momentum in parent rest frame" << std::endl;
678  break;
680  std::cout << " : Blatt-Weisskopf barrier factors use bachelor momentum in resonance rest frame" << std::endl;
681  break;
683  std::cout << " : Blatt-Weisskopf barrier factors use covariant expression" << std::endl;
684  break;
685  }
686 
687  summaryPrinted_ = kTRUE;
688  }
689 
690  // Loop over all possible resonance states we have defined in
691  // createResonanceVector() until we get a match with the name of the resonance
692 
693  LauResonanceInfo* resInfo(0);
694  for (std::vector<LauResonanceInfo*>::const_iterator iter=resInfo_.begin(); iter!=resInfo_.end(); ++iter) {
695 
696  if (resName == (*iter)->getName()) {
697  // We have recognised the resonance name.
698  std::cout<<"INFO in LauResonanceMaker::getResonance : Creating resonance: "<<resName<<std::endl;
699 
700  resInfo = (*iter);
701 
702  // stop looping
703  break;
704  }
705  }
706 
707  // if we couldn't find the right resonance then we should return a null pointer
708  if ( resInfo == 0 ) {
709  std::cout<<"ERROR in LauResonanceMaker::getResonance : Unable to locate resonance info for: "<<resName<<std::endl;
710  return 0;
711  }
712 
713  // Now construct the resonance using the specified type
714  LauAbsResonance* theResonance(0);
715  switch ( resType ) {
716 
717  case LauAbsResonance::BW :
718  // Simple non-relativistic Breit-Wigner
719  std::cout<<" : Using simple Breit-Wigner lineshape. "<<std::endl;
720  theResonance = new LauBreitWignerRes(resInfo, resPairAmpInt, daughters);
721  break;
722 
724  {
725  // Relativistic Breit-Wigner with Blatt-Weisskopf factors.
726  std::cout<<" : Using relativistic Breit-Wigner lineshape. "<<std::endl;
727  theResonance = new LauRelBreitWignerRes(resInfo, resPairAmpInt, daughters);
729  LauBlattWeisskopfFactor::BlattWeisskopfCategory resCategory = bwCategory;
730  if ( bwCategory == LauBlattWeisskopfFactor::Default ) {
731  resCategory = resInfo->getBWCategory();
732  }
733  LauBlattWeisskopfFactor* resBWFactor = this->getBWFactor( resCategory, resInfo );
734  LauBlattWeisskopfFactor* parBWFactor = this->getBWFactor( parCategory, resInfo );
735  theResonance->setBarrierRadii( resBWFactor, parBWFactor );
736  break;
737  }
738 
739  case LauAbsResonance::GS :
740  {
741  // Gounaris-Sakurai function to try and model the rho(770) better
742  std::cout<<" : Using Gounaris-Sakurai lineshape. "<<std::endl;
743  theResonance = new LauGounarisSakuraiRes(resInfo, resPairAmpInt, daughters);
745  LauBlattWeisskopfFactor::BlattWeisskopfCategory resCategory = bwCategory;
746  if ( bwCategory == LauBlattWeisskopfFactor::Default ) {
747  resCategory = resInfo->getBWCategory();
748  }
749  LauBlattWeisskopfFactor* resBWFactor = this->getBWFactor( resCategory, resInfo );
750  LauBlattWeisskopfFactor* parBWFactor = this->getBWFactor( parCategory, resInfo );
751  theResonance->setBarrierRadii( resBWFactor, parBWFactor );
752  break;
753  }
754 
756  // Flatte distribution - coupled channel Breit-Wigner
757  std::cout<<" : Using Flatte lineshape. "<<std::endl;
758  theResonance = new LauFlatteRes(resInfo, resPairAmpInt, daughters);
759  break;
760 
762  // Sigma model - should only be used for the pi-pi system
763  std::cout<<" : Using Sigma lineshape. "<<std::endl;
764  theResonance = new LauSigmaRes(resInfo, resPairAmpInt, daughters);
765  break;
766 
768  // Kappa model - should only be used for the K-pi system
769  std::cout<<" : Using Kappa lineshape. "<<std::endl;
770  theResonance = new LauKappaRes(resInfo, resPairAmpInt, daughters);
771  break;
772 
774  // Dabba model - should only be used for the D-pi system
775  std::cout<<" : Using Dabba lineshape. "<<std::endl;
776  theResonance = new LauDabbaRes(resInfo, resPairAmpInt, daughters);
777  break;
778 
779  case LauAbsResonance::LASS :
780  // LASS function to try and model the K-pi S-wave better
781  std::cout<<" : Using LASS lineshape. "<<std::endl;
782  theResonance = new LauLASSRes(resInfo, resPairAmpInt, daughters);
783  break;
784 
786  // LASS function to try and model the K-pi S-wave better
787  std::cout<<" : Using LASS lineshape (resonant part only). "<<std::endl;
788  theResonance = new LauLASSBWRes(resInfo, resPairAmpInt, daughters);
789  break;
790 
792  // LASS function to try and model the K-pi S-wave better
793  std::cout<<" : Using LASS lineshape (nonresonant part only). "<<std::endl;
794  theResonance = new LauLASSNRRes(resInfo, resPairAmpInt, daughters);
795  break;
796 
798  // EFKLLM form-factor description of the K-pi S-wave
799  std::cout<<" : Using EFKLLM lineshape. "<<std::endl;
800  theResonance = new LauEFKLLMRes(resInfo, resPairAmpInt, daughters);
801  break;
802 
804  // K-matrix description of S-wave
805  std::cerr<<"ERROR in LauResonanceMaker::getResonance : K-matrix type specified, which should be separately handled."<<std::endl;
806  break;
807 
809  // uniform NR amplitude - arguments are there to preserve the interface
810  std::cout<<" : Using uniform NR lineshape. "<<std::endl;
811  // we override resPairAmpInt here
812  theResonance = new LauFlatNR(resInfo, 0, daughters);
813  break;
814 
816  // NR amplitude model - arguments are there to preserve the interface
817  std::cout<<" : Using NR-model lineshape. "<<std::endl;
818  // we override resPairAmpInt here
819  theResonance = new LauNRAmplitude(resInfo, 0, daughters);
820  break;
821 
824  // Belle NR amplitude model - arguments are there to preserve the interface
825  std::cout<<" : Using Belle NR lineshape. "<<std::endl;
826  theResonance = new LauBelleNR(resInfo, resType, resPairAmpInt, daughters);
827  break;
828 
832  // Belle NR amplitude model - arguments are there to preserve the interface
833  std::cout<<" : Using Belle symmetric NR lineshape. "<<std::endl;
834  theResonance = new LauBelleSymNR(resInfo, resType, resPairAmpInt, daughters);
835  break;
836 
838  // Polynomial NR amplitude model - arguments are there to preserve the interface
839  std::cout<<" : Using polynomial NR lineshape. "<<std::endl;
840  theResonance = new LauPolNR(resInfo, resPairAmpInt, daughters);
841  break;
842 
844  // Model independent partial wave
845  std::cout<<" : Using model independent partial wave lineshape (magnitude and phase). "<<std::endl;
846  theResonance = new LauModIndPartWaveMagPhase(resInfo, resPairAmpInt, daughters);
847  break;
848 
850  // Model independent partial wave
851  std::cout<<" : Using model independent partial wave lineshape (real and imaginary part). "<<std::endl;
852  theResonance = new LauModIndPartWaveRealImag(resInfo, resPairAmpInt, daughters);
853  break;
854 
856  // Incoherent Gaussian
857  std::cout<<" : Using incoherent Gaussian lineshape. "<<std::endl;
858  theResonance = new LauGaussIncohRes(resInfo, resPairAmpInt, daughters);
859  break;
860 
865  // Rho-omega mass mixing model
866  std::cout<<" : Using rho-omega mass mixing lineshape. "<<std::endl;
867  theResonance = new LauRhoOmegaMix(resInfo, resType, resPairAmpInt, daughters);
869  LauBlattWeisskopfFactor::BlattWeisskopfCategory resCategory = bwCategory;
870  if ( bwCategory == LauBlattWeisskopfFactor::Default ) {
871  resCategory = resInfo->getBWCategory();
872  }
873  LauBlattWeisskopfFactor* resBWFactor = this->getBWFactor( resCategory, resInfo );
874  LauBlattWeisskopfFactor* parBWFactor = this->getBWFactor( parCategory, resInfo );
875  theResonance->setBarrierRadii( resBWFactor, parBWFactor );
876  break;
877 
878  }
879 
880  // Set the spin formalism choice
881  theResonance->setSpinType( spinFormalism_ );
882 
883  return theResonance;
884 }
885 
886 Int_t LauResonanceMaker::resTypeInt(const TString& name) const
887 {
888  // Internal function that returns the resonance integer, specified by the
889  // order of the resonance vector defined in createResonanceVector(),
890  // for a given resonance name.
891  Int_t resTypInt(-99);
892  Int_t i(0);
893 
894  for (std::vector<LauResonanceInfo*>::const_iterator iter=resInfo_.begin(); iter!=resInfo_.end(); ++iter) {
895 
896  if (name.BeginsWith((*iter)->getName(), TString::kExact) == kTRUE) {
897  // We have recognised the resonance from those that are available
898  resTypInt = i;
899  return resTypInt;
900  }
901  ++i;
902  }
903 
904  return resTypInt;
905 }
906 
907 void LauResonanceMaker::printAll( std::ostream& stream ) const
908 {
909  for ( std::vector<LauResonanceInfo*>::const_iterator iter = resInfo_.begin(); iter != resInfo_.end(); ++iter ) {
910  stream << (**iter) << std::endl;
911  }
912 }
913 
914 LauResonanceInfo* LauResonanceMaker::getResInfo(const TString& resName) const
915 {
916  LauResonanceInfo* resInfo(0);
917  for (std::vector<LauResonanceInfo*>::const_iterator iter=resInfo_.begin(); iter!=resInfo_.end(); ++iter) {
918 
919  if (resName == (*iter)->getName()) {
920  // We have recognised the resonance name.
921  resInfo = (*iter);
922  // stop looping
923  break;
924  }
925  }
926 
927  return resInfo;
928 
929 }
File containing declaration of LauNRAmplitude class.
Bool_t fixed() const
Check whether the parameter is fixed or floated.
LauBlattWeisskopfFactor * getBWFactor(const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory, const LauResonanceInfo *resInfo)
Retrieve Blatt-Weisskopf factor for the given category.
LauBlattWeisskopfFactor::BlattWeisskopfCategory getBWCategory() const
Retrieve the BW category of the resonant particle.
Class for defining the LASS resonance model.
Definition: LauLASSRes.hh:31
void setBarrierRadii(LauBlattWeisskopfFactor *resFactor, LauBlattWeisskopfFactor *parFactor)
Set the form factor model and parameters.
File containing declaration of LauResonanceInfo class.
void setSpinFormalism(const LauAbsResonance::LauSpinType spinType)
Set the spin formalism to be used for all resonances.
ClassImp(LauAbsCoeffSet)
LauBlattWeisskopfFactor::RestFrame bwRestFrame_
The rest frame in which the bachelor momentum used in the Blatt-Weisskopf factors should be calculate...
BWFactorCategoryMap bwFactors_
The Blatt-Weisskopf factor objects (and related information) for each category.
Class for defining the properties of a resonant particle.
const TString & name() const
The parameter name.
Class that defines the particular 3-body decay under study.
Definition: LauDaughters.hh:33
static LauResonanceMaker * resonanceMaker_
The singleton instance.
std::vector< LauResonanceInfo * > resInfo_
The known resonances.
File containing declaration of LauBelleNR class.
void setBWType(const LauBlattWeisskopfFactor::BarrierType bwType)
Set the type of BW factor (for all categories)
Class for defininf the Gounaris-Sakurai resonance model.
void fixBWRadius(const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory, const Bool_t fixRadius)
Fix or release the Blatt-Weisskopf barrier radius for the given category.
File containing declaration of LauDaughters class.
const LauParameter * getRadiusParameter() const
Retrieve the radius parameter.
File containing declaration of LauModIndPartWaveRealImag class.
LauBlattWeisskopfFactor * bwFactor_
The BW factor object.
File containing declaration of LauGounarisSakuraiRes class.
Double_t defaultRadius_
The default value for the radius in this category.
File containing declaration of LauBelleSymNR class.
File containing declaration of LauModIndPartWaveMagPhase class.
File containing declaration of LauBreitWignerRes class.
LauSpinType
Define the allowed spin formalisms.
Int_t resTypeInt(const TString &name) const
Retrieve the integer index for the specified resonance.
File containing declaration of LauRelBreitWignerRes class.
File containing declaration of LauLASSBWRes class.
void setDefaultBWRadius(const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory, const Double_t bwRadius)
Set the BW radius for the given category.
Class for defining the EFKLLM K-pi S-wave model.
Definition: LauEFKLLMRes.hh:36
Class for defining the non resonant part of the LASS model.
Definition: LauLASSNRRes.hh:31
File containing declaration of LauEFKLLMRes class.
Class for defining the Sigma resonance model.
Definition: LauSigmaRes.hh:31
Singleton factory class for creating resonances.
void setBWBachelorRestFrame(const LauBlattWeisskopfFactor::RestFrame restFrame)
Set the rest frame in which the bachelor momentum should be calculated (for all BW categories) ...
void createResonanceVector()
Create the list of known resonances.
UInt_t getSpin() const
Retrieve the spin of the resonant particle.
File containing declaration of LauResonanceMaker class.
File containing declaration of LauRhoOmegaMix class.
File containing declaration of LauSigmaRes class.
Class for defining the Belle nonresonant model.
Definition: LauBelleNR.hh:34
Class for defining the NR amplitude model.
Class for defining the simple Breit-Wigner resonance model.
Bool_t summaryPrinted_
Boolean flag to control printing a summary of the formalism to be used when the first resonance is cr...
Class for defining a model independent partial wave component where the amplitudes are parameterised ...
LauBlattWeisskopfFactor::BarrierType bwBarrierType_
The type of the Blatt-Weisskopf barrier to use for all resonances.
Class for defining the terms of Babar nonresonant model.
Definition: LauPolNR.hh:33
Class for defining the Kappa resonance model.
Definition: LauKappaRes.hh:32
File containing declaration of LauDabbaRes class.
void printAll(std::ostream &stream) const
Print the information records, one per line, to the requested stream.
Class for defining the fit parameter objects.
Definition: LauParameter.hh:35
Class for defining the relativistic Breit-Wigner resonance model.
Class for defining the rho-omega resonance mixing model.
std::vector< LauBlattWeisskopfFactor * > bwIndepFactors_
The Blatt-Weisskopf factor objects for resonances in the independent category.
File containing declaration of LauLASSRes class.
File containing declaration of LauLASSNRRes class.
virtual ~LauResonanceMaker()
Destructor.
LauResonanceModel
Define the allowed resonance types.
Class for defining the symmetric Belle Non Resonant model.
LauAbsResonance::LauSpinType spinFormalism_
The spin formalism that should be used for all resonances.
LauResonanceMaker()
Constructor.
BarrierType
Define the allowed types of barrier factors.
LauAbsResonance * getResonance(const LauDaughters *daughters, const TString &resName, const Int_t resPairAmpInt, const LauAbsResonance::LauResonanceModel resType, const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory=LauBlattWeisskopfFactor::Default)
Create a resonance.
Class for defining an incoherent resonance with a Gaussian mass dependence.
Bool_t radiusFixed_
Whether or not the radius value for this category should be fixed in the fit.
Class for defining the resonant part of the LASS model.
Definition: LauLASSBWRes.hh:31
Abstract class for defining type for resonance amplitude models (Breit-Wigner, Flatte etc...
Class for defining the Dabba resonance model.
Definition: LauDabbaRes.hh:31
Double_t initValue() const
The initial value of the parameter.
File containing declaration of LauAbsResonance class.
RestFrame
Define the rest frame in which the momentum should be calculated (only relevant for bachelor) ...
Class that implements the Blatt-Weisskopf barrier factor.
void setSpinType(const LauSpinType spinType)
Set the spin formalism to be used.
LauBlattWeisskopfFactor * createClone(const UInt_t newSpin)
Method to create a new factor with cloned radius parameter.
Double_t value() const
The value of the parameter.
File containing declaration of LauPolNR class.
BlattWeisskopfCategory
Define resonance categories that will share common barrier factor radii.
UInt_t nResDefMax_
The number of known resonances.
Class for defining a model independent partial wave component where the amplitudes are parameterised ...
static LauResonanceMaker & get()
Get the factory instance.
Class for defining a uniform nonresonant amplitude.
Definition: LauFlatNR.hh:30
LauResonanceInfo * createChargeConjugate()
Create the charge conjugate particle info record.
LauResonanceInfo * createSharedParameterRecord(const TString &name)
Create another record that will share parameters with this one.
Class for defining the Flatte resonance model.
Definition: LauFlatteRes.hh:31
Double_t genValue() const
The value generated for the parameter.
File containing declaration of LauGaussIncohRes class.
LauResonanceInfo * getResInfo(const TString &resName) const
Get the information for the given resonance name.
Data structure to store information on a given Blatt-Weisskopf category.
File containing declaration of LauKappaRes class.
File containing declaration of LauFlatteRes class.
File containing declaration of LauFlatNR class.