GooFit  v2.1.3
Functions | Variables
SigGen.cpp File Reference
#include <TFile.h>
#include <TTree.h>
#include <goofit/Application.h>
#include <goofit/PDFs/basic/PolynomialPdf.h>
#include <goofit/PDFs/combine/AddPdf.h>
#include <goofit/PDFs/physics/DP4Pdf.h>
#include <goofit/UnbinnedDataSet.h>
#include <goofit/Variable.h>
#include <thrust/count.h>

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Variables

const fptype _mD0 = 1.8645
 
const fptype piPlusMass = 0.13957018
 
const fptype KmMass = .493677
 

Function Documentation

◆ main()

int main ( int  argc,
char **  argv 
)

Definition at line 23 of file SigGen.cpp.

References _mD0, GooFit::DecayInfo4::amplitudes, constantOne, constantZero, eventNumber, GooFit::DPPdf::GenerateSig(), GOOFIT_PARSE, KmMass, m12, GooFit::M_12, GooFit::M_12_3, GooFit::M_12_4, GooFit::M_13, GooFit::M_13_2, GooFit::M_24, GooFit::M_24_1, GooFit::M_24_3, GooFit::M_34, GooFit::M_34_2, GooFit::DecayInfo4::meson_radius, GooFit::DecayInfo4::particle_masses, piPlusMass, GooFit::DPPdf::setGenerationOffset(), and weights.

23  {
24  GooFit::Application app("Signal Generator Example", argc, argv);
25 
26  GOOFIT_PARSE(app);
27 
28  DecayInfo4 DK3P_DI;
29  DK3P_DI.meson_radius = 1.5;
30  DK3P_DI.particle_masses.push_back(_mD0);
31  DK3P_DI.particle_masses.push_back(piPlusMass);
32  DK3P_DI.particle_masses.push_back(piPlusMass);
33  DK3P_DI.particle_masses.push_back(KmMass);
34  DK3P_DI.particle_masses.push_back(piPlusMass);
35 
36  Variable RhoMass("rho_mass", 0.77526, 0.01, 0.7, 0.8);
37  Variable RhoWidth("rho_width", 0.1478, 0.01, 0.1, 0.2);
38  Variable KstarM("KstarM", 0.89581, 0.01, 0.9, 0.1);
39  Variable KstarW("KstarW", 0.0474, 0.01, 0.1, 0.2);
40  Variable f600M("f600M", 0.519, 0.01, 0.75, 0.85);
41  Variable f600W("f600W", 0.454, 0.01, 0.75, 0.85);
42  Variable a1M("a1M", 1.23, 0.01, 1.2, 1.3);
43  Variable a1W("a1W", 0.42, 0.01, 0.37, 0.47);
44  Variable K1M("K1M", 1.272, 0.01, 1.2, 1.3);
45  Variable K1W("K1W", 0.09, 0.01, 0.08, 0.1);
46  Variable K1430M("K1430M", 1.414, 0.01, 1.4, 1.5);
47  Variable K1430W("K1430W", .29, 0.01, 0.25, 0.35);
48 
49  // Spin factors: we have two due to the bose symmetrization of the two pi+
50  std::vector<SpinFactor *> SFKRS;
51  SFKRS.push_back(new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_S, 0, 1, 2, 3));
52  SFKRS.push_back(new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_S, 3, 1, 2, 0));
53 
54  std::vector<SpinFactor *> SFKRP;
55  SFKRP.push_back(new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_P, 0, 1, 2, 3));
56  SFKRP.push_back(new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_P, 3, 1, 2, 0));
57 
58  std::vector<SpinFactor *> SFKRD;
59  SFKRD.push_back(new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_D, 0, 1, 2, 3));
60  SFKRD.push_back(new SpinFactor("SF", SF_4Body::DtoV1V2_V1toP1P2_V2toP3P4_D, 3, 1, 2, 0));
61 
62  std::vector<SpinFactor *> SFKF;
63  SFKF.push_back(new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4, 2, 3, 0, 1));
64  SFKF.push_back(new SpinFactor("SF", SF_4Body::DtoVS_VtoP1P2_StoP3P4, 2, 0, 3, 1));
65 
66  std::vector<SpinFactor *> SFKK;
67  SFKK.push_back(new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4, 0, 1, 3, 2));
68  SFKK.push_back(new SpinFactor("SF", SF_4Body::DtoAP1_AtoSP2_StoP3P4, 3, 1, 0, 2));
69 
70  std::vector<SpinFactor *> SFK1R;
71  SFK1R.push_back(new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4, 3, 2, 0, 1));
72  SFK1R.push_back(new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4, 0, 2, 3, 1));
73 
74  std::vector<SpinFactor *> SFA1R;
75  SFA1R.push_back(new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4, 2, 3, 0, 1));
76  SFA1R.push_back(new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2_VtoP3P4, 2, 0, 3, 1));
77 
78  std::vector<SpinFactor *> SFA1RD;
79  SFA1RD.push_back(new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2Dwave_VtoP3P4, 2, 3, 0, 1));
80  SFA1RD.push_back(new SpinFactor("SF", SF_4Body::DtoAP1_AtoVP2Dwave_VtoP3P4, 2, 0, 3, 1));
81 
82  // Lineshapes, also for both pi+ configurations
83  std::vector<Lineshape *> LSKRS;
84  LSKRS.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_12));
85  LSKRS.push_back(new Lineshapes::RBW("K*(892)bar", KstarM, KstarW, 1, M_34));
86  LSKRS.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_24));
87  LSKRS.push_back(new Lineshapes::RBW("K*(892)bar", KstarM, KstarW, 1, M_13));
88 
89  std::vector<Lineshape *> LSKRP;
90  LSKRP.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_12));
91  LSKRP.push_back(new Lineshapes::RBW("K*(892)bar", KstarM, KstarW, 1, M_34));
92  LSKRP.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_24));
93  LSKRP.push_back(new Lineshapes::RBW("K*(892)bar", KstarM, KstarW, 1, M_13));
94 
95  std::vector<Lineshape *> LSKRD;
96  LSKRD.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_12));
97  LSKRD.push_back(new Lineshapes::RBW("K*(892)bar", KstarM, KstarW, 1, M_34));
98  LSKRD.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_24));
99  LSKRD.push_back(new Lineshapes::RBW("K*(892)bar", KstarM, KstarW, 1, M_13));
100 
101  std::vector<Lineshape *> LSKF;
102  LSKF.push_back(new Lineshapes::RBW("K*(892)bar", KstarM, KstarW, 1, M_34));
103  LSKF.push_back(new Lineshapes::Bugg("f600", f600M, f600W, 0, M_12));
104  LSKF.push_back(new Lineshapes::RBW("K*(892)bar", KstarM, KstarW, 1, M_13));
105  LSKF.push_back(new Lineshapes::Bugg("f600", f600M, f600W, 0, M_24));
106 
107  std::vector<Lineshape *> LSKK;
108  LSKK.push_back(new Lineshapes::SBW("K(1)(1270)bar", K1M, K1W, 1, M_34_2));
109  LSKK.push_back(new Lineshapes::LASS("K(0)*(1430)bar", K1430M, K1430W, 0, M_34));
110  LSKK.push_back(new Lineshapes::SBW("K(1)(1270)bar2", K1M, K1W, 1, M_13_2));
111  LSKK.push_back(new Lineshapes::LASS("K(0)*(1430)bar2", K1430M, K1430W, 0, M_13));
112 
113  std::vector<Lineshape *> LSK1R;
114  LSK1R.push_back(new Lineshapes::SBW("K(1)(1270)bar", K1M, K1W, 0, M_12_3));
115  LSK1R.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_12));
116  LSK1R.push_back(new Lineshapes::SBW("K(1)(1270)bar", K1M, K1W, 0, M_24_3));
117  LSK1R.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_24));
118 
119  std::vector<Lineshape *> LSA1R;
120  LSA1R.push_back(new Lineshapes::SBW("a(1)(1260)+", a1M, a1W, 0, M_12_4));
121  LSA1R.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_12));
122  LSA1R.push_back(new Lineshapes::SBW("a(1)(1260)+", a1M, a1W, 0, M_24_1));
123  LSA1R.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_24));
124 
125  std::vector<Lineshape *> LSA1RD;
126  LSA1RD.push_back(new Lineshapes::SBW("a(1)(1260)+", a1M, a1W, 2, M_12_4));
127  LSA1RD.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_12));
128  LSA1RD.push_back(new Lineshapes::SBW("a(1)(1260)+", a1M, a1W, 2, M_24_1));
129  LSA1RD.push_back(new Lineshapes::RBW("rho(770)", RhoMass, RhoWidth, 1, M_24));
130 
131  // the very last parameter means that we have two permutations. so the first half of the Lineshapes
132  // and the first half of the spinfactors are amplitude 1, rest is amplitude 2
133  // This means that it is important for symmetrized amplitueds that the spinfactors and lineshapes are in the "right"
134  // order
135  Amplitude *Bose_symmetrized_AMP_S = new Amplitude(
136  "K*(892)rho(770)_S", Variable("amp_real1", -0.115177), Variable("amp_imag1", 0.153976), LSKRS, SFKRS, 2);
137  Amplitude *Bose_symmetrized_AMP_P = new Amplitude(
138  "K*(892)rho(770)_P", Variable("amp_real2", -0.0298697), Variable("amp_imag2", -0.0722874), LSKRP, SFKRP, 2);
139  Amplitude *Bose_symmetrized_AMP_D = new Amplitude(
140  "K*(892)rho(770)_D", Variable("amp_real3", -0.452212), Variable("amp_imag3", 0.426521), LSKRD, SFKRD, 2);
141 
142  Amplitude *Bose_symmetrized_KF
143  = new Amplitude("KF", Variable("amp_real3", 0.0120787), Variable("amp_imag3", -0.0332525), LSKF, SFKF, 2);
144  Amplitude *Bose_symmetrized_KK
145  = new Amplitude("LSKK", Variable("amp_real3", 0.0109033), Variable("amp_imag3", -0.00186219), LSKK, SFKK, 2);
146  Amplitude *Bose_symmetrized_K1R
147  = new Amplitude("LSK1R", Variable("amp_real3", -0.10728), Variable("amp_imag3", -0.130213), LSK1R, SFK1R, 2);
148  Amplitude *Bose_symmetrized_A1R
149  = new Amplitude("LSA1R", Variable("amp_real3", 1.0), Variable("amp_imag3", 0.0), LSA1R, SFA1R, 2);
150  Amplitude *Bose_symmetrized_A1RD
151  = new Amplitude("LSA1RD", Variable("amp_real3", -0.94921), Variable("amp_imag3", -1.73407), LSA1RD, SFA1RD, 2);
152 
153  DK3P_DI.amplitudes.push_back(Bose_symmetrized_KF);
154  DK3P_DI.amplitudes.push_back(Bose_symmetrized_AMP_S);
155  DK3P_DI.amplitudes.push_back(Bose_symmetrized_AMP_P);
156  DK3P_DI.amplitudes.push_back(Bose_symmetrized_AMP_D);
157  DK3P_DI.amplitudes.push_back(Bose_symmetrized_KK);
158  DK3P_DI.amplitudes.push_back(Bose_symmetrized_K1R);
159  DK3P_DI.amplitudes.push_back(Bose_symmetrized_A1R);
160  DK3P_DI.amplitudes.push_back(Bose_symmetrized_A1RD);
161 
162  Observable m12("m12", 0, 3);
163  Observable m34("m34", 0, 3);
164  Observable cos12("cos12", -1, 1);
165  Observable cos34("m12", -1, 1);
166  Observable phi("phi", -3.5, 3.5);
167  EventNumber eventNumber("eventNumber");
168  Variable constantOne("constantOne", 1);
169  Variable constantZero("constantZero", 0);
170 
171  vector<Observable> observables;
172  vector<Variable> coefficients;
173  vector<Variable> offsets;
174 
175  observables.push_back(m12);
176  observables.push_back(m34);
177  observables.push_back(cos12);
178  observables.push_back(cos34);
179  observables.push_back(phi);
180  observables.push_back(eventNumber);
181  offsets.push_back(constantZero);
182  offsets.push_back(constantZero);
183  coefficients.push_back(constantOne);
184 
185  PolynomialPdf *eff = new PolynomialPdf("constantEff", observables, coefficients, offsets, 0);
186  DPPdf *dp = new DPPdf("test", observables, DK3P_DI, eff, 5);
187 
188  TFile *file = new TFile("SigGen.root", "RECREATE");
189  TTree *tree = new TTree("events", "events");
190 
191  double tm12, tm34, tc12, tc34, tphi;
192  tree->Branch("m12", &tm12, "m12/D");
193  tree->Branch("m34", &tm34, "m34/D");
194  tree->Branch("c12", &tc12, "c12/D");
195  tree->Branch("c34", &tc34, "c34/D");
196  tree->Branch("phi", &tphi, "phi/D");
197 
198  for(int k = 0; k < 4; ++k) {
199  int numEvents = 1e6;
200  dp->setGenerationOffset(k * numEvents);
201  auto tuple = dp->GenerateSig(numEvents);
202 
203  auto particles = std::get<0>(tuple);
204  auto variables = std::get<1>(tuple);
205  auto weights = std::get<2>(tuple);
206  auto flags = std::get<3>(tuple);
207  long accepted = thrust::count_if(flags.begin(), flags.end(), thrust::identity<bool>());
208  fmt::print(
209  stderr, "Using accept-reject method would leave you with {} out of {} events\n", accepted, numEvents);
210 
211  for(int i = 0; i < weights.size(); ++i) {
212  if(flags[i] == 1) {
213  // printf("%.5g %.5g %.5g %.5g %.5g %.5g %.5g\n", (*(variables[0]))[i], (*(variables[1]))[i],
214  // (*(variables[2]))[i], (*(variables[3]))[i], (*(variables[4]))[i], weights[i], flags[i]);
215  tm12 = (*(variables[0]))[i];
216  tm34 = (*(variables[1]))[i];
217  tc12 = (*(variables[2]))[i];
218  tc34 = (*(variables[3]))[i];
219  tphi = (*(variables[4]))[i];
220  tree->Fill();
221  }
222  }
223 
224  delete variables[0];
225  delete variables[1];
226  delete variables[2];
227  delete variables[3];
228  delete variables[4];
229 
230  delete particles[0];
231  delete particles[1];
232  delete particles[2];
233  delete particles[3];
234  }
235 
236  tree->Write();
237  file->Close();
238  return 0;
239 }
KmMass
const fptype KmMass
Definition: SigGen.cpp:21
weights
std::vector< Variable > weights
Definition: pipipi0DPFit.cpp:151
constantZero
Variable constantZero("constantZero", 0)
GooFit::DPPdf::setGenerationOffset
__host__ void setGenerationOffset(int off)
Definition: DP4Pdf.h:44
m12
Observable * m12
Definition: pipipi0DPFit.cpp:69
GooFit::Observable
Special class for observables. Used in DataSets.
Definition: Variable.h:109
GooFit::DPPdf::GenerateSig
__host__ std::tuple< mcbooster::ParticlesSet_h, mcbooster::VariableSet_h, mcbooster::RealVector_h, mcbooster::RealVector_h > GenerateSig(unsigned int numEvents)
piPlusMass
const fptype piPlusMass
Definition: SigGen.cpp:20
constantOne
Variable constantOne("constantOne", 1)
eventNumber
EventNumber * eventNumber
Definition: pipipi0DPFit.cpp:71
GooFit::DecayInfo4::particle_masses
std::vector< fptype > particle_masses
Definition: DalitzPlotHelpers.h:98
GOOFIT_PARSE
#define GOOFIT_PARSE(app,...)
Definition: Application.h:11
GooFit::DecayInfo4::amplitudes
std::vector< Amplitude * > amplitudes
Definition: DalitzPlotHelpers.h:101
_mD0
const fptype _mD0
Definition: SigGen.cpp:19

Variable Documentation

◆ _mD0

const fptype _mD0 = 1.8645

Definition at line 19 of file SigGen.cpp.

Referenced by main().

◆ KmMass

const fptype KmMass = .493677

Definition at line 21 of file SigGen.cpp.

Referenced by main().

◆ piPlusMass

const fptype piPlusMass = 0.13957018

Definition at line 20 of file SigGen.cpp.

Referenced by main().

Generated on Sat Apr 21 2018 20:13:18 for GooFit by   doxygen 1.8.13