FreeSASA  2.0
Open source SASA calculations
View on GitHub
freesasa.pyx
1 ## @package freesasa
2 # @author Simon Mitternacht
3 # @copyright [MIT License](md_license.html)
4 #
5 # @brief FreeSASA Python interface.
6 #
7 # A minimal program would be something like
8 #
9 # ~~~{.py}
10 # structure = freesasa.Structure("1abc.pdb")
11 # result = freesasa.calc(structure)
12 # print "SASA = %f" % result.totalArea()
13 # ~~~
14 #
15 # See documentation of the classes and functions for how to customize behavior.
16 #
17 
18 from libc.stdio cimport FILE, fopen, fclose
19 from libc.stdlib cimport free, realloc, malloc
20 from libc.string cimport memcpy
21 from cfreesasa cimport *
22 
23 ## Used to specify the algorithm by Shrake & Rupley
24 ShrakeRupley = 'ShrakeRupley'
25 
26 ## Used to specify the algorithm by Lee & Richards
27 LeeRichards = 'LeeRichards'
28 
29 ## Used for classification
30 polar = 'Polar'
31 
32 ## Used for classification
33 apolar = 'Apolar'
34 
35 ## int: Suppress all warnings and errors (used by setVerbosity())
36 silent = FREESASA_V_SILENT
37 
38 ## int: Suppress all warnings but not errors (used by setVerbosity())
39 nowarnings = FREESASA_V_NOWARNINGS
40 
41 ## int: Normal verbosity (used by setVerbosity())
42 normal = FREESASA_V_NORMAL
43 
44 ## int: Print debug messages (used by setVerbosity())
45 debug = FREESASA_V_DEBUG
46 
47 ## The default values for calculation parameters
48 defaultParameters = {
49  'algorithm' : LeeRichards,
50  'probe-radius' : freesasa_default_parameters.probe_radius,
51  'n-points' : freesasa_default_parameters.shrake_rupley_n_points,
52  'n-slices' : freesasa_default_parameters.lee_richards_n_slices,
53  'n-threads' : freesasa_default_parameters.n_threads
54 }
55 
56 ## Stores parameter values to be used by calculation.
57 #
58 # Wraps the C struct freesasa_parameters.
59 cdef class Parameters:
60 
61  cdef freesasa_parameters _c_param
62  ## Initializes Parameters object.
63  #
64  # @param param (dict) optional argument to specify parameter-values, modeled after
65  # ::defaultParameters
66  # @exception AssertionError invalid parameter values supplied
67  # @see defaultParameters
68  def __init__(self,param=None):
69 
70  self._c_param = freesasa_default_parameters
71  if param != None:
72  if 'algorithm' in param: self.setAlgorithm(param['algorithm'])
73  if 'probe-radius' in param: self.setProbeRadius(param['probe-radius'])
74  if 'n-points' in param: self.setNPoints(param['n-points'])
75  if 'n-slices' in param: self.setNSlices(param['n-slices'])
76  if 'n-threads' in param: self.setNThreads(param['n-threads'])
77  unknownKeys = []
78  for key in param:
79  if not key in defaultParameters:
80  unknownKeys.append(key)
81  if len(unknownKeys) > 0:
82  raise AssertionError('Key(s): ',unknownKeys,', unknown');
83 
84  ## Set algorithm.
85  #
86  # @param alg (str) algorithm name, only allowed values are ::ShrakeRupley and ::LeeRichards
87  # @exception AssertionError unknown algorithm specified
88  def setAlgorithm(self,alg):
89  if alg == ShrakeRupley:
90  self._c_param.alg = FREESASA_SHRAKE_RUPLEY
91  elif alg == LeeRichards:
92  self._c_param.alg = FREESASA_LEE_RICHARDS
93  else:
94  raise AssertionError("Algorithm '%s' is unknown" % alg)
95 
96  ## Get algorithm.
97  #
98  # @return Name of algorithm
99  def algorithm(self):
100  if self._c_param.alg == FREESASA_SHRAKE_RUPLEY:
101  return ShrakeRupley
102  if self._c_param.alg == FREESASA_LEE_RICHARDS:
103  return LeeRichards
104  raise Exception("No algorithm specified, shouldn't be possible")
105 
106  ## Set probe radius.
107  # @param r probe radius in Å (>= 0)
108  # @exception AssertionError r < 0
109  def setProbeRadius(self,r):
110  assert(r >= 0)
111  self._c_param.probe_radius = r
112 
113  ## Get probe radius.
114  # @return Probe radius in Å
115  def probeRadius(self):
116  return self._c_param.probe_radius
117 
118  ## Set number of test points in Shrake & Rupley algorithm.
119  # @param n (int) Number of points (> 0).
120  # @exception AssertionError n <= 0.
121  def setNPoints(self,n):
122  assert(n > 0)
123  self._c_param.shrake_rupley_n_points = n
124 
125  ## Get number of test points in Shrake & Rupley algorithm.
126  # @return Number of points.
127  def nPoints(self):
128  return self._c_param.shrake_rupley_n_points
129 
130  ## Set the number of slices per atom in Lee & Richards algorithm.
131  # @param n (int) Number of slices (> 0)
132  # @exception AssertionError n <= 0
133  def setNSlices(self,n):
134  assert(n> 0)
135  self._c_param.lee_richards_n_slices = n
136 
137  ## Get the number of slices per atom in Lee & Richards algorithm.
138  # @return Number of slices.
139  def nSlices(self):
140  return self._c_param.lee_richards_n_slices
141 
142 
143  ## Set the number of threads to use in calculations.
144  # @param n (int) Number of points (> 0)
145  # @exception AssertionError n <= 0
146  def setNThreads(self,n):
147  assert(n>0)
148  self._c_param.n_threads = n
149 
150 
151  ## Get the number of threads to use in calculations.
152  # @return Number of threads.
153  def nThreads(self):
154  return self._c_param.n_threads
155 
156  # not pretty, but only way I've found to pass pointers around
157  def _get_address(self, size_t ptr2ptr):
158  cdef freesasa_parameters **p = <freesasa_parameters**> ptr2ptr
159  p[0] = &self._c_param
160 
161 ## Stores results from SASA calculation.
162 # The type of object returned by calc(), not intended to be used
163 # outside of that context.
164 #
165 # Wraps the C struct freesasa_result.
166 cdef class Result:
167  cdef freesasa_result* _c_result
168 
169  ## The constructor
170  def __cinit__ (self):
171  self._c_result = NULL
172 
173  ## The destructor
174  def __dealloc__(self):
175  if self._c_result is not NULL:
177 
178  ## Number of atoms in the results
179  # @return Number of atoms
180  def nAtoms(self):
181  if self._c_result is not NULL:
182  return self._c_result.n_atoms
183  return 0
184 
185 
186  ## Total SASA.
187  # @return The total area in Å^2.
188  # @exception AssertionError If no results have been associated
189  # with the object
190  def totalArea(self):
191  assert(self._c_result is not NULL)
192  return self._c_result.total
193 
194  ## SASA for a given atom.
195  # @param i (int) index of atom.
196  # @return SASA of atom i in Å^2.
197  # @exception AssertionError If no results have been associated
198  # with the object or if index is out of bounds
199  def atomArea(self,i):
200  assert(self._c_result is not NULL)
201  assert(i < self._c_result.n_atoms)
202  return self._c_result.sasa[i]
203 
204  def _get_address(self, size_t ptr2ptr):
205  cdef freesasa_result **p = <freesasa_result**> ptr2ptr
206  p[0] = self._c_result
207 
208 
209 ## Assigns class and radius to atom by residue and atom name.
210 #
211 # Subclasses derived from Classifier can be used to define custom
212 # atomic radii and/or classes. Can also be initialized from a @ref
213 # Config-file "configuration file" with a custom classifier.
214 #
215 # Wraps a C freesasa_classifier. If initialized without arguments the
216 # default classifier is used.
217 #
218 # Residue names should be of the format `"ALA"`,`"ARG"`, etc.
219 #
220 # Atom names should be of the format `"CA"`, `"N"`, etc.
221 #
222 cdef class Classifier:
223  cdef freesasa_classifier* _c_classifier
224 
225  ## Constructor.
226  #
227  # If no file is provided the default classifier is used.
228  #
229  # @see @ref Config-file.
230  #
231  # @param fileName Name of file with classifier configuration.
232  # @exception IOError Problem opening/reading file
233  # @exception Exception Problem parsing provided configuration or
234  # initializing defaults
235  def __cinit__ (self,fileName=None):
236  cdef FILE *config
237  self._c_classifier = NULL
238  if fileName is not None:
239  config = fopen(fileName,'r')
240  if config is NULL:
241  raise IOError("File '%s' could not be opened." % fileName)
243  fclose(config)
244  if self._c_classifier is NULL:
245  raise Exception("Error parsing configuration in '%s'." % fileName)
246  else:
247  self._c_classifier = &freesasa_default_classifier
248  ## The destructor
249  def __dealloc__(self):
250  if self._c_classifier is not &freesasa_default_classifier:
252 
253  ## Class of atom.
254  #
255  # Depending on the configuration these classes can be
256  # anything, but typically they will be 'Polar' and 'Apolar'.
257  # Unrecognized atoms will get the class 'Unknown'.
258  #
259  # @param residueName (str) Residue name (`"ALA"`,`"ARG"`,...).
260  # @param atomName (str) Atom name (`"CA"`,`"C"`,...).
261  # @return A string describing the class
262  def classify(self,residueName,atomName):
263  classIndex = freesasa_classifier_class(self._c_classifier, residueName, atomName)
264  return freesasa_classifier_class2str(classIndex)
265 
266  ## Radius of atom.
267  #
268  # This allows the classifier to be used to calculate the atomic
269  # radii used in calculations. Unknown atoms will get a negative
270  # radius.
271  #
272  # @param residueName (str) Residue name (`"ALA"`, `"ARG"`, ...).
273  # @param atomName (str) Atom name (`"CA"`, `"C"`, ...).
274  # @return The radius in Å.
275  def radius(self,residueName,atomName):
276  return freesasa_classifier_radius(self._c_classifier, residueName, atomName)
277 
278  def _get_address(self, size_t ptr2ptr):
279  cdef freesasa_classifier **p = <freesasa_classifier**> ptr2ptr
280  p[0] = self._c_classifier
281 
282 ## Represents a protein structure, including its atomic radii.
283 #
284 # Initialized from PDB-file. Calculates atomic radii using default
285 # classifier, or custom one provided as argument to initalizer
286 #
287 # Wraps the C struct freesasa_structure.
288 #
289 # Since it is intended to be a static structure the word 'get' is
290 # omitted in the getter-functions.
291 cdef class Structure:
292  cdef freesasa_structure* _c_structure
293  ## By default ignore HETATM, Hydrogens, only use first model. For unknown atoms
294  # try to guess the radius, if this fails, assign radius 0 (to
295  # allow changing the radius later).
296  defaultOptions = { 'hetatm' : False,
297  'hydrogen' : False,
298  'join-models' : False,
299  'skip-unknown' : False,
300  'halt-at-unknown' : False
301  }
302 
303  ## Constructor
304  #
305  # If PDB file is provided, the structure will be constructed
306  # based on the file. If not, this simply initializes an empty
307  # structure and the other arguments are ignored. In this case
308  # atoms will have to be added manually using addAtom().
309  #
310  # @param fileName (str) PDB file (if None empty structure generated).
311  # @param classifier An optional Classifier to calculate atomic
312  # radii, uses default if none provided
313  # @param options specify which atoms and models to include
314  # @exception IOError Problem opening/reading file.
315  # @exception Exception Problem parsing PDB file or calculating
316  # atomic radii.
317  # @exception Exception If option 'halt-at-unknown' selected and
318  # unknown atom encountered.
319  def __init__(self,fileName=None,classifier=None,
320  options = defaultOptions):
321 
322  self._c_structure = NULL
323 
324  if fileName is None:
326  return
327  cdef FILE *input
328  input = fopen(fileName,'r')
329  if input is NULL:
330  raise IOError("File '%s' could not be opened." % fileName)
331  structure_options = Structure._get_structure_options(options)
332  self._c_structure = freesasa_structure_from_pdb(input,NULL,structure_options)
333  fclose(input)
334 
335  if self._c_structure is NULL:
336  raise Exception("Error reading '%s'." % fileName)
337 
338  # this means we might be calculating the radii twice, the
339  # advantage of doing it this way is that we can define new
340  # classifiers using a Python interface
341  if (classifier is not None):
342  self.setRadiiWithClassifier(classifier)
343 
344 
345  ## Add atom to structure.
346  #
347  # This function is meant to be used if the structure was not
348  # initialized from a PDB. Default radii will be assigned to each
349  # atom. This can be overriden by calling
350  # setRadiiWithClassifier() afterwards.
351  #
352  # There are no restraints on string lengths for the arguments, but
353  # the atom won't be added if the default classifier doesn't
354  # recognize the atom and also cannot deduce its element from the
355  # atom name.
356  #
357  # @param atomName (str) atom name (e.g. `"CA"`)
358  # @param residueName (str) residue name (e.g. `"ALA"`)
359  # @param residueNumber (str or int) residue number (e.g. `'12'`)
360  # or integer. Some PDBs have residue-numbers that aren't
361  # regular numbers. Therefore treated as a string primarily.
362  # @param chainLabel (str) 1-character string with chain label (e.g. 'A')
363  # @param x,y,z (float) coordinates
364  #
365  # @exception Exception Residue-number invalid
366  def addAtom(self, atomName, residueName, residueNumber, chainLabel, x, y, z):
367  if (type(residueNumber) is str):
368  resnum = residueNumber
369  elif (type(residueNumber) is int):
370  resnum = "%d" % residueNumber
371  else:
372  raise Exception("Residue-number invalid, must be either string or number")
373  cdef const char *label = chainLabel
374  ret = freesasa_structure_add_atom(self._c_structure, atomName,
375  residueName, resnum, label[0],
376  x, y, z)
377  assert(ret != FREESASA_FAIL)
378 
379  ## Assign radii to atoms in structure using a classifier.
380  #
381  # @param classifier A classifier to use to calculate radii
382  # @exception AssertionError if structure not properly initialized
383  def setRadiiWithClassifier(self,classifier):
384  assert(self._c_structure is not NULL)
385  n = self.nAtoms()
386  r = []
387  for i in range(0,n):
388  r.append(classifier.radius(self.residueName(i), self.atomName(i)))
389  self.setRadii(r)
390 
391  ## Set atomic radii from an array
392  # @param radiusArray Array of atomic radii in Ångström, should
393  # have nAtoms() elements.
394  # @exception AssertionError if radiusArray has wrong dimension, structure
395  # not properly initialized, or if the array contains
396  # negative radii (not properly classified?)
397  def setRadii(self,radiusArray):
398  assert(self._c_structure is not NULL)
399  n = self.nAtoms()
400  assert len(radiusArray) == n
401  cdef double *r = <double *>malloc(sizeof(double)*n)
402  assert(r is not NULL)
403  for i in range(0,n):
404  r[i] = radiusArray[i]
405  assert(r[i] >= 0)
407 
408  ## Number of atoms.
409  #
410  # @return Number of atoms
411  # @exception AssertionError if not properly initialized
412  def nAtoms(self):
413  assert(self._c_structure is not NULL)
415 
416  ## Radius of atom.
417  # @param i (int) Index of atom.
418  # @return Radius in Å.
419  # @exception AssertionError if index out of bounds, object not properly initalized.
420  def radius(self,i):
421  assert(i >= 0 and i < self.nAtoms())
422  assert(self._c_structure is not NULL)
423  cdef const double *r = freesasa_structure_radius(self._c_structure)
424  assert(r is not NULL)
425  return r[i]
426 
427  ## Set radius for a given atom
428  # @param atomIndex Index of atom
429  # @param radius Value of radius
430  # @exception AssertionError if index out of bounds, radius
431  # negative, or structure not properly initialized
432  def setRadius(self, atomIndex, radius):
433  assert(self._c_structure is not NULL)
434  assert(atomIndex >= 0 and atomIndex < self.nAtoms())
435  assert(radius >= 0)
436  freesasa_structure_atom_set_radius(self._c_structure, atomIndex, radius);
437 
438  ## Get atom name
439  # @param i (int) Atom index.
440  # @return Atom name as 4-character string.
441  # @exception AssertionError: if index out of range or Structure not properly initialized.
442  def atomName(self,i):
443  assert(i >= 0 and i < self.nAtoms())
444  assert(self._c_structure is not NULL)
446 
447  ## Get residue name of given atom.
448  # @param i (int) Atom index.
449  # @return Residue name as 3-character string.
450  # @exception AssertionError if index out of range or Structure not properly initialized
451  def residueName(self,i):
452  assert(i >= 0 and i < self.nAtoms())
453  assert(self._c_structure is not NULL)
455 
456  ## Get residue number for given atom
457  # @param i (int) Atom index.
458  # @return Residue number as 4-character string
459  # @exception AssertionError if index out of range or Structure not properly initialized
460  def residueNumber(self,i):
461  assert(i >= 0 and i < self.nAtoms())
462  assert(self._c_structure is not NULL)
464 
465  ## Get chain label for given atom.
466  # @param i (int) Atom index.
467  # @return Chain label as 1-character string
468  # @exception AssertionError if index out of range or Structure not properly initialized
469  def chainLabel(self,i):
470  assert(i >= 0 and i < self.nAtoms())
471  assert(self._c_structure is not NULL)
472  cdef char label[2]
473  label[0] = freesasa_structure_atom_chain(self._c_structure,i);
474  label[1] = '\0'
475  return label
476 
477  ## Get coordinates of given atom.
478  # @param i (int) Atom index.
479  # @return array of x, y, and z coordinates
480  # @exception AssertionError if index out of range or Structure not properly initialized
481  def coord(self, i):
482  assert(i >= 0 and i < self.nAtoms())
483  assert(self._c_structure is not NULL)
484  cdef const double *coord = freesasa_structure_coord_array(self._c_structure);
485  return [coord[3*i], coord[3*i+1], coord[3*i+2]]
486 
487  @staticmethod
488  def _get_structure_options(param):
489  options = 0
490 
491  # check validity of options
492  knownOptions = {'hetatm','hydrogen','join-models','separate-models',
493  'separate-chains','skip-unknown','halt-at-unknown'}
494  unknownOptions = []
495  for key in param:
496  if not key in knownOptions:
497  unknownOptions.append(key)
498  if len(unknownOptions) > 0:
499  raise AssertionError("Option(s): ",unknownOptions," unknown.");
500 
501  # calculate bitfield
502  if 'hetatm' in param and param['hetatm']:
503  options |= FREESASA_INCLUDE_HETATM
504  if 'hydrogen' in param and param['hydrogen']:
505  options |= FREESASA_INCLUDE_HYDROGEN
506  if 'join-models' in param and param['join-models']:
507  options |= FREESASA_JOIN_MODELS
508  if 'separate-models' in param and param['separate-models']:
509  options |= FREESASA_SEPARATE_MODELS
510  if 'separate-chains' in param and param['separate-chains']:
511  options |= FREESASA_SEPARATE_CHAINS
512  if 'skip-unknown' in param and param['skip-unknown']:
513  options |= FREESASA_SKIP_UNKNOWN
514  if 'halt-at-unknown' in param and param['halt-at-unknown']:
515  options |= FREESASA_HALT_AT_UNKNOWN
516  return options
517 
518  def _get_address(self, size_t ptr2ptr):
519  cdef freesasa_structure **p = <freesasa_structure**> ptr2ptr
520  p[0] = self._c_structure
521 
522  def _set_address(self, size_t ptr2ptr):
523  cdef freesasa_structure **p = <freesasa_structure**> ptr2ptr
524  self._c_structure = p[0]
525 
526  ## The destructor
527  def __dealloc__(self):
528  if self._c_structure is not NULL:
530 
531 ## Default options for structureArray.
532 # Defined separately for Doxygen's sake.
533 defaultStructureArrayOptions = {
534  'hetatm' : False,
535  'hydrogen' : False,
536  'separate-chains' : True,
537  'separate-models' : False
538 }
539 
540 ## Create array of structures from PDB file.
541 #
542 # Split PDB file into several structures by either by treating
543 # chains separately, by treating each MODEL as a separate
544 # structure, or both.
545 #
546 # @param fileName (str) The PDB file.
547 # @param options (dict) Specification for how to read the PDB-file
548 # (see def value for options).
549 # @param classifier: Classifier to assign atoms radii, default is used
550 # if none specified.
551 # @exception AssertionError if fileName is None
552 # @exception AssertionError if an option value is not recognized
553 # @exception AssertionError if neither of the options 'separate-chains'
554 # and 'separate-models' are specified.
555 # @exception IOError if can't open file
556 # @exception Exception: if there are problems parsing the input
557 # @return An array of Structures
558 def structureArray(fileName,
559  options = defaultStructureArrayOptions,
560  classifier = None):
561  assert fileName is not None
562  # we need to have at least one of these
563  assert(('separate-chains' in options and options['separate-chains'] is True)
564  or ('separate-models' in options and options['separate-models'] is True))
565  structure_options = Structure._get_structure_options(options)
566  cdef FILE *input
567  input = fopen(fileName,'r')
568  if input is NULL:
569  raise IOError("File '%s' could not be opened." % fileName)
570  cdef int n;
571  cdef freesasa_structure** sArray = freesasa_structure_array(input,&n,NULL,structure_options)
572  fclose(input)
573  if sArray is NULL:
574  raise Exception("Problems reading structures in '%s'." % fileName)
575  structures = []
576  for i in range(0,n):
577  structures.append(Structure())
578  structures[-1]._set_address(<size_t> &sArray[i])
579  if classifier is not None:
580  structures[-1].setRadiiWithClassifier(classifier)
581  free(sArray)
582  return structures
583 
584 
585 ## Calculate SASA of Structure
586 # @param structure Structure to be used
587 # @param parameters Parameters to use (if not specified defaults are used)
588 # @return A Result object
589 # @exception Exception something went wrong in calculation (see C library error messages)
590 def calc(structure,parameters=None):
591  cdef const freesasa_parameters *p = NULL
592  cdef const freesasa_structure *s = NULL
593  if parameters is not None: parameters._get_address(<size_t>&p)
594  structure._get_address(<size_t>&s)
595  result = Result()
596  result._c_result = <freesasa_result*> freesasa_calc_structure(s,p)
597  if result._c_result is NULL:
598  raise Exception("Error calculating SASA.")
599  return result
600 
601 ## Break SASA result down into classes.
602 # @param result Result from sasa calculation.
603 # @param structure Structure used in calculation.
604 # @param classifier Classifier (if not specified default is used).
605 # @return Dictionary with names of classes as keys and their SASA values as values.
606 # @exception Exception: Problems with classification, see C library error messages
607 # (or Python exceptions if run with derived classifier).
608 def classifyResults(result,structure,classifier=None):
609  if classifier is None:
610  classifier = Classifier()
611  ret = dict()
612  for i in range(0,structure.nAtoms()):
613  name = classifier.classify(structure.residueName(i),structure.atomName(i))
614  if name not in ret:
615  ret[name] = 0
616  ret[name] += result.atomArea(i)
617  return ret
618 
619 ## Sum SASA result over a selection of atoms
620 # @param commands A list of commands with selections using Pymol
621 # syntax, e.g. `"s1, resn ala+arg"` or `"s2, chain A and resi 1-5"`
622 # (see @ref Selection).
623 # @param structure A Structure.
624 # @param result Result from sasa calculation on structure.
625 # @return Dictionary with names of selections (`"s1"`,`"s2"`,...) as
626 # keys, and the corresponding SASA values as values.
627 # @exception Exception: Parser failed (typically syntax error), see
628 # library error messages.
629 def selectArea(commands, structure, result):
630  cdef freesasa_structure *s
631  cdef freesasa_result *r
632  cdef freesasa_selection *selection
633  structure._get_address(<size_t> &s)
634  result._get_address(<size_t> &r)
635  value = dict()
636  for cmd in commands:
637  selection = freesasa_selection_new(cmd, s, r)
638  if selection == NULL:
639  raise Exception("Error parsing '%s'" % cmd)
640  value[freesasa_selection_name(selection)] = freesasa_selection_area(selection)
641  freesasa_selection_free(selection)
642  return value
643 
644 ## Set global verbosity
645 # @param verbosity Can have values freesasa.silent, freesasa.nowarnings or freesasa.normal
646 # @exception AssertionError if verbosity has illegal value
647 def setVerbosity(verbosity):
648  assert(verbosity in [silent, nowarnings, normal])
649  freesasa_set_verbosity(verbosity)
650 
651 
652 ## Get global verbosity
653 # @return Verbosity (freesasa.silent, freesasa.nowarnings or freesasa.normal)
655  return freesasa_get_verbosity()
656 
657 ## Create a freesasa structure from a Bio.PDB structure
658 #
659 # @remark Experimental, not thorougly tested yet
660 #
661 # @param bioPDBStructure a Bio.PDB structure
662 # @param classifier an optional classifier to specify atomic radii
663 # @param options Options supported are 'hetatm', 'skip-unknown' and 'halt-at-unknown'
664 # @return a freesasa.Structure
665 #
666 # @exception Exception if option 'halt-at-unknown' is selected and
667 # unknown atoms are encountered. Passes on exceptions from
668 # Structure.addAtom() and
669 # Structure.setRadiiWithClassifier().
670 def structureFromBioPDB(bioPDBStructure, classifier=None, options = Structure.defaultOptions):
671  structure = Structure()
672  if (classifier is None):
673  classifier = Classifier()
674  optbitfield = Structure._get_structure_options(options)
675 
676  atoms = bioPDBStructure.get_atoms()
677 
678  for a in atoms:
679  r = a.get_parent()
680  hetflag, resseq, icode = r.get_id()
681 
682  if (hetflag is not ' ' and not (optbitfield & FREESASA_INCLUDE_HETATM)):
683  continue
684 
685  c = r.get_parent()
686  v = a.get_vector()
687 
688  if (classifier.classify(r.get_resname(), a.get_fullname()) is 'Unknown'):
689  if (optbitfield & FREESASA_SKIP_UNKNOWN):
690  continue
691  if (optbitfield & FREESASA_HALT_AT_UNKNOWN):
692  raise Exception("Halting at unknown atom")
693 
694  structure.addAtom(a.get_fullname(), r.get_resname(), resseq, c.get_id(),
695  v[0], v[1], v[2])
696 
697  structure.setRadiiWithClassifier(classifier)
698  return structure
699 
700 ## Calc SASA from Bio.PDB structure
701 #
702 # Usage
703 #
704 # result, sasa_classes = calcBioPDB(structure, ...)
705 #
706 # @remark Experimental, not thorougly tested yet
707 #
708 # @param bioPDBStructure A Bio.PDB structure
709 # @param parameters A freesasa.Paremeters object
710 # @param classifier A freesasa.Classifier object
711 # @param options Options supported are 'hetatm', 'skip-unknown' and 'halt-at-unknown'
712 #
713 # @return A freesasa.Result object and a dictionary with classes
714 # defined by the classifier and associated areas
715 #
716 # @exception Exception if unknown atom is encountered and the option
717 # 'halt-at-unknown' is active. Passes on exceptions from
718 # calc(), classifyResults() and structureFromBioPDB().
719 def calcBioPDB(bioPDBStructure, parameters = Parameters(),
720  classifier = None, options = Structure.defaultOptions):
721  structure = structureFromBioPDB(bioPDBStructure, classifier, options)
722  result = calc(structure, parameters)
723  sasa_classes = classifyResults(result, structure, classifier)
724  return result, sasa_classes
725 
726 
freesasa_verbosity freesasa_get_verbosity(void)
Get the current verbosity level.
const char * freesasa_structure_atom_name(const freesasa_structure *structure, int i)
Get atom name.
def structureFromBioPDB
Create a freesasa structure from a Bio.PDB structure.
Definition: freesasa.pyx:670
int freesasa_structure_n(const freesasa_structure *structure)
Get number of atoms.
void freesasa_selection_free(freesasa_selection *selection)
Free selection.
def setProbeRadius(self, r)
Set probe radius.
Definition: freesasa.pyx:109
Assigns class and radius to atom by residue and atom name.
Definition: freesasa.pyx:222
Stores parameter values to be used by calculation.
Definition: freesasa.pyx:59
def calc
Calculate SASA of Structure.
Definition: freesasa.pyx:590
def __init__
Constructor.
Definition: freesasa.pyx:320
def probeRadius(self)
Get probe radius.
Definition: freesasa.pyx:115
def getVerbosity()
Get global verbosity.
Definition: freesasa.pyx:654
const char * freesasa_structure_atom_res_name(const freesasa_structure *structure, int i)
Get residue name.
def __dealloc__(self)
The destructor.
Definition: freesasa.pyx:527
void freesasa_structure_set_radius(freesasa_structure *structure, const double *radii)
Override the radii set when the structure was initialized.
int freesasa_set_verbosity(freesasa_verbosity v)
Set the global verbosity level.
def selectArea(commands, structure, result)
Sum SASA result over a selection of atoms.
Definition: freesasa.pyx:629
const char * freesasa_selection_name(const freesasa_selection *selection)
Name of the selection.
def atomArea(self, i)
SASA for a given atom.
Definition: freesasa.pyx:199
freesasa_selection * freesasa_selection_new(const char *command, const freesasa_structure *structure, const freesasa_result *result)
Get area of a selection.
freesasa_result * freesasa_calc_structure(const freesasa_structure *structure, const freesasa_parameters *parameters)
Calculates SASA based on a given structure.
def __init__
Initializes Parameters object.
Definition: freesasa.pyx:68
def nSlices(self)
Get the number of slices per atom in Lee & Richards algorithm.
Definition: freesasa.pyx:139
def structureArray
Create array of structures from PDB file.
Definition: freesasa.pyx:560
def radius(self, residueName, atomName)
Radius of atom.
Definition: freesasa.pyx:275
freesasa_atom_class freesasa_classifier_class(const freesasa_classifier *classifier, const char *res_name, const char *atom_name)
Use a classifier to determine the class of a given atom.
freesasa_structure * freesasa_structure_from_pdb(FILE *pdb, const freesasa_classifier *classifier, int options)
Init structure with coordinates from pdb-file.
Represents a protein structure, including its atomic radii.
Definition: freesasa.pyx:291
def addAtom(self, atomName, residueName, residueNumber, chainLabel, x, y, z)
Add atom to structure.
Definition: freesasa.pyx:366
double freesasa_selection_area(const freesasa_selection *selection)
Area of the selection.
void freesasa_classifier_free(freesasa_classifier *classifier)
Frees a classifier object.
freesasa_structure ** freesasa_structure_array(FILE *pdb, int *n, const freesasa_classifier *classifier, int options)
Init array of structures from PDB.
double freesasa_classifier_radius(const freesasa_classifier *classifier, const char *res_name, const char *atom_name)
Use a classifier to determine the radius of a given atom.
void freesasa_result_free(freesasa_result *result)
Frees a freesasa_result object.
int freesasa_structure_add_atom(freesasa_structure *structure, const char *atom_name, const char *residue_name, const char *residue_number, char chain_label, double x, double y, double z)
Add individual atom to structure using default behavior.
def setRadiiWithClassifier(self, classifier)
Assign radii to atoms in structure using a classifier.
Definition: freesasa.pyx:383
def radius(self, i)
Radius of atom.
Definition: freesasa.pyx:420
def nAtoms(self)
Number of atoms.
Definition: freesasa.pyx:412
def __cinit__(self)
The constructor.
Definition: freesasa.pyx:170
def calcBioPDB
Calc SASA from Bio.PDB structure.
Definition: freesasa.pyx:720
def atomName(self, i)
Get atom name.
Definition: freesasa.pyx:442
def setRadii(self, radiusArray)
Set atomic radii from an array.
Definition: freesasa.pyx:397
const double * freesasa_structure_radius(const freesasa_structure *structure)
Returns a pointer to an array of the radii of each atom.
def nThreads(self)
Get the number of threads to use in calculations.
Definition: freesasa.pyx:153
def classifyResults
Break SASA result down into classes.
Definition: freesasa.pyx:608
void freesasa_structure_free(freesasa_structure *structure)
Free structure.
def coord(self, i)
Get coordinates of given atom.
Definition: freesasa.pyx:481
def nAtoms(self)
Number of atoms in the results.
Definition: freesasa.pyx:180
def residueNumber(self, i)
Get residue number for given atom.
Definition: freesasa.pyx:460
def algorithm(self)
Get algorithm.
Definition: freesasa.pyx:99
const char * freesasa_classifier_class2str(freesasa_atom_class atom_class)
Names for freesasa_atom_class.
def setRadius(self, atomIndex, radius)
Set radius for a given atom.
Definition: freesasa.pyx:432
def setVerbosity(verbosity)
Set global verbosity.
Definition: freesasa.pyx:647
freesasa_classifier * freesasa_classifier_from_file(FILE *file)
Generate a classifier from a config-file.
def __cinit__
Constructor.
Definition: freesasa.pyx:235
char freesasa_structure_atom_chain(const freesasa_structure *structure, int i)
Get chain label.
def __dealloc__(self)
The destructor.
Definition: freesasa.pyx:249
def nPoints(self)
Get number of test points in Shrake & Rupley algorithm.
Definition: freesasa.pyx:127
def setNPoints(self, n)
Set number of test points in Shrake & Rupley algorithm.
Definition: freesasa.pyx:121
const double * freesasa_structure_coord_array(const freesasa_structure *structure)
Get array of coordinates.
freesasa_structure * freesasa_structure_new(void)
Allocate empty structure.
def setAlgorithm(self, alg)
Set algorithm.
Definition: freesasa.pyx:88
Stores results from SASA calculation.
Definition: freesasa.pyx:166
def residueName(self, i)
Get residue name of given atom.
Definition: freesasa.pyx:451
void freesasa_structure_atom_set_radius(freesasa_structure *structure, int i, double radius)
Set atom radius.
def setNThreads(self, n)
Set the number of threads to use in calculations.
Definition: freesasa.pyx:146
def classify(self, residueName, atomName)
Class of atom.
Definition: freesasa.pyx:262
def totalArea(self)
Total SASA.
Definition: freesasa.pyx:190
const char * freesasa_structure_atom_res_number(const freesasa_structure *structure, int i)
Get residue number.
def chainLabel(self, i)
Get chain label for given atom.
Definition: freesasa.pyx:469
def setNSlices(self, n)
Set the number of slices per atom in Lee & Richards algorithm.
Definition: freesasa.pyx:133
def __dealloc__(self)
The destructor.
Definition: freesasa.pyx:174