Merge lp:~maddevelopers/mg5amcnlo/unleashed_reweighting into lp:~maddevelopers/mg5amcnlo/2.3.1

[Olivier Mattelaer]

This contains all the modifications for the LO reweighting including
1) the possibility to change model at LO
2) the possibility to change the process definition at LO
3) the possibility to use Loop-induced definition for the process
4) the kamikaze reweighting

This introduces a new way to interface python and fortran which do not require to have a lot of running executables/ open files for communication.
This method is based on the f2py program which is part of numpy.
Therefore you need to check that you have numpy install (on mac easy-install numpy) in order to run this package.

Cheers,

Olivier

[Olivier Mattelaer (olivier-mattelaer)]

review: Abstain

[Valentin Hirschi (valentin-hirschi)]

So here are my comments on this really great functionality.

------------------------------------------
General comments:
------------------------------------------

A1) What does this question really imply when launching the process?

      5 Add weight to events based on coupling parameters reweight=OFF

I thought that you were multiplying/dividing by the matrix element of the original and new hypothesis, and not just the ratio of particular couplings (which would be faster but much less general). So is this some hybrid implementation? I see the xml tag '<mgrwt>' added when turning this on, what is it used for?

A2) What is the status of the parrallelization? In the code it says that it's unstable, but is it still the case with f2py?

A3) What about adding the possibility of an 'unweight <filename>' command to the reweight_card.dat to allow the user to decide to create a new event file with only the unweighted events corresponding to the new hypothesis?

A4) How is the indicative error computed?

------------------------------------------
Comments on my trials of the reweighting:
------------------------------------------

B1) Debug left at [common_run_interface.py at line 3089]

B2) INFO: storring files of Previous run -> INFO: Storing files of previous run
    INFO: Do remember that the reweighting -> INFO: Remember that the reweighting

B3) I tried to reweight h > g g from HEFT with loop-induced, using:

change model loop_sm
change process g g > h [virt=QCD]
launch
/Users/valentin/Documents/Work/MG5/unleashed_reweighting/TEST_LI_gg_hg/Cards/param_card2.dat

and it crashed in a weird way telling me that the param_card I used was identical to the older one.
raise self.InvalidCmd, 'original card and new card are identical'.
Notice that you don't want to crash if they are identical but the model/process changed of course, so make sure that this case is correctly handled by the code.
Even then, I would put a warning and let the code go on with the same card because it is a good check to make sure that with an identical card one recovers the same result, and users might want to do this to get confidence in the tool.

Anyway, I realized that the problem in this case was that I should have used:

change model loop_sm
change process g g > h [sqrvirt=QCD]
launch
/Users/valentin/Documents/Work/MG5/LoopInduced_newrefine/PROC_loop_sm_2/Cards/param_card.dat

which makes more sense. But even then it crashed because the file:
"template_files/loop_optimized/check_py.f"
was not added to the revision. So I'll continue the test of the loop functionality once this file has been added.
Anyway, I think that the fact that one should preferably use the syntax 'sqrvirt=QCD' and also that launch should come after the 'change ...' is a bit tricky to guess, so you should add an explicit example of loop-induced reweighting in the header of the card reweight_card.dat.

B4) I tried NLO reweighting for p p > e+ ve [QCD] by changing the irrelevant top Yukawa with
    set ymt 200.0
and it first crashed because my f2py is called f2py-2.7. The crash was not helpful to the lambda user:

MadGraph5Error : A compilation Error occurs when trying to compile...

review: Needs Fixing

[lp:~maddevelopers/mg5amcnlo/unleashed_reweighting]

327. By Olivier Mattelaer on 2015-07-21

Improve the code following the suggestion of Valentin

[Olivier Mattelaer (olivier-mattelaer)]

HI Valentin,

Thanks a lot for your comment.
That’s what I call a good review.

> A1) What does this question really imply when launching the process?
>
> 5 Add weight to events based on coupling parameters reweight=OFF
>
> I thought that you were multiplying/dividing by the matrix element of the original and new hypothesis, and not just the ratio of particular couplings (which would be faster but much less general). So is this some hybrid implementation?

You are correct this is multiplying/dividing by the matrix element of the original and new hypothesis.
This nomenclature is mainly to distinguish from SysCalc type of reweighing.

> I see the xml tag ‘<mgrwt>' added when turning this on, what is it used for?

je ne vois pas ce tag pour des events NLO. @LO, c’est le tag pour SysCalc. Donc rien a voir avec cette branche.

> A2) What is the status of the parrallelization? In the code it says that it’s unstable, but is it still the case with f2py?

The parrallelization needs to be handle in a completely different way for f2py. The first problem is the JIL that I’m not sure how he is going to handle f2py.
if it consider it as a standard python function, then we are basically in trouble. So it is not going to work.

> A3) What about adding the possibility of an ‘unweight <filename>' command to the reweight_card.dat to allow the user to decide to create a new event file with only the unweighted events corresponding to the new hypothesis?

The options exists to create a new weighted sample. (the command is “change output 2")

> A4) How is the indicative error computed?

This has two part:
1) the original statistical error (rescale by the ratio of the cross-section).
2) a component proportional to the variance of the weight factor (variance/math.sqrt(event_nb) * orig_cross)
I combine those linearly.

> B1) Debug left at [common_run_interface.py at line 3089]

removed

> B2) INFO: storring files of Previous run -> INFO: Storing files of previous run
> INFO: Do remember that the reweighting -> INFO: Remember that the reweighting

thanks.

> change model loop_sm
> change process g g > h [virt=QCD]
> launch
> /Users/valentin/Documents/Work/MG5/unleashed_reweighting/TEST_LI_gg_hg/Cards/param_card2.dat
>
> and it crashed in a weird way telling me that the param_card I used was identical to the older one.

for me it crash with the expected message:
Command "generate_events run_01" interrupted with error:
InvalidCmd : NLO processes can’t be reweight (for Loop induced reweighting use [sqrvirt =])

> But even then it crashed because the file:
> “template_files/loop_optimized/check_py.f"

file added

> Even then, I would put a warning and let the code go on with the same card because it is a good check to make sure that with an identical card one recovers the same result, and users might want to do this to get confidence in the tool.

Ok, good point.

> B4) I tried NLO reweighting for p p > e+ ve [QCD] by changing the irrelevant top Yukawa with
> set ymt 200.0
> and it first crashed because my f2py is called f2py-2.7. The crash was not helpful to the lambda user:

ok now I run with the first of those three:
f2py, f2p...

[lp:~maddevelopers/mg5amcnlo/unleashed_reweighting]

328. By Olivier Mattelaer on 2015-07-24

additional change after skype comment with Valentin

[Olivier Mattelaer (olivier-mattelaer)]

Hi Valentin,

I have added the two point that we speak about yesterday.

Cheers,

Olivier

[lp:~maddevelopers/mg5amcnlo/unleashed_reweighting]

329. By Olivier Mattelaer on 2015-07-24

merge with latest 2.3.2

330. By Olivier Mattelaer on 2015-07-27

remove border effect linked to f2py

331. By Olivier Mattelaer on 2015-07-31

fix small bugs

332. By Olivier Mattelaer on 2015-08-08

implement gridpack mode (similar to the MS one)

333. By Olivier Mattelaer on 2015-08-09

try to force the fortran compiler for f2py

334. By Olivier Mattelaer on 2015-08-10

add security for matrix3 if the simple replacement rule do not work

335. By Olivier Mattelaer on 2015-08-10

forcing the compilation to go trough

336. By Olivier Mattelaer on 2015-08-10

correct the default output type for [virt=QCD]

337. By Frederix <email address hidden> on 2015-08-11

changed the comments in the default reweight_card.dat

338. By Olivier Mattelaer on 2015-08-12

fixing the import of matrix3 on linux/change the line in the question interface

339. By Olivier Mattelaer on 2015-08-12

review of Rik

340. By Olivier Mattelaer on 2015-08-13

add launch automatically if found a path or 'set' command

341. By Olivier Mattelaer on 2015-08-14

allow f2py to be specify via the configuration file

342. By Olivier Mattelaer on 2015-08-14

make sure that if "-p -R" are use simulataneously it works

343. By Olivier Mattelaer on 2015-08-14

fix small bug in the question (option 5) was not supported anymore

344. By Olivier Mattelaer on 2015-08-14

fixing the output of ML standalone --problem with the compiler change

345. By Olivier Mattelaer on 2015-08-16

merge with latest 2.3.2

346. By Olivier Mattelaer on 2015-08-16

add entry in the update note

Unmerged revisions

346. By Olivier Mattelaer on 2015-08-16

add entry in the update note

345. By Olivier Mattelaer on 2015-08-16

merge with latest 2.3.2

344. By Olivier Mattelaer on 2015-08-14

fixing the output of ML standalone --problem with the compiler change

343. By Olivier Mattelaer on 2015-08-14

fix small bug in the question (option 5) was not supported anymore

342. By Olivier Mattelaer on 2015-08-14

make sure that if "-p -R" are use simulataneously it works

341. By Olivier Mattelaer on 2015-08-14

allow f2py to be specify via the configuration file

340. By Olivier Mattelaer on 2015-08-13

add launch automatically if found a path or 'set' command

339. By Olivier Mattelaer on 2015-08-12

review of Rik

338. By Olivier Mattelaer on 2015-08-12

fixing the import of matrix3 on linux/change the line in the question interface

337. By Frederix <email address hidden> on 2015-08-11

changed the comments in the default reweight_card.dat

=== renamed file 'Template/LO/Cards/reweight_card_default.dat' => 'Template/Common/Cards/reweight_card_default.dat'
--- Template/LO/Cards/reweight_card_default.dat	2013-11-15 12:29:09 +0000
+++ Template/Common/Cards/reweight_card_default.dat	2015-07-24 15:59:43 +0000
@@ -1,25 +1,40 @@
-#******************************************************************
-#                       Reweight Module                           *
-#******************************************************************
+#*************************************************************************
+#                       Reweight Module                                  *
+#*************************************************************************
+#
+# Reweighting by Leading order (Tree level / loop-induced)
+#   element
+# This is not (in general) NLO accurate
+#
+# Note: 
+#   1) the value of alphas will be used from the event
+#      so the value of the param_card is not taken into account.
+#   2) It is (in general) dangerous/wrong  to change a particle's mass.
+#
+# Possible options: 
+# You can enter one of the following line to customize the re-weighting 
+#        procedure
+#
+# change model NAME : use another model for the re-weighted matrix-element
+#        In that case you need to provide the path to a correct param_card
+# change process DEF [-add]: change the process by which you re-weighted. 
+#        All the initial/final state of the original process needs to be 
+#	 defined.
+# change helicity False: Make the reweighting by full matrix-element and 
+#        not by the one of the associate helicity
+#        
+# Example of (standard) code -computing two weights:
+#
+# launch                 ! tag to start one computation of weight
+#    set BLOCKNAME ID VALUE ! rule to modity the param_card
+#    set BLOCKNAME ID VALUE ! rule to modity the param_card
+# launch                 ! start to compute a second weight
+#    /home/Cards/param_card_2.dat ! you can also enter a path to a valid card
+#
+#*************************************************************************
+# ENTER YOUR COMMAND BELOW. 
+#*************************************************************************
+
 launch
-#* Use the set command to specify the new set of parameter 
-#* Or specify a path to a valid param_card/banner
-#* Example of valid command:
-#*     set aewm1 137
-#*     ~/param_card.dat
-#*
-#* Note: 
-#*   1) the value of alphas will be used from the event
-#*      so the value of the param_card is not taken into account.
-#*   2) It is dangerous to change a mass of any particle. 
-
-
-#* If you want to compute the weight for more than one hyppothesis
-#* you need first to uncomment the following line:
-# launch
-# and then use the set command to specify your parameter.
-# All modification will start from the ORIGINAL card not from the 
-# last define one.
-#* You can have as many weight as you want.
-
-
+#  SPECIFY A PATH OR USE THE SET COMMAND like
+# set sminputs 1 130 # modify 1/alpha_EW
\ No newline at end of file

=== modified file 'Template/NLO/MCatNLO/include/LHEFRead.h'
--- Template/NLO/MCatNLO/include/LHEFRead.h	2014-01-23 17:04:40 +0000
+++ Template/NLO/MCatNLO/include/LHEFRead.h	2015-07-24 15:59:43 +0000
@@ -16,7 +16,7 @@
   LHEF::Reader reader;
 
 
-  void lhef_read_wgtsinfo_(int &cwgtinfo_nn, char (cwgtinfo_weights_info[250][15])) {
+  void lhef_read_wgtsinfo_(int &cwgtinfo_nn, char (cwgtinfo_weights_info[350][15])) {
 
     // Read header of event file
     std::stringstream hss;
@@ -60,11 +60,27 @@
 	  ++cwgtinfo_nn;
 	}
       }
+
+      // Read the mg_reweighting block
+      if (hs.find("<weightgroup type='mg_reweighting'") != std::string::npos) {
+	while (true) {
+	  std::getline(hss,hs,'\n');
+	  if (hs.find("</weightgroup>") != std::string::npos) break;
+	  if (hs.find("<weight id") != std::string::npos) {
+	    std::string sRWGT = hs.substr(hs.find("weight id")+11,hs.length());
+	    sRWGT = sRWGT.substr(0,sRWGT.find("'>"));
+	    //store the reweight label
+	    sprintf(cwgtinfo_weights_info[cwgtinfo_nn], "%15s", sRWGT.c_str());
+	    ++cwgtinfo_nn;
+	  }
+	}
+      }
+
     }
   }
 
 
-  void lhef_read_wgts_(double (cwgt_ww[250])) {
+  void lhef_read_wgts_(double (cwgt_ww[350])) {
     
     // Read events
     if (reader.readEvent()) {
@@ -106,7 +122,8 @@
 	  std::getline(ss,s,'\n');
 	  if (s.find("</rwgt>") != std::string::npos) break;
 	  if (s.find("id=") != std::string::npos) {
-	    std::string sww = s.substr(s.find("id=")+11,s.length());
+	    int ioffs=s.find("'>")+2;
+	    std::string sww = s.substr(ioffs,s.length());
 	    sww = sww.substr(0,sww.find("</w")-1);
 	    cwgt_ww[iww] = atof(sww.c_str());
 	    ++iww;

=== modified file 'Template/NLO/MCatNLO/srcHerwig/madfks_hwlhin.f'
--- Template/NLO/MCatNLO/srcHerwig/madfks_hwlhin.f	2014-06-17 08:36:10 +0000
+++ Template/NLO/MCatNLO/srcHerwig/madfks_hwlhin.f	2015-07-24 15:59:43 +0000
@@ -24,6 +24,7 @@
       PARAMETER (IUNIT=61)
       CHARACTER*80 STRING
       CHARACTER*140 BUFF
+      character*20 cidwgt
       character*12 dummy12
       character*2 dummy2
       CHARACTER*9 CH1
@@ -35,9 +36,14 @@
       COMMON/CEVWGT_LH/EVWGT_LH
       include 'reweight0.inc'
       integer iww,max_weight
-      parameter (max_weight=maxscales*maxscales+maxpdfs+1)
+      integer maxRWGT
+      parameter (maxRWGT=100)
+      double precision wgtxsecRWGT(maxRWGT)
+      parameter (max_weight=maxscales*maxscales+maxpdfs+maxRWGT+1)
       double precision ww(max_weight)
       common/cww/ww
+      integer numRWGTS
+      common/cnrwgts/numRWGTS
 C
       IF (IERROR.NE.0) RETURN
 c
@@ -124,7 +130,8 @@
         elseif(jwgtinfo.eq.9)then
           if (numscales.eq.0 .and. numPDFpairs.eq.0) then
              write (*,*) 'event file not correct format'
-             call HWWARN('UPEVNT',502)
+             write(*,*)'FATAL ERROR #3 IN UPEVNT'
+             stop
           endif
           read(iunit,'(a)')string
           wgtref=XWGTUP/MQQ
@@ -160,6 +167,11 @@
              iww=iww+1
              ww(iww)=wgtxsecPDF(i)
           enddo
+          do i=1,numRWGTS
+             call read_rwgt_line_RWGT(iunit,cidwgt,wgtxsecRWGT(i))
+             iww=iww+1
+             ww(iww)=wgtxsecRWGT(i)
+          enddo
           if (numscales.eq.0) then
              wgtxsecmu(1,1)=wgtref
           endif
@@ -176,7 +188,8 @@
         buff_tlh=' '
       endif
       if(INDEX(STRING,'</event>').eq.0)then
-        CALL HWWARN('UPEVNT',501)
+        write(*,*)'FATAL ERROR #2 IN UPEVNT'
+        stop
       endif
 c Modify what follows to set scale of H or S events in a different way
 c$$$      IF(ISORH_LHE.EQ.2)THEN
@@ -230,10 +243,13 @@
       character*15 weights_info(max_weight)
       common/cwgtsinfo/weights_info
       double precision xmuR,xmuF
-      integer iPDF,i
+      integer iPDF,i,numRWGTS
+      common/cnrwgts/numRWGTS
+      character*20 sRWGT
 C
       numscales=0
       numPDFpairs=0
+      numRWGTS=0
       nwgt=1
       weights_info(nwgt)="central value  "
       IF (IERROR.NE.0) RETURN
@@ -277,6 +293,20 @@
             ENDDO
             nwgt=nwgt+numPDFpairs
             numPDFpairs=numPDFpairs/2
+         ELSEIF( INDEX(STRING,"<weightgroup type='mg_reweighting'").ne.0
+     $           .and.STRING(1:1).ne.'#') then
+            DO WHILE (.TRUE.)
+               READ(61,'(a)') STRING
+               if (INDEX(STRING,"<weight id").ne.0 .and.
+     $              STRING(1:1).ne.'#')then
+                  numRWGTS=numRWGTS+1
+                  read(string(index(string,"weight id")+11:index(string,"'>")-1),*)sRWGT
+                  write(weights_info(numscales**2+2*numPDFpairs+numRWGTS+1),113)sRWGT
+               endif
+               if (INDEX(STRING,"</weightgroup>").ne.0 .and.
+     $              STRING(1:1).ne.'#') exit
+            ENDDO
+            nwgt=nwgt+numRWGTS
          ELSEIF ( INDEX(STRING,'</header>').ne.0 .and.
      &        STRING(1:1).ne.'#' ) then
             EXIT
@@ -297,8 +327,10 @@
       enddo
  111  format(a4,f3.1,x,a4,f3.1)
  112  format(a4,i8,a3)
+ 113  format(a15)
       return
- 998  CALL HWWARN('UPINIT',500)
+ 998  write(*,*)'FATAL ERROR #2 IN UPINIT'
+      stop
  999  END
 
 
@@ -312,7 +344,10 @@
       DOUBLE PRECISION PP(5,*),P(5,20),P2(20),M2(20),SP(5),
      & TINY,FAC,ECM,DCM,EP,STEP,FRT,HWUSQR
       DATA TINY,NT/1D-9,20/
-      IF (NP.GT.20) CALL HWWARN('HWURSC',300+NP)
+      IF (NP.GT.20)THEN
+        write(*,*)'FATAL ERROR #1 IN HWURSC'
+        stop
+      ENDIF
 C--COMPUTE CM MOMENTUM
       CALL HWVZRO(4,SP)
       DO IP=1,NP
@@ -337,15 +372,21 @@
                DCM=DCM+P2(IP)/EP
             ENDIF
          ENDDO
-         IF (DCM.EQ.0D0) CALL HWWARN('HWURSC',390)
+         IF (DCM.EQ.0D0)THEN
+            write(*,*)'FATAL ERROR #2 IN HWURSC'
+            stop
+         ENDIF
          STEP=2D0*(ECM-SP(5))/DCM
          FAC=FAC-STEP
          IF (ABS(STEP).LT.TINY) GOTO 100
       ENDDO
 C--FAILED TO CONVERGE
-      CALL HWWARN('HWURSC',1)
+        write(*,*)'WARNING #1 IN HWURSC'
 C--CONVERGED: RESCALE 3-MOMENTA AND BOOST BACK 
- 100  IF (FAC.LT.0D0) CALL HWWARN('HWURSC',391)
+ 100  IF (FAC.LT.0D0)THEN
+         write(*,*)'FATAL ERROR #3 IN HWURSC'
+         stop
+      ENDIF
       FRT=SQRT(FAC)
       DO IP=1,NP
          CALL HWVSCA(3,FRT,P(1,IP),P(1,IP))
@@ -372,3 +413,26 @@
       read (buff(wgt_start:100),*) wgt
       return
       end
+
+
+      subroutine read_rwgt_line_RWGT(unit,cid,wgt)
+c read a line in the <rwgt> tag. The syntax should be
+c  <wgt id='1001'> 0.1234567e+01 </wgt>
+c The id should be exactly 4 digits long.
+      implicit none
+      integer unit,wgt_start,id_start,id_end
+      double precision wgt
+      character*100 buff
+      character*20 cid
+      read (unit,'(a)') buff
+c Use char() to make sure that the non-standard characters are compiler
+c independent (char(62)=">", char(61)="=", char(39)="'")
+      wgt_start=index(buff,CHAR(39)//CHAR(62))+2
+      id_start=index(buff,'id'//CHAR(61)//CHAR(39))+4
+      id_end=wgt_start-3
+
+      read (buff(id_start:id_end),*) cid
+      read (buff(wgt_start:100),*) wgt
+
+      return
+      end

=== modified file 'Template/NLO/MCatNLO/srcPythia/madfks_pylhin.f'
--- Template/NLO/MCatNLO/srcPythia/madfks_pylhin.f	2014-06-17 08:36:10 +0000
+++ Template/NLO/MCatNLO/srcPythia/madfks_pylhin.f	2015-07-24 15:59:43 +0000
@@ -23,6 +23,7 @@
       PARAMETER (IUNIT=61)
       CHARACTER*80 STRING
       CHARACTER*140 BUFF
+      character*20 cidwgt
       character*12 dummy12
       character*2 dummy2
       CHARACTER*9 CH1
@@ -34,9 +35,14 @@
       COMMON/CEVWGT_LH/EVWGT_LH
       include 'reweight0.inc'
       integer iww,max_weight
-      parameter (max_weight=maxscales*maxscales+maxpdfs+1)
+      integer maxRWGT
+      parameter (maxRWGT=100)
+      double precision wgtxsecRWGT(maxRWGT)
+      parameter (max_weight=maxscales*maxscales+maxpdfs+maxRWGT+1)
       double precision ww(max_weight)
       common/cww/ww
+      integer numRWGTS
+      common/cnrwgts/numRWGTS
 C
       COMMON/PYPARS/MSTP(200),PARP(200),MSTI(200),PARI(200)
       DOUBLE PRECISION EVWEIGHT
@@ -164,6 +170,11 @@
              iww=iww+1
              ww(iww)=wgtxsecPDF(i)
           enddo
+          do i=1,numRWGTS
+             call read_rwgt_line_RWGT(iunit,cidwgt,wgtxsecRWGT(i))
+             iww=iww+1
+             ww(iww)=wgtxsecRWGT(i)
+          enddo
           if (numscales.eq.0) then
              wgtxsecmu(1,1)=wgtref
           endif
@@ -235,10 +246,13 @@
       character*15 weights_info(max_weight)
       common/cwgtsinfo/weights_info
       double precision xmuR,xmuF
-      integer iPDF,i
+      integer iPDF,i,numRWGTS
+      common/cnrwgts/numRWGTS
+      character*20 sRWGT
 C
       numscales=0
       numPDFpairs=0
+      numRWGTS=0
       nwgt=1
       weights_info(nwgt)="central value  "
 C--SET UP INPUT FILES
@@ -281,6 +295,20 @@
             ENDDO
             nwgt=nwgt+numPDFpairs
             numPDFpairs=numPDFpairs/2
+         ELSEIF( INDEX(STRING,"<weightgroup type='mg_reweighting'").ne.0
+     $           .and.STRING(1:1).ne.'#') then
+            DO WHILE (.TRUE.)
+               READ(61,'(a)') STRING
+               if (INDEX(STRING,"<weight id").ne.0 .and.
+     $              STRING(1:1).ne.'#')then
+                  numRWGTS=numRWGTS+1
+                  read(string(index(string,"weight id")+11:index(string,"'>")-1),*)sRWGT
+                  write(weights_info(numscales**2+2*numPDFpairs+numRWGTS+1),113)sRWGT
+               endif
+               if (INDEX(STRING,"</weightgroup>").ne.0 .and.
+     $              STRING(1:1).ne.'#') exit
+            ENDDO
+            nwgt=nwgt+numRWGTS
          ELSEIF ( INDEX(STRING,'</header>').ne.0 .and.
      &        STRING(1:1).ne.'#' ) then
             EXIT
@@ -301,6 +329,7 @@
       enddo
  111  format(a4,f3.1,x,a4,f3.1)
  112  format(a4,i8,a3)
+ 113  format(a15)
       return
  998  write(*,*)'FATAL ERROR #2 IN UPINIT'
       stop
@@ -516,3 +545,26 @@
       read (buff(wgt_start:100),*) wgt
       return
       end
+
+
+      subroutine read_rwgt_line_RWGT(unit,cid,wgt)
+c read a line in the <rwgt> tag. The syntax should be
+c  <wgt id='1001'> 0.1234567e+01 </wgt>
+c The id should be exactly 4 digits long.
+      implicit none
+      integer unit,wgt_start,id_start,id_end
+      double precision wgt
+      character*100 buff
+      character*20 cid
+      read (unit,'(a)') buff
+c Use char() to make sure that the non-standard characters are compiler
+c independent (char(62)=">", char(61)="=", char(39)="'")
+      wgt_start=index(buff,CHAR(39)//CHAR(62))+2
+      id_start=index(buff,'id'//CHAR(61)//CHAR(39))+4
+      id_end=wgt_start-3
+
+      read (buff(id_start:id_end),*) cid
+      read (buff(wgt_start:100),*) wgt
+
+      return
+      end

=== modified file 'Template/NLO/MCatNLO/srcPythia8/Pythia8.cc'
--- Template/NLO/MCatNLO/srcPythia8/Pythia8.cc	2015-03-30 12:01:21 +0000
+++ Template/NLO/MCatNLO/srcPythia8/Pythia8.cc	2015-07-24 15:59:43 +0000
@@ -30,8 +30,8 @@
   Pythia pythia;
 
   int cwgtinfo_nn;
-  char cwgtinfo_weights_info[250][15];
-  double cwgt_ww[250];
+  char cwgtinfo_weights_info[350][15];
+  double cwgt_ww[350];
 
   string inputname="Pythia8.cmd",outputname="Pythia8.hep";
 

=== modified file 'Template/NLO/MCatNLO/srcPythia8/Pythia82.cc'
--- Template/NLO/MCatNLO/srcPythia8/Pythia82.cc	2015-03-30 12:01:21 +0000
+++ Template/NLO/MCatNLO/srcPythia8/Pythia82.cc	2015-07-24 15:59:43 +0000
@@ -31,8 +31,8 @@
   Pythia pythia;
 
   int cwgtinfo_nn;
-  char cwgtinfo_weights_info[250][15];
-  double cwgt_ww[250];
+  char cwgtinfo_weights_info[350][15];
+  double cwgt_ww[350];
 
   string inputname="Pythia8.cmd",outputname="Pythia8.hep";
 

=== modified file 'Template/loop_material/StandAlone/SubProcesses/makefile'
--- Template/loop_material/StandAlone/SubProcesses/makefile	2015-05-10 14:33:18 +0000
+++ Template/loop_material/StandAlone/SubProcesses/makefile	2015-07-24 15:59:43 +0000
@@ -11,6 +11,7 @@
 
 include $(ROOT)/Source/make_opts
 include $(ROOT)/SubProcesses/MadLoop_makefile_definitions
+SHELL = /bin/bash
 
 LIBDIR = $(ROOT)/lib/
 PROG   = check
@@ -80,4 +81,23 @@
 	mv libMadLoop.$(libext) $(MADLOOP_LIB)
 
 clean:
-	@rm -f *.o
+	@rm -f *.o *.so
+
+../$(OLP):
+	cd ..; make $(OLP); cd -
+	rm -f libMadLoop.dylib
+	ln -s ../libMadLoop.dylib
+
+libMadLoop.dylib: ../$(OLP)
+
+matrix2py.so: ../$(OLP)
+	touch __init__.py
+	@ if [ -f `which f2py` ]; then \
+		f2py -L. -lMadLoop -m matrix2py -c check_sa.f;\
+	  else \
+		if [ -f `which f2py-2.7` ]; then \
+			f2py-2.7 -L. -lMadLoop -m matrix2py -c check_sa.f;\
+		else \
+			f2py-2.6 -L. -lMadLoop -m matrix2py -c check_sa.f;\
+		fi; \
+	  fi
\ No newline at end of file

=== modified file 'UpdateNotes.txt'
--- UpdateNotes.txt	2015-06-28 14:44:22 +0000
+++ UpdateNotes.txt	2015-07-24 15:59:43 +0000
@@ -70,6 +70,15 @@
                 RF: For fNLO runs the virtuals were included twice in the setting of the integration grids. 
                     This was not leading to any bias in previous version of the code.
 
+2.2.3(XX/XX/XX) RF: If a NAN is found, the code now skips that PS point and continues.
+                OM: Fix a bug in MadWeight (correlated param_card was not creating the correct input file)
+                RF: For fNLO runs the virtuals were included twice in the setting of the integration grids. 
+		RF: When requiring more than 1M events for (N)LO+PS runs, do not go to higher precision than 0.001
+                    for the grids and cross section (can be overwritten with the req_acc run_card parameter).
+		RF: Make sure that reweight info (for PDF and scale uncertainties) also works for UNLOPS events.
+
+
+
 2.2.2(06/11/14) OM: Correct a bug in the integration grid (introduces in 2.1.2). This was biasing the cross-section of 
                     processes like a a > mu+ mu- in the Effective Photon Approximation by three order of magnitude.
                     For LHC processes no sizeable effect have been observe so far.

=== modified file 'madgraph/interface/amcatnlo_run_interface.py'
--- madgraph/interface/amcatnlo_run_interface.py	2015-06-24 16:56:55 +0000
+++ madgraph/interface/amcatnlo_run_interface.py	2015-07-24 15:59:43 +0000
@@ -1209,6 +1209,7 @@
 
         if not mode in ['LO', 'NLO']:
             assert evt_file == pjoin(self.me_dir,'Events', self.run_name, 'events.lhe'), '%s != %s' %(evt_file, pjoin(self.me_dir,'Events', self.run_name, 'events.lhe.gz'))
+            self.exec_cmd('reweight -from_cards', postcmd=False)
             self.exec_cmd('decay_events -from_cards', postcmd=False)
             evt_file = pjoin(self.me_dir,'Events', self.run_name, 'events.lhe')
         
@@ -1851,7 +1852,7 @@
                           '\n      Total cross-section: %(xsect)8.3e +- %(errt)6.1e pb' % \
                         self.cross_sect_dict
 
-                if int(self.run_card['nevents'])>=10000 and self.run_card['reweight_scale']:
+                if self.run_card['nevents']>=10000 and self.run_card['reweight_scale']:
                    message = message + \
                        ('\n      Ren. and fac. scale uncertainty: +%0.1f%% -%0.1f%%') % \
                        (scale_pdf_info['scale_upp'], scale_pdf_info['scale_low'])
@@ -2356,7 +2357,7 @@
         Event dir. Return the name of the event file created
         """
         scale_pdf_info={}
-        if (self.run_card['reweight_scale'] or self.run_card['reweight_PDF']):
+        if self.run_card['reweight_scale'] or self.run_card['reweight_PDF'] :
             scale_pdf_info = self.run_reweight(options['reweightonly'])
 
         self.update_status('Collecting events', level='parton', update_results=True)
@@ -4038,7 +4039,7 @@
         
         
         void = 'NOT INSTALLED'
-        switch_order = ['order', 'fixed_order', 'shower','madspin']
+        switch_order = ['order', 'fixed_order', 'shower','madspin', 'reweight']
         switch_default = {'order': 'NLO', 'fixed_order': 'OFF', 'shower': void,
                   'madspin': void}
         if not switch:
@@ -4050,16 +4051,19 @@
         allowed_switch_value = {'order': ['LO', 'NLO'],
                                 'fixed_order': default_switch,
                                 'shower': default_switch,
-                                'madspin': default_switch}
+                                'madspin': default_switch,
+                                'reweight': default_switch}
         
         description = {'order':  'Perturbative order of the calculation:',
                        'fixed_order': 'Fixed order (no event generation and no MC@[N]LO matching):',
                        'shower': 'Shower the generated events:',
-                       'madspin': 'Decay particles with the MadSpin module:' }
+                       'madspin': 'Decay particles with the MadSpin module:',
+                       'reweight': 'Add weight to events based on coupling parameters'}
 
         force_switch = {('shower', 'ON'): {'fixed_order': 'OFF'},
                        ('madspin', 'ON'): {'fixed_order':'OFF'},
-                       ('fixed_order', 'ON'): {'shower': 'OFF', 'madspin': 'OFF'}
+                       ('reweight', 'ON'): {'fixed_order':'OFF'},
+                       ('fixed_order', 'ON'): {'shower': 'OFF', 'madspin': 'OFF', 'reweight':'ON'}
                        }
         special_values = ['LO', 'NLO', 'aMC@NLO', 'aMC@LO', 'noshower', 'noshowerLO']
 
@@ -4085,11 +4089,25 @@
                 switch['madspin'] = 'ON'
             else:
                 switch['madspin'] = 'OFF'
+            if misc.which('f2py'):
+                if os.path.exists(pjoin(self.me_dir,'Cards','reweight_card.dat')):
+                    switch['reweight'] = 'ON'
+                else:
+                    switch['reweight'] = 'OFF'
+            else:
+                switch['reweight'] = 'Numpy python package not available.'
+
+        if not aMCatNLO or self.options['mg5_path']:
+            available_mode.append('5')
+            if os.path.exists(pjoin(self.me_dir,'Cards','reweight_card.dat')):
+                switch['reweight'] = 'ON'
+            else:
+                switch['reweight'] = 'OFF'
             
         answers = list(available_mode) + ['auto', 'done']
         alias = {}
         for id, key in enumerate(switch_order):
-            if switch[key] != void:
+            if switch[key] != void and switch[key] in allowed_switch_value[key]:
                 answers += ['%s=%s' % (key, s) for s in allowed_switch_value[key]]
                 #allow lower case for on/off
                 alias.update(dict(('%s=%s' % (key, s.lower()), '%s=%s' % (key, s))
@@ -4121,36 +4139,32 @@
                 return 
             elif answer in special_values:
                 logger.info('Enter mode value: Go to the related mode', '$MG:color:BLACK')
+                assign_switch('reweight', 'OFF')
+                assign_switch('madspin', 'OFF')
                 if answer == 'LO':
                     switch['order'] = 'LO'
                     switch['fixed_order'] = 'ON'
                     assign_switch('shower', 'OFF')
-                    assign_switch('madspin', 'OFF')
                 elif answer == 'NLO':
                     switch['order'] = 'NLO'
                     switch['fixed_order'] = 'ON'
                     assign_switch('shower', 'OFF')
-                    assign_switch('madspin', 'OFF')
                 elif answer == 'aMC@NLO':
                     switch['order'] = 'NLO'
                     switch['fixed_order'] = 'OFF'
                     assign_switch('shower', 'ON')
-                    assign_switch('madspin', 'OFF')
                 elif answer == 'aMC@LO':
                     switch['order'] = 'LO'
                     switch['fixed_order'] = 'OFF'
                     assign_switch('shower', 'ON')
-                    assign_switch('madspin', 'OFF')
                 elif answer == 'noshower':
                     switch['order'] = 'NLO'
                     switch['fixed_order'] = 'OFF'
-                    assign_switch('shower', 'OFF')
-                    assign_switch('madspin', 'OFF')                                                    
+                    assign_switch('shower', 'OFF')                                                  
                 elif answer == 'noshowerLO':
                     switch['order'] = 'LO'
                     switch['fixed_order'] = 'OFF'
                     assign_switch('shower', 'OFF')
-                    assign_switch('madspin', 'OFF')
                 if mode:
                     return
             return switch
@@ -4207,8 +4221,11 @@
             options['parton'] = True
             ignore = ['shower_card.dat', 'madspin_card.dat']
             cards.append('FO_analyse_card.dat')
-        elif switch['madspin'] == 'ON':
-            cards.append('madspin_card.dat')
+        else:
+            if switch['madspin'] == 'ON':
+                cards.append('madspin_card.dat')
+            if switch['reweight'] == 'ON':
+                cards.append('reweight_card.dat')
         if 'aMC@' in mode:
             cards.append('shower_card.dat')
         if mode == 'onlyshower':

=== modified file 'madgraph/interface/common_run_interface.py'
--- madgraph/interface/common_run_interface.py	2015-06-14 00:32:59 +0000
+++ madgraph/interface/common_run_interface.py	2015-07-24 15:59:43 +0000
@@ -774,11 +774,11 @@
            madweight_card.dat [MW]
         """
 
-        text = open(path).read(50000)
-        if text == '':
+        fulltext = open(path).read(50000)
+        if fulltext == '':
             logger.warning('File %s is empty' % path)
             return 'unknown'
-        text = re.findall('(<MGVersion>|ParticlePropagator|<mg5proccard>|CEN_max_tracker|#TRIGGER CARD|parameter set name|muon eta coverage|QES_over_ref|MSTP|b_stable|FO_ANALYSIS_FORMAT|MSTU|Begin Minpts|gridpack|ebeam1|block\s+mw_run|BLOCK|DECAY|launch|madspin|transfer_card\.dat|set)', text, re.I)
+        text = re.findall('(<MGVersion>|ParticlePropagator|<mg5proccard>|CEN_max_tracker|#TRIGGER CARD|parameter set name|muon eta coverage|QES_over_ref|MSTP|b_stable|FO_ANALYSIS_FORMAT|MSTU|Begin Minpts|gridpack|ebeam1|block\s+mw_run|BLOCK|DECAY|launch|madspin|transfer_card\.dat|set)', fulltext, re.I)
         text = [t.lower() for t in text]
         if '<mgversion>' in text or '<mg5proccard>' in text:
             return 'banner'
@@ -809,12 +809,19 @@
             return 'shower_card.dat'
         elif 'fo_analysis_format' in text:
             return 'FO_analyse_card.dat'
-        elif 'decay' in text and 'launch' in text and 'madspin' in text:
-            return 'madspin_card.dat'
-        elif 'launch' in text and 'set' in text:
-            return 'reweight_card.dat'
-        elif 'decay' in text and 'launch' in text:
-            return 'madspin_card.dat'
+        elif 'launch' in text:
+            # need to separate madspin/reweight.
+            # decay/set can be in both...
+            if 'madspin' in text:
+                return 'madspin_card.dat'
+            if 'decay' in text:
+                # need to check if this a line like "decay w+" or "set decay"
+                if re.search("(^|;)\s*decay", fulltext):
+                    return 'madspin_card.dat'
+                else:
+                    return 'reweight_card.dat'
+            else:
+                return 'reweight_card.dat'
         else:
             return 'unknown'
 
@@ -1010,6 +1017,54 @@
         """Dummy routine, to be overwritten by daughter classes"""
 
         pass
+    ############################################################################
+    def do_reweight(self, line):
+        """ Allow to reweight the events generated with a new choices of model
+            parameter.
+        """
+        
+        if '-from_cards' in line and not os.path.exists(pjoin(self.me_dir, 'Cards', 'reweight_card.dat')):
+            return
+        
+        # Check that MG5 directory is present .
+        if MADEVENT and not self.options['mg5_path']:
+            raise self.InvalidCmd, '''The module reweight requires that MG5 is installed on the system.
+            You can install it and set its path in ./Cards/me5_configuration.txt'''
+        elif MADEVENT:
+            sys.path.append(self.options['mg5_path'])
+        try:
+            import madgraph.interface.reweight_interface as reweight_interface
+        except ImportError:
+            raise self.ConfigurationError, '''Can\'t load Reweight module.
+            The variable mg5_path might not be correctly configured.'''
+        
+        self.to_store.append('event')
+        if not '-from_cards' in line:
+            self.keep_cards(['reweight_card.dat'])
+            self.ask_edit_cards(['reweight_card.dat'], 'fixed', plot=False)        
+
+        # forbid this function to create an empty item in results.
+        if self.results.current['cross'] == 0 and self.run_name:
+            self.results.delete_run(self.run_name, self.run_tag)
+
+        # load the name of the event file
+        args = self.split_arg(line) 
+        self.check_decay_events(args) 
+        # args now alway content the path to the valid files
+        reweight_cmd = reweight_interface.ReweightInterface(args[0])
+        reweight_cmd.mother = self
+        self.update_status('Running Reweight', level='madspin')
+        
+        
+        path = pjoin(self.me_dir, 'Cards', 'reweight_card.dat')
+        reweight_cmd.me_dir = self.me_dir
+        reweight_cmd.import_command_file(path)
+        
+        # re-define current run
+        try:
+            self.results.def_current(self.run_name, self.run_tag)
+        except Exception:
+            pass
 
     ############################################################################
     def do_pgs(self, line):
@@ -3022,7 +3077,7 @@
         # check if input is a file
         elif hasattr(self, 'do_%s' % args[0]):
             self.do_set(' '.join(args[1:]))
-        elif os.path.exists(line):
+        elif os.path.isfile(line):
             self.copy_file(line)
             self.value = 'repeat'
         elif os.path.exists(pjoin(self.me_dir, line)):
@@ -3114,6 +3169,7 @@
         return self.mother_interface.complete_compute_widths(*args,**opts)
 
 
+
     def help_asperge(self):
         """Help associated to the asperge command"""
         signal.alarm(0)
@@ -3191,7 +3247,7 @@
             logger.warning('Fail to determine the type of the file. Not copied')
         if card_name != 'banner':
             logger.info('copy %s as %s' % (path, card_name))
-            files.cp(path, pjoin(self.mother_interface.me_dir, 'Cards', card_name))
+            files.cp(path, pjoin(self.me_dir, 'Cards', card_name))
         elif card_name == 'banner':
             banner_mod.split_banner(path, self.mother_interface.me_dir, proc_card=False)
             logger.info('Splitting the banner in it\'s component')

=== modified file 'madgraph/interface/extended_cmd.py'
--- madgraph/interface/extended_cmd.py	2015-03-31 06:48:20 +0000
+++ madgraph/interface/extended_cmd.py	2015-07-24 15:59:43 +0000
@@ -705,7 +705,7 @@
             return self.check_answer_in_input_file(question_instance, default, path)
         elif path:
             line = os.path.expanduser(os.path.expandvars(line))
-            if os.path.exists(line):
+            if os.path.isfile(line):
                 return line
         # No valid answer provides
         if self.haspiping:

=== modified file 'madgraph/interface/loop_interface.py'
--- madgraph/interface/loop_interface.py	2015-05-22 17:46:32 +0000
+++ madgraph/interface/loop_interface.py	2015-07-24 15:59:43 +0000
@@ -25,6 +25,7 @@
 import madgraph
 from madgraph import MG4DIR, MG5DIR, MadGraph5Error
 import madgraph.interface.madgraph_interface as mg_interface
+import madgraph.interface.launch_ext_program as launch_ext
 import madgraph.core.base_objects as base_objects
 import madgraph.core.diagram_generation as diagram_generation
 import madgraph.loop.loop_diagram_generation as loop_diagram_generation
@@ -90,9 +91,10 @@
         
         mg_interface.MadGraphCmd.check_output(self,args)
         
-        if self._export_format not in ["standalone", "matchbox"]: 
-            self._export_format = 'standalone'
 
+        if self._export_format not in self.supported_ML_format:
+            raise self.InvalidCmd, "not supported format %s" % self._export_format
+        
     def check_launch(self, args, options):
         """ Further check that only valid options are given to the MadLoop
         default launcher."""
@@ -228,10 +230,10 @@
 #                                  "will come out with the next release.")                           
 
         if isinstance(proc, base_objects.ProcessDefinition) and mode.startswith('ML5'):
-            if proc.has_multiparticle_label():
-                raise self.InvalidCmd(
+            if proc.has_multiparticle_label():        
+                raise self.InvalidCmd, \
                   "When running ML5 standalone, multiparticle labels cannot be"+\
-                  " employed. Please use the FKS5 interface instead.")
+                  " employed."
         
         if proc['decay_chains']:
             raise self.InvalidCmd(
@@ -341,7 +343,9 @@
       " corrections with this model because it does not support Feynman gauge.")
 
 class LoopInterface(CheckLoop, CompleteLoop, HelpLoop, CommonLoopInterface):
-        
+          
+    supported_ML_format = ['standalone', 'standalone_rw', 'matchbox'] 
+    
     def __init__(self, mgme_dir = '', *completekey, **stdin):
         """ Special init tasks for the Loop Interface """
 
@@ -420,16 +424,16 @@
         aloha_original_quad_mode = aloha.mp_precision
         aloha.mp_precision = True
 
-        if self._export_format not in ['standalone', 'matchbox']:
-            raise self.InvalidCmd('ML5 only support standalone/matchbox as export format.')
+        if self._export_format not in self.supported_ML_format:
+            raise self.InvalidCmd('ML5 only support "%s" as export format.' % \
+                                  ''.join(self.supported_ML_format))
 
-        if not os.path.isdir(self._export_dir) and \
-           self._export_format in ['matrix']:
+        if not os.path.isdir(self._export_dir) and self._export_format in ['matrix']:
             raise self.InvalidCmd('Specified export directory %s does not exist.'\
                                                          %str(self._export_dir))
 
         if not force and not noclean and os.path.isdir(self._export_dir)\
-               and self._export_format in ['standalone']:
+               and self._export_format.startswith('standalone'):
             # Don't ask if user already specified force or noclean
             logger.info('INFO: directory %s already exists.' % self._export_dir)
             logger.info('If you continue this directory will be cleaned')
@@ -443,17 +447,22 @@
                     raise self.InvalidCmd('Could not remove directory %s.'\
                                                          %str(self._export_dir))
 
-        if self._export_format == 'standalone':
+        if self._export_format.startswith('standalone'):
             output_type = 'madloop'
         elif self._export_format == 'matchbox':
             output_type = 'madloop_matchbox'
-
+        
         self._curr_exporter = export_v4.ExportV4Factory(self, \
                                                  noclean, output_type=output_type)
 
         if self._export_format in ['standalone', 'matchbox']:
             self._curr_exporter.copy_v4template(modelname=self._curr_model.get('name'))
 
+        if self._export_format == "standalone_rw":
+            self._export_format = "standalone"
+            self._curr_exporter.copy_v4template(modelname=self._curr_model.get('name'))
+            self._export_format = "standalone_rw"
+
         # Reset _done_export, since we have new directory
         self._done_export = False
 
@@ -509,7 +518,7 @@
         calls = 0
 
         path = self._export_dir
-        if self._export_format in ['standalone','matchbox']:
+        if self._export_format in self.supported_ML_format:
             path = pjoin(path, 'SubProcesses')
             
         cpu_time1 = time.time()
@@ -519,8 +528,8 @@
                         self._curr_matrix_elements.get_matrix_elements()
 
         # Fortran MadGraph5_aMC@NLO Standalone
-        if self._export_format in ['standalone', 'matchbox']:
-            for me in matrix_elements:            
+        if self._export_format in self.supported_ML_format:
+            for me in matrix_elements:
                 calls = calls + \
                         self._curr_exporter.generate_subprocess_directory_v4(\
                             me, self._curr_fortran_model)
@@ -531,12 +540,16 @@
             if self.options['loop_optimized_output'] and len(matrix_elements)>1:
                 max_lwfspins = [m.get_max_loop_particle_spin() for m in \
                                                                 matrix_elements]
-                max_loop_vert_ranks = [me.get_max_loop_vertex_rank() for me in \
+                try:
+                    max_loop_vert_ranks = [me.get_max_loop_vertex_rank() for me in \
                                                                 matrix_elements]
-                if len(set(max_lwfspins))>1 or len(set(max_loop_vert_ranks))>1:
-                    self._curr_exporter.fix_coef_specs(max(max_lwfspins),\
+                except MadGraph5Error:
+                    pass
+                else:
+                    if len(set(max_lwfspins))>1 or len(set(max_loop_vert_ranks))>1:
+                        self._curr_exporter.fix_coef_specs(max(max_lwfspins),\
                                                        max(max_loop_vert_ranks))
-                    logger.warning('ML5 has just output processes which do not'+\
+                        logger.warning('ML5 has just output processes which do not'+\
                       ' share the same maximum loop wavefunction size or the '+\
                       ' same maximum loop vertex rank. This is potentially '+\
                       ' dangerous. Please prefer to output them separately.')
@@ -580,7 +593,7 @@
         """Copy necessary sources and output the ps representation of 
         the diagrams, if needed"""
 
-        if self._export_format in ['standalone','matchbox']:
+        if self._export_format in self.supported_ML_format:
             logger.info('Export UFO model to MG4 format')
             # wanted_lorentz are the lorentz structures which are
             # actually used in the wavefunctions and amplitudes in
@@ -599,7 +612,7 @@
                                            wanted_lorentz,
                                            wanted_couplings)
 
-        if self._export_format in ['standalone', 'matchbox']:
+        if self._export_format in self.supported_ML_format:
             self._curr_exporter.finalize_v4_directory( \
                                            self._curr_matrix_elements,
                                            self.history,
@@ -607,7 +620,7 @@
                                            online,
                                            self.options['fortran_compiler'])
 
-        if self._export_format in ['standalone','matchbox']:
+        if self._export_format in self.supported_ML_format:
             logger.info('Output to directory ' + self._export_dir + ' done.')
 
     def do_launch(self, line, *args,**opt):

=== modified file 'madgraph/interface/madevent_interface.py'
--- madgraph/interface/madevent_interface.py	2015-06-23 21:24:06 +0000
+++ madgraph/interface/madevent_interface.py	2015-07-24 15:59:43 +0000
@@ -3061,54 +3061,6 @@
         self.update_status('End Parton', level='parton', makehtml=False)
         devnull.close()
     
-    ############################################################################
-    def do_reweight(self, line):
-        """ Allow to reweight the events generated with a new choices of model
-            parameter.
-        """
-        
-        if '-from_cards' in line and not os.path.exists(pjoin(self.me_dir, 'Cards', 'reweight_card.dat')):
-            return
-        
-        # Check that MG5 directory is present .
-        if MADEVENT and not self.options['mg5_path']:
-            raise self.InvalidCmd, '''The module reweight requires that MG5 is installed on the system.
-            You can install it and set its path in ./Cards/me5_configuration.txt'''
-        elif MADEVENT:
-            sys.path.append(self.options['mg5_path'])
-        try:
-            import madgraph.interface.reweight_interface as reweight_interface
-        except ImportError:
-            raise self.ConfigurationError, '''Can\'t load Reweight module.
-            The variable mg5_path might not be correctly configured.'''
-        
-        self.to_store.append('event')
-        if not '-from_cards' in line:
-            self.keep_cards(['reweight_card.dat'])
-            self.ask_edit_cards(['reweight_card.dat'], 'fixed', plot=False)        
-
-        # forbid this function to create an empty item in results.
-        if self.results.current['cross'] == 0 and self.run_name:
-            self.results.delete_run(self.run_name, self.run_tag)
-
-        # load the name of the event file
-        args = self.split_arg(line) 
-        self.check_decay_events(args) 
-        # args now alway content the path to the valid files
-        reweight_cmd = reweight_interface.ReweightInterface(args[0])
-        reweight_cmd. mother = self
-        self.update_status('Running Reweight', level='madspin')
-        
-        
-        path = pjoin(self.me_dir, 'Cards', 'reweight_card.dat')
-        reweight_cmd.me_dir = self.me_dir
-        reweight_cmd.import_command_file(path)
-        
-        # re-define current run
-        try:
-            self.results.def_current(self.run_name, self.run_tag)
-        except Exception:
-            pass
         
     ############################################################################ 
     def do_create_gridpack(self, line):
@@ -4280,16 +4232,19 @@
                 switch['madspin'] = 'ON'
             else:
                 switch['madspin'] = 'OFF'
-            if os.path.exists(pjoin(self.me_dir,'Cards','reweight_card.dat')):
-                switch['reweight'] = 'ON'
-            else:
-                switch['reweight'] = 'OFF'
+            if misc.which('f2py'):
+                if os.path.exists(pjoin(self.me_dir,'Cards','reweight_card.dat')):
+                    switch['reweight'] = 'ON'
+                else:
+                    switch['reweight'] = 'OFF'
+            else: 
+                switch['reweight'] = 'Numpy python package not available.'
                  
 
 
         options = list(available_mode) + ['auto', 'done']
         for id, key in enumerate(switch_order):
-            if switch[key] != void:
+            if switch[key] not in [void, 'Numpy python package not available.']:
                 options += ['%s=%s' % (key, s) for s in ['ON','OFF']]
                 options.append(key)
         options.append('parton')    

=== modified file 'madgraph/interface/madgraph_interface.py'
--- madgraph/interface/madgraph_interface.py	2015-06-26 10:55:43 +0000
+++ madgraph/interface/madgraph_interface.py	2015-07-24 15:59:43 +0000
@@ -1057,7 +1057,8 @@
         if not args:
             if self._done_export:
                 mode = self.find_output_type(self._done_export[0])
-                if mode != self._done_export[1]:
+                
+                if not self._done_export[1].startswith(mode):
                     print mode, self._done_export[1]
                     raise self.InvalidCmd, \
                           '%s not valid directory for launch' % self._done_export[0]
@@ -1384,10 +1385,6 @@
             text = 'No processes generated. Please generate a process first.'
             raise self.InvalidCmd(text)
 
-
-
-
-
         if args and args[0][0] != '-':
             # This is a path
             path = args.pop(0)
@@ -1415,6 +1412,7 @@
                 self.get_default_path()
                 if '-noclean' not in args and os.path.exists(self._export_dir):
                     args.append('-noclean')
+                    
             else:
                 if self.options['pythia8_path']:
                     self._export_dir = self.options['pythia8_path']
@@ -1548,7 +1546,7 @@
                                     (self._curr_model['name'], i)
             auto_path = lambda i: pjoin(self.writing_dir,
                                                name_dir(i))
-        elif self._export_format == 'standalone':
+        elif self._export_format.startswith('standalone'):
             name_dir = lambda i: 'PROC_SA_%s_%s' % \
                                     (self._curr_model['name'], i)
             auto_path = lambda i: pjoin(self.writing_dir,

=== modified file 'madgraph/interface/reweight_interface.py'
--- madgraph/interface/reweight_interface.py	2015-03-25 01:19:54 +0000
+++ madgraph/interface/reweight_interface.py	2015-07-24 15:59:43 +0000
@@ -20,6 +20,7 @@
 import os
 import re
 import shutil
+import sys
 import tempfile
 import time
 import subprocess
@@ -38,6 +39,8 @@
 import madgraph.various.banner as banner
 import madgraph.various.lhe_parser as lhe_parser
 import madgraph.various.combine_plots as combine_plots
+import madgraph.various.cluster as cluster
+import madgraph.fks.fks_common as fks_common
 
 import models.import_ufo as import_ufo
 import models.check_param_card as check_param_card 
@@ -75,8 +78,12 @@
         
         self.events_file = None
         self.processes = {}
+        self.second_model = None
+        self.second_process = None
         self.mg5cmd = master_interface.MasterCmd()
         self.seed = None
+        self.output_type = "default"
+        self.helicity_reweighting = True
         
         if event_path:
             logger.info("Extracting the banner ...")
@@ -91,6 +98,10 @@
             
     def do_import(self, inputfile, allow_madspin=False):
         """import the event file"""
+
+        args = self.split_arg(inputfile)
+        if not args:
+            return self.InvalidCmd, 'import requires arguments'
         
         # change directory where to write the output
         self.options['curr_dir'] = os.path.realpath(os.path.dirname(inputfile))
@@ -121,6 +132,17 @@
             self.lhe_input.banner = open(value).read()
         self.banner = self.lhe_input.get_banner()
         
+        #get original cross-section/error
+        if 'init' not in self.banner:
+            raise self.InvalidCmd('Event file does not contain init information')
+        for line in self.banner['init'].split('\n'):
+                split = line.split()
+                if len(split) == 4:
+                    cross, error = float(split[0]), float(split[1])
+        self.orig_cross = (cross, error)
+        
+        
+        
         # Check the validity of the banner:
         if 'slha' not in self.banner:
             misc.sprint(self.banner)
@@ -137,8 +159,8 @@
         # load information
         process = self.banner.get_detail('proc_card', 'generate')
         if '[' in process:
-            msg = 'Reweighting options is valid only for LO events'
-            raise Exception, msg            
+            logger.warning("Remember that the reweighting is performed at Leading Order. NLO precision is not guarantee.")
+
         if not process:
             msg = 'Invalid proc_card information in the file (no generate line):\n %s' % self.banner['mg5proccard']
             raise Exception, msg
@@ -149,6 +171,49 @@
         logger.info("options: %s" % option)
 
 
+    def get_LO_definition_from_NLO(self, proc):
+        """return the LO definitions of the process corresponding to the born/real"""
+        
+        # split the line definition with the part before and after the NLO tag
+        process, order, final = re.split('\[\s*(.*)\s*\]', proc)
+        # add the part without any additional jet.
+        commandline="add process %s %s --no_warning=duplicate;" % (process, final)
+        if not order:
+            #NO NLO tag => nothing to do actually return input
+            return proc
+        elif not order.startswith(('virt=','loonly=')):
+            # OK this a standard NLO process
+            if '=' in order:
+                # get the type NLO QCD/QED/...
+                order = order.split('=',1)[1]
+            # define the list of particles that are needed for the radiation
+            pert = fks_common.find_pert_particles_interactions(self.model,
+                                           pert_order = order)['soft_particles']
+            commandline += "define pert_%s = %s;" % (order.replace(' ',''), ' '.join(map(str,pert)) )
+            
+            # check if we have to increase by one the born order
+            if '%s=' % order in process:
+                result=re.split(' ',process)
+                process=''
+                for r in result:
+                    if '%s=' % order in r:
+                        ior=re.split('=',r)
+                        r='QCD=%i' % (int(ior[1])+1)
+                    process=process+r+' '
+            #handle special tag $ | / @
+            result = re.split('([/$@]|\w+=\w+)', process, 1)                    
+            if len(result) ==3:
+                process, split, rest = result
+                commandline+="add process %s pert_%s %s%s %s --no_warning=duplicate;" % (process, order.replace(' ','') ,split, rest, final)
+            else:
+                commandline +='add process %s pert_%s %s --no_warning=duplicate;' % (process,order.replace(' ',''), final)
+        else:
+            #just return the input. since this Madloop.
+            return proc                                       
+        logger.info(commandline)
+        return commandline
+
+
     def check_events(self):
         """Check some basic property of the events file"""
         
@@ -202,35 +267,39 @@
                                     pjoin(*[a for a in args if \
                                                       a.endswith(os.path.sep)]))
     
-    def check_set(self, args):
-        """checking the validity of the set command"""
-        
-        if len(args) < 2:
-            raise self.InvalidCmd('set command requires at least two argument.')
-        
-        valid = ['max_weight','seed','curr_dir']
-        if args[0] not in self.options and args[0] not in valid:
-            raise self.InvalidCmd('Unknown options %s' % args[0])        
-    
-        if args[0] == 'max_weight':
-            try:
-                args[1] = float(args[1].replace('d','e'))
-            except ValueError:
-                raise self.InvalidCmd('second argument should be a real number.')
-        
-        elif args[0] == 'BW_effect':
-            if args[1] in [0, False,'.false.', 'F', 'f', 'False', 'no']:
-                args[1] = 0
-            elif args[1] in [1, True,'.true.', 'T', 't', 'True', 'yes']:
-                args[1] = 1
+    def help_change(self):
+        """help for change command"""
+    
+        print "change model X :use model X for the reweighting"
+        print "change process p p > e+ e-: use a new process for the reweighting"
+        print "change process p p > mu+ mu- --add : add one new process to existing ones"
+    
+    def do_change(self, line):
+        """allow to define a second model/processes"""
+        
+        args = self.split_arg(line)
+        if len(args)<2:
+            logger.critical("not enough argument (need at least two). Discard line")
+        if args[0] == "model":
+            self.second_model = " ".join(args[1:])
+            if self.has_standalone_dir:
+                raise Exception, "command %s can only be run before the first launch. please run multiple times the re-weighting module to achieve this. (please quit python between each run)"
+        elif args[0] == "process":
+            if self.has_standalone_dir:
+                raise Exception, "command %s can only be run before the first launch. please run multiple times the re-weighting module to achieve this.(please quit python between each run)"
+            if args[-1] == "--add":
+                self.second_process.append(" ".join(args[1:-1]))
             else:
-                raise self.InvalidCmd('second argument should be either T or F.')
+                self.second_process = [" ".join(args[1:])]
+        elif args[0] == "output":
+            if args[1] in ['default', '2.0', 'unweight']:
+                self.output_type = args[1]
+        elif args[0] == "helicity":
+            self.helicity_reweighting = banner.ConfigFile.format_variable(args[1], bool, "helicity")
+        else:
+            logger.critical("unknown option! %s.  Discard line." % args[0])
         
-        elif args[0] == 'curr_dir':
-            if not os.path.isdir(args[1]):
-                raise self.InvalidCmd('second argument should be a path to a existing directory')
-
-            
+             
     def check_launch(self, args):
         """check the validity of the launch command"""
         
@@ -258,23 +327,27 @@
         model_line = self.banner.get('proc_card', 'full_model_line')
 
         if not self.has_standalone_dir:
+            misc.sprint("recreate standalone directory")
             self.create_standalone_directory()        
+        
+        if self.second_model or self.second_process:
+            rw_dir = pjoin(self.me_dir, 'rw_me_second')
+        else:
+            rw_dir = pjoin(self.me_dir, 'rw_me')
             
-        ff = open(pjoin(self.me_dir, 'rw_me','Cards', 'param_card.dat'), 'w')
+        ff = open(pjoin(rw_dir,'Cards', 'param_card.dat'), 'w')
         ff.write(self.banner['slha'])
         ff.close()
         ff = open(pjoin(self.me_dir, 'rw_me','Cards', 'param_card_orig.dat'), 'w')
         ff.write(self.banner['slha'])
         ff.close()        
         cmd = common_run_interface.CommonRunCmd.ask_edit_card_static(cards=['param_card.dat'],
-                                   ask=self.ask, pwd=pjoin(self.me_dir,'rw_me'))
-
-        new_card = open(pjoin(self.me_dir, 'rw_me', 'Cards', 'param_card.dat')).read()        
+                                   ask=self.ask, pwd=rw_dir)
+        new_card = open(pjoin(rw_dir, 'Cards', 'param_card.dat')).read()        
         # check if "Auto" is present for a width parameter
         if "auto" in new_card.lower():            
             self.mother.check_param_card(pjoin(self.me_dir, 'rw_me', 'Cards', 'param_card.dat'))
-            new_card = open(pjoin(self.me_dir, 'rw_me', 'Cards', 'param_card.dat')).read() 
-        
+            new_card = open(pjoin(rw_dir, 'Cards', 'param_card.dat')).read() 
 
         # Find new tag in the banner and add information if needed
         if 'initrwgt' in self.banner:
@@ -282,7 +355,7 @@
                 blockpat = re.compile(r'''<weightgroup type=\'mg_reweighting\'\s*>(?P<text>.*?)</weightgroup>''', re.I+re.M+re.S)
                 before, content, after = blockpat.split(self.banner['initrwgt'])
                 header_rwgt_other = before + after
-                pattern = re.compile('<weight id=\'mg_reweight_(?P<id>\d+)\'>(?P<info>[^<>]*)</weight>', re.S+re.I+re.M)
+                pattern = re.compile('<weight id=\'mg_reweight_(?P<id>\d+)\'>(?P<info>[^<]*)</weight>', re.S+re.I+re.M)
                 mg_rwgt_info = pattern.findall(content)
                 maxid = 0
                 for i, diff in mg_rwgt_info:
@@ -299,21 +372,38 @@
             header_rwgt_other = ''
             mg_rwgt_info = []
             rewgtid = 1
-        
 
         
-        
         # add the reweighting in the banner information:
         #starts by computing the difference in the cards.
         s_orig = self.banner['slha']
         s_new = new_card
-        old_param = check_param_card.ParamCard(s_orig.splitlines())
-        new_param =  check_param_card.ParamCard(s_new.splitlines())
-        card_diff = old_param.create_diff(new_param)
-        if card_diff == '':
-            raise self.InvalidCmd, 'original card and new card are identical'
-        mg_rwgt_info.append((str(rewgtid), card_diff))
-        
+        if not self.second_model:
+            old_param = check_param_card.ParamCard(s_orig.splitlines())
+            new_param =  check_param_card.ParamCard(s_new.splitlines())
+            card_diff = old_param.create_diff(new_param)
+            if card_diff == '' and not self.second_process:
+                logger.warning(' REWEIGHTING: original card and new card are identical. Is this really the run that you expect?')
+                #raise self.InvalidCmd, 'original card and new card are identical'
+            try:
+                if old_param['sminputs'].get(3)- new_param['sminputs'].get(3) > 1e-3 * new_param['sminputs'].get(3):
+                    logger.warning("We found different value of alpha_s. Note that the value of alpha_s used is the one associate with the event and not the one from the cards.")
+            except Exception, error:
+                logger.debug("error in check of alphas: %s" % str(error))
+                pass #this is a security                
+            if not self.second_process:
+                mg_rwgt_info.append((str(rewgtid), card_diff))
+            else:
+                str_proc = "\n change process  ".join([""]+self.second_process)
+                mg_rwgt_info.append((str(rewgtid), str_proc + '\n'+ card_diff))
+        else:
+            str_info = "change model %s" % self.second_model
+            if self.second_process:
+                str_info += "\n change process  ".join([""]+self.second_process)
+            card_diff = str_info
+            str_info += '\n' + s_new
+            mg_rwgt_info.append((str(rewgtid), str_info))
+
         # re-create the banner.
         self.banner['initrwgt'] = header_rwgt_other
         self.banner['initrwgt'] += '\n<weightgroup type=\'mg_reweighting\'>\n'
@@ -335,15 +425,16 @@
         if self.mother:
             out_path = pjoin(self.mother.me_dir, 'Events', 'reweight.lhe')
             output2 = open(out_path, 'w')
+            lha_strategy = self.banner.get_lha_strategy() 
+            self.banner.set_lha_strategy(4*lha_strategy/abs(lha_strategy)) 
             self.banner.write(output2, close_tag=False)
+            self.banner.set_lha_strategy(lha_strategy)
             new_banner = banner.Banner(self.banner)
             if not hasattr(self, 'run_card'):
                 self.run_card = new_banner.charge_card('run_card')
             self.run_card['run_tag'] = 'reweight_%s' % rewgtid
             new_banner['slha'] = s_new   
             del new_banner['initrwgt']
-            #ensure that original banner is kept untouched
-            assert new_banner['slha'] != self.banner['slha']
             assert 'initrwgt' in self.banner 
             ff = open(pjoin(self.mother.me_dir,'Events',self.mother.run_name, '%s_%s_banner.txt' % \
                           (self.mother.run_name, self.run_card['run_tag'])),'w') 
@@ -354,29 +445,69 @@
         tag_name = 'mg_reweight_%s' % rewgtid
         start = time.time()
         cross = 0
+        ratio, ratio_square = 0, 0 # to compute the variance and associate error
         
         os.environ['GFORTRAN_UNBUFFERED_ALL'] = 'y'
         if self.lhe_input.closed:
+            misc.sprint("using", self.lhe_input.name)
             self.lhe_input = lhe_parser.EventFile(self.lhe_input.name)
+
+#        Multicore option not really stable -> not use it
+        nb_core = 1
+#        if nb_core >1:
+#            multicore = cluster.MultiCore(nb_core)
+
         self.lhe_input.seek(0)
         for event_nb,event in enumerate(self.lhe_input):
             #control logger
-            if (event_nb % max(int(10**int(math.log10(float(event_nb)+1))),1000)==0): 
+            if (event_nb % max(int(10**int(math.log10(float(event_nb)+1))),10)==0): 
                     running_time = misc.format_timer(time.time()-start)
                     logger.info('Event nb %s %s' % (event_nb, running_time))
             if (event_nb==10001): logger.info('reducing number of print status. Next status update in 10000 events')
 
-
-            weight = self.calculate_weight(event)
-            cross += weight
-            event.reweight_data[tag_name] = weight
-            #write this event with weight
-            output.write(str(event))
-            if self.mother:
-                event.wgt = weight
-                event.reweight_data = {}
-                output2.write(str(event))
-
+            if nb_core > 1:
+                #        Multicore option not really stable -> not use it
+                while 1:
+                    if multicore.queue.qsize() < 100 * nb_core:
+                        multicore.submit(self.write_reweighted_event, argument=[event, tag_name])
+                        break
+                    #else:
+                    #    time.sleep(0.001)
+                continue
+            else:
+                weight = self.calculate_weight(event)
+                cross += weight
+                ratio += weight/event.wgt
+                ratio_square += (weight/event.wgt)**2
+                if self.output_type == "default":
+                    event.reweight_data[tag_name] = weight
+                    #write this event with weight
+                    output.write(str(event))
+                    if self.mother:
+                        event.wgt = weight
+                        event.reweight_data = {}
+                        output2.write(str(event))
+    
+                else:
+                    event.wgt = weight
+                    event.reweight_data = {}
+                    if self.mother:
+                        output2.write(str(event))
+                    else:
+                        output.write(str(event))
+
+        if nb_core >1:
+            multicore.wait_time = 1
+            multicore.wait(rw_dir,lambda *x:misc.sprint(x))
+            for t in multicore.demons:
+                cross += t.local.cross
+                
+        # check normalisation of the events:
+        if 'event_norm' in self.run_card:
+            if self.run_card['event_norm'] == 'average':
+                cross /= event_nb+1
+
+                
         running_time = misc.format_timer(time.time()-start)
         logger.info('All event done  (nb_event: %s) %s' % (event_nb+1, running_time))        
         
@@ -393,7 +524,12 @@
                 results.add_run(run_name, self.run_card, current=True)
                 results.add_detail('nb_event', event_nb+1)
                 results.add_detail('cross', cross)
-                results.add_detail('error', 'nan')
+                event_nb +=1
+                variance = ratio_square/event_nb - (ratio/event_nb)**2
+                orig_cross, orig_error = self.orig_cross
+                error = variance/math.sqrt(event_nb) * orig_cross + ratio/event_nb * orig_error
+                
+                results.add_detail('error', error)
                 self.mother.create_plot(mode='reweight', event_path=output2.name,
                                         tag=self.run_card['run_tag'])
                 #modify the html output to add the original run
@@ -420,82 +556,177 @@
                 #self.run_card['run_tag'] = self.run_card['run_tag'][9:]
                 #self.mother.run_name = old_name
         self.lhe_input.close()
-        files.mv(output.name, self.lhe_input.name)
+        if not self.mother or self.output_type != "default" :
+            target = pjoin(self.mother.me_dir, 'Events', run_name, 'events.lhe')
+        else:
+            target = self.lhe_input.name
+        
+        if self.output_type == "default":
+            files.mv(output.name, target)
+        elif self.output_type == "unweight":
+            output2.close()
+            lhe = lhe_parser.EventFile(output2.name)
+            nb_event = lhe.unweight(target)
+            if self.mother and  hasattr(self.mother, 'results'):
+                results = self.mother.results
+                results.add_detail('nb_event', nb_event)
+                results.current.parton.append('lhe')
+                
+        else:
+            files.mv(output2.name, self.lhe_input.name)     
+            if self.mother and  hasattr(self.mother, 'results'):
+                results = self.mother.results
+                results.current.parton.append('lhe')   
+
         logger.info('Event %s have now the additional weight' % self.lhe_input.name)
-        logger.info('new cross-section is : %g pb' % cross)
+        logger.info('new cross-section is : %g pb (indicative error: %g pb)' % (cross,error))
         self.terminate_fortran_executables(new_card_only=True)
         #store result
-        self.all_cross_section[rewgtid] = cross
-        
-
-
-            
-    def calculate_weight(self, event):
+        self.all_cross_section[rewgtid] = (cross, error)
+        
+
+    def write_reweighted_event(self, event, tag_name, **opt):
+        """a function for running in multicore"""
+        
+        if not hasattr(opt['thread_space'], "calculator"):
+            opt['thread_space'].calculator = {}
+            opt['thread_space'].calculator_nbcall = {}
+            opt['thread_space'].cross = 0
+            opt['thread_space'].output = open( self.lhe_input.name +'rw.%s' % opt['thread_id'], 'w')
+            if self.mother:
+                out_path = pjoin(self.mother.me_dir, 'Events', 'reweight.lhe.%s' % opt['thread_id'])
+                opt['thread_space'].output2 = open(out_path, 'w')
+                
+        weight = self.calculate_weight(event, space=opt['thread_space'])
+        opt['thread_space'].cross += weight
+        if self.output_type == "default":
+            event.reweight_data[tag_name] = weight
+            #write this event with weight
+            opt['thread_space'].output.write(str(event))
+            if self.mother:
+                event.wgt = weight
+                event.reweight_data = {}
+                opt['thread_space'].output2.write(str(event))
+        else:
+            event.wgt = weight
+            event.reweight_data = {}
+            if self.mother:
+                opt['thread_space'].output2.write(str(event))
+            else:
+                opt['thread_space'].output.write(str(event))
+        
+        return 0
+    
+    def calculate_weight(self, event, space=None):
+        """space defines where to find the calculator (in multicore)"""
+        
+        if not space: 
+            space = self
         
         event.parse_reweight()
-        w_orig = self.calculate_matrix_element(event, 0)
-        w_new =  self.calculate_matrix_element(event, 1)
         
+        w_orig = self.calculate_matrix_element(event, 0, space)
+        w_new =  self.calculate_matrix_element(event, 1, space)
+        if w_orig == 0:
+            tag, order = event.get_tag_and_order()
+            orig_order, Pdir, hel_dict = self.id_to_path[tag]
+            misc.sprint(w_orig, w_new)
+            misc.sprint(event)
+            raise Exception, "Invalid matrix element for original computation (weight=0)"
         return w_new/w_orig*event.wgt
-    
-    
-    def calculate_matrix_element(self, event, hypp_id):
+     
+    @staticmethod   
+    def invert_momenta(p):
+        """ fortran/C-python do not order table in the same order"""
+        new_p = []
+        for i in range(len(p[0])):  new_p.append([0]*len(p))
+        for i, onep in enumerate(p):
+            for j, x in enumerate(onep):
+                new_p[j][i] = x
+        return new_p
+    
+    def calculate_matrix_element(self, event, hypp_id, space):
         """routine to return the matrix element"""
 
         tag, order = event.get_tag_and_order()
-        orig_order, Pdir = self.id_to_path[tag]
+        
+        if (not self.second_model and not self.second_process) or hypp_id==0:
+            orig_order, Pdir, hel_dict = self.id_to_path[tag]
+        else:
+            orig_order, Pdir, hel_dict = self.id_to_path_second[tag]
 
         run_id = (tag, hypp_id)
 
-        if run_id in self.calculator:
-            external = self.calculator[run_id]
-            self.calculator_nbcall[run_id] += 1
-        else:
-            # create the executable for this param_card            
-
-            tmpdir = pjoin(self.me_dir,'rw_me', 'SubProcesses', Pdir)
-            executable_prod="./check"
-            if not os.path.exists(pjoin(tmpdir, 'check')):
-                misc.compile( cwd=tmpdir)
-            external = Popen(executable_prod, stdout=PIPE, stdin=PIPE, 
-                                                      stderr=STDOUT, cwd=tmpdir)
-            self.calculator[run_id] = external 
-            self.calculator_nbcall[run_id] = 1       
-            # set the param_card
+
+
+        start = False
+        if run_id in space.calculator:
+            external = space.calculator[run_id]
+        elif (not self.second_model and not self.second_process) or hypp_id==0:
+            # create the executable for this param_card  
+            subdir = pjoin(self.me_dir,'rw_me ', 'SubProcesses')
+            if self.me_dir not in sys.path:
+                sys.path.insert(0,self.me_dir)
+            Pname = os.path.basename(Pdir)
+            if hypp_id == 0:         
+                misc.compile(['matrix2py.so'], cwd=Pdir)
+                mymod = __import__('rw_me.SubProcesses.%s.matrix2py' % Pname, globals(), locals(), [],-1)
+                S = mymod.SubProcesses
+                P = getattr(S, Pname)
+                mymod = P.matrix2py
+                with misc.chdir(Pdir):
+                    mymod.initialise('param_card_orig.dat')
             if hypp_id == 1:
-                external.stdin.write('param_card.dat\n')
-            elif hypp_id == 0:
-                external.stdin.write('param_card_orig.dat\n')
-        #import the value of alphas
-        external.stdin.write('%g\n' % event.aqcd)
-        stdin_text = event.get_momenta_str(orig_order)
-        external.stdin.write(stdin_text)
-        me_value = external.stdout.readline()
-        try: 
-            me_value = float(me_value)
-        except Exception:
-            print 'ZERO DETECTED'
-            print stdin_text
-            print me_value
-            os.system('lsof -p %s' % external.pid)
-            me_value = 0
-        
-        if len(self.calculator) > 100:
-            logger.debug('more than 100 calculator. Perform cleaning')
-            nb_calls = self.calculator_nbcall.values()
-            nb_calls.sort()
-            cut = max([nb_calls[len(nb_calls)//2], 0.001 * nb_calls[-1]])
-            for key, external in list(self.calculator.items()):
-                nb = self.calculator_nbcall[key]
-                if nb < cut:
-                    external.stdin.close()
-                    external.stdout.close()
-                    external.terminate()
-                    del self.calculator[key]
-                    del self.calculator_nbcall[key]
-                else:
-                    self.calculator_nbcall[key] = self.calculator_nbcall[key] //10
-        
+                if not os.path.exists(pjoin(Pdir, 'matrix3py.so')):
+                    open(pjoin(Pdir, 'matrix3py.so'),'w').write(open(pjoin(Pdir, 'matrix2py.so')
+                                        ).read().replace('matrix2py', 'matrix3py'))
+                mymod = __import__('rw_me.SubProcesses.%s.matrix3py' % Pname, globals(), locals(), [],-1)
+                S = mymod.SubProcesses
+                P = getattr(S, Pname)
+                mymod = P.matrix3py
+                with misc.chdir(Pdir):
+                    mymod.initialise('param_card.dat') 
+            space.calculator[run_id] = mymod.get_me
+            external = space.calculator[run_id]                      
+        else:
+            subdir = pjoin(self.me_dir,'rw_me_second', 'SubProcesses')
+            if self.me_dir not in sys.path:
+                sys.path.append(self.me_dir)
+
+            assert hypp_id == 1
+            Pname = os.path.basename(Pdir)
+            misc.compile(['matrix2py.so'], cwd=pjoin(subdir, Pdir))
+            with misc.chdir(Pdir):
+                mymod = __import__("rw_me_second.SubProcesses.%s.matrix2py" % Pname)
+                reload(mymod)
+                S = mymod.SubProcesses
+                P = getattr(S, Pname)
+                mymod = P.matrix2py
+                mymod.initialise('param_card.dat')              
+            space.calculator[run_id] = mymod.get_me
+            external = space.calculator[run_id]                
+        
+        
+        p = self.invert_momenta(event.get_momenta(orig_order))
+        # add helicity information
+        
+        hel_order = event.get_helicity(orig_order)
+        if self.helicity_reweighting and 9 not in hel_order:
+            nhel = hel_dict[tuple(hel_order)]
+        else:
+            nhel = 0
+
+        with misc.chdir(Pdir):
+            with misc.stdchannel_redirected(sys.stdout, os.devnull):
+                me_value = external(p,event.aqcd, nhel)
+        # for NLO we have also the stability status code
+        if isinstance(me_value, tuple):
+            me_value, code = me_value
+            #if code points unstability -> returns 0
+            hundred_value = (code % 1000) //100
+            if hundred_value in [1,4]:
+                me_value = 0.
+            
         return me_value
     
     def terminate_fortran_executables(self, new_card_only=False):
@@ -504,11 +735,7 @@
         for (mode, production) in dict(self.calculator):
             
             if new_card_only and production == 0:
-                continue
-            external = self.calculator[(mode, production)]
-            external.stdin.close()
-            external.stdout.close()
-            external.terminate()            
+                continue            
             del self.calculator[(mode, production)]
     
     def do_quit(self, line):
@@ -526,7 +753,7 @@
         keys = self.all_cross_section.keys()
         keys.sort()
         for key in keys:
-            logger.info('%s : %s' % (key,self.all_cross_section[key]))  
+            logger.info('%s : %s +- %s pb' % (key,self.all_cross_section[key][0],self.all_cross_section[key][1] ))  
         self.terminate_fortran_executables()
             
     def __del__(self):
@@ -542,7 +769,7 @@
         
     
     @misc.mute_logger()
-    def create_standalone_directory(self):
+    def create_standalone_directory(self, second=False):
         """generate the various directory for the weight evaluation"""
         
         # 0. clean previous run ------------------------------------------------
@@ -599,9 +826,9 @@
         commandline=''
         for proc in processes:
             if '[' not in proc:
-                commandline+="add process %s ;" % proc
+                commandline += "add process %s ;" % proc
             else:
-                raise self.InvalidCmd('NLO processes can\'t be reweight')
+                commandline += self.get_LO_definition_from_NLO(proc)
         
         commandline = commandline.replace('add process', 'generate',1)
         logger.info(commandline)
@@ -610,8 +837,8 @@
         mgcmd.exec_cmd(commandline, precmd=True)        
         logger.info('Done %.4g' % (time.time()-start))
         self.has_standalone_dir = True
+        
 
-        
         # 3. Store id to directory information ---------------------------------
         matrix_elements = mgcmd._curr_matrix_elements.get_matrix_elements()
         
@@ -638,8 +865,142 @@
                         raise self.InvalidCmd, '2 different process have the same final states. This module can not handle such situation'
                     else:
                         continue
-                self.id_to_path[tag] = [order, Pdir]
-
+                # build the helicity dictionary
+                hel_nb = 0
+                hel_dict = {9:0} # unknown helicity -> use full ME
+                for helicities in me.get_helicity_matrix():
+                    hel_nb +=1 #fortran starts at 1
+                    hel_dict[tuple(helicities)] = hel_nb
+
+                self.id_to_path[tag] = [order, Pdir, hel_dict]        
+        
+        # 3. If we need a new model/process-------------------------------------
+        if self.second_model or self.second_process:
+            self.create_second_standalone_directory()    
+
+    @misc.mute_logger()
+    def create_second_standalone_directory(self, second=False):
+        """generate the various directory for the weight evaluation"""
+        
+        # 0. clean previous run ------------------------------------------------
+        path_me = self.me_dir
+        try:
+            shutil.rmtree(pjoin(path_me,'rw_me_second'))
+        except Exception: 
+            pass
+        
+        mgcmd = self.mg5cmd
+        # 1. Load model---------------------------------------------------------  
+        if self.second_model:
+            use_mgdefault= True
+            complex_mass = False
+            if ' ' in self.second_model:
+                args = self.second_model
+                if '--modelname' in args:
+                    use_mgdefault = False
+                model_name = args[0]
+            else:
+                model_name = self.second_model
+            self.load_model(model_name, use_mgdefault, complex_mass)
+                
+            modelpath = self.model.get('modelpath')
+            if os.path.basename(modelpath) != mgcmd._curr_model['name']:
+                name, restrict = mgcmd._curr_model['name'].rsplit('-',1)
+                if os.path.exists(pjoin(os.path.dirname(modelpath),name, 'restrict_%s.dat' % restrict)):
+                    modelpath = pjoin(os.path.dirname(modelpath), mgcmd._curr_model['name'])
+                
+            commandline="import model %s " % modelpath
+            mgcmd.exec_cmd(commandline)
+            
+        # 2. compute the production matrix element -----------------------------
+        if self.second_process:
+            processes = self.second_process
+        else:
+            processes = [line[9:].strip() for line in self.banner.proc_card
+                         if line.startswith('generate')]
+            processes += [' '.join(line.split()[2:]) for line in self.banner.proc_card
+                          if re.search('^\s*add\s+process', line)]   
+        mgcmd.exec_cmd("set group_subprocesses False")
+
+        logger.info('generating the square matrix element for reweighting')
+        start = time.time()
+        commandline=''
+        for proc in processes:
+            if '[' not in proc:
+                commandline+="add process %s ;" % proc
+            elif 'sqrvirt' in proc:
+                commandline+="add process %s ;" % proc
+            else:
+                raise self.InvalidCmd('NLO processes can\'t be reweight (for Loop induced reweighting use [sqrvirt =])')
+        
+        commandline = commandline.replace('add process', 'generate',1)
+        logger.info(commandline)
+        mgcmd.exec_cmd(commandline, precmd=True)
+        
+        matrix_elements = mgcmd._curr_matrix_elements.get_matrix_elements()
+        
+        commandline = 'output standalone_rw %s' % pjoin(path_me,'rw_me_second')
+        mgcmd.exec_cmd(commandline, precmd=True)        
+        logger.info('Done %.4g' % (time.time()-start))
+        
+        # 3. Store id to directory information ---------------------------------
+        matrix_elements = mgcmd._curr_matrix_elements.get_matrix_elements()
+        
+        self.id_to_path_second = {}
+        to_check = [] # list of tag that do not have a Pdir at creation time. 
+        for me in matrix_elements:
+            for proc in me.get('processes'):
+                initial = []    #filled in the next line
+                final = [l.get('id') for l in proc.get('legs')\
+                      if l.get('state') or initial.append(l.get('id'))]
+                order = (initial, final)
+                tag = proc.get_initial_final_ids()
+                decay_finals = proc.get_final_ids_after_decay()
+
+                if tag[1] != decay_finals:
+                    order = (initial, list(decay_finals))
+                    decay_finals.sort()
+                    tag = (tag[0], tuple(decay_finals))
+                Pdir = pjoin(path_me, 'rw_me_second', 'SubProcesses', 
+                                  'P%s' % me.get('processes')[0].shell_string())
+                if not os.path.exists(Pdir):
+                    to_check.append(tag)
+                    continue
+                if tag in self.id_to_path_second:
+                    if not Pdir == self.id_to_path_second[tag][1]:
+                        misc.sprint(tag, Pdir, self.id_to_path_second[tag][1])
+                        raise self.InvalidCmd, '2 different process have the same final states. This module can not handle such situation'
+                    else:
+                        continue
+                
+                # build the helicity dictionary
+                hel_nb = 0
+                hel_dict = {9:0} # unknown helicity -> use full ME
+                for helicities in me.get_helicity_matrix():
+                    hel_nb +=1 #fortran starts at 1
+                    hel_dict[tuple(helicities)] = hel_nb    
+                self.id_to_path_second[tag] = [order, Pdir, hel_dict]
+
+        for tag in to_check:
+            if tag not in self.id_to_path:
+                logger.warning("no valid path for %s" % (tag,))
+                #raise self.InvalidCmd, "no valid path for %s" % (tag,)
+        
+        # 4. Check MadLoopParam
+        if os.path.exists(pjoin(path_me, 'rw_me_second', 'Cards', 'MadLoopParams.dat')):
+            MLCard = banner.MadLoopParam(pjoin(path_me, 'rw_me_second', 'Cards', 'MadLoopParams.dat'))
+            MLCard.set('WriteOutFilters', False)
+            MLCard.set('UseLoopFilter', False)
+            MLCard.set("DoubleCheckHelicityFilter", False)
+            MLCard.set("HelicityFilterLevel", 0)
+            MLCard.write(pjoin(path_me, 'rw_me_second', 'SubProcesses', 'MadLoopParams.dat'),
+                         pjoin(path_me, 'rw_me_second', 'Cards', 'MadLoopParams.dat'), 
+                         commentdefault=False)
+            
+            
+            
+        
+        
 
     def load_model(self, name, use_mg_default, complex_mass=False):
         """load the model"""

=== modified file 'madgraph/iolibs/export_v4.py'
--- madgraph/iolibs/export_v4.py	2015-07-01 22:02:22 +0000
+++ madgraph/iolibs/export_v4.py	2015-07-24 15:59:43 +0000
@@ -1776,9 +1776,6 @@
         if self.format == 'standalone':
             shutil.copy(pjoin(self.mgme_dir, 'madgraph', 'iolibs', 'template_files', 'check_sa.f'), 
                     pjoin(self.dir_path, 'SubProcesses', 'check_sa.f'))
-        elif self.format == 'standalone_rw':
-            shutil.copy(pjoin(self.mgme_dir, 'madgraph', 'iolibs', 'template_files', 'driver_reweight.f'), 
-                    pjoin(self.dir_path, 'SubProcesses', 'check_sa.f'))
                         
         # Add file in Source
         shutil.copy(pjoin(temp_dir, 'Source', 'make_opts'), 
@@ -1795,8 +1792,16 @@
 
         super(ProcessExporterFortranSA,self).export_model_files(model_path)
         # Add the routine update_as_param in v4 model 
-        # This is a function created in the UFO 
-        
+        # This is a function created in the UFO  
+        text="""
+        subroutine update_as_param()
+          call setpara('param_card.dat',.false.)
+          return
+        end
+        """
+        ff = open(os.path.join(self.dir_path, 'Source', 'MODEL', 'couplings.f'),'a')
+        ff.write(text)
+        ff.close()        
         
         text = open(pjoin(self.dir_path,'SubProcesses','check_sa.f')).read()
         text = text.replace('call setpara(\'param_card.dat\')', 'call setpara(\'param_card.dat\', .true.)')
@@ -1836,6 +1841,9 @@
             history.write(output_file)
         
         ProcessExporterFortran.finalize_v4_directory(self, matrix_elements, history, makejpg, online, compiler)
+        open(pjoin(self.dir_path,'__init__.py'),'w')
+        open(pjoin(self.dir_path,'SubProcesses','__init__.py'),'w')
+
 
     def compiler_choice(self, compiler):
         """ Different daughter classes might want different compilers.
@@ -2119,7 +2127,7 @@
             matrix_template = 'matrix_standalone_msP_v4.inc'
         elif self.opt['export_format']=='standalone_msF':
             matrix_template = 'matrix_standalone_msF_v4.inc'
-	elif self.opt['export_format']=='matchbox':
+        elif self.opt['export_format']=='matchbox':
             replace_dict["proc_prefix"] = 'MG5_%i_' % matrix_element.get('processes')[0].get('id')
             replace_dict["color_information"] = self.get_color_string_lines(matrix_element)
 
@@ -5995,7 +6003,9 @@
       'fortran_compiler':cmd.options['fortran_compiler'],
       'output_dependencies':cmd.options['output_dependencies'],
       'SubProc_prefix':'P',
-      'compute_color_flows':cmd.options['loop_color_flows']}
+      'compute_color_flows':cmd.options['loop_color_flows'],
+      'mode': 'reweight' if cmd._export_format == "standalone_rw" else ''
+      }
 
     if output_type.startswith('madloop'):
         import madgraph.loop.loop_exporters as loop_exporters

=== modified file 'madgraph/iolibs/file_writers.py'
--- madgraph/iolibs/file_writers.py	2014-11-25 02:14:42 +0000
+++ madgraph/iolibs/file_writers.py	2015-07-24 15:59:43 +0000
@@ -192,7 +192,7 @@
     # Private variables
     __indent = 0
     __keyword_list = []
-    __comment_pattern = re.compile(r"^(\s*#|c$|(c\s+([^=]|$)))", re.IGNORECASE)
+    __comment_pattern = re.compile(r"^(\s*#|c$|(c\s+([^=]|$))|cf2py)", re.IGNORECASE)
 
     def write_line(self, line):
         """Write a fortran line, with correct indent and line splits"""
@@ -299,6 +299,10 @@
         # write_comment_line must have a single line as argument
         assert(isinstance(line, str) and line.find('\n') == -1)
 
+        if line.startswith('F2PY'):
+            return ["C%s\n" % line.strip()]
+        
+
         res_lines = []
 
         # This is a comment

=== removed file 'madgraph/iolibs/template_files/driver_reweight.f'
--- madgraph/iolibs/template_files/driver_reweight.f	2013-06-25 06:14:17 +0000
+++ madgraph/iolibs/template_files/driver_reweight.f	1970-01-01 00:00:00 +0000
@@ -1,103 +0,0 @@
-      PROGRAM DRIVER
-C**************************************************************************
-C     THIS IS THE DRIVER FOR CHECKING THE STANDALONE MATRIX ELEMENT.
-C     IT USES A SIMPLE PHASE SPACE GENERATOR
-C     Fabio Maltoni - 3rd Febraury 2007
-C     Modified Version for the decay Package (Pierre Artoisenet)
-C**************************************************************************
-      IMPLICIT NONE
-C     
-C     CONSTANTS  
-C     
-      REAL*8 ZERO
-      PARAMETER (ZERO=0D0)
-C     
-C     INCLUDE FILES
-C     
-C---  the include file with the values of the parameters and masses	
-      INCLUDE "coupl.inc"
-C---  integer nexternal ! number particles (incoming+outgoing) in the me 
-      INCLUDE "nexternal.inc" 
-C---  particle masses
-      REAL*8 PMASS(NEXTERNAL)	
-C---  integer    n_max_cg
-      INCLUDE "ngraphs.inc"     !how many diagrams (could be useful to know...)
-
-      DOUBLE PRECISION AMP2(n_max_cg)
-      COMMON/TO_AMPS/  AMP2
-
-
-C     
-C     LOCAL
-C     
-      INTEGER I,J,K
-      REAL*8 P(0:3,NEXTERNAL)   ! four momenta. Energy is the zeroth component.
-      REAL*8 SQRTS,MATELEM           ! sqrt(s)= center of mass energy 
-      REAL*8 PIN(0:3), POUT(0:3)
-      CHARACTER*120 BUFF(NEXTERNAL)
-      CHARACTER*200 param_card
-      real*8 new_alphas
-      real*8 pi
-C     
-C     EXTERNAL
-C     
-      REAL*8 DOT
-      EXTERNAL DOT
-      
-      pi = 3.141592653589793d0
-C-----
-C     BEGIN CODE
-C-----
-C     
-C---  INITIALIZATION CALLS
-C     
-c---  Call to initialize the values of the couplings, masses and widths 
-c     used in the evaluation of the matrix element. The primary parameters of the
-c     models are read from Cards/param_card.dat. The secondary parameters are calculated
-c     in Source/MODEL/couplings.f. The values are stored in common blocks that are listed
-c     in coupl.inc .
-      read(*,*) param_card
-      call setpara(param_card)  !first call to setup the paramaters
-
-c      include "../parameters.inc"
-      call coup()
-      include "pmass.inc"             !set up masses
-
-c      read phase-space point (main loop)
-1     read(*,*) new_alphas
-      G = 2* DSQRT(new_alphas*pi)
-      call UPDATE_AS_PARAM()
-      do i=1,nexternal
-         read (*,*) P(0,i),P(1,i),P(2,i),P(3,i) 
-      enddo
-
-
-
-c     
-c     Now we can call the matrix element!
-c
-
-      CALL SMATRIX(P,MATELEM)
-c
-
-c      write (*,*) "Matrix element = ", MATELEM, " GeV^",-(2*nexternal-8)	
-c      write (*,*) "-----------------------------------------------------------------------------"
-
-       write(*,*) MATELEM
-       call flush()
-       goto 1
-
-      end
-	
-	  
-	  
-	  
-	   double precision function dot(p1,p2)
-C****************************************************************************
-C     4-Vector Dot product
-C****************************************************************************
-      implicit none
-      double precision p1(0:3),p2(0:3)
-      dot=p1(0)*p2(0)-p1(1)*p2(1)-p1(2)*p2(2)-p1(3)*p2(3)
-      end
-

=== modified file 'madgraph/iolibs/template_files/loop/loop_matrix_standalone.inc'
--- madgraph/iolibs/template_files/loop/loop_matrix_standalone.inc	2015-03-10 10:10:05 +0000
+++ madgraph/iolibs/template_files/loop/loop_matrix_standalone.inc	2015-07-24 15:59:43 +0000
@@ -336,7 +336,7 @@
 109  CONTINUE
 CLOSE(1)
 IF(BOOTANDSTOP) THEN
-  WRITE(*,*) 'Stopped by user request.'
+  WRITE(*,*) '##Stopped by user request.'
   STOP
 ENDIF
 
@@ -622,8 +622,8 @@
   ELSEIF (.NOT.HELDOUBLECHECKED)THEN
     IF ((.NOT.GOODHEL(HELPICKED)).AND.(.NOT.%(proc_prefix)sISZERO(ABS(ANS(1))+ABS(ANS(2))+ABS(ANS(3)),REF/DBLE(NCOMB),-1))) THEN
 	  write(*,*) '##W15 Helicity filter could not be successfully double checked.'
-	  write(*,*) 'One reason for this is that you have changed sensible parameters which affected what are the zero helicity configurations.'
-	  write(*,*) 'MadLoop will try to reset the Helicity filter with the next PS points it receives.'
+	  write(*,*) '##One reason for this is that you have changed sensible parameters which affected what are the zero helicity configurations.'
+	  write(*,*) '##MadLoop will try to reset the Helicity filter with the next PS points it receives.'
 	  NTRY=0
 	  OPEN(30,FILE=HelFilterFN,err=349)
 349   CONTINUE
@@ -721,19 +721,19 @@
 	  NEPS=NEPS+1
       CALL %(proc_prefix)sCOMPUTE_ACCURACY(DP_RES,N_DP_EVAL,TEMP1,TEMP)	  
       WRITE(*,*) '##W03 WARNING An unstable PS point was',       ' detected.'
-	  WRITE(*,*) '(DP,QP) accuracies : (',TEMP1,',',ACC,')'
-	  WRITE(*,*) 'Best estimate (fin,1eps,2eps) :',(ANS(I),I=1,3)
+	  WRITE(*,*) '##(DP,QP) accuracies : (',TEMP1,',',ACC,')'
+	  WRITE(*,*) '##Best estimate (fin,1eps,2eps) :',(ANS(I),I=1,3)
 	  IF(NEPS.LE.10) THEN
-	    WRITE(*,*) 'Double precision evaluations :',(DP_RES(1,I),I=1,N_DP_EVAL)
-	    WRITE(*,*) 'Quad   precision evaluations :',(QP_RES(1,I),I=1,N_QP_EVAL)		 
-	    WRITE(*,*) 'PS point specification :'
-	    WRITE(*,*) 'Renormalization scale MU_R=',MU_R	
+	    WRITE(*,*) '##Double precision evaluations :',(DP_RES(1,I),I=1,N_DP_EVAL)
+	    WRITE(*,*) '##Quad   precision evaluations :',(QP_RES(1,I),I=1,N_QP_EVAL)		 
+	    WRITE(*,*) '##PS point specification :'
+	    WRITE(*,*) '##Renormalization scale MU_R=',MU_R	
 	    DO I=1,NEXTERNAL
           WRITE (*,'(i2,1x,4e27.17)') i, P(0,i),P(1,i),P(2,i),P(3,i) 
         ENDDO
 	  ENDIF
 	  IF(NEPS.EQ.10) THEN
-	    WRITE(*,*) 'Further output of the details of these unstable PS points will now be suppressed.'
+	    WRITE(*,*) '##Further output of the details of these unstable PS points will now be suppressed.'
 	  ENDIF
     ENDIF
   ELSE
@@ -883,7 +883,7 @@
 	  ENDIF
 
 	  IF(N_DP_EVALS.GT.20.OR.N_QP_EVALS.GT.20) THEN
-	    WRITE(*,*) 'ERROR:: Increase hardcoded maxstabilitylength.'
+	    WRITE(*,*) '##ERROR:: Increase hardcoded maxstabilitylength.'
 		STOP
 	  ENDIF
 
@@ -927,7 +927,7 @@
 C ----------
 C BEGIN CODE
 C ----------
-      write(*,*) 'WARNING:: Ignored, the possibility of selecting specific squared order contributions is not available in the default mode.'
+      write(*,*) '##WARNING:: Ignored, the possibility of selecting specific squared order contributions is not available in the default mode.'
 
 	  END
 

=== added file 'madgraph/iolibs/template_files/loop_optimized/check_py.f'
--- madgraph/iolibs/template_files/loop_optimized/check_py.f	1970-01-01 00:00:00 +0000
+++ madgraph/iolibs/template_files/loop_optimized/check_py.f	2015-07-24 15:59:43 +0000
@@ -0,0 +1,112 @@
+      Subroutine Initialise(path)
+
+      CHARACTER(128) path
+CF2PY intent(in):path
+ 
+C     INCLUDE FILES
+C     
+C     the include file with the values of the parameters and masses   
+C        
+      INCLUDE 'coupl.inc'
+      CALL SETPARA(path)
+      return 
+      end
+
+      SUBROUTINE GET_ME(P, ALPHAS, NHEL , ANS,RETURNCODE)
+      IMPLICIT NONE
+C     
+C     CONSTANTS  
+C     
+      REAL*8 ZERO
+      PARAMETER (ZERO=0D0)
+
+C     integer nexternal C number particles (incoming+outgoing) in the
+C      me 
+      include 'nexternal.inc'
+
+C      CHARACTER(512) MADLOOPRESOURCEPATH
+C     
+C     INCLUDE FILES
+C     
+C     the include file with the values of the parameters and masses   
+C        
+      INCLUDE 'coupl.inc'
+C     particle masses
+      REAL*8 PMASS(NEXTERNAL)
+C     integer    n_max_cg
+      INCLUDE 'ngraphs.inc'
+      INCLUDE 'nsquaredSO.inc'
+
+C     LOCAL
+C     
+      INTEGER I
+C     four momenta. Energy is the zeroth component.
+      REAL*8 P(0:3,NEXTERNAL)
+      INTEGER MATELEM_ARRAY_DIM
+      REAL*8 , ALLOCATABLE :: MATELEM(:,:)
+      INTEGER RETURNCODE
+      INTEGER NSQUAREDSO_LOOP
+      REAL*8 , ALLOCATABLE :: PREC_FOUND(:)
+
+      DOUBLE PRECISION ANS
+      INTEGER NHEL
+      DOUBLE PRECISION ALPHAS
+CF2PY INTENT(OUT) :: ANS
+CF2PY INTENT(OUT) :: RETURNCODE
+CF2PY INTENT(IN) :: NHEL   
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL) 
+CF2PY INTENT(IN) :: ALPHAS
+
+
+C     
+C     GLOBAL VARIABLES
+C     
+C     This is from ML code for the list of split orders selected by
+C     the process definition
+C     
+      INTEGER NLOOPCHOSEN
+      CHARACTER*20 CHOSEN_LOOP_SO_INDICES(NSQUAREDSO)
+      LOGICAL CHOSEN_LOOP_SO_CONFIGS(NSQUAREDSO)
+      COMMON/%(proc_prefix)sCHOSEN_LOOP_SQSO/CHOSEN_LOOP_SO_CONFIGS
+C     
+C     BEGIN CODE
+C     
+      CALL %(proc_prefix)sGET_ANSWER_DIMENSION(MATELEM_ARRAY_DIM)
+      ALLOCATE(MATELEM(0:3,0:MATELEM_ARRAY_DIM))
+      CALL %(proc_prefix)sGET_NSQSO_LOOP(NSQUAREDSO_LOOP)
+      ALLOCATE(PREC_FOUND(0:NSQUAREDSO_LOOP))
+      INCLUDE 'pmass.inc'
+
+C     Start by initializing what is the squared split orders indices
+C      chosen
+      NLOOPCHOSEN=0
+      DO I=1,NSQUAREDSO
+        IF (CHOSEN_LOOP_SO_CONFIGS(I)) THEN
+          NLOOPCHOSEN=NLOOPCHOSEN+1
+          WRITE(CHOSEN_LOOP_SO_INDICES(NLOOPCHOSEN),'(I3,A2)') I,'L)'
+        ENDIF
+      ENDDO
+
+C       Update the couplings with the new MU_R
+        CALL UPDATE_AS_PARAM2(2*P(0,1)*P(0,2), ALPHAS)
+
+C       
+C       Now we can call the matrix element
+C       
+        if (NHEL.eq.0) then
+          CALL %(proc_prefix)sSLOOPMATRIX_THRES(P,MATELEM,-1.0D0, PREC_FOUND, RETURNCODE)
+       else
+       CALL %(proc_prefix)sSLOOPMATRIXHEL_THRES(P,NHEL, MATELEM,-1.0D0, PREC_FOUND, RETURNCODE)
+       endif
+       
+c      loop induce -> only finite part 
+       ANS =  MATELEM(1,0)
+
+      END
+
+
+
+
+
+
+

=== modified file 'madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc'
--- madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc	2015-06-13 23:01:58 +0000
+++ madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc	2015-07-24 15:59:43 +0000
@@ -436,15 +436,15 @@
   ENDIF 
 ## if(LoopInduced){
   IF(.NOT.LoopInitStartOver) THEN
-    WRITE(*,*) 'INFO: For loop-induced processes it is preferable to always set the parameter LoopInitStartOver to True, so it is hard-set here to True.'
+    WRITE(*,*) '##INFO: For loop-induced processes it is preferable to always set the parameter LoopInitStartOver to True, so it is hard-set here to True.'
     LoopInitStartOver=.TRUE.
   ENDIF
   IF(.NOT.HelInitStartOver) THEN
-    WRITE(*,*) "INFO: For loop-induced processes it is preferable to always set the parameter HelInitStartOver to True, so it is hard-set here to True.'
+    WRITE(*,*) "##INFO: For loop-induced processes it is preferable to always set the parameter HelInitStartOver to True, so it is hard-set here to True.'
     HelInitStartOver=.TRUE.
   ENDIF
   IF (CheckCycle.LT.5) THEN
-    WRITE(*,*) "INFO: Due to the dynamic setting of the reference scale for contributions comparisons, it is preferable to set the parameter CheckCycle to a value larger than 4, so it is hard-set here to 5.'
+    WRITE(*,*) "##INFO: Due to the dynamic setting of the reference scale for contributions comparisons, it is preferable to set the parameter CheckCycle to a value larger than 4, so it is hard-set here to 5.'
     CheckCycle=5
   ENDIF
 ## }
@@ -452,7 +452,7 @@
 C Make sure that NROTATIONS_QP and NROTATIONS_DP are set to zero if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
   if(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
     IF(NROTATIONS_DP.NE.0.or.NROTATIONS_QP.NE.0) THEN
-      WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disabled, so MadLoop automatically resets NROTATIONS_DP and NROTATIONS_QP to 0.'
+      WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING is disabled, so MadLoop automatically resets NROTATIONS_DP and NROTATIONS_QP to 0.'
       NROTATIONS_QP=0
       NROTATIONS_DP=0
     ENDIF
@@ -542,7 +542,7 @@
 	ENDDO
   ENDDO
   IF(BOOTANDSTOP) THEN
-    WRITE(*,*) 'Stopped by user request.'
+    WRITE(*,*) '##Stopped by user request.'
     stop
   ENDIF
 ENDIF
@@ -1304,8 +1304,8 @@
     IF (GOODHEL(HELPICKED).EQ.-HELOFFSET) THEN
 	  IF (.NOT.%(proc_prefix)sIsZero(DABS(HELSAVED(1,HELPICKED))+DABS(HELSAVED(2,HELPICKED))+DABS(HELSAVED(2,HELPICKED)),REF/DBLE(NCOMB),-1,-1)) THEN
 	    write(*,*) '##W15 Helicity filter could not be successfully double checked.'
-	    write(*,*) 'One reason for this is that you might have changed sensible parameters which affected what are the zero helicity configurations.'
-	    write(*,*) 'MadLoop will try to reset the Helicity filter with the next PS points it receives.'
+	    write(*,*) '##One reason for this is that you might have changed sensible parameters which affected what are the zero helicity configurations.'
+	    write(*,*) '##MadLoop will try to reset the Helicity filter with the next PS points it receives.'
 	    NTRY=0
 	    OPEN(29,FILE=HelFilterFN,err=348)
 348     CONTINUE
@@ -1315,8 +1315,8 @@
 	IF (GOODHEL(HELPICKED).LT.-HELOFFSET.AND.NTRY.NE.0) THEN
 	  IF(%(proc_prefix)sISSAME(HELSAVED(1,HELPICKED),HELSAVED(1,ABS(GOODHEL(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
 	    write(*,*) '##W15 Helicity filter could not be successfully double checked.'
-	    write(*,*) 'One reason for this is that you might have changed sensible parameters which affected the helicity dependance relations.'
-	    write(*,*) 'MadLoop will try to reset the Helicity filter with the next PS points it receives.'
+	    write(*,*) '##One reason for this is that you might have changed sensible parameters which affected the helicity dependance relations.'
+	    write(*,*) '##MadLoop will try to reset the Helicity filter with the next PS points it receives.'
 	    NTRY=0
 	    OPEN(30,FILE=HelFilterFN,err=349)
 349     CONTINUE
@@ -1532,12 +1532,12 @@
 	  IF(NEPS.LE.10) THEN	  
         WRITE(*,*) '##W03 WARNING An unstable PS point was',       ' detected.'
 		IF (NSQUAREDSO.NE.1) THEN
-		  WRITE(*,*) 'Accuracies for each split order, starting with the summed case'
-	      WRITE(*,*) 'DP accuracies (for each split order): ',(TEMP1(I),I=0,NSQUAREDSO)
-	      WRITE(*,*) 'QP accuracies (for each split order): ',(ACC(I),I=0,NSQUAREDSO)
+		  WRITE(*,*) '##Accuracies for each split order, starting with the summed case'
+	      WRITE(*,*) '##DP accuracies (for each split order): ',(TEMP1(I),I=0,NSQUAREDSO)
+	      WRITE(*,*) '##QP accuracies (for each split order): ',(ACC(I),I=0,NSQUAREDSO)
 		ELSE
-	      WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-	      WRITE(*,*) 'QP accuracy: ',ACC(1)
+	      WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+	      WRITE(*,*) '##QP accuracy: ',ACC(1)
 		ENDIF
 	    DO J=0,NSQUAREDSO
 		  IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== modified file 'madgraph/iolibs/template_files/makefile_sa_f_sp'
--- madgraph/iolibs/template_files/makefile_sa_f_sp	2013-10-30 08:15:19 +0000
+++ madgraph/iolibs/template_files/makefile_sa_f_sp	2015-07-24 15:59:43 +0000
@@ -1,5 +1,5 @@
 include ../../Source/make_opts
-
+SHELL = /bin/bash
 LIBDIR = ../../lib/
 PROG   = check
 PROG_SPLITORDERS = check_sa_born_splitOrders
@@ -15,4 +15,17 @@
 $(PROG_SPLITORDERS): $(PROCESS) $(CHECK_SA_SPLITORDERS) makefile $(LIBS)
 	$(FC) $(FFLAGS) -o $(PROG) $(PROCESS) $(CHECK_SA_SPLITORDERS) $(LINKLIBS)
 
-driver.f: nexternal.inc pmass.inc ngraphs.inc coupl.inc
\ No newline at end of file
+driver.f: nexternal.inc pmass.inc ngraphs.inc coupl.inc
+
+# For python linking (require f2py part of numpy)
+matrix2py.so: matrix.f makefile $(LIBS)
+	touch __init__.py
+	@ if [ -f `which f2py` ]; then \
+		f2py $(LINKLIBS) -c  matrix.f -m matrix2py;\
+	  else \
+		if [ -f `which f2py-2.7` ]; then \
+			f2py-2.7 $(LINKLIBS) -c  matrix.f -m matrix2py;\
+		else \
+			f2py-2.6 $(LINKLIBS) -c  matrix.f -m matrix2py;\
+		fi; \
+	  fi
\ No newline at end of file

=== modified file 'madgraph/iolibs/template_files/matrix_standalone_v4.inc'
--- madgraph/iolibs/template_files/matrix_standalone_v4.inc	2013-12-07 23:38:41 +0000
+++ madgraph/iolibs/template_files/matrix_standalone_v4.inc	2015-07-24 15:59:43 +0000
@@ -7,6 +7,10 @@
       PARAMETER (NEXTERNAL=%(nexternal)d)
       INTEGER                 NCOMB         
       PARAMETER (             NCOMB=%(ncomb)d)
+CF2PY INTENT(OUT) :: ANS                                                                                                                                                                                     
+CF2PY INTENT(IN) :: HEL                                                                                                                                                                                    
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL)  
+
 C  
 C ARGUMENTS 
 C  
@@ -52,6 +56,8 @@
 C ARGUMENTS 
 C  
       REAL*8 P(0:3,NEXTERNAL),ANS
+CF2PY INTENT(OUT) :: ANS
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL)
 C  
 C LOCAL VARIABLES 
 C  
@@ -76,15 +82,15 @@
 C ----------
 C BEGIN CODE
 C ----------
-      NTRY=NTRY+1
+      IF(USERHEL.EQ.-1) NTRY=NTRY+1
       DO IHEL=1,NEXTERNAL
          JC(IHEL) = +1
       ENDDO
       ANS = 0D0
           DO IHEL=1,NCOMB
 	         IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN		  
-                IF (GOODHEL(IHEL) .OR. NTRY .LT. 20) THEN
-                    T=MATRIX(P ,NHEL(1,IHEL),JC(1))            
+                IF (GOODHEL(IHEL) .OR. NTRY .LT. 20.OR.USERHEL.NE.-1) THEN
+                    T=%(proc_prefix)sMATRIX(P ,NHEL(1,IHEL),JC(1))            
                     ANS=ANS+T
                     IF (T .NE. 0D0 .AND. .NOT.    GOODHEL(IHEL)) THEN
                         GOODHEL(IHEL)=.TRUE.
@@ -99,7 +105,7 @@
       END
        
        
-      REAL*8 FUNCTION MATRIX(P,NHEL,IC)
+      REAL*8 FUNCTION %(proc_prefix)sMATRIX(P,NHEL,IC)
 C  
 %(info_lines)s
 C
@@ -162,3 +168,46 @@
       ENDDO
 %(amp2_lines)s
       END
+
+      SUBROUTINE %(proc_prefix)sGET_ME(P, ALPHAS, NHEL ,ANS)
+      IMPLICIT NONE   
+C
+C CONSTANT
+C
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=%(nexternal)d)
+C  
+C ARGUMENTS 
+C  
+      REAL*8 P(0:3,NEXTERNAL),ANS
+      INTEGER NHEL
+      DOUBLE PRECISION ALPHAS 
+      real*8 pi     
+CF2PY INTENT(OUT) :: ANS  
+CF2PY INTENT(IN) :: NHEL   
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL) 
+CF2PY INTENT(IN) :: ALPHAS
+C     ROUTINE FOR F2PY to read the benchmark point.    
+C     the include file with the values of the parameters and masses 
+      include "coupl.inc"
+      
+      pi = 3.141592653589793d0
+      G = 2* DSQRT(ALPHAS*pi)
+      call UPDATE_AS_PARAM()
+      if (NHEL.ne.0) then
+        CALL SMATRIXHEL(P, NHEL, ANS)
+      else
+        CALL SMATRIX(P, ANS)
+      endif
+      return 
+      end
+
+      SUBROUTINE %(proc_prefix)sINITIALISE(PATH)
+C     ROUTINE FOR F2PY to read the benchmark point.    
+      IMPLICIT NONE   
+      CHARACTER*180 PATH
+CF2PY INTENT(IN) :: PATH 
+      call setpara(PATH)  !first call to setup the paramaters    
+      return 
+      end      
+      
\ No newline at end of file

=== modified file 'madgraph/loop/loop_exporters.py'
--- madgraph/loop/loop_exporters.py	2015-06-27 00:16:19 +0000
+++ madgraph/loop/loop_exporters.py	2015-07-24 15:59:43 +0000
@@ -82,13 +82,15 @@
                         'fortran_compiler':'gfortran',
                         'SubProc_prefix': 'P',
                         'output_dependencies': 'external',
-                        'compute_color_flows': False}
+                        'compute_color_flows': False,
+                        'mode':''}
 
 
     def __init__(self, mgme_dir="", dir_path = "", opt=None):
         """Initiate the LoopExporterFortran with directory information on where
         to find all the loop-related source files, like CutTools"""
 
+
         self.opt = dict(self.default_opt)
         if opt:
             self.opt.update(opt)
@@ -977,9 +979,12 @@
                 replace_dict[key]=''
         if matrix_element.get('processes')[0].get('has_born'):
             file = open(os.path.join(self.template_dir,'check_sa.inc')).read()
+        elif self.opt['mode'] == 'reweight':
+            file = open(os.path.join(self.template_dir,\
+                                          'check_py.f')).read()            
         else:
             file = open(os.path.join(self.template_dir,\
-                                          'check_sa_loop_induced.inc')).read()            
+                                          'check_sa_loop_induced.inc')).read()
         file=file%replace_dict
         writer.writelines(file)
 

=== modified file 'madgraph/madevent/gen_crossxhtml.py'
--- madgraph/madevent/gen_crossxhtml.py	2015-05-22 16:18:38 +0000
+++ madgraph/madevent/gen_crossxhtml.py	2015-07-24 15:59:43 +0000
@@ -30,11 +30,13 @@
     import internal.save_load_object as save_load_object
     import internal.lhe_parser as lhe_parser
     import internal.misc as misc
+    import internal.banner as bannerlib
 else:
     import madgraph.iolibs.files as files
     import madgraph.iolibs.save_load_object as save_load_object
     import madgraph.various.lhe_parser as lhe_parser
     import madgraph.various.misc as misc
+    import madgraph.various.banner as bannerlib
 
 pjoin = os.path.join
 exists = os.path.exists
@@ -160,6 +162,7 @@
         # Check if some directory already exists and if so add them
         runs = [d for d in os.listdir(pjoin(path, 'Events')) if 
                                       os.path.isdir(pjoin(path, 'Events', d))]
+
         if runs:
             if recreateold:
                 for run in runs:
@@ -174,6 +177,7 @@
         
         event_path = pjoin(self.path, "Events", run_name, "unweighted_events.lhe")
         
+
         try:
             import madgraph
         except ImportError:
@@ -625,8 +629,8 @@
      
         # No pythia only a single run:
         if not tags:
-            self[0]['nb_event'] = nb_event
-            self[0]['cross'] = cross
+            self[-1]['nb_event'] = nb_event
+            self[-1]['cross'] = cross
           
         #Loop over pythia run
         for tag in tags:

=== modified file 'madgraph/madevent/sum_html.py'
--- madgraph/madevent/sum_html.py	2015-06-14 00:32:59 +0000
+++ madgraph/madevent/sum_html.py	2015-07-24 15:59:43 +0000
@@ -287,6 +287,7 @@
                         if m:
                             return float(m.group(1))*10**(float(m.group(2)))
                         return 
+
                 data = [secure_float(d) for d in line.split()]
                 self.axsec, self.xerru, self.xerrc, self.nevents, self.nw,\
                          self.maxit, self.nunwgt, self.luminosity, self.wgt, \

=== modified file 'madgraph/various/banner.py'
--- madgraph/various/banner.py	2015-06-18 16:49:50 +0000
+++ madgraph/various/banner.py	2015-07-24 15:59:43 +0000
@@ -292,7 +292,62 @@
                 if pid not in pid2label.keys(): 
                     block.remove((pid,))
 
-
+    def get_lha_strategy(self):
+        """get the lha_strategy: how the weight have to be handle by the shower"""
+        
+        if not self["init"]:
+            raise Exception, "No init block define"
+        
+        data = self["init"].split('\n')[0].split()
+        if len(data) != 10:
+            misc.sprint(len(data), self['init'])
+            raise Exception, "init block has a wrong format"
+        return int(float(data[-2]))
+        
+    def set_lha_strategy(self, value):
+        """set the lha_strategy: how the weight have to be handle by the shower"""
+        
+        if not (-4 <= int(value) <= 4):
+            raise Exception, "wrong value for lha_strategy", value
+        if not self["init"]:
+            raise Exception, "No init block define"
+        
+        all_lines = self["init"].split('\n')
+        data = all_lines[0].split()
+        if len(data) != 10:
+            misc.sprint(len(data), self['init'])
+            raise Exception, "init block has a wrong format"
+        data[-2] = '%s' % value
+        all_lines[0] = ' '.join(data)
+        self['init'] = '\n'.join(all_lines)
+
+    def modify_init_cross(self, cross):
+        """modify the init information with the associate cross-section"""
+
+        assert isinstance(cross, dict)
+#        assert "all" in cross
+        assert "init" in self
+        
+        all_lines = self["init"].split('\n')
+        new_data = []
+        new_data.append(all_lines[0])
+        for i in range(1, len(all_lines)):
+            line = all_lines[i]
+            split = line.split()
+            if len(split) == 4:
+                xsec, xerr, xmax, pid = split 
+            else:
+                new_data += all_lines[i:]
+                break
+            if int(pid) not in cross:
+                raise Exception
+            pid = int(pid)
+            ratio = cross[pid]/float(xsec)
+            line = "   %+13.7e %+13.7e %+13.7e %i" % \
+                (float(cross[pid]), ratio* float(xerr), ratio*float(xmax), pid)
+            new_data.append(line)
+        self['init'] = '\n'.join(new_data)
+                
     ############################################################################
     #  WRITE BANNER
     ############################################################################
@@ -374,6 +429,8 @@
                 tag = 'madspin'
             elif 'FO_analyse_card' in card_name:
                 tag = 'foanalyse'
+            elif 'reweight_card' in card_name:
+                tag='reweight_card'
             else:
                 raise Exception, 'Impossible to know the type of the card'
 
@@ -1878,9 +1935,17 @@
         if self['pdlabel'] not in possible_set:
             raise InvalidRunCard, 'Invalid PDF set (argument of pdlabel) possible choice are:\n %s' % ','.join(possible_set)
     
-        # check that we use lhapdf if reweighting is ON
-        if self['reweight_pdf'] and self['pdlabel'] != "lhapdf":
-            raise InvalidRunCard, 'Reweight PDF option requires to use pdf sets associated to lhapdf. Please either change the pdlabel or set reweight_pdf to False.'
+
+        # PDF reweighting check
+        if self['reweight_pdf']:
+            # check that we use lhapdf if reweighting is ON
+            if self['pdlabel'] != "lhapdf":
+                raise InvalidRunCard, 'Reweight PDF option requires to use pdf sets associated to lhapdf. Please either change the pdlabel or set reweight_pdf to False.'
+            
+            # check that the number of pdf set is coherent for the reweigting:    
+            if (self['pdf_set_max'] - self['pdf_set_min'] + 1) % 2:
+                raise InvalidRunCard, "The number of PDF error sets must be even" 
+        
 
     def write(self, output_file, template=None, python_template=False):
         """Write the run_card in output_file according to template 

=== modified file 'madgraph/various/combine_plots.py'
--- madgraph/various/combine_plots.py	2013-11-13 02:21:10 +0000
+++ madgraph/various/combine_plots.py	2015-07-24 15:59:43 +0000
@@ -50,13 +50,14 @@
         self.histo[tag] = {}
         self.histo[tag]["values"] = string_values
         old_max = self.max
+
         for line in string_values.split('\n'):
             split = line.split()
-            if len(split) ==3:
+            if len(split) in [2,3]:
                 self.max = max(self.max, float(line.split()[1]))
         if self.max != old_max:
             self.max_file = tag
-            
+        
     def get_histo(self, tag, norm):
         """return a string with the histogram values, and the normalization """
         if tag not in self.histo or self.histo[tag]["values"] == '':
@@ -100,6 +101,7 @@
         for index, tag_plot in enumerate(self.order_plots):
             norm1 = 1 
             norm2 = 1
+            
             if self.plots[tag_plot].max_file == file1:
                 color1, color2 = 'WHITE','BLUE'
                 histtype1, histtype2 = 'HIST SOLID\n', 'SET PATTERN .05 .07\n' 

=== modified file 'madgraph/various/lhe_parser.py'
--- madgraph/various/lhe_parser.py	2015-07-01 22:02:22 +0000
+++ madgraph/various/lhe_parser.py	2015-07-24 15:59:43 +0000
@@ -283,6 +283,7 @@
             def weight(event):
                 return event.wgt
             get_wgt  = weight
+            unwgt_name = "central weight"
         elif isinstance(get_wgt, str):
             unwgt_name =get_wgt 
             def get_wgt(event):
@@ -1262,50 +1263,59 @@
                       'comments': self.comment,
                       'reweight': reweight_str}
         return re.sub('[\n]+', '\n', out)
+
+    def get_momenta(self, get_order, allow_reversed=True):
+        """return the momenta vector in the order asked for"""
+        
+        #avoid to modify the input
+        order = [list(get_order[0]), list(get_order[1])] 
+        out = [''] *(len(order[0])+len(order[1]))
+        for i, part in enumerate(self):
+            if part.status == 1: #final
+                try:
+                    ind = order[1].index(part.pid)
+                except ValueError, error:
+                    if not allow_reversed:
+                        raise error
+                    else:
+                        order = [[-i for i in get_order[0]],[-i for i in get_order[1]]]
+                        try:
+                            return self.get_momenta_str(order, False)
+                        except ValueError:
+                            raise error     
+                position = len(order[0]) + ind
+                order[1][ind] = 0   
+            elif part.status == -1:
+                try:
+                    ind = order[0].index(part.pid)
+                except ValueError, error:
+                    if not allow_reversed:
+                        raise error
+                    else:
+                        order = [[-i for i in get_order[0]],[-i for i in get_order[1]]]
+                        try:
+                            return self.get_momenta_str(order, False)
+                        except ValueError:
+                            raise error
+                 
+                position =  ind
+                order[0][ind] = 0
+            else: #intermediate
+                continue
+
+            out[position] = (part.E, part.px, part.py, part.pz)
+            
+        return out
+
     
     def get_momenta_str(self, get_order, allow_reversed=True):
         """return the momenta str in the order asked for"""
         
-        
-        #avoid to modify the input
-        order = [list(get_order[0]), list(get_order[1])] 
-        out = [''] *(len(order[0])+len(order[1]))
-        for i, part in enumerate(self):
-            if part.status == 1: #final
-                try:
-                    ind = order[1].index(part.pid)
-                except ValueError, error:
-                    if not allow_reversed:
-                        raise error
-                    else:
-                        order = [[-i for i in get_order[0]],[-i for i in get_order[1]]]
-                        try:
-                            return self.get_momenta_str(order, False)
-                        except ValueError:
-                            raise error     
-                position = len(order[0]) + ind
-                order[1][ind] = 0   
-            elif part.status == -1:
-                try:
-                    ind = order[0].index(part.pid)
-                except ValueError, error:
-                    if not allow_reversed:
-                        raise error
-                    else:
-                        order = [[-i for i in get_order[0]],[-i for i in get_order[1]]]
-                        try:
-                            return self.get_momenta_str(order, False)
-                        except ValueError:
-                            raise error
-                 
-                position =  ind
-                order[0][ind] = 0
-            else: #intermediate
-                continue
-            format = '%.12f'
-            format_line = ' '.join([format]*4) + ' \n'
-            out[position] = format_line % (part.E, part.px, part.py, part.pz)
-            
+        out = self.get_momenta(get_order, allow_reversed)
+        #format
+        format = '%.12f'
+        format_line = ' '.join([format]*4) + ' \n'
+        out = [format_line % one for one in out]
         out = ''.join(out).replace('e','d')
         return out    
 

=== modified file 'madgraph/various/misc.py'
--- madgraph/various/misc.py	2015-06-24 12:54:39 +0000
+++ madgraph/various/misc.py	2015-07-24 15:59:43 +0000
@@ -15,6 +15,7 @@
 
 """A set of functions performing routine administrative I/O tasks."""
 
+import contextlib
 import logging
 import os
 import re
@@ -169,6 +170,20 @@
                 return exe_file
     return None
 
+def has_f2py():
+    has_f2py = False
+    if which('f2py'):
+        has_f2py = True
+    elif sys.version_info[1] == 6:
+        if which('f2py-2.6'):
+            has_f2py = True
+    else:
+        if which('f2py-2.7'):
+            has_f2py = True 
+    return has_f2py       
+        
+        
+        
 #===============================================================================
 # find a library
 #===============================================================================
@@ -447,6 +462,49 @@
             #for i, h in enumerate(my_logger.handlers):
             #    h.setLevel(cls.logger_saved_info[logname][2][i])
 
+nb_open =0
[e-mail]
+def stdchannel_redirected(stdchannel, dest_filename):
+    """                                                                                                                                                                                                     
+    A context manager to temporarily redirect stdout or stderr                                                                                                                                              
+                                                                                                                                                                                                            
+    e.g.:                                                                                                                                                                                                   
+                                                                                                                                                                                                            
+                                                                                                                                                                                                            
+    with stdchannel_redirected(sys.stderr, os.devnull):                                                                                                                                                     
+        if compiler.has_function('clock_gettime', libraries=['rt']):                                                                                                                                        
+            libraries.append('rt')                                                                                                                                                                          
+    """
+
+    try:
+        oldstdchannel = os.dup(stdchannel.fileno())
+        dest_file = open(dest_filename, 'w')
+        os.dup2(dest_file.fileno(), stdchannel.fileno())
+        yield
+    finally:
+        if oldstdchannel is not None:
+            os.dup2(oldstdchannel, stdchannel.fileno())
+            os.close(oldstdchannel)
+        if dest_file is not None:
+            dest_file.close()
+        
+def get_open_fds():
+    '''
+    return the number of open file descriptors for current process
+
+    .. warning: will only work on UNIX-like os-es.
+    '''
+    import subprocess
+    import os
+
+    pid = os.getpid()
+    procs = subprocess.check_output( 
+        [ "lsof", '-w', '-Ff', "-p", str( pid ) ] )
+    nprocs = filter( 
+            lambda s: s and s[ 0 ] == 'f' and s[1: ].isdigit(),
+            procs.split( '\n' ) )
+        
+    return nprocs
 
 def detect_current_compiler(path, compiler_type='fortran'):
     """find the current compiler for the current directory"""
@@ -1031,6 +1089,29 @@
     def __exit__(self, etype, value, traceback):
         os.chdir(self.savedPath)
 
+################################################################################
+# Timeout FUNCTION
+################################################################################
+
+def timeout(func, args=(), kwargs={}, timeout_duration=1, default=None):
+    '''This function will spwan a thread and run the given function using the args, kwargs and 
+    return the given default value if the timeout_duration is exceeded 
+    ''' 
+    import threading
+    class InterruptableThread(threading.Thread):
+        def __init__(self):
+            threading.Thread.__init__(self)
+            self.result = default
+        def run(self):
+            try:
+                self.result = func(*args, **kwargs)
+            except Exception,error:
+                print error
+                self.result = default
+    it = InterruptableThread()
+    it.start()
+    it.join(timeout_duration)
+    return it.result
 
 
 ################################################################################

=== modified file 'tests/acceptance_tests/test_cmd.py'
--- tests/acceptance_tests/test_cmd.py	2015-07-01 22:02:22 +0000
+++ tests/acceptance_tests/test_cmd.py	2015-07-24 15:59:43 +0000
@@ -637,7 +637,7 @@
                                                'lib', 'libmodel.a')))
         # Check that check_sa.f compiles
         subprocess.call(['make', 'check'],
-                        stdout=devnull, stderr=devnull, 
+#                        stdout=devnull, stderr=devnull, 
                         cwd=os.path.join(self.out_dir, 'SubProcesses',
                                          'P0_gg_hgg'))
         self.assertTrue(os.path.exists(os.path.join(self.out_dir,
@@ -1099,6 +1099,13 @@
     def test_madevent_subproc_group_symmetry(self):
         """Check that symmetry.f gives right output"""
 
+        def analyse(fsock):
+            data = []
+            for line in fsock:
+                if line.strip():
+                    data.append([int(i) for i in line.split()])
+            return data
+
         if os.path.isdir(self.out_dir):
             shutil.rmtree(self.out_dir)
 
@@ -1114,12 +1121,7 @@
                                                     'SubProcesses',
                                                     'P0_qq_gogo_go_qqn1_go_qqn1')))
         
-        # Check the contents of the symfact.dat file
-        self.assertEqual(open(os.path.join(self.out_dir,
-                                           'SubProcesses',
-                                           'P0_qq_gogo_go_qqn1_go_qqn1',
-                                           'symfact_orig.dat')).read().split('\n'),
-                         """ 1   1
+        target=""" 1   1
  2  -1
  3  -1
  4  -1
@@ -1131,7 +1133,13 @@
 10  -9
 11  -9
 12  -9
-""".split('\n'))
+"""
+        
+        self.assertEqual(analyse(target.split('\n')), 
+                         analyse(open(os.path.join(self.out_dir,
+                                           'SubProcesses',
+                                           'P0_qq_gogo_go_qqn1_go_qqn1',
+                                           'symfact_orig.dat'))))
 
         # Compile the Source directory
         status = subprocess.call(['make'],
@@ -1152,8 +1160,13 @@
         proc.communicate('100 4 0.1 .false.\n')
         self.assertEqual(proc.returncode, 0)
 
+<<<<<<< TREE
 
         target =                          """   1   1
+=======
+
+        target ="""   1   1
+>>>>>>> MERGE-SOURCE
    2  -1
    3  -1
    4  -1
@@ -1165,6 +1178,7 @@
   10  -9
   11  -9
   12  -9
+<<<<<<< TREE
 """
         data = [[int(i) for i in line.split()] for line in target.split('\n')]
         result = []
@@ -1177,6 +1191,18 @@
             self.assertEqual(info, data[ii])
             ii+=1
 
+=======
+"""
+            
+        # Check the new contents of the symfact.dat file
+        self.assertEqual(analyse(open(os.path.join(self.out_dir,
+                                           'SubProcesses',
+                                           'P0_qq_gogo_go_qqn1_go_qqn1',
+                                           'symfact.dat'))), 
+                         analyse(target.split('\n')))
+                         
+        
+>>>>>>> MERGE-SOURCE
     def test_madevent_subproc_group_decay_chain(self):
         """Test decay chain output using the SubProcess group functionality"""
 

=== added directory 'tests/input_files/IOTestsComparison/ExportV4IOTest'
=== added directory 'tests/input_files/IOTestsComparison/ExportV4IOTest/export_matrix_element_v4_standalone'
=== added file 'tests/input_files/IOTestsComparison/ExportV4IOTest/export_matrix_element_v4_standalone/matrix.f'
--- tests/input_files/IOTestsComparison/ExportV4IOTest/export_matrix_element_v4_standalone/matrix.f	1970-01-01 00:00:00 +0000
+++ tests/input_files/IOTestsComparison/ExportV4IOTest/export_matrix_element_v4_standalone/matrix.f	2015-07-24 15:59:43 +0000
@@ -0,0 +1,272 @@
+      SUBROUTINE SMATRIXHEL(P,HEL,ANS)
+      IMPLICIT NONE
+C     
+C     CONSTANT
+C     
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=5)
+      INTEGER                 NCOMB
+      PARAMETER (             NCOMB=32)
+CF2PY INTENT(OUT) :: ANS
+CF2PY INTENT(IN) :: HEL
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL)
+
+C     
+C     ARGUMENTS 
+C     
+      REAL*8 P(0:3,NEXTERNAL),ANS
+      INTEGER HEL
+C     
+C     GLOBAL VARIABLES
+C     
+      INTEGER USERHEL
+      COMMON/HELUSERCHOICE/USERHEL
+C     ----------
+C     BEGIN CODE
+C     ----------
+      USERHEL=HEL
+      CALL SMATRIX(P,ANS)
+      USERHEL=-1
+
+      END
+
+      SUBROUTINE SMATRIX(P,ANS)
+C     
+C     Generated by MadGraph5_aMC@NLO v. %(version)s, %(date)s
+C     By the MadGraph5_aMC@NLO Development Team
+C     Visit launchpad.net/madgraph5 and amcatnlo.web.cern.ch
+C     
+C     MadGraph5_aMC@NLO StandAlone Version
+C     
+C     Returns amplitude squared summed/avg over colors
+C     and helicities
+C     for the point in phase space P(0:3,NEXTERNAL)
+C     
+C     Process: e+ e- > a a a
+C     
+      IMPLICIT NONE
+C     
+C     CONSTANTS
+C     
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=5)
+      INTEGER                 NCOMB
+      PARAMETER (             NCOMB=32)
+      INTEGER HELAVGFACTOR
+      PARAMETER (HELAVGFACTOR=4)
+C     
+C     ARGUMENTS 
+C     
+      REAL*8 P(0:3,NEXTERNAL),ANS
+CF2PY INTENT(OUT) :: ANS
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL)
+C     
+C     LOCAL VARIABLES 
+C     
+      INTEGER NHEL(NEXTERNAL,NCOMB),NTRY
+      REAL*8 T
+      REAL*8 MATRIX
+      INTEGER IHEL,IDEN, I
+      INTEGER JC(NEXTERNAL)
+      LOGICAL GOODHEL(NCOMB)
+      DATA NTRY/0/
+      DATA GOODHEL/NCOMB*.FALSE./
+
+C     
+C     GLOBAL VARIABLES
+C     
+      INTEGER USERHEL
+      COMMON/HELUSERCHOICE/USERHEL
+      DATA USERHEL/-1/
+
+      DATA (NHEL(I,   1),I=1,5) /-1, 1,-1,-1,-1/
+      DATA (NHEL(I,   2),I=1,5) /-1, 1,-1,-1, 1/
+      DATA (NHEL(I,   3),I=1,5) /-1, 1,-1, 1,-1/
+      DATA (NHEL(I,   4),I=1,5) /-1, 1,-1, 1, 1/
+      DATA (NHEL(I,   5),I=1,5) /-1, 1, 1,-1,-1/
+      DATA (NHEL(I,   6),I=1,5) /-1, 1, 1,-1, 1/
+      DATA (NHEL(I,   7),I=1,5) /-1, 1, 1, 1,-1/
+      DATA (NHEL(I,   8),I=1,5) /-1, 1, 1, 1, 1/
+      DATA (NHEL(I,   9),I=1,5) /-1,-1,-1,-1,-1/
+      DATA (NHEL(I,  10),I=1,5) /-1,-1,-1,-1, 1/
+      DATA (NHEL(I,  11),I=1,5) /-1,-1,-1, 1,-1/
+      DATA (NHEL(I,  12),I=1,5) /-1,-1,-1, 1, 1/
+      DATA (NHEL(I,  13),I=1,5) /-1,-1, 1,-1,-1/
+      DATA (NHEL(I,  14),I=1,5) /-1,-1, 1,-1, 1/
+      DATA (NHEL(I,  15),I=1,5) /-1,-1, 1, 1,-1/
+      DATA (NHEL(I,  16),I=1,5) /-1,-1, 1, 1, 1/
+      DATA (NHEL(I,  17),I=1,5) / 1, 1,-1,-1,-1/
+      DATA (NHEL(I,  18),I=1,5) / 1, 1,-1,-1, 1/
+      DATA (NHEL(I,  19),I=1,5) / 1, 1,-1, 1,-1/
+      DATA (NHEL(I,  20),I=1,5) / 1, 1,-1, 1, 1/
+      DATA (NHEL(I,  21),I=1,5) / 1, 1, 1,-1,-1/
+      DATA (NHEL(I,  22),I=1,5) / 1, 1, 1,-1, 1/
+      DATA (NHEL(I,  23),I=1,5) / 1, 1, 1, 1,-1/
+      DATA (NHEL(I,  24),I=1,5) / 1, 1, 1, 1, 1/
+      DATA (NHEL(I,  25),I=1,5) / 1,-1,-1,-1,-1/
+      DATA (NHEL(I,  26),I=1,5) / 1,-1,-1,-1, 1/
+      DATA (NHEL(I,  27),I=1,5) / 1,-1,-1, 1,-1/
+      DATA (NHEL(I,  28),I=1,5) / 1,-1,-1, 1, 1/
+      DATA (NHEL(I,  29),I=1,5) / 1,-1, 1,-1,-1/
+      DATA (NHEL(I,  30),I=1,5) / 1,-1, 1,-1, 1/
+      DATA (NHEL(I,  31),I=1,5) / 1,-1, 1, 1,-1/
+      DATA (NHEL(I,  32),I=1,5) / 1,-1, 1, 1, 1/
+      DATA IDEN/24/
+C     ----------
+C     BEGIN CODE
+C     ----------
+      IF(USERHEL.EQ.-1) NTRY=NTRY+1
+      DO IHEL=1,NEXTERNAL
+        JC(IHEL) = +1
+      ENDDO
+      ANS = 0D0
+      DO IHEL=1,NCOMB
+        IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN
+          IF (GOODHEL(IHEL) .OR. NTRY .LT. 20.OR.USERHEL.NE.-1) THEN
+            T=MATRIX(P ,NHEL(1,IHEL),JC(1))
+            ANS=ANS+T
+            IF (T .NE. 0D0 .AND. .NOT.    GOODHEL(IHEL)) THEN
+              GOODHEL(IHEL)=.TRUE.
+            ENDIF
+          ENDIF
+        ENDIF
+      ENDDO
+      ANS=ANS/DBLE(IDEN)
+      IF(USERHEL.NE.-1) THEN
+        ANS=ANS*HELAVGFACTOR
+      ENDIF
+      END
+
+
+      REAL*8 FUNCTION MATRIX(P,NHEL,IC)
+C     
+C     Generated by MadGraph5_aMC@NLO v. %(version)s, %(date)s
+C     By the MadGraph5_aMC@NLO Development Team
+C     Visit launchpad.net/madgraph5 and amcatnlo.web.cern.ch
+C     
+C     Returns amplitude squared summed/avg over colors
+C     for the point with external lines W(0:6,NEXTERNAL)
+C     
+C     Process: e+ e- > a a a
+C     
+      IMPLICIT NONE
+C     
+C     CONSTANTS
+C     
+      INTEGER    NGRAPHS
+      PARAMETER (NGRAPHS=6)
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=5)
+      INTEGER    NWAVEFUNCS, NCOLOR
+      PARAMETER (NWAVEFUNCS=9, NCOLOR=1)
+      REAL*8     ZERO
+      PARAMETER (ZERO=0D0)
+      COMPLEX*16 IMAG1
+      PARAMETER (IMAG1=(0D0,1D0))
+C     
+C     ARGUMENTS 
+C     
+      REAL*8 P(0:3,NEXTERNAL)
+      INTEGER NHEL(NEXTERNAL), IC(NEXTERNAL)
+C     
+C     LOCAL VARIABLES 
+C     
+      INTEGER I,J
+      COMPLEX*16 ZTEMP
+      REAL*8 DENOM(NCOLOR), CF(NCOLOR,NCOLOR)
+      COMPLEX*16 AMP(NGRAPHS), JAMP(NCOLOR)
+      COMPLEX*16 W(18,NWAVEFUNCS)
+      COMPLEX*16 DUM0,DUM1
+      DATA DUM0, DUM1/(0D0, 0D0), (1D0, 0D0)/
+C     
+C     GLOBAL VARIABLES
+C     
+      INCLUDE 'coupl.inc'
+
+C     
+C     COLOR DATA
+C     
+      DATA DENOM(1)/1/
+      DATA (CF(I,  1),I=  1,  1) /    1/
+C     1 ColorOne()
+C     ----------
+C     BEGIN CODE
+C     ----------
+      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))
+      CALL VXXXXX(P(0,3),ZERO,NHEL(3),+1*IC(3),W(1,3))
+      CALL VXXXXX(P(0,4),ZERO,NHEL(4),+1*IC(4),W(1,4))
+      CALL VXXXXX(P(0,5),ZERO,NHEL(5),+1*IC(5),W(1,5))
+      CALL FVOXXX(W(1,1),W(1,3),MGVX12,ZERO,ZERO,W(1,6))
+      CALL FVIXXX(W(1,2),W(1,4),MGVX12,ZERO,ZERO,W(1,7))
+C     Amplitude(s) for diagram number 1
+      CALL IOVXXX(W(1,7),W(1,6),W(1,5),MGVX12,AMP(1))
+      CALL FVIXXX(W(1,2),W(1,5),MGVX12,ZERO,ZERO,W(1,8))
+C     Amplitude(s) for diagram number 2
+      CALL IOVXXX(W(1,8),W(1,6),W(1,4),MGVX12,AMP(2))
+      CALL FVOXXX(W(1,1),W(1,4),MGVX12,ZERO,ZERO,W(1,6))
+      CALL FVIXXX(W(1,2),W(1,3),MGVX12,ZERO,ZERO,W(1,9))
+C     Amplitude(s) for diagram number 3
+      CALL IOVXXX(W(1,9),W(1,6),W(1,5),MGVX12,AMP(3))
+C     Amplitude(s) for diagram number 4
+      CALL IOVXXX(W(1,8),W(1,6),W(1,3),MGVX12,AMP(4))
+      CALL FVOXXX(W(1,1),W(1,5),MGVX12,ZERO,ZERO,W(1,6))
+C     Amplitude(s) for diagram number 5
+      CALL IOVXXX(W(1,9),W(1,6),W(1,4),MGVX12,AMP(5))
+C     Amplitude(s) for diagram number 6
+      CALL IOVXXX(W(1,7),W(1,6),W(1,3),MGVX12,AMP(6))
+      JAMP(1)=-AMP(1)-AMP(2)-AMP(3)-AMP(4)-AMP(5)-AMP(6)
+
+      MATRIX = 0.D0
+      DO I = 1, NCOLOR
+        ZTEMP = (0.D0,0.D0)
+        DO J = 1, NCOLOR
+          ZTEMP = ZTEMP + CF(J,I)*JAMP(J)
+        ENDDO
+        MATRIX = MATRIX+ZTEMP*DCONJG(JAMP(I))/DENOM(I)
+      ENDDO
+
+      END
+
+      SUBROUTINE GET_ME(P, ALPHAS, NHEL ,ANS)
+      IMPLICIT NONE
+C     
+C     CONSTANT
+C     
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=5)
+C     
+C     ARGUMENTS 
+C     
+      REAL*8 P(0:3,NEXTERNAL),ANS
+      INTEGER NHEL
+      DOUBLE PRECISION ALPHAS
+      REAL*8 PI
+CF2PY INTENT(OUT) :: ANS
+CF2PY INTENT(IN) :: NHEL
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL)
+CF2PY INTENT(IN) :: ALPHAS
+C     ROUTINE FOR F2PY to read the benchmark point.    
+C     the include file with the values of the parameters and masses 
+      INCLUDE 'coupl.inc'
+
+      PI = 3.141592653589793D0
+      G = 2* DSQRT(ALPHAS*PI)
+      CALL UPDATE_AS_PARAM()
+      IF (NHEL.NE.0) THEN
+        CALL SMATRIXHEL(P, NHEL, ANS)
+      ELSE
+        CALL SMATRIX(P, ANS)
+      ENDIF
+      RETURN
+      END
+
+      SUBROUTINE INITIALISE(PATH)
+C     ROUTINE FOR F2PY to read the benchmark point.    
+      IMPLICIT NONE
+      CHARACTER*180 PATH
+CF2PY INTENT(IN) :: PATH
+      CALL SETPARA(PATH)  !first call to setup the paramaters    
+      RETURN
+      END
+

=== added directory 'tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_standalone'
=== added 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	1970-01-01 00:00:00 +0000
+++ tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_standalone/matrix.f	2015-07-24 15:59:43 +0000
@@ -0,0 +1,272 @@
+      SUBROUTINE SMATRIXHEL(P,HEL,ANS)
+      IMPLICIT NONE
+C     
+C     CONSTANT
+C     
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=5)
+      INTEGER                 NCOMB
+      PARAMETER (             NCOMB=32)
+CF2PY INTENT(OUT) :: ANS
+CF2PY INTENT(IN) :: HEL
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL)
+
+C     
+C     ARGUMENTS 
+C     
+      REAL*8 P(0:3,NEXTERNAL),ANS
+      INTEGER HEL
+C     
+C     GLOBAL VARIABLES
+C     
+      INTEGER USERHEL
+      COMMON/HELUSERCHOICE/USERHEL
+C     ----------
+C     BEGIN CODE
+C     ----------
+      USERHEL=HEL
+      CALL SMATRIX(P,ANS)
+      USERHEL=-1
+
+      END
+
+      SUBROUTINE SMATRIX(P,ANS)
+C     
+C     Generated by MadGraph5_aMC@NLO v. %(version)s, %(date)s
+C     By the MadGraph5_aMC@NLO Development Team
+C     Visit launchpad.net/madgraph5 and amcatnlo.web.cern.ch
+C     
+C     MadGraph5_aMC@NLO StandAlone Version
+C     
+C     Returns amplitude squared summed/avg over colors
+C     and helicities
+C     for the point in phase space P(0:3,NEXTERNAL)
+C     
+C     Process: e+ e- > a a a
+C     
+      IMPLICIT NONE
+C     
+C     CONSTANTS
+C     
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=5)
+      INTEGER                 NCOMB
+      PARAMETER (             NCOMB=32)
+      INTEGER HELAVGFACTOR
+      PARAMETER (HELAVGFACTOR=4)
+C     
+C     ARGUMENTS 
+C     
+      REAL*8 P(0:3,NEXTERNAL),ANS
+CF2PY INTENT(OUT) :: ANS
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL)
+C     
+C     LOCAL VARIABLES 
+C     
+      INTEGER NHEL(NEXTERNAL,NCOMB),NTRY
+      REAL*8 T
+      REAL*8 MATRIX
+      INTEGER IHEL,IDEN, I
+      INTEGER JC(NEXTERNAL)
+      LOGICAL GOODHEL(NCOMB)
+      DATA NTRY/0/
+      DATA GOODHEL/NCOMB*.FALSE./
+
+C     
+C     GLOBAL VARIABLES
+C     
+      INTEGER USERHEL
+      COMMON/HELUSERCHOICE/USERHEL
+      DATA USERHEL/-1/
+
+      DATA (NHEL(I,   1),I=1,5) /-1, 1,-1,-1,-1/
+      DATA (NHEL(I,   2),I=1,5) /-1, 1,-1,-1, 1/
+      DATA (NHEL(I,   3),I=1,5) /-1, 1,-1, 1,-1/
+      DATA (NHEL(I,   4),I=1,5) /-1, 1,-1, 1, 1/
+      DATA (NHEL(I,   5),I=1,5) /-1, 1, 1,-1,-1/
+      DATA (NHEL(I,   6),I=1,5) /-1, 1, 1,-1, 1/
+      DATA (NHEL(I,   7),I=1,5) /-1, 1, 1, 1,-1/
+      DATA (NHEL(I,   8),I=1,5) /-1, 1, 1, 1, 1/
+      DATA (NHEL(I,   9),I=1,5) /-1,-1,-1,-1,-1/
+      DATA (NHEL(I,  10),I=1,5) /-1,-1,-1,-1, 1/
+      DATA (NHEL(I,  11),I=1,5) /-1,-1,-1, 1,-1/
+      DATA (NHEL(I,  12),I=1,5) /-1,-1,-1, 1, 1/
+      DATA (NHEL(I,  13),I=1,5) /-1,-1, 1,-1,-1/
+      DATA (NHEL(I,  14),I=1,5) /-1,-1, 1,-1, 1/
+      DATA (NHEL(I,  15),I=1,5) /-1,-1, 1, 1,-1/
+      DATA (NHEL(I,  16),I=1,5) /-1,-1, 1, 1, 1/
+      DATA (NHEL(I,  17),I=1,5) / 1, 1,-1,-1,-1/
+      DATA (NHEL(I,  18),I=1,5) / 1, 1,-1,-1, 1/
+      DATA (NHEL(I,  19),I=1,5) / 1, 1,-1, 1,-1/
+      DATA (NHEL(I,  20),I=1,5) / 1, 1,-1, 1, 1/
+      DATA (NHEL(I,  21),I=1,5) / 1, 1, 1,-1,-1/
+      DATA (NHEL(I,  22),I=1,5) / 1, 1, 1,-1, 1/
+      DATA (NHEL(I,  23),I=1,5) / 1, 1, 1, 1,-1/
+      DATA (NHEL(I,  24),I=1,5) / 1, 1, 1, 1, 1/
+      DATA (NHEL(I,  25),I=1,5) / 1,-1,-1,-1,-1/
+      DATA (NHEL(I,  26),I=1,5) / 1,-1,-1,-1, 1/
+      DATA (NHEL(I,  27),I=1,5) / 1,-1,-1, 1,-1/
+      DATA (NHEL(I,  28),I=1,5) / 1,-1,-1, 1, 1/
+      DATA (NHEL(I,  29),I=1,5) / 1,-1, 1,-1,-1/
+      DATA (NHEL(I,  30),I=1,5) / 1,-1, 1,-1, 1/
+      DATA (NHEL(I,  31),I=1,5) / 1,-1, 1, 1,-1/
+      DATA (NHEL(I,  32),I=1,5) / 1,-1, 1, 1, 1/
+      DATA IDEN/24/
+C     ----------
+C     BEGIN CODE
+C     ----------
+      IF(USERHEL.EQ.-1) NTRY=NTRY+1
+      DO IHEL=1,NEXTERNAL
+        JC(IHEL) = +1
+      ENDDO
+      ANS = 0D0
+      DO IHEL=1,NCOMB
+        IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN
+          IF (GOODHEL(IHEL) .OR. NTRY .LT. 20.OR.USERHEL.NE.-1) THEN
+            T=MATRIX(P ,NHEL(1,IHEL),JC(1))
+            ANS=ANS+T
+            IF (T .NE. 0D0 .AND. .NOT.    GOODHEL(IHEL)) THEN
+              GOODHEL(IHEL)=.TRUE.
+            ENDIF
+          ENDIF
+        ENDIF
+      ENDDO
+      ANS=ANS/DBLE(IDEN)
+      IF(USERHEL.NE.-1) THEN
+        ANS=ANS*HELAVGFACTOR
+      ENDIF
+      END
+
+
+      REAL*8 FUNCTION MATRIX(P,NHEL,IC)
+C     
+C     Generated by MadGraph5_aMC@NLO v. %(version)s, %(date)s
+C     By the MadGraph5_aMC@NLO Development Team
+C     Visit launchpad.net/madgraph5 and amcatnlo.web.cern.ch
+C     
+C     Returns amplitude squared summed/avg over colors
+C     for the point with external lines W(0:6,NEXTERNAL)
+C     
+C     Process: e+ e- > a a a
+C     
+      IMPLICIT NONE
+C     
+C     CONSTANTS
+C     
+      INTEGER    NGRAPHS
+      PARAMETER (NGRAPHS=6)
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=5)
+      INTEGER    NWAVEFUNCS, NCOLOR
+      PARAMETER (NWAVEFUNCS=9, NCOLOR=1)
+      REAL*8     ZERO
+      PARAMETER (ZERO=0D0)
+      COMPLEX*16 IMAG1
+      PARAMETER (IMAG1=(0D0,1D0))
+C     
+C     ARGUMENTS 
+C     
+      REAL*8 P(0:3,NEXTERNAL)
+      INTEGER NHEL(NEXTERNAL), IC(NEXTERNAL)
+C     
+C     LOCAL VARIABLES 
+C     
+      INTEGER I,J
+      COMPLEX*16 ZTEMP
+      REAL*8 DENOM(NCOLOR), CF(NCOLOR,NCOLOR)
+      COMPLEX*16 AMP(NGRAPHS), JAMP(NCOLOR)
+      COMPLEX*16 W(18,NWAVEFUNCS)
+      COMPLEX*16 DUM0,DUM1
+      DATA DUM0, DUM1/(0D0, 0D0), (1D0, 0D0)/
+C     
+C     GLOBAL VARIABLES
+C     
+      INCLUDE 'coupl.inc'
+
+C     
+C     COLOR DATA
+C     
+      DATA DENOM(1)/1/
+      DATA (CF(I,  1),I=  1,  1) /    1/
+C     1 ColorOne()
+C     ----------
+C     BEGIN CODE
+C     ----------
+      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))
+      CALL VXXXXX(P(0,3),ZERO,NHEL(3),+1*IC(3),W(1,3))
+      CALL VXXXXX(P(0,4),ZERO,NHEL(4),+1*IC(4),W(1,4))
+      CALL VXXXXX(P(0,5),ZERO,NHEL(5),+1*IC(5),W(1,5))
+      CALL FVOXXX(W(1,1),W(1,3),MGVX12,ZERO,ZERO,W(1,6))
+      CALL FVIXXX(W(1,2),W(1,4),MGVX12,ZERO,ZERO,W(1,7))
+C     Amplitude(s) for diagram number 1
+      CALL IOVXXX(W(1,7),W(1,6),W(1,5),MGVX12,AMP(1))
+      CALL FVIXXX(W(1,2),W(1,5),MGVX12,ZERO,ZERO,W(1,8))
+C     Amplitude(s) for diagram number 2
+      CALL IOVXXX(W(1,8),W(1,6),W(1,4),MGVX12,AMP(2))
+      CALL FVOXXX(W(1,1),W(1,4),MGVX12,ZERO,ZERO,W(1,6))
+      CALL FVIXXX(W(1,2),W(1,3),MGVX12,ZERO,ZERO,W(1,9))
+C     Amplitude(s) for diagram number 3
+      CALL IOVXXX(W(1,9),W(1,6),W(1,5),MGVX12,AMP(3))
+C     Amplitude(s) for diagram number 4
+      CALL IOVXXX(W(1,8),W(1,6),W(1,3),MGVX12,AMP(4))
+      CALL FVOXXX(W(1,1),W(1,5),MGVX12,ZERO,ZERO,W(1,6))
+C     Amplitude(s) for diagram number 5
+      CALL IOVXXX(W(1,9),W(1,6),W(1,4),MGVX12,AMP(5))
+C     Amplitude(s) for diagram number 6
+      CALL IOVXXX(W(1,7),W(1,6),W(1,3),MGVX12,AMP(6))
+      JAMP(1)=-AMP(1)-AMP(2)-AMP(3)-AMP(4)-AMP(5)-AMP(6)
+
+      MATRIX = 0.D0
+      DO I = 1, NCOLOR
+        ZTEMP = (0.D0,0.D0)
+        DO J = 1, NCOLOR
+          ZTEMP = ZTEMP + CF(J,I)*JAMP(J)
+        ENDDO
+        MATRIX = MATRIX+ZTEMP*DCONJG(JAMP(I))/DENOM(I)
+      ENDDO
+
+      END
+
+      SUBROUTINE GET_ME(P, ALPHAS, NHEL ,ANS)
+      IMPLICIT NONE
+C     
+C     CONSTANT
+C     
+      INTEGER    NEXTERNAL
+      PARAMETER (NEXTERNAL=5)
+C     
+C     ARGUMENTS 
+C     
+      REAL*8 P(0:3,NEXTERNAL),ANS
+      INTEGER NHEL
+      DOUBLE PRECISION ALPHAS
+      REAL*8 PI
+CF2PY INTENT(OUT) :: ANS
+CF2PY INTENT(IN) :: NHEL
+CF2PY INTENT(IN) :: P(0:3,NEXTERNAL)
+CF2PY INTENT(IN) :: ALPHAS
+C     ROUTINE FOR F2PY to read the benchmark point.    
+C     the include file with the values of the parameters and masses 
+      INCLUDE 'coupl.inc'
+
+      PI = 3.141592653589793D0
+      G = 2* DSQRT(ALPHAS*PI)
+      CALL UPDATE_AS_PARAM()
+      IF (NHEL.NE.0) THEN
+        CALL SMATRIXHEL(P, NHEL, ANS)
+      ELSE
+        CALL SMATRIX(P, ANS)
+      ENDIF
+      RETURN
+      END
+
+      SUBROUTINE INITIALISE(PATH)
+C     ROUTINE FOR F2PY to read the benchmark point.    
+      IMPLICIT NONE
+      CHARACTER*180 PATH
+CF2PY INTENT(IN) :: PATH
+      CALL SETPARA(PATH)  !first call to setup the paramaters    
+      RETURN
+      END
+

=== 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	2015-06-30 16:20:47 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_eq_4.f	2015-07-24 15:59:43 +0000
@@ -420,9 +420,9 @@
 C        if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
         IF(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
           IF(NROTATIONS_DP.NE.0.OR.NROTATIONS_QP.NE.0) THEN
-            WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disable'
-     $       //'d, so MadLoop automatically resets NROTATIONS_DP an'
-     $       //'d NROTATIONS_QP to 0.'
+            WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING i'
+     $       //'s disabled, so MadLoop automatically resets NROTATIONS'
+     $       //'_DP and NROTATIONS_QP to 0.'
             NROTATIONS_QP=0
             NROTATIONS_DP=0
           ENDIF
@@ -518,10 +518,17 @@
             ENDDO
           ENDDO
         ENDDO
+<<<<<<< TREE
         IF(BOOTANDSTOP) THEN
           WRITE(*,*) 'Stopped by user request.'
           STOP
         ENDIF
+=======
+      ENDDO
+      IF(BOOTANDSTOP) THEN
+        WRITE(*,*) '##Stopped by user request.'
+        STOP
+>>>>>>> MERGE-SOURCE
       ENDIF
 
       IF(NTRY.EQ.0) THEN
@@ -2258,10 +2265,10 @@
      $       ,REF/DBLE(NCOMB),-1,-1)) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected wha'
      $         //'t are the zero helicity configurations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(29,FILE=HELFILTERFN,ERR=348)
@@ -2274,10 +2281,10 @@
      $       L(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected th'
      $         //'e helicity dependance relations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(30,FILE=HELFILTERFN,ERR=349)
@@ -2430,15 +2437,15 @@
                 WRITE(*,*) '##W03 WARNING An unstable PS point was'
      $           ,       ' detected.'
                 IF (NSQUAREDSO.NE.1) THEN
-                  WRITE(*,*) 'Accuracies for each split orde'
+                  WRITE(*,*) '##Accuracies for each split orde'
      $             //'r, starting with the summed case'
-                  WRITE(*,*) 'DP accuracies (for each split order): '
-     $             ,(TEMP1(I),I=0,NSQUAREDSO)
-                  WRITE(*,*) 'QP accuracies (for each split order): '
-     $             ,(ACC(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##DP accuracies (for each split orde'
+     $             //'r): ',(TEMP1(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##QP accuracies (for each split orde'
+     $             //'r): ',(ACC(I),I=0,NSQUAREDSO)
                 ELSE
-                  WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-                  WRITE(*,*) 'QP accuracy: ',ACC(1)
+                  WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+                  WRITE(*,*) '##QP accuracy: ',ACC(1)
                 ENDIF
                 DO J=0,NSQUAREDSO
                   IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== 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	2015-06-30 16:20:47 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_0_QEDAmpAndQEDsq_gt_2.f	2015-07-24 15:59:43 +0000
@@ -420,9 +420,9 @@
 C        if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
         IF(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
           IF(NROTATIONS_DP.NE.0.OR.NROTATIONS_QP.NE.0) THEN
-            WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disable'
-     $       //'d, so MadLoop automatically resets NROTATIONS_DP an'
-     $       //'d NROTATIONS_QP to 0.'
+            WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING i'
+     $       //'s disabled, so MadLoop automatically resets NROTATIONS'
+     $       //'_DP and NROTATIONS_QP to 0.'
             NROTATIONS_QP=0
             NROTATIONS_DP=0
           ENDIF
@@ -518,10 +518,17 @@
             ENDDO
           ENDDO
         ENDDO
+<<<<<<< TREE
         IF(BOOTANDSTOP) THEN
           WRITE(*,*) 'Stopped by user request.'
           STOP
         ENDIF
+=======
+      ENDDO
+      IF(BOOTANDSTOP) THEN
+        WRITE(*,*) '##Stopped by user request.'
+        STOP
+>>>>>>> MERGE-SOURCE
       ENDIF
 
       IF(NTRY.EQ.0) THEN
@@ -9141,10 +9148,10 @@
      $       ,REF/DBLE(NCOMB),-1,-1)) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected wha'
      $         //'t are the zero helicity configurations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(29,FILE=HELFILTERFN,ERR=348)
@@ -9157,10 +9164,10 @@
      $       L(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected th'
      $         //'e helicity dependance relations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(30,FILE=HELFILTERFN,ERR=349)
@@ -9313,15 +9320,15 @@
                 WRITE(*,*) '##W03 WARNING An unstable PS point was'
      $           ,       ' detected.'
                 IF (NSQUAREDSO.NE.1) THEN
-                  WRITE(*,*) 'Accuracies for each split orde'
+                  WRITE(*,*) '##Accuracies for each split orde'
      $             //'r, starting with the summed case'
-                  WRITE(*,*) 'DP accuracies (for each split order): '
-     $             ,(TEMP1(I),I=0,NSQUAREDSO)
-                  WRITE(*,*) 'QP accuracies (for each split order): '
-     $             ,(ACC(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##DP accuracies (for each split orde'
+     $             //'r): ',(TEMP1(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##QP accuracies (for each split orde'
+     $             //'r): ',(ACC(I),I=0,NSQUAREDSO)
                 ELSE
-                  WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-                  WRITE(*,*) 'QP accuracy: ',ACC(1)
+                  WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+                  WRITE(*,*) '##QP accuracy: ',ACC(1)
                 ENDIF
                 DO J=0,NSQUAREDSO
                   IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== 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	2015-06-30 16:20:47 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_4.f	2015-07-24 15:59:43 +0000
@@ -412,9 +412,9 @@
 C        if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
         IF(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
           IF(NROTATIONS_DP.NE.0.OR.NROTATIONS_QP.NE.0) THEN
-            WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disable'
-     $       //'d, so MadLoop automatically resets NROTATIONS_DP an'
-     $       //'d NROTATIONS_QP to 0.'
+            WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING i'
+     $       //'s disabled, so MadLoop automatically resets NROTATIONS'
+     $       //'_DP and NROTATIONS_QP to 0.'
             NROTATIONS_QP=0
             NROTATIONS_DP=0
           ENDIF
@@ -510,10 +510,17 @@
             ENDDO
           ENDDO
         ENDDO
+<<<<<<< TREE
         IF(BOOTANDSTOP) THEN
           WRITE(*,*) 'Stopped by user request.'
           STOP
         ENDIF
+=======
+      ENDDO
+      IF(BOOTANDSTOP) THEN
+        WRITE(*,*) '##Stopped by user request.'
+        STOP
+>>>>>>> MERGE-SOURCE
       ENDIF
 
       IF(NTRY.EQ.0) THEN
@@ -1304,10 +1311,10 @@
      $       ,REF/DBLE(NCOMB),-1,-1)) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected wha'
      $         //'t are the zero helicity configurations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(29,FILE=HELFILTERFN,ERR=348)
@@ -1320,10 +1327,10 @@
      $       L(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected th'
      $         //'e helicity dependance relations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(30,FILE=HELFILTERFN,ERR=349)
@@ -1476,15 +1483,15 @@
                 WRITE(*,*) '##W03 WARNING An unstable PS point was'
      $           ,       ' detected.'
                 IF (NSQUAREDSO.NE.1) THEN
-                  WRITE(*,*) 'Accuracies for each split orde'
+                  WRITE(*,*) '##Accuracies for each split orde'
      $             //'r, starting with the summed case'
-                  WRITE(*,*) 'DP accuracies (for each split order): '
-     $             ,(TEMP1(I),I=0,NSQUAREDSO)
-                  WRITE(*,*) 'QP accuracies (for each split order): '
-     $             ,(ACC(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##DP accuracies (for each split orde'
+     $             //'r): ',(TEMP1(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##QP accuracies (for each split orde'
+     $             //'r): ',(ACC(I),I=0,NSQUAREDSO)
                 ELSE
-                  WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-                  WRITE(*,*) 'QP accuracy: ',ACC(1)
+                  WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+                  WRITE(*,*) '##QP accuracy: ',ACC(1)
                 ENDIF
                 DO J=0,NSQUAREDSO
                   IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== 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	2015-06-30 16:20:47 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QEDsq_le_4.f	2015-07-24 15:59:43 +0000
@@ -420,9 +420,9 @@
 C        if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
         IF(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
           IF(NROTATIONS_DP.NE.0.OR.NROTATIONS_QP.NE.0) THEN
-            WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disable'
-     $       //'d, so MadLoop automatically resets NROTATIONS_DP an'
-     $       //'d NROTATIONS_QP to 0.'
+            WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING i'
+     $       //'s disabled, so MadLoop automatically resets NROTATIONS'
+     $       //'_DP and NROTATIONS_QP to 0.'
             NROTATIONS_QP=0
             NROTATIONS_DP=0
           ENDIF
@@ -518,10 +518,17 @@
             ENDDO
           ENDDO
         ENDDO
+<<<<<<< TREE
         IF(BOOTANDSTOP) THEN
           WRITE(*,*) 'Stopped by user request.'
           STOP
         ENDIF
+=======
+      ENDDO
+      IF(BOOTANDSTOP) THEN
+        WRITE(*,*) '##Stopped by user request.'
+        STOP
+>>>>>>> MERGE-SOURCE
       ENDIF
 
       IF(NTRY.EQ.0) THEN
@@ -9386,10 +9393,10 @@
      $       ,REF/DBLE(NCOMB),-1,-1)) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected wha'
      $         //'t are the zero helicity configurations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(29,FILE=HELFILTERFN,ERR=348)
@@ -9402,10 +9409,10 @@
      $       L(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected th'
      $         //'e helicity dependance relations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(30,FILE=HELFILTERFN,ERR=349)
@@ -9558,15 +9565,15 @@
                 WRITE(*,*) '##W03 WARNING An unstable PS point was'
      $           ,       ' detected.'
                 IF (NSQUAREDSO.NE.1) THEN
-                  WRITE(*,*) 'Accuracies for each split orde'
+                  WRITE(*,*) '##Accuracies for each split orde'
      $             //'r, starting with the summed case'
-                  WRITE(*,*) 'DP accuracies (for each split order): '
-     $             ,(TEMP1(I),I=0,NSQUAREDSO)
-                  WRITE(*,*) 'QP accuracies (for each split order): '
-     $             ,(ACC(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##DP accuracies (for each split orde'
+     $             //'r): ',(TEMP1(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##QP accuracies (for each split orde'
+     $             //'r): ',(ACC(I),I=0,NSQUAREDSO)
                 ELSE
-                  WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-                  WRITE(*,*) 'QP accuracy: ',ACC(1)
+                  WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+                  WRITE(*,*) '##QP accuracy: ',ACC(1)
                 ENDIF
                 DO J=0,NSQUAREDSO
                   IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== 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	2015-06-30 16:20:47 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_WGTsq_le_10_QEDAmpAndQEDsq_gt_2.f	2015-07-24 15:59:43 +0000
@@ -420,9 +420,9 @@
 C        if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
         IF(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
           IF(NROTATIONS_DP.NE.0.OR.NROTATIONS_QP.NE.0) THEN
-            WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disable'
-     $       //'d, so MadLoop automatically resets NROTATIONS_DP an'
-     $       //'d NROTATIONS_QP to 0.'
+            WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING i'
+     $       //'s disabled, so MadLoop automatically resets NROTATIONS'
+     $       //'_DP and NROTATIONS_QP to 0.'
             NROTATIONS_QP=0
             NROTATIONS_DP=0
           ENDIF
@@ -518,10 +518,17 @@
             ENDDO
           ENDDO
         ENDDO
+<<<<<<< TREE
         IF(BOOTANDSTOP) THEN
           WRITE(*,*) 'Stopped by user request.'
           STOP
         ENDIF
+=======
+      ENDDO
+      IF(BOOTANDSTOP) THEN
+        WRITE(*,*) '##Stopped by user request.'
+        STOP
+>>>>>>> MERGE-SOURCE
       ENDIF
 
       IF(NTRY.EQ.0) THEN
@@ -9141,10 +9148,10 @@
      $       ,REF/DBLE(NCOMB),-1,-1)) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected wha'
      $         //'t are the zero helicity configurations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(29,FILE=HELFILTERFN,ERR=348)
@@ -9157,10 +9164,10 @@
      $       L(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected th'
      $         //'e helicity dependance relations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(30,FILE=HELFILTERFN,ERR=349)
@@ -9313,15 +9320,15 @@
                 WRITE(*,*) '##W03 WARNING An unstable PS point was'
      $           ,       ' detected.'
                 IF (NSQUAREDSO.NE.1) THEN
-                  WRITE(*,*) 'Accuracies for each split orde'
+                  WRITE(*,*) '##Accuracies for each split orde'
      $             //'r, starting with the summed case'
-                  WRITE(*,*) 'DP accuracies (for each split order): '
-     $             ,(TEMP1(I),I=0,NSQUAREDSO)
-                  WRITE(*,*) 'QP accuracies (for each split order): '
-     $             ,(ACC(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##DP accuracies (for each split orde'
+     $             //'r): ',(TEMP1(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##QP accuracies (for each split orde'
+     $             //'r): ',(ACC(I),I=0,NSQUAREDSO)
                 ELSE
-                  WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-                  WRITE(*,*) 'QP accuracy: ',ACC(1)
+                  WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+                  WRITE(*,*) '##QP accuracy: ',ACC(1)
                 ENDIF
                 DO J=0,NSQUAREDSO
                   IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== 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	2015-06-30 16:20:47 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_default.f	2015-07-24 15:59:43 +0000
@@ -420,9 +420,9 @@
 C        if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
         IF(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
           IF(NROTATIONS_DP.NE.0.OR.NROTATIONS_QP.NE.0) THEN
-            WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disable'
-     $       //'d, so MadLoop automatically resets NROTATIONS_DP an'
-     $       //'d NROTATIONS_QP to 0.'
+            WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING i'
+     $       //'s disabled, so MadLoop automatically resets NROTATIONS'
+     $       //'_DP and NROTATIONS_QP to 0.'
             NROTATIONS_QP=0
             NROTATIONS_DP=0
           ENDIF
@@ -518,10 +518,17 @@
             ENDDO
           ENDDO
         ENDDO
+<<<<<<< TREE
         IF(BOOTANDSTOP) THEN
           WRITE(*,*) 'Stopped by user request.'
           STOP
         ENDIF
+=======
+      ENDDO
+      IF(BOOTANDSTOP) THEN
+        WRITE(*,*) '##Stopped by user request.'
+        STOP
+>>>>>>> MERGE-SOURCE
       ENDIF
 
       IF(NTRY.EQ.0) THEN
@@ -2258,10 +2265,10 @@
      $       ,REF/DBLE(NCOMB),-1,-1)) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected wha'
      $         //'t are the zero helicity configurations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(29,FILE=HELFILTERFN,ERR=348)
@@ -2274,10 +2281,10 @@
      $       L(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected th'
      $         //'e helicity dependance relations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(30,FILE=HELFILTERFN,ERR=349)
@@ -2430,15 +2437,15 @@
                 WRITE(*,*) '##W03 WARNING An unstable PS point was'
      $           ,       ' detected.'
                 IF (NSQUAREDSO.NE.1) THEN
-                  WRITE(*,*) 'Accuracies for each split orde'
+                  WRITE(*,*) '##Accuracies for each split orde'
      $             //'r, starting with the summed case'
-                  WRITE(*,*) 'DP accuracies (for each split order): '
-     $             ,(TEMP1(I),I=0,NSQUAREDSO)
-                  WRITE(*,*) 'QP accuracies (for each split order): '
-     $             ,(ACC(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##DP accuracies (for each split orde'
+     $             //'r): ',(TEMP1(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##QP accuracies (for each split orde'
+     $             //'r): ',(ACC(I),I=0,NSQUAREDSO)
                 ELSE
-                  WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-                  WRITE(*,*) 'QP accuracy: ',ACC(1)
+                  WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+                  WRITE(*,*) '##QP accuracy: ',ACC(1)
                 ENDIF
                 DO J=0,NSQUAREDSO
                   IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== 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	2015-06-30 16:20:47 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDpert_default.f	2015-07-24 15:59:43 +0000
@@ -412,9 +412,9 @@
 C        if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
         IF(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
           IF(NROTATIONS_DP.NE.0.OR.NROTATIONS_QP.NE.0) THEN
-            WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disable'
-     $       //'d, so MadLoop automatically resets NROTATIONS_DP an'
-     $       //'d NROTATIONS_QP to 0.'
+            WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING i'
+     $       //'s disabled, so MadLoop automatically resets NROTATIONS'
+     $       //'_DP and NROTATIONS_QP to 0.'
             NROTATIONS_QP=0
             NROTATIONS_DP=0
           ENDIF
@@ -510,10 +510,17 @@
             ENDDO
           ENDDO
         ENDDO
+<<<<<<< TREE
         IF(BOOTANDSTOP) THEN
           WRITE(*,*) 'Stopped by user request.'
           STOP
         ENDIF
+=======
+      ENDDO
+      IF(BOOTANDSTOP) THEN
+        WRITE(*,*) '##Stopped by user request.'
+        STOP
+>>>>>>> MERGE-SOURCE
       ENDIF
 
       IF(NTRY.EQ.0) THEN
@@ -1304,10 +1311,10 @@
      $       ,REF/DBLE(NCOMB),-1,-1)) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected wha'
      $         //'t are the zero helicity configurations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(29,FILE=HELFILTERFN,ERR=348)
@@ -1320,10 +1327,10 @@
      $       L(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected th'
      $         //'e helicity dependance relations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(30,FILE=HELFILTERFN,ERR=349)
@@ -1476,15 +1483,15 @@
                 WRITE(*,*) '##W03 WARNING An unstable PS point was'
      $           ,       ' detected.'
                 IF (NSQUAREDSO.NE.1) THEN
-                  WRITE(*,*) 'Accuracies for each split orde'
+                  WRITE(*,*) '##Accuracies for each split orde'
      $             //'r, starting with the summed case'
-                  WRITE(*,*) 'DP accuracies (for each split order): '
-     $             ,(TEMP1(I),I=0,NSQUAREDSO)
-                  WRITE(*,*) 'QP accuracies (for each split order): '
-     $             ,(ACC(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##DP accuracies (for each split orde'
+     $             //'r): ',(TEMP1(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##QP accuracies (for each split orde'
+     $             //'r): ',(ACC(I),I=0,NSQUAREDSO)
                 ELSE
-                  WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-                  WRITE(*,*) 'QP accuracy: ',ACC(1)
+                  WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+                  WRITE(*,*) '##QP accuracy: ',ACC(1)
                 ENDIF
                 DO J=0,NSQUAREDSO
                   IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== 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	2015-06-30 16:20:47 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QEDpert_default.f	2015-07-24 15:59:43 +0000
@@ -412,9 +412,9 @@
 C        if AUTOMATIC_TIR_CACHE_CLEARING is disabled.
         IF(.NOT.AUTOMATIC_TIR_CACHE_CLEARING) THEN
           IF(NROTATIONS_DP.NE.0.OR.NROTATIONS_QP.NE.0) THEN
-            WRITE(*,*) 'INFO: AUTOMATIC_TIR_CACHE_CLEARING is disable'
-     $       //'d, so MadLoop automatically resets NROTATIONS_DP an'
-     $       //'d NROTATIONS_QP to 0.'
+            WRITE(*,*) '##INFO: AUTOMATIC_TIR_CACHE_CLEARING i'
+     $       //'s disabled, so MadLoop automatically resets NROTATIONS'
+     $       //'_DP and NROTATIONS_QP to 0.'
             NROTATIONS_QP=0
             NROTATIONS_DP=0
           ENDIF
@@ -510,10 +510,17 @@
             ENDDO
           ENDDO
         ENDDO
+<<<<<<< TREE
         IF(BOOTANDSTOP) THEN
           WRITE(*,*) 'Stopped by user request.'
           STOP
         ENDIF
+=======
+      ENDDO
+      IF(BOOTANDSTOP) THEN
+        WRITE(*,*) '##Stopped by user request.'
+        STOP
+>>>>>>> MERGE-SOURCE
       ENDIF
 
       IF(NTRY.EQ.0) THEN
@@ -1443,10 +1450,10 @@
      $       ,REF/DBLE(NCOMB),-1,-1)) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected wha'
      $         //'t are the zero helicity configurations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(29,FILE=HELFILTERFN,ERR=348)
@@ -1459,10 +1466,10 @@
      $       L(HELPICKED)+HELOFFSET)),REF,.TRUE.).EQ.0) THEN
               WRITE(*,*) '##W15 Helicity filter could not be successfu'
      $         //'lly double checked.'
-              WRITE(*,*) 'One reason for this is that you might hav'
+              WRITE(*,*) '##One reason for this is that you might hav'
      $         //'e changed sensible parameters which affected th'
      $         //'e helicity dependance relations.'
-              WRITE(*,*) 'MadLoop will try to reset the Helicit'
+              WRITE(*,*) '##MadLoop will try to reset the Helicit'
      $         //'y filter with the next PS points it receives.'
               NTRY=0
               OPEN(30,FILE=HELFILTERFN,ERR=349)
@@ -1615,15 +1622,15 @@
                 WRITE(*,*) '##W03 WARNING An unstable PS point was'
      $           ,       ' detected.'
                 IF (NSQUAREDSO.NE.1) THEN
-                  WRITE(*,*) 'Accuracies for each split orde'
+                  WRITE(*,*) '##Accuracies for each split orde'
      $             //'r, starting with the summed case'
-                  WRITE(*,*) 'DP accuracies (for each split order): '
-     $             ,(TEMP1(I),I=0,NSQUAREDSO)
-                  WRITE(*,*) 'QP accuracies (for each split order): '
-     $             ,(ACC(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##DP accuracies (for each split orde'
+     $             //'r): ',(TEMP1(I),I=0,NSQUAREDSO)
+                  WRITE(*,*) '##QP accuracies (for each split orde'
+     $             //'r): ',(ACC(I),I=0,NSQUAREDSO)
                 ELSE
-                  WRITE(*,*) 'DP accuracy: ',TEMP1(1)
-                  WRITE(*,*) 'QP accuracy: ',ACC(1)
+                  WRITE(*,*) '##DP accuracy: ',TEMP1(1)
+                  WRITE(*,*) '##QP accuracy: ',ACC(1)
                 ENDIF
                 DO J=0,NSQUAREDSO
                   IF (NSQUAREDSO.NE.1.OR.J.NE.0) THEN

=== renamed file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%MadLoopCommons.f' => 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%MadLoopCommons.f.THIS'
=== removed file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%MadLoopParamReader.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%MadLoopParamReader.f	2015-01-31 06:28:49 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%MadLoopParamReader.f	1970-01-01 00:00:00 +0000
@@ -1,270 +0,0 @@
-      subroutine MadLoopParamReader(filename, printParam) 
-
-      implicit none
-
-      CHARACTER(512) fileName, buff, buff2, mode
-      CHARACTER*10 MLReductionLib_str,MLReductionLib_str_save
-      CHARACTER*2  MLReductionLib_char
-      INTEGER MLRed,i,j,k
-
-      include "MadLoopParams.inc"
-
-      logical printParam, couldRead, paramPrinted, find
-      data paramPrinted/.FALSE./
-      couldRead=.False.
-!     Default parameters
-
-      open(666, file=fileName, err=676, status='OLD', action='READ')
-      do
-         read(666,*,end=999) buff
-         if(index(buff,'#').eq.1) then
-
-           if (buff .eq. '#CTModeInit') then
-             read(666,*,end=999) CTModeInit
-             if (CTModeInit .lt. 0 .or.
-     &           CTModeInit .gt. 6 ) then
-               stop 'CTModeInit must be >= 0 and <=6.'
-             endif
-           
-           else if (buff .eq. '#CTModeRun') then
-             read(666,*,end=999) CTModeRun
-             if (CTModeRun .lt. -1 .or.
-     &           CTModeRun .gt. 6 ) then
-               stop 'CTModeRun must be >= -1 and <=6.'
-             endif 
-
-           else if (buff .eq. '#NRotations_DP') then
-             read(666,*,end=999) NRotations_DP
-             if (NRotations_DP .lt. 0 .or.
-     &           NRotations_DP .gt. 2 ) then
-               stop 'NRotations_DP must be >= 0 and <=2.'
-             endif
-
-           else if (buff .eq. '#NRotations_QP') then
-             read(666,*,end=999) NRotations_QP
-             if (NRotations_QP .lt. 0 .or.
-     &           NRotations_QP .gt. 2 ) then
-               stop 'NRotations_QP must be >= 0 and <=2.'
-             endif
-
-           else if (buff .eq. '#MLStabThres') then
-             read(666,*,end=999) MLStabThres
-             if (MLStabThres.lt.0.0d0) then
-                 stop 'MLStabThres must be >= 0'
-             endif
-          
-           else if (buff .eq. '#CTLoopLibrary') then
-             read(666,*,end=999) CTLoopLibrary
-             if (CTLoopLibrary.lt.2 .or.
-     &           CTLoopLibrary.gt.3) then
-                 stop 'CTLoopLibrary must be >= 2 and <=3.'
-             endif
-
-           else if (buff .eq. '#CTStabThres') then
-             read(666,*,end=999) CTStabThres
-             if (CTStabThres.le.0.0d0) then
-                 stop 'CTStabThres must be > 0'
-             endif
-
-           else if (buff .eq. '#ZeroThres') then
-             read(666,*,end=999) ZeroThres
-             if (ZeroThres.le.0.0d0) then
-                 stop 'ZeroThres must be > 0'
-             endif
-
-           else if (buff .eq. '#OSThres') then
-              read(666,*,end=999) OSThres
-              if (OSThres.le.0.0d0) then
-                 stop 'OSThres must be > 0'
-              endif
-
-           else if (buff .eq. '#CheckCycle') then
-             read(666,*,end=999) CheckCycle
-             if (CheckCycle.lt.1) then
-                 stop 'CheckCycle must be >= 1'
-             endif
-
-           else if (buff .eq. '#MaxAttempts') then
-             read(666,*,end=999) MaxAttempts
-             if (MaxAttempts.lt.1) then
-                 stop 'MaxAttempts must be >= 1'
-             endif
-
-          else if (buff .eq. '#UseLoopFilter') then
-             read(666,*,end=999) UseLoopFilter
-
-          else if (buff .eq. '#DoubleCheckHelicityFilter') then
-             read(666,*,end=999) DoubleCheckHelicityFilter
-
-          else if (buff .eq. '#LoopInitStartOver') then
-             read(666,*,end=999) LoopInitStartOver
-
-          else if (buff .eq. '#HelInitStartOver') then
-             read(666,*,end=999) HelInitStartOver
-
-          else if (buff .eq. '#WriteOutFilters') then
-             read(666,*,end=999) WriteOutFilters
-
-          else if (buff .eq. '#ImprovePSPoint') then
-             read(666,*,end=999) ImprovePSPoint
-             if (ImprovePSPoint .lt. -1 .or.
-     &           ImprovePSPoint .gt. 2 ) then
-               stop 'ImprovePSPoint must be >= -1 and <=2.'
-             endif
-
-          else if (buff .eq. '#HelicityFilterLevel') then
-             read(666,*,end=999) HelicityFilterLevel
-             if (HelicityFilterLevel .lt. 0 .or.
-     &           HelicityFilterLevel .gt. 2 ) then
-               stop 'HelicityFilterLevel must be >= 0 and <=2.'
-             endif
-
-          else if (buff .eq. '#MLReductionLib') then
-             read(666,*,end=999) MLReductionLib_str
-             MLReductionLib(1:5)=0
-             MLReductionLib_str_save=MLReductionLib_str
-             j=0
-             DO
-                i=index(MLReductionLib_str,'|')
-                IF(i.EQ.0)THEN
-                   MLReductionLib_char=MLReductionLib_str
-                ELSE
-                   MLReductionLib_char=MLReductionLib_str(:i-1)
-                ENDIF
-                IF(MLReductionLib_char.EQ.'1 ')THEN
-                   MLRed=1
-                ELSEIF(MLReductionLib_char.EQ.'2 ')THEN
-                   MLRed=2
-                ELSEIF(MLReductionLib_char.EQ.'3 ')THEN
-                   MLRed=3
-                ELSEIF(MLReductionLib_char.EQ.'4 ')THEN
-                   MLRed=4
-                ELSE
-                   PRINT *, 'MLReductionLib is wrong: '//
-     $                  TRIM(MLReductionLib_str_save)
-                   STOP
-                ENDIF
-                find=.FALSE.
-                DO k=1,j
-                   IF(MLReductionLib(k).EQ.MLRed)THEN
-                      find=.TRUE.
-                      EXIT
-                   ENDIF
-                ENDDO
-                IF(.NOT.find)THEN
-                   j=j+1
-                   MLReductionLib(j)=MLRed
-                ENDIF
-                IF(i.EQ.0)THEN
-                   EXIT
-                ELSE
-                   MLReductionLib_str=MLReductionLib_str(i+1:)
-                ENDIF
-             ENDDO
-          else if (buff .eq. '#IREGIRECY') then
-             read(666,*,end=999) IREGIRECY
-          else if (buff .eq. '#IREGIMODE') then
-             read(666,*,end=999) IREGIMODE
-             if (IREGIMODE .lt. 0 .or.
-     &            IREGIMODE .gt.2) then
-                stop 'IREGIMODE must be >=0 and <=2.'
-             endif
-          else
-             write(*,*) 'The parameter name ',buff(2:),
-     &' is not reckognized.'
-             stop
-           endif
-
-         endif
-      enddo
-  999 continue
-      couldRead=.True.
-      goto 998      
-
-  676 continue
-      write(*,*) '##E00 Error:: MadLoop parameter file ',fileName,
-     &' could not be found or is malformed. Please specify it.'
-      stop 
-C     Below is the code if one desires to let the code continue with
-C     a non existing or malformed parameter file
-      write(*,*) '##I01 INFO :: The file ',fileName,' could not be ',
-     & ' open or did not contain the necessary information. The ',
-     & ' default MadLoop parameters will be used.'
-      call DefaultParam()
-      goto 999
-
-  998 continue
-
-      if(printParam.and..not.paramPrinted) then
-      write(*,*)
-     & '==============================================================='
-      if (couldRead) then      
-        write(*,*) 'INFO: MadLoop read these parameters from '
-     &,filename
-      else
-        write(*,*) 'INFO: MadLoop used the default parameters.'
-      endif
-      write(*,*)
-     & '==============================================================='
-      write(*,*) ' > MLReductionLib            = '
-     $     //TRIM(MLReductionLib_str_save)
-      write(*,*) ' > IREGIMODE                 = ',IREGIMODE
-      write(*,*) ' > IREGIRECY                 = ',IREGIRECY
-      write(*,*) ' > CTModeRun                 = ',CTModeRun
-      write(*,*) ' > MLStabThres               = ',MLStabThres
-      write(*,*) ' > NRotations_DP             = ',NRotations_DP
-      write(*,*) ' > NRotations_QP             = ',NRotations_QP
-      write(*,*) ' > CTStabThres               = ',CTStabThres
-      write(*,*) ' > CTLoopLibrary             = ',CTLoopLibrary
-      write(*,*) ' > CTModeInit                = ',CTModeInit
-      write(*,*) ' > CheckCycle                = ',CheckCycle
-      write(*,*) ' > MaxAttempts               = ',MaxAttempts
-      write(*,*) ' > UseLoopFilter             = ',UseLoopFilter
-      write(*,*) ' > HelicityFilterLevel       = ',HelicityFilterLevel
-      write(*,*) ' > ImprovePSPoint            = ',ImprovePSPoint
-      write(*,*) ' > DoubleCheckHelicityFilter = ',
-     &DoubleCheckHelicityFilter
-      write(*,*) ' > LoopInitStartOver         = ',LoopInitStartOver
-      write(*,*) ' > HelInitStartOver          = ',HelInitStartOver
-      write(*,*) ' > ZeroThres                 = ',ZeroThres
-      write(*,*) ' > OSThres                   = ',OSThres
-      write(*,*) ' > WriteOutFilters           = ',WriteOutFilters      
-      write(*,*)
-     & '==============================================================='
-      paramPrinted=.TRUE.
-      endif
-
-      close(666)
-
-      end
-
-      subroutine DefaultParam() 
-
-      implicit none
-      
-      include "MadLoopParams.inc"
-
-      MLReductionLib(1)=1
-      MLReductionLib(2:5)=0
-      IREGIMODE=2
-      IREGIRECY=.TRUE.
-      CTModeInit=0
-      CTModeRun=-1
-      NRotations_DP=1
-      NRotations_QP=0
-      MLStabThres=1.0d-3
-      CTStabThres=1.0d-2
-      CTLoopLibrary=3
-      CheckCycle=3
-      MaxAttempts=10
-      HelicityFilterLevel=2
-      UseLoopFilter=.False.
-      DoubleCheckHelicityFilter=.True.
-      LoopInitStartOver=.False.
-      HelInitStartOver=.False.
-      WriteOutFilters=.True.
-      ZeroThres=1.0d-9
-      OSThres=1.0d-13
-      ImprovePSPoint=2
-
-      end

=== removed file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%CT_interface.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%CT_interface.f	2014-10-20 01:07:00 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%CT_interface.f	1970-01-01 00:00:00 +0000
@@ -1,770 +0,0 @@
-      SUBROUTINE ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,RES,STABLE)
-C     
-C     Generated by MadGraph5_aMC@NLO v. %(version)s, %(date)s
-C     By the MadGraph5_aMC@NLO Development Team
-C     Visit launchpad.net/madgraph5 and amcatnlo.web.cern.ch
-C     
-C     Interface between MG5 and CutTools.
-C     
-C     Process: g g > t t~ [ virt = QCD ]
-C     
-C     
-C     CONSTANTS 
-C     
-      INTEGER    NEXTERNAL
-      PARAMETER (NEXTERNAL=4)
-      LOGICAL CHECKPCONSERVATION
-      PARAMETER (CHECKPCONSERVATION=.TRUE.)
-      REAL*8 NORMALIZATION
-      PARAMETER (NORMALIZATION = 1.D0/(16.D0*3.14159265358979323846D0*
-     $ *2))
-C     
-C     ARGUMENTS 
-C     
-      INTEGER NLOOPLINE, RANK
-      REAL*8 PL(0:3,NLOOPLINE)
-      REAL*8 PCT(0:3,0:NLOOPLINE-1)
-      COMPLEX*16 M2L(NLOOPLINE)
-      COMPLEX*16 M2LCT(0:NLOOPLINE-1)
-      COMPLEX*16 RES(3)
-      LOGICAL STABLE
-C     
-C     LOCAL VARIABLES 
-C     
-      COMPLEX*16 R1, ACC
-      INTEGER I, J, K
-      LOGICAL CTINIT, TIRINIT, GOLEMINIT
-      COMMON/REDUCTIONCODEINIT/CTINIT,TIRINIT,GOLEMINIT
-C     
-C     EXTERNAL FUNCTIONS
-C     
-      EXTERNAL ML5_0_LOOPNUM
-      EXTERNAL ML5_0_MPLOOPNUM
-C     
-C     GLOBAL VARIABLES
-C     
-      INCLUDE 'coupl.inc'
-      INTEGER CTMODE
-      REAL*8 LSCALE
-      COMMON/ML5_0_CT/LSCALE,CTMODE
-
-      INTEGER ID,SQSOINDEX,R
-      COMMON/ML5_0_LOOP/ID,SQSOINDEX,R
-
-C     ----------
-C     BEGIN CODE
-C     ----------
-
-C     INITIALIZE CUTTOOLS IF NEEDED
-      IF (CTINIT) THEN
-        CTINIT=.FALSE.
-        CALL ML5_0_INITCT()
-      ENDIF
-
-C     YOU CAN FIND THE DETAILS ABOUT THE DIFFERENT CTMODE AT THE
-C      BEGINNING OF THE FILE CTS_CUTS.F90 IN THE CUTTOOLS DISTRIBUTION
-
-C     CONVERT THE MASSES TO BE COMPLEX
-      DO I=1,NLOOPLINE
-        M2LCT(I-1)=M2L(I)
-      ENDDO
-
-C     CONVERT THE MOMENTA FLOWING IN THE LOOP LINES TO CT CONVENTIONS
-      DO I=0,3
-        DO J=0,(NLOOPLINE-1)
-          PCT(I,J)=0.D0
-        ENDDO
-      ENDDO
-      DO I=0,3
-        DO J=1,NLOOPLINE
-          PCT(I,0)=PCT(I,0)+PL(I,J)
-        ENDDO
-      ENDDO
-      IF (CHECKPCONSERVATION) THEN
-        IF (PCT(0,0).GT.1.D-6) THEN
-          WRITE(*,*) 'energy is not conserved ',PCT(0,0)
-          STOP 'energy is not conserved'
-        ELSEIF (PCT(1,0).GT.1.D-6) THEN
-          WRITE(*,*) 'px is not conserved ',PCT(1,0)
-          STOP 'px is not conserved'
-        ELSEIF (PCT(2,0).GT.1.D-6) THEN
-          WRITE(*,*) 'py is not conserved ',PCT(2,0)
-          STOP 'py is not conserved'
-        ELSEIF (PCT(3,0).GT.1.D-6) THEN
-          WRITE(*,*) 'pz is not conserved ',PCT(3,0)
-          STOP 'pz is not conserved'
-        ENDIF
-      ENDIF
-      DO I=0,3
-        DO J=1,(NLOOPLINE-1)
-          DO K=1,J
-            PCT(I,J)=PCT(I,J)+PL(I,K)
-          ENDDO
-        ENDDO
-      ENDDO
-
-      CALL CTSXCUT(CTMODE,LSCALE,MU_R,NLOOPLINE,ML5_0_LOOPNUM
-     $ ,ML5_0_MPLOOPNUM,RANK,PCT,M2LCT,RES,ACC,R1,STABLE)
-      RES(1)=NORMALIZATION*2.0D0*DBLE(RES(1))
-      RES(2)=NORMALIZATION*2.0D0*DBLE(RES(2))
-      RES(3)=NORMALIZATION*2.0D0*DBLE(RES(3))
-C     WRITE(*,*) 'Loop ID',ID,' =',RES(1),RES(2),RES(3)
-      END
-
-      SUBROUTINE ML5_0_INITCT()
-C     
-C     INITIALISATION OF CUTTOOLS
-C     
-C     LOCAL VARIABLES 
-C     
-      REAL*8 THRS
-      LOGICAL EXT_NUM_FOR_R1
-C     
-C     GLOBAL VARIABLES 
-C     
-      INCLUDE 'MadLoopParams.inc'
-C     ----------
-C     BEGIN CODE
-C     ----------
-
-C     DEFAULT PARAMETERS FOR CUTTOOLS
-C     -------------------------------  
-C     THRS1 IS THE PRECISION LIMIT BELOW WHICH THE MP ROUTINES
-C      ACTIVATES
-      THRS=CTSTABTHRES
-C     LOOPLIB SET WHAT LIBRARY CT USES
-C     1 -> LOOPTOOLS
-C     2 -> AVH
-C     3 -> QCDLOOP
-      LOOPLIB=CTLOOPLIBRARY
-C     MADLOOP'S NUMERATOR IN THE OPEN LOOP IS MUCH FASTER THAN THE
-C      RECONSTRUCTED ONE IN CT. SO WE BETTER USE MADLOOP ONE IN THIS
-C      CASE.
-      EXT_NUM_FOR_R1=.TRUE.
-C     -------------------------------	  
-
-C     The initialization below is for CT v1.8.+
-      CALL CTSINIT(THRS,LOOPLIB,EXT_NUM_FOR_R1)
-C     The initialization below is for the older stable CT v1.7, still
-C      used for now in the beta release.
-C     CALL CTSINIT(THRS,LOOPLIB)
-
-      END
-
-      SUBROUTINE ML5_0_LOOP_2_3(P1, P2, W1, W2, W3, M1, M2,  RANK
-     $ , SQUAREDSOINDEX, LOOPNUM)
-      INTEGER    NEXTERNAL
-      PARAMETER (NEXTERNAL=4)
-      INTEGER    NLOOPLINE
-      PARAMETER (NLOOPLINE=2)
-      INTEGER    NWAVEFUNCS
-      PARAMETER (NWAVEFUNCS=10)
-      INTEGER    NLOOPGROUPS
-      PARAMETER (NLOOPGROUPS=26)
-      INTEGER    NCOMB
-      PARAMETER (NCOMB=16)
-C     These are constants related to the split orders
-      INTEGER    NSQUAREDSO
-      PARAMETER (NSQUAREDSO=1)
-C     
-C     ARGUMENTS 
-C     
-      INTEGER W1, W2, W3
-      COMPLEX*16 M1, M2
-      INTEGER P1, P2
-      INTEGER RANK, LSYMFACT
-      INTEGER LOOPNUM, SQUAREDSOINDEX
-C     
-C     LOCAL VARIABLES 
-C     
-      REAL*8 PL(0:3,NLOOPLINE)
-      COMPLEX*16 M2L(NLOOPLINE)
-      INTEGER PAIRING(NLOOPLINE),WE(3)
-      INTEGER I, J, K, TEMP,I_LIB
-      LOGICAL COMPLEX_MASS,DOING_QP
-C     
-C     GLOBAL VARIABLES
-C     
-      INCLUDE 'MadLoopParams.inc'
-      INTEGER ID,SQSOINDEX,R
-      COMMON/ML5_0_LOOP/ID,SQSOINDEX,R
-
-      LOGICAL CHECKPHASE, HELDOUBLECHECKED
-      COMMON/ML5_0_INIT/CHECKPHASE, HELDOUBLECHECKED
-
-      INTEGER HELOFFSET
-      INTEGER GOODHEL(NCOMB)
-      LOGICAL GOODAMP(NSQUAREDSO,NLOOPGROUPS)
-      COMMON/ML5_0_FILTERS/GOODAMP,GOODHEL,HELOFFSET
-
-      COMPLEX*16 LOOPRES(3,NSQUAREDSO,NLOOPGROUPS)
-      LOGICAL S(NSQUAREDSO,NLOOPGROUPS)
-      COMMON/ML5_0_LOOPRES/LOOPRES,S
-
-
-      COMPLEX*16 W(20,NWAVEFUNCS)
-      COMMON/ML5_0_W/W
-      REAL*8 LSCALE
-      INTEGER CTMODE
-      COMMON/ML5_0_CT/LSCALE,CTMODE
-      INTEGER LIBINDEX
-      COMMON/ML5_0_I_LIB/LIBINDEX
-
-C     ----------
-C     BEGIN CODE
-C     ----------
-
-      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
-     $ EX,LOOPNUM)) THEN
-        WE(1)=W1
-        WE(2)=W2
-        WE(3)=W3
-        M2L(1)=M2**2
-        M2L(2)=M1**2
-        PAIRING(1)=P1
-        PAIRING(2)=P2
-        R=RANK
-        ID=LOOPNUM
-        SQSOINDEX=SQUAREDSOINDEX
-        DO I=0,3
-          TEMP=1
-          DO J=1,NLOOPLINE
-            PL(I,J)=0.D0
-            DO K=TEMP,(TEMP+PAIRING(J)-1)
-              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
-            ENDDO
-            TEMP=TEMP+PAIRING(J)
-          ENDDO
-        ENDDO
-C       Determine whether the integral is with complex masses or not
-C       since some reduction libraries, e.g.PJFry++ and IREGI, are
-C        still
-C       not able to deal with complex masses
-        COMPLEX_MASS=.FALSE.
-        DO I=1,NLOOPLINE
-          IF(DIMAG(M2L(I)).EQ.0D0)CYCLE
-          IF(ABS(DIMAG(M2L(I)))/MAX(ABS(M2L(I)),1D-2).GT.1D-15)THEN
-            COMPLEX_MASS=.TRUE.
-            EXIT
-          ENDIF
-        ENDDO
-C       Determine it uses qp or not
-        DOING_QP=.FALSE.
-        IF(CTMODE.GE.4)DOING_QP=.TRUE.
-C       Choose the correct loop library
-        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
-     $   ,DOING_QP,I_LIB)
-        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
-C         USE CUTTOOLS
-          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
-     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
-        ELSE
-C         USE TIR
-          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
-     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
-     $     ,LOOPNUM))
-        ENDIF
-      ELSE
-        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
-        LOOPRES(2,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
-        LOOPRES(3,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
-        S(SQUAREDSOINDEX,LOOPNUM)=.TRUE.
-      ENDIF
-      END
-
-      SUBROUTINE ML5_0_LOOP_4(W1, W2, W3, W4, M1, M2, M3, M4,  RANK
-     $ , SQUAREDSOINDEX, LOOPNUM)
-      INTEGER    NEXTERNAL
-      PARAMETER (NEXTERNAL=4)
-      INTEGER    NLOOPLINE
-      PARAMETER (NLOOPLINE=4)
-      INTEGER    NWAVEFUNCS
-      PARAMETER (NWAVEFUNCS=10)
-      INTEGER    NLOOPGROUPS
-      PARAMETER (NLOOPGROUPS=26)
-      INTEGER    NCOMB
-      PARAMETER (NCOMB=16)
-C     These are constants related to the split orders
-      INTEGER    NSQUAREDSO
-      PARAMETER (NSQUAREDSO=1)
-C     
-C     ARGUMENTS 
-C     
-      INTEGER W1, W2, W3, W4
-      COMPLEX*16 M1, M2, M3, M4
-
-      INTEGER RANK, LSYMFACT
-      INTEGER LOOPNUM, SQUAREDSOINDEX
-C     
-C     LOCAL VARIABLES 
-C     
-      REAL*8 PL(0:3,NLOOPLINE)
-      COMPLEX*16 M2L(NLOOPLINE)
-      INTEGER PAIRING(NLOOPLINE),WE(4)
-      INTEGER I, J, K, TEMP,I_LIB
-      LOGICAL COMPLEX_MASS,DOING_QP
-C     
-C     GLOBAL VARIABLES
-C     
-      INCLUDE 'MadLoopParams.inc'
-      INTEGER ID,SQSOINDEX,R
-      COMMON/ML5_0_LOOP/ID,SQSOINDEX,R
-
-      LOGICAL CHECKPHASE, HELDOUBLECHECKED
-      COMMON/ML5_0_INIT/CHECKPHASE, HELDOUBLECHECKED
-
-      INTEGER HELOFFSET
-      INTEGER GOODHEL(NCOMB)
-      LOGICAL GOODAMP(NSQUAREDSO,NLOOPGROUPS)
-      COMMON/ML5_0_FILTERS/GOODAMP,GOODHEL,HELOFFSET
-
-      COMPLEX*16 LOOPRES(3,NSQUAREDSO,NLOOPGROUPS)
-      LOGICAL S(NSQUAREDSO,NLOOPGROUPS)
-      COMMON/ML5_0_LOOPRES/LOOPRES,S
-
-
-      COMPLEX*16 W(20,NWAVEFUNCS)
-      COMMON/ML5_0_W/W
-      REAL*8 LSCALE
-      INTEGER CTMODE
-      COMMON/ML5_0_CT/LSCALE,CTMODE
-      INTEGER LIBINDEX
-      COMMON/ML5_0_I_LIB/LIBINDEX
-
-C     ----------
-C     BEGIN CODE
-C     ----------
-
-      IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED).OR.GOODAMP(SQUAREDSOIND
-     $ EX,LOOPNUM)) THEN
-        WE(1)=W1
-        WE(2)=W2
-        WE(3)=W3
-        WE(4)=W4
-        M2L(1)=M4**2
-        M2L(2)=M1**2
-        M2L(3)=M2**2
-        M2L(4)=M3**2
-        DO I=1,NLOOPLINE
-          PAIRING(I)=1
-        ENDDO
-
-        R=RANK
-        ID=LOOPNUM
-        SQSOINDEX=SQUAREDSOINDEX
-        DO I=0,3
-          TEMP=1
-          DO J=1,NLOOPLINE
-            PL(I,J)=0.D0
-            DO K=TEMP,(TEMP+PAIRING(J)-1)
-              PL(I,J)=PL(I,J)-DBLE(W(1+I,WE(K)))
-            ENDDO
-            TEMP=TEMP+PAIRING(J)
-          ENDDO
-        ENDDO
-C       Determine whether the integral is with complex masses or not
-C       since some reduction libraries, e.g.PJFry++ and IREGI, are
-C        still
-C       not able to deal with complex masses
-        COMPLEX_MASS=.FALSE.
-        DO I=1,NLOOPLINE
-          IF(DIMAG(M2L(I)).EQ.0D0)CYCLE
-          IF(ABS(DIMAG(M2L(I)))/MAX(ABS(M2L(I)),1D-2).GT.1D-15)THEN
-            COMPLEX_MASS=.TRUE.
-            EXIT
-          ENDIF
-        ENDDO
-C       Determine it uses qp or not
-        DOING_QP=.FALSE.
-        IF(CTMODE.GE.4)DOING_QP=.TRUE.
-C       Choose the correct loop library
-        CALL ML5_0_CHOOSE_LOOPLIB(LIBINDEX,NLOOPLINE,RANK,COMPLEX_MASS
-     $   ,DOING_QP,I_LIB)
-        IF(MLREDUCTIONLIB(I_LIB).EQ.1)THEN
-C         USE CUTTOOLS
-          CALL ML5_0_CTLOOP(NLOOPLINE,PL,M2L,RANK,LOOPRES(1,SQUAREDSOIN
-     $     DEX,LOOPNUM),S(SQUAREDSOINDEX,LOOPNUM))
-        ELSE
-C         USE TIR
-          CALL ML5_0_TIRLOOP(SQUAREDSOINDEX,LOOPNUM,I_LIB,NLOOPLINE,PL
-     $     ,M2L,RANK,LOOPRES(1,SQUAREDSOINDEX,LOOPNUM),S(SQUAREDSOINDEX
-     $     ,LOOPNUM))
-        ENDIF
-      ELSE
-        LOOPRES(1,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
-        LOOPRES(2,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
-        LOOPRES(3,SQUAREDSOINDEX,LOOPNUM)=(0.0D0,0.0D0)
-        S(SQUAREDSOINDEX,LOOPNUM)=.TRUE.
-      ENDIF
-      END
-
-      SUBROUTINE ML5_0_LOOP_3_4(P1, P2, P3, W1, W2, W3, W4, M1, M2, M3
-     $ ,  RANK, SQUAREDSOINDEX, LOOPNUM)
-      INTEGER    NEXTERNAL
-      PARAMETER (NEXTERNAL=4)
-      INTEGER    NLOOPLINE
-      PARAMETER (NLOOPLINE=3)
-      INTEGER    NWAVEFUNCS
-      PARAMETER (NWAVEFUNCS=10)
-      INTEGER    NLOOPGROUPS
-      PARAMETER (NLOOPGROUPS=26)
-      INTEGER    NCOMB
-      PARAMETER (NCOMB=16)
-C     These are constants related to the split orders
-      INTEGER    NSQUAREDSO
-      PARAMETER (NSQUAREDSO=1)
-C     
-C     ARGUMENTS 
-C     
-      INTEGER W1, W2, W3, W4
-      COMPLEX*16 M1, M2, M3
-      INTEGER P1, P2, P3
-      INTEGER RANK, LSYMFACT
-      INTEGER LOOPNUM, SQUAREDSOINDEX
-C     
-C     LOCAL VARIABLES 
-C     
-      REAL*8 PL(0:3,NLOOPLINE)
-      COMPLEX*16 M2L(NLOOPLINE)
-      INTEGER PAIRING(NLOOPLINE),WE(4)
-      INTEGER I, J, K, TEMP,I_LIB
-      LOGICAL COMPLEX_MASS,DOING_QP
-C     
-C     GLOBAL VARIABLES
-C     
-      INCLUDE 'MadLoopParams.inc'
-      INTEGER ID,SQSOINDEX,R
-      COMMON/ML5_0_LOOP/ID,SQSOINDEX,R
-
-      LOGICAL CHECKPHASE, HELDOUBLECHECKED
-      COMMON/ML5_0_INIT/CHECKPHASE, HELDOUBLECHECKED
-
-      INTEGER HELOFFSET
-      INTEGER GOODHEL(NCOMB)
-      LOGICAL GOODAMP(NSQUAREDSO,NLOOPGROUPS)
-      COMMON/ML5_0_FILTERS/GOODAMP,GOODHEL,HELOFFSET
-
-      COMPLEX*16 LOOPRES(3,NSQUAREDSO,NLOOPGROUPS)
-      LOGICAL S(NSQUAREDSO,NLOOPGROUPS)
-      COMMON/ML5_0_LOOPRES/LOOPRES,S
-
-
-      COMPLEX*16 W(20,NWAVEFUNCS)
-      COMMON/ML5_0_W/W
-      REAL*8 LSCALE
-      INTEGER CTMODE