Syntax to calculate NLO EW corrections

Dear experts,

I am trying to calculate the NLO electroweak (aS^2*aEW^1) corrections
to top pair production at the LHC.

As far as I know, the syntax reads

   generate LO_process QCD_QED_order_of_LO_matrix_element [QED]

So naively I expect the following five statements are equivalent,
because, without the [QED] tag, they generate the same pure QCD LO
contribution (aS^2 * aEW^0).

# process1
> generate p p > t t~ [QED]
# process2
> generate p p > t t~ QED^2=0 QCD^2=4 [QED]
# process3
> generate p p > t t~ QED^2==0 QCD^2=4 [QED]
# process4
> generate p p > t t~ QED^2==0 [QED]
# process5
> generate p p > t t~ QED^2==0 QCD^2==4 [QED]

Processes 1 and 2 give the same result 4.535e+02 +- 1.0e+00 pb

Process 3 gives 6.374e+01 +- 7.5e-02 pb, which is much smaller than
the LO prediction, so it cannot be what I want.

Processes 4 and 5 are interrupted by the error:
  No Born diagrams found.
  NoDiagramException : No amplitudes generated

Here are my questions.
(1) What is the difference between processes 1~5?
    Does process 3 give any physical predictions?

(2) Do processes 1 and 2 calculate the NLO EW corrections (aS^2*aEW^1)
to the top pair production?

(3) How to calculate full LO + NLO cross sections?
LO = aS^2 * aEW^0 + aS^1 * aEW^1 + aS^1 * aEW^2
NLO = aS^3 * aEW^0 + ... + aS^0 * aEW^3

/******************************************************************/

As an exercise, I have also tried to reproduce the NLO EW corrections
(LO_1 + NLO_2) of ttZ production in table 2 of the NLO EW automation paper
1804.10017, where
  LO_1 ~ as^2 * aEW^1 (without subleading contributions ~ as^1*aEW^2 etc.)
  NLO_2 ~ as^2 * aEW^2

Here is the script I used
#######################################################
set complex_mass_scheme true
import model loop_qcd_qed_sm_Gmu
define p = g d d~ u u~ s s~ c c~ b b~ a
define j = p
generate p p > t t~ z QCD^2=4 QED^2=2 [QED]

output ttzNLOEW
launch
fixed_order = ON

set ebeam1 6500.0
set ebeam2 6500.0

set pdlabel lhapdf
set lhaid 82000 # LUXqed plus PDF4LHC15 nnlo 100

set fixed_ren_scale False
set fixed_fac_scale False
set dynamical_scale_choice 3 # half of the sum of transverse masses
# no kinematic cuts added
#######################################################

which gives
   LO = 0.5063 +- 0.001458 pb
   NLOEW = 0.5150 +- 5.1e-04 pb

and fails to reproduce the results in the NLO EW automation paper 1804.10017
   LO=0.50456,
   NLOEW=0.50033 (pb).

The LO predictions are consistent up to the Monte-Carlo error,
but the sign of the NLO EW corrections is wrong.

(4) Do you have any idea what is wrong with the script above?

Thanks.

Xiaomin