BmnRoot
Loading...
Searching...
No Matches
MpdDecayerPyt8.cxx
Go to the documentation of this file.
1// -------------------------------------------------------------------------
2// ----- MpdDecayerPyt8 header file -----
3// ----- Created 23/07/2020 -----
4// ----- A. Zinchenko -----
5// ----- External decayer for MPD using Pythia8 -----
6// ----- (adapted from TPythia8Decayer) -----
7// -------------------------------------------------------------------------
8
9#include "MpdDecayerPyt8.h"
10
11#include <Pythia8/Pythia.h>
12
13#include <TClonesArray.h>
14#include <TDatabasePDG.h>
15#include <TF1.h>
16#include <TLorentzVector.h>
17#include <TObjString.h>
18#include <TParticle.h>
19#include <TPythia8.h>
20#include <TString.h>
21#include <TSystem.h>
22#include <TVirtualMC.h>
23
24using namespace Pythia8;
25
26class MyRandomClass : public RndmEngine {
27 double flat() { return gRandom->Rndm(); }
28};
29
30MpdDecayerPyt8 *MpdDecayerPyt8::fgInstance = 0;
31
34
35MpdDecayerPyt8 *MpdDecayerPyt8::Instance()
36{
37 if (!fgInstance) fgInstance = new MpdDecayerPyt8;
38 return fgInstance;
39}
40
43
44MpdDecayerPyt8::MpdDecayerPyt8() : fPythia8(new TPythia8()), fDebug(0), fGlobalPolar(0), fLambda(NULL), fRandom(NULL)
45{
46
47 fPythia8->Pythia8()->readString("SoftQCD:elastic = on");
48
49 // Disable cascade decays in one go
50 fPythia8->Pythia8()->readString("ParticleDecays:limitRadius = on");
51 fPythia8->Pythia8()->readString("ParticleDecays:rMax = 0.");
52
53 fPythia8->Pythia8()->particleData.readString("59:m0 = 100.1396"); //??? - to avoid crash with Geant4
54
55 fPythia8->Pythia8()->particleData.list(223);
56
57 // Modify branching ratios
58 ChangeBranchings();
59
60// Random number generator
61#if PYTHIA_VERSION_INTEGER < 8310
62 RndmEngine *rndm = new MyRandomClass();
63#else
64 auto rndm = make_shared<MyRandomClass>();
65#endif
66 fPythia8->Pythia8()->setRndmEnginePtr(rndm);
67
68 fPythia8->Pythia8()->init();
69
70 Init();
71}
72
75
76void MpdDecayerPyt8::Init()
77{
78
79 static Bool_t init = kFALSE;
80 if (init) return;
81 init = kTRUE;
82 // ForceDecay();
83 fParticles = new TClonesArray("TParticle");
84
85 // Lambda branching to p+\pi-
86 fPythia8->Pythia8()->particleData.list(3122);
87 // ParticleDataEntry* part = fPythia8->Pythia8()->particleData.particleDataEntryPtr(3122); // Lambda
88 fLambda = fPythia8->Pythia8()->particleData.particleDataEntryPtr(3122); // Lambda
89 DecayChannel &channel = fLambda->channel(0);
90 fBranch = channel.bRatio();
91
92 // Random number generator
93 fRandom = new TRandom(0); // time-dependent seed
94 // fRandom = gRandom;
95
96 // fPythia8->Pythia8()->particleData.isResonance(113,kTRUE);
97 fPythia8->Pythia8()->particleData.list(113);
98}
99
102
103void MpdDecayerPyt8::Decay(Int_t pdg, TLorentzVector *p)
104{
105
106 //{ cout << " !!! Out !!! " << endl; exit(0); }
107 // Reset internal particle container
108 fParticles->Delete();
109 fSourceFlag = kPythia;
110 fMother.SetXYZT(0, 0, 0, 0);
111
112 if (!p) return;
113
114 TRandom *saveRandom = gRandom;
115 gRandom = fRandom;
116
117 TParticle *part = gMC->GetStack()->GetCurrentTrack();
118
119 if (fMothersPdg.find(pdg) == fMothersPdg.end()) {
120 // Particle is not defined
121 new ((*fParticles)[0]) TParticle(*part); // store mother particle
122 fSourceFlag = kCustom;
123 return;
124 }
125
126 TVector3 polar;
127 part->GetPolarisation(polar);
128 Bool_t polarFlag = kFALSE;
129
130 if (TMath::Abs(polar.X()) > 0.0001 || TMath::Abs(polar.Z()) > 0.0001) {
131 // Polarized lambda - simulate anysotropic decay only to p + \pi-
132 polarFlag = kTRUE;
133 if (gRandom->Rndm() < fBranch) {
134 // Anysotropic decay
135 // cout << " ----------------- " << fBranch << endl;
136 new ((*fParticles)[0]) TParticle(*part); // store mother particle
137 Gdecay(pdg, p);
138 gRandom = saveRandom;
139 return;
140 } else {
141 // cout << " ++++++++++ " << fBranch << endl;
142 // Force the other channels
143 DecayChannel &channel = fLambda->channel(0);
144 channel.bRatio(0.0);
145 }
146 }
147
148 ClearEvent();
149 AppendParticle(pdg, p);
150 Int_t idPart = fPythia8->Pythia8()->event[0].id();
151 fPythia8->Pythia8()->particleData.mayDecay(idPart, kTRUE);
152 fPythia8->Pythia8()->moreDecays();
153 if (fDebug > 0) fPythia8->EventListing();
154 if (polarFlag) {
155 DecayChannel &channel = fLambda->channel(0);
156 channel.bRatio(fBranch);
157 }
158
159 gRandom = saveRandom;
160}
161
164
165Int_t MpdDecayerPyt8::ImportParticles(TClonesArray *particles)
166{
167
168 Int_t npart = 0;
169
170 if (fSourceFlag == kPythia) {
171 fParticles->Delete();
172 npart = fPythia8->ImportParticles(fParticles, "All");
173 TLorentzVector lvDaugh;
174
175 for (Int_t j = 0; j < npart; ++j) {
176 TParticle *part = (TParticle *)fParticles->UncheckedAt(j);
177 part->ProductionVertex(lvDaugh);
178 lvDaugh += fMother;
179 part->SetProductionVertex(lvDaugh);
180 new ((*particles)[j]) TParticle(*part);
181 }
182 } else {
183 npart = fParticles->GetEntriesFast();
184 for (Int_t j = 0; j < npart; ++j) new ((*particles)[j]) TParticle(*((TParticle *)fParticles->UncheckedAt(j)));
185 }
186
187 // return (fPythia8->ImportParticles(particles, "All"));
188 return npart;
189}
190
193
194void MpdDecayerPyt8::SetForceDecay(Int_t /*type*/)
195{
196 printf("SetForceDecay not yet implemented !\n");
197}
200
201void MpdDecayerPyt8::ForceDecay()
202{
203 // printf("ForceDecay not yet implemented !\n");
204}
206
207Float_t MpdDecayerPyt8::GetPartialBranchingRatio(Int_t /*ipart*/)
208{
209 return 0.0;
210}
213
214Float_t MpdDecayerPyt8::GetLifetime(Int_t pdg)
215{
216 return (fPythia8->Pythia8()->particleData.tau0(pdg) * 3.3333e-12);
217}
218
221
222void MpdDecayerPyt8::ReadDecayTable() {}
223
226
227void MpdDecayerPyt8::AppendParticle(Int_t pdg, TLorentzVector *p)
228{
229 fPythia8->Pythia8()->event.append(pdg, 11, 0, 0, p->Px(), p->Py(), p->Pz(), p->E(), p->M());
230}
231
234
235void MpdDecayerPyt8::ClearEvent()
236{
237 fPythia8->Pythia8()->event.clear();
238}
239
242
243void MpdDecayerPyt8::Gdecay(Int_t idpart, TLorentzVector *p)
244{
245
246 Double_t pcm[2][4] = {{0}, {0}};
247 Double_t xm0 = p->M(), xm1 = TDatabasePDG::Instance()->GetParticle(2212)->Mass();
248 Double_t xm2 = TDatabasePDG::Instance()->GetParticle(-211)->Mass();
249
250 Gdeca2(xm0, xm1, xm2, pcm);
251
252 // First, second decay products.
253
254 TLorentzVector daughter1, daughter2;
255 daughter1.SetPxPyPzE(pcm[0][0], pcm[0][1], pcm[0][2], pcm[0][3]);
256 daughter2.SetPxPyPzE(pcm[1][0], pcm[1][1], pcm[1][2], pcm[1][3]);
257
258 // Boost to lab frame
259 // TVector3 boost = p->BoostVector(); //AZ - this is wrong !!!!
260 TVector3 boost(0, 0, p->P() / p->E());
261 daughter1.Boost(boost);
262 daughter2.Boost(boost);
263 daughter1.Transform(fRotation);
264 daughter2.Transform(fRotation);
265
266 TLorentzVector pos;
267 gMC->TrackPosition(pos);
268 // pos.Print();
269
270 Int_t npart = fParticles->GetEntriesFast();
271 new ((*fParticles)[npart++]) TParticle(2212, 1, 1, -1, 0, 0, daughter1, pos);
272 new ((*fParticles)[npart]) TParticle(-211, 1, 1, -1, 0, 0, daughter2, pos);
273
274 fSourceFlag = kCustom;
275
276 // Modify mother particle parameters
277 TParticle *mother = (TParticle *)fParticles->UncheckedAt(0);
278 mother->SetFirstMother(1);
279 mother->SetLastMother(-1);
280 mother->SetFirstDaughter(2);
281 mother->SetLastDaughter(3);
282 mother->SetStatusCode(11);
283}
284
287
288void MpdDecayerPyt8::Gdeca2(Double_t xm0, Double_t xm1, Double_t xm2, Double_t pcm[2][4])
289{
290
291 Double_t random[2], pvert[3], costh = 0.0, sinth = 0.0, phi = 0.0;
292
293 Double_t e1 = (xm0 * xm0 + xm1 * xm1 - xm2 * xm2) / (2. * xm0);
294 Double_t p1 = TMath::Sqrt(TMath::Abs((e1 - xm1) * (e1 + xm1)));
295
296 static TRandom *myRandom = NULL;
297 if (!myRandom) {
298 myRandom = new TRandom();
299 myRandom->SetSeed(0);
300 }
301 TRandom *oldRandom = gRandom;
302 gRandom = myRandom;
303
304 gRandom->RndmArray(2, random);
305
306 // Sanity check - should not happen (legacy code)
307 TParticle *part = gMC->GetStack()->GetCurrentTrack();
308 if (part->GetPdgCode() != 3122) {
309 cout << " ??? Not Lambda - exit " << part->GetPdgCode() << endl;
310 exit(0);
311 }
312
313 TLorentzVector lorV;
314 part->Momentum(lorV);
315 lorV.Vect().GetXYZ(pvert);
316
317 TVector3 polar;
318 part->GetPolarisation(polar);
319 if (fGlobalPolar) polar.SetMag(part->GetWeight()); // particle polarization value is passed as its weight
320
321 if (TMath::Abs(polar.X()) < 0.0001 && TMath::Abs(polar.Z()) < 0.0001 && TMath::Abs(polar.Y()) < 0.0001) {
322 costh = 2. * random[0] - 1.0;
323 phi = TMath::TwoPi() * random[1];
324 } else {
325 // Anisotropic decay angular distribution (0.5*(1+alpha*P*cos)).
326 // Anisotropy (pvert, random, polar.X(), phi, costh);
327 Anisotropy(pvert, random, polar, phi, costh);
328 }
329
330 if (TMath::Abs(costh) >= 1.0)
331 costh = TMath::Sign(1.0, costh);
332 else
333 sinth = TMath::Sqrt((1. - costh) * (1. + costh));
334
335 // Polar co-ordinates to momentum components.
336
337 pcm[0][0] = p1 * sinth * TMath::Cos(phi);
338 pcm[0][1] = p1 * sinth * TMath::Sin(phi);
339 pcm[0][2] = p1 * costh;
340 pcm[0][3] = e1;
341
342 // Generate second decay product.
343
344 pcm[1][0] = -pcm[0][0];
345 pcm[1][1] = -pcm[0][1];
346 pcm[1][2] = -pcm[0][2];
347 pcm[1][3] = TMath::Sqrt(pcm[1][0] * pcm[1][0] + pcm[1][1] * pcm[1][1] + pcm[1][2] * pcm[1][2] + xm2 * xm2);
348
349 gRandom = oldRandom;
350}
351
354
355// void MpdDecayerPyt8::Anisotropy (Double_t *pvert, Double_t *rndm, Double_t polar, Double_t &phi, Double_t &costh)
356void MpdDecayerPyt8::Anisotropy(Double_t *pvert, Double_t *rndm, TVector3 &polar3, Double_t &phi, Double_t &costh)
357{
358 // Simulate anisotropic (due to polarization) decay of lambda
359
360 // exit(0);
361 // std::ofstream outfile ("costh.txt");
362 // freopen ("costh.txt", "w", stdout);
363
364 Double_t polar = 0.0;
365 if (fGlobalPolar)
366 polar = polar3.Mag(); // global polarization - CHECK THIS
367 else
368 polar = polar3.X(); // transverse polar.
369
370 // const Double_t alpha = 0.642;
371 // static TF1 f("f","1+0.642*[0]*x",-1,1);
372 static TF1 f("f", "1+0.732*[0]*x", -1, 1); // AZ 10.06.21 - new value
373 f.SetParameter(0, polar);
374 f.SetNpx(1000);
375
376 Double_t costhe = f.GetRandom();
377
378 Double_t sinth = 0.0;
379 if (TMath::Abs(costhe) >= 1.0)
380 costhe = TMath::Sign(1.0, costhe);
381 else
382 sinth = TMath::Sqrt(1.0 - costhe * costhe);
383 phi = TMath::TwoPi() * rndm[1];
384 // std::cout << " Cos: " << costhe << " " << phi << " " << pvert[0] << " " << pvert[1] << " " << pvert[2] <<
385 // std::endl;
386
387 // Compute normal vector
388 TVector3 beam(0.0, 0.0, 1.0);
389 TVector3 lambda(pvert[0], pvert[1], pvert[2]);
390 TVector3 norm;
391 if (fGlobalPolar)
392 norm = polar3.Unit(); // global polarization
393 else
394 norm = beam.Cross(lambda).Unit(); // transverse polarization
395
396 // Unit vector with theta and phi w.r.t. normal
397 TVector3 unit(sinth * TMath::Cos(phi), sinth * TMath::Sin(phi), costhe);
398
399 // Coordinate system transformation
400 // (from lambda to lab)
401 TVector3 lambU = lambda.Unit(), lambZ(0, 0, 1), lambX(1, 0, 0), lambY(0, 1, 0);
402 lambZ.RotateUz(lambU);
403 lambX.RotateUz(lambU);
404 lambY.RotateUz(lambU);
405 TRotation rotL;
406 rotL.RotateAxes(lambX, lambY, lambZ); // transformation to lab. system
407 fRotation = rotL;
408 rotL.Invert(); // from lab. to lambda
409 unit.RotateUz(norm); // to lab. system
410 unit.Transform(rotL); // to lambda system
411 /*
412 TVector3 normZ(0,0,1), normX(1,0,0), normY(0,1,0);
413 normZ.RotateUz(norm);
414 normX.RotateUz(norm);
415 normY.RotateUz(norm);
416 TRotation rotN;
417 rotN.RotateAxes(normX,normY,normZ); // transformation to lab. system
418 unit.Transform(rotN); // to lab. system
419 unit.Transform(rotL); // to lambda system
420 */
421 costh = TMath::Cos(unit.Theta());
422 phi = unit.Phi();
423 // std::cout << costh << " " << phi << std::endl;
424 // outfile << i << " " << costh << endl;
425 // outfile.close();
426}
427
430
431void MpdDecayerPyt8::AddMotherPdg(Int_t pdg)
432{
433
434 if (fMothersPdg.find(pdg) != fMothersPdg.end()) return;
435 fMothersPdg.insert(pdg);
436 gMC->SetUserDecay(pdg);
437}
438
441
442void MpdDecayerPyt8::Decay(TParticle *part)
443{
444 // Decay function
445
446 //{ cout << " !!! Out !!! " << endl; exit(0); }
447 // Reset internal particle container
448 fParticles->Delete();
449 fSourceFlag = kPythia;
450
451 /*
452 if (!p) return;
453 p->ProductionVertex(fMother);
454 TPythia6 *pythia = TPythia6::Instance();
455
456 pythia->Py1ent(0, p->GetPdgCode(), p->Energy(), p->Theta(), p->Phi());
457 pythia->GetPrimaries();
458 */
459 if (!part) return;
460 // part->Print();
461
462 TRandom *saveRandom = gRandom;
463 gRandom = fRandom;
464
465 part->ProductionVertex(fMother);
466
467 Int_t pdg = part->GetPdgCode();
468 TLorentzVector p;
469 part->Momentum(p);
470
471 if (part->GetFirstMother() < 0) {
472 // Primary particle - select mass
473 Double_t mass = fPythia8->Pythia8()->particleData.mSel(pdg);
474 // cout << "mmm " << mass << " " << part->GetFirstMother() << endl;
475 p.SetE(TMath::Sqrt(p.P() * p.P() + mass * mass));
476 }
477 // p.Print();
478
479 ClearEvent();
480 AppendParticle(pdg, &p);
481 Int_t idPart = fPythia8->Pythia8()->event[0].id();
482 fPythia8->Pythia8()->particleData.mayDecay(idPart, kTRUE);
483 fPythia8->Pythia8()->moreDecays();
484 if (fDebug > 0) fPythia8->EventListing();
485
486 gRandom = saveRandom;
487}
488
489//__________________________________________________________________________
490
491void MpdDecayerPyt8::ChangeBranchings()
492{
493 // Change branching ratios of some decay channels
494 // (using description file)
495
496 TString path = gSystem->ExpandPathName("$VMCWORKDIR");
497 path += "/gconfig/MpdDecayConfig.txt";
498 ifstream fin(path.Data());
499
500 TString chline;
501 TObjArray *tokens = NULL;
502
503 while (1) {
504 chline.ReadLine(fin);
505 if (fin.eof()) break;
506 cout << chline << endl;
507 if (chline.Contains("#")) continue; // comment
508 // Tokenize line
509 // AZ-161123 tokens = chline.Tokenize(";"); // ";" - channel terminator / separator
510 tokens = chline.Tokenize(" "); // " " - channel separator
511 // if (tokens) cout << " aaa " << tokens->GetEntriesFast() << endl;
512 Bool_t exclusive = kTRUE;
513 if (chline.Contains("*")) exclusive = kFALSE;
514 ChangeParticleBr(tokens, exclusive);
515 if (tokens) {
516 tokens->Delete();
517 delete tokens;
518 tokens = NULL;
519 }
520 }
521}
522
523//__________________________________________________________________________
524
525void MpdDecayerPyt8::ChangeParticleBr(TObjArray *channels, Bool_t exclusive)
526{
527 // Change branching ratios of selected channels for single particle
528
529 Int_t nchan = channels->GetEntriesFast(); // number of channels to be changed
530 TString last = ((TObjString *)channels->Last())->String(); // AZ-121123 - safety
531 if (!last.Contains(":")) --nchan; // AZ-121123 - safety
532
533 Int_t pdg = ((TObjString *)channels->UncheckedAt(0))->String().Atoi();
534 // cout << pdg << endl;
535 ParticleDataEntryPtr part = fPythia8->Pythia8()->particleData.particleDataEntryPtr(pdg); // mother particle
536 Int_t ndecays = part->sizeChannels();
537
538 TObjArray *tokens = NULL;
539 Double_t changedBrs = 0.0, unchangedBrs = 1.0;
540 vector<Int_t> decayInds;
541 vector<Double_t> bRatios;
542 set<int> activeChan; // AZ-211123
543
544 for (Int_t ich = 0; ich < nchan; ++ich) {
545 TString chline = ((TObjString *)channels->UncheckedAt(ich))->String();
546 cout << chline << endl;
547 if (chline.Contains("*")) continue; // flag for inclusive decays
548 // Tokenize channel
549 tokens = chline.Tokenize(":");
550
551 // Get daughters' PDGs
552 Int_t nprongs = tokens->GetEntriesFast() - 2;
553 vector<Int_t> pdgs;
554 for (Int_t id = 0; id < nprongs; ++id) {
555 pdgs.push_back(((TObjString *)tokens->UncheckedAt(1 + id))->String().Atoi());
556 }
557 TString scale = ((TObjString *)tokens->Last())->String();
558 scale.Remove(0, 1);
559 Double_t scaleBr = scale.Atof();
560
561 // Pick channel requested
562 for (Int_t idecay = 0; idecay < ndecays; ++idecay) {
563 DecayChannel &channel = part->channel(idecay);
564 if (channel.multiplicity() != nprongs) continue;
565 // Check daughters
566 Int_t ok = 0;
567 for (Int_t id = 0; id < nprongs; ++id) {
568 if (channel.product(id) != pdgs[id]) break;
569 ++ok;
570 }
571 if (ok != nprongs) continue;
572 changedBrs += channel.bRatio() * scaleBr;
573 unchangedBrs -= channel.bRatio();
574 decayInds.push_back(idecay);
575 bRatios.push_back(channel.bRatio() * scaleBr);
576 activeChan.insert(idecay); // AZ-211123
577 // cout << channel.bRatio() << endl;
578 }
579 // cout << changedBrs << endl;
580
581 if (tokens) {
582 tokens->Delete();
583 delete tokens;
584 tokens = NULL;
585 }
586 } // for (Int_t ich = 0; ich < nchan;
587
588 part->rescaleBR((1.0 - changedBrs) / unchangedBrs);
589 if (!exclusive) --nchan;
590
591 for (Int_t ich = 0; ich < nchan; ++ich) {
592 DecayChannel &channel = part->channel(decayInds[ich]);
593 channel.bRatio(bRatios[ich]);
594 }
595
596 if (exclusive) {
597 // Disable inactive decay channels
598
599 for (Int_t idecay = 0; idecay < ndecays; ++idecay) {
600 if (activeChan.find(idecay) == activeChan.end()) {
601 DecayChannel &channel = part->channel(idecay);
602 channel.bRatio(0.0);
603 }
604 }
605 }
606
607 // Check
608 Double_t totalBr = 0;
609 for (Int_t idecay = 0; idecay < ndecays; ++idecay) {
610 DecayChannel &channel = part->channel(idecay);
611 totalBr += channel.bRatio();
612 }
613 cout << " Total branching after rescaling: " << totalBr << endl;
614}
f
float f
Definition P4_F32vec4.h:21
MpdDecayerPyt8::ClearEvent
void ClearEvent()
Clear the event stack.
Definition MpdDecayerPyt8.cxx:235
MpdDecayerPyt8::AppendParticle
void AppendParticle(Int_t pdg, TLorentzVector *p)
Append a particle to the stack.
Definition MpdDecayerPyt8.cxx:227
MpdDecayerPyt8::Instance
static MpdDecayerPyt8 * Instance()
Get the singleton object.
Definition MpdDecayerPyt8.cxx:35
MpdDecayerPyt8::SetForceDecay
virtual void SetForceDecay(Int_t type)
Set forced decay mode.
Definition MpdDecayerPyt8.cxx:194
MpdDecayerPyt8::GetPartialBranchingRatio
virtual Float_t GetPartialBranchingRatio(Int_t ipart)
Definition MpdDecayerPyt8.cxx:207
MpdDecayerPyt8::Init
virtual void Init()
Initialize the decayer.
Definition MpdDecayerPyt8.cxx:76
MpdDecayerPyt8::MpdDecayerPyt8
MpdDecayerPyt8()
constructor
Definition MpdDecayerPyt8.cxx:44
MpdDecayerPyt8::ReadDecayTable
virtual void ReadDecayTable()
to read a decay table (not yet implemented)
Definition MpdDecayerPyt8.cxx:222
MpdDecayerPyt8::GetLifetime
virtual Float_t GetLifetime(Int_t kf)
return lifetime in seconds of teh particle with PDG number pdg
Definition MpdDecayerPyt8.cxx:214
MpdDecayerPyt8::ForceDecay
virtual void ForceDecay()
ForceDecay not yet implemented.
Definition MpdDecayerPyt8.cxx:201
MpdDecayerPyt8::ImportParticles
virtual Int_t ImportParticles(TClonesArray *particles)
import the decay products into particles array
Definition MpdDecayerPyt8.cxx:165
MpdDecayerPyt8::Decay
virtual void Decay(Int_t pdg, TLorentzVector *p)
Decay a single particle.
Definition MpdDecayerPyt8.cxx:103
MpdDecayerPyt8::AddMotherPdg
void AddMotherPdg(Int_t pdg)
Add PDG of particle to be decayed by this package.
Definition MpdDecayerPyt8.cxx:431
ParticleDataEntryPtr
Pythia8::ParticleDataEntry * ParticleDataEntryPtr
Definition MpdDecayerPyt8.h:26