wiki:FAQ-General-15

How can I remove some unwanted diagram?

Warning

Be carefull about gauge invariance (a non gauge invariant subset will also break lorentz invariance) Additionally, this method is only implemented for leading order processes

method:

  1. You need to add a file PLUGIN/user_filter.py this file should define a function named "remove_diag" (see example below). The function should return True if you want to discard the diagram and False if you want to keep it.
  2. When running the code, you have to add the flag "--diagram_filter"
    generate p p > w+ w+ j j --diagram_filter
    
    If you remove any diagram you should see on the screen:
    INFO: Trying process: c c > w+ w+ s s WEIGHTED<=6 @1
    WARNING: Diagram filter is ON and removed 12 diagrams for this subprocess.
    

example of user_filter.py for version of the code for 2.7.3 and later (and 3.0.2 and later)

example using the drawing module

import madgraph.core.drawing as drawing

def remove_diag(diag, model):
    """remove all diagram with quark in T-channel"""

    #diag is a Direct Accyclic Graph of the Feynman diagram
    #convert it to a full graph (the one we used for the plotting)
    #In that representation each vertex is associated to all the legs/particle and we can move quite freely inside the diagram
    # The structure representation is simpler and should help to apply complex filtering

    draw = drawing.FeynmanDiagram(diag, model)
    draw.load_diagram()
    print(draw._debug_load_diagram())
    
    # Diagram content  can be loop over three list:
    #  - all vertex of diagram (note that for each external particle we have a vertex with a single particle attached to it)
    #    for v in draw.vertexList:
    #  - all vertex corresponding to the initial/final state particles
    #    for v in draw.initial_vertex:
    # - all the particles (initial states / final states and propagator) of the diagram
    #    for p in draw.lineList
    #
    # All vertex can be assigned to a level by running
    draw.define_level()
    #      0: correspond to the initial state vertex
    #      1: correspond to the T-channel level
    #   you can use draw._debug_level() to text representation in that case
    #   For a single line the "begin" level will always be one level lower than the "end" level
    #    BUT for t-channel where both are equal and set to 1
    #print(draw._debug_level())
    #
    # All vertex are of type VertexPoint
    #   They have only two relevant attributes
    #      self.lines : all the particle connecting to this vertex
    #      self.level : level of this vertex (see above)
    #
    # All particles are of type FeynmanLine
    #   They have the following relevant attributes
    #       self.begin: vertex associated to the beginning of the particles
    #       self.end:  vertex associated to the end of the particles
    #       self.get('id'): pdg code
    #       self.get('number'): number associated to the original DAG
    #       other attributes are ['polarization', 'number', 'onshell', 'state', 'model', 'from_group', 'loop_line']

    # remove all 4 point interaction
    for v in draw.vertexList:
        if len(v.lines) > 3:
            return True
        
    # check if initial/final  states quark are emitting a gluon
    for external in draw.initial_vertex: # vertex with only one leg for initial/final state particles
        # external is VertexPoint object
        p = external.lines[0] # FeynmanLine object
        if abs(p.id) < 7: # PDG code less than 7 means a SM quark
            other_vertex = p.end # VertexPoint object
            if other_vertex == external:
                other_vertex = p.begin 
            if any(p2.id ==21 for p2 in other_vertex.lines):
                return True      
    return False

Additionally all the example below are still working but requires one more argument "model" "

example of user_filter.py for version of the code before 2.7.3 (not included).

def remove_diag(diag, model=None):
    """remove all diagram with quark in T-channel"""
    for vertex in diag['vertices']: #last
         if vertex['id'] == 0: #special final vertex
             continue
         if vertex['legs'][-1]['number'] < 3: #this means T-channel
             if abs(vertex['legs'][-1]['id']) <6:
                  return True
    return False
def remove_diag(diag,model=None):
    """remove all diagram which are not pure QCD or pure QED: Fake example"""

    if diag.get('orders')['QCD'] >0 and diag.get('orders')['QED']>0:
        return True

    return False
import madgraph.various.misc as misc

def remove_diag(diag,model=None):
    """force to have two s-channel W decaying into a specific way. 
       Designed for interference studies where decay chain are not allowed
       This is designed for keeping the following diagram generate p p > j j w- w-, w- > e- ve~ 
       from the following syntax:  p p > j j e- ve~ e- ve~
"""

    found_w = 0
    for vertex in diag['vertices']:
        #misc.sprint(vertex.nice_string())
        if vertex.get('id') == 0: #special final vertex
             continue
        if len(vertex['legs'])!=3:
            continue

        for i, leg in enumerate(list(vertex['legs'])):
            #misc.sprint(i, leg)
            if abs(leg['id']) != 24:
                continue
            else:
                break
        else:
            continue
        # we have a vertex with a W
        if i == 2: #This means a vertex A B > W or W > A B
            if leg['number'] > 2: # means that we have W > A B
                #misc.sprint(vertex.nice_string())
                # check that A B are what we want
                if not all(oleg['state'] for oleg in vertex['legs'][:2]) : # this means not final state
                    continue
                #misc.sprint([oleg['id'] for oleg in vertex['legs'][:2] ])
                if set([oleg['id'] for oleg in vertex['legs'][:2] ]) == set([11,-12]):
                    found_w +=1
        else: # this means A W > B
            if vertex['legs'][-1].get('number') < 3:
                # T-channel -> continue
                continue
            else:
                #check that they are both final state
                if not vertex['legs'][-1]['state']:
                    continue
                if not vertex['legs'][1 if i==0 else 0]['state']:
                    continue

                other_id = [vertex['legs'][-1].get('id'), vertex['legs'][1 if i==0 else 0]['id']]
                if set(other_id) == set([11,-12]):
                    found_w += 1

    misc.sprint(found_w)
    if found_w == 2:
        return True
    else:
        return False

Last modified 5 years ago Last modified on Jun 20, 2020, 8:45:58 PM
Note: See TracWiki for help on using the wiki.