forked from pierrepo/PBxplore
-
Notifications
You must be signed in to change notification settings - Fork 0
/
PBassign.py
executable file
·157 lines (134 loc) · 5.08 KB
/
PBassign.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
#! /usr/bin/env python
# -*- coding: utf-8 -*-
"""
Read PDB structures and assign protein blocs (PBs).
2013 - P. Poulain, A. G. de Brevern
"""
## Use print as a function for python 3 compatibility
from __future__ import print_function
## standard modules
import os
import sys
import glob
import argparse
## local module
import PBlib as PB
import PDBlib as PDB
# MDAnalysis is an optional requirement
try:
import MDAnalysis
except ImportError:
IS_MDANALYSIS = False
else:
IS_MDANALYSIS = True
## Python2/Python3 compatibility
# The range function in python 3 behaves as the range function in python 2
# and returns a generator rather than a list. To produce a list in python 3,
# one should use list(range). Here we change range to behave the same in
# python 2 and in python 3. In both cases, range will return a generator.
try:
range = xrange
except NameError:
pass
def user_inputs():
"""
Handle the user parameter from the command line
Parse the command line parameter and build the list of input files.
"""
parser = argparse.ArgumentParser(
description='Read PDB structures and assign protein blocs (PBs).')
# arguments
parser.add_argument("-p", action="append",
help=("name of a pdb file "
"or name of a directory containing pdb files"))
parser.add_argument("-o", action="store", required=True,
help="name for results")
# arguments for MDanalysis
group = parser.add_argument_group(
title='other options [if MDanalysis module is available]')
group.add_argument("-x", action="store",
help="name of xtc file (Gromacs)")
group.add_argument("-g", action="store",
help="name of gro file (Gromacs)")
# optional arguments
parser.add_argument("--phipsi", action="store_true", default=False,
help="writes phi and psi angle")
parser.add_argument("--flat", action="store_true", default=False,
help="writes one PBs sequence per line")
parser.add_argument('-v', '--version', action='version',
version='%(prog)s 1.0')
# get all arguments
options = parser.parse_args()
# check options
if not options.p:
if not IS_MDANALYSIS:
parser.error("MDAnalysis is not installed; options -x and -d are not available")
if not options.x:
parser.print_help()
parser.error("use at least option -p or -x")
elif not options.g:
parser.print_help()
parser.error("option -g is mandatory, with use of option -x")
# check files
pdb_name_lst = []
if options.p:
for name in options.p:
# input is a file: store file name
if os.path.isfile(name):
pdb_name_lst.append(name)
# input is a directory: list and store all PDB and PDBx/mmCIF files
elif os.path.isdir(name):
for extension in (PDB.PDB_EXTENSIONS + PDB.PDBx_EXTENSIONS):
pdb_name_lst += glob.glob(os.path.join(name, "*" + extension))
# input is not a file neither a directory: say it
elif (not os.path.isfile(name) or not os.path.isdir(name)):
parser.error("{0}: not a valid file or directory".format(name))
else:
if not os.path.isfile(options.x):
sys.exit("{0}: not a valid file".format(options.x))
elif not os.path.isfile(options.g):
sys.exit("{0}: not a valid file".format(options.g))
return options, pdb_name_lst
def pbassign_cli():
"""
PBassign command line.
"""
options, pdb_name_lst = user_inputs()
if options.p:
if pdb_name_lst:
print("{} PDB file(s) to process".format(len(pdb_name_lst)))
else:
print('Nothing to do. Good bye.')
return
# PB assignement of PDB structures
chains = PDB.chains_from_files(pdb_name_lst)
else:
# PB assignement of a Gromacs trajectory
chains = PDB.chains_from_trajectory(options.x, options.g)
all_comments = []
all_sequences = []
all_dihedrals = []
for comment, chain in chains:
dihedrals = chain.get_phi_psi_angles()
sequence = PB.assign(dihedrals)
all_comments.append(comment)
all_dihedrals.append(dihedrals)
all_sequences.append(sequence)
fasta_name = options.o + ".PB.fasta"
with open(fasta_name, 'w') as outfile:
PB.write_fasta(outfile, all_sequences, all_comments)
if options.flat:
flat_name = options.o + ".PB.flat"
with open(flat_name, 'w') as outfile:
PB.write_flat(outfile, all_sequences)
if options.phipsi:
phipsi_name = options.o + ".PB.phipsi"
with open(phipsi_name, 'w') as outfile:
PB.write_phipsi(outfile, all_dihedrals, all_comments)
print("wrote {0}".format(fasta_name))
if options.flat:
print("wrote {0}".format(flat_name))
if options.phipsi:
print("wrote {0}".format(phipsi_name))
if __name__ == '__main__':
pbassign_cli()