Bug in Pythia Interface for semileptonic decays

[Bijan Chokoufe Nejad]

Here is a reduced test case to reproduce the problem that Mikael reported.

process ww_sl0muq = e1, E1 => u:c:U:C, d:s:b:D:S:B,e2:E2,n1:n2:n3:N1:N2:N3

iterations = 1:200
seed = 814995

! OPAL tune:
$ps_PYTHIA_PYGIVE = "MSTJ(28)=0; PMAS(25,1)=120.; PMAS(25,2)=0.3605E-02; MSTJ(41)=2; MSTU(22)=2000; PARJ(21)=0.40000; PARJ(41)=0.11000; PARJ(42)=0.52000; PARJ(81)=0.25000; PARJ(82)=1.90000; MSTJ(11)=3; PARJ(54)=-0.03100; PARJ(55)=-0.00200; PARJ(1)=0.08500; PARJ(3)=0.45000; PARJ(4)=0.02500; PARJ(2)=0.31000; PARJ(11)=0.60000; PARJ(12)=0.40000; PARJ(13)=0.72000; PARJ(14)=0.43000; PARJ(15)=0.08000; PARJ(16)=0.08000; PARJ(17)=0.17000;   MSTP(125)=2; MSTP(3)=1;MSTP(71)=1"
?ps_fsr_active = true
$shower_method = "PYTHIA6"
?ps_taudec_active = true
n_events = 1
sqrts = 500 GeV

simulate ( ww_sl0muq )

Run it with whizard *sin --debug shower and you directly see that the sum of outgoing momenta is not balanced. Interestingly, if I disable the $ps_PYTHIA_PYGIVE command, this does not happen in 100 events.

[Juergen Reuter]

There is some more information by Mikael Berggren coming in 2017-05-24:

I made a possibly helpful observation about the problem with quarks flagged as stable in the output in qqll events. The thing is
that this problem only shows up if one sets MSTP(125)=2 in Pythia (Pythia manual: "= 2 : list complete documentation of intermediate
steps of parton-shower evolution.") rather than keeping it at it's default value (=1). To study the performance of jet-clustering,
it is very useful to have the full history, so we definitely want this setting to work.

In addition, with  MSTP(125)=2, it is more up-front what the difference is in colour-connection in Pythia between Whizard 2.x.x and
Whizard 1.95+ Tim's user routine: In 2.x.x, one sees that all four quarks are joined into a single "CMshower" object (Pythia-code
94), while with Tim's interface, the two colour-connected pairs goes into two different CMshowers. Also, in qqll events, one sees
that all for fermions goes into a single CMshower in 2.x.x, while only the quarks do that with Tim's code. All this is seen if
MSTP(125)=2; if it is kept at the default, the CMshowers are deleted from the record.

I've enclosed a file with listings of two cscs events, one from the DBD production, one from whizard-2.5.0-head-rev7960, to
illustrate this. Apart from the CMshower, note that Tim makes a copy of the whizard event (isr-gammas + the four fermions, where
only the order is changed to be (q-qbar)[colour=501] (q-qpar)[colour=502], while nothing like that is done in the Whizard-2.5.0 case
(actually the order is qbar-qbar-q-q, and the adjacent entries do not have the same colour-codes).

[Bijan Chokoufe Nejad]

File that Mikael attached

[Juergen Reuter]

This most likely be automagically resolved by #814.

[Bijan Chokoufe Nejad]

#814 is about to be merged and the other issue was the usage of MSTP(125)=2, which if I understood Mikael correctly, they will stop using in the future and rely on the history instead.