laura is hosted by Hepforge, IPPP Durham
Laura++  v3r0p1
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 "LauFlatteRes.hh"
24 #include "LauFlatNR.hh"
25 #include "LauModIndPartWave.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"
32 #include "LauNRAmplitude.hh"
33 #include "LauPolNR.hh"
34 #include "LauRelBreitWignerRes.hh"
35 #include "LauResonanceInfo.hh"
36 #include "LauResonanceMaker.hh"
37 #include "LauSigmaRes.hh"
38 
40 
41 
43 
44 
46  nResDefMax_(0)
47 {
48  this->createResonanceVector();
50 }
51 
53 {
54  for ( std::vector<LauBlattWeisskopfFactor*>::iterator iter = bwIndepFactors_.begin(); iter != bwIndepFactors_.end(); ++iter ) {
55  delete *iter;
56  }
57  bwIndepFactors_.clear();
58  for ( BWFactorCategoryMap::iterator iter = bwFactors_.begin(); iter != bwFactors_.end(); ++iter ) {
59  delete iter->second;
60  }
61  bwFactors_.clear();
62 }
63 
65 {
66  if ( resonanceMaker_ == 0 ) {
68  }
69 
70  return *resonanceMaker_;
71 }
72 
74 {
75  // Function to create all possible resonances that this class supports.
76  // Also add in the sigma and kappa - but a special paramterisation is used
77  // instead of the PDG "pole mass and width" values.
78 
79  std::cout << "INFO in LauResonanceMaker::createResonanceVector : Setting up possible resonance states..." << std::endl;
80 
81  LauResonanceInfo* neutral(0);
82  LauResonanceInfo* positve(0);
83  LauResonanceInfo* negatve(0);
84 
85  // Define the resonance names and store them in the array
86  resInfo_.clear();
87  resInfo_.reserve(100);
88  // rho resonances name, mass, width, spin, charge, default BW category, BW radius parameter (defaults to 4.0)
89  // rho(770)
90  neutral = new LauResonanceInfo("rho0(770)", 0.77526, 0.1478, 1, 0, LauBlattWeisskopfFactor::Light, 5.3);
91  positve = new LauResonanceInfo("rho+(770)", 0.77511, 0.1491, 1, 1, LauBlattWeisskopfFactor::Light, 5.3);
92  negatve = positve->createChargeConjugate();
93  resInfo_.push_back( neutral );
94  resInfo_.push_back( positve );
95  resInfo_.push_back( negatve );
96  // rho(1450)
97  neutral = new LauResonanceInfo("rho0(1450)", 1.465, 0.400, 1, 0, LauBlattWeisskopfFactor::Light );
98  positve = new LauResonanceInfo("rho+(1450)", 1.465, 0.400, 1, 1, LauBlattWeisskopfFactor::Light );
99  negatve = positve->createChargeConjugate();
100  resInfo_.push_back( neutral );
101  resInfo_.push_back( positve );
102  resInfo_.push_back( negatve );
103  // rho(1700)
104  neutral = new LauResonanceInfo("rho0(1700)", 1.720, 0.250, 1, 0, LauBlattWeisskopfFactor::Light );
105  positve = new LauResonanceInfo("rho+(1700)", 1.720, 0.250, 1, 1, LauBlattWeisskopfFactor::Light );
106  negatve = positve->createChargeConjugate();
107  resInfo_.push_back( neutral );
108  resInfo_.push_back( positve );
109  resInfo_.push_back( negatve );
110 
111  // K* resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
112  // K*(892)
113  neutral = new LauResonanceInfo("K*0(892)", 0.89581, 0.0474, 1, 0, LauBlattWeisskopfFactor::Kstar, 3.0);
114  positve = new LauResonanceInfo("K*+(892)", 0.89166, 0.0508, 1, 1, LauBlattWeisskopfFactor::Kstar, 3.0);
115  negatve = positve->createChargeConjugate();
116  resInfo_.push_back( neutral );
117  resInfo_.push_back( positve );
118  resInfo_.push_back( negatve );
119  // K*(1410)
120  neutral = new LauResonanceInfo("K*0(1410)", 1.414, 0.232, 1, 0, LauBlattWeisskopfFactor::Kstar );
121  positve = new LauResonanceInfo("K*+(1410)", 1.414, 0.232, 1, 1, LauBlattWeisskopfFactor::Kstar );
122  negatve = positve->createChargeConjugate();
123  resInfo_.push_back( neutral );
124  resInfo_.push_back( positve );
125  resInfo_.push_back( negatve );
126  // K*_0(1430)
127  neutral = new LauResonanceInfo("K*0_0(1430)", 1.425, 0.270, 0, 0, LauBlattWeisskopfFactor::Kstar );
128  positve = new LauResonanceInfo("K*+_0(1430)", 1.425, 0.270, 0, 1, LauBlattWeisskopfFactor::Kstar );
129  negatve = positve->createChargeConjugate();
130  resInfo_.push_back( neutral );
131  resInfo_.push_back( positve );
132  resInfo_.push_back( negatve );
133  // LASS nonresonant model
134  neutral = neutral->createSharedParameterRecord("LASSNR0");
135  positve = positve->createSharedParameterRecord("LASSNR+");
136  negatve = negatve->createSharedParameterRecord("LASSNR-");
137  resInfo_.push_back( neutral );
138  // K*_2(1430)
139  neutral = new LauResonanceInfo("K*0_2(1430)", 1.4324, 0.109, 2, 0, LauBlattWeisskopfFactor::Kstar );
140  positve = new LauResonanceInfo("K*+_2(1430)", 1.4256, 0.0985, 2, 1, LauBlattWeisskopfFactor::Kstar );
141  negatve = positve->createChargeConjugate();
142  resInfo_.push_back( neutral );
143  resInfo_.push_back( positve );
144  resInfo_.push_back( negatve );
145  // K*(1680)
146  neutral = new LauResonanceInfo("K*0(1680)", 1.717, 0.322, 1, 0, LauBlattWeisskopfFactor::Kstar );
147  positve = new LauResonanceInfo("K*+(1680)", 1.717, 0.322, 1, 1, LauBlattWeisskopfFactor::Kstar );
148  negatve = positve->createChargeConjugate();
149  resInfo_.push_back( neutral );
150  resInfo_.push_back( positve );
151  resInfo_.push_back( negatve );
152 
153  // phi resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
154  // phi(1020)
155  neutral = new LauResonanceInfo("phi(1020)", 1.019461, 0.004266, 1, 0, LauBlattWeisskopfFactor::Light );
156  resInfo_.push_back( neutral );
157  // phi(1680)
158  neutral = new LauResonanceInfo("phi(1680)", 1.680, 0.150, 1, 0, LauBlattWeisskopfFactor::Light );
159  resInfo_.push_back( neutral );
160 
161  // f resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
162  // f_0(980)
163  neutral = new LauResonanceInfo("f_0(980)", 0.990, 0.070, 0, 0, LauBlattWeisskopfFactor::Light );
164  resInfo_.push_back( neutral );
165  // f_2(1270)
166  neutral = new LauResonanceInfo("f_2(1270)", 1.2751, 0.1851, 2, 0, LauBlattWeisskopfFactor::Light );
167  resInfo_.push_back( neutral );
168  // f_0(1370)
169  neutral = new LauResonanceInfo("f_0(1370)", 1.370, 0.350, 0, 0, LauBlattWeisskopfFactor::Light );
170  resInfo_.push_back( neutral );
171  // f'_0(1300), from Belle's Kspipi paper
172  neutral = new LauResonanceInfo("f'_0(1300)", 1.449, 0.126, 0, 0, LauBlattWeisskopfFactor::Light );
173  resInfo_.push_back( neutral );
174  // f_0(1500)
175  neutral = new LauResonanceInfo("f_0(1500)", 1.505, 0.109, 0, 0, LauBlattWeisskopfFactor::Light );
176  resInfo_.push_back( neutral );
177  // f'_2(1525)
178  neutral = new LauResonanceInfo("f'_2(1525)", 1.525, 0.073, 2, 0, LauBlattWeisskopfFactor::Light );
179  resInfo_.push_back( neutral );
180  // f_0(1710)
181  neutral = new LauResonanceInfo("f_0(1710)", 1.722, 0.135, 0, 0, LauBlattWeisskopfFactor::Light );
182  resInfo_.push_back( neutral );
183  // f_2(2010)
184  neutral = new LauResonanceInfo("f_2(2010)", 2.011, 0.202, 2, 0, LauBlattWeisskopfFactor::Light );
185  resInfo_.push_back( neutral );
186 
187  // omega resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
188  // omega(782)
189  neutral = new LauResonanceInfo("omega(782)", 0.78265, 0.00849, 1, 0, LauBlattWeisskopfFactor::Light );
190  resInfo_.push_back( neutral );
191 
192  // a resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
193  // a_0(980)
194  neutral = new LauResonanceInfo("a0_0(980)", 0.980, 0.092, 0, 0, LauBlattWeisskopfFactor::Light );
195  positve = new LauResonanceInfo("a+_0(980)", 0.980, 0.092, 0, 1, LauBlattWeisskopfFactor::Light );
196  negatve = positve->createChargeConjugate();
197  resInfo_.push_back( neutral );
198  resInfo_.push_back( positve );
199  resInfo_.push_back( negatve );
200  // a_0(1450)
201  neutral = new LauResonanceInfo("a0_0(1450)", 1.474, 0.265, 0, 0, LauBlattWeisskopfFactor::Light );
202  positve = new LauResonanceInfo("a+_0(1450)", 1.474, 0.265, 0, 1, LauBlattWeisskopfFactor::Light );
203  negatve = positve->createChargeConjugate();
204  resInfo_.push_back( neutral );
205  resInfo_.push_back( positve );
206  resInfo_.push_back( negatve );
207  // a_2(1320)
208  neutral = new LauResonanceInfo("a0_2(1320)", 1.3190, 0.1050, 2, 0, LauBlattWeisskopfFactor::Light );
209  positve = new LauResonanceInfo("a+_2(1320)", 1.3190, 0.1050, 2, 1, LauBlattWeisskopfFactor::Light );
210  negatve = positve->createChargeConjugate();
211  resInfo_.push_back( neutral );
212  resInfo_.push_back( positve );
213  resInfo_.push_back( negatve );
214 
215  // charmonium resonances name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
216  // chi_c0
217  neutral = new LauResonanceInfo("chi_c0", 3.41475, 0.0105, 0, 0, LauBlattWeisskopfFactor::Charmonium );
218  resInfo_.push_back( neutral );
219  // chi_c1
220  neutral = new LauResonanceInfo("chi_c1", 3.51066, 0.00084, 0, 0, LauBlattWeisskopfFactor::Charmonium );
221  resInfo_.push_back( neutral );
222  // chi_c2
223  neutral = new LauResonanceInfo("chi_c2", 3.55620, 0.00193, 2, 0, LauBlattWeisskopfFactor::Charmonium );
224  resInfo_.push_back( neutral );
225  // X(3872)
226  neutral = new LauResonanceInfo("X(3872)", 3.87169, 0.0012, 1, 0, LauBlattWeisskopfFactor::Charmonium );
227  resInfo_.push_back( neutral );
228 
229  // unknown scalars name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
230  // sigma
231  neutral = new LauResonanceInfo("sigma0", 0.475, 0.550, 0, 0, LauBlattWeisskopfFactor::Light );
232  positve = new LauResonanceInfo("sigma+", 0.475, 0.550, 0, 1, LauBlattWeisskopfFactor::Light );
233  negatve = positve->createChargeConjugate();
234  resInfo_.push_back( neutral );
235  resInfo_.push_back( positve );
236  resInfo_.push_back( negatve );
237  // kappa
238  neutral = new LauResonanceInfo("kappa0", 0.682, 0.547, 0, 0, LauBlattWeisskopfFactor::Kstar );
239  positve = new LauResonanceInfo("kappa+", 0.682, 0.547, 0, 1, LauBlattWeisskopfFactor::Kstar );
240  negatve = positve->createChargeConjugate();
241  resInfo_.push_back( neutral );
242  resInfo_.push_back( positve );
243  resInfo_.push_back( negatve );
244  // dabba
245  neutral = new LauResonanceInfo("dabba0", 2.098, 0.520, 0, 0, LauBlattWeisskopfFactor::Charm );
246  positve = new LauResonanceInfo("dabba+", 2.098, 0.520, 0, 1, LauBlattWeisskopfFactor::Charm );
247  negatve = positve->createChargeConjugate();
248  resInfo_.push_back( neutral );
249  resInfo_.push_back( positve );
250  resInfo_.push_back( negatve );
251 
252  // excited charm states name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
253  // D*
254  neutral = new LauResonanceInfo("D*0", 2.00696, 0.0021, 1, 0, LauBlattWeisskopfFactor::Charm );
255  positve = new LauResonanceInfo("D*+", 2.01026, 83.4e-6, 1, 1, LauBlattWeisskopfFactor::Charm );
256  negatve = positve->createChargeConjugate();
257  resInfo_.push_back( neutral );
258  resInfo_.push_back( positve );
259  resInfo_.push_back( negatve );
260  // D*_0
261  neutral = new LauResonanceInfo("D*0_0", 2.318, 0.267, 0, 0, LauBlattWeisskopfFactor::Charm );
262  positve = new LauResonanceInfo("D*+_0", 2.403, 0.283, 0, 1, LauBlattWeisskopfFactor::Charm );
263  negatve = positve->createChargeConjugate();
264  resInfo_.push_back( neutral );
265  resInfo_.push_back( positve );
266  resInfo_.push_back( negatve );
267  // D*_2
268  //AVERAGE--neutral = new LauResonanceInfo("D*0_2", 2.4618, 0.049, 2, 0 );
269  neutral = new LauResonanceInfo("D*0_2", 2.4626, 0.049, 2, 0, LauBlattWeisskopfFactor::Charm );
270  positve = new LauResonanceInfo("D*+_2", 2.4643, 0.037, 2, 1, LauBlattWeisskopfFactor::Charm );
271  negatve = positve->createChargeConjugate();
272  resInfo_.push_back( neutral );
273  resInfo_.push_back( positve );
274  resInfo_.push_back( negatve );
275  // D1(2420)
276  neutral = new LauResonanceInfo("D0_1(2420)", 2.4214, 0.0274, 1, 0, LauBlattWeisskopfFactor::Charm );
277  positve = new LauResonanceInfo("D+_1(2420)", 2.4232, 0.025, 1, 1, LauBlattWeisskopfFactor::Charm );
278  negatve = positve->createChargeConjugate();
279  resInfo_.push_back( neutral );
280  resInfo_.push_back( positve );
281  resInfo_.push_back( negatve );
282  // D(2600)
283  //OLD--neutral = new LauResonanceInfo("D0(2600)", 2.6087, 0.093, 0, 0 );
284  //OLD--positve = new LauResonanceInfo("D+(2600)", 2.6213, 0.093, 0, 1 );
285  neutral = new LauResonanceInfo("D0(2600)", 2.612, 0.093, 0, 0, LauBlattWeisskopfFactor::Charm );
286  positve = new LauResonanceInfo("D+(2600)", 2.612, 0.093, 0, 1, LauBlattWeisskopfFactor::Charm );
287  negatve = positve->createChargeConjugate();
288  resInfo_.push_back( neutral );
289  resInfo_.push_back( positve );
290  resInfo_.push_back( negatve );
291  // D(2760)
292  //OLD-- neutral = new LauResonanceInfo("D0(2760)", 2.7633, 0.061, 1, 0 );
293  //OLD-- positve = new LauResonanceInfo("D+(2760)", 2.7697, 0.061, 1, 1 );
294  neutral = new LauResonanceInfo("D0(2760)", 2.761, 0.063, 1, 0, LauBlattWeisskopfFactor::Charm );
295  positve = new LauResonanceInfo("D+(2760)", 2.761, 0.063, 1, 1, LauBlattWeisskopfFactor::Charm );
296  negatve = positve->createChargeConjugate();
297  resInfo_.push_back( neutral );
298  resInfo_.push_back( positve );
299  resInfo_.push_back( negatve );
300  // D(2900)
301  neutral = new LauResonanceInfo("D0(3000)", 3.0, 0.15, 0, 0, LauBlattWeisskopfFactor::Charm );
302  resInfo_.push_back( neutral );
303  // D(3400)
304  neutral = new LauResonanceInfo("D0(3400)", 3.4, 0.15, 0, 0, LauBlattWeisskopfFactor::Charm );
305  resInfo_.push_back( neutral );
306 
307  // excited strange charm name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
308  // Ds*
309  positve = new LauResonanceInfo("Ds*+", 2.1121, 0.0019, 1, 1, LauBlattWeisskopfFactor::StrangeCharm );
310  negatve = positve->createChargeConjugate();
311  resInfo_.push_back( positve );
312  resInfo_.push_back( negatve );
313  // Ds0*(2317)
314  positve = new LauResonanceInfo("Ds*+_0(2317)", 2.3177, 0.0038, 0, 1, LauBlattWeisskopfFactor::StrangeCharm );
315  negatve = positve->createChargeConjugate();
316  resInfo_.push_back( positve );
317  resInfo_.push_back( negatve );
318  // Ds2*(2573)
319  positve = new LauResonanceInfo("Ds*+_2(2573)", 2.5719, 0.017, 2, 1, LauBlattWeisskopfFactor::StrangeCharm );
320  negatve = positve->createChargeConjugate();
321  resInfo_.push_back( positve );
322  resInfo_.push_back( negatve );
323  // Ds1*(2700)
324  positve = new LauResonanceInfo("Ds*+_1(2700)", 2.709, 0.117, 1, 1, LauBlattWeisskopfFactor::StrangeCharm );
325  negatve = positve->createChargeConjugate();
326  resInfo_.push_back( positve );
327  resInfo_.push_back( negatve );
328 
329  // excited bottom states name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
330  // B*
331  neutral = new LauResonanceInfo("B*0", 5.3252, 0.00, 1, 0, LauBlattWeisskopfFactor::Beauty, 6.0);
332  positve = new LauResonanceInfo("B*+", 5.3252, 0.00, 1, 1, LauBlattWeisskopfFactor::Beauty, 6.0);
333  negatve = positve->createChargeConjugate();
334  resInfo_.push_back( neutral );
335  resInfo_.push_back( positve );
336  resInfo_.push_back( negatve );
337 
338  // excited strange bottom name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
339  // Bs*
340  neutral = new LauResonanceInfo("Bs*0", 5.4154, 0.00, 1, 0, LauBlattWeisskopfFactor::StrangeBeauty, 6.0);
341  resInfo_.push_back( neutral );
342 
343  // nonresonant models name, mass, width, spin, charge, BW category, BW radius parameter (defaults to 4.0)
344  // Phase-space nonresonant model
345  neutral = new LauResonanceInfo("NonReson", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
346  resInfo_.push_back( neutral );
347  // Theory-based nonresonant model
348  neutral = new LauResonanceInfo("NRModel", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
349  resInfo_.push_back( neutral );
350  // Belle nonresonant models
351  neutral = new LauResonanceInfo("BelleSymNR", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
352  resInfo_.push_back( neutral );
353  neutral = new LauResonanceInfo("BelleNR", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
354  resInfo_.push_back( neutral );
355  positve = new LauResonanceInfo("BelleNR+", 0.0, 0.0, 0, 1, LauBlattWeisskopfFactor::Light );
356  negatve = positve->createChargeConjugate();
357  resInfo_.push_back( positve );
358  resInfo_.push_back( negatve );
359  neutral = new LauResonanceInfo("BelleNR_Swave", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
360  resInfo_.push_back( neutral );
361  positve = new LauResonanceInfo("BelleNR_Swave+",0.0, 0.0, 0, 1, LauBlattWeisskopfFactor::Light );
362  negatve = positve->createChargeConjugate();
363  resInfo_.push_back( positve );
364  resInfo_.push_back( negatve );
365  neutral = new LauResonanceInfo("BelleNR_Pwave", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light );
366  resInfo_.push_back( neutral );
367  positve = new LauResonanceInfo("BelleNR_Pwave+",0.0, 0.0, 1, 1, LauBlattWeisskopfFactor::Light );
368  negatve = positve->createChargeConjugate();
369  resInfo_.push_back( positve );
370  resInfo_.push_back( negatve );
371  neutral = new LauResonanceInfo("BelleNR_Dwave", 0.0, 0.0, 2, 0, LauBlattWeisskopfFactor::Light );
372  resInfo_.push_back( neutral );
373  positve = new LauResonanceInfo("BelleNR_Dwave+",0.0, 0.0, 2, 1, LauBlattWeisskopfFactor::Light );
374  negatve = positve->createChargeConjugate();
375  resInfo_.push_back( positve );
376  resInfo_.push_back( negatve );
377  // Taylor expansion nonresonant model
378  neutral = new LauResonanceInfo("NRTaylor", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
379  resInfo_.push_back( neutral );
380  // Polynomial nonresonant models
381  neutral = new LauResonanceInfo("PolNR_S0", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
382  resInfo_.push_back( neutral );
383  neutral = new LauResonanceInfo("PolNR_S1", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
384  resInfo_.push_back( neutral );
385  neutral = new LauResonanceInfo("PolNR_S2", 0.0, 0.0, 0, 0, LauBlattWeisskopfFactor::Light );
386  resInfo_.push_back( neutral );
387  neutral = new LauResonanceInfo("PolNR_P0", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light );
388  resInfo_.push_back( neutral );
389  neutral = new LauResonanceInfo("PolNR_P1", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light );
390  resInfo_.push_back( neutral );
391  neutral = new LauResonanceInfo("PolNR_P2", 0.0, 0.0, 1, 0, LauBlattWeisskopfFactor::Light );
392  resInfo_.push_back( neutral );
393 
394  nResDefMax_ = resInfo_.size();
395 }
396 
398 {
399  if ( bwCategory == LauBlattWeisskopfFactor::Default || bwCategory == LauBlattWeisskopfFactor::Indep ) {
400  std::cerr << "WARNING in LauResonanceMaker::setDefaultBWRadius : cannot set radius values for Default or Indep categories" << std::endl;
401  return;
402  }
403 
404  // Check if a LauBlattWeisskopfFactor object has been created for this category
405  // If it has then we can set it directly
406  // If not then we just store the value to be used when it is created later
407  BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory );
408  if ( factor_iter != bwFactors_.end() ) {
409  LauBlattWeisskopfFactor* bwFactor = factor_iter->second;
410  LauParameter* radius = bwFactor->getRadiusParameter();
411  radius->value(bwRadius);
412  radius->initValue(bwRadius);
413  radius->genValue(bwRadius);
414  }
415 
416  bwDefaultRadii_[bwCategory] = bwRadius;
417 }
418 
420 {
421  if ( bwCategory == LauBlattWeisskopfFactor::Default || bwCategory == LauBlattWeisskopfFactor::Indep ) {
422  std::cerr << "WARNING in LauResonanceMaker::fixBWRadius : cannot fix/float radius values for Default or Indep categories" << std::endl;
423  return;
424  }
425 
426  // Check if a LauBlattWeisskopfFactor object has been created for this category
427  // If it has then we can set it directly
428  // If not then we just store the value to be used when it is created later
429  BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory );
430  if ( factor_iter != bwFactors_.end() ) {
431  LauBlattWeisskopfFactor* bwFactor = factor_iter->second;
432  LauParameter* radius = bwFactor->getRadiusParameter();
433  radius->fixed(fixRadius);
434  }
435 
436  bwFixRadii_[bwCategory] = fixRadius;
437 }
438 
440 {
441  LauBlattWeisskopfFactor* bwFactor(0);
442 
443  BWFactorCategoryMap::iterator factor_iter = bwFactors_.find( bwCategory );
444  BWRadiusCategoryMap::const_iterator radius_iter = bwDefaultRadii_.find( bwCategory );
445 
446  if ( bwCategory == LauBlattWeisskopfFactor::Indep ) {
447  bwFactor = new LauBlattWeisskopfFactor( *resInfo, bwType, bwCategory );
448  bwIndepFactors_.push_back(bwFactor);
449  } else if ( factor_iter == bwFactors_.end() ) {
450  if ( radius_iter == bwDefaultRadii_.end() ) {
451  bwFactor = new LauBlattWeisskopfFactor( *resInfo, bwType, bwCategory );
452  } else {
453  bwFactor = new LauBlattWeisskopfFactor( *resInfo, radius_iter->second, bwType, bwCategory );
454  }
455  bwFactors_[bwCategory] = bwFactor;
456  } else {
457  const UInt_t resSpin = resInfo->getSpin();
458  bwFactor = factor_iter->second->createClone( resSpin );
459  if ( bwFactor->getBarrierType() != bwType ) {
460  std::cerr << "WARNING in LauResonanceMaker::getBWFactor : A barrier factor already exists for the specified category but does not have the specified barrier type.\n";
461  std::cerr << " : If you want to use that type you will need to use a different category for this resonance." << std::endl;
462  }
463  }
464 
465  BWRadiusFixedCategoryMap::const_iterator radius_fixed_iter = bwFixRadii_.find( bwCategory );
466 
467  if ( radius_fixed_iter != bwFixRadii_.end() ) {
468  const Bool_t fixed = radius_fixed_iter->second;
469 
470  LauParameter* radius = bwFactor->getRadiusParameter();
471  radius->fixed( fixed );
472  }
473 
474  return bwFactor;
475 }
476 
477 LauAbsResonance* LauResonanceMaker::getResonance(const LauDaughters* daughters, const TString& resName, const Int_t resPairAmpInt, const LauAbsResonance::LauResonanceModel resType, const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory, const LauBlattWeisskopfFactor::BarrierType bwType)
478 {
479  // Routine to return the appropriate LauAbsResonance object given the resonance
480  // name (resName), which daughter is the bachelor track (resPairAmpInt = 1,2 or 3),
481  // and the resonance type ("BW" = Breit-Wigner, "Flatte" = Flatte distribution).
482 
483  // Loop over all possible resonance states we have defined in
484  // createResonanceVector() until we get a match with the name of the resonance
485 
486  LauResonanceInfo* resInfo(0);
487  for (std::vector<LauResonanceInfo*>::const_iterator iter=resInfo_.begin(); iter!=resInfo_.end(); ++iter) {
488 
489  if (resName == (*iter)->getName()) {
490  // We have recognised the resonance name.
491  std::cout<<"INFO in LauResonanceMaker::getResonance : Creating resonance: "<<resName<<std::endl;
492 
493  resInfo = (*iter);
494 
495  // stop looping
496  break;
497  }
498  }
499 
500  // if we couldn't find the right resonance then we should return a null pointer
501  if ( resInfo == 0 ) {
502  std::cout<<"WARNING in LauResonanceMaker::getResonance : Creating resonance: "<<resName<<std::endl;
503  return 0;
504  }
505 
506  LauAbsResonance* theResonance(0);
507 
508  // Now construct the resonnace using the right type.
509  // If we don't recognise the resonance model name, just use a simple Breit-Wigner.
510 
511  if ( resType == LauAbsResonance::Flatte ) {
512 
513  // Flatte distribution - coupled channel Breit-Wigner
514  std::cout<<" : Using Flatte lineshape. "<<std::endl;
515  theResonance =
516  new LauFlatteRes(resInfo, resPairAmpInt, daughters);
517 
518  } else if ( resType == LauAbsResonance::RelBW ) {
519 
520  // Relativistic Breit-Wigner with Blatt-Weisskopf factors.
521  std::cout<<" : Using relativistic Breit-Wigner lineshape. "<<std::endl;
522  theResonance =
523  new LauRelBreitWignerRes(resInfo, resPairAmpInt, daughters);
524 
526  LauBlattWeisskopfFactor::BlattWeisskopfCategory resCategory = bwCategory;
527  if ( bwCategory == LauBlattWeisskopfFactor::Default ) {
528  resCategory = resInfo->getBWCategory();
529  }
530  LauBlattWeisskopfFactor* resBWFactor = this->getBWFactor( resCategory, resInfo, bwType );
531  LauBlattWeisskopfFactor* parBWFactor = this->getBWFactor( parCategory, resInfo, bwType );
532  theResonance->setBarrierRadii( resBWFactor, parBWFactor );
533 
534  } else if ( resType == LauAbsResonance::Dabba ) {
535 
536  // Dabba model - should only be used for the D-pi system
537  std::cout<<" : Using Dabba lineshape. "<<std::endl;
538  theResonance =
539  new LauDabbaRes(resInfo, resPairAmpInt, daughters);
540 
541  } else if ( resType == LauAbsResonance::Kappa ) {
542 
543  // Kappa model - should only be used for the K-pi system
544  std::cout<<" : Using Kappa lineshape. "<<std::endl;
545  theResonance =
546  new LauKappaRes(resInfo, resPairAmpInt, daughters);
547 
548  } else if ( resType == LauAbsResonance::Sigma ) {
549 
550  // Sigma model - should only be used for the pi-pi system
551  std::cout<<" : Using Sigma lineshape. "<<std::endl;
552  theResonance =
553  new LauSigmaRes(resInfo, resPairAmpInt, daughters);
554 
555  } else if ( resType == LauAbsResonance::LASS_BW ) {
556 
557  // LASS function to try and model the K-pi S-wave better
558  std::cout<<" : Using LASS lineshape (resonant part only). "<<std::endl;
559  theResonance =
560  new LauLASSBWRes(resInfo, resPairAmpInt, daughters);
561 
562  } else if ( resType == LauAbsResonance::LASS_NR ) {
563 
564  // LASS function to try and model the K-pi S-wave better
565  std::cout<<" : Using LASS lineshape (nonresonant part only). "<<std::endl;
566  theResonance =
567  new LauLASSNRRes(resInfo, resPairAmpInt, daughters);
568 
569  } else if ( resType == LauAbsResonance::LASS ) {
570 
571  // LASS function to try and model the K-pi S-wave better
572  std::cout<<" : Using LASS lineshape. "<<std::endl;
573  theResonance =
574  new LauLASSRes(resInfo, resPairAmpInt, daughters);
575 
576  } else if ( resType == LauAbsResonance::GS ) {
577 
578  // Gounaris-Sakurai function to try and model the rho(770) better
579  std::cout<<" : Using Gounaris-Sakurai lineshape. "<<std::endl;
580  theResonance =
581  new LauGounarisSakuraiRes(resInfo, resPairAmpInt, daughters);
582 
584  LauBlattWeisskopfFactor::BlattWeisskopfCategory resCategory = bwCategory;
585  if ( bwCategory == LauBlattWeisskopfFactor::Default ) {
586  resCategory = resInfo->getBWCategory();
587  }
588  LauBlattWeisskopfFactor* resBWFactor = this->getBWFactor( resCategory, resInfo, bwType );
589  LauBlattWeisskopfFactor* parBWFactor = this->getBWFactor( parCategory, resInfo, bwType );
590  theResonance->setBarrierRadii( resBWFactor, parBWFactor );
591 
592  } else if ( resType == LauAbsResonance::FlatNR ) {
593 
594  // uniform NR amplitude - arguments are there to preserve the interface
595  std::cout<<" : Using uniform NR lineshape. "<<std::endl;
596  // we override resPairAmpInt here
597  theResonance =
598  new LauFlatNR(resInfo, 0, daughters);
599 
600  } else if ( resType == LauAbsResonance::NRModel ) {
601 
602  // NR amplitude model - arguments are there to preserve the interface
603  std::cout<<" : Using NR-model lineshape. "<<std::endl;
604  // we override resPairAmpInt here
605  theResonance =
606  new LauNRAmplitude(resInfo, 0, daughters);
607 
608  } else if ( resType == LauAbsResonance::BelleSymNR || resType == LauAbsResonance::TaylorNR ) {
609 
610  // Belle NR amplitude model - arguments are there to preserve the interface
611  std::cout<<" : Using Belle symmetric NR lineshape. "<<std::endl;
612  theResonance =
613  new LauBelleSymNR(resInfo, resType, resPairAmpInt, daughters);
614 
615  } else if ( resType == LauAbsResonance::BelleNR || resType == LauAbsResonance::PowerLawNR ) {
616 
617  // Belle NR amplitude model - arguments are there to preserve the interface
618  std::cout<<" : Using Belle NR lineshape. "<<std::endl;
619  theResonance =
620  new LauBelleNR(resInfo, resType, resPairAmpInt, daughters);
621 
622  } else if ( resType == LauAbsResonance::PolNR ) {
623 
624  // Polynomial NR amplitude model - arguments are there to preserve the interface
625  std::cout<<" : Using polynomial NR lineshape. "<<std::endl;
626  theResonance =
627  new LauPolNR(resInfo, resPairAmpInt, daughters);
628 
629  } else if ( resType == LauAbsResonance::MIPW ) {
630 
631  // Model independent partial wave
632  std::cout<<" : Using model independent partial wave lineshape. "<<std::endl;
633  theResonance =
634  new LauModIndPartWave(resInfo, resPairAmpInt, daughters);
635 
636  } else if ( resType == LauAbsResonance::GaussIncoh ) {
637 
638  // Incoherent Gaussian
639  std::cout<<" : Using incoherent Gaussian lineshape. "<<std::endl;
640  theResonance =
641  new LauGaussIncohRes(resInfo, resPairAmpInt, daughters);
642 
643  } else if ( resType == LauAbsResonance::BW ) {
644 
645  // Simple non-relativistic Breit-Wigner
646  std::cout<<" : Using simple Breit-Wigner lineshape. "<<std::endl;
647  theResonance =
648  new LauBreitWignerRes(resInfo, resPairAmpInt, daughters);
649 
650  } else {
651  std::cerr << "ERROR in LauResonanceMaker::getResonance : Could not match resonace type \"" << resType << "\"." << std::endl;
652  return 0;
653  }
654 
655  return theResonance;
656 
657 }
658 
659 Int_t LauResonanceMaker::resTypeInt(const TString& name) const
660 {
661  // Internal function that returns the resonance integer, specified by the
662  // order of the resonance vector defined in createResonanceVector(),
663  // for a given resonance name.
664  Int_t resTypInt(-99);
665  Int_t i(0);
666 
667  for (std::vector<LauResonanceInfo*>::const_iterator iter=resInfo_.begin(); iter!=resInfo_.end(); ++iter) {
668 
669  if (name.BeginsWith((*iter)->getName(), TString::kExact) == kTRUE) {
670  // We have recognised the resonance from those that are available
671  resTypInt = i;
672  return resTypInt;
673  }
674  ++i;
675  }
676 
677  return resTypInt;
678 }
679 
680 void LauResonanceMaker::printAll( std::ostream& stream ) const
681 {
682  for ( std::vector<LauResonanceInfo*>::const_iterator iter = resInfo_.begin(); iter != resInfo_.end(); ++iter ) {
683  stream << (**iter) << std::endl;
684  }
685 }
686 
Class for defining a model independent partial wave component.
File containing declaration of LauNRAmplitude class.
Bool_t fixed() const
Check whether the parameter is fixed or floated.
LauBlattWeisskopfFactor::BlattWeisskopfCategory getBWCategory() const
Retrieve the BW category of the resonant particle.
Class for defining the LASS resonance model.
Definition: LauLASSRes.hh:31
BWRadiusFixedCategoryMap bwFixRadii_
void setBarrierRadii(LauBlattWeisskopfFactor *resFactor, LauBlattWeisskopfFactor *parFactor)
Set the form factor model and parameters.
BWRadiusCategoryMap bwDefaultRadii_
File containing declaration of LauResonanceInfo class.
ClassImp(LauAbsCoeffSet)
BWFactorCategoryMap bwFactors_
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.
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.
LauAbsResonance * getResonance(const LauDaughters *daughters, const TString &resName, const Int_t resPairAmpInt, const LauAbsResonance::LauResonanceModel resType, const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory=LauBlattWeisskopfFactor::Default, const LauBlattWeisskopfFactor::BarrierType bwType=LauBlattWeisskopfFactor::BWPrimeBarrier)
Create a resonance.
const LauParameter * getRadiusParameter() const
Retrieve the radius parameter.
File containing declaration of LauGounarisSakuraiRes class.
File containing declaration of LauBelleSymNR class.
File containing declaration of LauBreitWignerRes class.
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 non resonant part of the LASS model.
Definition: LauLASSNRRes.hh:31
Class for defining the Sigma resonance model.
Definition: LauSigmaRes.hh:31
Singleton factory class for creating resonances.
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 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.
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:34
Class for defining the relativistic Breit-Wigner resonance model.
std::vector< LauBlattWeisskopfFactor * > bwIndepFactors_
The Blatt-Weisskopf factor objects for resonances in the independent category.
BarrierType getBarrierType() const
Retrieve the barrier type.
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.
LauResonanceMaker()
Constructor.
BarrierType
Define the allowed types of barrier factors.
LauBlattWeisskopfFactor * getBWFactor(const LauBlattWeisskopfFactor::BlattWeisskopfCategory bwCategory, const LauResonanceInfo *resInfo, const LauBlattWeisskopfFactor::BarrierType bwType)
Retrieve Blatt-Weisskopf factor for the given category.
Class for defining a uniform nonresonant amplitude.
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.
Class that implements the Blatt-Weisskopf barrier factor.
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.
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.
File containing declaration of LauModIndPartWave class.
File containing declaration of LauKappaRes class.
File containing declaration of LauFlatteRes class.
File containing declaration of LauFlatNR class.