Merge lp:~mg5hpcteam/mg5amcnlo/mg5hpc into lp:~mg5core2/mg5amcnlo/2.5.5

[Olivier Mattelaer]

Modification to the core of MG5 to increase modularity of the code in a way which allows
the MPI plugin to work.

[Valentin Hirschi (valentin-hirschi)]

Hi Olivier,

There are quite a few things I don't understand well. So could you please elaborate on those in reply to this message or directly on skype if you think this would be more convenient?

a) In 'madgraph/interface/common_run_interface.py':

415 - if amcatnlo:
416 - run_card = banner_mod.RunCardNLO(opt['run_card'])
417 - else:
418 - run_card = banner_mod.RunCard(opt['run_card'])
419 + run_card = banner_mod.RunCard(opt['run_card'])

How is this giving you the same behavior as before given that I don't see banner_mod.RunCard modifications that would turn it into a factory. Or was that feature already there before?

b)

+ self.run_card = run_card

What is the reason for adding the attribute self.run_card to the inteface?

c)
I dont get this:

+ if self.configured > time_mod and hasattr(self, 'random') and hasattr(self, 'run_card'):
553 + #just ensure that cluster specific are correctly handled
554 + self.cluster.modify_interface(self)

what are the cases where this has an effect, and why discriminating against the random and run_card attribute?

d) this addition here:

+ IF(CONFSUB(I,SYMCONF(J)).NE.0) THEN

in super_auto_dsig is probably for you to better control what processes are selected. Did you check that it couldn't have any side effect in the normal workflow of madevent?

e) I see:

         self.add_param('survey_splitting', -1, hidden=True, include=False)
841 + self.add_param('survey_nchannel_per_job', 2, hidden=True, include=False)

Isn't it a bit weird to have a new banner parameter that is only relevant in the presence of a plugin? Also, I imagine that having both active at the same time would crash the code?

f) I don't understand this:

+ def modify_interface(self, run_interface):
856 + """routine which allow to modify the run_card/mg5cmd object to change the
857 + default behavior of the runs.
858 + This is called at the time of the compilation of the run_card.
859 + Note that this function can be called multiple times by run.
860 + """
861 + #run_card = run_interface.run_card
862 + return

and especially why it is in the cluster.py file?

g) I am a bit against this change:

- logger.error("Plugin %s is not supported by this version of MG5aMC." % name)
876 - plugin_support[name] = False
877 + if __debug__:
878 + logger.error("Plugin %s seems not supported by this version of MG5aMC. Keep it active (please update status)" % name)
879 + plugin_support[name] = True
880 + else:
881 + logger.error("Plugin %s is not supported by this version of MG5aMC." % name)
882 + plugin_support[name] = False

I like the idea, but please change the message in non-debug mode so as to let the user know that he can still try to run the plugin in the '--debug' mode, at his own risk of course.

review: Needs Information

[Olivier Mattelaer (olivier-mattelaer)]

Hi Valentin,

Thanks a lot for this.

> How is this giving you the same behavior as before given that I don't see banner_mod.RunCard modifications that would turn it into a factory. Or was that feature already there before?

Yes the RunCard class is and has always been a factory class.
The funny part is since this is the motherless of RunCardLO and RunCardNLO they are both factory class as well.
So in principle if you call RunCardLO(path) the instance created can be RunCardNLO. Not sure it really work though.

> b)
>
> + self.run_card = run_card
>
> What is the reason for adding the attribute self.run_card to the inteface?

This is always useful to keep such information directly available rather than recompute it systematically.

> c)
> I dont get this:
>
> + if self.configured > time_mod and hasattr(self, 'random') and hasattr(self, 'run_card'):
> 553 + #just ensure that cluster specific are correctly handled
> 554 + self.cluster.modify_interface(self)

This is more security than anything else.

> d) this addition here:
>
> + IF(CONFSUB(I,SYMCONF(J)).NE.0) THEN
>
> in super_auto_dsig is probably for you to better control what processes are selected. Did you check that it couldn't have any side effect in the normal workflow of madevent?

This line is not needed in madevent since those will never be selected. In MPI, since I run multiple integral in a row, it is better/required to double check that those are not selected.

> e) I see:
>
> self.add_param('survey_splitting', -1, hidden=True, include=False)
> 841 + self.add_param('survey_nchannel_per_job', 2, hidden=True, include=False)
>
> Isn't it a bit weird to have a new banner parameter that is only relevant in the presence of a plugin? Also, I imagine that having both active at the same time would crash the code?

This parameter has some impact for normal survey. It was just hardcoded before and I allow like this more flexibility.
Having both active, will not make the code to crash since they depend of different job strategy.

> f) I don't understand this:
>
> + def modify_interface(self, run_interface):
> 856 + """routine which allow to modify the run_card/mg5cmd object to change the
> 857 + default behavior of the runs.
> 858 + This is called at the time of the compilation of the run_card.
> 859 + Note that this function can be called multiple times by run.
> 860 + """
> 861 + #run_card = run_interface.run_card
> 862 + return
>
> and especially why it is in the cluster.py file?

The idea is that a plugin might want to change the default job strategy and therefore influence the way the job are created.
(like setting the survey_nchannel_per_job higher or lower or survey_splitting higher or lower)

In MPI, I'm using it for the following:
        run_card.set('survey_splitting',0, changeifuserset=False)
        run_card.set('survey_nchannel_per_job', 50, changeifuserset=False)
        run_card.set('gridpack', True, changeifuserset=True)

        if run_card['nhel']:
            logger.warning('MPI mode does not support MC over helicity')
            run_card.se...

[Valentin Hirschi (valentin-hirschi)]

Hi Olivier,

On Fri, May 26, 2017 at 10:57 AM, Olivier Mattelaer <
<email address hidden>> wrote:

> Hi Valentin,
>
> Thanks a lot for this.
>
> > How is this giving you the same behavior as before given that I don't
> see banner_mod.RunCard modifications that would turn it into a factory. Or
> was that feature already there before?
>
> Yes the RunCard class is and has always been a factory class.
> The funny part is since this is the motherless of RunCardLO and RunCardNLO
> they are both factory class as well.
> So in principle if you call RunCardLO(path) the instance created can be
> RunCardNLO. Not sure it really work though.

Well, then it's clearly much better as you put it now, so perfect.

>
> > b)
> >
> > + self.run_card = run_card
> >
> > What is the reason for adding the attribute self.run_card to the
> inteface?
>
>
> This is always useful to keep such information directly available rather
> than recompute it systematically.
>

Well I'm just afraid of side-effects of such "optimization" that can happen
when the actual run_card and the attribute become de-synchronized.
But I assume that you take care of updating it whenever the code asks to
edit the card.
Also, I assume you make sure that this attribute is assigned whenever
loading the interface on an existing directory.
If so, then fine by me.

> > c)
> > I dont get this:
> >
> > + if self.configured > time_mod and hasattr(self, 'random') and
> hasattr(self, 'run_card'):
> > 553 + #just ensure that cluster specific are correctly
> handled
> > 554 + self.cluster.modify_interface(self)
>
> This is more security than anything else.

You mean that you only want this to happen for certain "type" of interface.
Then shouldn't you define a new class attribute that specifies if it
qualifies for such behaviors?
That maybe an overkill here, but would anyway be useful for future
interface-type of plugins no?

>
> > d) this addition here:
> >
> > + IF(CONFSUB(I,SYMCONF(J)).NE.0) THEN
> >
> > in super_auto_dsig is probably for you to better control what processes
> are selected. Did you check that it couldn't have any side effect in the
> normal workflow of madevent?
>
> This line is not needed in madevent since those will never be selected. In
> MPI, since I run multiple integral in a row, it is better/required to
> double check that those are not selected.

You mean multiple integral *in parallel*, so that you needed this to avoid
concurrency issues? Ok, I'm not sure how these could arise, but if you
confirm that this line cannot possibly impact the normal madevent workflow,
then I'm fine with it.

>
> > e) I see:
> >
> > self.add_param('survey_splitting', -1, hidden=True,
> include=False)
> > 841 + self.add_param('survey_nchannel_per_job', 2,
> hidden=True, include=False)
> >
> > Isn't it a bit weird to have a new banner parameter that is only
> relevant in the presence of a plugin? Also, I imagine that having both
> active at the same time would crash the code?
>
> This parameter has some impact for normal survey. It was just hardcoded
> before and I allow like this more flexibility.
> Having both active, will ...

[Valentin Hirschi (valentin-hirschi)]

You answered all my queries and as a result nothing needs immediate modification. I therefore approve this merge.

review: Approve

[lp:~mg5hpcteam/mg5amcnlo/mg5hpc]

352. By olivier-mattelaer on 2017-07-13

merge with latest 2.5.6

=== modified file 'MadSpin/decay.py'
--- MadSpin/decay.py	2017-05-22 08:05:00 +0000
+++ MadSpin/decay.py	2017-07-13 10:23:41 +0000
@@ -207,6 +207,7 @@
         line+="\n"
         return line
 
+
     def reshuffle_resonances(self,mother):
         """ reset the momentum of each resonance in the production event
                 to the sum of momenta of the daughters
@@ -222,9 +223,9 @@
                 if self.resonance[index]["mothup1"]==mother:
                     daughters.append(index)
 
-        if len(daughters)!=2:
-            logger.info("Got more than 2 (%s) daughters for one particles" % len(daughters))
-            logger.info("in one production event (before decay)")
+#        if len(daughters)!=2:
+#            logger.info("Got more than 2 (%s) daughters for one particles" % len(daughters))
+#            logger.info("in one production event (before decay)")
 
 
         if daughters[0]>0:
@@ -1633,11 +1634,15 @@
                 self.ask_edit_cards(['param_card'],[], plot=False)
         
         
+        commandline = 'import model %s' % model.get('modelpath+restriction') 
+        if not model.mg5_name:
+            commandline += ' --modelname'
+        cmd.exec_cmd(commandline)
+
         line = 'compute_widths %s %s' % \
                 (' '.join([str(i) for i in opts['particles']]),
                  ' '.join('--%s=%s' % (key,value) for (key,value) in opts.items()
                         if key not in ['model', 'force', 'particles'] and value))
-        cmd.exec_cmd('import model %s' % model.get('modelpath+restriction'))
 
         #pattern for checking complex mass scheme.
         has_cms = re.compile(r'''set\s+complex_mass_scheme\s*(True|T|1|true|$|;)''', re.M)

=== modified file 'Template/LO/Source/.make_opts'
--- Template/LO/Source/.make_opts	2016-12-02 21:56:34 +0000
+++ Template/LO/Source/.make_opts	2017-07-13 10:23:41 +0000
@@ -37,9 +37,6 @@
 CFLAGS= -O $(STDLIB_FLAG) $(MACFLAG)
 endif
 
-# Increase the number of allowed charcters in a Fortran line
-FFLAGS+= -ffixed-line-length-132
-
 # Set FC unless it's defined by an environment variable
 ifeq ($(origin FC),default)
 FC=$(DEFAULT_F_COMPILER)
@@ -48,6 +45,15 @@
 F2PY=$(DEFAULT_F2PY_COMPILER)
 endif
 
+# Increase the number of allowed charcters in a Fortran line
+ifeq ($(FC), ftn)
+FFLAGS+= -extend-source # for ifort type of compiler
+else
+FFLAGS+= -ffixed-line-length-132
+endif
+
+
+
 UNAME := $(shell uname -s)
 ifeq ($(origin LDFLAGS), undefined)
 LDFLAGS=$(STDLIB) $(MACFLAG)
@@ -91,7 +97,7 @@
 ifneq ($(lhapdf),)
 CXXFLAGS += $(shell $(lhapdf) --cppflags)
 alfas_functions=alfas_functions_lhapdf
-llhapdf+= -lLHAPDF
+llhapdf+= $(shell $(lhapdf) --libs) -lLHAPDF
 else
 alfas_functions=alfas_functions
 llhapdf=

=== modified file 'Template/LO/Source/dsample.f'
--- Template/LO/Source/dsample.f	2015-10-01 10:08:45 +0000
+++ Template/LO/Source/dsample.f	2017-07-13 10:23:41 +0000
@@ -997,6 +997,43 @@
 
       end subroutine write_discrete_grids
 
+      subroutine write_grid(name)
+c************************************************************************
+c     Write out the grid
+c************************************************************************
+      implicit none
+
+      character*(*) name
+
+      include 'genps.inc'
+
+      double precision tmean, trmean, tsigma
+      integer             dim, events, itm, kn, cur_it, invar, configs
+      common /sample_common/
+     .     tmean, trmean, tsigma, dim, events, itm, kn, cur_it, invar, configs
+
+      double precision twgt, maxwgt,swgt(maxevents)
+      integer                             lun, nw, itmin
+      common/to_unwgt/twgt, maxwgt, swgt, lun, nw, itmin
+
+      double precision    grid(2, ng, 0:maxinvar)
+      common /data_grid/ grid
+
+      double precision force_max_wgt
+      common/unwgt_secondary_max/force_max_wgt
+
+      integer i,j
+
+      open(26, file='ftn26',status='unknown')
+      write(26,fmt='(4f21.17)') ((grid(2,i,j),i=1,ng),j=1,invar)
+      write(26,*) twgt, force_max_wgt
+c                  write(26,fmt='(4f21.16)') (alpha(i),i=1,maxconfigs)
+      call write_discrete_grids(26,'ref')
+      close(26)
+      return
+      end
+
+
       subroutine read_discrete_grids(stream_id)
 c************************************************************************
 c     Write out the grid using the DiscreteSampler module
@@ -1714,7 +1751,7 @@
          if (kn .ge. max_events .and. non_zero .le. 5) then
             call none_pass(max_events)
          endif
-         if (non_zero .eq. events .or. (kn .gt. 200*events .and.
+         if (non_zero .ge. events .or. (kn .gt. 200*events .and.
      $        non_zero .gt. 5)) then
 
 c          # special mode where we store information to combine them
@@ -2179,13 +2216,7 @@
      .              /23X,11HCross sec =,e12.4/
      .              13X,21HChi**2 per DoF.     =,f12.4/1X,79(1H-))
                if (use_cut .ne. 0) then
-               open(26, file='ftn26',status='unknown')
-               write(26,fmt='(4f21.17)')
-     $              ((grid(2,i,j),i=1,ng),j=1,invar)
-               write(26,*) twgt, force_max_wgt
-               call write_discrete_grids(26,'ref')
-c               write(26,fmt='(4f21.17)') (alpha(i),i=1,maxconfigs)
-               close(26)
+                   call write_grid('ftn26')
                endif
                call sample_writehtm()
 c               open(unit=22,file=result_file,status='old',

=== modified file 'Template/LO/SubProcesses/genps.f'
--- Template/LO/SubProcesses/genps.f	2016-09-06 15:01:18 +0000
+++ Template/LO/SubProcesses/genps.f	2017-07-13 10:23:41 +0000
@@ -124,7 +124,6 @@
       integer sprop(maxsproc,-max_branch:-1,lmaxconfigs)
       integer tprid(-max_branch:-1,lmaxconfigs)
       common/to_sprop/sprop,tprid
-      integer          lwgt(0:maxconfigs,maxinvar)
       logical firsttime
 
       double precision xprop(3,nexternal),tprop(3,nexternal)
@@ -132,7 +131,7 @@
       integer imatch
       save maxwgt
 
-      integer ninvar, nconfigs
+      integer ninvar
       
 c
 c     External
@@ -142,6 +141,9 @@
 c
 c     Global
 c
+      integer          lwgt(0:maxconfigs,maxinvar)
+      common/to_lwgt/lwgt
+
       double precision pmass(nexternal)
       common/to_mass/  pmass
 
@@ -172,7 +174,6 @@
       double precision stot,m1,m2
       common/to_stot/stot,m1,m2
 
-      save lwgt
       save ndim,nfinal,nbranch,nparticles
 
       integer jfig,k
@@ -200,40 +201,7 @@
 c      write(*,*) 'using iconfig',iconfig
       if (firsttime) then
          firsttime=.false.
-         do i=1,nexternal
-            m(i)=pmass(i)
-         enddo
-c        Set stot
-         if (nincoming.eq.1) then
-            stot=m(1)**2
-         else
-            m1=m(1)
-            m2=m(2)
-            if (abs(lpp(1)) .eq. 1 .or. abs(lpp(1)) .eq. 2) m1 = 0.938d0
-            if (abs(lpp(2)) .eq. 1 .or. abs(lpp(2)) .eq. 2) m2 = 0.938d0
-            if (abs(lpp(1)) .eq. 3) m1 = 0.000511d0
-            if (abs(lpp(2)) .eq. 3) m2 = 0.000511d0
-            if(ebeam(1).lt.m1) ebeam(1)=m1
-            if(ebeam(2).lt.m2) ebeam(2)=m2
-            pi1(0)=ebeam(1)
-            pi1(3)=sqrt(max(ebeam(1)**2-m1**2, 0d0))
-            pi2(0)=ebeam(2)
-            pi2(3)=-sqrt(max(ebeam(2)**2-m2**2, 0d0))
-            stot=m1**2+m2**2+2*(pi1(0)*pi2(0)-pi1(3)*pi2(3))
-         endif
-         write(*,'(x,a,f13.2)') 'Set CM energy to ',sqrt(stot)
-c        Start graph mapping
-         do i=1,mapconfig(0)
-            if (mapconfig(i) .eq. iconfig) this_config=i
-         enddo
-         write(*,*) 'Mapping Graph',iconfig,' to config',this_config
-         iconfig = this_config
-         nconfigs = 1
-         mincfig=iconfig
-         maxcfig=iconfig
-         call map_invarients(minvar,nconfigs,ninvar,mincfig,maxcfig,nexternal,nincoming)
-         maxwgt=0d0
-c         write(*,'(a,12i4)') 'Summing configs',(isym(i),i=1,isym(0))
+         call configure_integral(this_config,mincfig, maxcfig, invar,maxwgt)
          nparticles   = nexternal
          nfinal       = nparticles-nincoming
          nbranch      = nparticles-2
@@ -241,51 +209,10 @@
          if (ndim .lt. 0) ndim = 0   !For 2->1 processes  tjs 5/24/2010
          if (abs(lpp(1)) .ge. 1) ndim=ndim+1
          if (abs(lpp(2)) .ge. 1) ndim=ndim+1
-         call set_peaks
-         if (.false. ) then
-            call find_matches(iconfig,isym(0))
-            write(*,'(a,12i4)') 'Summing configs',(isym(i),i=1,isym(0))
-         endif
-         if (.false.) then
-            i=1
-            do while (mapconfig(i) .ne. iconfig
-     $          .and. i .lt. mapconfig(0))
-               i=i+1
-            enddo
-         endif
-
+         do i=1,nexternal
+            m(i)=pmass(i)
+         enddo
          write(*,'(a,12e10.3)') ' Masses:',(m(i),i=1,nparticles)
-         do j=1,invar
-            lwgt(0,j)=0
-         enddo
-c
-c     Here we set up which diagrams contribute to each variable
-c     in principle more than 1 diagram can contribute to a variable
-c     if we believe they will have identical structure.
-c
-c         do i=1,mapconfig(0)
-         do i=mincfig,maxcfig
-c         do k=1,isym(0)
-c            i = isym(k)
-            write(*,'(15i4)') i,(minvar(j,i),j=1,ndim)
-            do j=1,ndim
-               ipole = minvar(j,i)
-               if (ipole .ne. 0) then
-                  n = lwgt(0,ipole)+1
-                  lwgt(n,ipole)=mapconfig(i)  
-                  lwgt(0,ipole)=n
-               endif
-            enddo
-         enddo
-
-c     Initialize dsig (needed for subprocess group running mode)
-         dum=dsig(0,0,1)
-
-      else
-         do i=1,11
-c            swidth(i)=-5d0         !tells us to use the same point over again
-         enddo
-c         swidth(10)=0d0
       endif                          !First_time
 
       if (.false.) then
@@ -548,6 +475,138 @@
       end
 
 
+      subroutine configure_integral(iconfig,mincfig,maxcfig,invar,maxwgt)
+c**************************************************************************
+c     inputs  iconfig   == Current configuration working on
+c     output  m
+c**************************************************************************
+
+      implicit none
+
+      include 'genps.inc'
+      include 'maxconfigs.inc'
+      include 'nexternal.inc'
+      include 'maxamps.inc'
+      include 'run.inc'
+
+c     local
+      double precision pi1(0:3),pi2(0:3),p0,p3
+      double precision dum
+      integer i,j,ipole,n
+      integer nbranch,ndim,nconfigs
+      integer ninvar
+      integer nparticles,nfinal
+
+
+c
+c     Arguments
+c
+      integer iconfig,mincfig,maxcfig,invar
+      double precision maxwgt
+c
+c     External
+c
+      double precision lambda,dot,dsig
+      logical passcuts
+
+
+      logical firsttime
+      data firsttime/.true./
+      save firsttime
+c
+c     global
+c
+      double precision M(-max_branch:max_particles)
+
+      double precision pmass(nexternal)
+      common/to_mass/  pmass
+
+      double precision stot,m1,m2
+      common/to_stot/stot,m1,m2
+
+      integer            mapconfig(0:lmaxconfigs), this_config
+      common/to_mconfigs/mapconfig, this_config
+
+      integer           Minvar(maxdim,lmaxconfigs)
+      common /to_invar/ Minvar
+
+      integer          lwgt(0:maxconfigs,maxinvar)
+      common/to_lwgt/lwgt
+
+      if (firsttime)then
+         firsttime=.false.
+         do i=1,nexternal
+            m(i)=pmass(i)
+         enddo
+c        Set stot
+         if (nincoming.eq.1) then
+            stot=m(1)**2
+         else
+            m1=m(1)
+            m2=m(2)
+            if (abs(lpp(1)) .eq. 1 .or. abs(lpp(1)) .eq. 2) m1 = 0.938d0
+            if (abs(lpp(2)) .eq. 1 .or. abs(lpp(2)) .eq. 2) m2 = 0.938d0
+            if (abs(lpp(1)) .eq. 3) m1 = 0.000511d0
+            if (abs(lpp(2)) .eq. 3) m2 = 0.000511d0
+            if(ebeam(1).lt.m1) ebeam(1)=m1
+            if(ebeam(2).lt.m2) ebeam(2)=m2
+            pi1(0)=ebeam(1)
+            pi1(3)=sqrt(max(ebeam(1)**2-m1**2, 0d0))
+            pi2(0)=ebeam(2)
+            pi2(3)=-sqrt(max(ebeam(2)**2-m2**2, 0d0))
+            stot=m1**2+m2**2+2*(pi1(0)*pi2(0)-pi1(3)*pi2(3))
+         endif
+         write(*,'(x,a,f13.2)') 'Set CM energy to ',sqrt(stot)
+         endif
+c        Start graph mapping
+         do i=1,mapconfig(0)
+            if (mapconfig(i) .eq. iconfig) this_config=i
+         enddo
+         write(*,*) 'Mapping Graph',iconfig,' to config',this_config
+         iconfig = this_config
+         nconfigs = 1
+         mincfig=iconfig
+         maxcfig=iconfig
+         call map_invarients(minvar,nconfigs,ninvar,mincfig,maxcfig,nexternal,nincoming)
+         maxwgt=0d0
+c         write(*,'(a,12i4)') 'Summing configs',(isym(i),i=1,isym(0))
+         nparticles   = nexternal
+         nfinal       = nparticles-nincoming
+         nbranch      = nparticles-2
+         ndim         = 3*nfinal-4
+         if (ndim .lt. 0) ndim = 0   !For 2->1 processes  tjs 5/24/2010
+         if (abs(lpp(1)) .ge. 1) ndim=ndim+1
+         if (abs(lpp(2)) .ge. 1) ndim=ndim+1
+         call set_peaks
+         do j=1,invar
+            lwgt(0,j)=0
+         enddo
+c
+c     Here we set up which diagrams contribute to each variable
+c     in principle more than 1 diagram can contribute to a variable
+c     if we believe they will have identical structure.
+c
+c         do i=1,mapconfig(0)
+         do i=mincfig,maxcfig
+c         do k=1,isym(0)
+c            i = isym(k)
+            write(*,'(15i4)') i,(minvar(j,i),j=1,ndim)
+            do j=1,ndim
+               ipole = minvar(j,i)
+               if (ipole .ne. 0) then
+                  n = lwgt(0,ipole)+1
+                  lwgt(n,ipole)=mapconfig(i)
+                  lwgt(0,ipole)=n
+               endif
+            enddo
+         enddo
+
+c     Initialize dsig (needed for subprocess group running mode)
+         dum=dsig(0,0,1)
+
+      return
+      end
+
       subroutine one_tree(itree,iconfig,nbranch,P,M,S,X,jac,pswgt)
 c************************************************************************
 c     Calculates the momentum for everything below in the tree until

=== modified file 'Template/LO/SubProcesses/refine.sh'
--- Template/LO/SubProcesses/refine.sh	2016-06-12 00:01:23 +0000
+++ Template/LO/SubProcesses/refine.sh	2017-07-13 10:23:41 +0000
@@ -89,7 +89,7 @@
      if [ "$keeplog" = true ] ; then
 	 echo "" >> $k; echo "ls status:" >> $k; ls >> $k	 
      else
-	 rm ftn26 > /dev/null
+	 rm ftn26 &> /dev/null
      fi
 
 

=== modified file 'Template/LO/SubProcesses/unwgt.f'
--- Template/LO/SubProcesses/unwgt.f	2016-09-03 13:26:00 +0000
+++ Template/LO/SubProcesses/unwgt.f	2017-07-13 10:23:41 +0000
@@ -603,7 +603,7 @@
       if (nincoming.eq.2) then
          if (xbk(1) .gt. 0d0 .and. xbk(1) .le. 1d0 .and.
      $       xbk(2) .gt. 0d0 .and. xbk(2) .le. 1d0) then
-            if(xbk(1).eq.1d0.or.pmass(1).eq.0d0) then
+            if(lpp(2).ne.0.and.(xbk(1).eq.1d0.or.pmass(1).eq.0d0)) then
                ! construct the beam momenta in each frame and compute the related (z)boost
                ebi(0) = p(0,1)/xbk(1) ! this assumes that particle 1 is massless or mass equal to beam
                ebi(1) = 0

=== modified file 'Template/NLO/FixedOrderAnalysis/HwU.f'
--- Template/NLO/FixedOrderAnalysis/HwU.f	2017-02-03 13:37:36 +0000
+++ Template/NLO/FixedOrderAnalysis/HwU.f	2017-07-13 10:23:41 +0000
@@ -1,7 +1,7 @@
 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
 C                                                                          C
 C                  HwU: Histograms with Uncertainties                      C
-C                   By Rikkert Frederix, Dec. 2014                         C
+C                 By Rikkert Frederix, 12-2014--05-2017                    C
 C                                                                          C
 C     Book, fill and write out histograms. Particularly suited for NLO     C
 C     computations with correlations between points (ie. event and         C
@@ -9,26 +9,41 @@
 C     and PDF uncertainties through reweighting).                          C
 C                                                                          C
 CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
+
+c The module contains effectively the common block with allocatable
+c variables (something not possible in old fortran version)
+      module HwU_variables
+         implicit none
+         integer :: max_plots,max_points,max_bins,nwgts,np
+         integer :: error_estimation=3
+         logical, allocatable :: booked(:)
+         integer, allocatable :: nbin(:),histi(:,:),p_bin(:),p_label(:)
+         character(len=50), allocatable :: title(:)
+         character(len=50), allocatable :: wgts_info(:)
+         double precision, allocatable :: histy(:,:,:),histy_acc(:,:,:)
+     $        ,histy2(:,:),histy_err(:,:),histxl(:,:),histxm(:,:)
+     $        ,step(:),p_wgts(:,:)
+         save
+      end module HwU_variables
+
       
 c To be called once at the start of each run. Initialises the packages
 c and sets the number of weights that need to be included for each point.
       subroutine HwU_inithist(nweights,wgt_info)
+      use HwU_variables
       implicit none
-      include "HwU.inc"
       integer i,nweights
       character*(*) wgt_info(*)
-      do i=1,max_plots
-         booked(i)=.false.
-      enddo
+      call HwU_deallocate_all
+      max_plots=0
+      max_points=0
+      max_bins=0
+      np=0
 c     Number of weights associated to each point. Note that the first
 c     weight should always be the 'central value' and it should not be
 c     zero if any of the other weights are non-zero.
       nwgts=nweights
-      if (nwgts.gt.max_wgts) then
-         write (*,*) 'ERROR: increase max_wgts in HwU histogramming'
-     $        ,max_wgts,nwgts
-         stop 1
-      endif
+      allocate(wgts_info(nwgts))
       do i=1,nwgts
          wgts_info(i)=wgt_info(i)
       enddo
@@ -48,10 +63,9 @@
 c       square root of the sum of the squares. Perform a weighted average
 c       iteration-by-iteration
 c     input 3: Same as input 2, but weighted average is same as from MINT
+      use HwU_variables
       implicit none
       integer input
-      integer error_estimation
-      common /HwU_common2/ error_estimation
       if (input.ge.0 .and. input.le.3) then
          error_estimation=input
       else
@@ -60,39 +74,20 @@
       endif
       return
       end
-      
-      block data HwU
-c set the default for the error estimation method. To reduce the size of
-c the executable, put the error_estimation variable in a separate common
-c block. If we would have included the 'HwU.inc' file here, that
-c complete common block seems to be included in the size executable
-c (approx. 110 MB).
-      integer error_estimation
-      common /HwU_common2/ error_estimation
-      data error_estimation /3/
-      end
-      
+            
 c Book the histograms at the start of the run. Give a 'label' (an
 c integer) that identifies the plot when filling it and a title
 c ('title_l') for each plot. Also the number of bins ('nbin_l') and the
 c plot range (from 'xmin' to 'xmax') should be given.
       subroutine HwU_book(label,title_l,nbin_l,xmin,xmax)
+      use HwU_variables
       implicit none
-      include "HwU.inc"
       integer label,nbin_l,i,j
       character*(*) title_l
       double precision xmin,xmax
-c     Check that label and number of bins are reasonable
-      if (label.gt.max_plots) then
-         write (*,*) 'ERROR: increase max_plots in HwU histogramming'
-     $        ,max_plots, label
-         stop 1
-      endif
-      if (nbin_l.gt.max_bins) then
-         write (*,*) 'ERROR: increase max_bins in HwU histogramming'
-     $        ,max_bins,nbin_l
-         stop 1
-      endif
+c     Allocate space for new histograms if needed      
+      call HwU_allocate_histo(label,nbin_l)
+c     Setup the histogram
       booked(label)=.true.
       title(label)=title_l
       nbin(label)=nbin_l
@@ -111,7 +106,6 @@
          histy2(label,i)=0d0
          histy_err(label,i)=0d0
       enddo
-      np=0
       return
       end
       
@@ -121,8 +115,8 @@
 c the 'HwU_inithist' subroutine. That means that each point should have
 c the same number of weights.
       subroutine HwU_fill(label,x,wgts)
+      use HwU_variables
       implicit none
-      include "HwU.inc"
       integer label,i,j,bin
       double precision x, wgts(*)
 c     If central weight is zero do not add this point.
@@ -147,11 +141,7 @@
       enddo
 c     If a new bin, add it to the list of points
       np=np+1
-      if (np.gt.max_points) then
-         write (*,*) 'ERROR: increase max_points in HwU histogramming'
-     $        ,max_points
-         stop 1
-      endif
+      call HwU_allocate_p
       p_label(np)=label
       p_bin(np)=bin
       do j=1,nwgts
@@ -159,7 +149,7 @@
       enddo
       return
       end
-
+      
 c Call after all correlated contributions for a give phase-space
 c point. I.e., every time you get a new set of random numbers from
 c MINT/VEGAS. It adds the current list of points to the histograms. Add
@@ -168,8 +158,8 @@
 c this way, correlations between events and counter-events can be
 c correctly taken into account.
       subroutine HwU_add_points
+      use HwU_variables
       implicit none
-      include "HwU.inc"
       integer i,j
       do i=1,np
          do j=1,nwgts
@@ -192,12 +182,11 @@
 c the current iteration so that they can be filled with the next
 c iteration.
       subroutine HwU_accum_iter(inclde,nPSpoints,values)
+      use HwU_variables
       implicit none
-      include "HwU.inc"
       logical inclde
       integer nPSpoints,label,i,j
-      double precision nPSinv,etot,vtot(max_wgts),niter,y_squared
-     $     ,values(2)
+      double precision nPSinv,etot,niter,y_squared,values(2)
       data niter /0d0/
       nPSinv = 1d0/dble(nPSpoints)
       if (inclde) niter = niter+1d0
@@ -226,8 +215,8 @@
 c intermediate stages this function can be called (together with
 c HwU_output) to write intermediate plots to disk.
       subroutine finalize_histograms(nPSpoints)
+      use HwU_variables
       implicit none
-      include "HwU.inc"
       integer label,nPSpoints,i,j
       double precision nPSinv,niter,dummy(2)
       nPSinv=1d0/dble(nPSpoints)
@@ -257,13 +246,13 @@
 c histograms the central value should not be zero if any of the other
 c weights are non-zero.
       subroutine accumulate_results(label,nPSinv,niter,values)
+      use HwU_variables
       implicit none
-      include "HwU.inc"
       integer label,i,j
-      double precision nPSinv,etot,vtot(max_wgts),niter,y_squared
+      double precision nPSinv,etot,niter,y_squared
      $     ,values(2),a1,a2
-      integer error_estimation
-      common /HwU_common2/ error_estimation
+      double precision,allocatable :: vtot(:)
+      if (.not. allocated(vtot)) allocate(vtot(nwgts))
       if (error_estimation.eq.2) then
 c     Use the weighted average bin-by-bin. This is not really justified
 c     for fNLO computations, because for bins with low statistics, the
@@ -395,14 +384,15 @@
 c Write the histograms to disk at the end of the run, multiplying the
 c output by 'xnorm'
       subroutine HwU_output(unit,xnorm)
+      use HwU_variables
       implicit none
-      include "HwU.inc"
       integer unit,i,j,label
       integer max_length
-      parameter (max_length=(max_wgts+3)*17)
-      character*(max_length) buffer
+      character(len=:), allocatable :: buffer
       character*4 str_nbin
       double precision xnorm
+      if (.not. allocated(buffer))
+     &     allocate(character(len=(nwgts+3)*17) :: buffer)
 c     column info: x_min, x_max, y (central value), dy, {extra
 c     weights}.
       write (unit,'(a$)') '##& xmin'
@@ -440,9 +430,156 @@
          write (unit,'(a)') ''
          write (unit,'(a)') ''         
       enddo
-      return
-      end
-
+      deallocate(buffer)
+      return
+      end
+
+c Clean all the allocatable variables:
+      subroutine HwU_deallocate_all
+      use HwU_variables
+      implicit none
+      if (allocated(wgts_info)) deallocate(wgts_info)
+      if (allocated(booked)) deallocate(booked)
+      if (allocated(title)) deallocate(title)
+      if (allocated(nbin)) deallocate(nbin)
+      if (allocated(step)) deallocate(step)
+      if (allocated(histxl)) deallocate(histxl)
+      if (allocated(histxm)) deallocate(histxm)
+      if (allocated(histy)) deallocate(histy)
+      if (allocated(histy_acc)) deallocate(histy_acc)
+      if (allocated(histi)) deallocate(histi)
+      if (allocated(histy2)) deallocate(histy2)
+      if (allocated(histy_err)) deallocate(histy_err)
+      if (allocated(p_bin)) deallocate(p_bin)
+      if (allocated(p_label)) deallocate(p_label)
+      if (allocated(p_wgts)) deallocate(p_wgts)
+      return
+      end
+
+
+      subroutine HwU_allocate_p
+      use HwU_variables
+      implicit none
+      integer,allocatable :: itemp1(:)
+      double precision, allocatable :: temp2(:,:)
+      if (.not. allocated(p_bin)) then
+         allocate(p_bin(max_plots))
+         allocate(p_label(max_plots))
+         allocate(p_wgts(nwgts,max_plots))
+         max_points=max_plots
+      else
+         if (np.gt.max_points) then
+c p_bin
+            allocate(itemp1(np+max_plots))
+            itemp1(1:max_points)=p_bin
+            call move_alloc(itemp1,p_bin)
+      
+c p_label
+            allocate(itemp1(np+max_plots))
+            itemp1(1:max_points)=p_label
+            call move_alloc(itemp1,p_label)
+c p_wgts
+            allocate(temp2(nwgts,np+max_plots))
+            temp2(1:nwgts,1:max_points)=p_wgts
+            call move_alloc(temp2,p_wgts)
+            max_points=np+max_plots
+         endif
+      endif
+      return
+      end
+
+      subroutine HwU_allocate_histo(label,nbin_l)
+      use HwU_variables
+      implicit none
+      logical,allocatable :: ltemp(:)
+      integer,allocatable :: itemp1(:),itemp2(:,:)
+      character(len=50),allocatable :: ctemp(:)
+      double precision, allocatable :: temp1(:),temp2(:,:),temp3(:,:,:)
+      integer label,i,nbin_l,label_max,nbin_max
+      logical debug
+      parameter (debug=.true.)
+c Check if variables are already allocated. If not, simply allocate a
+c single histogram
+      if (.not. allocated(booked)) then
+         allocate(booked(1))
+         booked(1)=.false.
+         allocate(title(1))
+         allocate(nbin(1))
+         allocate(step(1))
+         allocate(histxl(1,nbin_l))
+         allocate(histxm(1,nbin_l))
+         allocate(histy(nwgts,1,nbin_l))
+         allocate(histy_acc(nwgts,1,nbin_l))
+         allocate(histi(1,nbin_l))
+         allocate(histy2(1,nbin_l))
+         allocate(histy_err(1,nbin_l))
+         max_plots=1
+         max_bins=nbin_l
+      endif
+c If current label is greater than the plots already allocated, increase
+c the size of the allocated arrays. This is kind of slow, but shouldn't
+c really matter since it's only done at the start of a run.
+      if (label.gt.max_plots .or. nbin_l.gt.max_bins) then
+         label_max=max(label,max_plots)
+         nbin_max=max(nbin_l,max_bins)
+c booked
+         allocate(ltemp(label_max))
+         ltemp(1:max_plots)=booked
+         call move_alloc(ltemp,booked)
+         do i=max_plots+1,label_max
+            booked(i)=.false. ! histos have not yet been setup
+         enddo
+c title
+         allocate(ctemp(label_max))
+         ctemp(1:max_plots)=title
+         call move_alloc(ctemp,title)
+c nbin
+         allocate(itemp1(label_max))
+         itemp1(1:max_plots)=nbin
+         call move_alloc(itemp1,nbin)
+c step
+         allocate(temp1(label_max))
+         temp1(1:max_plots)=step
+         call move_alloc(temp1,step)
+c histxl
+         allocate(temp2(label_max,nbin_max))
+         temp2(1:max_plots,1:max_bins)=histxl
+         call move_alloc(temp2,histxl)
+c histxm
+         allocate(temp2(label_max,nbin_max))
+         temp2(1:max_plots,1:max_bins)=histxm
+         call move_alloc(temp2,histxm)
+c histy
+         allocate(temp3(nwgts,label_max,nbin_max))
+         temp3(1:nwgts,1:max_plots,1:max_bins)=histy
+         call move_alloc(temp3,histy)
+c histy_acc
+         allocate(temp3(nwgts,label_max,nbin_max))
+         temp3(1:nwgts,1:max_plots,1:max_bins)=histy_acc
+         call move_alloc(temp3,histy_acc)
+c histi
+         allocate(itemp2(label_max,nbin_max))
+         itemp2(1:max_plots,1:max_bins)=histi
+         call move_alloc(itemp2,histi)
+c histy2
+         allocate(temp2(label_max,nbin_max))
+         temp2(1:max_plots,1:max_bins)=histy2
+         call move_alloc(temp2,histy2)
+c histy_err
+         allocate(temp2(label_max,nbin_max))
+         temp2(1:max_plots,1:max_bins)=histy_err
+         call move_alloc(temp2,histy_err)
+c Update maximums
+         max_plots=label_max
+         max_bins=nbin_max
+      elseif (booked(label)) then
+         write (*,*) 'ERROR in HwU.f: histogram already booked',label
+         stop
+      endif
+      return
+      end
+      
+      
 c dummy subroutine
       subroutine accum(idummy)
       integer idummy
@@ -451,3 +588,6 @@
       subroutine addfil(string)
       character*(*) string
       end
+
+
+

=== removed file 'Template/NLO/FixedOrderAnalysis/HwU.inc'
--- Template/NLO/FixedOrderAnalysis/HwU.inc	2016-02-18 14:05:45 +0000
+++ Template/NLO/FixedOrderAnalysis/HwU.inc	1970-01-01 00:00:00 +0000
@@ -1,22 +0,0 @@
-*     -*-fortran-*-
-
-      integer max_plots,max_bins,max_wgts,max_points
-      parameter (max_plots=200)
-      parameter (max_bins=100)
-      parameter (max_wgts=1024)
-      parameter (max_points=max_plots*40)
-
-      logical booked(max_plots)
-      integer nbin(max_plots),nwgts,np,p_bin(max_points)
-     &     ,p_label(max_points),histi(max_plots,max_bins)
-      character*50 title(max_plots)
-      character*50 wgts_info(max_wgts)
-      double precision histy(max_wgts,max_plots,max_bins)
-     $     ,histy_acc(max_wgts,max_plots,max_bins),histy2(max_plots
-     $     ,max_bins),histy_err(max_plots,max_bins),histxl(max_plots
-     $     ,max_bins),histxm(max_plots,max_bins),step(max_plots)
-     $     ,p_wgts(max_wgts,max_points)
-
-      common/HwU_common/histy,histy_acc,histy2,histy_err,histxl,histxm
-     &     ,p_wgts,step,histi,nbin,p_bin,p_label,np,nwgts
-     &     ,booked,title,wgts_info

=== modified file 'Template/NLO/FixedOrderAnalysis/analysis_HwU_template.f'
--- Template/NLO/FixedOrderAnalysis/analysis_HwU_template.f	2014-12-03 11:51:14 +0000
+++ Template/NLO/FixedOrderAnalysis/analysis_HwU_template.f	2017-07-13 10:23:41 +0000
@@ -17,12 +17,12 @@
 c     o) The first argument is an integer that labels the histogram. In
 c     the analysis_end and analysis_fill subroutines this label is used
 c     to keep track of the histogram. The label should be a number
-c     between 1 and max_plots=200 (can be increased in HwU.inc).
+c     starting at 1 and be increased for each plot.
 c     o) The second argument is a string that will apear above the
 c     histogram. Do not use brackets "(" or ")" inside this string.
 c     o) The third, forth and fifth arguments are the number of bis, the
 c     lower edge of the first bin and the upper edge of the last
-c     bin. There is a maximum of 100 bins per histogram.
+c     bin.
 c     o) When including scale and/or PDF uncertainties, declare a
 c     histogram for each weight, and compute the uncertainties from the
 c     final set of histograms

=== removed symlink 'Template/NLO/MCatNLO/include/HwU.inc'
=== target was u'../../FixedOrderAnalysis/HwU.inc'
=== modified file 'Template/NLO/SubProcesses/check_poles.f'
--- Template/NLO/SubProcesses/check_poles.f	2017-04-13 09:51:06 +0000
+++ Template/NLO/SubProcesses/check_poles.f	2017-07-13 10:23:41 +0000
@@ -183,6 +183,7 @@
             enddo
           enddo
 
+          CALL UPDATE_AS_PARAM()
           call sborn(p_born, born)
           call sloopmatrix_thres(p_born,virt_wgts,tolerance,
      1 accuracies,return_code)

=== modified file 'Template/NLO/SubProcesses/collect_events.f'
--- Template/NLO/SubProcesses/collect_events.f	2016-02-24 12:45:15 +0000
+++ Template/NLO/SubProcesses/collect_events.f	2017-07-13 10:23:41 +0000
@@ -1,6 +1,6 @@
       program collect_events
       implicit none
-      character*120 string120,eventfile
+      character*512 string512,eventfile
       character*19 basicfile,nextbasicfile
       character*15 outputfile
       integer istep,i,numoffiles,nbunches,nevents,ievents,junit(80)
@@ -65,9 +65,9 @@
       nevents=0
       xtotal=0.d0
       do while (.true.)
-         read(10,'(120a)',err=2,end=2) string120
-         eventfile=string120(2:index(string120,'   '))
-         read(string120(index(string120,'   '):120),*)
+         read(10,'(512a)',err=2,end=2) string512
+         eventfile=string512(2:index(string512,'   '))
+         read(string512(index(string512,'   '):512),*)
      $    ievents,absxsec,xsecfrac
          if (ievents.eq.0) cycle
          nevents=nevents+ievents

=== modified file 'Template/NLO/SubProcesses/driver_mintFO.f'
--- Template/NLO/SubProcesses/driver_mintFO.f	2017-05-09 11:54:11 +0000
+++ Template/NLO/SubProcesses/driver_mintFO.f	2017-07-13 10:23:41 +0000
@@ -505,6 +505,7 @@
 
       subroutine update_fks_dir(nFKS)
       implicit none
+      include 'run.inc'
       integer nFKS
       integer              nFKSprocess
       common/c_nFKSprocess/nFKSprocess
@@ -513,6 +514,7 @@
       call leshouche_inc_chooser()
       call setcuts
       call setfksfactor(.false.)
+      if (ickkw.eq.3) call configs_and_props_inc_chooser()
       return
       end
       

=== modified file 'Template/NLO/SubProcesses/driver_mintMC.f'
--- Template/NLO/SubProcesses/driver_mintMC.f	2017-05-09 12:41:29 +0000
+++ Template/NLO/SubProcesses/driver_mintMC.f	2017-07-13 10:23:41 +0000
@@ -719,14 +719,7 @@
 c
 c     Here I want to set up with B.W. we map and which we don't
 c
-      dconfig = dconfig-iconfig
-      if (dconfig .eq. 0) then
-         write(*,*) 'Not subdividing B.W.'
-         lbw(0)=0
-      else
-         write(*,*) 'Error BW setting: not supported at NLO'
-         stop 1
-      endif
+      lbw(0)=0
  10   format( a)
  12   format( a,i4)
       end
@@ -1207,6 +1200,7 @@
 
       subroutine update_fks_dir(nFKS)
       implicit none
+      include 'run.inc'
       integer nFKS
       integer              nFKSprocess
       common/c_nFKSprocess/nFKSprocess
@@ -1215,6 +1209,7 @@
       call leshouche_inc_chooser()
       call setcuts
       call setfksfactor(.true.)
+      if (ickkw.eq.3) call configs_and_props_inc_chooser()
       return
       end
 

=== modified file 'Template/RWGTNLO/makefile'
--- Template/RWGTNLO/makefile	2016-05-26 23:36:59 +0000
+++ Template/RWGTNLO/makefile	2017-07-13 10:23:41 +0000
@@ -10,7 +10,7 @@
 rwgt2py.so: rwgt.f makefile $(LIBS)
 	touch __init__.py	
 	touch ../__init__.py
-	$(F2PY) $(LINKLIBS) $(LDFLAGS) -c  $(RWGT) -m rwgt2py --fcompiler=$(FC)
+	$(F2PY) $(LINKLIBS) $(LDFLAGS) -c  $(RWGT) -m rwgt2py --f77exec=$(FC)
 
 $(LIBS):
 	cd PDF; make

=== modified file 'UpdateNotes.txt'
--- UpdateNotes.txt	2017-05-24 22:17:29 +0000
+++ UpdateNotes.txt	2017-07-13 10:23:41 +0000
@@ -1,6 +1,14 @@
 Update notes for MadGraph5_aMC@NLO (in reverse time order)
 
-2.5.5(25/05/17)
+2.5.6 (XX/XX/XX)
+      OM: Fixing an issue, leading to NAN for some of the channel of integration for complicated processes.
+      RF: Refactor of the HwU.f fortran code. Gives more flexibility and potentially lower
+          memory requirements.
+      RF: Fixed bug #1694548 (problem with NLO for QCD-charged heavy vector bosons)
+      RF: Another fix (adding on a fix in 2.5.5) related to FxFx merging in case there are
+          diagrams with 1->3 decays
+
+2.5.5(26/05/17)
       OM: Fixing bug in the creation of the LO gridpack introduced in 2.4.3. Since 2.4.3 the generated 
           gridpack were lacking to include the generated grid for each channel. This does not lead to 
           bias but to a significant slow down of the associated gridpack.

=== modified file 'VERSION'
--- VERSION	2017-05-24 22:17:29 +0000
+++ VERSION	2017-07-13 10:23:41 +0000
@@ -1,5 +1,6 @@
-version = 2.5.5
-date = 2017-05-25
+version = 2.5.6.beta
+date = 2017-05-26
+
 
 
 

=== modified file 'aloha/create_aloha.py'
--- aloha/create_aloha.py	2017-02-23 11:28:51 +0000
+++ aloha/create_aloha.py	2017-07-13 10:23:41 +0000
@@ -407,19 +407,22 @@
     @staticmethod
     def mod_propagator_expression(tag, text):
         """Change the index of the propagator to match the current need"""
-        data = re.split(r'(\b[a-zA-Z]\w*?)\(([\'\w,\s]*?)\)',text)
 
+        data = re.split(r'(\b[a-zA-Z]\w*?)\(([\'\w,\s\"\+\-]*?)\)',text)
+        to_change = {}
+        for old, new in tag.items():
+            if isinstance(new, str):
+                new='\'%s\'' % new
+            else:
+                new = str(new)
+            to_change[r'%s' % old] = new
         pos=-2
         while pos +3 < len(data):
             pos = pos+3
             ltype = data[pos]
             if ltype != 'complex':
-                for old, new in tag.items():
-                    if isinstance(new, str):
-                        new='\'%s\'' % new
-                    else:
-                        new = str(new)
-                    data[pos+1] = re.sub(r'\b%s\b' % old, new, data[pos+1])
+                data[pos+1] = re.sub(r'\b(?<!-)(%s)\b' % '|'.join(to_change),
+                                     lambda x: to_change[x.group()], data[pos+1])
             data[pos+1] = '(%s)' % data[pos+1]
         text=''.join(data)
         return text

=== modified file 'madgraph/interface/amcatnlo_run_interface.py'
--- madgraph/interface/amcatnlo_run_interface.py	2017-05-15 09:41:27 +0000
+++ madgraph/interface/amcatnlo_run_interface.py	2017-07-13 10:23:41 +0000
@@ -1265,7 +1265,8 @@
         if self.param_card_iterator:
             param_card_iterator = self.param_card_iterator
             self.param_card_iterator = [] #avoid to next generate go trough here
-            param_card_iterator.store_entry(self.run_name, self.results.current['cross'])
+            param_card_iterator.store_entry(self.run_name, self.results.current['cross'],
+                                            error=self.results.current['error'])
             orig_name = self.run_name
             #go trough the scal
             with misc.TMP_variable(self, 'allow_notification_center', False):
@@ -1282,7 +1283,8 @@
                         argss[0] = mode
                     self.do_launch("", options=options, argss=argss, switch=switch)
                     #self.exec_cmd("launch -f ",precmd=True, postcmd=True,errorhandling=False)
-                    param_card_iterator.store_entry(self.run_name, self.results.current['cross'])
+                    param_card_iterator.store_entry(self.run_name, self.results.current['cross'],
+                                                    error=self.results.current['error'])
             #restore original param_card
             param_card_iterator.write(pjoin(self.me_dir,'Cards','param_card.dat'))
             name = misc.get_scan_name(orig_name, self.run_name)

=== modified file 'madgraph/interface/common_run_interface.py'
--- madgraph/interface/common_run_interface.py	2017-05-24 22:17:29 +0000
+++ madgraph/interface/common_run_interface.py	2017-07-13 10:23:41 +0000
@@ -110,13 +110,11 @@
         logger.info("      by condor cluster (since condor has it's own way to prevent it).")
 
     def help_plot(self):
-        logger.info("syntax: help [RUN] [%s] [-f]" % '|'.join(self._plot_mode))
+        logger.info("syntax: plot [RUN] [%s] [-f]" % '|'.join(self._plot_mode))
         logger.info("-- create the plot for the RUN (current run by default)")
         logger.info("     at the different stage of the event generation")
         logger.info("     Note than more than one mode can be specified in the same command.")
-        logger.info("   This require to have MadAnalysis and td require. By default")
-        logger.info("     if those programs are installed correctly, the creation")
-        logger.info("     will be performed automaticaly during the event generation.")
+        logger.info("   This requires to have MadAnalysis and td installed.")
         logger.info("   -f options: answer all question by default.")
 
     def help_compute_widths(self):
@@ -645,9 +643,7 @@
                 fsock = open(pjoin(me_dir,'RunWeb'),'w')
                 fsock.write(`pid`)
                 fsock.close()
-    
-                misc.Popen([os.path.relpath(pjoin(self.dirbin, 'gen_cardhtml-pl'), me_dir)],
-                            cwd=me_dir)
+                self.gen_card_html()
 
         self.to_store = []
         self.run_name = None
@@ -749,6 +745,16 @@
     def do_treatcards(self, line, amcatnlo=False):
         """Advanced commands: create .inc files from param_card.dat/run_card.dat"""
 
+
+        #ensure that the cluster/card are consistent
+        if hasattr(self, 'run_card'):
+            self.cluster.modify_interface(self)
+        else:   
+            try:
+                self.cluster.modify_interface(self)
+            except Exception, error:
+                misc.sprint(str(error))
+                
         keepwidth = False
         if '--keepwidth' in line:
             keepwidth = True
@@ -758,10 +764,7 @@
 
         if mode in ['run', 'all']:
             if not hasattr(self, 'run_card'):
-                if amcatnlo:
-                    run_card = banner_mod.RunCardNLO(opt['run_card'])
-                else:
-                    run_card = banner_mod.RunCard(opt['run_card'])
+                run_card = banner_mod.RunCard(opt['run_card'])
             else:
                 run_card = self.run_card
 
@@ -1164,6 +1167,17 @@
                 auto_open=False)
         return True
     
+    def gen_card_html(self):
+        """ """
+        devnull = open(os.devnull, 'w')        
+        try:
+            misc.call(['./bin/internal/gen_cardhtml-pl'], cwd=self.me_dir,
+                        stdout=devnull, stderr=devnull)
+        except Exception:
+            pass
+        devnull.close()
+            
+    
     def create_plot(self, mode='parton', event_path=None, output=None, tag=None):
         """create the plot"""
 
@@ -2216,7 +2230,8 @@
     #===============================================================================
     @staticmethod
     def get_MadAnalysis5_interpreter(mg5_path, ma5_path, mg5_interface=None, 
-                    logstream = sys.stdout, loglevel =logging.INFO, forced = True):
+                    logstream = sys.stdout, loglevel =logging.INFO, forced = True,
+                    compilation=False):
         """ Makes sure to correctly setup paths and constructs and return an MA5 path"""
         
         MA5path = os.path.normpath(pjoin(mg5_path,ma5_path)) 
@@ -2240,7 +2255,8 @@
             with misc.stdchannel_redirected(sys.stdout, os.devnull):
                 with misc.stdchannel_redirected(sys.stderr, os.devnull):
                     MA5_interpreter = MA5Interpreter(MA5path, LoggerLevel=loglevel,
-                                                     LoggerStream=logstream,forced=forced, no_compilation=True)
+                                                     LoggerStream=logstream,forced=forced, 
+                                                     no_compilation=not compilation)
         except Exception as e:
             logger.warning('MadAnalysis5 failed to start so that MA5 analysis will be skipped.')
             error=StringIO.StringIO()
@@ -2579,7 +2595,8 @@
                 self.options['madanalysis5_path'],
                 logstream=sys.stdout,
                 loglevel=100,
-                forced=True)
+                forced=True,
+                compilation=True)
 
 
         # If failed to start MA5, then just leave
@@ -3069,9 +3086,9 @@
                                 logger.info("cluster handling will be done with PLUGIN: %s" % plug,'$MG:color:BLACK')
                                 self.cluster = plugin.new_cluster[cluster_name](**opt)
                                 break
-                        else:
-                            continue
-                        break
+                    else:
+                        continue
+                    break
                 else:
                     raise self.InvalidCmd, "%s is not recognized as a supported cluster format." % cluster_name
                 
@@ -3186,14 +3203,8 @@
             self.update_status('', level=None)
         except Exception, error:
             pass
-        devnull = open(os.devnull, 'w')
-        try:
-            misc.call(['./bin/internal/gen_cardhtml-pl'], cwd=self.me_dir,
-                        stdout=devnull, stderr=devnull)
-        except Exception:
-            pass
-        devnull.close()
 
+        self.gen_card_html()
         return super(CommonRunCmd, self).do_quit(line)
 
     # Aliases
@@ -5224,9 +5235,9 @@
                 if extrapaths:
                     self.do_set('shower_card extrapaths %s ' % ' '.join(extrapaths))
                 else:
-                    self.do_set('shower_card extrapaths None ')   
-                    
-                        
+                    self.do_set('shower_card extrapaths None ') 
+
+
     def reask(self, *args, **opt):
         
         cmd.OneLinePathCompletion.reask(self,*args, **opt)
@@ -5550,7 +5561,7 @@
             fsock.write(text) 
         self.reload_card(path)
 
-        
+    
 
     def do_compute_widths(self, line):
         signal.alarm(0) # avoid timer if any

=== modified file 'madgraph/interface/extended_cmd.py'
--- madgraph/interface/extended_cmd.py	2017-04-06 21:08:27 +0000
+++ madgraph/interface/extended_cmd.py	2017-07-13 10:23:41 +0000
@@ -1150,14 +1150,15 @@
             raise 'No such file or directory %s' % args[1]
         
     
-    def check_answer_in_input_file(self, question_instance, default, path=False):
+    def check_answer_in_input_file(self, question_instance, default, path=False, line=None):
         """Questions can have answer in output file (or not)"""
 
         if not self.inputfile:
             return None# interactive mode
 
-        line = self.get_stored_line()
-        # line define if a previous answer was not answer correctly 
+        if line is None:
+            line = self.get_stored_line()
+            # line define if a previous answer was not answer correctly 
         if not line:
             try:
                 line = self.inputfile.next()
@@ -1198,6 +1199,12 @@
             possibility = [opt for opt in options if line.lower()==opt.lower()]
             if len (possibility)==1:
                 return possibility[0]
+        if '=' in line and ' ' in line.strip():
+            line,n = re.subn('\s*=\s*','=', line)
+            if n:
+                return self.check_answer_in_input_file(question_instance, default, path=path, line=line)
+            
+            
             
         # No valid answer provides
         if self.haspiping:
@@ -1205,8 +1212,8 @@
             return None # print the question and use the pipe
         else:
             logger.info(question_instance.question)
-            logger.warning('The answer to the previous question is not set in your input file')
-            logger.warning('Use %s value' % default)
+            logger.warning('found line : %s' % line)
+            logger.warning('This answer is not valid for current question. Keep it for next question and use here default: %s', default) 
             self.store_line(line)
             return str(default)
 
@@ -2135,7 +2142,12 @@
                 return self.reask()
             elif len(self.allow_arg)==0:
                 return True
-            elif not self.casesensitive:
+            elif ' ' in line.strip() and '=' in self.value:
+                line,n = re.subn(r'\s*=\s*', '=', line)
+                if n:
+                    self.default(line)
+                    return self.postcmd(stop, line)
+            if not self.casesensitive:
                 for ans in self.allow_arg:
                     if ans.lower() == self.value.lower():
                         self.value = ans
@@ -2143,6 +2155,8 @@
                         break
                 else:
                     raise Exception
+
+                
             else: 
                 raise Exception
         except Exception,error:

=== modified file 'madgraph/interface/loop_interface.py'
--- madgraph/interface/loop_interface.py	2017-04-06 21:08:27 +0000
+++ madgraph/interface/loop_interface.py	2017-07-13 10:23:41 +0000
@@ -890,7 +890,8 @@
     local_installer = ['ninja', 'collier']
     required = ['cuttools', 'iregi']
     order = ['cuttools', 'iregi', 'ninja', 'collier', 'golem', 'pjfry']
-    
+    bypassed = ['pjfry']
+
     @property
     def answer(self):
         return self.code
@@ -971,6 +972,8 @@
                   'required': 'will be installed (required)'}
         
         for i,key in enumerate(self.order,1):
+            if key in self.bypassed and self.code[key] == 'off':
+                continue
             if os.path.sep not in self.code[key]:
                 question += '%s. %%(start_blue)s%-9s %-5s %-13s%%(stop)s : %s%s\n' % \
                    tuple([i,]+descript[key]+[status[self.code[key]],]+\

=== modified file 'madgraph/interface/madevent_interface.py'
--- madgraph/interface/madevent_interface.py	2017-05-03 09:21:01 +0000
+++ madgraph/interface/madevent_interface.py	2017-07-13 10:23:41 +0000
@@ -425,6 +425,11 @@
         self.run_options_help([("--" + key,value[-1]) for (key,value) in \
                                self._survey_options.items()])
      
+     
+    def help_restart_gridpack(self):
+        logger.info("syntax: restart_gridpack --precision= --restart_zero")
+
+
     def help_launch(self):
         """exec generate_events for 2>N and calculate_width for 1>N"""
         logger.info("syntax: launch [run_name] [options])")
@@ -2014,6 +2019,57 @@
         return
 
     ############################################################################
+    
+    ############################################################################
+    def do_restart_gridpack(self, line):
+        """ syntax restart_gridpack --precision=1.0 --restart_zero
+        collect the result of the current run and relaunch each channel
+        not completed or optionally a completed one with a precision worse than 
+        a threshold (and/or the zero result channel)"""
+        
+    
+        args = self.split_arg(line)
+        # Check argument's validity
+        self.check_survey(args)
+    
+        # initialize / remove lhapdf mode
+        #self.run_card = banner_mod.RunCard(pjoin(self.me_dir, 'Cards', 'run_card.dat'))
+        #self.configure_directory()
+        
+        gensym = gen_ximprove.gensym(self)
+        
+        min_precision = 1.0
+        resubmit_zero=False
+        if '--precision=' in line:
+            s = line.index('--precision=') + len('--precision=')
+            arg=line[s:].split(1)[0]
+            min_precision = float(arg)
+        
+        if '--restart_zero' in line:
+            resubmit_zero = True
+            
+            
+        gensym.resubmit(min_precision, resubmit_zero)
+        self.monitor(run_type='All jobs submitted for gridpack', html=True)
+
+                        #will be done during the refine (more precisely in gen_ximprove)
+        cross, error = sum_html.make_all_html_results(self)
+        self.results.add_detail('cross', cross)
+        self.results.add_detail('error', error)  
+        self.exec_cmd("print_results %s" % self.run_name,
+                       errorhandling=False, printcmd=False, precmd=False, postcmd=False)      
+        
+        self.results.add_detail('run_statistics', dict(gensym.run_statistics))
+
+        
+        #self.exec_cmd('combine_events', postcmd=False)
+        #self.exec_cmd('store_events', postcmd=False)
+        self.exec_cmd('decay_events -from_cards', postcmd=False)
+        self.exec_cmd('create_gridpack', postcmd=False)
+        
+    
+
+    ############################################################################    
 
     ############################################################################
     def do_generate_events(self, line):
@@ -2633,12 +2689,14 @@
                 run_card = banner_mod.RunCard(opt['run_card'])
             else:
                 run_card = self.run_card
+            self.run_card = run_card
+            self.cluster.modify_interface(self)
             if self.ninitial == 1:
                 run_card['lpp1'] =  0
                 run_card['lpp2'] =  0
                 run_card['ebeam1'] = 0
                 run_card['ebeam2'] = 0
-            
+                
             # Ensure that the bias parameters has all the required input from the
             # run_card
             if run_card['bias_module'].lower() not in ['dummy','none']:
@@ -2887,7 +2945,6 @@
             self.pass_in_difficult_integration_mode()
         
         jobs, P_zero_result = ajobcreator.launch()
-        
         # Check if all or only some fails
         if P_zero_result:
             if len(P_zero_result) == len(subproc):
@@ -3279,10 +3336,15 @@
         bannerfile = lhe_parser.EventFile(pjoin(self.me_dir, 'Events', self.run_name, '.banner.tmp.gz'),'w')
         banner = banner_mod.Banner(lhe.banner)
         banner.modify_init_cross(cross)
+        banner.set_lha_strategy(-4)
         banner.write(bannerfile, close_tag=False)
         bannerfile.close()
         # replace the lhe file by the new one
-        os.system('cat %s %s > %s' %(bannerfile.name, output.name, lhe.name))
+        if lhe.name.endswith('.gz'):
+            os.system('cat %s %s > %s' %(bannerfile.name, output.name, lhe.name))
+        else:
+            os.system('cat %s %s > %s.gz' %(bannerfile.name, output.name, lhe.name))
+            os.remove(lhe.name)
         os.remove(bannerfile.name)
         os.remove(output.name)
         
@@ -3365,9 +3427,8 @@
                         os.remove(pjoin(G_path, 'ftn25'))
 
         # 3) Update the index.html
-        misc.call(['%s/gen_cardhtml-pl' % self.dirbin],
-                            cwd=pjoin(self.me_dir))
-        
+        self.gen_card_html()
+
         
         # 4) Move the Files present in Events directory
         E_path = pjoin(self.me_dir, 'Events')
@@ -4772,7 +4833,7 @@
         self.check_plot(args)
         logger.info('plot for run %s' % self.run_name)
         if not self.force:
-            self.ask_edit_cards([], args, plot=True)
+            self.ask_edit_cards(['plot_card.dat'], args, plot=True)
                 
         if any([arg in ['all','parton'] for arg in args]):
             filename = pjoin(self.me_dir, 'Events', self.run_name, 'unweighted_events.lhe')
@@ -5148,7 +5209,9 @@
         #see when the last file was modified
         time_mod = max([os.path.getctime(pjoin(self.me_dir,'Cards','run_card.dat')),
                         os.path.getctime(pjoin(self.me_dir,'Cards','param_card.dat'))])
-        if self.configured > time_mod and hasattr(self, 'random'):
+        if self.configured > time_mod and hasattr(self, 'random') and hasattr(self, 'run_card'):
+            #just ensure that cluster specific are correctly handled
+            self.cluster.modify_interface(self)
             return
         else:
             self.configured = time.time()
@@ -5257,6 +5320,7 @@
         for nb_proc,subdir in enumerate(subproc):
             Pdir = pjoin(self.me_dir, 'SubProcesses',subdir.strip())
             self.compile(['clean'], cwd=Pdir)
+        self.configured = time.time()
 
     ############################################################################
     ##  HELPING ROUTINE

=== modified file 'madgraph/interface/madgraph_interface.py'
--- madgraph/interface/madgraph_interface.py	2017-04-06 21:08:27 +0000
+++ madgraph/interface/madgraph_interface.py	2017-07-13 10:23:41 +0000
@@ -1193,7 +1193,7 @@
         if not args:
             if self._done_export:
                 mode = self.find_output_type(self._done_export[0])
-                if (self._done_export[1] == 'plugin' and mode not in self._export_formats):
+                if (self._done_export[1] == 'plugin' and mode in self._export_formats):
                     args.append(mode)
                     args.append(self._done_export[0])
                 elif self._done_export[1].startswith(mode):
@@ -5057,6 +5057,8 @@
                 try:
                     self._curr_model = import_ufo.import_model(args[1], prefix=prefix,
                         complex_mass_scheme=self.options['complex_mass_scheme'])
+                    if os.path.sep in args[1] and "import" in self.history[-1]:
+                        self.history[-1] = 'import model %s' % self._curr_model.get('modelpath+restriction')
                 except import_ufo.UFOImportError, error:
                     if 'not a valid UFO model' in str(error):
                         logger_stderr.warning('WARNING: %s' % error)

=== modified file 'madgraph/iolibs/export_fks.py'
--- madgraph/iolibs/export_fks.py	2017-05-09 13:42:19 +0000
+++ madgraph/iolibs/export_fks.py	2017-07-13 10:23:41 +0000
@@ -1009,6 +1009,13 @@
                         (ifstring, part.get_pdg_code(), part.get_anti_pdg_code())
             iflines_width += 'get_width_from_id=abs(%s)\n' % part.get('width')
 
+        # Make sure it compiles with an if-statement if the above lists are empty
+        if len(mass_particles)==0:
+            iflines_mass = 'if (.True.) then\n'
+
+        if len(width_particles)==0:
+            iflines_width = 'if (.True.) then\n'
+
         replace_dict = {'iflines_mass' : iflines_mass,
                         'iflines_width' : iflines_width}
 
@@ -2749,10 +2756,13 @@
         info_list = fks_born.get_fks_info_list()
         lines = []
         if info_list:
-            # if the reals have been generated, fill with the corresponding value of ij
+            # if the reals have been generated, fill with the corresponding value of ij if
+            # ij is massless, or with 0 if ij is massive (no collinear singularity)
+            ij_list = [info['fks_info']['ij']if \
+                    fks_born.born_matrix_element['processes'][0]['legs'][info['fks_info']['ij']-1]['massless'] \
+                    else 0 for info in info_list]
             lines.append('INTEGER IJ_VALUES(%d)' % len(info_list))
-            lines.append('DATA IJ_VALUES /' + \
-                         ', '.join(['%d' % info['fks_info']['ij'] for info in info_list]) + '/')
+            lines.append('DATA IJ_VALUES /' + ', '.join(['%d' % ij for ij in ij_list]) + '/')
         else:
             #otherwise just put the first leg
             lines.append('INTEGER IJ_VALUES(1)')

=== modified file 'madgraph/iolibs/export_v4.py'
--- madgraph/iolibs/export_v4.py	2017-05-24 22:17:29 +0000
+++ madgraph/iolibs/export_v4.py	2017-07-13 10:23:41 +0000
@@ -639,10 +639,17 @@
         else:
             model_line='''$(LIBDIR)libmodel.$(libext): MODEL param_card.inc\n\tcd MODEL; make    
 param_card.inc: ../Cards/param_card.dat\n\t../bin/madevent treatcards param\n'''
-        text = open(path).read() % {'libraries': set_of_lib, 'model':model_line} 
-        writer.write(text)
-        
-        return True
+        
+        replace_dict= {'libraries': set_of_lib, 
+                       'model':model_line,
+                       'additional_dsample': '',
+                       'additional_dependencies':''} 
+        
+        if writer:
+            text = open(path).read() % replace_dict
+            writer.write(text)
+            
+        return replace_dict
 
     #===========================================================================
     # write_nexternal_madspin
@@ -2149,10 +2156,18 @@
         path = pjoin(_file_path,'iolibs','template_files','madevent_makefile_source')
         set_of_lib = '$(LIBDIR)libdhelas.$(libext) $(LIBDIR)libmodel.$(libext)'
         model_line='''$(LIBDIR)libmodel.$(libext): MODEL\n\t cd MODEL; make\n'''
-        text = open(path).read() % {'libraries': set_of_lib, 'model':model_line} 
-        writer.write(text)
-        
-        return True
+
+        replace_dict= {'libraries': set_of_lib, 
+                       'model':model_line,
+                       'additional_dsample': '',
+                       'additional_dependencies':''} 
+
+        text = open(path).read() % replace_dict
+        
+        if writer:
+            writer.write(text)
+        
+        return replace_dict
 
     #===========================================================================
     # write_matrix_element_v4

=== modified file 'madgraph/iolibs/file_writers.py'
--- madgraph/iolibs/file_writers.py	2016-08-23 23:06:25 +0000
+++ madgraph/iolibs/file_writers.py	2017-07-13 10:23:41 +0000
@@ -19,6 +19,12 @@
 
 import re
 import collections
+try:
+    import madgraph
+except ImportError:
+        import internal.misc
+else:
+    import madgraph.various.misc as misc
 
 class FileWriter(file):
     """Generic Writer class. All writers should inherit from this class."""
@@ -392,6 +398,52 @@
 #===============================================================================
 # CPPWriter
 #===============================================================================
+
+
+    def remove_routine(self, text, fct_names, formatting=True):
+        """write the incoming text but fully removing the associate routine/function
+           text can be a path to a file, an iterator, a string
+           fct_names should be a list of functions to remove
+        """
+
+        f77_type = ['real*8', 'integer', 'double precision']
+        pattern = re.compile('^\s+(?:SUBROUTINE|(?:%(type)s)\s+function)\s+([a-zA-Z]\w*)' \
+                             % {'type':'|'.join(f77_type)}, re.I)
+        
+        removed = []
+        if isinstance(text, str):   
+            if '\n' in text:
+                text = text.split('\n')
+            else:
+                text = open(text)
+        if isinstance(fct_names, str):
+            fct_names = [fct_names]
+        
+        to_write=True     
+        for line in text:
+            fct = pattern.findall(line)
+            if fct:
+                if fct[0] in fct_names:
+                    to_write = False
+                else:
+                    to_write = True
+
+            if to_write:
+                if formatting:
+                    if line.endswith('\n'):
+                        line = line[:-1]
+                    self.writelines(line)
+                else:
+                    if not line.endswith('\n'):
+                        line = '%s\n' % line
+                    file.writelines(self, line)
+            else:
+                removed.append(line)
+                
+        return removed
+        
+
+
 class CPPWriter(FileWriter):
     """Routines for writing C++ lines. Keeps track of brackets,
     spaces, indentation and splitting of long lines"""

=== modified file 'madgraph/iolibs/files.py'
--- madgraph/iolibs/files.py	2016-01-22 15:02:25 +0000
+++ madgraph/iolibs/files.py	2017-07-13 10:23:41 +0000
@@ -202,8 +202,12 @@
             starting_dir = os.path.join(cwd, starting_dir)        
 
     # Remove existing link if necessary
-    if os.path.exists(os.path.join(starting_dir, name)):
-        os.remove(os.path.join(starting_dir, name))
+    path = os.path.join(starting_dir, name)
+    if os.path.exists(path):
+        if os.path.realpath(path) != os.path.realpath(file_pos):
+            os.remove(os.path.join(starting_dir, name))
+        else:
+            return
 
     if not abspath:
         target = os.path.relpath(file_pos, starting_dir)

=== modified file 'madgraph/iolibs/template_files/born_fks.inc'
--- madgraph/iolibs/template_files/born_fks.inc	2017-05-09 13:46:36 +0000
+++ madgraph/iolibs/template_files/born_fks.inc	2017-07-13 10:23:41 +0000
@@ -58,6 +58,7 @@
       common/c_nfksprocess/nfksprocess
       double precision       wgt_ME_born,wgt_ME_real
       common /c_wgt_ME_tree/ wgt_ME_born,wgt_ME_real
+      logical cond_ij
 C ----------
 C BEGIN CODE
 C ----------
@@ -65,11 +66,15 @@
       glu_ij = ij_values(nfksprocess)
       NTRY(nFKSprocess)=NTRY(nFKSprocess)+1
       if (NTRY(nFKSprocess).lt.2) then
+        if (glu_ij.eq.0) then
+      	  skip(nfksprocess)=0
+       else
          skip(nFKSprocess)=1			 
          do while(nhel(glu_ij ,skip(nFKSprocess)).ne.-NHEL(GLU_IJ ,1))
             skip(nFKSprocess)=skip(nFKSprocess)+1
          enddo
          skip(nFKSprocess)=skip(nFKSprocess)-1
+       endif
       endif
       DO JJ=1,NGRAPHS
           amp2(jj)=0d0
@@ -101,7 +106,11 @@
       ANS(2) = 0D0
       hel_fac=1d0
       DO IHEL=1,NCOMB
-        if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1)) then
+        ! the following lines are to avoid segfaults when glu_ij=0
+        cond_ij=skip(nfksprocess).eq.0
+        if (.not.cond_ij) cond_ij=cond_ij.or.nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1)
+        !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        if (cond_ij) then
           IF ((GOODHEL(IHEL,nFKSprocess) .OR. GOODHEL(IHEL+SKIP(nFKSprocess),nFKSprocess) .OR. NTRY(nFKSprocess) .LT. 2) ) THEN
             ANS(1)=ANS(1)+BORN(P1,NHEL(1,IHEL),IHEL,BORNTILDE,borns)
             ANS(2)=ANS(2)+BORNTILDE
@@ -174,6 +183,7 @@
       integer nfksprocess
       common/c_nfksprocess/nfksprocess
       integer step_hel
+      logical cond_ij
 %(ij_lines)s
 C  
 C COLOR DATA
@@ -188,20 +198,33 @@
       back_hel = nhel(glu_ij)
       borns(1) = 0d0
       borns(2) = 0d0
-      if (back_hel.ne.0) then
-        step_hel=-2*back_hel
+      if (glu_ij.ne.0) then
+        back_hel = nhel(glu_ij)
+        if (back_hel.ne.0) then
+          step_hel=-2*back_hel
+        else
+          step_hel=1
+        endif
       else
+        back_hel=0
         step_hel=1
       endif
       DO IHEL=back_hel,-back_hel,step_hel
-         IF (IHEL.EQ.back_hel.OR.NHEL(GLU_IJ).NE.0) THEN
-        if (nhel(glu_ij).ne.0) nhel(glu_ij) = ihel
+        if (glu_ij.ne.0) then
+          cond_ij=IHEL.EQ.back_hel.OR.NHEL(GLU_IJ).NE.0
+        else
+          cond_ij=IHEL.EQ.back_hel 
+        endif
+        IF (cond_ij) THEN
+        if (glu_ij.ne.0) then
+          if (nhel(glu_ij).ne.0) nhel(glu_ij) = ihel
+        endif
         if (.not. calculatedBorn) then
 %(helas_calls)s
         do i=1,ngraphs
-          if(ihel.eq.BACK_HEL)then
+          if(ihel.eq.back_hel)then
             saveamp(i,hell)=amp(i)
-          elseif(ihel.eq.-BACK_HEL)then
+          elseif(ihel.eq.-back_hel)then
             saveamp(i,hell+skip(nFKSprocess))=amp(i)
           else
             write(*,*) "ERROR #1 in born.f"
@@ -210,9 +233,9 @@
         enddo
         elseif (calculatedBorn) then
         do i=1,ngraphs
-          if(ihel.eq.BACK_HEL)then
+          if(ihel.eq.back_hel)then
             amp(i)=saveamp(i,hell)
-          elseif(ihel.eq.-BACK_HEL)then
+          elseif(ihel.eq.-back_hel)then
             amp(i)=saveamp(i,hell+skip(nFKSprocess))
           else
             write(*,*) "ERROR #1 in born.f"
@@ -245,7 +268,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      nhel(glu_ij) = back_hel
+      if (glu_ij.ne.0) nhel(glu_ij) = back_hel
       END
        
 

=== modified file 'madgraph/iolibs/template_files/loop_optimized/helas_calls_split.inc'
--- madgraph/iolibs/template_files/loop_optimized/helas_calls_split.inc	2016-04-05 23:57:17 +0000
+++ madgraph/iolibs/template_files/loop_optimized/helas_calls_split.inc	2017-07-13 10:23:41 +0000
@@ -76,8 +76,8 @@
 	  %(complex_dp_format)s W(20,NWAVEFUNCS)
 	  common/%(proc_prefix)sW/W
 
-	  %(complex_dp_format)s WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,0:NLOOPWAVEFUNCS)
-	  %(complex_dp_format)s PL(0:3,0:NLOOPWAVEFUNCS)
+	  %(complex_dp_format)s WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,-1:NLOOPWAVEFUNCS)
+	  %(complex_dp_format)s PL(0:3,-1:NLOOPWAVEFUNCS)
 	  common/%(proc_prefix)sWL/WL,PL
 
 ## if(not AmplitudeReduction){

=== modified file 'madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc'
--- madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc	2017-03-27 09:01:13 +0000
+++ madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc	2017-07-13 10:23:41 +0000
@@ -367,8 +367,8 @@
       %(complex_mp_format)s MPW(20,NWAVEFUNCS)
 	  common/%(proc_prefix)sMP_W/MPW
 
-	  %(complex_dp_format)s WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,0:NLOOPWAVEFUNCS)
-	  %(complex_dp_format)s PL(0:3,0:NLOOPWAVEFUNCS)
+	  %(complex_dp_format)s WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,-1:NLOOPWAVEFUNCS)
+	  %(complex_dp_format)s PL(0:3,-1:NLOOPWAVEFUNCS)
 	  common/%(proc_prefix)sWL/WL,PL
 
 ## if(not AmplitudeReduction){
@@ -597,12 +597,15 @@
 
 C SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C The index -1 is for the charge-conjugated fermions with flipped fermion flow.
   DO I=0,3
+    PL(I,-1)=DCMPLX(0.0d0,0.0d0)
     PL(I,0)=DCMPLX(0.0d0,0.0d0)
   ENDDO
   DO I=1,MAXLWFSIZE
     DO J=0,LOOPMAXCOEFS-1
 	  DO K=1,MAXLWFSIZE
+	    WL(I,J,K,-1)=(0.0d0,0.0d0)
 	    IF(I.EQ.K.AND.J.EQ.0) then
         	WL(I,J,K,0)=(1.0d0,0.0d0)
 		ELSE
@@ -617,6 +620,13 @@
   ENDIF
 ENDIF
 
+C This is the chare conjugate version of the unit 4-currents in the canonical cartesian basis.
+C This, for now, is only defined for 4-fermionic currents.
+WL(1,0,2,-1) = DCMPLX(-1.0d0,0.0d0)
+WL(2,0,1,-1) = DCMPLX(1.0d0,0.0d0)
+WL(3,0,4,-1) = DCMPLX(1.0d0,0.0d0)
+WL(4,0,3,-1) = DCMPLX(-1.0d0,0.0d0)
+
 C Make sure that lorentz rotation tests are not used if there is external loop wavefunction of spin 2 and that one specific helicity is asked
 NROTATIONS_DP_BU = NROTATIONS_DP
 NROTATIONS_QP_BU = NROTATIONS_QP

=== modified file 'madgraph/iolibs/template_files/loop_optimized/mp_compute_loop_coefs.inc'
--- madgraph/iolibs/template_files/loop_optimized/mp_compute_loop_coefs.inc	2016-09-08 23:15:34 +0000
+++ madgraph/iolibs/template_files/loop_optimized/mp_compute_loop_coefs.inc	2017-07-13 10:23:41 +0000
@@ -130,12 +130,12 @@
 	  %(complex_dp_format)s DPW(20,NWAVEFUNCS)
 	  common/%(proc_prefix)sW/DPW
 
-	  %(complex_mp_format)s WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,0:NLOOPWAVEFUNCS)
-	  %(complex_mp_format)s PL(0:3,0:NLOOPWAVEFUNCS)
+	  %(complex_mp_format)s WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,-1:NLOOPWAVEFUNCS)
+	  %(complex_mp_format)s PL(0:3,-1:NLOOPWAVEFUNCS)
 	  common/%(proc_prefix)sMP_WL/WL,PL
 
-      %(complex_dp_format)s DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,0:NLOOPWAVEFUNCS)
-	  %(complex_dp_format)s DP_PL(0:3,0:NLOOPWAVEFUNCS)
+      %(complex_dp_format)s DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,-1:NLOOPWAVEFUNCS)
+	  %(complex_dp_format)s DP_PL(0:3,-1:NLOOPWAVEFUNCS)
 	  common/%(proc_prefix)sWL/DP_WL,DP_PL
 
 ## if(not AmplitudeReduction){
@@ -222,8 +222,10 @@
 ENDDO
 
 DO I=0,3
+  PL(I,-1)=CMPLX(ZERO,ZERO,KIND=16)
   PL(I,0)=CMPLX(ZERO,ZERO,KIND=16)
   IF (.NOT.COMPUTE_INTEGRAND_IN_QP) THEN
+    DP_PL(I,-1)=DCMPLX(0.0d0,0.0d0)
     DP_PL(I,0)=DCMPLX(0.0d0,0.0d0)
   ENDIF
 ENDDO
@@ -234,6 +236,8 @@
 DO I=1,MAXLWFSIZE
   DO J=0,LOOPMAXCOEFS-1
     DO K=1,MAXLWFSIZE
+	  WL(I,J,K,-1)=(ZERO,ZERO)
+	  DP_WL(I,J,K,-1)=(0.0d0,0.0d0)
       IF (I.EQ.K.AND.J.EQ.0) THEN
         WL(I,J,K,0)=(1.0e0_16,ZERO)
       ELSE
@@ -250,6 +254,18 @@
   ENDDO
 ENDDO
 
+C This is the chare conjugate version of the unit 4-currents in the canonical cartesian basis.
+C This, for now, is only defined for 4-fermionic currents.
+WL(1,0,2,-1) = (-1.0e0_16,ZERO)
+WL(2,0,1,-1) = (1.0e0_16,ZERO)
+WL(3,0,4,-1) = (1.0e0_16,ZERO)
+WL(4,0,3,-1) = (-1.0e0_16,ZERO)
+DP_WL(1,0,2,-1) = DCMPLX(-1.0d0,0.0d0)
+DP_WL(2,0,1,-1) = DCMPLX(1.0d0,0.0d0)
+DP_WL(3,0,4,-1) = DCMPLX(1.0d0,0.0d0)
+DP_WL(4,0,3,-1) = DCMPLX(-1.0d0,0.0d0)
+
+
 DO K=1, 3
 ## if(not AmplitudeReduction){
   DO I=1,NCTAMPS

=== modified file 'madgraph/iolibs/template_files/loop_optimized/mp_helas_calls_split.inc'
--- madgraph/iolibs/template_files/loop_optimized/mp_helas_calls_split.inc	2016-04-05 23:57:17 +0000
+++ madgraph/iolibs/template_files/loop_optimized/mp_helas_calls_split.inc	2017-07-13 10:23:41 +0000
@@ -60,8 +60,8 @@
 	  %(complex_mp_format)s W(20,NWAVEFUNCS)
 	  common/%(proc_prefix)sMP_W/W
 
-	  %(complex_mp_format)s WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,0:NLOOPWAVEFUNCS)
-	  %(complex_mp_format)s PL(0:3,0:NLOOPWAVEFUNCS)
+	  %(complex_mp_format)s WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,-1:NLOOPWAVEFUNCS)
+	  %(complex_mp_format)s PL(0:3,-1:NLOOPWAVEFUNCS)
 	  common/%(proc_prefix)sMP_WL/WL,PL
 
 ## if(not AmplitudeReduction){

=== modified file 'madgraph/iolibs/template_files/madevent_combine_events.f'
--- madgraph/iolibs/template_files/madevent_combine_events.f	2015-04-05 23:25:07 +0000
+++ madgraph/iolibs/template_files/madevent_combine_events.f	2017-07-13 10:23:41 +0000
@@ -65,6 +65,9 @@
       data s_buff/7*''/
       data jseed/-1/
       data buffclus/max_particles*' '/
+      double precision bias_weight
+      logical impact_xsec
+      common/bias/bias_weight,impact_xsec
 c-----
 c  Begin Code
 c-----
@@ -118,7 +121,7 @@
 c
       I4 = 4
       R8 = 8
-      record_length = 4*I4+maxexternal*I4*7+maxexternal*5*R8+3*R8+
+      record_length = 4*I4+maxexternal*I4*7+maxexternal*5*R8+4*R8+
      &   1000+7*s_bufflen+max_particles*clus_bufflen
 C $B$ scratch_name $B$ !this is tag for automatic modification by MW
       filename='scratch'
@@ -164,7 +167,9 @@
             read(sfnum,rec=iarray(i)) wgt,n,
      &           ((ic(m,j),j=1,maxexternal),m=1,7),ievent,
      &           ((p(m,j),m=0,4),j=1,maxexternal),sscale,aqcd,aqed,
-     &           buff,(s_buff(j),j=1,7),(buffclus(j),j=1,max_particles)
+     &           buff,(s_buff(j),j=1,7),(buffclus(j),j=1,max_particles),
+     &           bias_weight
+            if(bias_weight.ne.1d0) impact_xsec=.false.
 c     Systematics info on/off
          if(s_buff(1)(1:7).eq.'<mgrwt>') then
             u_syst=.true.
@@ -205,7 +210,7 @@
             read(sfnum,rec=iarray(i)) wgt,n,
      &           ((ic(m,j),j=1,maxexternal),m=1,7),ievent,
      &           ((p(m,j),m=0,4),j=1,maxexternal),sscale,aqcd,aqed,
-     &        buff
+     &        buff       
             if (dabs(wgt) .gt. goal_wgt*xran1(jseed)) then
                keep(i) = .true.
                if (wgt.lt.0d0) has_negative = .true.
@@ -260,7 +265,8 @@
             read(sfnum,rec=iarray(i)) wgt,n,
      &           ((ic(m,j),j=1,maxexternal),m=1,7),ievent,
      &           ((p(m,j),m=0,4),j=1,maxexternal),sscale,aqcd,aqed,
-     &           buff,(s_buff(j),j=1,7),(buffclus(j),j=1,max_particles)
+     &           buff,(s_buff(j),j=1,7),(buffclus(j),j=1,max_particles),
+     &           bias_weight
             wgt=dsign(xsec/nreq,wgt)
 c     Systematics info on/off
             if(s_buff(1)(1:7).eq.'<mgrwt>') then
@@ -638,6 +644,10 @@
       character*(clus_bufflen) buffclus(max_particles)
       data buffclus/max_particles*' '/
 c
+      double precision bias_weight
+      logical impact_xsec
+      common/bias/bias_weight,impact_xsec
+c
 c     Les Houches init block (for the <init> info)
 c
       integer maxpup
@@ -696,7 +706,8 @@
                write(sfnum,rec=kevent) wgt,n,
      &           ((ic(m,j),j=1,maxexternal),m=1,7),ievent,
      &           ((p(m,j),m=0,4),j=1,maxexternal),sscale,aqcd,aqed,
-     &           buff,(s_buff(j),j=1,7),(buffclus(j),j=1,max_particles)
+     &           buff,(s_buff(j),j=1,7),(buffclus(j),j=1,max_particles),
+     &           bias_weight
                sum=sum+dabs(wgt)
                found=.false.
                do i=1,nprup

=== modified file 'madgraph/iolibs/template_files/madevent_makefile_source'
--- madgraph/iolibs/template_files/madevent_makefile_source	2016-11-18 09:46:17 +0000
+++ madgraph/iolibs/template_files/madevent_makefile_source	2017-07-13 10:23:41 +0000
@@ -11,7 +11,7 @@
 # Source files
 
 PROCESS= hfill.o matrix.o myamp.o 
-DSAMPLE = dsample.o ranmar.o DiscreteSampler.o StringCast.o
+DSAMPLE = dsample.o ranmar.o DiscreteSampler.o StringCast.o %(additional_dsample)s
 HBOOK = hfill.o hcurve.o hbook1.o hbook2.o
 GENERIC = $(alfas_functions).o transpole.o invarients.o hfill.o pawgraphs.o ran1.o \
           rw_events.o rw_routines.o kin_functions.o open_file.o basecode.o setrun.o \
@@ -71,7 +71,7 @@
 select_events.o: select_events.f run_config.inc
 setrun.o: setrun.f nexternal.inc leshouche.inc run_card.inc run_config.inc
 rw_events.o: rw_events.f run_config.inc
-
+%(additional_dependencies)s
 run_card.inc: ../Cards/run_card.dat
 	../bin/madevent treatcards run
 

=== modified file 'madgraph/iolibs/template_files/makefile_sa_f_sp'
--- madgraph/iolibs/template_files/makefile_sa_f_sp	2016-11-13 13:28:20 +0000
+++ madgraph/iolibs/template_files/makefile_sa_f_sp	2017-07-13 10:23:41 +0000
@@ -24,4 +24,4 @@
 
 matrix$(MENUM)py.so: matrix.f makefile $(LIBS)
 	touch __init__.py
-	$(F2PY) $(LINKLIBS) -c  matrix.f -m matrix$(MENUM)py
\ No newline at end of file
+	$(F2PY) $(LINKLIBS) -c  matrix.f -m matrix$(MENUM)py --f77exec=$(FC)
\ No newline at end of file

=== modified file 'madgraph/iolibs/template_files/matrix_standalone_splitOrders_v4.inc'
--- madgraph/iolibs/template_files/matrix_standalone_splitOrders_v4.inc	2016-03-22 20:44:25 +0000
+++ madgraph/iolibs/template_files/matrix_standalone_splitOrders_v4.inc	2017-07-13 10:23:41 +0000
@@ -116,8 +116,9 @@
 	  COMMON/%(proc_prefix)sHELUSERCHOICE/USERHEL
       
 	  INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/%(proc_prefix)sBORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/%(proc_prefix)sBORN_BEAM_POL/POLARIZATIONS
+
 C
 C     FUNCTIONS
 C

=== modified file 'madgraph/iolibs/template_files/matrix_standalone_v4.inc'
--- madgraph/iolibs/template_files/matrix_standalone_v4.inc	2016-03-22 20:44:25 +0000
+++ madgraph/iolibs/template_files/matrix_standalone_v4.inc	2017-07-13 10:23:41 +0000
@@ -88,8 +88,9 @@
 %(den_factor_line)s
 
 	  INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/%(proc_prefix)sBORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/%(proc_prefix)sBORN_BEAM_POL/POLARIZATIONS
+
 C
 C     FUNCTIONS
 C

=== modified file 'madgraph/iolibs/template_files/super_auto_dsig_group_v4.inc'
--- madgraph/iolibs/template_files/super_auto_dsig_group_v4.inc	2015-10-12 19:04:37 +0000
+++ madgraph/iolibs/template_files/super_auto_dsig_group_v4.inc	2017-07-13 10:23:41 +0000
@@ -281,6 +281,7 @@
       TOTWGT=0D0
       DO J=1,SYMCONF(0)
         DO I=1,MAXSPROC
+	  IF(CONFSUB(I,SYMCONF(J)).NE.0) THEN
           DO K=1,2
             TOTWGT=TOTWGT+SELPROC(K,I,J)
             IF(R.LT.TOTWGT)THEN
@@ -290,6 +291,7 @@
                 GOTO 50
              ENDIF
           ENDDO
+	  ENDIF
         ENDDO
       ENDDO
  50   CONTINUE

=== modified file 'madgraph/loop/loop_helas_objects.py'
--- madgraph/loop/loop_helas_objects.py	2017-05-24 17:46:52 +0000
+++ madgraph/loop/loop_helas_objects.py	2017-07-13 10:23:41 +0000
@@ -2119,8 +2119,13 @@
                 loopwfnumber=loopwfnumber+1
             for loopamp in loopdiag.get_loop_amplitudes():
                 # Set the number of the starting wavefunction (common to all 
-                # diagrams) to one.
-                loopamp.get_starting_loop_wavefunction().set('number',0)
+                # diagrams) to 0 or -1 if it requires complex conjugation. 
+                start_loop_wf = loopamp.get_starting_loop_wavefunction()
+                if start_loop_wf.get('fermionflow')==1:
+                    start_loop_wf.set('number',0)
+                else:
+                    # External loop WF for flipped fermionflow.
+                    start_loop_wf.set('number',-1)
                 for amp in loopamp['amplitudes']:
                     amp.set('number',loop_ampnumber)
                     loop_ampnumber=loop_ampnumber+1

=== modified file 'madgraph/madevent/combine_runs.py'
--- madgraph/madevent/combine_runs.py	2016-09-03 13:26:00 +0000
+++ madgraph/madevent/combine_runs.py	2017-07-13 10:23:41 +0000
@@ -117,7 +117,7 @@
         #back out with the appropriate scaled weight
         fsock = open(pjoin(channel, 'events.lhe'), 'w')
         wgt = results.axsec / results.nunwgt
-        tot_nevents=0
+        tot_nevents, nb_file = 0, 0
         for result in results:  
             i = result.name
             if channel.endswith(os.path.pathsep):
@@ -125,8 +125,9 @@
             else:
                 path = channel + i
             nw = self.copy_events(fsock, pjoin(path,'events.lhe'), wgt)
-            #tot_events += nw
-        #logger.debug("Combined %s events to %s " % (tot_events, channel))
+            tot_nevents += nw
+            nb_file += 1
+        logger.debug("Combined %s file generating %s events for %s " , nb_file, tot_nevents, channel)
 
     @staticmethod
     def get_fortran_str(nb):
@@ -143,8 +144,10 @@
 
         new_wgt = self.get_fortran_str(new_wgt)
         old_line = ""
+        nb_evt =0 
         for line in open(input):
             if old_line.startswith("<event>"):
+                nb_evt+=1
                 data = line.split()
                 if not len(data) == 6:
                     raise MadGraph5Error, "Line after <event> should have 6 entries"
@@ -156,7 +159,7 @@
                                            new_wgt, '  '.join(data[3:]))
             fsock.write(line)
             old_line = line
-           
+        return nb_evt
     def get_channels(self, proc_path):
         """Opens file symfact.dat to determine all channels"""
         sympath = os.path.join(proc_path, 'symfact.dat')

=== modified file 'madgraph/madevent/gen_ximprove.py'
--- madgraph/madevent/gen_ximprove.py	2017-03-27 08:37:47 +0000
+++ madgraph/madevent/gen_ximprove.py	2017-07-13 10:23:41 +0000
@@ -104,13 +104,15 @@
         #if the user defines it in the run_card:
         if self.run_card['survey_splitting'] != -1:
             self.splitted_grid = self.run_card['survey_splitting']
+        if self.run_card['survey_nchannel_per_job'] != -1:
+            self.combining_job = self.run_card['survey_nchannel_per_job']        
         
         self.splitted_Pdir = {}
         self.splitted_for_dir = lambda x,y: self.splitted_grid
         self.combining_job_for_Pdir = lambda x: self.combining_job
         self.lastoffset = {}
     
-    def launch(self):
+    def launch(self, to_submit=True, clean=True):
         """ """
         
         self.subproc = [l.strip() for l in open(pjoin(self.me_dir,'SubProcesses', 
@@ -120,8 +122,8 @@
         P_zero_result = [] # check the number of times where they are no phase-space
         
         nb_tot_proc = len(subproc)
+        job_list = {}        
         for nb_proc,subdir in enumerate(subproc):
-            job_list = {}
             self.cmd.update_status('Compiling for process %s/%s. <br> (previous processes already running)' % \
                                (nb_proc+1,nb_tot_proc), level=None)
 
@@ -130,15 +132,16 @@
             logger.info('    %s ' % subdir)
             
             # clean previous run
-            for match in misc.glob('*ajob*', Pdir):
-                if os.path.basename(match)[:4] in ['ajob', 'wait', 'run.', 'done']:
-                    os.remove(match)
-            for match in misc.glob('G*', Pdir):
-                if os.path.exists(pjoin(match,'results.dat')):
-                    os.remove(pjoin(match, 'results.dat')) 
-                if os.path.exists(pjoin(match, 'ftn25')):
-                    os.remove(pjoin(match, 'ftn25')) 
-                    
+            if clean:
+                for match in misc.glob('*ajob*', Pdir):
+                    if os.path.basename(match)[:4] in ['ajob', 'wait', 'run.', 'done']:
+                        os.remove(match)
+                for match in misc.glob('G*', Pdir):
+                    if os.path.exists(pjoin(match,'results.dat')):
+                        os.remove(pjoin(match, 'results.dat')) 
+                    if os.path.exists(pjoin(match, 'ftn25')):
+                        os.remove(pjoin(match, 'ftn25')) 
+                        
             #compile gensym
             self.cmd.compile(['gensym'], cwd=Pdir)
             if not os.path.exists(pjoin(Pdir, 'gensym')):
@@ -179,9 +182,55 @@
                     raise Exception, 'Parsing error in gensym: %s' % stdout
                      
             self.cmd.compile(['madevent'], cwd=Pdir)
-            self.submit_to_cluster(job_list)
+            if to_submit:
+                self.submit_to_cluster(job_list)
+                job_list = {}
+                
         return job_list, P_zero_result
             
+    def resubmit(self, min_precision=1.0, resubmit_zero=False):
+        """collect the result of the current run and relaunch each channel
+        not completed or optionally a completed one with a precision worse than 
+        a threshold (and/or the zero result channel)"""
+        
+        job_list, P_zero_result = self.launch(to_submit=False, clean=False)
+        
+        for P , jobs in dict(job_list).items():
+            misc.sprint(jobs)
+            to_resub = []
+            for job in jobs:
+                if os.path.exists(pjoin(P, 'G%s' % job)) and os.path.exists(pjoin(P, 'G%s' % job, 'results.dat')):
+                    one_result = sum_html.OneResult(job)
+                    try:
+                        one_result.read_results(pjoin(P, 'G%s' % job, 'results.dat'))
+                    except:
+                        to_resub.append(job)
+                    if one_result.xsec == 0:
+                        if resubmit_zero:
+                            to_resub.append(job)
+                    elif max(one_result.xerru, one_result.xerrc)/one_result.xsec > min_precision:
+                        to_resub.append(job)
+                else:
+                    to_resub.append(job)   
+            if to_resub:
+                for G in to_resub:
+                    try:
+                        shutil.rmtree(pjoin(P, 'G%s' % G))
+                    except Exception, error:
+                        misc.sprint(error)
+                        pass
+            misc.sprint(to_resub) 
+            self.submit_to_cluster({P: to_resub})
+                    
+                    
+                    
+                    
+                
+                
+                
+        
+        
+           
             
     def submit_to_cluster(self, job_list):
         """ """

=== modified file 'madgraph/madevent/sum_html.py'
--- madgraph/madevent/sum_html.py	2017-02-14 10:10:56 +0000
+++ madgraph/madevent/sum_html.py	2017-07-13 10:23:41 +0000
@@ -436,16 +436,21 @@
         if update_statistics:
             self.run_statistics.aggregate_statistics([_.run_statistics for _ in self])
 
-    def compute_average(self):
+    def compute_average(self, error=None):
         """compute the value associate to this combination"""
 
         nbjobs = len(self)
         if not nbjobs:
             return
+        max_xsec = max(one.xsec for one in self)
+        min_xsec = min(one.xsec for one in self)
         self.axsec = sum([one.axsec for one in self]) / nbjobs
         self.xsec = sum([one.xsec for one in self]) /nbjobs
         self.xerrc = sum([one.xerrc for one in self]) /nbjobs
         self.xerru = math.sqrt(sum([one.xerru**2 for one in self])) /nbjobs
+        if error:
+            self.xerrc = error
+            self.xerru = error
 
         self.nevents = sum([one.nevents for one in self])
         self.nw = 0#sum([one.nw for one in self])
@@ -464,6 +469,20 @@
             self.eff_iter += result.eff_iter
             self.maxwgt_iter += result.maxwgt_iter
 
+        #check full consistency
+        onefail = False
+        for one in list(self):
+            if one.xsec < (self.xsec - 25* one.xerru):
+                if not onefail:
+                    logger.debug('multi run are inconsistent: %s < %s - 25* %s: assign error %s', one.xsec, self.xsec, one.xerru, error if error else max_xsec-min_xsec)
+                onefail = True
+                self.remove(one)
+        if onefail:
+            if error:
+                return self.compute_average(error)
+            else:
+                return self.compute_average((max_xsec-min_xsec)/2.)
+            
 
     
     def compute_iterations(self):

=== modified file 'madgraph/various/banner.py'
--- madgraph/various/banner.py	2017-03-03 10:01:29 +0000
+++ madgraph/various/banner.py	2017-07-13 10:23:41 +0000
@@ -149,7 +149,7 @@
                     self[tag] = text
                     text = ''
                     store = False
-            if store:
+            if store and not line.startswith(('<![CDATA[',']]>')):
                 if line.endswith('\n'):
                     text += line
                 else:
@@ -382,18 +382,19 @@
         ff.write(header % { 'version':float(self.lhe_version)})
 
 
-        for tag in [t for t in self.ordered_items if t in self.keys()]:
-            if tag in exclude: 
-                continue
-            capitalized_tag = self.capitalized_items[tag] if tag in self.capitalized_items else tag
-            ff.write('<%(tag)s>\n%(text)s\n</%(tag)s>\n' % \
-                     {'tag':capitalized_tag, 'text':self[tag].strip()})
-        for tag in [t for t in self.keys() if t not in self.ordered_items]:
-            if tag in ['init'] or tag in exclude:
-                continue
-            capitalized_tag = self.capitalized_items[tag] if tag in self.capitalized_items else tag
-            ff.write('<%(tag)s>\n%(text)s\n</%(tag)s>\n' % \
-                     {'tag':capitalized_tag, 'text':self[tag].strip()})
+        for tag in [t for t in self.ordered_items if t in self.keys()]+ \
+            [t for t in self.keys() if t not in self.ordered_items]:
+            if tag in ['init'] or tag in exclude: 
+                continue
+            capitalized_tag = self.capitalized_items[tag] if tag in self.capitalized_items else tag
+            start_data, stop_data = '', ''
+            if '<' in self[tag] or '@' in self[tag]:
+                start_data = '\n<![CDATA['
+                stop_data = ']]>\n'
+            ff.write('<%(tag)s>%(start_data)s\n%(text)s\n%(stop_data)s</%(tag)s>\n' % \
+                     {'tag':capitalized_tag, 'text':self[tag].strip(),
+                      'start_data': start_data, 'stop_data':stop_data})
+        
         
         if not '/header' in exclude:
             ff.write('</header>\n')    
@@ -2618,7 +2619,9 @@
         #job handling of the survey/ refine
         self.add_param('job_strategy', 0, hidden=True, include=False)
         self.add_param('survey_splitting', -1, hidden=True, include=False)
+        self.add_param('survey_nchannel_per_job', 2, hidden=True, include=False)
         self.add_param('refine_evt_by_job', -1, hidden=True, include=False)
+        
         # Specify what particle IDs to use for the CKKWL merging cut ktdurham
         
     def check_validity(self):

=== modified file 'madgraph/various/cluster.py'
--- madgraph/various/cluster.py	2017-03-13 11:26:51 +0000
+++ madgraph/various/cluster.py	2017-07-13 10:23:41 +0000
@@ -257,11 +257,11 @@
             target = path.rsplit('/SubProcesses',1)[0]
         elif 'MCatNLO' in path:
             target = path.rsplit('/MCatNLO',1)[0]
+        elif 'PY8_parallelization' in path:
+            target = path.rsplit('/PY8_parallelization',1)[0]
         elif second_path:
             target=path
             logger.warning("cluster.get_job_identifier runs unexpectedly. This should be fine but report this message if you have problem.")
-        elif 'PY8_parallelization' in path:
-            target = path.rsplit('/PY8_parallelization',1)[0]
         else:
             target = path
             
@@ -537,6 +537,15 @@
         logger.warning("""This cluster didn't support metajob submit.""")
         return 0
 
+    def modify_interface(self, run_interface):
+        """routine which allow to modify the run_card/mg5cmd object to change the
+           default behavior of the runs.
+           This is called at the time of the compilation of the run_card. 
+           Note that this function can be called multiple times by run.
+           """
+        #run_card = run_interface.run_card
+        return 
+
 class Packet(object):
     """ an object for handling packet of job, it is designed to be thread safe
     """

=== modified file 'madgraph/various/lhe_parser.py'
--- madgraph/various/lhe_parser.py	2017-05-03 09:21:01 +0000
+++ madgraph/various/lhe_parser.py	2017-07-13 10:23:41 +0000
@@ -174,6 +174,10 @@
     def __init__(self, path, mode='r', *args, **opt):
         """open file and read the banner [if in read mode]"""
         
+        if path.endswith('.gz') and mode == 'w' and\
+                                              isinstance(self, EventFileNoGzip):
+            path = path[:-3]
+        
         self.parsing = True # check if/when we need to parse the event.
         self.eventgroup  = False
         try:

=== modified file 'madgraph/various/misc.py'
--- madgraph/various/misc.py	2017-05-22 09:11:33 +0000
+++ madgraph/various/misc.py	2017-07-13 10:23:41 +0000
@@ -1755,8 +1755,12 @@
 It has been validated for the last time with version: %s""",
                                         name, '.'.join(str(i) for i in val_ver))
     else:
-        logger.error("Plugin %s is not supported by this version of MG5aMC." % name)
-        plugin_support[name] = False
+        if __debug__:
+            logger.error("Plugin %s seems not supported by this version of MG5aMC. Keep it active (please update status)" % name)
+            plugin_support[name] = True            
+        else:
+            logger.error("Plugin %s is not supported by this version of MG5aMC." % name)
+            plugin_support[name] = False
     return plugin_support[name]
     
 

=== modified file 'mg5decay/decay_objects.py'
--- mg5decay/decay_objects.py	2017-02-14 15:08:19 +0000
+++ mg5decay/decay_objects.py	2017-07-13 10:23:41 +0000
@@ -1992,6 +1992,9 @@
     def running_externals(self, q, loopnum=2):
         """ Recalculate external parameters at the given scale. """
         
+        if q < 0.5:
+            return
+
         # Raise error for wrong type of q
         if not isinstance(q, int) and not isinstance(q, long) and \
                 not isinstance(q, float):

=== modified file 'models/check_param_card.py'
--- models/check_param_card.py	2017-03-03 23:36:24 +0000
+++ models/check_param_card.py	2017-07-13 10:23:41 +0000
@@ -900,15 +900,17 @@
             yield param_card
         
     
-    def store_entry(self, run_name, cross):
+    def store_entry(self, run_name, cross, error=None):
         """store the value of the cross-section"""
         if isinstance(cross, dict):
             info = dict(cross)
             info.update({'bench' : self.itertag, 'run_name': run_name})
             self.cross.append(info)
         else:
-            self.cross.append({'bench' : self.itertag, 'run_name': run_name, 'cross(pb)':cross})
-        
+            if error is None:
+                self.cross.append({'bench' : self.itertag, 'run_name': run_name, 'cross(pb)':cross})
+            else:
+                self.cross.append({'bench' : self.itertag, 'run_name': run_name, 'cross(pb)':cross, 'error(pb)':error})        
 
     def write_summary(self, path, order=None, lastline=False, nbcol=20):
         """ """

=== modified file 'models/import_ufo.py'
--- models/import_ufo.py	2017-02-16 15:16:38 +0000
+++ models/import_ufo.py	2017-07-13 10:23:41 +0000
@@ -78,10 +78,9 @@
                     last_model_path = os.path.join(MG5DIR, p, model_name)
                 return os.path.join(MG5DIR, p, model_name)
     if os.path.isdir(model_name):
-        if last_model_path != os.path.join(MG5DIR, model_name):
-            logger.info("model loaded from: %s", os.path.join(os.getcwd(), model_name))
-            last_model_path = os.path.join(MG5DIR, model_name)
-        return model_name   
+        logger.warning('No model %s found in default path. Did you mean \'import model ./%s\'',
+                       model_name, model_name)
+        raise UFOImportError("Path %s is not a valid pathname" % model_name)    
     else:
         raise UFOImportError("Path %s is not a valid pathname" % model_name)    
     

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_gg_ttx%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_gg_ttx%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_gg_ttx%born.f	2017-07-13 10:23:41 +0000
@@ -77,6 +77,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -84,11 +85,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -121,7 +126,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -201,6 +211,7 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(1)
       DATA IJ_VALUES / 1 /
 C     
@@ -221,14 +232,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL VXXXXX(P(0,1),ZERO,NHEL(1),-1*IC(1),W(1,1))
             CALL VXXXXX(P(0,2),ZERO,NHEL(2),-1*IC(2),W(1,2))
@@ -292,7 +316,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_uux_ttx%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_uux_ttx%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_uux_ttx%born.f	2017-07-13 10:23:41 +0000
@@ -80,6 +80,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -87,11 +88,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -124,7 +129,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -207,6 +217,7 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(1)
       DATA IJ_VALUES / 1 /
 C     
@@ -227,14 +238,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL IXXXXX(P(0,1),ZERO,NHEL(1),+1*IC(1),W(1,1))
             CALL OXXXXX(P(0,2),ZERO,NHEL(2),-1*IC(2),W(1,2))
@@ -292,7 +316,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_uxu_ttx%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_uxu_ttx%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fks_loonly/%SubProcesses%P0_uxu_ttx%born.f	2017-07-13 10:23:41 +0000
@@ -80,6 +80,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -87,11 +88,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -124,7 +129,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -207,6 +217,7 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(1)
       DATA IJ_VALUES / 1 /
 C     
@@ -227,14 +238,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL OXXXXX(P(0,1),ZERO,NHEL(1),-1*IC(1),W(1,1))
             CALL IXXXXX(P(0,2),ZERO,NHEL(2),+1*IC(2),W(1,2))
@@ -292,7 +316,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_gg_ttx%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_gg_ttx%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_gg_ttx%born.f	2017-07-13 10:23:41 +0000
@@ -55,7 +55,7 @@
       INTEGER IDEN_VALUES(8)
       DATA IDEN_VALUES /256, 256, 256, 256, 256, 256, 256, 256/
       INTEGER IJ_VALUES(8)
-      DATA IJ_VALUES /1, 2, 3, 4, 1, 1, 2, 2/
+      DATA IJ_VALUES /1, 2, 0, 0, 1, 1, 2, 2/
 C     
 C     GLOBAL VARIABLES
 C     
@@ -77,6 +77,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -84,11 +85,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -121,7 +126,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -201,8 +211,9 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(8)
-      DATA IJ_VALUES /1, 2, 3, 4, 1, 1, 2, 2/
+      DATA IJ_VALUES /1, 2, 0, 0, 1, 1, 2, 2/
 C     
 C     COLOR DATA
 C     
@@ -221,14 +232,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL VXXXXX(P(0,1),ZERO,NHEL(1),-1*IC(1),W(1,1))
             CALL VXXXXX(P(0,2),ZERO,NHEL(2),-1*IC(2),W(1,2))
@@ -292,7 +316,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_uux_ttx%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_uux_ttx%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_uux_ttx%born.f	2017-07-13 10:23:41 +0000
@@ -58,7 +58,7 @@
       INTEGER IDEN_VALUES(6)
       DATA IDEN_VALUES /36, 36, 36, 36, 36, 36/
       INTEGER IJ_VALUES(6)
-      DATA IJ_VALUES /1, 2, 3, 4, 1, 2/
+      DATA IJ_VALUES /1, 2, 0, 0, 1, 2/
 C     
 C     GLOBAL VARIABLES
 C     
@@ -80,6 +80,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -87,11 +88,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -124,7 +129,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -207,8 +217,9 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(6)
-      DATA IJ_VALUES /1, 2, 3, 4, 1, 2/
+      DATA IJ_VALUES /1, 2, 0, 0, 1, 2/
 C     
 C     COLOR DATA
 C     
@@ -227,14 +238,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL IXXXXX(P(0,1),ZERO,NHEL(1),+1*IC(1),W(1,1))
             CALL OXXXXX(P(0,2),ZERO,NHEL(2),-1*IC(2),W(1,2))
@@ -292,7 +316,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_uxu_ttx%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_uxu_ttx%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_pptt_fksreal/%SubProcesses%P0_uxu_ttx%born.f	2017-07-13 10:23:41 +0000
@@ -58,7 +58,7 @@
       INTEGER IDEN_VALUES(6)
       DATA IDEN_VALUES /36, 36, 36, 36, 36, 36/
       INTEGER IJ_VALUES(6)
-      DATA IJ_VALUES /1, 2, 3, 4, 1, 2/
+      DATA IJ_VALUES /1, 2, 0, 0, 1, 2/
 C     
 C     GLOBAL VARIABLES
 C     
@@ -80,6 +80,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -87,11 +88,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -124,7 +129,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -207,8 +217,9 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(6)
-      DATA IJ_VALUES /1, 2, 3, 4, 1, 2/
+      DATA IJ_VALUES /1, 2, 0, 0, 1, 2/
 C     
 C     COLOR DATA
 C     
@@ -227,14 +238,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL OXXXXX(P(0,1),ZERO,NHEL(1),-1*IC(1),W(1,1))
             CALL IXXXXX(P(0,2),ZERO,NHEL(2),+1*IC(2),W(1,2))
@@ -292,7 +316,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%born_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%born_matrix.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%born_matrix.f	2017-07-13 10:23:41 +0000
@@ -130,8 +130,9 @@
       COMMON/HELUSERCHOICE/USERHEL
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%coef_construction_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_CT_calls_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_CT_calls_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_CT_calls_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f	2017-03-28 07:06:02 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f	2017-07-13 10:23:41 +0000
@@ -353,9 +353,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -583,12 +583,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -603,6 +607,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_coef_construction_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_compute_loop_coefs.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_compute_loop_coefs.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_compute_loop_coefs.f	2017-07-13 10:23:41 +0000
@@ -119,14 +119,14 @@
       COMPLEX*16 DPW(20,NWAVEFUNCS)
       COMMON/W/DPW
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MP_WL/WL,PL
 
-      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 DP_PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/DP_WL,DP_PL
 
       COMPLEX*32 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -191,8 +191,10 @@
       ENDDO
 
       DO I=0,3
+        PL(I,-1)=CMPLX(ZERO,ZERO,KIND=16)
         PL(I,0)=CMPLX(ZERO,ZERO,KIND=16)
         IF (.NOT.COMPUTE_INTEGRAND_IN_QP) THEN
+          DP_PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           DP_PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDIF
       ENDDO
@@ -200,6 +202,8 @@
       DO I=1,MAXLWFSIZE
         DO J=0,LOOPMAXCOEFS-1
           DO K=1,MAXLWFSIZE
+            WL(I,J,K,-1)=(ZERO,ZERO)
+            DP_WL(I,J,K,-1)=(0.0D0,0.0D0)
             IF (I.EQ.K.AND.J.EQ.0) THEN
               WL(I,J,K,0)=(1.0E0_16,ZERO)
             ELSE
@@ -216,6 +220,19 @@
         ENDDO
       ENDDO
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = (-1.0E0_16,ZERO)
+      WL(2,0,1,-1) = (1.0E0_16,ZERO)
+      WL(3,0,4,-1) = (1.0E0_16,ZERO)
+      WL(4,0,3,-1) = (-1.0E0_16,ZERO)
+      DP_WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      DP_WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      DP_WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      DP_WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
+
       DO K=1, 3
         DO I=1,NCTAMPS
           AMPL(K,I)=(ZERO,ZERO)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%mp_helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%born.f	2017-07-13 10:23:41 +0000
@@ -74,6 +74,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -81,11 +82,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -118,7 +123,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -199,6 +209,7 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(4)
       DATA IJ_VALUES /1, 2, 1, 2/
 C     
@@ -216,14 +227,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL OXXXXX(P(0,1),ZERO,NHEL(1),-1*IC(1),W(1,1))
             CALL IXXXXX(P(0,2),ZERO,NHEL(2),+1*IC(2),W(1,2))
@@ -278,7 +302,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%born_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%born_matrix.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%born_matrix.f	2017-07-13 10:23:41 +0000
@@ -130,8 +130,9 @@
       COMMON/HELUSERCHOICE/USERHEL
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%coef_construction_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_CT_calls_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_CT_calls_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_CT_calls_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_matrix.f	2017-03-28 07:06:02 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_matrix.f	2017-07-13 10:23:41 +0000
@@ -353,9 +353,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -583,12 +583,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -603,6 +607,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_coef_construction_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_compute_loop_coefs.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_compute_loop_coefs.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_compute_loop_coefs.f	2017-07-13 10:23:41 +0000
@@ -119,14 +119,14 @@
       COMPLEX*16 DPW(20,NWAVEFUNCS)
       COMMON/W/DPW
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MP_WL/WL,PL
 
-      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 DP_PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/WL/DP_WL,DP_PL
 
       COMPLEX*32 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -191,8 +191,10 @@
       ENDDO
 
       DO I=0,3
+        PL(I,-1)=CMPLX(ZERO,ZERO,KIND=16)
         PL(I,0)=CMPLX(ZERO,ZERO,KIND=16)
         IF (.NOT.COMPUTE_INTEGRAND_IN_QP) THEN
+          DP_PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           DP_PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDIF
       ENDDO
@@ -200,6 +202,8 @@
       DO I=1,MAXLWFSIZE
         DO J=0,LOOPMAXCOEFS-1
           DO K=1,MAXLWFSIZE
+            WL(I,J,K,-1)=(ZERO,ZERO)
+            DP_WL(I,J,K,-1)=(0.0D0,0.0D0)
             IF (I.EQ.K.AND.J.EQ.0) THEN
               WL(I,J,K,0)=(1.0E0_16,ZERO)
             ELSE
@@ -216,6 +220,19 @@
         ENDDO
       ENDDO
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = (-1.0E0_16,ZERO)
+      WL(2,0,1,-1) = (1.0E0_16,ZERO)
+      WL(3,0,4,-1) = (1.0E0_16,ZERO)
+      WL(4,0,3,-1) = (-1.0E0_16,ZERO)
+      DP_WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      DP_WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      DP_WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      DP_WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
+
       DO K=1, 3
         DO I=1,NCTAMPS
           AMPL(K,I)=(ZERO,ZERO)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%mp_helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%born.f	2017-07-13 10:23:41 +0000
@@ -74,6 +74,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -81,11 +82,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -118,7 +123,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -199,6 +209,7 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(4)
       DATA IJ_VALUES /1, 2, 1, 2/
 C     
@@ -216,14 +227,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL IXXXXX(P(0,1),ZERO,NHEL(1),+1*IC(1),W(1,1))
             CALL OXXXXX(P(0,2),ZERO,NHEL(2),-1*IC(2),W(1,2))
@@ -278,7 +302,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_tdecay_fksreal/%SubProcesses%P0_t_budx%born.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_tdecay_fksreal/%SubProcesses%P0_t_budx%born.f	2017-05-09 14:31:51 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_tdecay_fksreal/%SubProcesses%P0_t_budx%born.f	2017-07-13 10:23:41 +0000
@@ -56,7 +56,7 @@
       INTEGER IDEN_VALUES(3)
       DATA IDEN_VALUES /6, 6, 6/
       INTEGER IJ_VALUES(3)
-      DATA IJ_VALUES /2, 3, 4/
+      DATA IJ_VALUES /0, 3, 4/
 C     
 C     GLOBAL VARIABLES
 C     
@@ -78,6 +78,7 @@
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       DOUBLE PRECISION       WGT_ME_BORN,WGT_ME_REAL
       COMMON /C_WGT_ME_TREE/ WGT_ME_BORN,WGT_ME_REAL
+      LOGICAL COND_IJ
 C     ----------
 C     BEGIN CODE
 C     ----------
@@ -85,11 +86,15 @@
       GLU_IJ = IJ_VALUES(NFKSPROCESS)
       NTRY(NFKSPROCESS)=NTRY(NFKSPROCESS)+1
       IF (NTRY(NFKSPROCESS).LT.2) THEN
-        SKIP(NFKSPROCESS)=1
-        DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
-          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
-        ENDDO
-        SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        IF (GLU_IJ.EQ.0) THEN
+          SKIP(NFKSPROCESS)=0
+        ELSE
+          SKIP(NFKSPROCESS)=1
+          DO WHILE(NHEL(GLU_IJ ,SKIP(NFKSPROCESS)).NE.-NHEL(GLU_IJ ,1))
+            SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)+1
+          ENDDO
+          SKIP(NFKSPROCESS)=SKIP(NFKSPROCESS)-1
+        ENDIF
       ENDIF
       DO JJ=1,NGRAPHS
         AMP2(JJ)=0D0
@@ -122,7 +127,12 @@
       ANS(2) = 0D0
       HEL_FAC=1D0
       DO IHEL=1,NCOMB
-        IF (NHEL(GLU_IJ,IHEL).EQ.NHEL(GLU_IJ,1)) THEN
+          ! the following lines are to avoid segfaults when glu_ij=0
+        COND_IJ=SKIP(NFKSPROCESS).EQ.0
+        IF (.NOT.COND_IJ) COND_IJ=COND_IJ.OR.NHEL(GLU_IJ,IHEL)
+     $   .EQ.NHEL(GLU_IJ,1)
+          !if (nhel(glu_ij,ihel).EQ.NHEL(GLU_IJ,1).or.skip(nfksprocess).eq.0) then
+        IF (COND_IJ) THEN
           IF ((GOODHEL(IHEL,NFKSPROCESS) .OR. GOODHEL(IHEL
      $     +SKIP(NFKSPROCESS),NFKSPROCESS) .OR. NTRY(NFKSPROCESS) .LT.
      $      2) ) THEN
@@ -203,8 +213,9 @@
       INTEGER NFKSPROCESS
       COMMON/C_NFKSPROCESS/NFKSPROCESS
       INTEGER STEP_HEL
+      LOGICAL COND_IJ
       INTEGER IJ_VALUES(3)
-      DATA IJ_VALUES /2, 3, 4/
+      DATA IJ_VALUES /0, 3, 4/
 C     
 C     COLOR DATA
 C     
@@ -220,14 +231,27 @@
       BACK_HEL = NHEL(GLU_IJ)
       BORNS(1) = 0D0
       BORNS(2) = 0D0
-      IF (BACK_HEL.NE.0) THEN
-        STEP_HEL=-2*BACK_HEL
+      IF (GLU_IJ.NE.0) THEN
+        BACK_HEL = NHEL(GLU_IJ)
+        IF (BACK_HEL.NE.0) THEN
+          STEP_HEL=-2*BACK_HEL
+        ELSE
+          STEP_HEL=1
+        ENDIF
       ELSE
+        BACK_HEL=0
         STEP_HEL=1
       ENDIF
       DO IHEL=BACK_HEL,-BACK_HEL,STEP_HEL
-        IF (IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0) THEN
-          IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+        IF (GLU_IJ.NE.0) THEN
+          COND_IJ=IHEL.EQ.BACK_HEL.OR.NHEL(GLU_IJ).NE.0
+        ELSE
+          COND_IJ=IHEL.EQ.BACK_HEL
+        ENDIF
+        IF (COND_IJ) THEN
+          IF (GLU_IJ.NE.0) THEN
+            IF (NHEL(GLU_IJ).NE.0) NHEL(GLU_IJ) = IHEL
+          ENDIF
           IF (.NOT. CALCULATEDBORN) THEN
             CALL IXXXXX(P(0,1),MDL_MT,NHEL(1),+1*IC(1),W(1,1))
             CALL OXXXXX(P(0,2),MDL_MB,NHEL(2),+1*IC(2),W(1,2))
@@ -284,7 +308,7 @@
         ENDDO
         BORNTILDE = BORNTILDE + ZTEMP*DCONJG(JAMPH(1,I))/DENOM(I)
       ENDDO
-      NHEL(GLU_IJ) = BACK_HEL
+      IF (GLU_IJ.NE.0) NHEL(GLU_IJ) = BACK_HEL
       END
 
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_madevent_group/super_auto_dsig.f'
--- tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_madevent_group/super_auto_dsig.f	2016-02-27 14:33:01 +0000
+++ tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_madevent_group/super_auto_dsig.f	2017-07-13 10:23:41 +0000
@@ -303,15 +303,17 @@
         TOTWGT=0D0
         DO J=1,SYMCONF(0)
           DO I=1,MAXSPROC
-            DO K=1,2
-              TOTWGT=TOTWGT+SELPROC(K,I,J)
-              IF(R.LT.TOTWGT)THEN
-                IPROC=I
-                ICONF=J
-                IMIRROR=K
-                GOTO 50
-              ENDIF
-            ENDDO
+            IF(CONFSUB(I,SYMCONF(J)).NE.0) THEN
+              DO K=1,2
+                TOTWGT=TOTWGT+SELPROC(K,I,J)
+                IF(R.LT.TOTWGT)THEN
+                  IPROC=I
+                  ICONF=J
+                  IMIRROR=K
+                  GOTO 50
+                ENDIF
+              ENDDO
+            ENDIF
           ENDDO
         ENDDO
  50     CONTINUE

=== modified file 'tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_standalone/matrix.f'
--- tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_standalone/matrix.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_standalone/matrix.f	2017-07-13 10:23:41 +0000
@@ -121,8 +121,9 @@
       DATA IDEN/24/
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_eq_4.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_eq_4.f	2017-04-06 21:08:27 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_eq_4.f	2017-07-13 10:23:41 +0000
@@ -358,9 +358,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -591,12 +591,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -611,6 +615,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_0_QEDAmpAndQEDsq_gt_2.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_0_QEDAmpAndQEDsq_gt_2.f	2017-04-06 21:08:27 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_0_QEDAmpAndQEDsq_gt_2.f	2017-07-13 10:23:41 +0000
@@ -358,9 +358,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -591,12 +591,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -611,6 +615,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_4.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_4.f	2017-04-06 21:08:27 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_4.f	2017-07-13 10:23:41 +0000
@@ -352,9 +352,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -583,12 +583,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -603,6 +607,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QEDsq_le_4.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QEDsq_le_4.f	2017-04-06 21:08:27 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QEDsq_le_4.f	2017-07-13 10:23:41 +0000
@@ -358,9 +358,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -591,12 +591,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -611,6 +615,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_WGTsq_le_10_QEDAmpAndQEDsq_gt_2.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_WGTsq_le_10_QEDAmpAndQEDsq_gt_2.f	2017-04-06 21:08:27 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_WGTsq_le_10_QEDAmpAndQEDsq_gt_2.f	2017-07-13 10:23:41 +0000
@@ -358,9 +358,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -591,12 +591,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -611,6 +615,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_default.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_default.f	2017-04-06 21:08:27 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_default.f	2017-07-13 10:23:41 +0000
@@ -358,9 +358,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -591,12 +591,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -611,6 +615,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDpert_default.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDpert_default.f	2017-04-06 21:08:27 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDpert_default.f	2017-07-13 10:23:41 +0000
@@ -352,9 +352,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -583,12 +583,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -603,6 +607,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QEDpert_default.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QEDpert_default.f	2017-04-06 21:08:27 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QEDpert_default.f	2017-07-13 10:23:41 +0000
@@ -352,9 +352,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -583,12 +583,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -603,6 +607,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%born_matrix.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%born_matrix.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%born_matrix.f	2017-07-13 10:23:41 +0000
@@ -133,8 +133,9 @@
       COMMON/ML5_0_HELUSERCHOICE/USERHEL
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/ML5_0_BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/ML5_0_BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%coef_construction_1.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/ML5_0_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/ML5_0_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%helas_calls_uvct_1.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%helas_calls_uvct_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%helas_calls_uvct_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/ML5_0_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_CT_calls_1.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_CT_calls_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_CT_calls_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/ML5_0_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_matrix.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_matrix.f	2017-03-28 14:32:20 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_matrix.f	2017-07-13 10:23:41 +0000
@@ -352,9 +352,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -583,12 +583,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -603,6 +607,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_coef_construction_1.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_compute_loop_coefs.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_compute_loop_coefs.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_compute_loop_coefs.f	2017-07-13 10:23:41 +0000
@@ -118,14 +118,14 @@
       COMPLEX*16 DPW(20,NWAVEFUNCS)
       COMMON/ML5_0_W/DPW
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_MP_WL/WL,PL
 
-      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 DP_PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/DP_WL,DP_PL
 
       COMPLEX*32 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -190,8 +190,10 @@
       ENDDO
 
       DO I=0,3
+        PL(I,-1)=CMPLX(ZERO,ZERO,KIND=16)
         PL(I,0)=CMPLX(ZERO,ZERO,KIND=16)
         IF (.NOT.COMPUTE_INTEGRAND_IN_QP) THEN
+          DP_PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           DP_PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDIF
       ENDDO
@@ -199,6 +201,8 @@
       DO I=1,MAXLWFSIZE
         DO J=0,LOOPMAXCOEFS-1
           DO K=1,MAXLWFSIZE
+            WL(I,J,K,-1)=(ZERO,ZERO)
+            DP_WL(I,J,K,-1)=(0.0D0,0.0D0)
             IF (I.EQ.K.AND.J.EQ.0) THEN
               WL(I,J,K,0)=(1.0E0_16,ZERO)
             ELSE
@@ -215,6 +219,19 @@
         ENDDO
       ENDDO
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = (-1.0E0_16,ZERO)
+      WL(2,0,1,-1) = (1.0E0_16,ZERO)
+      WL(3,0,4,-1) = (1.0E0_16,ZERO)
+      WL(4,0,3,-1) = (-1.0E0_16,ZERO)
+      DP_WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      DP_WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      DP_WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      DP_WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
+
       DO K=1, 3
         DO I=1,NCTAMPS
           AMPL(K,I)=(ZERO,ZERO)

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_helas_calls_uvct_1.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_helas_calls_uvct_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%mp_helas_calls_uvct_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_NoSQSO.f'
--- tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_NoSQSO.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_NoSQSO.f	2017-07-13 10:23:41 +0000
@@ -153,8 +153,9 @@
       DATA IDEN/144/
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_QCDsq_le_6.f'
--- tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_QCDsq_le_6.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_QCDsq_le_6.f	2017-07-13 10:23:41 +0000
@@ -181,8 +181,9 @@
       COMMON/HELUSERCHOICE/USERHEL
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_ampOrderQED2_eq_2_WGTsq_le_14_QCDsq_gt_4.f'
--- tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_ampOrderQED2_eq_2_WGTsq_le_14_QCDsq_gt_4.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_ampOrderQED2_eq_2_WGTsq_le_14_QCDsq_gt_4.f	2017-07-13 10:23:41 +0000
@@ -181,8 +181,9 @@
       COMMON/HELUSERCHOICE/USERHEL
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%coef_construction_1.f'
--- tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/MG5_1_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MG5_1_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/MG5_1_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MG5_1_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_CT_calls_1.f'
--- tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_CT_calls_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_CT_calls_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/MG5_1_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MG5_1_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_matrix.f'
--- tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_matrix.f	2017-03-28 07:06:02 +0000
+++ tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_matrix.f	2017-07-13 10:23:41 +0000
@@ -352,9 +352,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/MG5_1_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MG5_1_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -583,12 +583,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -603,6 +607,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_coef_construction_1.f'
--- tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/MG5_1_MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MG5_1_MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_compute_loop_coefs.f'
--- tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_compute_loop_coefs.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_compute_loop_coefs.f	2017-07-13 10:23:41 +0000
@@ -118,14 +118,14 @@
       COMPLEX*16 DPW(20,NWAVEFUNCS)
       COMMON/MG5_1_W/DPW
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MG5_1_MP_WL/WL,PL
 
-      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 DP_PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MG5_1_WL/DP_WL,DP_PL
 
       COMPLEX*32 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -190,8 +190,10 @@
       ENDDO
 
       DO I=0,3
+        PL(I,-1)=CMPLX(ZERO,ZERO,KIND=16)
         PL(I,0)=CMPLX(ZERO,ZERO,KIND=16)
         IF (.NOT.COMPUTE_INTEGRAND_IN_QP) THEN
+          DP_PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           DP_PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDIF
       ENDDO
@@ -199,6 +201,8 @@
       DO I=1,MAXLWFSIZE
         DO J=0,LOOPMAXCOEFS-1
           DO K=1,MAXLWFSIZE
+            WL(I,J,K,-1)=(ZERO,ZERO)
+            DP_WL(I,J,K,-1)=(0.0D0,0.0D0)
             IF (I.EQ.K.AND.J.EQ.0) THEN
               WL(I,J,K,0)=(1.0E0_16,ZERO)
             ELSE
@@ -215,6 +219,19 @@
         ENDDO
       ENDDO
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = (-1.0E0_16,ZERO)
+      WL(2,0,1,-1) = (1.0E0_16,ZERO)
+      WL(3,0,4,-1) = (1.0E0_16,ZERO)
+      WL(4,0,3,-1) = (-1.0E0_16,ZERO)
+      DP_WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      DP_WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      DP_WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      DP_WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
+
       DO K=1, 3
         DO I=1,NCTAMPS
           AMPL(K,I)=(ZERO,ZERO)

=== modified file 'tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%mp_helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/MG5_1_MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/MG5_1_MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/born_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/born_matrix.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/born_matrix.f	2017-07-13 10:23:41 +0000
@@ -121,8 +121,9 @@
       DATA IDEN/36/
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/ML5_0_BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/ML5_0_BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/born_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/born_matrix.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/born_matrix.f	2017-07-13 10:23:41 +0000
@@ -137,8 +137,9 @@
       DATA IDEN/256/
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/ML5_0_BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/ML5_0_BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/born_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/born_matrix.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/born_matrix.f	2017-07-13 10:23:41 +0000
@@ -121,8 +121,9 @@
       DATA IDEN/36/
 
       INTEGER POLARIZATIONS(0:NEXTERNAL,0:5)
+      COMMON/ML5_0_BORN_BEAM_POL/POLARIZATIONS
       DATA ((POLARIZATIONS(I,J),I=0,NEXTERNAL),J=0,5)/NPOLENTRIES*-1/
-      COMMON/ML5_0_BORN_BEAM_POL/POLARIZATIONS
+
 C     
 C     FUNCTIONS
 C     

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/coef_construction_1.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/ML5_0_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/helas_calls_ampb_1.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/helas_calls_ampb_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/helas_calls_ampb_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/ML5_0_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/helas_calls_uvct_1.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/helas_calls_uvct_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/helas_calls_uvct_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/ML5_0_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_CT_calls_1.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_CT_calls_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_CT_calls_1.f	2017-07-13 10:23:41 +0000
@@ -76,9 +76,9 @@
       COMPLEX*16 W(20,NWAVEFUNCS)
       COMMON/ML5_0_W/W
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_matrix.f	2017-03-28 07:06:02 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_matrix.f	2017-07-13 10:23:41 +0000
@@ -352,9 +352,9 @@
       COMPLEX*32 MPW(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/MPW
 
-      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/WL,PL
 
       COMPLEX*16 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -583,12 +583,16 @@
 
 C       SETUP OF THE COMMON STARTING EXTERNAL LOOP WAVEFUNCTION
 C       IT IS ALSO PS POINT INDEPENDENT, SO IT CAN BE DONE HERE.
+C       The index -1 is for the charge-conjugated fermions with
+C        flipped fermion flow.
         DO I=0,3
+          PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDDO
         DO I=1,MAXLWFSIZE
           DO J=0,LOOPMAXCOEFS-1
             DO K=1,MAXLWFSIZE
+              WL(I,J,K,-1)=(0.0D0,0.0D0)
               IF(I.EQ.K.AND.J.EQ.0) THEN
                 WL(I,J,K,0)=(1.0D0,0.0D0)
               ELSE
@@ -603,6 +607,14 @@
         ENDIF
       ENDIF
 
+C     This is the chare conjugate version of the unit 4-currents in
+C      the canonical cartesian basis.
+C     This, for now, is only defined for 4-fermionic currents.
+      WL(1,0,2,-1) = DCMPLX(-1.0D0,0.0D0)
+      WL(2,0,1,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(3,0,4,-1) = DCMPLX(1.0D0,0.0D0)
+      WL(4,0,3,-1) = DCMPLX(-1.0D0,0.0D0)
+
 C     Make sure that lorentz rotation tests are not used if there is
 C      external loop wavefunction of spin 2 and that one specific
 C      helicity is asked

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/mp_coef_construction_1.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/mp_coef_construction_1.f	2016-04-06 18:45:52 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/mp_coef_construction_1.f	2017-07-13 10:23:41 +0000
@@ -60,9 +60,9 @@
       COMPLEX*32 W(20,NWAVEFUNCS)
       COMMON/ML5_0_MP_W/W
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_MP_WL/WL,PL
 
       COMPLEX*32 AMPL(3,NCTAMPS)

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/mp_compute_loop_coefs.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/mp_compute_loop_coefs.f	2016-09-15 20:20:30 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/mp_compute_loop_coefs.f	2017-07-13 10:23:41 +0000
@@ -118,14 +118,14 @@
       COMPLEX*16 DPW(20,NWAVEFUNCS)
       COMMON/ML5_0_W/DPW
 
-      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*32 PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*32 WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*32 PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_MP_WL/WL,PL
 
-      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE
-     $ ,0:NLOOPWAVEFUNCS)
-      COMPLEX*16 DP_PL(0:3,0:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_WL(MAXLWFSIZE,0:LOOPMAXCOEFS-1,MAXLWFSIZE,
+     $ -1:NLOOPWAVEFUNCS)
+      COMPLEX*16 DP_PL(0:3,-1:NLOOPWAVEFUNCS)
       COMMON/ML5_0_WL/DP_WL,DP_PL
 
       COMPLEX*32 LOOPCOEFS(0:LOOPMAXCOEFS-1,NSQUAREDSO,NLOOPGROUPS)
@@ -190,8 +190,10 @@
       ENDDO
 
       DO I=0,3
+        PL(I,-1)=CMPLX(ZERO,ZERO,KIND=16)
         PL(I,0)=CMPLX(ZERO,ZERO,KIND=16)
         IF (.NOT.COMPUTE_INTEGRAND_IN_QP) THEN
+          DP_PL(I,-1)=DCMPLX(0.0D0,0.0D0)
           DP_PL(I,0)=DCMPLX(0.0D0,0.0D0)
         ENDIF
       ENDDO
@@ -199,6 +201,8 @@
       DO I=1,MAXLWFSIZE
         DO J=0,LOOPMAXCOEFS-1
           DO K=1,MAXLWFSIZE
+            WL(I,J,K,-1)=(ZERO,ZERO)