import pkg_resources
if pkg_resources.get_distribution('CGRtools').version.split('.')[:2] != ['4', '0']:
    print('WARNING. Tutorial was tested on 4.0 version of CGRtools')
else:
    print('Welcome!')

# load data for tutorial
from pickle import load
from traceback import format_exc

with open('reactions.dat', 'rb') as f:
    reactions = load(f) # list of ReactionContainer objects

r1 = reactions[0] # reaction
m6 = r1.reactants[1]
m6copy = m6.copy()
m6copy.atom(5)._Core__isotope = 13

"""
Explanation: 6. Reactor

(c) 2019, 2020 Dr. Ramil Nugmanov;
(c) 2019 Dr. Timur Madzhidov; Ravil Mukhametgaleev

Installation instructions of CGRtools package information and tutorial's files see on https://github.com/stsouko/CGRtools
NOTE: Tutorial should be performed sequentially from the start. Random cell running will lead to unexpected results.
End of explanation
"""


from CGRtools import CGRReactor, Reactor   # import of Reactor
from CGRtools.containers import * # import of required objects
from CGRtools.containers.bonds import DynamicBond

"""
Explanation: Reactor objects stores single transformation and can apply it to molecules or CGRs.
Transformations is ReactionContainer object which in reactant side consist of query for matching group and in product side patch for updating matched atoms and bonds
End of explanation
"""


acid = QueryContainer()                       # this query matches acids. Use construction possibilities.
acid.add_atom('C', neighbors=3)   # add carboxyl carbon. Hybridization is irrelevant here
acid.add_atom('O', neighbors=1)   # add hydroxyl oxygen. Hybridization is irrelevant here 
acid.add_atom('O')                                # add carbonyl oxygen. Number of neighbors is irrelevant here.
acid.add_bond(1, 2, 1) # create single bond between carbon and hydroxyl oxygen
acid.add_bond(1, 3, 2) # create double bond
print(acid)
acid.clean2d()
acid

methyl_ester = QueryContainer()  # create patch - how carboxyl group should be changed. We write methylated group
methyl_ester.add_atom('C', 1) # second argument is predefined atom mapping. Notice that mapping corresponds...  
methyl_ester.add_atom('O', 2) # ... to order in already created acid group. Atom 2 is released water.
methyl_ester.add_atom('O', 4)
methyl_ester.add_atom('O', 3)
methyl_ester.add_atom('C', 5)
methyl_ester.add_bond(1, 4, 1)
methyl_ester.add_bond(1, 3, 2)
methyl_ester.add_bond(4, 5, 1)
# No bond between atom 1 and atom 2. This bond will be broken. 
methyl_ester.clean2d()
methyl_ester

m6 # acid

template = ReactionContainer([acid], [methyl_ester]) # merge query and patch in template, which is ReactionContainer
reactor = CGRReactor(template)                        # Reactor is initialized
reacted_acid = next(reactor(m6))                            # application of Reactor to molecule

reacted_acid.clean2d() # calculate coordinates
reacted_acid       # desired methylated ester have been generated

"""
Explanation: 6.1. Products generation
Reactor works similar to ChemAxon Reactions enumeration.
Example here presents application of it to create esters from acids.
First we need to construct carboxy group matcher query. Then, ether group need to be specified. 
Atom numbers in query and patch should be mapped to each other. The same atoms should have same numbers.
End of explanation
"""


second_stage = next(reactor(reacted_acid)) # apply transformation on product of previous transformation
second_stage.clean2d() #  recalculate coordinates for correct drawing
second_stage

"""
Explanation: One can notice presence of separate oxygen (water) and ester group.
The second group can substituted by calling reactor on observed product.
End of explanation
"""


products = second_stage.split() # split product into individual molecules
react = ReactionContainer([m6], products) # unite reagent and product into reaction. 
react

"""
Explanation: second_stage has 3 components in a single MoleculeContainer object. We can split it into individual molecules and place all molecules into ReactionContainer object. Since in CGRtools atom-to-atom mapping corresponds to numbering of atoms in molecules, the resulting product has AAM according to the rule applied. Thus, reaction has correct AAM and nothing special should be made to keep or find it.
End of explanation
"""


m6copy

enums = set()              # the set enums is used to select structurally diverse products
for m in reactor(m6copy): # limit=0 is enumeration of all possible products by reactor
    print(m)                       # print signatures for observed molecules. Notice presence of water as component of product
    m.clean2d()         # recalculate coordinates
    enums.update(m.split())        # split product into separate molecules
enums = list(enums)                # set of all resulting molecules

"""
Explanation: For multicomponent reactions one can merge molecules of reactants into single MoleculeContainer object and apply reactor on it.
It is possible to generate all available products in case that molecule has several groups matching the query.
End of explanation
"""


enums[0]

enums[1]

enums[2]

"""
Explanation: Let's have a look at molecules in set.
Note to lost isotope mark.
End of explanation
"""


reactions[1] # reaction under study

cgr = ~reactions[1] # generate reaction CGR
print(cgr)
cgr.clean2d()
cgr

cgr.centers_list # reaction has two reaction centers. [10,11,12] - pseudo reaction appeared due to AAM error

"""
Explanation: 6.2. MetaReactions (reactions on CGRs).
Reactor could be applied to CGR to introduce changes into reaction. 
6.2.1. Example of atom-to-atom mapping fixing.
End of explanation
"""


nitro = QueryCGRContainer() # construct query for invalid reaction center - CGR of wrongly mapped nitro-group
nitro.add_atom('N', charge=1, p_charge=1) # atom 1
nitro.add_atom('O', charge=0, p_charge=-1) # atom 2. Notice that due to AAM error charge was changed
nitro.add_atom('O', charge=-1, p_charge=0) # atom 3. Notice that due to AAM error charge was changed
nitro.add_atom('C') # atom 4

nitro.add_bond(1, 2, DynamicBond(2, 1)) # bond between atoms 1 and 2. Due to AAM error bond is dynamic ('2>1' type) 
nitro.add_bond(1, 3, DynamicBond(1, 2)) # bond between atoms 1 and 3. Due to AAM error bond is dynamic ('1>2' type) 
nitro.add_bond(1, 4, 1) # ordinary bond
print(nitro)
# this query matches reaction center in CGR appeared due to AAM error.
nitro.clean2d()
nitro

nitro < cgr # query matches CGR of reaction with error.

valid_nitro = QueryCGRContainer() # construct nitro group without dynamic atoms. Notice that atom order should correspond object nitro
valid_nitro.add_atom('N', charge=1, p_charge=1) # ordinary N atom
valid_nitro.add_atom('O', charge=-1, p_charge=-1) # ordinary negatively charged oxygen atom
valid_nitro.add_atom('O')                            # ordinary oxygen atom

valid_nitro.add_bond(1, 2, 1) # ordinary single bond
valid_nitro.add_bond(1, 3, 2) # ordinary double bond
print(valid_nitro)
valid_nitro.clean2d()
valid_nitro

"""
Explanation: Reaction has AAM error in nitro-group
Lets try to use Reactor for AAM fixing
End of explanation
"""


template = ReactionContainer([nitro], [valid_nitro]) # template shows how wrong part of CGR is transformed into correct one.
print(template) # notice complex structure of query: CGR signature is given in braces, then >> and molecule signature
template

"""
Explanation: Now time to prepare and apply Template to CGR based on reaction with incorrect AAM.
Template is Reaction container with query in reactants and patch in products
End of explanation
"""


reactor = CGRReactor(template)

"""
Explanation: Reactor class accept single template. Existence of dynamic bond in it is not a problem.
End of explanation
"""


fixed = next(reactor(cgr)) # fix CGR

"""
Explanation: Reactor object is callable and accept as argument molecule or CGR.
NOTE: fixed is new CGR object
End of explanation
"""


print(fixed)
fixed

"""
Explanation: CGRreactor returns None if template could not be applied, otherwise patched structure is returned.
End of explanation
"""


fixed.centers_list # reaction center appeared due to AAM error before does not exist. Only 1 reaction center is found

"""
Explanation: One can see that nitro group has no dynamic bonds any more. CGR corresponds only to substitution.
End of explanation
"""


from CGRtools.files import MRVRead
from io import StringIO

e2 = next(MRVRead('e2.mrv')) # read E2 reaction from ChemAxon MRV file
e2

# create CGR query for E2 reaction side
e2query = QueryCGRContainer() 
e2query.add_atom('C', 1) # create carbon with mapping number 1
e2query.add_atom('C', 2) # create carbon with mapping number 2
# addition of iodine atom
e2query.add_atom('I', 3, neighbors=1, p_neighbors=0, charge=0, p_charge=-1)
# addition of OH- or RO- groups
e2query.add_atom('O', 4, neighbors=[0, 1], p_neighbors=[0, 1], charge=-1, p_charge=0)

e2query.add_bond(1, 2, DynamicBond(1, 2)) # bond between two carbon corresponds to formation of double from single
e2query.add_bond(1, 3, DynamicBond(1)) # bond between carbon and halogen breaks in E2 reaction
print(e2query) # it is CGR of E2 reaction center
e2query.clean2d()
e2query

e2_cgr = ~e2 # compose reaction into CGR
e2_cgr

e2query < e2_cgr # E2 CGR pattern works!

# create patch creating SN2 reaction. Notice that ordering of atoms correspond to that of E2 CGR query
sn2patch = QueryCGRContainer()
sn2patch.add_atom('C', 1) # save atom unchanged.
sn2patch.add_atom('C', 2) # it is central atom.
sn2patch.add_atom('I', 3, charge=0, p_charge=-1)
sn2patch.add_atom('O', 4, charge=-1, p_charge=0)

sn2patch.add_bond(1, 2, 1) # set carbon - carbon single bond that is unchanged in SN2 reaction
sn2patch.add_bond(1, 3, DynamicBond(1)) # this bond is broken in SN2 reaction
sn2patch.add_bond(1, 4, DynamicBond(None, 1)) # it corresponds to formation of bond O(S)-C bond in SN2 reaction
print(sn2patch)
sn2patch.clean2d()
sn2patch

reactor = CGRReactor(ReactionContainer([e2query], [sn2patch])) # create template and pass it to Reactor
sn2_cgr = next(reactor(e2_cgr)) # apply Reactor on E2 reaction

print(sn2_cgr)
sn2_cgr

# decompose CGR into reaction
sn2 = ReactionContainer.from_cgr(sn2_cgr)
sn2.clean2d()
sn2

"""
Explanation: 6.2.2 Reaction transformation
Example of E2 to SN2 transformation.
E2 and SN2 are concurrent reactions.
We can easily change reaction center of E2 reaction to SN2. It could be achieved by substitution of reaction center corresponding to double bond formation in E2 reaction by the one corresponding to formation of new single bond with base as in SN2.
End of explanation
"""


