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