From 13e24162e34202808788b015b8e4df9b68335ba9 Mon Sep 17 00:00:00 2001
From: Eduardo Rodrigues <eduardo.rodrigues@cern.ch>
Date: Wed, 31 Jan 2018 22:25:49 +0000
Subject: [PATCH] Merge branch 'decianm-TurboSPTeslaMCConfiguration' into
 '2017-patches'

Trying to get TurboSP and MC with truthmatching to work

See merge request lhcb/DaVinci!101

(cherry picked from commit bf52f8bf967c80c3990315aeabb2d91b2affd876)

87e70041 First commit of Configuration file, trying to get TurboSP MC to run with truthmatching
8e646c9a fix some typos
6fc9415d Commit latest try for neutral association in TurboSP
d6b201cf Add working neutral-matching logic, tidy up.
67d09999 Working version.
c4797d9e Remove temporary comments.
---
 Phys/Tesla/python/Tesla/Configuration.py | 212 +++++++++++++++++++++--
 1 file changed, 195 insertions(+), 17 deletions(-)

diff --git a/Phys/Tesla/python/Tesla/Configuration.py b/Phys/Tesla/python/Tesla/Configuration.py
index 933a3468f..4e0aef1f1 100644
--- a/Phys/Tesla/python/Tesla/Configuration.py
+++ b/Phys/Tesla/python/Tesla/Configuration.py
@@ -21,6 +21,8 @@ from Configurables import (
     DecodeRawEvent,
     DstConf,
     EventSelector,
+    Gaudi__DataCopy as DataCopy,
+    Gaudi__DataLink as DataLink,
     GaudiSequencer,
     HltPackedDataDecoder,
     HltRoutingBitsFilter,
@@ -68,6 +70,36 @@ from LHCbKernel.Configuration import *
 from RawEventCompat.Configuration import ReverseDict as RawEventLocationsToBanks
 
 
+def _strip_root(root, locs):
+    """Return list of paths `loc` stripped of the `root` prefix.
+
+    Examples
+    --------
+
+        >>> _strip_root('/Event/Turbo', ['/Event/Turbo/Foo'])
+        'Foo'
+        >>> _strip_root('/Event', ['/Event/Turbo/Foo'])
+        'Turbo/Foo'
+        >>> # Path is already relative
+        >>> _strip_root('/Event', ['Turbo/Foo'])
+        'Turbo/Foo'
+        >>> # Complains if the prefix doesn't match the root
+        >>> _strip_root('/Event/Ultra', ['/Event/Turbo/Foo'])
+        AssertionError
+    """
+    root = root.rstrip('/')
+    ret = []
+    for l in locs:
+        # Only need to modify paths that are already absolute
+        if l.startswith('/'):
+            # Strip off the root the prefixes the string
+            l = l.replace(root + '/', '')
+            assert not l.startswith('/'), \
+                   'TES path still prefixed: {0}'.format(l)
+        ret.append(l)
+    return ret
+
+
 class Tesla(LHCbConfigurableUser):
     __used_configurables__ = [ LHCbApp, LumiAlgsConf, RawEventJuggler, DecodeRawEvent, RawEventFormatConf, DstConf, RecombineRawEvent, PackParticlesAndVertices, TrackAssociator, ChargedPP2MC, TrackSys, TurboConf]
 
@@ -172,6 +204,7 @@ class Tesla(LHCbConfigurableUser):
     writerName = "DstWriter"
     teslaSeq = GaudiSequencer("TeslaSequence")
     base = "Turbo/"
+    neutral_pp2mc_loc = "Relations/Turbo/NeutralPP2MC"
 
 
     def _safeSet(self,conf,param):
@@ -211,7 +244,7 @@ class Tesla(LHCbConfigurableUser):
         SequencerTimerTool().OutputLevel = WARNING
 
     def _configureTrackTruth(self,assocpp,trackLoc) :
-        assoctr = TrackAssociator("TurboAssocTr"+trackLoc)
+        assoctr = TrackAssociator("TurboAssocTr"+trackLoc.replace("/", ""))
         assoctr.OutputLevel = self.getProp('OutputLevel')
         assoctr.TracksInContainer = trackLoc
         assocpp.TrackLocations += [ trackLoc ]
@@ -222,7 +255,7 @@ class Tesla(LHCbConfigurableUser):
         assocdigits = CaloDigit2MCLinks2Table("TurboDigitAssoc")
         return
 
-    def _configureClustersAndProtosTruth(self,digits,clusters,protos) :
+    def _configureClustersAndProtosTruth(self,digits,clusters,protos,protoTabLoc):
         retSeq = GaudiSequencer("NeutralTruth")
         clusterTabLoc = self.base + "Relations/CaloClusters"
 
@@ -232,11 +265,6 @@ class Tesla(LHCbConfigurableUser):
         assoccluster.Output = clusterTabLoc
         assoccluster.Clusters+=clusters
 
-        # When filtering MC, the relations table cloners will copy the tables
-        # to /Event/Turbo for us. Otherwise we create them there directly
-        protoTabLoc = "Relations/Turbo/NeutralPP2MC"
-        if not self.getProp("FilterMC"):
-            protoTabLoc = os.path.join(self.base, protoTabLoc)
         assocneutral = NeutralPP2MC("TurboNeutralPP2MC")
         assocneutral.InputData += protos
         assocneutral.OutputLevel = self.getProp('OutputLevel')
@@ -244,7 +272,7 @@ class Tesla(LHCbConfigurableUser):
         assocneutral.MCCaloTable = clusterTabLoc
 
         retSeq.Members += [assoccluster,assocneutral]
-        return protoTabLoc, retSeq
+        return retSeq
 
     def _unpackMC(self):
         DataOnDemandSvc().NodeMap['/Event/MC']   = 'DataObject'
@@ -267,10 +295,8 @@ class Tesla(LHCbConfigurableUser):
         microDST algorithms, using the relations tables Tesla makes to find the
         MC particles that have been associated to the reconstruction.
         """
-        if not self.getProp('FilterMC'):
-            return
-
-        output_prefix = self.base.replace('/', '')
+        output_prefix = self.base.rstrip('/')
+        output_level = self.getProp('OutputLevel')
 
         # Algorithm to clone the existing relations tables and the MC particles
         # and vertices these tables point to from their original location to
@@ -281,6 +307,7 @@ class Tesla(LHCbConfigurableUser):
         # Don't use the assumed copy of the input ProtoParticle objects, as
         # Tesla doesn't need to copy them (they're already under /Event/Turbo)
         copy_pp2mcp.UseOriginalFrom = True
+        copy_pp2mcp.OutputLevel = output_level
 
         # Algorithm to clone all MC particles and vertices that are associated
         # to the simulated signal process (using LHCb::MCParticle::fromSignal)
@@ -288,6 +315,7 @@ class Tesla(LHCbConfigurableUser):
         copy_signal_mcp.OutputPrefix = output_prefix
         # Don't copy the associated reconstruction, as Tesla already saves it
         copy_signal_mcp.SaveAssociatedRecoInfo = False
+        copy_signal_mcp.OutputLevel = output_level
 
         algorithms = [copy_pp2mcp, copy_signal_mcp]
         seq = GaudiSequencer('TurboMCFiltering', Members=algorithms)
@@ -494,7 +522,10 @@ class Tesla(LHCbConfigurableUser):
 
             # Configure
             self._configureDigitsTruth()
-            protoneutral, retSeq = self._configureClustersAndProtosTruth(outputDigiLoc,neutralClustersTot,neutralProtosTot)
+            protoneutral = self.neutral_pp2mc_loc
+            if not self.getProp('FilterMC'):
+                protoneutral = os.path.join(tesROOT, protoneutral)
+            retSeq = self._configureClustersAndProtosTruth(outputDigiLoc,neutralClustersTot,neutralProtosTot,protoneutral)
             NeutralProtoSeq.Members+=[retSeq]
 
             # Add standard Turbo locations if not packing
@@ -591,6 +622,113 @@ class Tesla(LHCbConfigurableUser):
         # Add the output sequence
         TeslaReportAlgoSeq.Members += [TeslaOutputStreamsSequence]
 
+    def _configureTruthMatching(self, name, packing):
+        output_level = self.getProp('OutputLevel')
+
+        # Locations of the objects we want to truth-match
+        track_locs = [packing.outputs['Hlt2LongTracks'],
+                      packing.outputs['Hlt2DownstreamTracks']]
+        proto_locs = [packing.outputs['Hlt2LongProtos'],
+                      packing.outputs['Hlt2DownstreamProtos']]
+        cluster_locs = [packing.outputs['Hlt2CaloClusters']]
+        neutral_locs = [packing.outputs['Hlt2NeutralProtos']]
+
+        assocpp = ChargedPP2MC("TurboProtoAssocPP")
+        # Configure Track -> MCParticle relations table maker
+        trackTruthSeq = GaudiSequencer('TurboTrackTruthSequencer')
+        trackTruthSeq.Members = [self._configureTrackTruth(assocpp, loc)
+                                 for loc in track_locs]
+
+        # Configure ProtoParticle -> MCParticle relations table maker
+        assocpp.OutputLevel = output_level
+        assocpp.VetoEmpty = True
+        assocpp.InputData += _strip_root('/Event', proto_locs)
+
+        assocDigits = CaloDigit2MCLinks2Table("TurboDigitAssoc")
+        # The digits relations table maker will try to find the linker table
+        # for digits at `Raw/{det}/Digits` at `Link/Raw/{det}/Digits`. Our
+        # digits are at `Turbo/Raw/{det}/Digits`, but we can't remake a linker
+        # table at `Link/Turbo/Raw/{det}/Digits because the information needed
+        # to make one is lost after Boole.
+        # Instead, exploit the fact that the digits under `Turbo/` are a subset
+        # of the originals, and just symlink the original linker table to the
+        # location expected for the Turbo digits
+        assocDigits.RespectInputDigitLocation = False
+        digit_locs = _strip_root('/Event', assocDigits.getDefaultProperty('Inputs'))
+        turbo_digit_locs = [os.path.join(self.base, l) for l in digit_locs]
+        linkingSeq = GaudiSequencer('TurboLinkTableLinkers')
+        for digit_loc in digit_locs:
+            link_loc = os.path.join('Link', digit_loc)
+            turbo_link_loc = os.path.join('Link', self.base, digit_loc)
+            name = '{0}Linker'.format(link_loc.replace('/', ''))
+            dl = DataLink(name, What=link_loc, Target=turbo_link_loc)
+            linkingSeq.Members.append(dl)
+            # DataLinks don't create the ancestor TES structure, so do it
+            # manually
+            path = ''
+            for node in turbo_link_loc.split(os.path.sep)[:-1]:
+                path = os.path.join(path, node)
+                DataOnDemandSvc().NodeMap[path] = 'DataObject'
+
+        # Configure CaloDigit -> MCParticle relations table maker
+        assocDigits.OutputLevel = output_level
+        # Can use a magic string here as this table is only temporary
+        assocDigits.Output = "Relations/CaloDigits"
+        assocDigits.Inputs = turbo_digit_locs
+
+        # Configure CaloCluster -> MCParticle and neutral ProtoParticle ->
+        # MCParticle relations table makers
+        neutral_rels_loc = self.neutral_pp2mc_loc
+        assocClustersNeutrals = self._configureClustersAndProtosTruth(
+            assocDigits.Output,
+            cluster_locs,
+            neutral_locs,
+            neutral_rels_loc
+        )
+
+        chargedProtoSeq = GaudiSequencer("ChargedTruthSequencer")
+        chargedProtoSeq.Members = [trackTruthSeq, assocpp]
+        neutralProtoSeq = GaudiSequencer("NeutralTruthSequencer")
+        neutralProtoSeq.Members = [linkingSeq, assocDigits, assocClustersNeutrals]
+        truthSeq = GaudiSequencer("TurboTruthSequencer")
+        truthSeq.Members = [chargedProtoSeq, neutralProtoSeq]
+
+        # We have now created the relations tables under /Event/Relations,
+        # and want them under /Event/Turbo/Relations. We have two strategies:
+        #   1. Use the CopyProtoParticle2MCRelations algorithm, which will copy
+        #   the tables *and* the MC particles they reference to under
+        #   /Event/Turbo, updating the copied tables to point to the copied MC
+        #   particles; or
+        #   2. Use Gaudi::DataCopy to copy the tables verbatim, so that they
+        #   still refer to MC objects in /Event/MC.
+        # The first case we call 'filtering', as in principle the original MC
+        # objects can now be deleted because we then only depend on those in
+        # /Event/Turbo/MC. The second case we choose to do when not filtering
+        # so that the paths to the final relations tables are the same as those
+        # in the first case.
+        relationsLocations = [os.path.join('Relations', path)
+                              for path in assocpp.InputData]
+        relationsLocations.append(neutral_rels_loc)
+        if self.getProp('FilterMC'):
+            filterMCSeq = self._filterMCParticlesSequence(relationsLocations)
+            truthSeq.Members += [filterMCSeq]
+        else:
+            base = self.base.rstrip('/')
+            truthSeq.Members.append(
+                GaudiSequencer('CopyRelationsTables', Members=[
+                    DataCopy('Copy{0}'.format(loc.replace('/', '')),
+                             What=loc,
+                             Target=os.path.join(base, loc),
+                             NeverFail=False,
+                             OutputLevel=output_level)
+                    for loc in relationsLocations
+                    if self.base in loc
+                ], IgnoreFilterPassed=True)
+            )
+
+        return truthSeq
+
+
     def _configureOutputStream(self, stream, lines_dict):
         '''
         The lines dictionary should have the following form:
@@ -830,7 +968,7 @@ class Tesla(LHCbConfigurableUser):
             )]
 
         # Kill any links the output file might have had to the input file
-        if online:
+        if online and not simulation:
             address_killer = [AddressKillerAlg()]
         else:
             address_killer = []
@@ -869,6 +1007,9 @@ class Tesla(LHCbConfigurableUser):
         pack = self.getProp('Pack')
         write_fsr = self.getProp('WriteFSR')
         simulation = self.getProp('Simulation')
+        filter_mc = self.getProp('FilterMC')
+        if filter_mc:
+            assert simulation, 'Expected Simulation = True as FilterMC = True'
 
         stream_seq = GaudiSequencer(namer('TeslaStreamSequence'))
 
@@ -901,6 +1042,10 @@ class Tesla(LHCbConfigurableUser):
             ]
         ))
 
+        if simulation:
+            self._unpackMC()
+            stream_seq.Members.append(self._configureTruthMatching(name, packing))
+
         # No need to clone if the output prefix is empty (everything stays
         # under /Event/Turbo)
         if output_prefix:
@@ -1015,6 +1160,26 @@ class Tesla(LHCbConfigurableUser):
                 )
                 packers.append(psandvs_packer)
 
+            # Pack the filtered MC particles and vertices
+            if filter_mc:
+                mcp_packer = PackMCParticle(
+                    'TurboPackMCParticle',
+                    AlwaysCreateOutput=True,
+                    DeleteInput=False,
+                    RootInTES=turbo_base
+                )
+                mcv_packer = PackMCVertex(
+                    'TurboPackMCVertex',
+                    AlwaysCreateOutput=True,
+                    DeleteInput=False,
+                    RootInTES=turbo_base
+                )
+
+                packers.append(mcp_packer)
+                packers.append(mcv_packer)
+
+                optional_output_locations += [os.path.join(turbo_base, 'pSim#*')]
+
             for packer in packers:
                 packer.OutputLevel = min(self.getProp('OutputLevel'),
                                          self.getProp('PackerOutputLevel'))
@@ -1031,6 +1196,12 @@ class Tesla(LHCbConfigurableUser):
                 os.path.join(turbo_base, 'pRec#*'),
                 os.path.join(turbo_base, 'Hlt2/pRec#*')
             ]
+
+            # Need to save the *unpacked* relations tables when we aren't
+            # streaming because we don't run PackParticlesAndVertices in that
+            # case
+            if simulation and not output_prefix:
+                optional_output_locations += [os.path.join(turbo_base, 'Relations#*')]
         else:
             # Save everything under the stream prefix
             optional_output_locations += [
@@ -1046,7 +1217,7 @@ class Tesla(LHCbConfigurableUser):
         else:
             fname = self.getProp('outputFile')
 
-        if online:
+        if online and not simulation:
             writer = OutputStream(namer(self.writerName))
             writer.AcceptAlgs += ['LumiSeq', 'PhysFilter']
             # In online mode we perform the raw event juggling (otherwise
@@ -1056,8 +1227,10 @@ class Tesla(LHCbConfigurableUser):
             # In offline mode, e.g. after Brunel, propagate everything from the
             # input file to the output
             writer = InputCopyStream(namer(self.writerName))
+
         writer.ItemList = required_output_locations
         writer.OptItemList = optional_output_locations
+
         # Keep the 2016 behaviour by not writing out anything if the
         # sel reports are missing
         writer.RequireAlgs = [self._selReportsCheck(), stream_seq]
@@ -1153,9 +1326,10 @@ class Tesla(LHCbConfigurableUser):
 
         ApplicationMgr().AppName="Tesla, utilising DaVinci"
         #
-        if self.getProp('Mode') is "Online":
+        if self.getProp('Mode') == "Online":
             self.setProp('WriteFSR',True)
-            self._configureLumi()
+            if not self.getProp('Simulation'):
+                self._configureLumi()
         else:
             DecodeRawEvent().DataOnDemand=True
             RecombineRawEvent(Version=self.getProp('SplitRawEventInput'))
@@ -1166,6 +1340,10 @@ class Tesla(LHCbConfigurableUser):
                 # /Event/Turbo, else it will be packed
                 TurboConf().setProp("RootInTES", "/Event")
 
+        # The MC logic requires the DataOnDemandSvc to be enabled
+        if self.getProp('Simulation'):
+            self._configureDataOnDemand()
+
         if self.getProp('DataType') == '2017' and not self.getProp('CandidatesFromSelReports'):
             # Unpack the DstData raw bank
             TurboConf().setProp("RunPackedDataDecoder", True)
-- 
GitLab