/usr/share/pyshared/Bio/PDB/PDBParser.py is in python-biopython 1.58-1.
This file is owned by root:root, with mode 0o644.
The actual contents of the file can be viewed below.
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 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 | # Copyright (C) 2002, Thomas Hamelryck (thamelry@binf.ku.dk)
# This code is part of the Biopython distribution and governed by its
# license. Please see the LICENSE file that should have been included
# as part of this package.
"""Parser for PDB files."""
import warnings
import numpy
from Bio.PDB.PDBExceptions import \
PDBConstructionException, PDBConstructionWarning
from Bio.PDB.StructureBuilder import StructureBuilder
from Bio.PDB.parse_pdb_header import _parse_pdb_header_list
# If PDB spec says "COLUMNS 18-20" this means line[17:20]
class PDBParser(object):
"""
Parse a PDB file and return a Structure object.
"""
def __init__(self, PERMISSIVE=True, get_header=False,
structure_builder=None, QUIET=False):
"""
The PDB parser call a number of standard methods in an aggregated
StructureBuilder object. Normally this object is instanciated by the
PDBParser object itself, but if the user provides his own StructureBuilder
object, the latter is used instead.
Arguments:
o PERMISSIVE - Evaluated as a Boolean. If false, exceptions in
constructing the SMCRA data structure are fatal. If true (DEFAULT),
the exceptions are caught, but some residues or atoms will be missing.
THESE EXCEPTIONS ARE DUE TO PROBLEMS IN THE PDB FILE!.
o structure_builder - an optional user implemented StructureBuilder class.
o QUIET - Evaluated as a Boolean. If true, warnings issued in constructing
the SMCRA data will be supressed. If false (DEFAULT), they will be shown.
These warnings might be indicative of problems in the PDB file!
"""
if structure_builder!=None:
self.structure_builder=structure_builder
else:
self.structure_builder=StructureBuilder()
self.header=None
self.trailer=None
self.line_counter=0
self.PERMISSIVE=bool(PERMISSIVE)
self.QUIET=bool(QUIET)
# Public methods
def get_structure(self, id, file):
"""Return the structure.
Arguments:
o id - string, the id that will be used for the structure
o file - name of the PDB file OR an open filehandle
"""
if self.QUIET:
warning_list = warnings.filters[:]
warnings.filterwarnings('ignore', category=PDBConstructionWarning)
self.header=None
self.trailer=None
# Make a StructureBuilder instance (pass id of structure as parameter)
self.structure_builder.init_structure(id)
handle_close = False
if isinstance(file, basestring):
file=open(file)
handle_close = True
self._parse(file.readlines())
self.structure_builder.set_header(self.header)
# Return the Structure instance
structure = self.structure_builder.get_structure()
if handle_close:
file.close()
if self.QUIET:
warnings.filters = warning_list
return structure
def get_header(self):
"Return the header."
return self.header
def get_trailer(self):
"Return the trailer."
return self.trailer
# Private methods
def _parse(self, header_coords_trailer):
"Parse the PDB file."
# Extract the header; return the rest of the file
self.header, coords_trailer=self._get_header(header_coords_trailer)
# Parse the atomic data; return the PDB file trailer
self.trailer=self._parse_coordinates(coords_trailer)
def _get_header(self, header_coords_trailer):
"Get the header of the PDB file, return the rest."
structure_builder=self.structure_builder
i = 0
for i in range(0, len(header_coords_trailer)):
structure_builder.set_line_counter(i+1)
line=header_coords_trailer[i]
record_type=line[0:6]
if(record_type=='ATOM ' or record_type=='HETATM' or record_type=='MODEL '):
break
header=header_coords_trailer[0:i]
# Return the rest of the coords+trailer for further processing
self.line_counter=i
coords_trailer=header_coords_trailer[i:]
header_dict=_parse_pdb_header_list(header)
return header_dict, coords_trailer
def _parse_coordinates(self, coords_trailer):
"Parse the atomic data in the PDB file."
local_line_counter=0
structure_builder=self.structure_builder
current_model_id=0
# Flag we have an open model
model_open=0
current_chain_id=None
current_segid=None
current_residue_id=None
current_resname=None
for i in range(0, len(coords_trailer)):
line=coords_trailer[i]
record_type=line[0:6]
global_line_counter=self.line_counter+local_line_counter+1
structure_builder.set_line_counter(global_line_counter)
if(record_type=='ATOM ' or record_type=='HETATM'):
# Initialize the Model - there was no explicit MODEL record
if not model_open:
structure_builder.init_model(current_model_id)
current_model_id+=1
model_open=1
fullname=line[12:16]
# get rid of whitespace in atom names
split_list=fullname.split()
if len(split_list)!=1:
# atom name has internal spaces, e.g. " N B ", so
# we do not strip spaces
name=fullname
else:
# atom name is like " CA ", so we can strip spaces
name=split_list[0]
altloc=line[16:17]
resname=line[17:20]
chainid=line[21:22]
try:
serial_number=int(line[6:11])
except:
serial_number=0
resseq=int(line[22:26].split()[0]) # sequence identifier
icode=line[26:27] # insertion code
if record_type=='HETATM': # hetero atom flag
if resname=="HOH" or resname=="WAT":
hetero_flag="W"
else:
hetero_flag="H"
else:
hetero_flag=" "
residue_id=(hetero_flag, resseq, icode)
# atomic coordinates
try:
x=float(line[30:38])
y=float(line[38:46])
z=float(line[46:54])
except:
#Should we allow parsing to continue in permissive mode?
#If so what coordindates should we default to? Easier to abort!
raise PDBConstructionException(\
"Invalid or missing coordinate(s) at line %i." \
% global_line_counter)
coord=numpy.array((x, y, z), 'f')
# occupancy & B factor
try:
occupancy=float(line[54:60])
except:
self._handle_PDB_exception("Invalid or missing occupancy",
global_line_counter)
occupancy = 0.0 #Is one or zero a good default?
try:
bfactor=float(line[60:66])
except:
self._handle_PDB_exception("Invalid or missing B factor",
global_line_counter)
bfactor = 0.0 #The PDB use a default of zero if the data is missing
segid=line[72:76]
element=line[76:78].strip()
if current_segid!=segid:
current_segid=segid
structure_builder.init_seg(current_segid)
if current_chain_id!=chainid:
current_chain_id=chainid
structure_builder.init_chain(current_chain_id)
current_residue_id=residue_id
current_resname=resname
try:
structure_builder.init_residue(resname, hetero_flag, resseq, icode)
except PDBConstructionException, message:
self._handle_PDB_exception(message, global_line_counter)
elif current_residue_id!=residue_id or current_resname!=resname:
current_residue_id=residue_id
current_resname=resname
try:
structure_builder.init_residue(resname, hetero_flag, resseq, icode)
except PDBConstructionException, message:
self._handle_PDB_exception(message, global_line_counter)
# init atom
try:
structure_builder.init_atom(name, coord, bfactor, occupancy, altloc,
fullname, serial_number, element)
except PDBConstructionException, message:
self._handle_PDB_exception(message, global_line_counter)
elif(record_type=='ANISOU'):
anisou=map(float, (line[28:35], line[35:42], line[43:49], line[49:56], line[56:63], line[63:70]))
# U's are scaled by 10^4
anisou_array=(numpy.array(anisou, 'f')/10000.0).astype('f')
structure_builder.set_anisou(anisou_array)
elif(record_type=='MODEL '):
try:
serial_num=int(line[10:14])
except:
self._handle_PDB_exception("Invalid or missing model serial number",
global_line_counter)
serial_num=0
structure_builder.init_model(current_model_id,serial_num)
current_model_id+=1
model_open=1
current_chain_id=None
current_residue_id=None
elif(record_type=='END ' or record_type=='CONECT'):
# End of atomic data, return the trailer
self.line_counter=self.line_counter+local_line_counter
return coords_trailer[local_line_counter:]
elif(record_type=='ENDMDL'):
model_open=0
current_chain_id=None
current_residue_id=None
elif(record_type=='SIGUIJ'):
# standard deviation of anisotropic B factor
siguij=map(float, (line[28:35], line[35:42], line[42:49], line[49:56], line[56:63], line[63:70]))
# U sigma's are scaled by 10^4
siguij_array=(numpy.array(siguij, 'f')/10000.0).astype('f')
structure_builder.set_siguij(siguij_array)
elif(record_type=='SIGATM'):
# standard deviation of atomic positions
sigatm=map(float, (line[30:38], line[38:45], line[46:54], line[54:60], line[60:66]))
sigatm_array=numpy.array(sigatm, 'f')
structure_builder.set_sigatm(sigatm_array)
local_line_counter=local_line_counter+1
# EOF (does not end in END or CONECT)
self.line_counter=self.line_counter+local_line_counter
return []
def _handle_PDB_exception(self, message, line_counter):
"""
This method catches an exception that occurs in the StructureBuilder
object (if PERMISSIVE), or raises it again, this time adding the
PDB line number to the error message.
"""
message="%s at line %i." % (message, line_counter)
if self.PERMISSIVE:
# just print a warning - some residues/atoms may be missing
warnings.warn("PDBConstructionException: %s\n"
"Exception ignored.\n"
"Some atoms or residues may be missing in the data structure."
% message, PDBConstructionWarning)
else:
# exceptions are fatal - raise again with new message (including line nr)
raise PDBConstructionException(message)
if __name__=="__main__":
import sys
p=PDBParser(PERMISSIVE=True)
filename = sys.argv[1]
s=p.get_structure("scr", filename)
for m in s:
p=m.get_parent()
assert(p is s)
for c in m:
p=c.get_parent()
assert(p is m)
for r in c:
print r
p=r.get_parent()
assert(p is c)
for a in r:
p=a.get_parent()
if not p is r:
print p, r
|