Merge lp:~maddevelopers/mg5amcnlo/dim6_eft into lp:~mg5core1/mg5amcnlo/2.6.6

[Olivier Mattelaer]

change in the reweighting for the EFT validation effort

[lp:~maddevelopers/mg5amcnlo/dim6_eft]

297. By olivier-mattelaer on 2019-05-16

merge with latest 2.6.6

298. By olivier-mattelaer on 2019-05-16

some cleaning + adapting test for new allowed input in RAMBO

[marco zaro (marco-zaro)]

Hi Olivier,
can you elaborate a bit more the relevant changes of this merge proposal?
thanks!
cheers,

Marco

On Thu, May 16, 2019 at 11:38 AM Olivier Mattelaer <
<email address hidden>> wrote:

> Olivier Mattelaer has proposed merging
> lp:~maddevelopers/mg5amcnlo/dim6_eft into lp:~mg5core1/mg5amcnlo/2.6.6.
>
> Commit message:
> change in the reweighting for the EFT validation effort
>
> Requested reviews:
> MadTeam (madteam)
>
> For more details, see:
> https://code.launchpad.net/~maddevelopers/mg5amcnlo/dim6_eft/+merge/367512
> --
> Your team mg5core1 is subscribed to branch lp:~mg5core1/mg5amcnlo/2.6.6.
>

[Olivier Mattelaer (olivier-mattelaer)]

The main point for this branch is to do some core change related to the
EFT validation effort (lead by Gauthier/Fabio). This is implemented as a plugin but the accord required some change in the re-weighting module (mainly for the formatting).

Those changes adds:
1) The following command: import process PATH is now supported. (allowing to provide a file with all your process definition whithout any generate/add process
2) Allow the re-weighting to distinguish process based on the process event written in the lhe file / @X [does not seems to work for loop]
3) More xml like formatting for the weight definition (including model version/...)
4) the support for NHEL=-2 in reweighting, which returns the sum over polarization (like -1) but clear all filter associated to the zero-helicity to avoid any possible bias)
5) a (not documented) flag --no_warning=nodiagram which prevent the code to crash if you do not have any diagram generated for a given process
6) possibility to call the reweighting interface as a standalone executable
7) improve the templating of the standalone output to support model with more than 1k parameters (and reduce RAM use for model with less parameters)

[lp:~maddevelopers/mg5amcnlo/dim6_eft]

299. By olivier-mattelaer on 2019-05-23

fix latest acceptance test

Unmerged revisions

299. By olivier-mattelaer on 2019-05-23

fix latest acceptance test

298. By olivier-mattelaer on 2019-05-16

some cleaning + adapting test for new allowed input in RAMBO

297. By olivier-mattelaer on 2019-05-16

merge with latest 2.6.6

296. By olivier-mattelaer on 2019-04-26

add reweight_card

295. By olivier-mattelaer on 2019-04-19

add more model information

294. By olivier-mattelaer on 2019-04-18

adding support for nincoming in standalone fortran as py lib

293. By olivier-mattelaer on 2019-04-18

merge with last 2.6.6

292. By olivier-mattelaer on 2019-01-17

avoid any issue with parameter with many input parameter in standalone mode

291. By olivier-mattelaer on 2019-01-13

fix issue with multiple launch in re=wgt outside of mg

290. By olivier-mattelaer on 2019-01-11

fix an issue related to NLO model for import process command

=== modified file 'madgraph/core/base_objects.py'
--- madgraph/core/base_objects.py	2018-12-03 22:10:49 +0000
+++ madgraph/core/base_objects.py	2019-05-23 07:08:37 +0000
@@ -1043,6 +1043,8 @@
         self['expansion_order'] = None
         self['version_tag'] = None # position of the directory (for security)
         self['gauge'] = [0, 1]
+        import collections
+        self['model_info'] = collections.defaultdict(str)
         self['case_sensitive'] = True
         # attribute which might be define if needed
         #self['name2pdg'] = {'name': pdg}
@@ -1184,7 +1186,7 @@
 
         if (name == 'order_hierarchy') and not self[name]:
             if self['interactions']:
-                self['order_hierarchy'] = self.get_order_hierarchy()    
+                self['order_hierarchy'] = self.get_order_hierarchy()  
 
         if (name == 'expansion_order') and self[name] == None:
             if self['interactions']:

=== modified file 'madgraph/interface/common_run_interface.py'
--- madgraph/interface/common_run_interface.py	2019-05-03 07:18:01 +0000
+++ madgraph/interface/common_run_interface.py	2019-05-23 07:08:37 +0000
@@ -5211,7 +5211,6 @@
     def do_set(self, line):
         """ edit the value of one parameter in the card"""
         
-        
         args = self.split_arg(line)
         
         
@@ -6696,6 +6695,7 @@
 
         if card_name == 'unknown':
             logger.warning('Fail to determine the type of the file. Not copied')
+            return
         if card_name != 'banner':
             logger.info('copy %s as %s' % (pathname, card_name))
             files.cp(path, self.paths[card_name.rsplit('_',1)[0]])

=== modified file 'madgraph/interface/extended_cmd.py'
--- madgraph/interface/extended_cmd.py	2019-04-24 11:02:09 +0000
+++ madgraph/interface/extended_cmd.py	2019-05-23 07:08:37 +0000
@@ -1652,7 +1652,6 @@
         # filepath can be overwritten during the run (leading to weird results)
         # Note also that we need a generator and not a list.
         for line in self.inputfile:
-            
             #remove pointless spaces and \n
             line = line.replace('\n', '').strip()
             # execute the line

=== modified file 'madgraph/interface/madgraph_interface.py'
--- madgraph/interface/madgraph_interface.py	2019-05-02 22:31:26 +0000
+++ madgraph/interface/madgraph_interface.py	2019-05-23 07:08:37 +0000
@@ -339,6 +339,9 @@
         logger.info("")
         logger.info("   import banner PATH  [--no_launch]:",'$MG:BOLD')
         logger.info("      Rerun the exact same run define in the valid banner.")
+        logger.info("")
+        logger.info("   import process PATH", '$MG:BOLD')
+        logger.info("      generate the list of processes defined in the file")
 
     def help_install(self):
         logger.info("syntax: install " + "|".join(self._install_opts),'$MG:color:BLUE')
@@ -2786,7 +2789,7 @@
     _switch_opts = ['mg5','aMC@NLO','ML5']
     _check_opts = ['full', 'timing', 'stability', 'profile', 'permutation',
                    'gauge','lorentz', 'brs', 'cms']
-    _import_formats = ['model_v4', 'model', 'proc_v4', 'command', 'banner']
+    _import_formats = ['model_v4', 'model', 'proc_v4', 'command', 'banner', 'process']
     _install_opts = ['Delphes', 'MadAnalysis4', 'ExRootAnalysis',
                      'update', 'Golem95', 'PJFry', 'QCDLoop', 'maddm', 'maddump',
                      'looptools']
@@ -3008,8 +3011,13 @@
         standalone_only = False
         if '--standalone' in args:
             standalone_only = True
-            args.remove('--standalone')            
-
+            args.remove('--standalone')   
+            
+        allow_no_diagram = False         
+        if '--no_warning=nodiagram' in args:
+            allow_no_diagram = True
+            args.remove('--no_warning=nodiagram')
+            
         # Check the validity of the arguments
         self.check_add(args)
 
@@ -3106,13 +3114,18 @@
                            self.options['ignore_six_quark_processes'] if \
                            "ignore_six_quark_processes" in self.options \
                            else []
-
-            myproc = diagram_generation.MultiProcess(myprocdef,
+            try:
+                myproc = diagram_generation.MultiProcess(myprocdef,
                                      collect_mirror_procs = collect_mirror_procs,
                                      ignore_six_quark_processes = ignore_six_quark_processes,
                                      optimize=optimize, diagram_filter=diagram_filter)
-
-
+            except diagram_generation.NoDiagramException:
+                if not allow_no_diagram:
+                    raise
+                else:
+                    logger.warning('no diagram generated. No error raised due to flag')
+                    return 
+                
             for amp in myproc.get('amplitudes'):
                 if amp not in self._curr_amps:
                     self._curr_amps.append(amp)
@@ -4732,9 +4745,10 @@
             split_orders=list(set(perturbation_couplings_list+squared_orders.keys()))
             try:
                 split_orders.sort(key=lambda elem: 0 if elem=='WEIGHTED' else
-                                       self._curr_model['order_hierarchy']
+                                       self._curr_model.get('order_hierarchy')
                                        [elem if not elem.endswith('.sqrt') else elem[:-5]])
             except KeyError:
+                misc.sprint(self._curr_model['order_hierarchy'], self._curr_model.get('modelpath'))
                 raise self.InvalidCmd, "The loaded model does not defined a "+\
                     " coupling order hierarchy for these couplings: %s"%\
                       str([so for so in split_orders if so!='WEIGHTED' and so not 
@@ -5234,6 +5248,7 @@
                 if os.path.sep in args[1] and "import" in self.history[-1]:
                     self.history[-1] = 'import model %s' % self._curr_model.get('modelpath+restriction')
 
+
                 if self.options['gauge']=='unitary':
                     if not force and isinstance(self._curr_model,\
                                               loop_base_objects.LoopModel) and \
@@ -5285,7 +5300,23 @@
                 self.check_for_export_dir(args[1])
                 # Execute the card
                 self.import_command_file(args[1])
-
+        elif args[0] == 'process':
+            if not os.path.isfile(args[1]):
+                raise self.InvalidCmd("Path %s is not a valid pathname" % args[1])
+            else:   
+                first = True         
+                for line in open(args[1]):
+                    line = line.split('#')[0].strip()
+                    if not line:
+                        continue
+                    if first:
+                        self.exec_cmd('generate %s' % line, precmd=True)
+                        first = False
+                    else:
+                        self.exec_cmd('add process %s' % line, precmd=True)
+                        
+                        
+                        
         elif args[0] == 'banner':
             type = madevent_interface.MadEventCmd.detect_card_type(args[1])
             if type != 'banner':
@@ -5456,8 +5487,6 @@
                                         (line.split()[0],why))
                 if self.history[-1] == 'define %s' % line.strip():
                     self.history.pop(-1)
-                else:
-                    misc.sprint([self.history[-1], 'define %s' % line.strip()])
 
         scheme = "old"
         for qcd_container in ['p', 'j']:
@@ -6011,7 +6040,7 @@
                 else:
                     othersource = 'uiuc'
                 # try with the mirror
-                misc.sprint('try other mirror', othersource, ' '.join(args))
+                logger.debug('try other mirror', othersource, ' '.join(args))
                 return self.do_install('%s --source=%s' % (' '.join(args), othersource), 
                                        paths, additional_options) 
             else:
@@ -6859,7 +6888,6 @@
         self.check_launch(args, options)
         options = options.__dict__
         # args is now MODE PATH
-
         if args[0].startswith('standalone'):
             if os.path.isfile(os.path.join(os.getcwd(),args[1],'Cards',\
               'MadLoopParams.dat')) and not os.path.isfile(os.path.join(\

=== modified file 'madgraph/interface/master_interface.py'
--- madgraph/interface/master_interface.py	2018-04-02 14:02:36 +0000
+++ madgraph/interface/master_interface.py	2019-05-23 07:08:37 +0000
@@ -272,6 +272,10 @@
 
     def do_import(self, *args, **opts):
         self.cmd.do_import(self, *args, **opts)
+        line = args[0]
+        import_mode = args[0].split()[0]
+        if import_mode in self._import_formats and not import_mode.startswith('model'):
+            return
         if self._curr_model:
             if isinstance(self._curr_model, loop_base_objects.LoopModel) and \
                self._curr_model['perturbation_couplings']!=[] and \
@@ -610,6 +614,7 @@
         
     def change_principal_cmd(self, name):
 
+
         old_cmd=self.current_interface
         if name in self.interface_names.keys():
             self.prompt= self.interface_names[name][0]+'>'

=== modified file 'madgraph/interface/reweight_interface.py'
--- madgraph/interface/reweight_interface.py	2019-05-09 10:05:05 +0000
+++ madgraph/interface/reweight_interface.py	2019-05-23 07:08:37 +0000
@@ -1,3 +1,4 @@
+#! /usr/bin/env python
 ################################################################################
 #
 # Copyright (c) 2009 The MadGraph5_aMC@NLO Development team and Contributors
@@ -26,8 +27,11 @@
 import subprocess
 from subprocess import Popen, PIPE, STDOUT
 
-
 pjoin = os.path.join
+if "__main__" == __name__:
+    local_path = os.path.dirname(os.path.realpath( __file__ ))
+    root_path = pjoin(local_path, '..','..')
+    sys.path.insert(0, root_path)
 
 import madgraph.interface.extended_cmd as extended_cmd
 import madgraph.interface.madgraph_interface as mg_interface
@@ -66,13 +70,12 @@
     
     prompt = 'Reweight>'
     debug_output = 'Reweight_debug'
-    rwgt_dir_possibility = ['rw_me','rw_me_second','rw_mevirt','rw_mevirt_second']
+    #rwgt_dir_possibility = ['rw_me','rw_me_second','rw_mevirt','rw_mevirt_second']
     
     @misc.mute_logger()
     def __init__(self, event_path=None, allow_madspin=False, mother=None, *completekey, **stdin):
         """initialize the interface with potentially an event_path"""
         
-        
         self.me_dir = os.getcwd()
         if not event_path:
             cmd_logger.info('************************************************************')
@@ -95,6 +98,7 @@
         self.f2pylib = {}
         self.second_model = None
         self.second_process = None
+        self.nb_library = 1
         self.dedicated_path = {}
         self.soft_threshold = None
         self.systematics = False # allow to run systematics in ouput2.0 mode
@@ -154,6 +158,7 @@
         if not self.lhe_input.banner:
             value = self.ask("What is the path to banner", 0, [0], "please enter a path", timeout=0)
             self.lhe_input.banner = open(value).read()
+
         self.banner = self.lhe_input.get_banner()
         
         #get original cross-section/error
@@ -361,6 +366,7 @@
             if self.has_standalone_dir:
                 self.terminate_fortran_executables()
                 self.has_standalone_dir = False
+
             if args[-1] == "--add":
                 self.second_process.append(" ".join(args[1:-1]))
             else:
@@ -465,6 +471,16 @@
             self.options['rwgt_name'] = opts['rwgt_name']
 
         model_line = self.banner.get('proc_card', 'full_model_line')
+        if not self.mother and not hasattr(self, 'output'):
+            name, ext = self.lhe_input.name.rsplit('.',1)
+            target = '%s_out.%s' % (name, ext)
+            
+            if os.path.exists(target) and os.path.getmtime(target) > os.path.getmtime(self.lhe_input.name):
+                nb_launch = len([1 for l in self.history if l.strip().startswith('launch')])
+                if nb_launch > 1:
+                    files.mv(target, self.lhe_input.name)
+                    self.exec_cmd('import %s' % self.lhe_input.name)
+
 
         if not self.has_standalone_dir:                           
             if self.rwgt_dir and os.path.exists(pjoin(self.rwgt_dir,'rw_me','rwgt.pkl')):
@@ -506,7 +522,7 @@
             path_me = self.me_dir 
             
         if self.second_model or self.second_process or self.dedicated_path:
-            rw_dir = pjoin(path_me, 'rw_me_second')
+            rw_dir = pjoin(path_me, 'rw_me_%s' % self.nb_library)
         else:
             rw_dir = pjoin(path_me, 'rw_me')
                 
@@ -566,7 +582,7 @@
             if self.output_type == "default":
                 for name in weight:
                     if 'orig' in name:
-                        continue             
+                        continue   
                     event.reweight_data['%s%s' % (tag_name,name)] = weight[name]
                     #write this event with weight
                 output.write(str(event))
@@ -583,7 +599,7 @@
                     output[(tag_name,name)].write(str(new_evt))
 
         # check normalisation of the events:
-        if 'event_norm' in self.run_card:
+        if self.run_card and 'event_norm' in self.run_card:
             if self.run_card['event_norm'] in ['average','bias']:
                 for key, value in cross.items():
                     cross[key] = value / (event_nb+1)
@@ -632,14 +648,16 @@
             import madgraph.interface.madevent_interface as ME_interface
 
         self.lhe_input.close()
-        if not self.mother:
+        if hasattr(self, 'output_path'):
+            target = self.output_path
+        elif not self.mother:
             name, ext = self.lhe_input.name.rsplit('.',1)
             target = '%s_out.%s' % (name, ext)            
         elif 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)
             logger.info('Event %s have now the additional weight' % self.lhe_input.name)
@@ -694,7 +712,7 @@
             path_me = self.me_dir 
             
         if self.second_model or self.second_process or self.dedicated_path:
-            rw_dir = pjoin(path_me, 'rw_me_second')
+            rw_dir = pjoin(path_me, 'rw_me_%s' % self.nb_library)
         else:
             rw_dir = pjoin(path_me, 'rw_me')
         
@@ -713,6 +731,8 @@
                 files.ln(ff.name, starting_dir=pjoin(path_me, 'rw_mevirt', 'Cards'))
             cmd = common_run_interface.CommonRunCmd.ask_edit_card_static(cards=['param_card.dat'],
                                    ask=self.ask, pwd=rw_dir, first_cmd=self.stored_line)
+
+            
             self.stored_line = None
         
         # check for potential scan in the new card 
@@ -751,9 +771,8 @@
                 blockpat = re.compile(r'''<weightgroup name=\'mg_reweighting\'\s*weight_name_strategy=\'includeIdInWeightName\'>(?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=\'(?:rwgt_(?P<id>\d+)|(?P<id2>[_\w]+))(?P<rwgttype>\s*|_\w+)\'>(?P<info>.*?)</weight>', re.S+re.I+re.M)
+                pattern = re.compile('<weight id=\'(?:rwgt_(?P<id>\d+)|(?P<id2>[_\w\-]+))(?P<rwgttype>\s*|_\w+)\'>(?P<info>.*?)</weight>', re.S+re.I+re.M)
                 mg_rwgt_info = pattern.findall(content)
-                
                 maxid = 0
                 for k,(i, fulltag, nlotype, diff) in enumerate(mg_rwgt_info):
                     if i:
@@ -796,12 +815,21 @@
         else:
             tag = str(rewgtid)
         
+        version = misc.get_pkg_info()['version']
+        try:
+            model_version = self.mg5cmd._curr_model['model_info']['version'] 
+        except TypeError:
+            model_version = '0.0.0'
+        str_version = "<weight_generator>MG5aMC_v%s</weight_generator>\n<weight_model_version>%s</weight_model_version>" %\
+            (version, model_version)
+        
         if not self.second_model and not self.dedicated_path:
             old_param = check_param_card.ParamCard(s_orig.splitlines())
             new_param =  self.new_param_card
             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.')
+            card_diff = "%s\n<weight_slha>%s</weight_slha>" % (str_version,card_diff)
             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.")
@@ -813,6 +841,7 @@
                     mg_rwgt_info.append((tag, name, card_diff))
             else:
                 str_proc = "\n change process  ".join([""]+self.second_process)
+                str_proc = "\n<weight_card>\n%s\n</weight_card>\n"
                 for name in type_rwgt:
                     mg_rwgt_info.append((tag, name, str_proc + '\n'+ card_diff))
         else:
@@ -825,8 +854,8 @@
             if self.dedicated_path:
                 for k,v in self.dedicated_path.items():
                     str_info += "\n change %s %s" % (k,v)
-            card_diff = str_info
-            str_info += '\n' + s_new
+            card_diff = '<weight_slha>\n%s\n</weight_slha>' % str_info
+            str_info = '\n<weight_card>\n%s\n</weight_card>\n%s\n<weight_slha>\n%s</weight_slha>\n' % (str_info,str_version,s_new)
             for name in type_rwgt:
                 mg_rwgt_info.append((tag, name, str_info))
         # re-create the banner.
@@ -846,8 +875,11 @@
         self.banner['initrwgt'] = self.banner['initrwgt'].replace('\n\n', '\n')
 
         logger.info('starts to compute weight for events with the following modification to the param_card:')
-        logger.info(card_diff.replace('\n','\nKEEP:'))
-        self.run_card = banner.Banner(self.banner).charge_card('run_card')
+        #logger.info(card_diff.replace('\n','\nKEEP:'))
+        try:
+            self.run_card = banner.Banner(self.banner).charge_card('run_card')
+        except Exception:
+            self.run_card = None
 
         if self.options['rwgt_name']:
             tag_name = self.options['rwgt_name']
@@ -858,7 +890,7 @@
         for (path,tag), module in self.f2pylib.items():
             with misc.chdir(pjoin(os.path.dirname(rw_dir), path)):
                 with misc.stdchannel_redirected(sys.stdout, os.devnull):
-                    if 'second' in path or tag == 3:
+                    if 'rw_me_' in path or tag == 3:
                         module.initialise(pjoin(rw_dir, 'Cards', 'param_card.dat'))
                     else:
                         module.initialise(pjoin(path_me, 'rw_me', 'Cards', 'param_card_orig.dat'))
@@ -986,6 +1018,7 @@
             misc.sprint(nhel, Pdir, hel_dict)                        
             raise Exception, "Invalid matrix element for original computation (weight=0)"
 
+
         return {'orig': orig_wgt, '': w_new/w_orig*orig_wgt*jac}
      
     def calculate_nlo_weight(self, event):
@@ -1158,14 +1191,20 @@
         if (not self.second_model and not self.second_process and not self.dedicated_path) 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] 
+            if tag in self.id_to_path_second:
+                orig_order, Pdir, hel_dict = self.id_to_path_second[tag]
+            else:
+                return 0.0 
             
         base = os.path.basename(os.path.dirname(Pdir))
-        if '_second' in base:
+        if base == 'rw_me':
+            moduletag = (base, 2+hypp_id)
+        else:
             moduletag = (base, 2)
-        else:
-            moduletag = (base, 2+hypp_id)
+    
+            
         
+    
         module = self.f2pylib[moduletag]
 
         p = event.get_momenta(orig_order)
@@ -1175,7 +1214,7 @@
         if self.helicity_reweighting and 9 not in hel_order:
             nhel = hel_dict[tuple(hel_order)]                
         else:
-            nhel = -1
+            nhel = -2
             
         # For 2>N pass in the center of mass frame
         #   - required for helicity by helicity re-weighitng
@@ -1193,8 +1232,12 @@
 
         with misc.chdir(Pdir):
             with misc.stdchannel_redirected(sys.stdout, os.devnull):
-                me_value = module.smatrixhel(pdg,p, event.aqcd, scale2, nhel)
-                                
+                if hasattr(event, 'ievent'):
+                    ievent = event.ievent
+                else:
+                    ievent = -1
+                me_value = module.smatrixhel(pdg, ievent, p, event.aqcd, scale2, nhel)
+               
         # for loop we have also the stability status code
         if isinstance(me_value, tuple):
             me_value, code = me_value
@@ -1202,7 +1245,7 @@
             hundred_value = (code % 1000) //100
             if hundred_value in [4]:
                 me_value = 0.
-            
+
         return me_value
     
     def terminate_fortran_executables(self, new_card_only=False):
@@ -1308,10 +1351,13 @@
                     nlo_order = nlo_order.replace('noborn', 'virt')
                 commandline += "add process %s [%s] %s;" % (base,nlo_order,post)
             commandline = commandline.replace('add process', 'generate',1)
-            logger.info("RETRY with %s", commandline)
-            mgcmd.exec_cmd(commandline, precmd=True)
-            has_nlo = False
+            if commandline:
+                logger.info("RETRY with %s", commandline)
+                mgcmd.exec_cmd(commandline, precmd=True)
+                has_nlo = False
         except Exception, error:
+            misc.sprint(mgcmd._curr_model.get('modelpath'))
+            misc.sprint(type(error))
             raise
         
         commandline = 'output standalone_rw %s --prefix=int' % pjoin(path_me,data['paths'][0])
@@ -1504,7 +1550,13 @@
             #self.id_to_path = {}
             #data['id2path'] = self.id_to_path
         else:
-            data['paths'] = ['rw_me_second', 'rw_mevirt_second']
+            for key in self.f2pylib.keys():
+                if 'rw_me_%s' % self.nb_library in key[0]:
+                    del self.f2pylib[key]
+                
+            
+            self.nb_library += 1
+            data['paths'] = ['rw_me_%s' % self.nb_library, 'rw_mevirt_%s' % self.nb_library]
             # model
             if self.second_model:
                 data['mg_names'] = True
@@ -1536,14 +1588,19 @@
         else:
             path_me = self.rwgt_dir
         data['path'] = path_me
-        try:
-            shutil.rmtree(pjoin(path_me,data['paths'][0]))
-        except Exception: 
-            pass
-        try:
-            shutil.rmtree(pjoin(path_me, data['paths'][1]))
-        except Exception: 
-            pass
+        for i in range(2):
+            pdir = pjoin(path_me,data['paths'][i])
+            if os.path.exists(pdir):
+                try:
+                    shutil.rmtree(pjoin(path_me,data['paths'][0]))
+                except Exception, error:
+                    misc.sprint(error) 
+                    pass                
+                
+                
+                    
+        
+
 
         # 1. prepare the interface----------------------------------------------
         mgcmd = self.mg5cmd
@@ -1691,12 +1748,16 @@
             path_me = self.me_dir
         else:
             path_me = self.rwgt_dir
-        for onedir in self.rwgt_dir_possibility:
+            
+        rwgt_dir_possibility =   ['rw_me','rw_me_%s' % self.nb_library,'rw_mevirt','rw_mevirt_%s' % self.nb_library]
+        for onedir in rwgt_dir_possibility:
             if not os.path.isdir(pjoin(path_me,onedir)):
                 continue
             pdir = pjoin(path_me, onedir, 'SubProcesses')
             if self.mother:
                 nb_core = self.mother.options['nb_core'] if self.mother.options['run_mode'] !=0 else 1
+            elif hasattr(self, 'nb_core'):
+                nb_core = self.nb_core
             else:
                 nb_core = 1
             os.environ['MENUM'] = '2'
@@ -1715,7 +1776,8 @@
         
         self.id_to_path = {}
         self.id_to_path_second = {}
-        for onedir in self.rwgt_dir_possibility:
+        rwgt_dir_possibility =   ['rw_me','rw_me_%s' % self.nb_library,'rw_mevirt','rw_mevirt_%s' % self.nb_library]
+        for onedir in rwgt_dir_possibility:
             if not os.path.exists(pjoin(path_me,onedir)):
                 continue 
             pdir = pjoin(path_me, onedir, 'SubProcesses')
@@ -1724,14 +1786,18 @@
                     mod_name = '%s.SubProcesses.allmatrix%spy' % (onedir, tag)
                     #mymod = __import__('%s.SubProcesses.allmatrix%spy' % (onedir, tag), globals(), locals(), [],-1)
                     if mod_name in sys.modules.keys():
+#                        if 'rw_me_' in mod_name:
+#                            raise Exception
                         del sys.modules[mod_name]
                         tmp_mod_name = mod_name
                         while '.' in tmp_mod_name:
                             tmp_mod_name = tmp_mod_name.rsplit('.',1)[0]
                             del sys.modules[tmp_mod_name]
-                        mymod = __import__(mod_name, globals(), locals(), [],-1)  
+                        mymod = __import__(mod_name, globals(), locals(), [],-1) 
+                        reload(mymod) 
                     else:
                         mymod = __import__(mod_name, globals(), locals(), [],-1)    
+                            
                     
                     S = mymod.SubProcesses
                     mymod = getattr(S, 'allmatrix%spy' % tag)
@@ -1744,12 +1810,16 @@
                     break
 
             data = self.id_to_path
-            if '_second' in onedir:
+            if onedir not in ["rw_me",  "rw_mevirt"]:
                 data = self.id_to_path_second
 
+            
+
             # get all the information
             all_pdgs = mymod.get_pdg_order()
-            all_pdgs = [[pdg for pdg in pdgs if pdg!=0] for pdgs in  mymod.get_pdg_order()]
+            allids, all_pids = mymod.get_pdg_order()
+            all_pdgs = [[pdg for pdg in pdgs if pdg!=0] for pdgs in  allids]
+            #all_pids = [pid for (pdgs, pid) in  allids]
             all_prefix = [''.join(j).strip().lower() for j in mymod.get_prefix()]
             prefix_set = set(all_prefix)
 
@@ -1772,7 +1842,7 @@
                     misc.sprint(os.path.exists(pjoin(path_me,onedir,'SubProcesses','MadLoop5_resources', '%sHelConfigs.dat' % prefix.upper())))
                     continue
 
-            for i,pdg in enumerate(all_pdgs):
+            for i,(pdg,pid) in enumerate(zip(all_pdgs,all_pids)):
                 if self.is_decay:
                     incoming = [pdg[0]]
                     outgoing = pdg[1:]
@@ -1785,7 +1855,7 @@
                     outgoing.sort()
                 tag = (tuple(incoming), tuple(outgoing))
                 if 'virt' in onedir:
-                    tag = (tag, 'V')
+                    tag = (tag, 'V')  
                 prefix = all_prefix[i]
                 hel = hel_dict[prefix]
                 if tag in data:
@@ -1907,10 +1977,111 @@
                     
         
         
-        
-
-
-
-
-
-        
+if "__main__" == __name__:
+    # Check if optimize mode is (and should be) activated
+    import optparse
+    # Write out nice usage message if called with -h or --help
+    usage = "usage: %prog [options] [FILE] "
+    parser = optparse.OptionParser(usage=usage)
+    parser.add_option("-l", "--logging", default='INFO',
+                  help="logging level (DEBUG|INFO|WARNING|ERROR|CRITICAL) [%default]")
+    parser.add_option("-f", "--file", default='',
+                    help="Use script file FILE")
+    parser.add_option("-d", "--mgme_dir", default='', dest = 'mgme_dir',
+                  help="Use MG_ME directory MGME_DIR")
+    parser.add_option("","--debug", action="store_true", default=False, dest='debug', \
+                 help='force to launch debug mode')
+    (options, args) = parser.parse_args()
+    if len(args) == 0:
+        args = ''
+
+    if __debug__ and not options.debug and \
+    (not os.path.exists(os.path.join(root_path, 'bin','create_release.py'))):
+        subprocess.call([sys.executable] + ['-O'] + sys.argv)
+        sys.exit()
+
+
+    import logging
+    import logging.config
+    import madgraph.interface.coloring_logging
+
+    try:
+        import readline
+    except ImportError:
+        try:
+            import pyreadline as readline
+        except:
+            print "For tab completion and history, install module readline."
+    else:
+        import rlcompleter
+
+        if 'r261:67515' in sys.version and  'GCC 4.2.1 (Apple Inc. build 5646)' in sys.version:
+            readline.parse_and_bind("bind ^I rl_complete")
+            readline.__doc__ = 'libedit'
+
+        elif hasattr(readline, '__doc__'):
+            if 'libedit' not in readline.__doc__:
+                readline.parse_and_bind("tab: complete")
+            else:
+                readline.parse_and_bind("bind ^I rl_complete")
+        else:
+            readline.__doc__ = 'GNU'
+            readline.parse_and_bind("tab: complete")
+
+    # charge history file
+    try:
+        history_file = os.path.join(os.environ['HOME'], '.mg5', 'rwgthistory')
+        readline.read_history_file(history_file)
+    except:
+        pass
+
+    try:
+        import psyco
+        psyco.full()
+    except:
+        pass
+
+    try:
+        if __debug__ and options.logging == 'INFO':
+            options.logging = 'DEBUG'
+        logging.config.fileConfig(os.path.join(root_path, 'madgraph', 'interface', '.mg5_logging.conf'))
+        logging.root.setLevel(eval('logging.' + options.logging))
+        logging.getLogger('madgraph').setLevel(eval('logging.' + options.logging))
+        logging.getLogger('madevent').setLevel(eval('logging.' + options.logging))
+    except:
+        pass
+
+    # Call the cmd interface main loop
+    try:
+        if options.file or args:
+            # They are an input file
+            if args:
+                input_file = os.path.realpath(args[0])
+            else:
+                input_file = os.path.realpath(options.file)
+            print "using input+file", input_file
+            cmd_line = ReweightInterface()
+            cmd_line.use_rawinput = False
+            cmd_line.haspiping = False
+            cmd_line.import_command_file(input_file)
+            cmd_line.run_cmd('quit')
+        else:
+            # Interactive mode
+            try:
+                cmd_line = ReweightInterface()
+                cmd_line.use_rawinput = True
+                cmd_line.cmdloop()
+            except:
+                pass
+        try:
+            cmd_line.exec_cmd('quit all', printcmd=False)
+            readline.set_history_length(100)
+            if not os.path.exists(os.path.join(os.environ['HOME'], '.mg5')):
+                os.mkdir(os.path.join(os.environ['HOME'], '.mg5'))
+            readline.write_history_file(history_file)
+        except Exception, error:
+            pass
+    except KeyboardInterrupt:
+        print 'writting history and quit on KeyboardInterrupt'
+        pass
+

=== modified file 'madgraph/iolibs/export_v4.py'
--- madgraph/iolibs/export_v4.py	2018-12-22 20:33:20 +0000
+++ madgraph/iolibs/export_v4.py	2019-05-23 07:08:37 +0000
@@ -2066,17 +2066,18 @@
         
         template = """
 %(python_information)s
-  subroutine smatrixhel(pdgs, npdg, p, ALPHAS, SCALE2, nhel, ANS)
+  subroutine smatrixhel(pdgs, procid, npdg, p, ALPHAS, SCALE2, nhel, ANS)
   IMPLICIT NONE
 
 CF2PY double precision, intent(in), dimension(0:3,npdg) :: p
 CF2PY integer, intent(in), dimension(npdg) :: pdgs
+CF2PY integer, intent(in):: procid
 CF2PY integer, intent(in) :: npdg
 CF2PY double precision, intent(out) :: ANS
 CF2PY double precision, intent(in) :: ALPHAS
 CF2PY double precision, intent(in) :: SCALE2
   integer pdgs(*)
-  integer npdg, nhel
+  integer npdg, nhel, procid
   double precision p(*)
   double precision ANS, ALPHAS, PI,SCALE2
   include 'coupl.inc'
@@ -2084,7 +2085,7 @@
   PI = 3.141592653589793D0
   G = 2* DSQRT(ALPHAS*PI)
   CALL UPDATE_AS_PARAM()
-  if (scale2.ne.0d0) stop 1
+ 1  if (scale2.ne.0d0) stop 1
 
 %(smatrixhel)s
 
@@ -2100,12 +2101,31 @@
       RETURN
       END
 
-    subroutine get_pdg_order(PDG)
+    subroutine get_nincoming(pdgs, procid, npdg, N)
+    IMPLICIT NONE
+CF2PY integer, intent(in), dimension(npdg) :: pdgs
+CF2PY integer, intent(in):: procid
+CF2PY integer, intent(in) :: npdg
+CF2PY INTEGER, intent(out) :: N
+  integer pdgs(*)
+  integer npdg, procid
+    INTEGER N
+    %(nincoming)s
+    return
+    end
+
+    subroutine get_pdg_order(PDG, PID)
   IMPLICIT NONE
 CF2PY INTEGER, intent(out) :: PDG(%(nb_me)i,%(maxpart)i)  
-  INTEGER PDG(%(nb_me)i,%(maxpart)i), PDGS(%(nb_me)i,%(maxpart)i)
+CF2PY INTEGER, intent(out) :: PID(%(nb_me)i)
+  
+    INTEGER PDG(%(nb_me)i,%(maxpart)i), PDGS(%(nb_me)i,%(maxpart)i)
+    INTEGER PID(%(nb_me)i),PIDs(%(nb_me)i)
   DATA PDGS/ %(pdgs)s /
+  DATA PIDS/ %(pids)s /
+  
   PDG = PDGS
+  PID = PIDS
   RETURN
   END 
 
@@ -2123,43 +2143,58 @@
          
         allids = self.prefix_info.keys()
         allprefix = [self.prefix_info[key][0] for key in allids]
-        min_nexternal = min([len(ids) for ids in allids])
-        max_nexternal = max([len(ids) for ids in allids])
+        min_nexternal = min([len(ids[0]) for ids in allids])
+        max_nexternal = max([len(ids[0]) for ids in allids])
 
         info = []
-        for key, (prefix, tag) in self.prefix_info.items():
-            info.append('#PY %s : %s # %s' % (tag, key, prefix))
-            
+        for (key, pid), (prefix, tag) in self.prefix_info.items():
+            info.append('#PY %s : %s # %s %s' % (tag, key, prefix, pid))
 
         text = []
+        text_incomming = []
         for n_ext in range(min_nexternal, max_nexternal+1):
-            current = [ids for ids in allids if len(ids)==n_ext]
-            if not current:
+            current_id = [ids[0] for ids in allids if len(ids[0])==n_ext]
+            current_pid = [ids[1] for ids in allids if len(ids[0])==n_ext]
+            if not current_id:
                 continue
             if min_nexternal != max_nexternal:
                 if n_ext == min_nexternal:
                     text.append('       if (npdg.eq.%i)then' % n_ext)
+                    text_incomming.append('       if (npdg.eq.%i)then' % n_ext)
                 else:
                     text.append('       else if (npdg.eq.%i)then' % n_ext)
-            for ii,pdgs in enumerate(current):
+                    text_incomming.append('       else if (npdg.eq.%i)then' % n_ext)
+            for ii,pdgs in enumerate(current_id):
+                pid = current_pid[ii] 
+                #misc.sprint([(i,pdg) for i,pdg in pdgs])
                 condition = '.and.'.join(['%i.eq.pdgs(%i)' %(pdg, i+1) for i, pdg in enumerate(pdgs)])
                 if ii==0:
-                    text.append( ' if(%s) then ! %i' % (condition, i))
+                    text.append( ' if(%s.and.(procid.le.0.or.procid.eq.%d)) then ! %i' % (condition, pid, i))
+                    text_incomming.append( ' if(%s.and.(procid.le.0.or.procid.eq.%d)) then ! %i' % (condition, pid, i))
                 else:
-                    text.append( ' else if(%s) then ! %i' % (condition,i))
-                text.append(' call %ssmatrixhel(p, nhel, ans)' % self.prefix_info[pdgs][0])
+                    text.append( ' else if(%s.and.(procid.le.0.or.procid.eq.%d)) then ! %i' % (condition,pid,i))
+                    text_incomming.append( ' else if(%s.and.(procid.le.0.or.procid.eq.%d)) then ! %i' % (condition,pid,i))
+                text.append(' call %ssmatrixhel(p, nhel, ans)' % self.prefix_info[(pdgs,pid)][0])
+                text_incomming.append('     N=%i' % len(tag[0]))
+            text.append( ' else if(procid.gt.0) then !')
+            text.append( ' procid = -1' )
+            text.append( ' goto 1' )
             text.append(' endif')
+            text_incomming.append(' endif')
         #close the function
         if min_nexternal != max_nexternal:
             text.append('endif')
+            text_incomming.append(' endif')
 
         formatting = {'python_information':'\n'.join(info), 
                           'smatrixhel': '\n'.join(text),
                           'maxpart': max_nexternal,
                           'nb_me': len(allids),
                           'pdgs': ','.join(str(pdg[i]) if i<len(pdg) else '0' 
-                                             for i in range(max_nexternal) for pdg in allids),
-                          'prefix':'\',\''.join(allprefix)
+                                             for i in range(max_nexternal) for (pdg,pid) in allids),
+                          'prefix':'\',\''.join(allprefix),
+                          'pids': ','.join(str(pid) for (pdg,pid) in allids),
+                          'nincoming':'\n'.join(text_incomming)
                           }
         formatting['lenprefix'] = len(formatting['prefix'])
         text = template % formatting
@@ -2257,7 +2292,7 @@
                 raise Exception, '--prefix options supports only \'int\' and \'proc\''
             for proc in matrix_element.get('processes'):
                 ids = [l.get('id') for l in proc.get('legs_with_decays')]
-                self.prefix_info[tuple(ids)] = [proc_prefix, proc.get_tag()] 
+                self.prefix_info[(tuple(ids), proc.get('id'))] = [proc_prefix, proc.get_tag()] 
                 
         calls = self.write_matrix_element_v4(
             writers.FortranWriter(filename),
@@ -2509,6 +2544,7 @@
             else:
                 matrix_template = "matrix_standalone_splitOrders_v4.inc"
 
+
         replace_dict['template_file'] = pjoin(_file_path, 'iolibs', 'template_files', matrix_template)
         replace_dict['template_file2'] = pjoin(_file_path, \
                                    'iolibs/template_files/split_orders_helping_functions.inc')
@@ -5688,8 +5724,14 @@
         else:
             load_card = ''
             lha_read_filename='lha_read.f'
-        cp( MG5DIR + '/models/template_files/fortran/' + lha_read_filename, \
-                                       os.path.join(self.dir_path,'lha_read.f'))
+        
+        nb_para= len(self.params_dep + self.params_indep + self.params_ext)
+        template = pjoin(MG5DIR, 'models', 'template_files', 'fortran', lha_read_filename)
+        open(os.path.join(self.dir_path,'lha_read.f'),'w').write(
+            open(template).read() % {'maxpara':nb_para}) 
+        
+        #cp( MG5DIR + '/models/template_files/fortran/' + lha_read_filename, \
+        #                               os.path.join(self.dir_path,'lha_read.f'))
         
         file=file%{'includes':'\n      '.join(includes),
                    'load_card':load_card}

=== modified file 'madgraph/iolibs/template_files/loop/loop_matrix_standalone.inc'
--- madgraph/iolibs/template_files/loop/loop_matrix_standalone.inc	2017-04-25 14:47:44 +0000
+++ madgraph/iolibs/template_files/loop/loop_matrix_standalone.inc	2019-05-23 07:08:37 +0000
@@ -448,7 +448,7 @@
 C Make sure that lorentz rotation tests are not used if there is external loop wavefunction of spin 2 and that one specific helicity is asked
 NROTATIONS_DP_BU = NROTATIONS_DP
 NROTATIONS_QP_BU = NROTATIONS_QP
-if(MAX_SPIN_EXTERNAL_PARTICLE.gt.3.AND.USERHEL.NE.-1) THEN
+if(MAX_SPIN_EXTERNAL_PARTICLE.gt.3.AND.USERHEL.GT.-1) THEN
    if(.NOT.WARNED_LORENTZ_STAB_TEST_OFF) THEN
      WRITE(*,*) '##WARNING: Evaluation of a specific helicity was asked for this PS point, and there is a spin-2 (or higher) particle in the external states.'
 	 WRITE(*,*) '##WARNING: As a result, MadLoop disabled the Lorentz rotation test for this phase-space point only.'
@@ -481,7 +481,7 @@
   ENDDO
 101  CONTINUE
 CLOSE(1)
-IF (HelicityFilterLevel.eq.0) then
+IF (HelicityFilterLevel.eq.0.or.userhel.eq.-2) then
   FOUNDHELFILTER=.TRUE.
   DO J=1,NCOMB
 	GOODHEL(J)=.TRUE.
@@ -576,7 +576,7 @@
   HELPICKED=1
   CTMODE=CTMODEINIT
 ELSE
-  IF (USERHEL.ne.-1) THEN
+  IF (USERHEL.gt.0) THEN
     IF(.NOT.GOODHEL(USERHEL)) THEN
       ANS(1)=0.0d0
       ANS(2)=0.0d0
@@ -584,7 +584,7 @@
       goto 9999
     ENDIF
   ENDIF
-  HELPICKED=USERHEL      
+  HELPICKED=MAX(USERHEL,-1)      
   IF (CTMODERUN.GT.-1) THEN
     CTMODE=CTMODERUN
   ELSE
@@ -694,7 +694,7 @@
 IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
 C  Update of NEXTREF, will be used for loop induced only.
   NEXTREF = NEXTREF + ANS(1) + ANS(2) + ANS(3)
-  IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+  IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
   BUFFR(1)=BUFFR(1)+ANS(1)
   BUFFR(2)=BUFFR(2)+ANS(2)
   BUFFR(3)=BUFFR(3)+ANS(3)
@@ -785,7 +785,7 @@
 
 DO K=1,3
   ANS(K)=ANS(K)/DBLE(IDEN)
-  IF (USERHEL.NE.-1) THEN
+  IF (USERHEL.GT.-1) THEN
     ANS(K)=ANS(K)*HELAVGFACTOR
   ELSE
 	DO J=1,NINITIAL

=== modified file 'madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc'
--- madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc	2017-07-31 22:08:40 +0000
+++ madgraph/iolibs/template_files/loop_optimized/loop_matrix_standalone.inc	2019-05-23 07:08:37 +0000
@@ -642,6 +642,10 @@
    CALL %(proc_prefix)sSET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
 ENDIF
 
+IF(USERHEL.eq.-2) THEN
+		  NTRY = 0
+ENDIF
+
 IF(NTRY.EQ.0) THEN
   HELDOUBLECHECKED=(.NOT.DoubleCheckHelicityFilter).OR.(HelicityFilterLevel.eq.0)
 OPEN(1, FILE=LoopFilterFN, err=100, status='OLD',           action='READ')
@@ -659,7 +663,7 @@
 101  CONTINUE
 CLOSE(1)
 
-IF (HelicityFilterLevel.eq.0) then
+IF (HelicityFilterLevel.eq.0.or.userhel.eq.-2) then
   FOUNDHELFILTER=.TRUE.
   DO J=1,NCOMB
 	GOODHEL(J)=1
@@ -816,7 +820,7 @@
   HELPICKED=1
   CTMODE=CTMODEINIT
 ELSE
-  IF (USERHEL.ne.-1) then
+  IF (USERHEL.gt.0) then
     IF(GOODHEL(USERHEL).eq.-HELOFFSET) THEN
 	  DO I=0,NSQUAREDSO
         ANS(1,I)=0.0d0
@@ -826,7 +830,7 @@
       goto 9999
     ENDIF
   ENDIF
-  HELPICKED=USERHEL
+  HELPICKED=MAX(USERHEL,-1)
   IF (CTMODERUN.NE.-1) THEN
     CTMODE=CTMODERUN
   ELSE
@@ -891,7 +895,7 @@
 C MadLoop jumps to this label during stability checks when it recomputes the same PS point with a different CTMode
 300 CONTINUE
 C Of course the trick of reusing coefficients when reducing at the amplitude level only works when computing one helicity at a time
-if (USERHEL.ne.-1) THEN
+if (USERHEL.gt.-1) THEN
   SKIP_LOOPNUM_COEFS_CONSTRUCTION = .TRUE.
 ENDIF
 308 CONTINUE
@@ -1147,7 +1151,7 @@
 
 ## if(ComputeColorFlows){
 C We should compute the color flow either if it contributes to the final result (i.e. not used just for the filtering), or if the computation is only done from the color flows
-IF (((.NOT.DIRECT_ME_COMPUTATION).AND.ME_COMPUTATION_FROM_JAMP).OR.((H.EQ.USERHEL.OR.USERHEL.EQ.-1).and.(POLARIZATIONS(0,0).eq.-1.or.%(proc_prefix)sIS_HEL_SELECTED(H)))) THEN
+IF (((.NOT.DIRECT_ME_COMPUTATION).AND.ME_COMPUTATION_FROM_JAMP).OR.((H.EQ.USERHEL.OR.USERHEL.LE.-1).and.(POLARIZATIONS(0,0).eq.-1.or.%(proc_prefix)sIS_HEL_SELECTED(H)))) THEN
 C The cumulative quantities must only be computed if that helicity contributes according to user request (second argument of the subroutine below).
   CALL %(proc_prefix)sCOMPUTE_COLOR_FLOWS(HEL_MULT)
   IF(ME_COMPUTATION_FROM_JAMP) THEN
@@ -1172,7 +1176,7 @@
         ENDDO
       ENDIF
 C   The if statement below is not strictly necessary but makes it clear when it is executed.
-    ELSEIF(H.EQ.USERHEL.OR.USERHEL.EQ.-1) THEN
+    ELSEIF(H.EQ.USERHEL.OR.USERHEL.LE.-1) THEN
 C     Make sure that that no polarization constraint filters out this helicity
       IF (POLARIZATIONS(0,0).eq.-1.or.%(proc_prefix)sIS_HEL_SELECTED(H)) THEN
 C       If both computational method is used, then we must just update RES_FROM_JAMP
@@ -1183,7 +1187,7 @@
         ENDDO
       ENDIF
     ENDIF
-	IF (H.EQ.USERHEL.OR.USERHEL.EQ.-1) THEN
+	IF (H.EQ.USERHEL.OR.USERHEL.LE.-1) THEN
 C     Make sure that that no polarization constraint filters out this helicity
       IF (POLARIZATIONS(0,0).eq.-1.or.%(proc_prefix)sIS_HEL_SELECTED(H)) THEN
 	    CALL %(proc_prefix)sCOMPUTE_COLOR_FLOWS_DERIVED_QUANTITIES(HEL_MULT)
@@ -1315,7 +1319,7 @@
 1226 CONTINUE
 
 IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-  IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+  IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C   Make sure that that no polarization constraint filters out this helicity
     IF (POLARIZATIONS(0,0).eq.-1.or.%(proc_prefix)sIS_HEL_SELECTED(HELPICKED)) THEN
 C     TO KEEP TRACK OF THE FINAL ANSWER TO BE RETURNED DURING CHECK PHASE
@@ -1527,7 +1531,7 @@
 ## }
   DO I=0,NSQUAREDSO
     ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-    IF (USERHEL.NE.-1) THEN
+    IF (USERHEL.GT.-1) THEN
       ANS(K,I)=ANS(K,I)*HELAVGFACTOR
 	ELSE
 	  DO J=1,NINITIAL
@@ -1560,7 +1564,7 @@
     DO K=0,3
 ## }
       RES_FROM_JAMP(K,I)=RES_FROM_JAMP(K,I)/DBLE(IDEN)
-      IF (USERHEL.NE.-1) THEN
+      IF (USERHEL.GT.-1) THEN
         RES_FROM_JAMP(K,I)=RES_FROM_JAMP(K,I)*HELAVGFACTOR
 	  ELSE
 	    DO L=1,NINITIAL

=== modified file 'madgraph/iolibs/template_files/matrix_standalone_splitOrders_v4.inc'
--- madgraph/iolibs/template_files/matrix_standalone_splitOrders_v4.inc	2018-05-17 12:03:54 +0000
+++ madgraph/iolibs/template_files/matrix_standalone_splitOrders_v4.inc	2019-05-23 07:08:37 +0000
@@ -127,24 +127,27 @@
 C ----------
 C BEGIN CODE
 C ----------
-      NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-         JC(IHEL) = +1
-      ENDDO
-	  DO I=1,NSQAMPSO
-      	ANS(I) = 0D0	  
-	  ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points have back-to-back initial states
 C     for which some of the spin-2 helicity configurations are zero. But they are no longer zero
 C     if the point is boosted on the z-axis. Remember that HELAS helicity amplitudes are no longer
 C     lorentz invariant with expternal spin-2 particles (only the helicity sum is).
 C     For this reason, we simply remove the filterin when there is only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+	  IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
 		  GOODHEL(IHEL)=.TRUE.
 		ENDDO
+		IF(USERHEL.EQ.-2) USERHEL = -1
 	  ENDIF
+
+	  NTRY=NTRY+1
+      DO IHEL=1,NEXTERNAL
+         JC(IHEL) = +1
+      ENDDO
+	  DO I=1,NSQAMPSO
+      	ANS(I) = 0D0
+	  ENDDO
+
       DO IHEL=1,NCOMB
 	     IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN
            IF (GOODHEL(IHEL) .OR. NTRY .LT. 2 .OR.USERHEL.NE.-1) THEN

=== modified file 'madgraph/iolibs/template_files/matrix_standalone_v4.inc'
--- madgraph/iolibs/template_files/matrix_standalone_v4.inc	2018-05-17 12:03:54 +0000
+++ madgraph/iolibs/template_files/matrix_standalone_v4.inc	2019-05-23 07:08:37 +0000
@@ -78,7 +78,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C
 C GLOBAL VARIABLES
 C
@@ -101,21 +100,22 @@
 C ----------
 C BEGIN CODE
 C ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-         JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points have back-to-back initial states
 C     for which some of the spin-2 helicity configurations are zero. But they are no longer zero
 C     if the point is boosted on the z-axis. Remember that HELAS helicity amplitudes are no longer
 C     lorentz invariant with expternal spin-2 particles (only the helicity sum is).
 C     For this reason, we simply remove the filterin when there is only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
-        DO IHEL=1,NCOMB
-		  GOODHEL(IHEL)=.TRUE.
-		ENDDO
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
+	      DO IHEL=1,NCOMB
+		      GOODHEL(IHEL)=.TRUE.
+		  ENDDO
+		  IF(USERHEL.EQ.-2)  USERHEL=-1
 	  ENDIF
+	  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		  

=== modified file 'madgraph/loop/loop_exporters.py'
--- madgraph/loop/loop_exporters.py	2018-06-03 12:13:47 +0000
+++ madgraph/loop/loop_exporters.py	2019-05-23 07:08:37 +0000
@@ -304,11 +304,12 @@
         CALL SETMADLOOPPATH(PATH)
       END
 
-  subroutine smatrixhel(pdgs, npdg, p, ALPHAS, SCALES2, nhel, ANS, RETURNCODE)
+  subroutine smatrixhel(pdgs, procid, npdg, p, ALPHAS, SCALES2, nhel, ANS, RETURNCODE)
   IMPLICIT NONE
 
 CF2PY double precision, intent(in), dimension(0:3,npdg) :: p
 CF2PY integer, intent(in), dimension(npdg) :: pdgs
+CF2PY integer, intent(in):: procid
 CF2PY integer, intent(in) :: npdg
 CF2PY double precision, intent(out) :: ANS
 CF2PY integer, intent(out) :: RETURNCODE
@@ -316,22 +317,26 @@
 CF2PY double precision, intent(in) :: SCALES2
 
   integer pdgs(*)
-  integer npdg, nhel, RETURNCODE
+  integer npdg, nhel, RETURNCODE, procid
   double precision p(*)
   double precision ANS, ALPHAS, PI,SCALES2
-
+ 1 continue
 %(smatrixhel)s
 
       return
       end
   
-  subroutine get_pdg_order(OUT)
+  subroutine get_pdg_order(PDG,PID)
   IMPLICIT NONE
-CF2PY INTEGER, intent(out) :: OUT(%(nb_me)i,%(maxpart)i)  
+CF2PY INTEGER, intent(out) :: PDG(%(nb_me)i,%(maxpart)i)  
+CF2PY INTEGER, intent(out) :: PID(%(nb_me)i)
   
-  INTEGER OUT(%(nb_me)i,%(maxpart)i), PDGS(%(nb_me)i,%(maxpart)i)
+  INTEGER PDG(%(nb_me)i,%(maxpart)i), PDGS(%(nb_me)i,%(maxpart)i)
+  INTEGER PID(%(nb_me)i),PIDs(%(nb_me)i)
   DATA PDGS/ %(pdgs)s /
-  OUT=PDGS
+  DATA PIDS/ %(pids)s /  
+  PDG=PDGS
+  PID=PIDS
   RETURN
   END
   
@@ -348,31 +353,39 @@
          
         allids = self.prefix_info.keys()
         allprefix = [self.prefix_info[key][0] for key in allids]
-        min_nexternal = min([len(ids) for ids in allids])
-        max_nexternal = max([len(ids) for ids in allids])
+        min_nexternal = min([len(ids[0]) for ids in allids])
+        max_nexternal = max([len(ids[0]) for ids in allids])
 
         info = []
-        for key, (prefix, tag) in self.prefix_info.items():
-            info.append('#PY %s : %s # %s' % (tag, key, prefix))
+        for (key, pid), (prefix, tag) in self.prefix_info.items():
+            info.append('#PY %s : %s # %s %s' % (tag, key, prefix, pid))
             
 
         text = []
         for n_ext in range(min_nexternal, max_nexternal+1):
-            current = [ids for ids in allids if len(ids)==n_ext]
-            if not current:
+            current_id = [ids[0] for ids in allids if len(ids[0])==n_ext]
+            current_pid = [ids[1] for ids in allids if len(ids[0])==n_ext]
+            if not current_id:
                 continue
             if min_nexternal != max_nexternal:
                 if n_ext == min_nexternal:
                     text.append('       if (npdg.eq.%i)then' % n_ext)
                 else:
                     text.append('       else if (npdg.eq.%i)then' % n_ext)
-            for ii,pdgs in enumerate(current):
+            for ii,pdgs in enumerate(current_id):
+                pid = current_pid[ii]
                 condition = '.and.'.join(['%i.eq.pdgs(%i)' %(pdg, i+1) for i, pdg in enumerate(pdgs)])
                 if ii==0:
-                    text.append( ' if(%s) then ! %i' % (condition, i))
+                    text.append( ' if(%s.and.(procid.le.0.or.procid.eq.%d)) then ! %i' % (condition, pid, i))
+                    #text.append( ' if(%s) then ! %i' % (condition, i))
                 else:
-                    text.append( ' else if(%s) then ! %i' % (condition,i))
-                text.append(' call %sget_me(p, ALPHAS, DSQRT(SCALES2), NHEL, ANS, RETURNCODE)' % self.prefix_info[pdgs][0])
+                    text.append( ' else if(%s.and.(procid.le.0.or.procid.eq.%d)) then ! %i' % (condition,pid,i))
+                    #text.append( ' else if(%s) then ! %i' % (condition,i))
+                text.append(' call %sget_me(p, ALPHAS, DSQRT(SCALES2), NHEL, ANS, RETURNCODE)' % self.prefix_info[(pdgs,pid)][0])
+            text.append( ' else if(procid.gt.0) then !')
+            text.append( ' procid = -1' )
+            text.append( ' goto 1' )
+            #text.append( ' call smatrixhel(pdgs, -1, npdg, p, ALPHAS, SCALES2, nhel, ANS, RETURNCODE)')
             text.append(' endif')
         #close the function
         if min_nexternal != max_nexternal:
@@ -384,7 +397,8 @@
                           'nb_me': len(allids),
                           'pdgs': ','.join([str(pdg[i]) if i<len(pdg) else '0' 
                                              for i in range(max_nexternal) \
-                                             for pdg in allids]),
+                                             for (pdg,pid) in allids]),
+                          'pids': ','.join(str(pid) for (pdg,pid) in allids),
                       'prefix':'\',\''.join(allprefix)
                       }
     
@@ -393,6 +407,7 @@
         fsock = writers.FortranWriter(pjoin(self.dir_path, 'SubProcesses', 'all_matrix.f'),'w')
         fsock.writelines(text)
         fsock.close()
+
         
     
     
@@ -984,7 +999,7 @@
         if 'prefix' in self.cmd_options and self.cmd_options['prefix'] in ['int','proc']:
             for proc in matrix_element.get('processes'):
                 ids = [l.get('id') for l in proc.get('legs_with_decays')]
-                self.prefix_info[tuple(ids)] = [dict['proc_prefix'], proc.get_tag()]
+                self.prefix_info[(tuple(ids), proc.get('id'))] = [dict['proc_prefix'], proc.get_tag()]
 
         # The proc_id is used for MadEvent grouping, so none of our concern here
         # and it is simply set to an empty string.        

=== modified file 'madgraph/various/banner.py'
--- madgraph/various/banner.py	2019-05-07 13:01:35 +0000
+++ madgraph/various/banner.py	2019-05-23 07:08:37 +0000
@@ -384,9 +384,12 @@
             header = open(pjoin(MG5DIR,'Template', 'LO', 'Source', 'banner_header.txt')).read()
             
         if not self.lhe_version:
-            self.lhe_version = self.get('run_card', 'lhe_version', default=1.0)
-            if float(self.lhe_version) < 3:
-                self.lhe_version = 1.0
+            if 'mgruncard' in self:
+                self.lhe_version = self.get('run_card', 'lhe_version', default=1.0)
+                if float(self.lhe_version) < 3:
+                    self.lhe_version = 1.0
+            else:
+                self.lhe_version = 3.0
         
         ff.write(header % { 'version':float(self.lhe_version)})
 

=== modified file 'madgraph/various/rambo.py'
--- madgraph/various/rambo.py	2015-10-01 16:00:08 +0000
+++ madgraph/various/rambo.py	2019-05-23 07:08:37 +0000
@@ -1,6 +1,7 @@
 from __future__ import division
 import math
 import random
+import madgraph.various.misc as misc
 
 class FortranList(list):
     
@@ -91,6 +92,12 @@
     V= FortranList(N)
     IWARN = [0,0]
 # Check input object
+    if isinstance(XM, list) and not isinstance(XM, FortranList):
+        tmp = FortranList(1, N)
+        for i,val in enumerate(XM):
+            tmp[i+1] = val
+        XM = tmp
+
     assert isinstance(XM, FortranList)
     assert XM.min == 1
     assert XM.max == N+1 

=== modified file 'madgraph/various/shower_card.py'
--- madgraph/various/shower_card.py	2019-05-03 07:40:21 +0000
+++ madgraph/various/shower_card.py	2019-05-23 07:08:37 +0000
@@ -117,7 +117,6 @@
                 if not i in list_dm: 
                     self['dm_'+str(i)] = ''
 
-        misc.sprint(self['nsplit_jobs'])
         self.text=content
 
 

=== modified file 'models/import_ufo.py'
--- models/import_ufo.py	2019-02-14 13:39:51 +0000
+++ models/import_ufo.py	2019-05-23 07:08:37 +0000
@@ -585,7 +585,7 @@
             else:
                 pass
         else:
-            self.model.set('order_hierarchy', hierarchy)            
+            self.model.set('order_hierarchy', hierarchy)     
         
         # Also set expansion_order, i.e., maximum coupling order per process
         expansion_order={}
@@ -602,11 +602,19 @@
                 pass
         else:
             self.model.set('expansion_order', expansion_order)
-            self.model.set('expansion_order', expansion_order)            
+
             
         #clean memory
         del self.checked_lor
 
+        for name in dir(self.ufomodel):
+            if name.startswith('__') and name.endswith('__'):
+                value = getattr(self.ufomodel, name)
+                if not isinstance(value, (str, int)):
+                    continue
+                self.model['model_info'][name[2:-2]] = value
+
+
         return self.model
     
     def optimise_interaction(self, interaction):

=== modified file 'models/template_files/fortran/lha_read.f'
--- models/template_files/fortran/lha_read.f	2016-10-27 15:40:27 +0000
+++ models/template_files/fortran/lha_read.f	2019-05-23 07:08:37 +0000
@@ -179,7 +179,7 @@
       implicit none
 
       integer maxpara
-      parameter (maxpara=1000)
+      parameter (maxpara=%(maxpara)i)
       character*20 param(maxpara),value(maxpara),val,par
       character*20 blockname
       integer npara
@@ -280,7 +280,7 @@
 c     parameters
 c
       integer maxpara
-      parameter (maxpara=1000)
+      parameter (maxpara=%(maxpara)i)
 c
 c     arguments
 c

=== modified file 'models/template_files/fortran/lha_read_mp.f'
--- models/template_files/fortran/lha_read_mp.f	2016-10-27 15:40:27 +0000
+++ models/template_files/fortran/lha_read_mp.f	2019-05-23 07:08:37 +0000
@@ -199,7 +199,7 @@
       implicit none
 
       integer maxpara
-      parameter (maxpara=1000)
+      parameter (maxpara=%(maxpara)i)
       character*20 param(maxpara),value(maxpara),val,par
       character*20 blockname
       integer npara
@@ -300,7 +300,7 @@
 c     parameters
 c
       integer maxpara
-      parameter (maxpara=1000)
+      parameter (maxpara=%(maxpara)i)
 c
 c     arguments
 c
@@ -356,7 +356,7 @@
 c     parameters
 c
       integer maxpara
-      parameter (maxpara=1000)
+      parameter (maxpara=%(maxpara)i)
 c
 c     arguments
 c

=== modified file 'tests/acceptance_tests/test_cmd_reweight.py'
--- tests/acceptance_tests/test_cmd_reweight.py	2018-06-01 12:14:52 +0000
+++ tests/acceptance_tests/test_cmd_reweight.py	2019-05-23 07:08:37 +0000
@@ -207,7 +207,7 @@
         ff.write(cmd_lines)
         ff.close()
         
-        with misc.stdchannel_redirected(sys.stdout, os.devnull):
+        if 1:#with misc.stdchannel_redirected(sys.stdout, os.devnull):
             me_cmd.run_cmd('reweight run_01 --from_cards')
         
         lhe = lhe_parser.EventFile(pjoin(self.run_dir,'Events','run_01', 'events.lhe.gz'))
@@ -281,7 +281,7 @@
         ff.write(cmd_lines)
         ff.close()
         
-        with misc.stdchannel_redirected(sys.stdout, os.devnull):
+        if 1 :#with misc.stdchannel_redirected(sys.stdout, os.devnull):
             me_cmd.run_cmd('reweight run_01 --from_cards')
         
         lhe = lhe_parser.EventFile(pjoin(self.run_dir,'Events','run_01', 'events.lhe.gz'))
@@ -312,7 +312,8 @@
         ff.write(cmd_lines)
         ff.close()
         
-        with misc.stdchannel_redirected(sys.stdout, os.devnull):
+        #with misc.stdchannel_redirected(sys.stdout, os.devnull):
+        if 1:
             me_cmd.run_cmd('reweight run_01 --from_cards')
             
 

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%born_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%born_matrix.f	2019-05-23 07:08:37 +0000
@@ -141,13 +141,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
-      DO I=1,NSQAMPSO
-        ANS(I) = 0D0
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -160,11 +153,21 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2) USERHEL = -1
       ENDIF
+
+      NTRY=NTRY+1
+      DO IHEL=1,NEXTERNAL
+        JC(IHEL) = +1
+      ENDDO
+      DO I=1,NSQAMPSO
+        ANS(I) = 0D0
+      ENDDO
+
       DO IHEL=1,NCOMB
         IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN
           IF (GOODHEL(IHEL) .OR. NTRY .LT. 2 .OR.USERHEL.NE.-1) THEN

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f	2018-11-08 22:24:49 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_dxu_wp%V0_dxu_wp%loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1069,7 +1073,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.IS_HEL_SELECTED(HELPICKED))
@@ -1287,7 +1291,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%born_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%born_matrix.f	2019-05-23 07:08:37 +0000
@@ -141,13 +141,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
-      DO I=1,NSQAMPSO
-        ANS(I) = 0D0
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -160,11 +153,21 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2) USERHEL = -1
       ENDIF
+
+      NTRY=NTRY+1
+      DO IHEL=1,NEXTERNAL
+        JC(IHEL) = +1
+      ENDDO
+      DO I=1,NSQAMPSO
+        ANS(I) = 0D0
+      ENDDO
+
       DO IHEL=1,NCOMB
         IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN
           IF (GOODHEL(IHEL) .OR. NTRY .LT. 2 .OR.USERHEL.NE.-1) THEN

=== modified file 'tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_matrix.f'
--- tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_matrix.f	2018-11-08 22:24:49 +0000
+++ tests/input_files/IOTestsComparison/IOExportFKSTest/test_ppw_fksall/%SubProcesses%P0_udx_wp%V0_udx_wp%loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1069,7 +1073,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.IS_HEL_SELECTED(HELPICKED))
@@ -1287,7 +1291,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_standalone/matrix.f'
--- tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_standalone/matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/IOExportV4IOTest/export_matrix_element_v4_standalone/matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -134,10 +133,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -150,11 +145,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/IOTest_Histogram/gnuplot_histo_output/HistoOut.gnuplot'
--- tests/input_files/IOTestsComparison/IOTest_Histogram/gnuplot_histo_output/HistoOut.gnuplot	2016-06-09 00:59:23 +0000
+++ tests/input_files/IOTestsComparison/IOTest_Histogram/gnuplot_histo_output/HistoOut.gnuplot	2019-05-23 07:08:37 +0000
@@ -26,62 +26,61 @@
 # color-blind individuals with either protanopia or deuteranopia. Bang
 # Wong [2011] Nature Methods 8, 441.
 
-set style line  1 lt 1 lc rgb "#009e73" lw 1.3
+set style line   1 lt 1 lc rgb "#009e73" lw 1.3
 set style line 101 lt 1 lc rgb "#009e73" lw 1.3 dt (6,3)
-set style line 11 lt 2 lc rgb "#009e73" lw 1.3 dt (6,3)
-set style line 21 lt 4 lc rgb "#009e73" lw 1.3 dt (3,2)
-set style line 31 lt 6 lc rgb "#009e73" lw 1.3 dt (2,1)
-set style line 41 lt 8 lc rgb "#009e73" lw 1.3 dt (4,3)
-
-set style line  2 lt 1 lc rgb "#0072b2" lw 1.3
-set style line  102 lt 1 lc rgb "#0072b2" lw 1.3 dt (6,3)
-set style line 12 lt 2 lc rgb "#0072b2" lw 1.3 dt (6,3)
-set style line 22 lt 4 lc rgb "#0072b2" lw 1.3 dt (3,2)
-set style line 32 lt 6 lc rgb "#0072b2" lw 1.3 dt (2,1)
-set style line 42 lt 8 lc rgb "#0072b2" lw 1.3 dt (4,3)
-
-
-set style line  3 lt 1 lc rgb "#d55e00" lw 1.3
+set style line  11 lt 2 lc rgb "#009e73" lw 1.3 dt (6,3)
+set style line  21 lt 4 lc rgb "#009e73" lw 1.3 dt (3,2)
+set style line  31 lt 6 lc rgb "#009e73" lw 1.3 dt (2,1)
+set style line  41 lt 8 lc rgb "#009e73" lw 1.3 dt (4,3)
+
+set style line   2 lt 1 lc rgb "#0072b2" lw 1.3
+set style line 102 lt 1 lc rgb "#0072b2" lw 1.3 dt (6,3)
+set style line  12 lt 2 lc rgb "#0072b2" lw 1.3 dt (6,3)
+set style line  22 lt 4 lc rgb "#0072b2" lw 1.3 dt (3,2)
+set style line  32 lt 6 lc rgb "#0072b2" lw 1.3 dt (2,1)
+set style line  42 lt 8 lc rgb "#0072b2" lw 1.3 dt (4,3)
+
+set style line   3 lt 1 lc rgb "#d55e00" lw 1.3
 set style line 103 lt 1 lc rgb "#d55e00" lw 1.3 dt (6,3)
-set style line 13 lt 2 lc rgb "#d55e00" lw 1.3 dt (6,3)
-set style line 23 lt 4 lc rgb "#d55e00" lw 1.3 dt (3,2)
-set style line 33 lt 6 lc rgb "#d55e00" lw 1.3 dt (2,1)
-set style line 43 lt 8 lc rgb "#d55e00" lw 1.3 dt (4,3)
-
-set style line  4 lt 1 lc rgb "#f0e442" lw 1.3
-set style line  104 lt 1 lc rgb "#f0e442" lw 1.3 dt (6,3)
-set style line 14 lt 2 lc rgb "#f0e442" lw 1.3 dt (6,3)
-set style line 24 lt 4 lc rgb "#f0e442" lw 1.3 dt (3,2)
-set style line 34 lt 6 lc rgb "#f0e442" lw 1.3 dt (2,1)
-set style line 44 lt 8 lc rgb "#f0e442" lw 1.3 dt (4,3)
-
-set style line  5 lt 1 lc rgb "#56b4e9" lw 1.3
-set style line  105 lt 1 lc rgb "#56b4e9" lw 1.3 dt (6,3)
-set style line 15 lt 2 lc rgb "#56b4e9" lw 1.3 dt (6,3)
-set style line 25 lt 4 lc rgb "#56b4e9" lw 1.3 dt (3,2)
-set style line 35 lt 6 lc rgb "#56b4e9" lw 1.3 dt (2,1)
-set style line 45 lt 8 lc rgb "#56b4e9" lw 1.3 dt (4,3)
-
-set style line  6 lt 1 lc rgb "#cc79a7" lw 1.3
-set style line  106 lt 1 lc rgb "#cc79a7" lw 1.3 dt (6,3)
-set style line 16 lt 2 lc rgb "#cc79a7" lw 1.3 dt (6,3)
-set style line 26 lt 4 lc rgb "#cc79a7" lw 1.3 dt (3,2)
-set style line 36 lt 6 lc rgb "#cc79a7" lw 1.3 dt (2,1)
-set style line 46 lt 8 lc rgb "#cc79a7" lw 1.3 dt (4,3)
-
-set style line  7 lt 1 lc rgb "#e69f00" lw 1.3
-set style line  107 lt 1 lc rgb "#e69f00" lw 1.3 dt (6,3)
-set style line 17 lt 2 lc rgb "#e69f00" lw 1.3 dt (6,3)
-set style line 27 lt 4 lc rgb "#e69f00" lw 1.3 dt (3,2)
-set style line 37 lt 6 lc rgb "#e69f00" lw 1.3 dt (2,1)
-set style line 47 lt 8 lc rgb "#e69f00" lw 1.3 dt (4,3)
-
-set style line  8 lt 1 lc rgb "black" lw 1.3
-set style line  108 lt 1 lc rgb "black" lw 1.3 dt (6,3)
-set style line 18 lt 2 lc rgb "black" lw 1.3 dt (6,3)
-set style line 28 lt 4 lc rgb "black" lw 1.3 dt (3,2)
-set style line 38 lt 6 lc rgb "black" lw 1.3 dt (2,1)
-set style line 48 lt 8 lc rgb "black" lw 1.3 dt (4,3)
+set style line  13 lt 2 lc rgb "#d55e00" lw 1.3 dt (6,3)
+set style line  23 lt 4 lc rgb "#d55e00" lw 1.3 dt (3,2)
+set style line  33 lt 6 lc rgb "#d55e00" lw 1.3 dt (2,1)
+set style line  43 lt 8 lc rgb "#d55e00" lw 1.3 dt (4,3)
+
+set style line   4 lt 1 lc rgb "#f0e442" lw 1.3
+set style line 104 lt 1 lc rgb "#f0e442" lw 1.3 dt (6,3)
+set style line  14 lt 2 lc rgb "#f0e442" lw 1.3 dt (6,3)
+set style line  24 lt 4 lc rgb "#f0e442" lw 1.3 dt (3,2)
+set style line  34 lt 6 lc rgb "#f0e442" lw 1.3 dt (2,1)
+set style line  44 lt 8 lc rgb "#f0e442" lw 1.3 dt (4,3)
+
+set style line   5 lt 1 lc rgb "#56b4e9" lw 1.3
+set style line 105 lt 1 lc rgb "#56b4e9" lw 1.3 dt (6,3)
+set style line  15 lt 2 lc rgb "#56b4e9" lw 1.3 dt (6,3)
+set style line  25 lt 4 lc rgb "#56b4e9" lw 1.3 dt (3,2)
+set style line  35 lt 6 lc rgb "#56b4e9" lw 1.3 dt (2,1)
+set style line  45 lt 8 lc rgb "#56b4e9" lw 1.3 dt (4,3)
+
+set style line   6 lt 1 lc rgb "#cc79a7" lw 1.3
+set style line 106 lt 1 lc rgb "#cc79a7" lw 1.3 dt (6,3)
+set style line  16 lt 2 lc rgb "#cc79a7" lw 1.3 dt (6,3)
+set style line  26 lt 4 lc rgb "#cc79a7" lw 1.3 dt (3,2)
+set style line  36 lt 6 lc rgb "#cc79a7" lw 1.3 dt (2,1)
+set style line  46 lt 8 lc rgb "#cc79a7" lw 1.3 dt (4,3)
+
+set style line   7 lt 1 lc rgb "#e69f00" lw 1.3
+set style line 107 lt 1 lc rgb "#e69f00" lw 1.3 dt (6,3)
+set style line  17 lt 2 lc rgb "#e69f00" lw 1.3 dt (6,3)
+set style line  27 lt 4 lc rgb "#e69f00" lw 1.3 dt (3,2)
+set style line  37 lt 6 lc rgb "#e69f00" lw 1.3 dt (2,1)
+set style line  47 lt 8 lc rgb "#e69f00" lw 1.3 dt (4,3)
+
+set style line   8 lt 1 lc rgb "black" lw 1.3
+set style line 108 lt 1 lc rgb "black" lw 1.3 dt (6,3)
+set style line  18 lt 2 lc rgb "black" lw 1.3 dt (6,3)
+set style line  28 lt 4 lc rgb "black" lw 1.3 dt (3,2)
+set style line  38 lt 6 lc rgb "black" lw 1.3 dt (2,1)
+set style line  48 lt 8 lc rgb "black" lw 1.3 dt (4,3)
 
 
 set style line 999 lt 1 lc rgb "gray" lw 1.3

=== 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	2018-05-25 09:13:24 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_eq_4.f	2019-05-23 07:08:37 +0000
@@ -644,6 +644,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -663,7 +667,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -804,7 +808,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -814,7 +818,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1099,7 +1103,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1317,7 +1321,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== 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	2018-05-25 09:13:24 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_0_QEDAmpAndQEDsq_gt_2.f	2019-05-23 07:08:37 +0000
@@ -644,6 +644,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -663,7 +667,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -804,7 +808,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -814,7 +818,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1100,7 +1104,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1318,7 +1322,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== 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	2018-05-25 09:13:24 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QCDsq_gt_4.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1069,7 +1073,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1287,7 +1291,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== 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	2018-05-25 09:13:24 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_QEDsq_le_4.f	2019-05-23 07:08:37 +0000
@@ -644,6 +644,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -663,7 +667,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -804,7 +808,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -814,7 +818,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1100,7 +1104,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1318,7 +1322,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== 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	2018-05-25 09:13:24 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_WGTsq_le_10_QEDAmpAndQEDsq_gt_2.f	2019-05-23 07:08:37 +0000
@@ -644,6 +644,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -663,7 +667,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -804,7 +808,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -814,7 +818,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1100,7 +1104,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1318,7 +1322,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== 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	2018-05-25 09:13:24 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDQEDpert_default.f	2019-05-23 07:08:37 +0000
@@ -644,6 +644,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -663,7 +667,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -804,7 +808,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -814,7 +818,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1099,7 +1103,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1317,7 +1321,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== 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	2018-05-25 09:13:24 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QCDpert_default.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1069,7 +1073,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1287,7 +1291,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QEDpert_default.f'
--- tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QEDpert_default.f	2017-08-15 19:34:40 +0000
+++ tests/input_files/IOTestsComparison/LoopSquaredOrder_IOTest/Loop_sqso_uux_ddx/loop_matrix_QEDpert_default.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1069,7 +1073,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1287,7 +1291,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%born_matrix.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%born_matrix.f	2019-05-23 07:08:37 +0000
@@ -144,13 +144,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
-      DO I=1,NSQAMPSO
-        ANS(I) = 0D0
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -163,11 +156,21 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2) USERHEL = -1
       ENDIF
+
+      NTRY=NTRY+1
+      DO IHEL=1,NEXTERNAL
+        JC(IHEL) = +1
+      ENDDO
+      DO I=1,NSQAMPSO
+        ANS(I) = 0D0
+      ENDDO
+
       DO IHEL=1,NCOMB
         IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN
           IF (GOODHEL(IHEL) .OR. NTRY .LT. 2 .OR.USERHEL.NE.-1) THEN

=== modified file 'tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_matrix.f'
--- tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_matrix.f	2018-11-08 22:24:49 +0000
+++ tests/input_files/IOTestsComparison/MadLoop_output_from_the_interface/TIR_output/%ggttx_IOTest%SubProcesses%P0_gg_ttx%loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1069,7 +1073,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1287,7 +1291,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_NoSQSO.f'
--- tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_NoSQSO.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_NoSQSO.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -166,10 +165,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -182,11 +177,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_QCDsq_le_6.f'
--- tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_QCDsq_le_6.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_QCDsq_le_6.f	2019-05-23 07:08:37 +0000
@@ -192,13 +192,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
-      DO I=1,NSQAMPSO
-        ANS(I) = 0D0
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -211,11 +204,21 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2) USERHEL = -1
       ENDIF
+
+      NTRY=NTRY+1
+      DO IHEL=1,NEXTERNAL
+        JC(IHEL) = +1
+      ENDDO
+      DO I=1,NSQAMPSO
+        ANS(I) = 0D0
+      ENDDO
+
       DO IHEL=1,NCOMB
         IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN
           IF (GOODHEL(IHEL) .OR. NTRY .LT. 2 .OR.USERHEL.NE.-1) THEN

=== modified file 'tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_ampOrderQED2_eq_2_WGTsq_le_14_QCDsq_gt_4.f'
--- tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_ampOrderQED2_eq_2_WGTsq_le_14_QCDsq_gt_4.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/SquaredOrder_IOTest/sqso_uux_uuxuuxx/matrix_ampOrderQED2_eq_2_WGTsq_le_14_QCDsq_gt_4.f	2019-05-23 07:08:37 +0000
@@ -192,13 +192,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
-      DO I=1,NSQAMPSO
-        ANS(I) = 0D0
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -211,11 +204,21 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2) USERHEL = -1
       ENDIF
+
+      NTRY=NTRY+1
+      DO IHEL=1,NEXTERNAL
+        JC(IHEL) = +1
+      ENDDO
+      DO I=1,NSQAMPSO
+        ANS(I) = 0D0
+      ENDDO
+
       DO IHEL=1,NCOMB
         IF (USERHEL.EQ.-1.OR.USERHEL.EQ.IHEL) THEN
           IF (GOODHEL(IHEL) .OR. NTRY .LT. 2 .OR.USERHEL.NE.-1) THEN

=== modified file 'tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_matrix.f'
--- tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_matrix.f	2018-11-08 22:24:49 +0000
+++ tests/input_files/IOTestsComparison/TestCmdMatchBox/MatchBoxOutput/%TEST%SubProcesses%P1_uux_uux%loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL MG5_1_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1069,7 +1073,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.MG5_1_IS_HEL_SELECTED(HELPICK
@@ -1287,7 +1291,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/born_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/born_matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -134,10 +133,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -150,11 +145,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/loop_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/loop_matrix.f	2017-08-14 07:04:15 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_default/dux_mumvmxg/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -492,7 +492,7 @@
 C      helicity is asked
       NROTATIONS_DP_BU = NROTATIONS_DP
       NROTATIONS_QP_BU = NROTATIONS_QP
-      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.NE.-1) THEN
+      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.GT.-1) THEN
         IF(.NOT.WARNED_LORENTZ_STAB_TEST_OFF) THEN
           WRITE(*,*) '##WARNING: Evaluation of a specific helicity was'
      $     //' asked for this PS point, and there is a spin-2 (or'
@@ -531,7 +531,7 @@
         ENDDO
  101    CONTINUE
         CLOSE(1)
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=.TRUE.
@@ -636,7 +636,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(.NOT.GOODHEL(USERHEL)) THEN
             ANS(1)=0.0D0
             ANS(2)=0.0D0
@@ -644,7 +644,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.GT.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -898,7 +898,7 @@
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
 C       Update of NEXTREF, will be used for loop induced only.
         NEXTREF = NEXTREF + ANS(1) + ANS(2) + ANS(3)
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
           BUFFR(1)=BUFFR(1)+ANS(1)
           BUFFR(2)=BUFFR(2)+ANS(2)
           BUFFR(3)=BUFFR(3)+ANS(3)
@@ -1010,7 +1010,7 @@
 
       DO K=1,3
         ANS(K)=ANS(K)/DBLE(IDEN)
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.-1) THEN
           ANS(K)=ANS(K)*HELAVGFACTOR
         ELSE
           DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/born_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/born_matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -150,10 +149,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -166,11 +161,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/loop_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/loop_matrix.f	2017-08-15 19:34:40 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_default/gg_wmtbx/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -492,7 +492,7 @@
 C      helicity is asked
       NROTATIONS_DP_BU = NROTATIONS_DP
       NROTATIONS_QP_BU = NROTATIONS_QP
-      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.NE.-1) THEN
+      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.GT.-1) THEN
         IF(.NOT.WARNED_LORENTZ_STAB_TEST_OFF) THEN
           WRITE(*,*) '##WARNING: Evaluation of a specific helicity was'
      $     //' asked for this PS point, and there is a spin-2 (or'
@@ -531,7 +531,7 @@
         ENDDO
  101    CONTINUE
         CLOSE(1)
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=.TRUE.
@@ -636,7 +636,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(.NOT.GOODHEL(USERHEL)) THEN
             ANS(1)=0.0D0
             ANS(2)=0.0D0
@@ -644,7 +644,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.GT.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -2014,7 +2014,7 @@
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
 C       Update of NEXTREF, will be used for loop induced only.
         NEXTREF = NEXTREF + ANS(1) + ANS(2) + ANS(3)
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
           BUFFR(1)=BUFFR(1)+ANS(1)
           BUFFR(2)=BUFFR(2)+ANS(2)
           BUFFR(3)=BUFFR(3)+ANS(3)
@@ -2126,7 +2126,7 @@
 
       DO K=1,3
         ANS(K)=ANS(K)/DBLE(IDEN)
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.-1) THEN
           ANS(K)=ANS(K)*HELAVGFACTOR
         ELSE
           DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/born_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/born_matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -134,10 +133,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -150,11 +145,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_matrix.f	2017-08-15 19:34:40 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/dux_mumvmxg/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1063,7 +1067,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1281,7 +1285,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/gg_wmtbx/born_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/gg_wmtbx/born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/gg_wmtbx/born_matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -150,10 +149,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -166,11 +161,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/gg_wmtbx/loop_matrix.f'
--- tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/gg_wmtbx/loop_matrix.f	2017-08-14 07:04:15 +0000
+++ tests/input_files/IOTestsComparison/long_ML_SMQCD_optimized/gg_wmtbx/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1063,7 +1067,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1281,7 +1285,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_LoopInduced/gg_hh/loop_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_LoopInduced/gg_hh/loop_matrix.f	2017-08-15 19:34:40 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_LoopInduced/gg_hh/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -490,7 +490,7 @@
 C      helicity is asked
       NROTATIONS_DP_BU = NROTATIONS_DP
       NROTATIONS_QP_BU = NROTATIONS_QP
-      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.NE.-1) THEN
+      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.GT.-1) THEN
         IF(.NOT.WARNED_LORENTZ_STAB_TEST_OFF) THEN
           WRITE(*,*) '##WARNING: Evaluation of a specific helicity was'
      $     //' asked for this PS point, and there is a spin-2 (or'
@@ -529,7 +529,7 @@
         ENDDO
  101    CONTINUE
         CLOSE(1)
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=.TRUE.
@@ -634,7 +634,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(.NOT.GOODHEL(USERHEL)) THEN
             ANS(1)=0.0D0
             ANS(2)=0.0D0
@@ -642,7 +642,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.GT.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -922,7 +922,7 @@
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
 C       Update of NEXTREF, will be used for loop induced only.
         NEXTREF = NEXTREF + ANS(1) + ANS(2) + ANS(3)
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
           BUFFR(1)=BUFFR(1)+ANS(1)
           BUFFR(2)=BUFFR(2)+ANS(2)
           BUFFR(3)=BUFFR(3)+ANS(3)
@@ -1034,7 +1034,7 @@
 
       DO K=1,3
         ANS(K)=ANS(K)/DBLE(IDEN)
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.-1) THEN
           ANS(K)=ANS(K)*HELAVGFACTOR
         ELSE
           DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_default/ddx_ttx/born_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_default/ddx_ttx/born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_default/ddx_ttx/born_matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -118,10 +117,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -134,11 +129,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_default/ddx_ttx/loop_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_default/ddx_ttx/loop_matrix.f	2017-08-15 19:34:40 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_default/ddx_ttx/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -492,7 +492,7 @@
 C      helicity is asked
       NROTATIONS_DP_BU = NROTATIONS_DP
       NROTATIONS_QP_BU = NROTATIONS_QP
-      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.NE.-1) THEN
+      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.GT.-1) THEN
         IF(.NOT.WARNED_LORENTZ_STAB_TEST_OFF) THEN
           WRITE(*,*) '##WARNING: Evaluation of a specific helicity was'
      $     //' asked for this PS point, and there is a spin-2 (or'
@@ -531,7 +531,7 @@
         ENDDO
  101    CONTINUE
         CLOSE(1)
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=.TRUE.
@@ -636,7 +636,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(.NOT.GOODHEL(USERHEL)) THEN
             ANS(1)=0.0D0
             ANS(2)=0.0D0
@@ -644,7 +644,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.GT.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -888,7 +888,7 @@
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
 C       Update of NEXTREF, will be used for loop induced only.
         NEXTREF = NEXTREF + ANS(1) + ANS(2) + ANS(3)
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
           BUFFR(1)=BUFFR(1)+ANS(1)
           BUFFR(2)=BUFFR(2)+ANS(2)
           BUFFR(3)=BUFFR(3)+ANS(3)
@@ -1000,7 +1000,7 @@
 
       DO K=1,3
         ANS(K)=ANS(K)/DBLE(IDEN)
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.-1) THEN
           ANS(K)=ANS(K)*HELAVGFACTOR
         ELSE
           DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_default/gg_ttx/born_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_default/gg_ttx/born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_default/gg_ttx/born_matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -118,10 +117,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -134,11 +129,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_default/gg_ttx/loop_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_default/gg_ttx/loop_matrix.f	2017-08-14 07:04:15 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_default/gg_ttx/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -492,7 +492,7 @@
 C      helicity is asked
       NROTATIONS_DP_BU = NROTATIONS_DP
       NROTATIONS_QP_BU = NROTATIONS_QP
-      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.NE.-1) THEN
+      IF(MAX_SPIN_EXTERNAL_PARTICLE.GT.3.AND.USERHEL.GT.-1) THEN
         IF(.NOT.WARNED_LORENTZ_STAB_TEST_OFF) THEN
           WRITE(*,*) '##WARNING: Evaluation of a specific helicity was'
      $     //' asked for this PS point, and there is a spin-2 (or'
@@ -531,7 +531,7 @@
         ENDDO
  101    CONTINUE
         CLOSE(1)
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=.TRUE.
@@ -636,7 +636,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(.NOT.GOODHEL(USERHEL)) THEN
             ANS(1)=0.0D0
             ANS(2)=0.0D0
@@ -644,7 +644,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.GT.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1151,7 +1151,7 @@
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
 C       Update of NEXTREF, will be used for loop induced only.
         NEXTREF = NEXTREF + ANS(1) + ANS(2) + ANS(3)
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
           BUFFR(1)=BUFFR(1)+ANS(1)
           BUFFR(2)=BUFFR(2)+ANS(2)
           BUFFR(3)=BUFFR(3)+ANS(3)
@@ -1263,7 +1263,7 @@
 
       DO K=1,3
         ANS(K)=ANS(K)/DBLE(IDEN)
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.-1) THEN
           ANS(K)=ANS(K)*HELAVGFACTOR
         ELSE
           DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/ddx_ttx/born_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/ddx_ttx/born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/ddx_ttx/born_matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -118,10 +117,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -134,11 +129,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/ddx_ttx/loop_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/ddx_ttx/loop_matrix.f	2017-08-14 07:04:15 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/ddx_ttx/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1063,7 +1067,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1281,7 +1285,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/gg_ttx/born_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/gg_ttx/born_matrix.f	2018-05-18 12:22:47 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/gg_ttx/born_matrix.f	2019-05-23 07:08:37 +0000
@@ -80,7 +80,6 @@
       LOGICAL GOODHEL(NCOMB)
       DATA NTRY/0/
       DATA GOODHEL/NCOMB*.FALSE./
-
 C     
 C     GLOBAL VARIABLES
 C     
@@ -118,10 +117,6 @@
 C     ----------
 C     BEGIN CODE
 C     ----------
-      IF(USERHEL.EQ.-1) NTRY=NTRY+1
-      DO IHEL=1,NEXTERNAL
-        JC(IHEL) = +1
-      ENDDO
 C     When spin-2 particles are involved, the Helicity filtering is
 C      dangerous for the 2->1 topology.
 C     This is because depending on the MC setup the initial PS points
@@ -134,11 +129,16 @@
 C      helicity sum is).
 C     For this reason, we simply remove the filterin when there is
 C      only three external particles.
-      IF (NEXTERNAL.LE.3) THEN
+      IF(USERHEL.EQ.-2.OR.NEXTERNAL.LE.3) THEN
         DO IHEL=1,NCOMB
           GOODHEL(IHEL)=.TRUE.
         ENDDO
+        IF(USERHEL.EQ.-2)  USERHEL=-1
       ENDIF
+      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

=== modified file 'tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/gg_ttx/loop_matrix.f'
--- tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/gg_ttx/loop_matrix.f	2017-08-14 07:04:15 +0000
+++ tests/input_files/IOTestsComparison/short_ML_SMQCD_optimized/gg_ttx/loop_matrix.f	2019-05-23 07:08:37 +0000
@@ -636,6 +636,10 @@
         CALL ML5_0_SET_N_EVALS(N_DP_EVAL,N_QP_EVAL)
       ENDIF
 
+      IF(USERHEL.EQ.-2) THEN
+        NTRY = 0
+      ENDIF
+
       IF(NTRY.EQ.0) THEN
         HELDOUBLECHECKED=(.NOT.DOUBLECHECKHELICITYFILTER)
      $   .OR.(HELICITYFILTERLEVEL.EQ.0)
@@ -655,7 +659,7 @@
  101    CONTINUE
         CLOSE(1)
 
-        IF (HELICITYFILTERLEVEL.EQ.0) THEN
+        IF (HELICITYFILTERLEVEL.EQ.0.OR.USERHEL.EQ.-2) THEN
           FOUNDHELFILTER=.TRUE.
           DO J=1,NCOMB
             GOODHEL(J)=1
@@ -796,7 +800,7 @@
         HELPICKED=1
         CTMODE=CTMODEINIT
       ELSE
-        IF (USERHEL.NE.-1) THEN
+        IF (USERHEL.GT.0) THEN
           IF(GOODHEL(USERHEL).EQ.-HELOFFSET) THEN
             DO I=0,NSQUAREDSO
               ANS(1,I)=0.0D0
@@ -806,7 +810,7 @@
             GOTO 9999
           ENDIF
         ENDIF
-        HELPICKED=USERHEL
+        HELPICKED=MAX(USERHEL,-1)
         IF (CTMODERUN.NE.-1) THEN
           CTMODE=CTMODERUN
         ELSE
@@ -1063,7 +1067,7 @@
  1226 CONTINUE
 
       IF (CHECKPHASE.OR.(.NOT.HELDOUBLECHECKED)) THEN
-        IF((USERHEL.EQ.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
+        IF((USERHEL.LE.-1).OR.(USERHEL.EQ.HELPICKED)) THEN
 C         Make sure that that no polarization constraint filters out
 C          this helicity
           IF (POLARIZATIONS(0,0).EQ.-1.OR.ML5_0_IS_HEL_SELECTED(HELPICK
@@ -1281,7 +1285,7 @@
       DO K=1,3
         DO I=0,NSQUAREDSO
           ANS(K,I)=ANS(K,I)/DBLE(IDEN)
-          IF (USERHEL.NE.-1) THEN
+          IF (USERHEL.GT.-1) THEN
             ANS(K,I)=ANS(K,I)*HELAVGFACTOR
           ELSE
             DO J=1,NINITIAL

=== modified file 'tests/unit_tests/various/test_rambo.py'
--- tests/unit_tests/various/test_rambo.py	2015-10-01 16:00:08 +0000
+++ tests/unit_tests/various/test_rambo.py	2019-05-23 07:08:37 +0000
@@ -34,8 +34,6 @@
         # not enough energy
         self.assertRaises(rambo.RAMBOError, rambo.RAMBO, 2,150,self.m2)
         
-        # not valid mass
-        self.assertRaises(AssertionError, rambo.RAMBO, 2,1500,[1,2])
         
         # not valid mass
         self.assertRaises(AssertionError, rambo.RAMBO, 2,1500,self.m1)
@@ -57,6 +55,13 @@
         for i in range(1,3):
             self.assertAlmostEqual(self.m2[i]**2, P[(4,i)]**2 - P[(1,i)]**2 - P[(2,i)]**2 - P[(3,i)]**2)
         
+        mass = [100,20]
+        P, wgt = rambo.RAMBO(2,500, [100,20])
+        for i in range(1,3):
+            m = mass[i-1]
+            self.assertAlmostEqual(m**2, P[(4,i)]**2 - P[(1,i)]**2 - P[(2,i)]**2 - P[(3,i)]**2)
+
+
 
 class test_wavefunctions(unittest.TestCase):
     """ check that the wavefunctions TXXXX, IXXXXX, ... are correctly define"""

=== modified file 'tests/unit_tests/various/test_shower_card.py'
--- tests/unit_tests/various/test_shower_card.py	2019-04-24 11:02:30 +0000
+++ tests/unit_tests/various/test_shower_card.py	2019-05-23 07:08:37 +0000
@@ -527,7 +527,6 @@
         new_text = self.card.text
         old_card = shower_card.ShowerCard(new_text, testing=True)
         new_text = old_card.set_param('ue_enabled', '.false.', write_to=True)
-        misc.sprint(new_text)
         new_card = shower_card.ShowerCard(new_text, testing=True)
         self.assertFalse(new_card['ue_enabled'])