"""
Modeling of point mutations with MODELLER.

Based on mutate_model.py from the MODELLER wiki (written by Ben Webb?):

http://salilab.org/modeller/wiki/Mutate%20model

TODO

* Fix modeller's treatment of residue numbering for outputted models so that reference model can start with standard numbering.
* Allow use of multiple reference models to build all mutants.
* Add treatment of double mutants.

"""

#===============================================================================
# IMPORTS
#===============================================================================

import sys
import os

#===============================================================================
# MUTATIONS TO MAKE
#===============================================================================

# List of mutations to make.
mutation_list = [
    # Wild-type.
    '',
    # Point mutants.
    'V281C',
    'M314H', 'M314Y', 'M314L',
    'V323C', 'V323S', 'V323T',
    'L325Y',
    'T338M', 'T338N', 'T338Q', 'T338V', 'T338I', 'T338F', 'T338L', 
    'Y340F', 'Y340H', 'Y340K', 
    'M341C', 'M341H',
    'A403C', 'A403S', 'A403T', 'A403V',
    # Double mutants.
    'V323C A403T', # (mimics EGFR ATP-binding site)
    'V323S A403T', # (mimics Her2 ATP-binding site)
    'T338M A403T', # (mimics PKA ATP-binding site)
    'T338M A403C', #
    'T338M A403S', #
    # Triple mutants.
    'M314L L325Y T338M', # kinases that bind bosutinib tightly but do not look like Src have these three mutations
    ]

#===============================================================================
# MAIN
#===============================================================================

if __name__ == '__main__':

    reference_model_filename = 'src-reference-model.pdb'
    chain = 'B'

    import mutate
    import amber
    import re
    import os.path

    directory = 'src-reference-mutants'
    os.mkdir(directory)

    for mutation in mutation_list:
        # Make mutation.
        print "Building mutation %s..." % (mutation)
        receptor_pdb_filename = '%s/Src %s.pdb' % (directory, mutation)
        mutate.mutate(reference_model_filename, mutation, chain=chain, verbose=False, outfile=receptor_pdb_filename)

        # Create AMBER system.
        print "Creating AMBER system..."
        crystallographic_waters_pdb_filename = '../modeller/closewater.pdb'
        ligand_frcmod_filename = '../antechamber/oechem/bosutinib.gaff.frcmod'
        ligand_mol2_filename = '../antechamber/oechem/bosutinib.gaff.mol2'
        prefix = '%s/Src %s' % (directory, mutation)
        ncations = 7
        nanions = 0        
        amber.create_system(receptor_pdb_filename, ligand_frcmod_filename, ligand_mol2_filename, prefix, ncations=ncations, nanions=nanions, otherpdb=crystallographic_waters_pdb_filename)