# $Id: LEP_higgs.sin 2293 2010-04-11 23:57:50Z jr_reuter $
# Simple complete physics example: Higgs search at LEP
########################################################################
#
# Copyright (C) 1999-2015 by 
#     Wolfgang Kilian <kilian@physik.uni-siegen.de>
#     Thorsten Ohl <ohl@physik.uni-wuerzburg.de>
#     Juergen Reuter <juergen.reuter@desy.de>
#     with contributions from
#     Christian Speckner <cnspeckn@googlemail.com>
#
# WHIZARD is free software; you can redistribute it and/or modify it
# under the terms of the GNU General Public License as published by 
# the Free Software Foundation; either version 2, or (at your option)
# any later version.
#
# WHIZARD is distributed in the hope that it will be useful, but
# WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the 
# GNU General Public License for more details.
#
# You should have received a copy of the GNU General Public License
# along with this program; if not, write to the Free Software
# Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
#
########################################################################

model = SM_CKM

# Alias definition
# include ( "%%ALIASES%%" )

# ww_h0cuxx ww_h0uubd ww_h0uusd ww_h0ccbs ww_h0ccds 
# #
#       WW: more than 2 flavours
#
# ww_h0cuxx       e1,E1 cq,down_type_q,uq,down_type_q                       omega w:c,c
# ww_h0uubd       e1,E1 uq,bq,uq,dq                                         omega w:c,c
# ww_h0uusd       e1,E1 uq,sq,uq,dq                                         omega w:c,c
# ww_h0ccbs       e1,E1 cq,bq,cq,sq                                         omega w:c,c
# ww_h0ccds       e1,E1 cq,dq,cq,sq                                         omega w:c,c
#
#yyuyyu_o	e1,E1	y,Y,u,Y,y,U                     omega	w:c,c
#yyuyyc_o	e1,E1	y,Y,u,Y,y,C                     omega	w:c,c
#yycyyu_o	e1,E1	y,Y,c,Y,y,U                     omega	w:c,c
#yycyyc_o	e1,E1	y,Y,c,Y,y,C                     omega	w:c,c
alias x = u:c
alias y = d:s:b
alias X = U:C
alias Y = D:S:B


   alias uq = u:U
   alias dq = d:D
   alias cq = c:C
   alias sq = s:S
   alias bq = b:B
#
   alias up_type_q = u:c:U:C
   alias down_type_q = d:s:b:D:S:B
   alias q = u:d:s:c:b
   alias Q = U:D:S:C:B

   alias pelectron = e1:E1
   alias pmu = e2:E2
   alias ptau = e3:E3
#
   alias pnu_e = n1:N1
   alias pnu_mu = n2:N2
   alias pnu_tau = n3:N3
#
   alias quark = u:d:s:c:b:U:D:S:C:B
 
#
#   alias up_type_q = u:c:U:C
#   alias down_type_q = d:s:b:D:S:B
#
   alias lepton = e1:e2:e3:E1:E2:E3:n1:n2:n3:N1:N2:N3
   alias clepton = e1:e2:e3:E1:E2:E3
   alias not_electron = e2:e3:E2:E3
#
   alias neutrino = n1:n2:n3:N1:N2:N3   
   alias not_nu_e = n2:n3:N2:N3
   alias n = n1:n2:n3
   alias N = N1:N2:N3   

   alias f = u:d:s:c:b:e1:e2:e3:n1:n2:n3
   alias F = U:D:S:C:B:E1:E2:E3:N1:N2:N3
#me = 0 GeV
#mmu = 0 GeV
#mmu = 0.10566 GeV
#   mu = 0 GeV
#   md = 0 GeV
  ms = 0 GeV
  mc = 0 GeV
  mb = 0 GeV 
  alphas = 0.000001
#
#      e+e- --> mu+mu- and tau+tau-
#
polarized e1,E1,e2,E2,e3,E3,n1,N1,n2,N2,n3,N3,u,U,d,D,s,S,c,C,b,B,t,T, "W+", "W-", Z

#process_num_id = 401010
process P6f_llxyyx = e1, E1 => e2+e3, E2+E3, x, Y, y, X { process_num_id = 401010 }
#process P6f_llxyyx = e1, E1 => (e2, E2, x, Y, y, X) + (e3, E3, x, Y, y, X) { process_num_id = 401010 }

compile

#
#logical ?is_higgs_process = true
logical ?is_higgs_process = false

#include ("%%PYTHIA6_PARAMETERS%%" )
include ( "/home/ilc/tianjp/generator/PostDBD/whizard2/sindarin/pythia6-parameters.sin" )

sqrts = 250 GeV

#beams = e1, E1 => beam_events => isr
beams = e1, E1 => circe2 => isr,isr
#beams = e1, E1 => isr
#beams = e1, E1 
#beams = e1, E1 => beam_events
#$beam_events_file = "%%BEAM_EVENTS_FILE%%"
#$beam_events_file ="/group/ilc/soft/gcc481/GP/beamdata/lumiee500-hs/lumiee500-run1to20.dat"
#$circe2_file = "/afs/desy.de/user/b/berggren/public/energy_spread/250_SetA_ee024.circe.gz"
$circe2_file = "/home/ilc/tianjp/generator/PostDBD/whizard2/250_SetA_ee024.circe"
beams_pol_density = @(+1), @(-1)
?circe2_polarized = false
beams_pol_fraction = 1.0,1.0
?isr_recoil = true

isr_alpha = 0.0072993
isr_mass = 0.000511 GeV

real default_Mqq_cut = 10 GeV
real default_Mll_cut = 4 GeV
#real default_Mll_cut = 10 GeV
real default_Mla_cut = 4 GeV
#cuts = all M > default_Mqq_cut [ quark,  quark ] and
#       (all M > default_Mll_cut [ lepton, lepton ] and all M < -default_Mll_cut [ incoming lepton, lepton ]) and
#       (all M > default_Mla_cut [ lepton, photon ] and all M < -default_Mla_cut [ incoming lepton, photon ]) and
#       all E > 10 GeV [photon]
cuts =  all M > default_Mqq_cut [ quark,  quark ] and
        all M > default_Mll_cut [ clepton, clepton ] and 
        all M < -default_Mll_cut [ incoming clepton, clepton ]

iterations = 1:500000, 10:500000, 1:1500000

relative_error_goal = 0.01

$integration_method = "vamp2"
$rng_method = "rng_stream"

integrate ( P6f_llxyyx) { ?polarized_events = true }

show (results)


?keep_remnants = true
?keep_beams = true
?hadronization_active = true
sample_split_n_kbytes = 500000

sample_format = lcio
?write_raw = false
$sample = "E250-SetA.P6f_llxyyx.Gwhizard-2_8_2.eR.pL.I401010"

#luminosity = 10
n_events = 20000
#n_events = 50
#simulate ( P6f_llxyyx) { ?polarized_events = true }
