diff --git a/Trigger/TrigTools/TrigByteStreamTools/CMakeLists.txt b/Trigger/TrigTools/TrigByteStreamTools/CMakeLists.txt
index c9139c7ad17255cda70897b5e76f372d135a5abb..eea9de68f9c7e83064c5e001afb44a38f3531d7a 100644
--- a/Trigger/TrigTools/TrigByteStreamTools/CMakeLists.txt
+++ b/Trigger/TrigTools/TrigByteStreamTools/CMakeLists.txt
@@ -25,11 +25,24 @@ atlas_install_data( share/*.ref )
 # Unit tests
 set( trigbs_testFileRun2 "/cvmfs/atlas-nightlies.cern.ch/repo/data/data-art/TrigP1Test/data18_13TeV.00360026.physics_EnhancedBias.merge.RAW._lb0151._SFO-1._0001.1")
 set( trigbs_testFileRun3 "/cvmfs/atlas-nightlies.cern.ch/repo/data/data-art/TrigP1Test/data20_test/data_test.00374861.physics_Main.daq.RAW._lb0004._SFO-1._0001.data")
+
 atlas_add_test( dumpHLTContentInBS_run3
                 SCRIPT trigbs_dumpHLTContentInBS_run3.py -n 5 -s 105 --l1 --hlt --hltres --stag --sizes --sizeSummary ${trigbs_testFileRun3}
                 LOG_IGNORE_PATTERN "^----" )
+
 atlas_add_test( updateMetadata
                 SCRIPT rm -f newMeta.*.data
                        && trigbs_updateBSMetadata.py -n 5 --copyFrom ${trigbs_testFileRun2} --outputName newMeta ${trigbs_testFileRun3}
                        && meta-reader -m full newMeta._0001.data
                 LOG_IGNORE_PATTERN "file_guid" )
+
+atlas_add_test( mergeBSfiles
+                SCRIPT rm -rf mergedBS.* test_mergeBSfiles_input*
+                       && mkdir test_mergeBSfiles_input1 test_mergeBSfiles_input2
+                       && ln -s ${trigbs_testFileRun2} test_mergeBSfiles_input1/file1.data
+                       && ln -s ${trigbs_testFileRun3} test_mergeBSfiles_input1/file2.data
+                       && ln -s ${trigbs_testFileRun2} test_mergeBSfiles_input2/file1.data
+                       && ln -s ${trigbs_testFileRun3} test_mergeBSfiles_input2/file2.data
+                       && trigbs_mergeBSfiles.py -v -r -N 5 -n 10 mergedBS test_mergeBSfiles_input1 test_mergeBSfiles_input2
+                       && meta-reader mergedBS.*
+                LOG_IGNORE_PATTERN "file_guid" )
diff --git a/Trigger/TrigTools/TrigByteStreamTools/bin/trigbs_mergeBSfiles.py b/Trigger/TrigTools/TrigByteStreamTools/bin/trigbs_mergeBSfiles.py
index 0ac008f0acc45da303f56c1c0732de355c9ccc3b..1592cdec4212de0036c7abd21b7f69026b925053 100755
--- a/Trigger/TrigTools/TrigByteStreamTools/bin/trigbs_mergeBSfiles.py
+++ b/Trigger/TrigTools/TrigByteStreamTools/bin/trigbs_mergeBSfiles.py
@@ -2,95 +2,102 @@
 
 # Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
 
-from __future__ import print_function
+'''
+Merges events from multiple raw files from one or more directories.
+This is mainly for use on small samples
+'''
+
 import eformat
 import libpyevent_storage as EventStorage
 import sys
-import getopt
-import subprocess as sp
-
-eventsPerFile=500
-verbose=False
-renumber=False
-
-def pexec(args):
-    return sp.Popen(args, stdout=sp.PIPE,stderr=sp.STDOUT).communicate()[0].strip()
-
-
-def usage():
-  print(" usage: %s [options] <outbasename> <dir1> [<dir2> ..]" % sys.argv[0])
-  print("    Merges events from multiple raw files from one or more directories")
-  print("    This is mainly to used on small samples")
-  print()
-  print("    Options:")
-  print("       -n <number of events per file> : defaults to 500")
-  print("       -r     : renumber event ids from 1")
-  print("       -v     : verbose output")  
-  sys.exit(1)
-  
-try:
-  opts, args = getopt.getopt(sys.argv[1:], "vn:")
-  if len(args)<2:
-      usage()
-except getopt:
-  usage()
-
-for opt,arg in opts:
-  if opt == '-n':
-    eventsPerFile=int(arg)
-  if opt == '-v':
-    verbose=True
-  if opt == '-r':
-    renumber=True
-
-
-basename=args[0]
-dirs = args[1:]
-files=[]
-
-for indir in dirs:
-  cmd='ls'
-  if indir.startswith('/castor'):
-    cmd='nsls'
-  fileNames=pexec([cmd,indir]).split('\n')
-  for name in fileNames:
-    files+=[indir+'/'+name]
-
-print('Merging',len(files),'files')
-
-ifile = 0
-id = 1
-outnum=1
-dr=EventStorage.pickDataReader(files[0])
-numFiles=0
-output = eformat.ostream(core_name=basename,
-                         run_number=dr.runNumber(),
-                         trigger_type=dr.triggerType(),
-                         detector_mask=dr.detectorMask(),
-                         beam_type=dr.beamType(),
-                         beam_energy=dr.beamEnergy())
-output.max_events_per_file(eventsPerFile)
-
-for input_file in files:
-  ifile += 1
-  print("Opening file %d: %s" % (ifile, input_file))
-
-  input = eformat.istream([input_file])
-
-  for event in input:
-    if verbose:
-      print("Reading event with Lvl1 Id = %ld and rewritting with %ld" % (event.lvl1_id(),id))
-    new_event = eformat.write.FullEventFragment(event)
-    if renumber:
-        new_event.lvl1_id(id)
-        new_event.global_id(id)
-        for rob in new_event.children():
-            rob.rod_lvl1_id(id)
-            rob.checksum_type(0)  #remove checksum
-            rob.rod_lvl1_trigger_type(event.lvl1_trigger_type())
-
-    output.write(new_event)
-
-    id+=1
-
-print('Wrote %d events to %d files' % ( id-1, 1+((id-1)//eventsPerFile) ))
+import os
+import argparse
+import logging
+
+
+def get_parser():
+    parser = argparse.ArgumentParser(usage='%(prog)s [options] outBaseName dir1 [dir2 ...]',
+                                     description=__doc__,
+                                     formatter_class=lambda prog: argparse.HelpFormatter(prog, max_help_position=36, width=100))
+    parser.add_argument('outBaseName',
+                        help='Output file name base')
+    parser.add_argument('dirs',
+                        metavar='dir',
+                        nargs='+',
+                        help='Input directory (can be multiple)')
+    parser.add_argument('--copyFrom', '-c',
+                        metavar='FILE', action='store',
+                        help='Copy metadata from other ByteStream file')
+    parser.add_argument('--nepfOutput', '-n',
+                        metavar='NUM', action='store', type=int, default=500,
+                        help='Number of events per output file (default=%(default)s)')
+    parser.add_argument('--nepfInput', '-N',
+                        metavar='NUM', action='store', type=int, default=0,
+                        help='Number of events per input file, 0 means all (default=%(default)s)')
+    parser.add_argument('--renumber', '-r',
+                        action='store_true',
+                        help='Renumber event IDs from 1')
+    parser.add_argument('--verbose', '-v',
+                        action='store_true',
+                        help='Increase output verbosity')
+    return parser
+
+
+def get_filenames(dirs):
+    filenames = []
+    for dirname in sorted(dirs):
+        files = os.listdir(dirname)
+        logging.debug('Found %d files in %s', len(files), dirname)
+        for f in sorted(files):
+            filenames.append(dirname + '/' + f)
+    return filenames
+
+
+def main():
+    args = get_parser().parse_args()
+    logging.basicConfig(stream=sys.stdout,
+                        format='%(levelname)-8s %(message)s',
+                        level=logging.DEBUG if args.verbose else logging.INFO)
+
+    filenames = get_filenames(args.dirs)
+    logging.info('Merging %d files', len(filenames))
+    logging.info('filenames: %s', filenames)
+
+    ievt = 1
+    dr = EventStorage.pickDataReader(filenames[0])
+    output = eformat.ostream(core_name=args.outBaseName,
+                             run_number=dr.runNumber(),
+                             trigger_type=dr.triggerType(),
+                             detector_mask=dr.detectorMask(),
+                             beam_type=dr.beamType(),
+                             beam_energy=dr.beamEnergy())
+    output.max_events_per_file(args.nepfOutput)
+
+    for ifile, input_file in enumerate(filenames):
+        logging.info("Opening file %d: %s", ifile+1, input_file)
+        inputFile = eformat.istream([input_file])
+
+        for ievt_in_file, event in enumerate(inputFile):
+            if args.nepfInput > 0 and ievt_in_file >= args.nepfInput:
+                break
+            logging.debug("Reading event with Lvl1 Id = %ld", event.lvl1_id())
+            new_event = eformat.write.FullEventFragment(event)
+            if args.renumber:
+                logging.debug("Rewriting the event with Lvl1 Id = %ld", ievt)
+                new_event.lvl1_id(ievt)
+                new_event.global_id(ievt)
+                for rob in new_event.children():
+                    rob.rod_lvl1_id(ievt)
+                    rob.checksum_type(0)  # remove checksum
+                    rob.rod_lvl1_trigger_type(event.lvl1_trigger_type())
+            else:
+                logging.debug("Rewriting the event with unchanged Lvl1 Id")
+
+            output.write(new_event)
+            ievt += 1
+
+    logging.info('Wrote %d events to %d files', ievt-1, 1+((ievt-2)//args.nepfOutput))
+
+
+if "__main__" in __name__:
+    sys.exit(main())
diff --git a/Trigger/TrigTools/TrigByteStreamTools/share/mergeBSfiles.ref b/Trigger/TrigTools/TrigByteStreamTools/share/mergeBSfiles.ref
new file mode 100644
index 0000000000000000000000000000000000000000..224ff3986aa8f28f572b826a054c26779e46e5ea
--- /dev/null
+++ b/Trigger/TrigTools/TrigByteStreamTools/share/mergeBSfiles.ref
@@ -0,0 +1,86 @@
+DEBUG    Found 2 files in test_mergeBSfiles_input1
+DEBUG    Found 2 files in test_mergeBSfiles_input2
+INFO     Merging 4 files
+INFO     filenames: ['test_mergeBSfiles_input1/file1.data', 'test_mergeBSfiles_input1/file2.data', 'test_mergeBSfiles_input2/file1.data', 'test_mergeBSfiles_input2/file2.data']
+INFO     Opening file 1: test_mergeBSfiles_input1/file1.data
+DEBUG    Reading event with Lvl1 Id = 1745102069
+DEBUG    Rewriting the event with Lvl1 Id = 1
+DEBUG    Reading event with Lvl1 Id = 1745110878
+DEBUG    Rewriting the event with Lvl1 Id = 2
+DEBUG    Reading event with Lvl1 Id = 1745129028
+DEBUG    Rewriting the event with Lvl1 Id = 3
+DEBUG    Reading event with Lvl1 Id = 1745120007
+DEBUG    Rewriting the event with Lvl1 Id = 4
+DEBUG    Reading event with Lvl1 Id = 1745115291
+DEBUG    Rewriting the event with Lvl1 Id = 5
+INFO     Opening file 2: test_mergeBSfiles_input1/file2.data
+DEBUG    Reading event with Lvl1 Id = 1762019442
+DEBUG    Rewriting the event with Lvl1 Id = 6
+DEBUG    Reading event with Lvl1 Id = 1778785142
+DEBUG    Rewriting the event with Lvl1 Id = 7
+DEBUG    Reading event with Lvl1 Id = 1795394603
+DEBUG    Rewriting the event with Lvl1 Id = 8
+DEBUG    Reading event with Lvl1 Id = 1795315619
+DEBUG    Rewriting the event with Lvl1 Id = 9
+DEBUG    Reading event with Lvl1 Id = 1795331759
+DEBUG    Rewriting the event with Lvl1 Id = 10
+INFO     Opening file 3: test_mergeBSfiles_input2/file1.data
+DEBUG    Reading event with Lvl1 Id = 1745102069
+DEBUG    Rewriting the event with Lvl1 Id = 11
+DEBUG    Reading event with Lvl1 Id = 1745110878
+DEBUG    Rewriting the event with Lvl1 Id = 12
+DEBUG    Reading event with Lvl1 Id = 1745129028
+DEBUG    Rewriting the event with Lvl1 Id = 13
+DEBUG    Reading event with Lvl1 Id = 1745120007
+DEBUG    Rewriting the event with Lvl1 Id = 14
+DEBUG    Reading event with Lvl1 Id = 1745115291
+DEBUG    Rewriting the event with Lvl1 Id = 15
+INFO     Opening file 4: test_mergeBSfiles_input2/file2.data
+DEBUG    Reading event with Lvl1 Id = 1762019442
+DEBUG    Rewriting the event with Lvl1 Id = 16
+DEBUG    Reading event with Lvl1 Id = 1778785142
+DEBUG    Rewriting the event with Lvl1 Id = 17
+DEBUG    Reading event with Lvl1 Id = 1795394603
+DEBUG    Rewriting the event with Lvl1 Id = 18
+DEBUG    Reading event with Lvl1 Id = 1795315619
+DEBUG    Rewriting the event with Lvl1 Id = 19
+DEBUG    Reading event with Lvl1 Id = 1795331759
+DEBUG    Rewriting the event with Lvl1 Id = 20
+INFO     Wrote 20 events to 2 files
+Py:MetaReader        INFO Current mode used: lite
+Py:MetaReader        INFO Current filenames: ['mergedBS._0001.data', 'mergedBS._0002.data']
+|-- mergedBS._0001.data: 
+|   |-- GeoAtlas          : None
+|   |-- beam_energy       : 6500
+|   |-- beam_type         : collisions
+|   |-- conditions_tag    : None
+|   |-- eventTypes        : ['IS_DATA', 'IS_ATLAS', 'Unknown']
+|   |-- evt_number        : [1]
+|   |-- file_guid         : 8EC3DF82-9113-EB11-8086-0242AC110002
+|   |-- file_size         : 10513140
+|   |-- file_type         : BS
+|   |-- lumiBlockNumbers  : [0]
+|   |-- nentries          : 10
+|   |-- processingTags    : ['EnhancedBias']
+|   |-- project_name      : 
+|   |-- runNumbers        : [360026]
+|   |-- run_type          : ['PHYSICS']
+|   `-- stream            : 
+`-- mergedBS._0002.data: 
+    |-- GeoAtlas          : None
+    |-- beam_energy       : 6500
+    |-- beam_type         : collisions
+    |-- conditions_tag    : None
+    |-- eventTypes        : ['IS_DATA', 'IS_ATLAS', 'Unknown']
+    |-- evt_number        : [11]
+    |-- file_guid         : 606A7083-9113-EB11-8086-0242AC110002
+    |-- file_size         : 10513172
+    |-- file_type         : BS
+    |-- lumiBlockNumbers  : [0]
+    |-- nentries          : 10
+    |-- processingTags    : ['EnhancedBias']
+    |-- project_name      : 
+    |-- runNumbers        : [360026]
+    |-- run_type          : ['PHYSICS']
+    `-- stream            : 
+