#!/usr/bin/env __XPLOR_DIR__/bin/pyXplor


from sys import argv, exit
usage=r"""usage:
    %s [-protein | -dna | -rna] [option] <sequence filename>
          
     where option is one or more of
       -segname <segment name>
       -startresid <lowest residue number>
       -cisPeptide resid
       -disulfide_bond resid1:resid2
       -disulfide_bridge resid1:resid2
       -protonateHIS
       -amidate_cterm
       -custom_rename

      If the sequence type is omitted, a guess is made- this will not work
      correctly for rna sequences- you must specify -rna.
    
      generates a psf file from protein or nucleic acid sequence
      the psf file is written in a filename derived the input filename.
   
      multiple chains (not covalently linked) can be obtained by using the
      special residue name '*'.
    
      use -disulfide_bond to make actual bonds, and -disulfide_bridge
      to remove the cysteine HG proton- to form disulfide with an
      NOE restraint.

      use -amidate_cterm to terminate the c-terminal end with a CONH2 group.

      by default the HIS HD1 atom is deleted. If -protonateHIS is specified,
      this proton is not deleted.

      If -custom_rename is specified, certain convenient atom renamings are 
      made for nucleic acids:
            ADE H61 --> HN'   
            ADE H62 --> HN''  
            GUA H21 --> HN'   
            GUA H22 --> HN''  
            CYT H41 --> HN'   
            CYT H42 --> HN''  
            THY C5A --> CM    

""" % argv[0]


seqType='auto'
segname='    '
seqFiles=[]
startResid=1
disulfide_bonds=[]
disulfide_bridges=[]
amidate_cterm=0
custom_rename=0
deprotonateHIS=True
cisPeptides=[]

cnt=1
while cnt<len(argv):
    if argv[cnt] == '-protein':
        seqType='prot'
    elif argv[cnt] == '-dna':
        seqType='dna'
    elif argv[cnt] == '-rna':
        seqType='rna'
    elif argv[cnt] == '-segname':
        segname=argv[cnt+1]
        cnt += 1
    elif argv[cnt] == '-protonateHIS':
	deprotonateHIS=False
    elif argv[cnt] == '-startresid':
        startResid=int(argv[cnt+1])
        cnt += 1
    elif argv[cnt] == '-disulfide_bond':
        cnt += 1
        disulfide_bonds.append( map(int,argv[cnt].split(':')) )
    elif argv[cnt] == '-disulfide_bridge':
        cnt += 1
        disulfide_bridges.append( map(int,argv[cnt].split(':')) )
    elif argv[cnt] == '-amidate_cterm':
        amidate_cterm=1
    elif argv[cnt] == '-custom_rename':
        custom_rename=1
    elif argv[cnt].startswith('-h') or argv[cnt].startswith('--h'):
        print usage
        exit(0)
    elif argv[cnt] == '-cisPeptide':
        cisPeptides.append( int(argv[cnt+1]) )
        cnt += 1
    elif argv[cnt][0] == '-':
        print "invalid argument:", argv[cnt]
        print usage
        exit(1)
    else:
        seqFiles.append(argv[cnt])
        pass
    cnt += 1
    pass

if not seqFiles: 
	print usage
	exit(1)

seq = open(seqFiles[0]).read()

for file in seqFiles[1:]:
    seq += ' * '
    seq += open(file).read()
    pass

seqFile = seqFiles[0]
dotIndex=seqFile.rfind('.')
outFile=seqFile
if dotIndex>=0:
    outFile=seqFile[:dotIndex]
    pass
outFile += ".psf"

seqs = seq.split('*')

from psfGen import seqToPSF, cisPeptide
for seq in seqs:
    seqToPSF(seq,
             startResid=startResid,
             seqType=seqType,
             segName=segname,
	     deprotonateHIS=deprotonateHIS,
             disulfide_bonds=disulfide_bonds,
             disulfide_bridges=disulfide_bridges,
             amidate_cterm=amidate_cterm,
             customRename=custom_rename         )
    startResid += len(seq.split())
    pass

for resid in cisPeptides:
    cisPeptide(resid)
    pass


xplorCommand="remark Created with the command:\n"
cmd='  '
for arg in argv:
    cmd += arg + ' '
    if len(cmd)>70:
        xplorCommand += "remark %s\n" % cmd
        cmd='    '
        pass
    pass
if len(cmd)>4:
    xplorCommand += "remark %s\n" % cmd
    pass

xplor.command(xplorCommand)
xplor.command('write structure output = "%s" end' % outFile)