! File generated automatically by O'Mega
!
!   /home/pestix/physik/svn/whizard/trunk/build/../install/bin/omega_SM.opt -o test.f90 -target:whizard -target:parameter_module parameters_SM -target:module test -target:md5sum A0AF22C6335BA7C3F095CF9CD11C69A3 -fusion:progress -scatter "e+ e- -> e+ e-"
!
! with all scattering amplitudes for the process(es)
!
!   flavor combinations:
!
!       1: e+ e- -> e+ e-
!
!   color flows:
!
!       1: (  0,  0) (  0,  0) -> (  0,  0) (  0,  0)
!
!     NB: i.g. not all color flows contribute to all flavor
!     combinations.  Consult the array FLV_COL_IS_ALLOWED
!     below for the allowed combinations.
!
!   Color Factors:
!
!     (  1,  1): + 1
!
!   vanishing or redundant flavor combinations:
!
!
! in minimal electroweak standard model in unitarity gauge
!
module test

  use kinds
  use omega95
  use omega_color, OCF => omega_color_factor

  use parameters_SM

  implicit none
  private

  public :: md5sum
  public :: number_particles_in, number_particles_out
  public :: number_spin_states, spin_states
  public :: number_flavor_states, flavor_states
  public :: number_color_flows, color_flows
  public :: number_color_indices, number_color_factors, color_factors, &
    color_sum
  public :: init, final, update_alpha_s
  public :: reset_helicity_selection

  public :: new_event, is_allowed, get_amplitude

  ! DON'T EVEN THINK of removing the following!
  ! If the compiler complains about undeclared
  ! or undefined variables, you are compiling
  ! against an incompatible omega95 module!
  integer, dimension(7), parameter, private :: require = &
    (/ omega_spinors_2010_01_A, omega_spinor_cpls_2010_01_A, &
       omega_vectors_2010_01_A, omega_polarizations_2010_01_A, &
       omega_couplings_2010_01_A, omega_color_2010_01_A, &
       omega_utils_2010_01_A /)

  integer, parameter, private :: n_prt = 4
  integer, parameter, private :: n_in = 2
  integer, parameter, private :: n_out = 2
  integer, parameter, private :: n_cflow = 1
  integer, parameter, private :: n_cindex = 2
  integer, parameter, private :: n_flv = 1
  integer, parameter, private :: n_hel = 16

  ! NB: you MUST NOT change the value of N_ here!!!
  !     It is defined here for convenience only and must be
  !     compatible with hardcoded values in the amplitude!
  real(kind=default), parameter, private :: N_ = 3
  logical, parameter, private :: F = .false.
  logical, parameter, private :: T = .true.

  integer, dimension(n_prt), parameter, private :: s0001 = (/ -1, -1, -1, -1 /)
  integer, dimension(n_prt), parameter, private :: s0002 = (/ -1, -1, -1,  1 /)
  integer, dimension(n_prt), parameter, private :: s0003 = (/ -1, -1,  1, -1 /)
  integer, dimension(n_prt), parameter, private :: s0004 = (/ -1, -1,  1,  1 /)
  integer, dimension(n_prt), parameter, private :: s0005 = (/ -1,  1, -1, -1 /)
  integer, dimension(n_prt), parameter, private :: s0006 = (/ -1,  1, -1,  1 /)
  integer, dimension(n_prt), parameter, private :: s0007 = (/ -1,  1,  1, -1 /)
  integer, dimension(n_prt), parameter, private :: s0008 = (/ -1,  1,  1,  1 /)
  integer, dimension(n_prt), parameter, private :: s0009 = (/  1, -1, -1, -1 /)
  integer, dimension(n_prt), parameter, private :: s0010 = (/  1, -1, -1,  1 /)
  integer, dimension(n_prt), parameter, private :: s0011 = (/  1, -1,  1, -1 /)
  integer, dimension(n_prt), parameter, private :: s0012 = (/  1, -1,  1,  1 /)
  integer, dimension(n_prt), parameter, private :: s0013 = (/  1,  1, -1, -1 /)
  integer, dimension(n_prt), parameter, private :: s0014 = (/  1,  1, -1,  1 /)
  integer, dimension(n_prt), parameter, private :: s0015 = (/  1,  1,  1, -1 /)
  integer, dimension(n_prt), parameter, private :: s0016 = (/  1,  1,  1,  1 /)
  integer, dimension(n_prt,n_hel), parameter, private :: table_spin_states = &
    reshape ( (/ s0001, s0002, s0003, s0004, s0005, s0006, s0007, s0008, s0009, &
    s0010, s0011, s0012, s0013, s0014, s0015, s0016 /), (/ n_prt, n_hel /) )

  integer, dimension(n_prt), parameter, private :: &
    f0001 = (/ -11,  11, -11,  11 /) ! e+ e- e+ e-
  integer, dimension(n_prt,n_flv), parameter, private :: table_flavor_states = &
    reshape ( (/ f0001 /), (/ n_prt, n_flv /) )

  integer, dimension(n_cindex, n_prt), parameter, private :: &
    c0001 = reshape ( (/  0,0,  0,0,  0,0,  0,0 /), (/ n_cindex, n_prt /) )
  integer, dimension(n_cindex, n_prt, n_cflow), parameter, private :: &
    table_color_flows = reshape ( (/ c0001 /), (/ n_cindex, n_prt, n_cflow /) )

  logical, dimension(n_prt), parameter, private :: g0001 = (/  F,  F,  F,  F /)
  logical, dimension(n_prt, n_cflow), parameter, private :: table_ghost_flags = &
    reshape ( (/ g0001 /), (/ n_prt, n_cflow /) )

  integer, parameter, private :: n_cfactors = 1
  type(OCF), dimension(n_cfactors), parameter, private :: &
    table_color_factors = (/  OCF(1,1,+one) /)

  logical, dimension(n_flv), parameter, private :: a0001 = (/ T /)
  logical, dimension(n_flv, n_cflow), parameter, private :: &
    flv_col_is_allowed = reshape ( (/ a0001 /), (/ n_flv, n_cflow /) )

  complex(kind=default), dimension(n_flv, n_hel, n_cflow), private, save :: amp

  logical, dimension(n_hel), private, save :: hel_is_allowed = T
  real(kind=default), dimension(n_hel), private, save :: hel_max_abs = 0
  real(kind=default), private, save :: hel_sum_abs = 0, hel_threshold = 1E10
  integer, private, save :: hel_count = 0, hel_cutoff = 100

contains

  pure function md5sum ()
    character(len=32) :: md5sum
    ! DON'T EVEN THINK of modifying the following line!
    md5sum = "A0AF22C6335BA7C3F095CF9CD11C69A3"
  end function md5sum

  subroutine init (par)
    real(default), dimension(*), intent(in) :: par
    call import_from_whizard (par)
  end subroutine init

  subroutine final ()
    use diagnostics, only: logfile_unit
    integer :: unit
    if (hel_threshold .gt. 0) then
      unit = logfile_unit (logfile=.true.)
      if (unit .ge. 0) then
        call omega_report_helicity_selection (hel_is_allowed, &
                table_spin_states, hel_threshold, unit)
      end if
    end if
  end subroutine final

  subroutine update_alpha_s (alpha_s)
    real(default), intent(in) :: alpha_s
    call model_update_alpha_s (alpha_s)
  end subroutine update_alpha_s

  pure function number_particles_in () result (n)
    integer :: n
    n = n_in
  end function number_particles_in

  pure function number_particles_out () result (n)
    integer :: n
    n = n_out
  end function number_particles_out

  pure function number_spin_states () result (n)
    integer :: n
    n = size (table_spin_states, dim=2)
  end function number_spin_states

  pure subroutine spin_states (a)
    integer, dimension(:,:), intent(out) :: a
    a = table_spin_states
  end subroutine spin_states

  pure function number_flavor_states () result (n)
    integer :: n
    n = size (table_flavor_states, dim=2)
  end function number_flavor_states

  pure subroutine flavor_states (a)
    integer, dimension(:,:), intent(out) :: a
    a = table_flavor_states
  end subroutine flavor_states

  pure function number_color_indices () result (n)
    integer :: n
    n = size (table_color_flows, dim=1)
  end function number_color_indices

  pure function number_color_flows () result (n)
    integer :: n
    n = size (table_color_flows, dim=3)
  end function number_color_flows

  pure subroutine color_flows (a, g)
    integer, dimension(:,:,:), intent(out) :: a
    logical, dimension(:,:), intent(out) :: g
    a = table_color_flows
    g = table_ghost_flags
  end subroutine color_flows

  pure function number_color_factors () result (n)
    integer :: n
    n = size (table_color_factors)
  end function number_color_factors

  pure subroutine color_factors (cf)
    type(OCF), dimension(:), intent(out) :: cf
    cf = table_color_factors
  end subroutine color_factors

  pure function color_sum (flv, hel) result (amp2)
    integer, intent(in) :: flv, hel
    real(kind=default) :: amp2
    amp2 = real (omega_color_sum (flv, hel, amp, table_color_factors))
  end function color_sum

  subroutine new_event (p)
    real(kind=default), dimension(0:3,*), intent(in) :: p
    call calculate_amplitudes (amp, p, hel_is_allowed)
    if ((hel_threshold .gt. 0) .and. (hel_count .le. hel_cutoff)) then
      call omega_update_helicity_selection (hel_count, amp, hel_max_abs, &
              hel_sum_abs, hel_is_allowed, hel_threshold, hel_cutoff)
    end if
  end subroutine new_event

  subroutine reset_helicity_selection (threshold, cutoff)
    real(kind=default), intent(in) :: threshold
    integer, intent(in) :: cutoff
    hel_is_allowed = T
    hel_max_abs = 0
    hel_sum_abs = 0
    hel_count = 0
    hel_threshold = threshold
    hel_cutoff = cutoff
  end subroutine reset_helicity_selection

  pure function is_allowed (flv, hel, col) result (yorn)
    logical :: yorn
    integer, intent(in) :: flv, hel, col
    yorn = hel_is_allowed(hel) .and. flv_col_is_allowed(flv,col)
  end function is_allowed

  pure function get_amplitude (flv, hel, col) result (amp_result)
    complex(kind=default) :: amp_result
    integer, intent(in) :: flv, hel, col
    amp_result = amp(flv, hel, col)
  end function get_amplitude

  pure subroutine calculate_amplitudes (amp, k, mask)
    complex(kind=default), dimension(:,:,:), intent(out) :: amp
    real(kind=default), dimension(0:3,*), intent(in) :: k
    logical, dimension(:), intent(in) :: mask
    integer, dimension(n_prt) :: s
    integer :: h
    type(momentum) :: p1, p2, p3, p4
    type(spinor) :: l1_3, l1_2
    type(conjspinor) :: l1b_4, l1b_1
    type(momentum) :: p12, p13
    type(vector) :: a_13, a_12, z_13, z_12
    complex(kind=default) :: l1bl1l1bl1
    p1 = - k(:,1) ! incoming
    p2 = - k(:,2) ! incoming
    p3 =   k(:,3) ! outgoing
    p4 =   k(:,4) ! outgoing
    p12 = p1 + p2
    p13 = p1 + p3
    amp = 1
  end subroutine calculate_amplitudes

end module test
! O'Mega revision control information:
!    Modellib.SM:
!      minimal electroweak standard model in unitarity gauge
!      Source: /trunk/src/omega/src/modellib_SM.ml 
!      revision: 2219  checked in by ohl  at 2010-04-04 18:05:44 +0200 (So, 04. Apr 2010) 
!    Targets.Make_Fortran():
!      Interface for Whizard 2.X
!      NB: non-gauge vector couplings are not available yet
!      Source: /trunk/src/omega/src/targets.ml 
!      revision: 2219  checked in by ohl  at 2010-04-04 18:05:44 +0200 (So, 04. Apr 2010) 
!    Targets.Fortran_Fermions():
!      generates Fortran95 code for Dirac fermions
!      using revision 2000_10_A of module omega95
!      Source: /trunk/src/omega/src/targets.ml 
!      revision: 2219  checked in by ohl  at 2010-04-04 18:05:44 +0200 (So, 04. Apr 2010) 
!    DAG.Graded():
!      Graded directed Acyclical Graph 
!      representing binary or n-ary trees
!      Source: /trunk/src/omega/src/dAG.ml 
!      revision: 2219  checked in by ohl  at 2010-04-04 18:05:44 +0200 (So, 04. Apr 2010) 
!    Topology.Mixed23:
!      phi**3 + phi**4 topology
!      Source: /trunk/src/omega/src/topology.ml 
!      revision: 759  checked in by ohl  at 2009-06-10 11:38:07 +0200 (Mi, 10. Jun 2009) 
!    Momentum.Bits():
!      Finite disjoint sums of momenta
!      using bitfields as representation.
!      Source: /trunk/src/omega/src/momentum.ml 
!      revision: 2000  checked in by ohl  at 2010-03-05 14:57:05 +0100 (Fr, 05. Mär 2010) 
!    Colorize.It():
!      Colorizing Generic Monochrome Models
!      Source: /trunk/src/omega/src/colorize.ml 
!      revision: 2219  checked in by ohl  at 2010-04-04 18:05:44 +0200 (So, 04. Apr 2010) 
!    Fusion.Make():
!      Fusions for arbitrary topologies
!      Source: /trunk/src/omega/src/fusion.ml 
!      revision: 2219  checked in by ohl  at 2010-04-04 18:05:44 +0200 (So, 04. Apr 2010)