cmd=r'''
topo
! patch to turn serine into phosphoserine
! this should probably be turned into a distinct residue
! CDS 2005/7/6
!
presidue TP
     group
     modify atom 1cb charge=0.0 end 
     modify atom 1og1  charge=0.0 end 
     delete atom 1hg1  charge=0.0 end
     group
     modify atom 2P charge=0.0 end
     modify atom 2O1P charge=0.0 end
     modify atom 2O2P charge=0.0 end
     modify atom 2O3P charge=0.0 end 
   add bond 1og1 2P
   add angle 1cb 1og1 2p
   add angle 1og1   2P 2O1P 
   add angle 1og1   2P 2O2P
   add angle 1og1   2P 2O3P 
end
end
'''

seq="""ALA THR ALA"""

import protocol
protocol.initTopology("pho.top")
protocol.initParams("pho.par")

import psfGen
psfGen.seqToPSF(seq)
psfGen.seqToPSF("PHO",startResid=500,seqType='metal')

xplor.command(cmd)
xplor.command("patch TP reference=1=(resid 2) reference=2=(resid 500) end ")

xplor.command("write structure output=tp.psf end")

try:
    protocol.genExtendedStructure()
except protocol.CovalentViolation:
    pass


protocol.writePDB("extended.pdb")