C++ Interface to Tauola
Tauola.cxx
1#include <fstream>
2#include <cstring>
3#include <vector>
4#include "Log.h"
5#include "Tauola.h"
6#include "TauolaEvent.h"
7
8namespace Tauolapp
9{
10bool Tauola::m_is_initialized = false;
11int Tauola::m_pdg_id = 15;
12int Tauola::m_firstDecayMode = 0;
13int Tauola::m_secondDecayMode = 0;
14bool Tauola::m_rad = true;
15double Tauola::m_rad_cut_off = 0.001;
16double Tauola::m_iniphy = 0.1;
17double Tauola::m_higgs_scalar_pseudoscalar_mix = M_PI/4;
18int Tauola::m_higgs_scalar_pseudoscalar_pdg = 35;
19int Tauola::m_helPlus = 0;
20int Tauola::m_helMinus = 0;
21double Tauola::m_wtEW = 0.0;
22double Tauola::m_wtEW0 = 0.0;
23double Tauola::table11A[NS1][NCOS][4][4] = {{{{0.0}}}};
24double Tauola::table1A [NS1][NCOS][4][4] = {{{{0.0}}}};
25double Tauola::table2A [NS1][NCOS][4][4] = {{{{0.0}}}};
26double Tauola::wtable11A[NS1][NCOS] = {{0.0}};
27double Tauola::wtable1A [NS1][NCOS] = {{0.0}};
28double Tauola::wtable2A [NS1][NCOS] = {{0.0}};
29double Tauola::w0table11A[NS1][NCOS] = {{0.0}};
30double Tauola::w0table1A [NS1][NCOS] = {{0.0}};
31double Tauola::w0table2A [NS1][NCOS] = {{0.0}};
32
33double Tauola::table11B[NS2][NCOS][4][4] = {{{{0.0}}}};
34double Tauola::table1B [NS2][NCOS][4][4] = {{{{0.0}}}};
35double Tauola::table2B [NS2][NCOS][4][4] = {{{{0.0}}}};
36double Tauola::wtable11B[NS2][NCOS] = {{0.0}};
37double Tauola::wtable1B [NS2][NCOS] = {{0.0}};
38double Tauola::wtable2B [NS2][NCOS] = {{0.0}};
39double Tauola::w0table11B[NS2][NCOS] = {{0.0}};
40double Tauola::w0table1B [NS2][NCOS] = {{0.0}};
41double Tauola::w0table2B [NS2][NCOS] = {{0.0}};
42
43double Tauola::table11C[NS3][NCOS][4][4] = {{{{0.0}}}};
44double Tauola::table1C [NS3][NCOS][4][4] = {{{{0.0}}}};
45double Tauola::table2C [NS3][NCOS][4][4] = {{{{0.0}}}};
46double Tauola::wtable11C[NS3][NCOS] = {{0.0}};
47double Tauola::wtable1C [NS3][NCOS] = {{0.0}};
48double Tauola::wtable2C [NS3][NCOS] = {{0.0}};
49double Tauola::w0table11C[NS3][NCOS] = {{0.0}};
50double Tauola::w0table1C [NS3][NCOS] = {{0.0}};
51double Tauola::w0table2C [NS3][NCOS] = {{0.0}};
52
53double Tauola::sminA = 0;
54double Tauola::smaxA = 0;
55
56double Tauola::sminB = 0;
57double Tauola::smaxB = 0;
58
59double Tauola::sminC = 0;
60double Tauola::smaxC = 0;
61
62int Tauola::ion[3] = {0};
63double Tauola::tau_lifetime = .08711;
64double Tauola::momentum_conservation_threshold = 0.1;
65
66Tauola::Particles Tauola::spin_correlation;
67
68Tauola::MomentumUnits Tauola::momentumUnit = Tauola::DEFAULT_MOMENTUM;
69Tauola::LengthUnits Tauola::lengthUnit = Tauola::DEFAULT_LENGTH;
70
74
75int Tauola::buf_incoming_pdg_id = 0;
76int Tauola::buf_outgoing_pdg_id = 0;
77double Tauola::buf_invariant_mass_squared = -1.;
78double Tauola::buf_cosTheta = 0.;
79
80double Tauola::buf_R[4][4] = {{0.0}}; //density matrix
81
82double (*Tauola::randomDouble)() = Tauola::defaultRandomGenerator;
83void (*Tauola::redefineTauPlusProperties)(TauolaParticle *) = defaultRedPlus;
84void (*Tauola::redefineTauMinusProperties)(TauolaParticle *) = defaultRedMinus;
85
86/**************************************************************/
88{
89 inirchl_(&mode);
90}
91
93 return rand()*1./RAND_MAX;
94}
95
96void Tauola::setRandomGenerator(double (*gen)()){
97 if(gen==NULL) randomDouble = defaultRandomGenerator;
98 else randomDouble = gen;
99}
100
101void Tauola::defaultRedPlus(TauolaParticle *tau) {}
102void Tauola::defaultRedMinus(TauolaParticle *tau) {}
103
104void Tauola::setRedefineTauMinus( void (*fun)(TauolaParticle *) ){
105 redefineTauMinusProperties=fun;
106}
107
108void Tauola::setRedefineTauPlus ( void (*fun)(TauolaParticle *) ){
109 redefineTauPlusProperties=fun;
110}
111
112void Tauola::getBornKinematics(int *incoming_pdg_id, int *outgoing_pdg_id, double *invariant_mass_squared,double *cosTheta){
113 *incoming_pdg_id = buf_incoming_pdg_id;
114 *outgoing_pdg_id = buf_outgoing_pdg_id;
115 *invariant_mass_squared = buf_invariant_mass_squared;
116 *cosTheta = buf_cosTheta;
117 // m_R[0][0] to be added in next step;
118}
119
120void Tauola::setUnits(MomentumUnits m, LengthUnits l){
121 Tauola::momentumUnit = m;
122 Tauola::lengthUnit = l;
123}
124
126 tau_lifetime = t;
127}
128
130 printf("\n");
131 printf(" *************************************\n");
132 printf(" * TAUOLA C++ Interface v1.1.8 *\n");
133 printf(" *-----------------------------------*\n");
134 printf(" * *\n");
135 printf(" * (c) Nadia Davidson, (1,2) *\n");
136 printf(" * (c) Gizo Nanava, (3) *\n");
137 printf(" * Tomasz Przedzinski,(4) *\n");
138 printf(" * Elzbieta Richter-Was,(2,4) *\n");
139 printf(" * Zbigniew Was (2,5) *\n");
140 printf(" * *\n");
141 printf(" * 1) Unimelb, Melbourne, Australia *\n");
142 printf(" * 2) INP, Krakow, Poland *\n");
143 printf(" * 3) University Bonn, Germany *\n");
144 printf(" * 4) UJ, Krakow, Poland *\n");
145 printf(" * 5) CERN, Geneva, Switzerland *\n");
146 printf(" *************************************\n");
147
148 if(m_is_initialized) {
149 printf(" * *\n");
150 printf(" * *\n");
151 printf(" *WARNING: Tauola already initialized*\n");
152 printf(" * reinitialization skipped *\n");
153 printf(" * *\n");
154 printf(" * *\n");
155 printf(" *************************************\n");
156 return;
157 }
158
159 m_is_initialized = true;
160
161 // Turn on all spin correlations
162 spin_correlation.setAll(true);
163
164 // Ininitalize tauola-fortran
165 f_interface_tauolaInitialize(m_pdg_id,m_firstDecayMode,
166 m_secondDecayMode,m_rad,
167 m_rad_cut_off, m_iniphy);
168
169 //---------------------------------------------------------------------------
170 // Initialize SANC tables
171 //---------------------------------------------------------------------------
172 cout<<"Reading SANC input files."<<endl;
173
174 ifstream f("table1-1.txt");
175
176 if(!f.is_open()){
177 cout<<"File 'table1-1.txt' missing... skipped."<<endl;
178 }
179 else{
180 string buf;
181
182 cout<<"Reading file 'table1-1.txt'..."<<endl;
183
184 int dbuf1,dbuf2,dbuf3,dbufcos;
185 f>>buf>>dbuf1>>dbuf2>>dbuf3>>dbufcos;
186
187 // Check table sizes
188 if(dbuf1!=NS1 || dbuf2!=NS2 || dbuf3!=NS3 || dbufcos!=NCOS) {
189 cout<<"mismatched NS1= "<<dbuf1 <<" <--> "<<NS1<<endl;
190 cout<<" NS2= "<<dbuf2 <<" <--> "<<NS2<<endl;
191 cout<<" NS3= "<<dbuf3 <<" <--> "<<NS3<<endl;
192 cout<<" NCOS= "<<dbufcos<<" <--> "<<NCOS<<endl;
193 return;
194 }
195
196 double buf1,buf2,buf3,buf4,buf5,buf6;
197 f>>buf>>buf1>>buf2>>buf3>>buf4>>buf5>>buf6;
198
199 // Set ranges
200 if(sminA==0.0){
201 sminA=buf1;
202 smaxA=buf2;
203 sminB=buf3;
204 smaxB=buf4;
205 sminC=buf5;
206 smaxC=buf6;
207 }
208
209 // Check ranges
210 if(buf1!=sminA || buf2!=smaxA || buf3!=sminB || buf4!=smaxB || buf5!=sminC || buf6!=smaxC) {
211 cout<<"mismatched sminA= "<<buf1<<" <--> "<<sminA<<endl;
212 cout<<" smaxA= "<<buf2<<" <--> "<<smaxA<<endl;
213 cout<<" sminB= "<<buf3<<" <--> "<<sminB<<endl;
214 cout<<" smaxB= "<<buf4<<" <--> "<<smaxB<<endl;
215 cout<<" sminC= "<<buf5<<" <--> "<<sminC<<endl;
216 cout<<" smaxC= "<<buf6<<" <--> "<<smaxC<<endl;
217 return;
218 }
219
220 // Print out header
221 while(!f.eof()){
222 char head[255];
223 f.getline(head,255);
224 if(strcmp(head,"BeginRange1")==0) break;
225 cout<<head<<endl;
226 }
227
228 // Read table
229 for (int i=0;i<NS1;i++){
230 for (int j=0;j<NCOS;j++){
231 for (int k=0;k<4;k++){
232 for (int l=0;l<4;l++){
233 f>>table1A[i][j][k][l];
234 } // for(l)
235 } // for(k)
236 f>>wtable1A[i][j];
237 f>>w0table1A[i][j];
238 } // for(j)
239 } // for(i)
240
241 // Find 2nd range
242 while(!f.eof()){
243 f>>buf;
244 if(strcmp(buf.c_str(),"BeginRange2")==0) break;
245 }
246
247 // Read table
248 for (int i=0;i<NS2;i++){
249 for (int j=0;j<NCOS;j++){
250 for (int k=0;k<4;k++){
251 for (int l=0;l<4;l++){
252 f>>table1B[i][j][k][l];
253 } // for(l)
254 } // for(k)
255 f>>wtable1B[i][j];
256 f>>w0table1B[i][j];
257 } // for(j)
258 } // for(i)
259
260 // Find 3rd range
261 while(!f.eof()){
262 f>>buf;
263 if(strcmp(buf.c_str(),"BeginRange3")==0) break;
264 }
265
266 // Read table
267 for (int i=0;i<NS3;i++){
268 for (int j=0;j<NCOS;j++){
269 for (int k=0;k<4;k++){
270 for (int l=0;l<4;l++){
271 f>>table1C[i][j][k][l];
272 } // for(l)
273 } // for(k)
274 f>>wtable1C[i][j];
275 f>>w0table1C[i][j];
276 } // for(j)
277 } // for(i)
278
279 // Check for proper file end
280 f>>buf;
281 if(buf.size() == 0 || strcmp(buf.c_str(),"End") != 0){
282 cout<<"...incorrect file version or file incomplete/damaged!"<<endl;
283
284 // In case of the error - do not use tables
285 table1A[0][0][0][0] = table1B[0][0][0][0] = table1C[0][0][0][0] = 0.0;
286 }
287 } // if (file is open)
288
289 f.close();
290 f.open("table2-2.txt");
291
292 if(!f.is_open()){
293 cout<<"File 'table2-2.txt' missing... skipped."<<endl;
294 }
295 else{
296 string buf;
297
298 cout<<"Reading file 'table2-2.txt'..."<<endl;
299
300 int dbuf1,dbuf2,dbuf3,dbufcos;
301 f>>buf>>dbuf1>>dbuf2>>dbuf3>>dbufcos;
302
303 // Check table sizes
304 if(dbuf1!=NS1 || dbuf2!=NS2 || dbuf3!=NS3 || dbufcos!=NCOS) {
305 cout<<"mismatched NS1= "<<dbuf1<<" <--> "<<NS1<<endl;
306 cout<<" NS2= "<<dbuf2<<" <--> "<<NS2<<endl;
307 cout<<" NS3= "<<dbuf3<<" <--> "<<NS3<<endl;
308 cout<<" NCOS= "<<dbufcos<<" <--> "<<NCOS<<endl;
309 return;
310 }
311
312 double buf1,buf2,buf3,buf4,buf5,buf6;
313 f>>buf>>buf1>>buf2>>buf3>>buf4>>buf5>>buf6;
314
315 // Set ranges
316 if(sminA==0.0)
317 {
318 sminA=buf1;
319 smaxA=buf2;
320 sminB=buf3;
321 smaxB=buf4;
322 sminC=buf5;
323 smaxC=buf6;
324 }
325
326 // Check ranges
327 if(buf1!=sminA || buf2!=smaxA || buf3!=sminB || buf4!=smaxB || buf5!=sminC || buf6!=smaxC) {
328 cout<<"mismatched sminA= "<<buf1<<" <--> "<<sminA<<endl;
329 cout<<" smaxA= "<<buf2<<" <--> "<<smaxA<<endl;
330 cout<<" sminB= "<<buf3<<" <--> "<<sminB<<endl;
331 cout<<" smaxB= "<<buf4<<" <--> "<<smaxB<<endl;
332 cout<<" sminC= "<<buf5<<" <--> "<<sminC<<endl;
333 cout<<" smaxC= "<<buf6<<" <--> "<<smaxC<<endl;
334 return;
335 }
336
337 // Print out header
338 while(!f.eof()){
339 char head[255];
340 f.getline(head,255);
341 if(strcmp(head,"BeginRange1")==0) break;
342 cout<<head<<endl;
343 }
344
345 // Read table
346 for (int i=0;i<NS1;i++){
347 for (int j=0;j<NCOS;j++){
348 for (int k=0;k<4;k++){
349 for (int l=0;l<4;l++){
350 f>>table2A[i][j][k][l];
351 } // for(l)
352 } // for(k)
353 f>>wtable2A[i][j];
354 f>>w0table2A[i][j];
355 } // for(j)
356 } // for(i)
357
358 // Find 2nd range
359 while(!f.eof()){
360 f>>buf;
361 if(strcmp(buf.c_str(),"BeginRange2")==0) break;
362 }
363
364 // Read table
365 for (int i=0;i<NS2;i++){
366 for (int j=0;j<NCOS;j++){
367 for (int k=0;k<4;k++){
368 for (int l=0;l<4;l++){
369 f>>table2B[i][j][k][l];
370 } // for(l)
371 } // for(k)
372 f>>wtable2B[i][j];
373 f>>w0table2B[i][j];
374 } // for(j)
375 } // for(i)
376
377 // Find 3rd range
378 while(!f.eof()){
379 f>>buf;
380 if(strcmp(buf.c_str(),"BeginRange3")==0) break;
381 }
382
383 // Read table
384 for (int i=0;i<NS3;i++){
385 for (int j=0;j<NCOS;j++){
386 for (int k=0;k<4;k++){
387 for (int l=0;l<4;l++){
388 f>>table2C[i][j][k][l];
389 } // for(l)
390 } // for(k)
391 f>>wtable2C[i][j];
392 f>>w0table2C[i][j];
393 } // for(j)
394 } // for(i)
395
396 // Check for proper file end
397 f>>buf;
398 if(buf.size()==0 || strcmp(buf.c_str(),"End")!=0){
399 cout<<"...incorrect file version or file incomplete/damaged!"<<endl;
400
401 // In case of the error - do not use tables
402 table2A[0][0][0][0] = table2B[0][0][0][0] = table2C[0][0][0][0] = 0.0;
403 }
404 } // if (file is open)
405
406 f.close();
407 f.open("table11-11.txt");
408
409 if(!f.is_open()){
410 cout<<"File 'table11-11.txt' missing... skipped."<<endl;
411 }
412 else{
413 string buf;
414
415 cout<<"Reading file 'table11-11.txt'..."<<endl;
416
417 int dbuf1,dbuf2,dbuf3,dbufcos;
418 f>>buf>>dbuf1>>dbuf2>>dbuf3>>dbufcos;
419
420 // Check table sizes
421 if(dbuf1!=NS1 || dbuf2!=NS2 || dbuf3!=NS3 || dbufcos!=NCOS) {
422 cout<<"mismatched NS1= "<<dbuf1<<" <--> "<<NS1<<endl;
423 cout<<" NS2= "<<dbuf2<<" <--> "<<NS2<<endl;
424 cout<<" NS3= "<<dbuf3<<" <--> "<<NS3<<endl;
425 cout<<" NCOS= "<<dbufcos<<" <--> "<<NCOS<<endl;
426 return;
427 }
428
429 double buf1,buf2,buf3,buf4,buf5,buf6;
430 f>>buf>>buf1>>buf2>>buf3>>buf4>>buf5>>buf6;
431
432 // Set ranges
433 if(sminA==0.0)
434 {
435 sminA=buf1;
436 smaxA=buf2;
437 sminB=buf3;
438 smaxB=buf4;
439 sminC=buf5;
440 smaxC=buf6;
441 }
442
443 // Check ranges
444 if(buf1!=sminA || buf2!=smaxA || buf3!=sminB || buf4!=smaxB || buf5!=sminC || buf6!=smaxC) {
445 cout<<"mismatched sminA= "<<buf1<<" <--> "<<sminA<<endl;
446 cout<<" smaxA= "<<buf2<<" <--> "<<smaxA<<endl;
447 cout<<" sminB= "<<buf3<<" <--> "<<sminB<<endl;
448 cout<<" smaxB= "<<buf4<<" <--> "<<smaxB<<endl;
449 cout<<" sminC= "<<buf5<<" <--> "<<sminC<<endl;
450 cout<<" smaxC= "<<buf6<<" <--> "<<smaxC<<endl;
451 return;
452 }
453
454 // Print out header
455 while(!f.eof()){
456 char head[255];
457 f.getline(head,255);
458 if(strcmp(head,"BeginRange1")==0) break;
459 cout<<head<<endl;
460 }
461
462 // Read table
463 for (int i=0;i<NS1;i++){
464 for (int j=0;j<NCOS;j++){
465 for (int k=0;k<4;k++){
466 for (int l=0;l<4;l++){
467 f>>table11A[i][j][k][l];
468 } // for(l)
469 } // for(k)
470 f>>wtable11A[i][j];
471 f>>w0table11A[i][j];
472 } // for(j)
473 } // for(i)
474
475 // Find 2nd range
476 while(!f.eof()){
477 f>>buf;
478 if(strcmp(buf.c_str(),"BeginRange2")==0) break;
479 }
480
481 // Read table
482 for (int i=0;i<NS2;i++){
483 for (int j=0;j<NCOS;j++){
484 for (int k=0;k<4;k++){
485 for (int l=0;l<4;l++){
486 f>>table11B[i][j][k][l];
487 } // for(l)
488 } // for(k)
489 f>>wtable11B[i][j];
490 f>>w0table11B[i][j];
491 } // for(j)
492 } // for(i)
493
494 // Find 3rd range
495 while(!f.eof()){
496 f>>buf;
497 if(strcmp(buf.c_str(),"BeginRange3")==0) break;
498 }
499
500 // Read table
501 for (int i=0;i<NS3;i++){
502 for (int j=0;j<NCOS;j++){
503 for (int k=0;k<4;k++){
504 for (int l=0;l<4;l++){
505 f>>table11C[i][j][k][l];
506 } // for(l)
507 } // for(k)
508 f>>wtable11C[i][j];
509 f>>w0table11C[i][j];
510 } // for(j)
511 } // for(i)
512
513 f>>buf;
514 if(buf.size()==0 || strcmp(buf.c_str(),"End")!=0){
515 cout<<"...incorrect file version or file incomplete/damaged!"<<endl;
516
517 // In case of the error - do not use tables
518 table11A[0][0][0][0] = table11B[0][0][0][0] = table11C[0][0][0][0] = 0.0;
519 }
520 } // if (file is open)
521
522 f.close();
523 cout<<endl;
524}
525
526void Tauola::decayOne(TauolaParticle *tau, bool undecay, double polx, double poly, double polz)
527{
528 if(!tau) return;
529
530 if(polx*polx+poly*poly+polz*polz>1)
531 {
532 Log::Warning()<<"decayOne(): ignoring wrong polarization vector: "<<polx<<" "<<poly<<" "<<polz<<endl;
533 polx=poly=polz=0;
534 }
535
536 // Let the interface know that we work in the decayOne mode
538
539 m_decay_one_polarization[0] = polx;
540 m_decay_one_polarization[1] = poly;
541 m_decay_one_polarization[2] = polz;
542
543 // Undecay if needed
544 if(tau->hasDaughters())
545 {
546 if(undecay) tau->undecay();
547 else
548 {
549 m_is_using_decay_one = false;
550 return;
551 }
552 }
553
554 std::vector<TauolaParticle *> list;
555 list.push_back(tau);
556
557 // Decay single tau
558 TauolaParticlePair t_pair(list);
559 t_pair.decayTauPair();
561
562 // Revert to normal mode
563 m_is_using_decay_one = false;
564}
565
567
568 Log::Warning() <<"Deprecated routine 'Tauola::initialise'"<<endl;
569 Log::Warning(0)<<"Use 'Tauola::initialize' instead."<<endl;
570
571 initialize();
572
573 // Deprecated routines: initialise, setInitialisePhy,
574 // f_interface_tauolaInitialise
575}
577{
578 return m_is_initialized;
579}
580
585
587{
588 return (bool) m_decay_one_boost_routine;
589}
590
595
597{
598 m_decay_one_boost_routine(mother,target);
599}
600
602{
604}
605
607 m_pdg_id=pdg_id;
608}
609
611 return abs(m_pdg_id);
612}
613
614void Tauola::setSameParticleDecayMode(int firstDecayMode){
615 m_firstDecayMode=firstDecayMode;
616}
617
618void Tauola::setOppositeParticleDecayMode(int secondDecayMode){
619 m_secondDecayMode=secondDecayMode;
620}
621
622void Tauola::setRadiation(bool rad){
623 m_rad=rad;
624}
625
626void Tauola::setRadiationCutOff(double rad_cut_off){
627 m_rad_cut_off=rad_cut_off;
628}
629
630
631void Tauola::setInitializePhy(double iniphy){
632 m_iniphy=iniphy;
633}
634
635void Tauola::setInitialisePhy(double iniphy){
636
637 Log::Warning() <<"Deprecated routine 'Tauola::setInitialisePhy'"<<endl;
638 Log::Warning(0)<<"Use 'Tauola::setInitializePhy' instead."<<endl;
639
640 setInitializePhy(iniphy);
641
642 // Deprecated routines: initialise, setInitialisePhy,
643 // f_interface_tauolaInitialise
644}
645
646void Tauola::setTauBr(int i, double value)
647{
648 if(taubra_.nchan==0)
649 Log::Warning()<<"setTauBr(): run Tauola::initialize() first."<<endl;
650 else if(i<1 || i>taubra_.nchan || value<0.)
651 Log::Warning()<<"setTauBr(): Invalid input. Value must be >= 0 and 0 < i <= "<<taubra_.nchan<<endl;
652 else taubra_.gamprt[i-1]=(float)value;
653}
654
655void Tauola::setTaukle(double bra1,double brk0, double brk0b, double brks)
656{
657 if(bra1<0 || bra1>1 || brk0<0 ||brk0>1 || brk0b<0 || brk0b>1 || brks<0 ||brks>1)
658 {
659 Log::Warning()<<"setTaukle(): variables must be in range [0,1]. Ignored."<<endl;
660 return;
661 }
662 taukle_.bra1 =(float)bra1;
663 taukle_.brk0 =(float)brk0;
664 taukle_.brk0b=(float)brk0b;
665 taukle_.brks =(float)brks;
666}
667
669 return f_getTauMass();
670}
671
672double Tauola::getHiggsScalarPseudoscalarMixingAngle(){
673 return m_higgs_scalar_pseudoscalar_mix;
674}
675
677 return m_higgs_scalar_pseudoscalar_pdg;
678}
679
680/** set the mixing angle. coupling: tau~(cos(phi)+isin(phi)gamma5)tau */
682 m_higgs_scalar_pseudoscalar_mix=angle;
683}
684
685/** set the pdg code of the higgs particle which tauola should
686 treat as a scalar-pseudoscalar mix */
688
689 if (particleCharge(pdg_code)!=0.0){
690 Log::Warning()<<"You want to use spin correlations of Higgs for particle of PDGID= "<<pdg_code<<endl
691 <<"This particle has charge="<<particleCharge(pdg_code)<<endl;
692 Log::Fatal("This choice is not appropriate.",0);
693 }
694 m_higgs_scalar_pseudoscalar_pdg=pdg_code;
695}
696
697int Tauola::getHelPlus(){
698 return m_helPlus;
699}
700int Tauola::getHelMinus(){
701 return m_helMinus;
702}
703double Tauola::getEWwt(){
704 return m_wtEW;
705}
706double Tauola::getEWwt0(){
707 return m_wtEW0;
708}
709void Tauola::setEWwt(double wt, double wt0)
710{
711 m_wtEW = wt;
712 m_wtEW0 = wt0;
713}
714void Tauola::setHelicities(int minus, int plus)
715{
716 m_helMinus = minus;
717 m_helPlus = plus;
718}
719void Tauola::setEtaK0sPi(int eta, int k, int pi)
720{
721 ion[0] = pi;
722 ion[1] = k;
723 ion[2] = eta;
724}
725
726void Tauola::summary()
727{
728 int sign_type=100;
729 double pol[4] = {0};
730
731 Log::Info() <<"Tauola::summary(): We use old TAUOLA FORTRAN printout."<<endl;
732 Log::Info(false)<<"As a consequence, there is a mismatch in printed TAUOLA version number."<<endl<<endl;
733
734 // Print summary taken from FORTRAN TAUOLA
735 dekay_(&sign_type,pol);
736}
737
738void Tauola::fill_val(int beg, int end, double* array, double value)
739{
740 for (int i = beg; i < end; i++)
741 array[i] = value;
742}
743
744double Tauola::particleCharge(int idhep)
745{
746 static double CHARGE[101] = { 0 };
747 static int j=0;
748 //--
749 //-- Array 'SPIN' contains the spin of the first 100 particles accor-
750 //-- ding to the PDG particle code...
751
752 if(j==0) // initialization
753 {
754 j=1;
755 fill_val(0 , 1, CHARGE, 0.0 );
756 fill_val(1 , 2, CHARGE,-0.3333333333);
757 fill_val(2 , 3, CHARGE, 0.6666666667);
758 fill_val(3 , 4, CHARGE,-0.3333333333);
759 fill_val(4 , 5, CHARGE, 0.6666666667);
760 fill_val(5 , 6, CHARGE,-0.3333333333);
761 fill_val(6 , 7, CHARGE, 0.6666666667);
762 fill_val(7 , 8, CHARGE,-0.3333333333);
763 fill_val(8 , 9, CHARGE, 0.6666666667);
764 fill_val(9 , 11, CHARGE, 0.0 );
765 fill_val(11 ,12, CHARGE,-1.0 );
766 fill_val(12 ,13, CHARGE, 0.0 );
767 fill_val(13 ,14, CHARGE,-1.0 );
768 fill_val(14, 15, CHARGE, 0.0 );
769 fill_val(15 ,16, CHARGE,-1.0 );
770 fill_val(16, 17, CHARGE, 0.0 );
771 fill_val(17 ,18, CHARGE,-1.0 );
772 fill_val(18, 24, CHARGE, 0.0 );
773 fill_val(24, 25, CHARGE, 1.0 );
774 fill_val(25, 37, CHARGE, 0.0 );
775 fill_val(37, 38, CHARGE, 1.0 );
776 fill_val(38,101, CHARGE, 0.0 );
777 }
778
779 int idabs=abs(idhep);
780 double phoch=0.0;
781
782 //--
783 //-- Charge of quark, lepton, boson etc....
784 if (idabs<=100) phoch=CHARGE[idabs];
785 else {
786 int Q3= idabs/1000 % 10;
787 int Q2= idabs/100 % 10;
788 int Q1= idabs/10 % 10;
789 if (Q3==0){
790 //--
791 //-- ...meson...
792 if(Q2 % 2==0) phoch=CHARGE[Q2]-CHARGE[Q1];
793 else phoch=CHARGE[Q1]-CHARGE[Q2];
794 }
795 else{
796 //--
797 //-- ...diquarks or baryon.
798 phoch=CHARGE[Q1]+CHARGE[Q2]+CHARGE[Q3];
799 }
800 }
801 //--
802 //-- Find the sign of the charge...
803 if (idhep<0.0) phoch=-phoch;
804 if (phoch*phoch<0.000001) phoch=0.0;
805
806 return phoch;
807}
808
809} // namespace Tauolapp
Abstract base class for particle in the event. This class also handles boosting.
static void Fatal(string text, unsigned short int code=0)
static void setSameParticleDecayMode(int firstDecayMode)
Definition Tauola.cxx:614
static void setTauBr(int i, double value)
Definition Tauola.cxx:646
static void setUnits(MomentumUnits m, LengthUnits l)
Definition Tauola.cxx:120
static void setInitializePhy(double iniphy)
Definition Tauola.cxx:631
static void initialise()
Definition Tauola.cxx:566
static void initialize()
Definition Tauola.cxx:129
static double defaultRandomGenerator()
Definition Tauola.cxx:92
static void setInitialisePhy(double iniphy)
Definition Tauola.cxx:635
static void setDecayingParticle(int pdg_id)
Definition Tauola.cxx:606
static double getTauMass()
Definition Tauola.cxx:668
static void setRandomGenerator(double(*gen)())
Definition Tauola.cxx:96
static void fill_val(int beg, int end, double *array, double value)
Definition Tauola.cxx:738
static void setBoostRoutine(void(*boost)(TauolaParticle *, TauolaParticle *))
Definition Tauola.cxx:591
static int getHiggsScalarPseudoscalarPDG()
Definition Tauola.cxx:676
static bool isUsingDecayOneBoost()
Definition Tauola.cxx:586
static void setNewCurrents(int mode)
Definition Tauola.cxx:87
static double m_decay_one_polarization[3]
static void setTauLifetime(double t)
Definition Tauola.cxx:125
static int getDecayingParticle()
Definition Tauola.cxx:610
static bool getIsTauolaIni()
Definition Tauola.cxx:576
static void decayOne(TauolaParticle *tau, bool undecay=false, double polx=0, double poly=0, double polz=0)
Definition Tauola.cxx:526
static bool isUsingDecayOne()
Definition Tauola.cxx:581
static void setHiggsScalarPseudoscalarMixingAngle(double angle)
Definition Tauola.cxx:681
static void setHiggsScalarPseudoscalarPDG(int pdg_id)
Definition Tauola.cxx:687
static const double * getDecayOnePolarization()
Definition Tauola.cxx:601
static void decayOneBoost(TauolaParticle *mother, TauolaParticle *target)
Definition Tauola.cxx:596
static void setOppositeParticleDecayMode(int secondDecayMode)
Definition Tauola.cxx:618
static void setRadiationCutOff(double rad_cut_off)
Definition Tauola.cxx:626
static void(* m_decay_one_boost_routine)(TauolaParticle *, TauolaParticle *)
static void setRadiation(bool rad)
Definition Tauola.cxx:622
static double particleCharge(int idhep)
Definition Tauola.cxx:744