diff --git a/decaylanguage/__init__.py b/decaylanguage/__init__.py
index 4e193864..6b1df889 100644
--- a/decaylanguage/__init__.py
+++ b/decaylanguage/__init__.py
@@ -1,6 +1,11 @@
+from __future__ import absolute_import
+
 
 __version__ = '0.2.0'
 
 version = __version__
 version_info = __version__.split('.')
 
+
+# Direct access to decay file parsing tools
+from .dec import DecFileParser
diff --git a/decaylanguage/data/decfile.lark b/decaylanguage/data/decfile.lark
index 618f378a..c191bdac 100644
--- a/decaylanguage/data/decfile.lark
+++ b/decaylanguage/data/decfile.lark
@@ -44,7 +44,7 @@ model_options : (value | LABEL)+
 _NEWLINE: ( /\r?\n[\t ]*/ | COMMENT )+
 
 // We must set priorities here to use lalr - match model name above label, and label above something else
-MODEL_NAME.2 : "BaryonPCR"|"BTO3PI_CP"|"BTOSLLALI"|"BTOSLLBALL"|"BTOXSGAMMA"|"BTOXSLL"|"CB3PI-MPP"|"CB3PI-P00"|"D_DALITZ"|"ETA_DALITZ"|"GOITY_ROBERTS"|"HELAMP"|"HQET"|"ISGW2"|"OMEGA_DALITZ"|"PARTWAVE"|"PHSP"|"PI0_DALITZ"|"PYTHIA"|"SLN"|"STS"|"SVP_HELAMP"|"SVS"|"SVV_HELAMP"|"TAUHADNU"|"TAULNUNU"|"TAUSCALARNU"|"TAUVECTORNU"|"TSS"|"TVS_PWAVE"|"VLL"|"VSP_PWAVE"|"VSS"|"VSS_BMIX"|"VUB"|"VVP"|"VVPIPI"|"VVS_PWAVE"|"SSD_CP"
+MODEL_NAME.2 : "BaryonPCR"|"BTO3PI_CP"|"BTOSLLALI"|"BTOSLLBALL"|"BTOXSGAMMA"|"BTOXSLL"|"CB3PI-MPP"|"CB3PI-P00"|"D_DALITZ"|"ETA_DALITZ"|"GOITY_ROBERTS"|"HELAMP"|"HQET"|"ISGW2"|"LbAmpGen"|"OMEGA_DALITZ"|"PARTWAVE"|"PHSP"|"PI0_DALITZ"|"PYTHIA"|"SLN"|"STS"|"SVP_HELAMP"|"SVS"|"SVV_HELAMP"|"TAUHADNU"|"TAULNUNU"|"TAUSCALARNU"|"TAUVECTORNU"|"TSS"|"TVS_PWAVE"|"VLL"|"VSP_PWAVE"|"VSS"|"VSS_BMIX"|"VUB"|"VVP"|"VVPIPI"|"VVS_PWAVE"|"SSD_CP"
 
 LABEL : /[a-zA-Z0-9\/\-+*_()']+/
 COMMENT : /[;#][^\n]*/
diff --git a/decaylanguage/dec/__init__.py b/decaylanguage/dec/__init__.py
index e69de29b..10e7de86 100644
--- a/decaylanguage/dec/__init__.py
+++ b/decaylanguage/dec/__init__.py
@@ -0,0 +1,3 @@
+from __future__ import absolute_import
+
+from .dec import DecFileParser
diff --git a/decaylanguage/dec/dec.py b/decaylanguage/dec/dec.py
index 1b72d88d..0d2861c2 100644
--- a/decaylanguage/dec/dec.py
+++ b/decaylanguage/dec/dec.py
@@ -1,8 +1,40 @@
+"""
+Submodule with classes and utilities to deal with and parse .dec decay files.
+
+Basic assumptions
+-----------------
+
+1) For standard (non-alias) particle (names):
+    - Decay modes defined via a 'Decay' statement.
+    - Related antiparticle decay modes either defined via a 'CDecay' statement
+      or via a 'Decay' statement. The latter option is often used if CP matters.
+2) For particle names defined via aliases:
+    - Particle decay modes defined as above.
+    - Related antiparticle decay modes defined with either options above,
+      *but* there needs to be a 'ChargeConj' statement specifying the
+      particle-antiparticle match. Typically::
+
+        Alias MyP+ P+
+        Alias MyP- P-
+        ChargeConj MyP+ MyP-
+        Decay MyP+
+        ...
+        Enddecay
+        CDecay MyP-
+
+3) As a consequence, particles that are self-conjugate should not be used
+   in 'CDecay' statements, obviously.
+"""
+
 from __future__ import absolute_import, division, print_function
 
 import os
 import warnings
-from lark import Lark, Transformer, Tree
+
+from lark import Lark
+from lark import Tree, Transformer, Visitor
+
+from particle import Particle, ParticleNotFound
 
 from ..data import open_text
 from .. import data
@@ -31,62 +63,103 @@ class DecFileParser(object):
 
     Example
     -------
-    >>> parser = DecFileParser('my-decay-file.dec')
-    >>> parser.parse()
+    >>> parsed_file = DecFileParser.from_file('my-decay-file.dec')
     """
 
     __slots__ = ("_grammar",
                  "_grammar_info",
-                 "_decay_file",
+                 "_dec_file_name",
                  "_parsed_dec_file",
-                 "_parsed_decays")
+                 "_parsed_decays",
+                 "_include_ccdecays")
 
-    def __init__(self, decay_file):
+    def __init__(self, filename):
         """
-        Parser constructor.
+        Default constructor.
 
         Parameters
         ----------
-        decay_file: str
-            Input .dec decay file.
+        filename: str
+            Input .dec decay file name.
         """
-        # Conversion to handle pathlib on Python < 3.6:
-        decay_file = str(decay_file)
-
-        # Check input file
-        if not os.path.exists(decay_file):
-            raise FileNotFoundError("'{0}'!".format(decay_file))
-
         self._grammar = None       # Loaded Lark grammar definition file
         self._grammar_info = None  # Name of Lark grammar definition file
 
-        self._decay_file = decay_file  # Name of the input decay file
-        self._parsed_dec_file = None   # Parsed decay file
+        # Conversion to handle pathlib on Python < 3.6:
+        self._dec_file_name = str(filename)  # Name of the input decay file
+        self._parsed_dec_file = None      # Parsed decay file
 
         self._parsed_decays = None  # Particle decays found in the decay file
 
-    def parse(self):
+        # By default, consider charge-conjugate decays when parsing
+        self._include_ccdecays = True
+
+    @classmethod
+    def from_file(cls, filename):
         """
-        Parse the given .dec decay file according to the default Lark options.
+        Parse a .dec decay file.
 
-        See method 'load_grammar' for how to explicitly set the Lark parser
-        and the parsing options.
+        Parameters
+        ----------
+        filename: str
+            Input .dec decay file name.
         """
-        # Has a fine been parsed already?
+        # Conversion to handle pathlib on Python < 3.6:
+        filename = str(filename)
+
+        # Check input file
+        if not os.path.exists(filename):
+            raise FileNotFoundError("'{0}'!".format(filename))
+
+        return cls(filename)
+
+    def parse(self, include_ccdecays=True):
+        """
+        Parse the given .dec decay file according to default Lark parser
+        and specified options, i.e.,
+            parser = 'lalr',
+            lexer = 'standard'.
+
+        See method 'load_grammar' for how to explicitly define the grammar
+        and set the Lark parsing options.
+
+        Parameters
+        ----------
+        include_ccdecays: boolean, optional, default=True
+            Choose whether or not to consider charge-conjugate decays,
+            which are specified via "CDecay <MOTHER>".
+            Make sure you understand the consequences of ignoring
+            charge conjugate decays - you won't have a complete picture!
+        """
+        # Has a file been parsed already?
         if self._parsed_decays is not None:
             warnings.warn("Input file being re-parsed ...")
 
-        decay_file = open(self._decay_file).read()
-        parser = Lark(self.grammar, parser='lalr', lexer='standard')
+        # Override the parsing settings for charge conjugate decays
+        self._include_ccdecays = include_ccdecays if include_ccdecays else False
 
-        self._parsed_dec_file = parser.parse(decay_file)
+        # Retrieve all info on the default Lark grammar and its default options,
+        # effectively loading it
+        opts = self.grammar_info()
+        extraopts = dict((k, v) for k, v in opts.items()
+                         if k not in ('lark_file', 'parser','lexer'))
 
+        # Instantiate the Lark parser according to chosen settings
+        parser = Lark(self.grammar(), parser=opts['parser'], lexer=opts['lexer'])
+
+        dec_file = open(self._dec_file_name).read()
+        self._parsed_dec_file = parser.parse(dec_file)
+
+        # At last, find all particle decays defined in the .dec decay file ...
         self._find_parsed_decays()
 
-    @property
+        # ... and create on the fly the charge conjugate decays, if requested
+        if self._include_ccdecays:
+            self._add_charge_conjugate_decays()
+
     def grammar(self):
         """
-        This accesses the internal Lark grammar definition file,
+        Access the internal Lark grammar definition file,
         loading it from the default location if needed.
 
         Returns
@@ -99,9 +172,25 @@ def grammar(self):
 
         return self._grammar
 
+    def grammar_info(self):
+        """
+        Access the internal Lark grammar definition file name and
+        parser options, loading the grammar from the default location if needed.
+
+        Returns
+        -------
+        out: dict
+            The Lark grammar definition file name and parser options.
+        """
+        if not self.grammar_loaded:
+            self.load_grammar()
+
+        return self._grammar_info
+
     def load_grammar(self, filename=None, parser='lalr', lexer='standard', **options):
         """
-        Load a Lark grammar definition file.
+        Load a Lark grammar definition file, either the default one,
+        or a user-specified one, optionally setting Lark parsing options.
 
         Parameters
         ----------
@@ -113,6 +202,7 @@ def load_grammar(self, filename=None, parser='lalr', lexer='standard', **options
             The Lark parser lexer mode to use.
         options: keyword arguments, optional
             Extra options to pass on to the parsing algorithm.
+
         See Lark's Lark class for a description of available options
         for parser, lexer and options.
         """
@@ -134,15 +224,174 @@ def grammar_loaded(self):
         """Check to see if the Lark grammar definition file is loaded."""
         return self._grammar is not None
 
+    def dict_definitions(self):
+        """
+        Return a dictionary of all definitions in the input parsed file,
+        of the form "Define <NAME> <VALUE>",
+        as {'NAME1': VALUE1, 'NAME2': VALUE2, ...}.
+        """
+        return get_definitions(self._parsed_dec_file)
+
+    def dict_aliases(self):
+        """
+        Return a dictionary of all alias definitions in the input parsed file,
+        of the form "Alias <NAME> <ALIAS>",
+        as {'NAME1': ALIAS1, 'NAME2': ALIAS2, ...}.
+        """
+        return get_aliases(self._parsed_dec_file)
+
+    def dict_charge_conjugates(self):
+        """
+        Return a dictionary of all charge conjugate definitions
+        in the input parsed file, of the form
+        "ChargeConj <PARTICLE> <CC_PARTICLE>", as
+        {'PARTICLE1': CC_PARTICLE1, 'PARTICLE2': CC_PARTICLE2, ...}.
+        """
+        return get_charge_conjugate_defs(self._parsed_dec_file)
+
+    def dict_pythia_definitions(self):
+        """
+        Return a dictionary of all Pythia definitions in the input parsed file,
+        of the form
+        "PythiaBothParam <NAME>=<LABEL>"
+        or
+        "PythiaBothParam <NAME>=<NUMBER>",
+        as {'NAME1': 'LABEL1', 'NAME2': VALUE2, ...}.
+        """
+        return get_pythia_definitions(self._parsed_dec_file)
+
+    def list_lineshape_definitions(self):
+        """
+        Return a dictionary of all SetLineshapePW definitions in the input parsed file,
+        of the form
+        "SetLineshapePW <MOTHER> <DAUGHTER1> <DAUGHTER2> <VALUE>",
+        as
+        [(['MOTHER1', 'DAUGHTER1-1', 'DAUGHTER1-2'], VALUE1),
+        (['MOTHER2', 'DAUGHTER2-1', 'DAUGHTER2-2'], VALUE2),
+        ...]
+        """
+        return get_lineshape_definitions(self._parsed_dec_file)
+
+    def global_photos_flag(self):
+        """
+        Return a boolean-like PhotosEnum enum specifying whether or not PHOTOS
+        has been enabled.
+
+        Note: PHOTOS is turned on(off) for all decays with the global flag
+        yesPhotos(noPhotos).
+
+        Returns
+        -------
+        out: PhotosEnum, default=PhotosEnum.no
+            PhotosEnum.yes / PhotosEnum.no if PHOTOS enabled / disabled
+        """
+        return get_global_photos_flag(self._parsed_dec_file)
+
+    def list_charge_conjugate_decays(self):
+        """
+        Return a (sorted) list of all charge conjugate decay definitions
+        in the input parsed file, of the form "CDecay <MOTHER>", as
+        ['MOTHER1', 'MOTHER2', ...].
+        """
+        return get_charge_conjugate_decays(self._parsed_dec_file)
+
     def _find_parsed_decays(self):
-        """Find all Tree instances of Tree.data='decay'."""
+        """
+        Return a tuple of all decay definitions in the input parsed file,
+        of the form
+        "Decay <MOTHER>",
+        as a tuple of Lark Tree instances with Tree.data=='decay', i.e.,
+        (Tree(decay, [Tree(particle, [Token(LABEL, <MOTHER1>]), ...),
+        Tree(decay, [Tree(particle, [Token(LABEL, <MOTHER2>]), ...)).
+
+        Duplicate definitions (a bug, of course) are removed, issuing a warning.
+
+        Note
+        ----
+        1) Method not meant to be used directly!
+        2) charge conjugates need to be dealt with differently,
+        see 'self._add_charge_conjugate_decays()'.
+        """
         if self._parsed_dec_file is not None:
-            self._parsed_decays = tuple(self._parsed_dec_file.find_data('decay'))
+            self._parsed_decays = get_decays(self._parsed_dec_file)
 
         # Check for duplicates - should be considered a bug in the .dec file!
         self._check_parsed_decays()
 
-        return self._parsed_decays
+    def _add_charge_conjugate_decays(self):
+        """
+        If requested (see the 'self._include_ccdecays' class attribute),
+        create the Lark Tree instances describing the charge conjugate decays
+        specified in the input parsed file via the statements of the form
+        "CDecay <MOTHER>".
+        These are added to the internal list of decays stored in the class
+        in variable 'self._parsed_decays', performing a charge conjugate (CC)
+        transformation on each CC-related decay, which is cloned.
+
+        Note
+        ----
+        1) If a decay file only defines 'Decay' decays and no 'CDecay',
+        then no charge conjugate decays will be created!
+        This seems correct given the "instructions" in the decay file:
+        - There is no 'CDecay' statement related to a 'Decay' statement
+          for a self-conjugate decaying particle such as the pi0.
+        - Else the decay file should be considered incomplete, hence buggy.
+        2) Method not meant to be used directly!
+        """
+
+        # Do not add any charge conjugate decays if the input parsed file
+        # does not define any!
+        mother_names_ccdecays = self.list_charge_conjugate_decays()
+        if len(mother_names_ccdecays) == 0:
+            return
+
+        # Cross-check - make sure charge conjugate decays are not defined
+        # with both 'Decay' and 'CDecay' statements!
+        mother_names_decays = [get_decay_mother_name(tree)
+                               for tree in self._parsed_decays]
+
+        duplicates = [n for n in mother_names_ccdecays if n in mother_names_decays]
+        if len(duplicates) > 0:
+            msg = """The following particles are defined in the input .dec file with both 'Decay' and 'CDecay': {0}!
+The 'CDecay' definition(s) will be ignored ...""".format(', '.join(d for d in duplicates))
+            warnings.warn(msg)
+
+        # If that's the case, proceed using the decay definitions specified
+        # via the 'Decay' statement, hence discard/remove the definition
+        # via the 'CDecay' statement.
+        for d in duplicates:
+            mother_names_ccdecays.remove(d)
+
+        # We're done if there are no more 'CDecay' decays to treat!
+        if len(mother_names_ccdecays) == 0:
+            return
+
+        # At last, create the charge conjugate decays:
+        # First, make a (deep) copy of the list of relevant Tree instances.
+        # Example:
+        # if mother_names_ccdecays = ['anti-M1', 'anti-M2'],
+        # the relevant Trees are the ones describing the decays of ['M1', 'M2'].
+        dict_cc_names = self.dict_charge_conjugates()
+
+        i = 0
+        trees_to_conjugate = []
+        for ccname in mother_names_ccdecays:
+            name = find_charge_conjugate_match(ccname, dict_cc_names)
+            i += 1
+            match = list(self._parsed_dec_file.find_pred(
+                lambda t: t.data=='decay' and t.children[0].children[0].value == name))
+            trees_to_conjugate.extend(match)
+
+        cdecays = [ tree.__deepcopy__(None) for tree in trees_to_conjugate]
+
+        # Finally, perform all particle -> anti(particle) replacements,
+        # taking care of charge conjugate decays defined via aliases,
+        # passing them as charge conjugates to be processed manually.
+        [ChargeConjugateReplacement(charge_conj_defs=dict_cc_names).visit(t)
+        for t in cdecays]
+
+        # ... and add all these charge-conjugate decays to the list of decays!
+        self._parsed_decays.extend(cdecays)
 
     def _check_parsing(self):
         """Has the .parse() method been called already?"""
@@ -150,9 +399,35 @@ def _check_parsing(self):
             raise DecFileNotParsed("Hint: call 'parse()'!")
 
     def _check_parsed_decays(self):
+        """
+        Is the number of decays parsed consistent with the number of
+        decay mother names? An inconsistency can arise if decays are redefined.
+
+        Duplicates are removed, starting from the second occurrence.
+        """
+        # Issue a helpful warning if duplicates are found
         lmn = self.list_decay_mother_names()
+        duplicates = []
         if self.number_of_decays != len(set(lmn)):
-            warnings.warn("Input .dec file redefines decays for particle(s) {0}!".format(set([n for n in lmn if lmn.count(n)>1])))
+            duplicates = set([n for n in lmn if lmn.count(n)>1])
+            msg = """The following particle(s) is(are) redefined in the input .dec file with 'Decay': {0}!
+All but the first occurence will be discarded/removed ...""".format(', '.join(d for d in duplicates))
+            warnings.warn(msg)
+
+        # Create a list with all occurrences to remove
+        # (duplications means multiple instances to remove)
+        duplicates_to_remove = []
+        for item in duplicates:
+            c = lmn.count(item)
+            if c>1:
+                duplicates_to_remove.extend([item]*(c-1))
+
+        # Actually remove all but the first occurence of duplicate decays
+        for tree in reversed(self._parsed_decays):
+            val = tree.children[0].children[0].value
+            if val in duplicates_to_remove:
+                duplicates_to_remove.remove(val)
+                self._parsed_decays.remove(tree)
 
     @property
     def number_of_decays(self):
@@ -162,11 +437,23 @@ def number_of_decays(self):
         return len(self._parsed_decays)
 
     def list_decay_mother_names(self):
+        """
+        Return a list of all decay mother names found in the parsed decay file.
+        """
         self._check_parsing()
 
         return [get_decay_mother_name(d) for d in self._parsed_decays]
 
     def _find_decay_modes(self, mother):
+        """
+        Return a tuple of Lark Tree instances describing all the decay modes
+        of the input mother particle as defined in the parsed .dec file.
+
+        Parameters
+        ----------
+        mother: str
+            Input mother particle name.
+        """
         self._check_parsing()
 
         for decay_Tree in self._parsed_decays:
@@ -291,16 +578,109 @@ def build_decay_chain(self, mother, stable_particles=[]):
     def __repr__(self):
         if self._parsed_dec_file is not None:
             return "<{self.__class__.__name__}: decfile='{decfile}', n_decays={n_decays}>".format(
-                self=self, decfile=self._decay_file, n_decays=self.number_of_decays)
+                self=self, decfile=self._dec_file_name, n_decays=self.number_of_decays)
         else:
             return "<{self.__class__.__name__}: decfile='{decfile}'>"\
-                    .format(self=self, decfile=self._decay_file)
+                    .format(self=self, decfile=self._dec_file_name)
 
     def __str__(self):
         return repr(self)
 
 
+class ChargeConjugateReplacement(Visitor):
+    """
+    Lark Visitor implementing the replacement of all particle names
+    with their charge conjugate particle names
+    in a Lark Tree of name 'particle' (Tree.data == 'particle').
+
+    Note
+    ----
+    1) There is no check of whether the mother particle is self-conjugate or not.
+    It is the responsibility of the caller to make sure the operation
+    is not trivial, meaning returning the same (self-conjugate) decay!
+    2) If a particle name (say, 'UNKOWN') is not found or known,
+    (search done via the Particle class in the particle package),
+    its charge conjugate name is denoted as 'ChargeConj(UNKOWN)'.
+
+    Parameters
+    ----------
+    charge_conj_defs: dict, optional, default={}
+        Dictionary with the charge conjugate particle definitions
+        in the parsed file. Argument to be passed to the class constructor.
+
+    Examples
+    --------
+    >>> from lark import Tree, Token
+    >>> ...
+    >>> t = Tree('decay', [Tree('particle', [Token('LABEL', 'D0')]), Tree('decayline', [Tree
+    ... ('value', [Token('SIGNED_NUMBER', '1.0')]), Tree('particle', [Token('LABEL', 'K-')])
+    ... , Tree('particle', [Token('LABEL', 'pi+')]), Tree('model', [Token('MODEL_NAME', 'PHS
+    ... P')])])])
+    >>> ChargeConjugateReplacement().visit(t)
+    Tree(decay, [Tree(particle, [Token(LABEL, 'D~0')]), Tree(decayline,
+    [Tree(value, [Token(SIGNED_NUMBER, '1.0')]), Tree(particle, [Token(LABEL, 'K+')]),
+    Tree(particle, [Token(LABEL, 'pi-')]), Tree(model, [Token(MODEL_NAME, 'PHSP')])])])
+    """
+    def __init__(self, charge_conj_defs=dict()):
+        self.charge_conj_defs = charge_conj_defs
+
+    # TODO: this can be improved!
+    def _last_chance_matching(self, pname):
+        try:
+            return Particle.from_dec(pname).invert().name
+        except ParticleNotFound:
+            return 'ChargeConj({0})'.format(pname)
+
+    def particle(self, tree):
+        """
+        Method for the rule (here, a replacement) we wish to implement.
+        """
+        assert tree.data == 'particle'
+        pname = tree.children[0].value
+        if len(self.charge_conj_defs) > 0:
+            for p, ccp in self.charge_conj_defs.items():
+                if ccp == pname:
+                    tree.children[0].value = p
+                # Yes, both 'ChargeConj P CCP' and 'ChargeConj CCP P' are relevant
+                elif p == pname:
+                    tree.children[0].value = ccp
+                else:
+                    tree.children[0].value = self._last_chance_matching(tree.children[0].value)
+        else:
+            tree.children[0].value = self._last_chance_matching(tree.children[0].value)
+
+
+def find_charge_conjugate_match(ccname, dict_cc_names=dict()):
+    # TODO: this can be improved! To be revisited ...
+    if len(dict_cc_names) > 0:
+        for p, ccp in dict_cc_names.items():
+            if ccp == ccname:
+                return p
+            # Yes, both 'ChargeConj P CCP' and 'ChargeConj CCP P' are relevant
+            elif p == ccname:
+                return ccp
+        else:
+            try:
+                return Particle.from_dec(ccname).invert().name
+            except ParticleNotFound:
+                return 'ChargeConj({0})'.format(ccname)
+    else:
+        try:
+            return Particle.from_dec(ccname).invert().name
+        except ParticleNotFound:
+            return 'ChargeConj({0})'.format(ccname)
+
+
 def get_decay_mother_name(decay_tree):
+    """
+    Return the mother particle name for the decay mode defined
+    in the input Tree of name 'decay'.
+
+    Parameters
+    ----------
+    decay_tree: Lark Tree instance
+        Input Tree satisfying Tree.data=='decay'.
+    """
     if not isinstance(decay_tree, Tree) or decay_tree.data != 'decay':
         raise RuntimeError("Input not an instance of a 'decay' Tree!")
 
@@ -310,6 +690,15 @@ def get_decay_mother_name(decay_tree):
 
 
 def get_branching_fraction(decay_mode):
+    """
+    Return the branching fraction (float) for the decay mode defined
+    in the input Tree of name 'decayline'.
+
+    Parameters
+    ----------
+    decay_mode: Lark Tree instance
+        Input Tree satisfying Tree.data=='decayline'.
+    """
     if not isinstance(decay_mode, Tree) or decay_mode.data != 'decayline':
         raise RuntimeError("Check your input, not an instance of a 'decayline' Tree!")
 
@@ -322,6 +711,25 @@ def get_branching_fraction(decay_mode):
 
 
 def get_final_state_particles(decay_mode):
+    """
+    Return a list of Lark Tree instances describing the final-state particles
+    for the decay mode defined in the input Tree of name 'decayline'.
+
+    Parameters
+    ----------
+    decay_mode: Lark Tree instance
+        Input Tree satisfying Tree.data=='decayline'.
+
+    Examples
+    --------
+    For a decay
+        Decay MyB_s0
+            1.000     K+     K-     SSD_CP 20.e12 0.1 1.0 0.04 9.6 -0.8 8.4 -0.6;
+        Enddecay
+    the list
+        [Tree(particle, [Token(LABEL, 'K+')]), Tree(particle, [Token(LABEL, 'K-')])]
+    will be returned.
+    """
     if not isinstance(decay_mode, Tree) or decay_mode.data != 'decayline':
         raise RuntimeError("Input not an instance of a 'decayline' Tree!")
 
@@ -330,6 +738,23 @@ def get_final_state_particles(decay_mode):
 
 
 def get_final_state_particle_names(decay_mode):
+    """
+    Return a list of final-state particle names for the decay mode defined
+    in the input Tree of name 'decayline'.
+
+    Parameters
+    ----------
+    decay_mode: Lark Tree instance
+        Input Tree satisfying Tree.data=='decayline'.
+
+    Examples
+    --------
+    For a decay
+        Decay MyB_s0
+            1.000     K+     K-     SSD_CP 20.e12 0.1 1.0 0.04 9.6 -0.8 8.4 -0.6;
+        Enddecay
+    the list ['K+', 'K-'] will be returned.
+    """
     if not isinstance(decay_mode, Tree) or decay_mode.data != 'decayline':
         raise RuntimeError("Input not an instance of a 'decayline' Tree!")
 
@@ -339,6 +764,22 @@ def get_final_state_particle_names(decay_mode):
 
 
 def get_model_name(decay_mode):
+    """
+    Return the decay model name in a Tree of name 'decayline'.
+
+    Parameters
+    ----------
+    decay_mode: Lark Tree instance
+        Input Tree satisfying Tree.data=='decayline'.
+
+    Examples
+    --------
+    For a decay
+        Decay MyB_s0
+            1.000     K+     K-     SSD_CP 20.e12 0.1 1.0 0.04 9.6 -0.8 8.4 -0.6;
+        Enddecay
+    the string 'SSD_CP' will be returned.
+    """
     if not isinstance(decay_mode, Tree) or decay_mode.data != 'decayline':
         raise RuntimeError("Input not an instance of a 'decayline' Tree!")
 
@@ -347,6 +788,25 @@ def get_model_name(decay_mode):
 
 
 def get_model_parameters(decay_mode):
+    """
+    Return a list of decay model parameters in a Tree of name 'decayline',
+    if defined, else an empty string.
+
+    Parameters
+    ----------
+    decay_mode: Lark Tree instance
+        Input Tree satisfying Tree.data=='decayline'.
+
+    Examples
+    --------
+    For a decay
+        Decay MyB_s0
+            1.000     K+     K-     SSD_CP 20.e12 0.1 1.0 0.04 9.6 -0.8 8.4 -0.6;
+        Enddecay
+    the list
+        [20000000000000.0, 0.1, 1.0, 0.04, 9.6, -0.8, 8.4, -0.6]
+    will be returned.
+    """
     if not isinstance(decay_mode, Tree) or decay_mode.data != 'decayline':
         raise RuntimeError("Input not an instance of a 'decayline' Tree!")
 
@@ -354,10 +814,53 @@ def get_model_parameters(decay_mode):
     return [float(tree.children[0].value) for tree in lmo[0].children] if len(lmo) == 1 else ''
 
 
+def get_decays(parsed_file):
+    """
+    Return a list of all decay definitions in the input parsed file,
+    of the form
+    "Decay <MOTHER>",
+    as a tuple of Lark Tree instances with Tree.data=='decay', i.e.,
+    [Tree(decay, [Tree(particle, [Token(LABEL, <MOTHER1>]), ...),
+     Tree(decay, [Tree(particle, [Token(LABEL, <MOTHER2>]), ...)].
+
+    Parameters
+    ----------
+    parsed_file: Lark Tree instance
+        Input parsed file.
+    """
+    if not isinstance(parsed_file, Tree) :
+        raise RuntimeError("Input not an instance of a Tree!")
+
+    try:
+        return list(parsed_file.find_data('decay'))
+    except:
+        RuntimeError("Input parsed file does not seem to have the expected structure.")
+
+
+def get_charge_conjugate_decays(parsed_file):
+    """
+    Return a (sorted) list of all charge conjugate decay definitions
+    in the input parsed file, of the form "CDecay <MOTHER>", as
+    ['MOTHER1', 'MOTHER2', ...].
+
+    Parameters
+    ----------
+    parsed_file: Lark Tree instance
+        Input parsed file.
+    """
+    if not isinstance(parsed_file, Tree) :
+        raise RuntimeError("Input not an instance of a Tree!")
+
+    try:
+        return sorted([tree.children[0].children[0].value for tree in parsed_file.find_data('cdecay')])
+    except:
+        RuntimeError("Input parsed file does not seem to have the expected structure.")
+
+
 def get_definitions(parsed_file):
     """
     Return a dictionary of all definitions in the input parsed file, of the form
-    "Define <NAME> <VALUE>" as {'NAME1': VALUE1, 'NAME2': VALUE2, ...}.
+    "Define <NAME> <VALUE>", as {'NAME1': VALUE1, 'NAME2': VALUE2, ...}.
 
     Parameters
     ----------
@@ -371,13 +874,13 @@ def get_definitions(parsed_file):
         return {tree.children[0].children[0].value:float(tree.children[1].children[0].value)
             for tree in parsed_file.find_data('define')}
     except:
-        warnings.error("Input parsed file does not seem to have the expected structure.")
+        RuntimeError("Input parsed file does not seem to have the expected structure.")
 
 
 def get_aliases(parsed_file):
     """
     Return a dictionary of all aliases in the input parsed file, of the form
-    "Alias <NAME> <ALIAS>" as {'NAME1': ALIAS1, 'NAME2': ALIAS2, ...}.
+    "Alias <NAME> <ALIAS>", as {'NAME1': ALIAS1, 'NAME2': ALIAS2, ...}.
 
     Parameters
     ----------
@@ -391,14 +894,14 @@ def get_aliases(parsed_file):
         return {tree.children[0].children[0].value:tree.children[1].children[0].value
             for tree in parsed_file.find_data('alias')}
     except:
-        warnings.error("Input parsed file does not seem to have the expected structure.")
+        RuntimeError("Input parsed file does not seem to have the expected structure.")
 
 
 def get_charge_conjugate_defs(parsed_file):
     """
-    Return a dictionary of all aliases in the input parsed file, of the form
-    "ChargeConj <PARTICLE> <CC_PARTICLE>" as
-    {'PARTICLE1': CC_PARTICLE1, 'PARTICLE2': CC_PARTICLE2, ...}.
+    Return a dictionary of all charge conjugate definitions
+    in the input parsed file, of the form "ChargeConj <PARTICLE> <CC_PARTICLE>",
+    as {'PARTICLE1': CC_PARTICLE1, 'PARTICLE2': CC_PARTICLE2, ...}.
 
     Parameters
     ----------
@@ -412,7 +915,104 @@ def get_charge_conjugate_defs(parsed_file):
         return {tree.children[0].children[0].value:tree.children[1].children[0].value
             for tree in parsed_file.find_data('chargeconj')}
     except:
-        warnings.error("Input parsed file does not seem to have the expected structure.")
+        RuntimeError("Input parsed file does not seem to have the expected structure.")
+
+
+def get_pythia_definitions(parsed_file):
+    """
+    Return a dictionary of all Pythia definitions in the input parsed file,
+    of the form
+    "PythiaBothParam <NAME>=<LABEL>"
+    or
+    "PythiaBothParam <NAME>=<NUMBER>",
+    as {'NAME1': 'LABEL1', 'NAME2': VALUE2, ...}.
+
+    Parameters
+    ----------
+    parsed_file: Lark Tree instance
+        Input parsed file.
+    """
+    if not isinstance(parsed_file, Tree) :
+        raise RuntimeError("Input not an instance of a Tree!")
+
+    def str_or_float(arg):
+        try:
+            return float(arg)
+        except:
+            return arg
+
+    try:
+        return {'{0}:{1}'.format(tree.children[0].value, tree.children[1].value):str_or_float(tree.children[2].value)
+            for tree in parsed_file.find_data('pythia_def')}
+    except:
+        RuntimeError("Input parsed file does not seem to have the expected structure.")
+
+
+def get_lineshape_definitions(parsed_file):
+    """
+    Return a dictionary of all SetLineshapePW definitions in the input parsed file,
+    of the form
+    "SetLineshapePW <MOTHER> <DAUGHTER1> <DAUGHTER2> <VALUE>",
+    as
+    [(['MOTHER1', 'DAUGHTER1-1', 'DAUGHTER1-2'], VALUE1),
+     (['MOTHER2', 'DAUGHTER2-1', 'DAUGHTER2-2'], VALUE2),
+     ...]
+
+    Parameters
+    ----------
+    parsed_file: Lark Tree instance
+        Input parsed file.
+    """
+    if not isinstance(parsed_file, Tree) :
+        raise RuntimeError("Input not an instance of a Tree!")
+
+    try:
+        d = list()
+        for tree in parsed_file.find_data('setlspw'):
+            particles = [p.children[0].value for p in tree.children[:-1]]
+            val = int(tree.children[3].children[0].value)
+            d.append((particles, val))
+        return d
+    except:
+        RuntimeError("Input parsed file does not seem to have the expected structure.")
+
+
+def get_global_photos_flag(parsed_file):
+    """
+    Return a boolean-like PhotosEnum enum specifying whether or not PHOTOS
+    has been enabled.
+
+    Note: PHOTOS is turned on(off) for all decays with the global flag
+    yesPhotos(noPhotos).
+
+    Parameters
+    ----------
+    parsed_file: Lark Tree instance
+        Input parsed file.
+
+    Returns
+    -------
+    out: PhotosEnum, default=PhotosEnum.no
+        PhotosEnum.yes / PhotosEnum.no if PHOTOS enabled / disabled
+    """
+    if not isinstance(parsed_file, Tree) :
+        raise RuntimeError("Input not an instance of a Tree!")
+
+    # Check if the flag is not set more than once, just in case ...
+    tree = tuple(parsed_file.find_data('global_photos'))
+    print('TREE:', tree)
+    if len(tree) == 0:
+        return PhotosEnum.no
+    elif len(tree) > 1:
+            warnings.warn("PHOTOS flag re-set! Using flag set in last ...")
+
+    tree = tree[-1]
+
+    try:
+        val = tree.children[0].data
+        return PhotosEnum.yes if val=='yes' else PhotosEnum.no
+    except:
+        RuntimeError("Input parsed file does not seem to have the expected structure.")
 
 
 class TreeToDec(Transformer):
diff --git a/tests/data/defs-aliases-chargeconj.dec b/tests/data/defs-aliases-chargeconj.dec
new file mode 100644
index 00000000..c3bacb3d
--- /dev/null
+++ b/tests/data/defs-aliases-chargeconj.dec
@@ -0,0 +1,271 @@
+# File for testing purposes.
+# Contains definitions and aliases extracted from a
+# typical DECAY.DEC file from LHCb.
+
+Define qoverp_incohMix_B_s0 1.0
+Define dm_incohMix_B_s0 17.8e12
+Define qoverp_incohMix_B0 1.0
+Define dm_incohMix_B0 0.507e12
+# Old definition of dm still in some decay models
+Define dm 0.507e12
+# Activate incoherent Mixing
+# define the values of the CKM angles  (alpha=70, beta=40)
+Define alpha 1.365
+Define beta  0.39
+Define gamma 1.387
+Define twoBetaPlusGamma  2.167
+Define betaPlusHalfGamma 1.0835
+Define minusGamma -1.387
+Define minusTwoBeta -0.78
+
+# New definitions for psiKstar modes (Lange, July 26, 2000)
+Define PKHplus  0.159
+Define PKHzero  0.775
+Define PKHminus 0.612
+Define PKphHplus  1.563
+Define PKphHzero  0.0
+Define PKphHminus 2.712
+
+Define Aplus  0.490
+Define Azero  1.10
+Define Aminus 0.4
+
+Define phAplus  2.5
+Define phAzero  0.0
+Define phAminus -0.17
+
+# These particle aliases are used in for CP violating decays
+# in which the decay distributions depends on how the K* decayed.
+#
+Alias K*L       K*0
+Alias K*S       K*0
+Alias K*BL       anti-K*0
+Alias K*BS       anti-K*0
+Alias K*0T      K*0
+Alias anti-K*0T anti-K*0
+Alias K*BR  anti-K*0
+Alias K*0R  K*0
+Alias anti-K_0*0N anti-K_0*0
+Alias K_0*0N K_0*0
+
+ChargeConj K*L    K*BL
+ChargeConj K*S    K*BS
+ChargeConj K*0T    anti-K*0T
+ChargeConj K_0*0N  anti-K_0*0N
+ChargeConj K*0R    K*BR
+#
+#  Aliases for signal productions in LHCb
+#
+Alias      B0sig             B0
+Alias      anti-B0sig        anti-B0
+ChargeConj B0sig             anti-B0sig
+Alias      B+sig             B+
+Alias      B-sig             B-
+ChargeConj B+sig             B-sig
+Alias      B_s0sig           B_s0
+Alias      anti-B_s0sig      anti-B_s0
+ChargeConj B_s0sig           anti-B_s0sig
+Alias      B_c+sig           B_c+
+Alias      B_c-sig           B_c-
+ChargeConj B_c+sig           B_c-sig
+Alias      eta_bsig          eta_b
+ChargeConj eta_bsig          eta_bsig
+Alias      h_bsig            h_b
+ChargeConj h_bsig            h_bsig
+Alias      Sigma_b-sig       Sigma_b-
+Alias      anti-Sigma_b+sig  anti-Sigma_b+
+ChargeConj Sigma_b-sig       anti-Sigma_b+sig
+Alias      Lambda_b0sig      Lambda_b0
+Alias      anti-Lambda_b0sig anti-Lambda_b0
+ChargeConj Lambda_b0sig      anti-Lambda_b0sig
+Alias      Omega_b-sig       Omega_b-
+Alias      anti-Omega_b+sig  anti-Omega_b+
+ChargeConj Omega_b-sig       anti-Omega_b+sig
+Alias      Xi_b-sig          Xi_b-
+Alias      anti-Xi_b+sig     anti-Xi_b+
+ChargeConj Xi_b-sig          anti-Xi_b+sig
+Alias      Xi_b0sig          Xi_b0
+Alias      anti-Xi_b0sig     anti-Xi_b0
+ChargeConj Xi_b0sig          anti-Xi_b0sig
+Alias      eta_csig          eta_c
+ChargeConj eta_csig          eta_csig
+Alias      J/psisig          J/psi
+ChargeConj J/psisig          J/psisig
+Alias      chi_c0sig         chi_c0
+ChargeConj chi_c0sig         chi_c0sig
+Alias      chi_c1sig         chi_c1
+ChargeConj chi_c1sig         chi_c1sig
+Alias      chi_c2sig         chi_c2
+ChargeConj chi_c2sig         chi_c2sig
+Alias      psi(2S)sig        psi(2S)
+ChargeConj psi(2S)sig        psi(2S)sig
+Alias      psi(3770)sig      psi(3770)
+ChargeConj psi(3770)sig      psi(3770)sig
+Alias      D*+sig            D*+
+Alias      D*-sig            D*-
+ChargeConj D*+sig            D*-sig
+Alias      D*0sig            D*0
+Alias      anti-D*0sig       anti-D*0
+ChargeConj D*0sig            anti-D*0sig
+Alias      D0sig             D0
+Alias      anti-D0sig        anti-D0
+ChargeConj D0sig             anti-D0sig
+Alias      D+sig             D+
+Alias      D-sig             D-
+ChargeConj D+sig             D-sig
+Alias      D_s+sig           D_s+
+Alias      D_s-sig           D_s-
+ChargeConj D_s+sig           D_s-sig
+Alias      Lambda_c+sig      Lambda_c+
+Alias      anti-Lambda_c-sig anti-Lambda_c-
+ChargeConj Lambda_c+sig      anti-Lambda_c-sig
+Alias      tau+sig           tau+
+Alias      tau-sig           tau-
+ChargeConj tau+sig           tau-sig
+Alias      Upsilonsig        Upsilon
+ChargeConj Upsilonsig        Upsilonsig
+Alias      Upsilon(2S)sig    Upsilon(2S)
+ChargeConj Upsilon(2S)sig    Upsilon(2S)sig
+Alias      Upsilon(3S)sig    Upsilon(3S)
+ChargeConj Upsilon(3S)sig    Upsilon(3S)sig
+Alias      Upsilon(4S)sig    Upsilon(4S)
+ChargeConj Upsilon(4S)sig    Upsilon(4S)sig
+Alias      Upsilon(5S)sig    Upsilon(5S)
+ChargeConj Upsilon(5S)sig    Upsilon(5S)sig
+Alias      X_1(3872)sig      X_1(3872)
+ChargeConj X_1(3872)sig      X_1(3872)sig
+Alias      h_csig            h_c
+ChargeConj h_csig            h_csig
+Alias      Sigma+sig         Sigma+
+Alias      anti-Sigma-sig    anti-Sigma-
+ChargeConj Sigma+sig         anti-Sigma-sig
+Alias      K+sig             K+
+Alias      K-sig             K-
+ChargeConj K+sig             K-sig
+Alias      Lambda0sig        Lambda0
+Alias      anti-Lambda0sig   anti-Lambda0
+ChargeConj Lambda0sig        anti-Lambda0sig
+Alias      B_10sig           B_10
+Alias      anti-B_10sig      anti-B_10
+ChargeConj B_10sig           anti-B_10sig
+Alias      B_2*0sig          B_2*0
+Alias      anti-B_2*0sig     anti-B_2*0
+ChargeConj B_2*0sig          anti-B_2*0sig
+Alias      B_s10sig          B_s10
+Alias      anti-B_s10sig     anti-B_s10
+ChargeConj B_s10sig          anti-B_s10sig
+Alias      B_s2*0sig         B_s2*0
+Alias      anti-B_s2*0sig    anti-B_s2*0
+ChargeConj B_s2*0sig         anti-B_s2*0sig
+Alias      B_1+sig           B_1+
+Alias      B_1-sig           B_1-
+ChargeConj B_1+sig           B_1-sig
+Alias      B_2*+sig          B_2*+
+Alias      B_2*-sig          B_2*-
+ChargeConj B_2*+sig          B_2*-sig
+Alias      K_S0sig           K_S0
+ChargeConj K_S0sig           K_S0sig
+Alias      K_L0sig           K_L0
+ChargeConj K_L0sig           K_L0sig
+Alias      phisig            phi
+ChargeConj phisig            phisig
+Alias	   Xi_cc+sig        Xi_cc+
+Alias	   anti-Xi_cc-sig   anti-Xi_cc-
+ChargeConj Xi_cc+sig		anti-Xi_cc-sig
+Alias      D_s0*+sig           D_s0*+
+Alias      D_s0*-sig           D_s0*-
+ChargeConj D_s0*+sig           D_s0*-sig
+Alias      D'_s1+sig           D'_s1+
+Alias      D'_s1-sig           D'_s1-
+ChargeConj D'_s1+sig           D'_s1-sig
+Alias      D_s1+sig           D_s1+
+Alias      D_s1-sig           D_s1-
+ChargeConj D_s1+sig           D_s1-sig
+Alias      chi_b0sig         chi_b0
+ChargeConj chi_b0sig         chi_b0sig
+Alias      chi_b1sig         chi_b1
+ChargeConj chi_b1sig         chi_b1sig
+Alias      chi_b2sig         chi_b2
+ChargeConj chi_b2sig         chi_b2sig
+Alias      Sigma_c0sig          Sigma_c0
+Alias      anti-Sigma_c0sig     anti-Sigma_c0
+ChargeConj Sigma_c0sig          anti-Sigma_c0sig
+Alias      Sigma_c+sig          Sigma_c+
+Alias      anti-Sigma_c-sig     anti-Sigma_c-
+ChargeConj Sigma_c+sig          anti-Sigma_c-sig
+Alias      Sigma_c++sig         Sigma_c++
+Alias      anti-Sigma_c--sig    anti-Sigma_c--
+ChargeConj Sigma_c++sig         anti-Sigma_c--sig
+Alias      Sigma_c*0sig         Sigma_c*0
+Alias      anti-Sigma_c*0sig    anti-Sigma_c*0
+ChargeConj Sigma_c*0sig         anti-Sigma_c*0sig
+Alias      Sigma_c*+sig         Sigma_c*+
+Alias      anti-Sigma_c*-sig    anti-Sigma_c*-
+ChargeConj Sigma_c*+sig         anti-Sigma_c*-sig
+Alias      Sigma_c*++sig        Sigma_c*++
+Alias      anti-Sigma_c*--sig   anti-Sigma_c*--
+ChargeConj Sigma_c*++sig        anti-Sigma_c*--sig
+Alias      Xi_c0sig             Xi_c0
+Alias      anti-Xi_c0sig        anti-Xi_c0
+ChargeConj Xi_c0sig             anti-Xi_c0sig
+Alias      Xi_c+sig             Xi_c+
+Alias      anti-Xi_c-sig        anti-Xi_c-
+ChargeConj Xi_c+sig             anti-Xi_c-sig
+Alias      D_2*0sig             D_2*0
+Alias      anti-D_2*0sig        anti-D_2*0
+ChargeConj D_2*0sig             anti-D_2*0sig
+Alias      D_2*-sig             D_2*-
+Alias      D_2*+sig             D_2*+
+ChargeConj D_2*+sig             D_2*-sig
+Alias      D_s2*-sig            D_s2*-
+Alias      D_s2*+sig            D_s2*+
+ChargeConj D_s2*+sig            D_s2*-sig
+Alias      D_s*-sig             D_s*-
+Alias      D_s*+sig             D_s*+
+ChargeConj D_s*+sig             D_s*-sig
+Alias      Omega-sig            Omega-
+Alias      anti-Omega+sig       anti-Omega+
+ChargeConj Omega-sig            anti-Omega+sig
+Alias      Sigma_b0sig          Sigma_b0
+Alias      anti-Sigma_b0sig     anti-Sigma_b0
+ChargeConj Sigma_b0sig          anti-Sigma_b0sig
+Alias      Sigma_b+sig          Sigma_b+
+Alias      anti-Sigma_b-sig     anti-Sigma_b-
+ChargeConj Sigma_b+sig          anti-Sigma_b-sig
+Alias      Xi_bc0sig            Xi_bc0
+Alias      anti-Xi_bc0sig       anti-Xi_bc0
+ChargeConj Xi_bc0sig            anti-Xi_bc0sig
+Alias      Xi_bc+sig            Xi_bc+
+Alias      anti-Xi_bc-sig       anti-Xi_bc-
+ChargeConj Xi_bc+sig            anti-Xi_bc-sig
+Alias	     Xi_cc++sig           Xi_cc++
+Alias	     anti-Xi_cc--sig      anti-Xi_cc--
+ChargeConj Xi_cc++sig           anti-Xi_cc--sig
+Alias      Xi-sig               Xi-
+Alias      anti-Xi+sig          anti-Xi+
+ChargeConj Xi-sig               anti-Xi+sig
+Alias      Xi0sig               Xi0
+Alias      anti-Xi0sig          anti-Xi0
+ChargeConj Xi0sig               anti-Xi0sig
+Alias      Omega_cc+sig         Omega_cc+
+Alias      anti-Omega_cc-sig    anti-Omega_cc-
+ChargeConj Omega_cc+sig         anti-Omega_cc-sig
+
+#(Very important that there are no blank spaces in the parameter string!)
+#Turn off B0-B0B mixing in Pythia
+PythiaBothParam ParticleDecays:mixB=off
+PythiaBothParam Init:showChangedSettings=off
+PythiaBothParam Init:showChangedParticleData=off
+PythiaBothParam Next:numberShowEvent=0
+
+# Turn on PHOTOS for all decays
+yesPhotos
+
+SetLineshapePW D_1+ D*+ pi0 2
+SetLineshapePW D_1+ D*0 pi+ 2
+SetLineshapePW D_1- D*- pi0 2
+SetLineshapePW D_1- anti-D*0 pi- 2
+SetLineshapePW D_10 D*0 pi0 2
+SetLineshapePW D_10 D*+ pi- 2
+SetLineshapePW anti-D_10 anti-D*0 pi0 2
+SetLineshapePW anti-D_10 D*- pi+ 2
diff --git a/tests/data/duplicate-decays.dec b/tests/data/duplicate-decays.dec
new file mode 100644
index 00000000..4142a310
--- /dev/null
+++ b/tests/data/duplicate-decays.dec
@@ -0,0 +1,40 @@
+# File for testing purposes.
+# Example bug taken from an old DECAY.DEC file ;-) :
+# a particle decay is duplicated and the CP conjugate definition with
+# 'CDecay anti-Sigma(1775)0' is effectively a different way of defining what is
+# already defined via 'Decay anti-Sigma(1775)0' !
+
+Decay Sigma(1775)0
+0.215  p+             K-              PHSP;
+0.215  n0             anti-K0         PHSP;
+0.20   Lambda0        pi0             PHSP;
+0.02   Sigma+         pi-             PHSP;
+0.02   Sigma-         pi+             PHSP;
+0.055  Sigma*+        pi-             PHSP;
+0.055  Sigma*-        pi+             PHSP;
+0.22   Lambda(1520)0  pi0             PHSP;
+Enddecay
+
+Decay anti-Sigma(1775)0
+0.215  anti-p-             K+         PHSP;
+0.215  anti-n0             K0         PHSP;
+0.20   anti-Lambda0        pi0        PHSP;
+0.02   anti-Sigma+         pi-        PHSP;
+0.02   anti-Sigma-         pi+        PHSP;
+0.055  anti-Sigma*+        pi-        PHSP;
+0.055  anti-Sigma*-        pi+        PHSP;
+0.22   anti-Lambda(1520)0  pi0        PHSP;
+Enddecay
+
+Decay Sigma(1775)0    # PDG 3216
+  0.23   Lambda0 pi0 PHSP;
+  0.23   Lambda(1405)0 pi0 PHSP;
+  0.2    p+       K-        PHSP;
+  0.1    n0       K_S0      PHSP;
+  0.1    n0       K_L0      PHSP;
+  0.05   Sigma*+  pi-       PHSP;
+  0.05   Sigma*-  pi+       PHSP;
+  0.02   Sigma+   pi-       PHSP;
+  0.02   Sigma-   pi+       PHSP;
+Enddecay
+CDecay anti-Sigma(1775)0
diff --git a/tests/data/test_example_Dst.dec b/tests/data/test_example_Dst.dec
index d7993c2e..938c4f20 100644
--- a/tests/data/test_example_Dst.dec
+++ b/tests/data/test_example_Dst.dec
@@ -1,13 +1,6 @@
 # Example decay chain for testing purposes
-
-Alias      D0sig             D0
-Alias      anti-D0sig        anti-D0
-ChargeConj D0sig             anti-D0sig
-
-Alias      D*+sig            D*+
-Alias      D*-sig            D*-
-ChargeConj D*+sig            D*-sig
-
+# Considered by itself, this file in in fact incomplete,
+# as there are no instructions on how to decay the anti-D0 and the D-!
 
 Decay D*+
 0.6770    D0  pi+                    VSS;
diff --git a/tests/test_dec.py b/tests/test_dec.py
index 9aa026a5..20b04b58 100644
--- a/tests/test_dec.py
+++ b/tests/test_dec.py
@@ -5,9 +5,15 @@
 except ImportError:
     from pathlib import Path
 
+from lark import Tree, Token
+
 from decaylanguage import data
+
 from decaylanguage.dec.dec import DecFileParser
 from decaylanguage.dec.dec import DecFileNotParsed, DecayNotFound
+from decaylanguage.dec.dec import ChargeConjugateReplacement
+from decaylanguage.dec.dec import get_decay_mother_name
+from decaylanguage.dec.dec import get_final_state_particle_names
 
 # New in Python 3
 try:
@@ -19,14 +25,24 @@
 DIR = Path(__file__).parent.resolve()
 
 
+def test_default_constructor():
+    p = DecFileParser(DIR / 'data/test_example_Dst.dec')
+    assert p is not None
+
+
+def test_from_file():
+    p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
+    assert p is not None
+
+
 def test_unknown_decfile():
     with pytest.raises(FileNotFoundError):
-        p = DecFileParser('non-existent.dec')
+        p = DecFileParser.from_file('non-existent.dec')
 
 
 def test_non_parsed_decfile():
     with pytest.raises(DecFileNotParsed):
-        p = DecFileParser(DIR / 'data/test_example_Dst.dec')
+        p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
         p.list_decay_mother_names()
 
 
@@ -38,20 +54,19 @@ def test_non_existent_decay():
 
 
 def test_default_grammar_loading():
-    p = DecFileParser(DIR / 'data/test_example_Dst.dec')
-
+    p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
     assert p.grammar is not None
 
 
 def test_explicit_grammar_loading():
-    p = DecFileParser(DIR / 'data/test_example_Dst.dec')
+    p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
     p.load_grammar(DIR / '../decaylanguage/data/decfile.lark')
 
     assert p.grammar_loaded is True
 
 
 def test_string_representation():
-        p = DecFileParser(DIR / 'data/test_example_Dst.dec')
+        p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
 
         assert "n_decays" not in p.__str__()
 
@@ -59,8 +74,72 @@ def test_string_representation():
         assert "n_decays=5" in p.__str__()
 
 
+def test_definitions_parsing():
+    p = DecFileParser.from_file(DIR / 'data/defs-aliases-chargeconj.dec')
+    p.parse()
+
+    assert len(p.dict_definitions()) == 24
+
+
+def test_aliases_parsing():
+    p = DecFileParser.from_file(DIR / 'data/defs-aliases-chargeconj.dec')
+    p.parse()
+
+    assert len(p.dict_aliases()) == 132
+
+
+def test_charge_conjugates_parsing():
+    p = DecFileParser.from_file(DIR / 'data/defs-aliases-chargeconj.dec')
+    p.parse()
+
+    assert len(p.dict_charge_conjugates()) == 77
+
+
+def test_pythia_definitions_parsing():
+    p = DecFileParser.from_file(DIR / 'data/defs-aliases-chargeconj.dec')
+    p.parse()
+
+    assert p.dict_pythia_definitions() == {'ParticleDecays:mixB': 'off',
+                                           'Init:showChangedSettings': 'off',
+                                           'Init:showChangedParticleData': 'off',
+                                           'Next:numberShowEvent': 0.0}
+
+
+def test_list_lineshape_definitions():
+    p = DecFileParser.from_file(DIR / 'data/defs-aliases-chargeconj.dec')
+    p.parse()
+
+    assert p.list_lineshape_definitions() == [(['D_1+', 'D*+', 'pi0'], 2),
+                                              (['D_1+', 'D*0', 'pi+'], 2),
+                                              (['D_1-', 'D*-', 'pi0'], 2),
+                                              (['D_1-', 'anti-D*0', 'pi-'], 2),
+                                              (['D_10', 'D*0', 'pi0'], 2),
+                                              (['D_10', 'D*+', 'pi-'], 2),
+                                              (['anti-D_10', 'anti-D*0', 'pi0'], 2),
+                                              (['anti-D_10', 'D*-', 'pi+'], 2)]
+
+
+def test_global_photos_flag():
+    p = DecFileParser.from_file(DIR / 'data/defs-aliases-chargeconj.dec')
+    p.parse()
+
+    assert p.global_photos_flag() == True
+
+def test_missing_global_photos_flag():
+    p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
+    p.parse()
+
+    assert p.global_photos_flag() == False
+
+def test_list_charge_conjugate_decays():
+    p = DecFileParser.from_file(DIR / 'data/test_Bd2DmTauNu_Dm23PiPi0_Tau2MuNu.dec')
+    p.parse()
+
+    assert p.list_charge_conjugate_decays() == ['MyD+', 'MyTau+', 'Mya_1-', 'anti-B0sig']
+
+
 def test_simple_dec():
-    p = DecFileParser(DIR / 'data/test_example_Dst.dec')
+    p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
     p.parse()
 
     assert p.list_decay_mother_names() == ['D*+', 'D*-', 'D0', 'D+', 'pi0']
@@ -69,7 +148,7 @@ def test_simple_dec():
 
 
 def test_decay_mode_details():
-    p = DecFileParser(DIR / 'data/test_example_Dst.dec')
+    p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
     p.parse()
 
     tree_Dp = p._find_decay_modes('D+')[0]
@@ -77,9 +156,89 @@ def test_decay_mode_details():
     assert p._decay_mode_details(tree_Dp) == output
 
 
+def test_duplicate_decay_definitions():
+    p = DecFileParser(DIR / 'data/duplicate-decays.dec')
+    p.parse()
+
+    assert p.number_of_decays == 2
+
+    assert p.list_decay_mother_names() == ['Sigma(1775)0', 'anti-Sigma(1775)0']
+
+
 def test_build_decay_chain():
-    p = DecFileParser(DIR / 'data/test_example_Dst.dec')
+    p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
     p.parse()
 
     output = {'D+': [{'bf': 1.0, 'fs': ['K-', 'pi+', 'pi+', 'pi0'], 'm': 'PHSP', 'mp': ''}]}
     assert p.build_decay_chain('D+', stable_particles=['pi0']) == output
+
+
+def test_Lark_ChargeConjugateReplacement_Visitor():
+    """
+    A simple example usage of the ChargeConjugateReplacement implementation
+    of a Lark's Visitor, here replacing all particles in a 'decay' Tree
+    by their antiparticles.
+    """
+    t = Tree('decay', [Tree('particle', [Token('LABEL', 'D0')]),
+            Tree('decayline', [Tree('value', [Token('SIGNED_NUMBER', '1.0')]),
+            Tree('particle', [Token('LABEL', 'K-')]),
+            Tree('particle', [Token('LABEL', 'pi+')]),
+            Tree('model', [Token('MODEL_NAME', 'PHSP')])])])
+
+    ChargeConjugateReplacement().visit(t)
+
+    assert get_decay_mother_name(t) == 'D~0'
+    assert get_final_state_particle_names(t.children[1]) == ['K+', 'pi-']
+
+
+def test_Lark_ChargeConjugateReplacement_Visitor_with_aliases():
+    """
+    Example with a D0 decay specified via an alias (MyD0).
+    As such, it is necessary to state what the particle-antiparticle match is,
+    which in decay files would mean the following lines:
+       Alias       MyD0        D0
+       Alias       MyAnti-D0   anti-D0
+       ChargeConj  MyD0        MyAnti-D0
+    A dictionary of matches should be passed to the Lark Visitor instance.
+    """
+    t = Tree('decay', [Tree('particle', [Token('LABEL', 'MyD0')]),
+            Tree('decayline', [Tree('value', [Token('SIGNED_NUMBER', '1.0')]),
+            Tree('particle', [Token('LABEL', 'K-')]),
+            Tree('particle', [Token('LABEL', 'pi+')]),
+            Tree('model', [Token('MODEL_NAME', 'PHSP')])])])
+
+    dict_ChargeConj_defs = {'MyD0': 'MyAnti-D0'}
+
+    ChargeConjugateReplacement(charge_conj_defs=dict_ChargeConj_defs).visit(t)
+
+    assert get_decay_mother_name(t) == 'MyAnti-D0'
+    assert get_final_state_particle_names(t.children[1]) == ['K+', 'pi-']
+
+
+def test_creation_charge_conjugate_decays_in_decfile_with_aliases():
+    """
+    Decay file contains 5 particle decays defined via a 'Decay' statement
+    and the 5 charge-conjugate decays defined via a 'CDecay' statement.
+    The decay modes for the latter 5 should be created on the fly,
+    hence providing in total 10 sets of particle decays parsed.
+    """
+    p = DecFileParser.from_file(DIR / 'data/test_Bd2DstDst.dec')
+    p.parse()
+
+    assert p.number_of_decays == 10
+
+
+def test_creation_charge_conjugate_decays_in_decfile_without_CDecay_defs():
+    """
+    Decay file contains 5 particle decays defined via a 'Decay' statement,
+    but no related charge-conjugate (CC) decays defined via a 'CDecay' statement
+    since the CC particle decays are also defined via a 'Decay' statement
+    for all cases except self-conjugate (mother) particles, obviously ;-).
+    This being said, this decay file is in fact incomplete by itself,
+    as there are no instructions on how to decay the anti-D0 and the D-!
+    In short, there should only be 5 sets of decay modes parsed.
+    """
+    p = DecFileParser.from_file(DIR / 'data/test_example_Dst.dec')
+    p.parse()
+
+    assert p.number_of_decays == 5