FreeSASA  2.0.1
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 # Derived classifiers should set the member purePython to True
219 #
220 # Residue names should be of the format `"ALA"`,`"ARG"`, etc.
221 #
222 # Atom names should be of the format `"CA"`, `"N"`, etc.
223 #
224 cdef class Classifier:
225  cdef freesasa_classifier* _c_classifier
226  purePython = False
227 
228  ## Constructor.
229  #
230  # If no file is provided the default classifier is used.
231  #
232  # @see @ref Config-file.
233  #
234  # @param fileName Name of file with classifier configuration.
235  # @exception IOError Problem opening/reading file
236  # @exception Exception Problem parsing provided configuration or
237  # initializing defaults
238  def __cinit__ (self,fileName=None):
239  cdef FILE *config
240  self._c_classifier = NULL
241  if fileName is not None:
242  config = fopen(fileName,'r')
243  if config is NULL:
244  raise IOError("File '%s' could not be opened." % fileName)
246  fclose(config)
247  if self._c_classifier is NULL:
248  raise Exception("Error parsing configuration in '%s'." % fileName)
249  else:
250  self._c_classifier = &freesasa_default_classifier
251  ## The destructor
252  def __dealloc__(self):
253  if self._c_classifier is not &freesasa_default_classifier:
255 
256 
257  # This is used internally to determine if a Classifier wraps a C
258  # classifier or not (necessary when generating structures)
259  # @return Boolean
260  def _isCClassifier(self):
261  return not self.purePython
262 
263  ## Class of atom.
264  #
265  # Depending on the configuration these classes can be
266  # anything, but typically they will be 'Polar' and 'Apolar'.
267  # Unrecognized atoms will get the class 'Unknown'.
268  #
269  # @param residueName (str) Residue name (`"ALA"`,`"ARG"`,...).
270  # @param atomName (str) Atom name (`"CA"`,`"C"`,...).
271  # @return A string describing the class
272  def classify(self,residueName,atomName):
273  classIndex = freesasa_classifier_class(self._c_classifier, residueName, atomName)
274  return freesasa_classifier_class2str(classIndex)
275 
276  ## Radius of atom.
277  #
278  # This allows the classifier to be used to calculate the atomic
279  # radii used in calculations. Unknown atoms will get a negative
280  # radius.
281  #
282  # @param residueName (str) Residue name (`"ALA"`, `"ARG"`, ...).
283  # @param atomName (str) Atom name (`"CA"`, `"C"`, ...).
284  # @return The radius in Å.
285  def radius(self,residueName,atomName):
286  return freesasa_classifier_radius(self._c_classifier, residueName, atomName)
287 
288  def _get_address(self, size_t ptr2ptr):
289  cdef freesasa_classifier **p = <freesasa_classifier**> ptr2ptr
290  p[0] = self._c_classifier
291 
292 ## Represents a protein structure, including its atomic radii.
293 #
294 # Initialized from PDB-file. Calculates atomic radii using default
295 # classifier, or custom one provided as argument to initalizer
296 #
297 # Wraps the C struct freesasa_structure.
298 #
299 # Since it is intended to be a static structure the word 'get' is
300 # omitted in the getter-functions.
301 cdef class Structure:
302  cdef freesasa_structure* _c_structure
303  ## By default ignore HETATM, Hydrogens, only use first model. For unknown atoms
304  # try to guess the radius, if this fails, assign radius 0 (to
305  # allow changing the radius later).
306  defaultOptions = { 'hetatm' : False,
307  'hydrogen' : False,
308  'join-models' : False,
309  'skip-unknown' : False,
310  'halt-at-unknown' : False
311  }
312 
313  ## Constructor
314  #
315  # If PDB file is provided, the structure will be constructed
316  # based on the file. If not, this simply initializes an empty
317  # structure and the other arguments are ignored. In this case
318  # atoms will have to be added manually using addAtom().
319  #
320  # @param fileName (str) PDB file (if None empty structure generated).
321  # @param classifier An optional Classifier to calculate atomic
322  # radii, uses default if none provided
323  # @param options specify which atoms and models to include
324  # @exception IOError Problem opening/reading file.
325  # @exception Exception Problem parsing PDB file or calculating
326  # atomic radii.
327  # @exception Exception If option 'halt-at-unknown' selected and
328  # unknown atom encountered.
329  def __init__(self,fileName=None,classifier=None,
330  options = defaultOptions):
331 
332  self._c_structure = NULL
333  cdef freesasa_classifier *c = NULL
334  if classifier is None:
335  classifier = Classifier()
336  if classifier._isCClassifier():
337  classifier._get_address(<size_t>&c)
338 
339  if fileName is None:
341  return
342  cdef FILE *input
343  input = fopen(fileName,'r')
344  if input is NULL:
345  raise IOError("File '%s' could not be opened." % fileName)
346  structure_options = Structure._get_structure_options(options)
347 
348  if not classifier._isCClassifier(): # supress warnings
349  setVerbosity(silent)
350 
351  self._c_structure = freesasa_structure_from_pdb(input, c, structure_options)
352 
353  if not classifier._isCClassifier():
354  setVerbosity(normal)
355 
356  fclose(input)
357 
358  if self._c_structure is NULL:
359  raise Exception("Error reading '%s'." % fileName)
360 
361  # for pure Python classifiers we use the default
362  # classifier above to initialize the structure and then
363  # reassign radii using the provided classifier here
364  if (not classifier._isCClassifier()):
365  self.setRadiiWithClassifier(classifier)
366 
367 
368  ## Add atom to structure.
369  #
370  # This function is meant to be used if the structure was not
371  # initialized from a PDB. Default radii will be assigned to each
372  # atom. This can be overriden by calling
373  # setRadiiWithClassifier() afterwards.
374  #
375  # There are no restraints on string lengths for the arguments, but
376  # the atom won't be added if the default classifier doesn't
377  # recognize the atom and also cannot deduce its element from the
378  # atom name.
379  #
380  # @param atomName (str) atom name (e.g. `"CA"`)
381  # @param residueName (str) residue name (e.g. `"ALA"`)
382  # @param residueNumber (str or int) residue number (e.g. `'12'`)
383  # or integer. Some PDBs have residue-numbers that aren't
384  # regular numbers. Therefore treated as a string primarily.
385  # @param chainLabel (str) 1-character string with chain label (e.g. 'A')
386  # @param x,y,z (float) coordinates
387  #
388  # @exception Exception Residue-number invalid
389  def addAtom(self, atomName, residueName, residueNumber, chainLabel, x, y, z):
390  if (type(residueNumber) is str):
391  resnum = residueNumber
392  elif (type(residueNumber) is int):
393  resnum = "%d" % residueNumber
394  else:
395  raise Exception("Residue-number invalid, must be either string or number")
396  cdef const char *label = chainLabel
397  ret = freesasa_structure_add_atom(self._c_structure, atomName,
398  residueName, resnum, label[0],
399  x, y, z)
400  assert(ret != FREESASA_FAIL)
401 
402  ## Assign radii to atoms in structure using a classifier.
403  #
404  # @param classifier A classifier to use to calculate radii
405  # @exception AssertionError if structure not properly initialized
406  def setRadiiWithClassifier(self,classifier):
407  assert(self._c_structure is not NULL)
408  n = self.nAtoms()
409  r = []
410  for i in range(0,n):
411  r.append(classifier.radius(self.residueName(i), self.atomName(i)))
412  self.setRadii(r)
413 
414  ## Set atomic radii from an array
415  # @param radiusArray Array of atomic radii in Ångström, should
416  # have nAtoms() elements.
417  # @exception AssertionError if radiusArray has wrong dimension, structure
418  # not properly initialized, or if the array contains
419  # negative radii (not properly classified?)
420  def setRadii(self,radiusArray):
421  assert(self._c_structure is not NULL)
422  n = self.nAtoms()
423  assert len(radiusArray) == n
424  cdef double *r = <double *>malloc(sizeof(double)*n)
425  assert(r is not NULL)
426  for i in range(0,n):
427  r[i] = radiusArray[i]
428  assert(r[i] >= 0)
430 
431  ## Number of atoms.
432  #
433  # @return Number of atoms
434  # @exception AssertionError if not properly initialized
435  def nAtoms(self):
436  assert(self._c_structure is not NULL)
438 
439  ## Radius of atom.
440  # @param i (int) Index of atom.
441  # @return Radius in Å.
442  # @exception AssertionError if index out of bounds, object not properly initalized.
443  def radius(self,i):
444  assert(i >= 0 and i < self.nAtoms())
445  assert(self._c_structure is not NULL)
446  cdef const double *r = freesasa_structure_radius(self._c_structure)
447  assert(r is not NULL)
448  return r[i]
449 
450  ## Set radius for a given atom
451  # @param atomIndex Index of atom
452  # @param radius Value of radius
453  # @exception AssertionError if index out of bounds, radius
454  # negative, or structure not properly initialized
455  def setRadius(self, atomIndex, radius):
456  assert(self._c_structure is not NULL)
457  assert(atomIndex >= 0 and atomIndex < self.nAtoms())
458  assert(radius >= 0)
459  freesasa_structure_atom_set_radius(self._c_structure, atomIndex, radius);
460 
461  ## Get atom name
462  # @param i (int) Atom index.
463  # @return Atom name as 4-character string.
464  # @exception AssertionError: if index out of range or Structure not properly initialized.
465  def atomName(self,i):
466  assert(i >= 0 and i < self.nAtoms())
467  assert(self._c_structure is not NULL)
469 
470  ## Get residue name of given atom.
471  # @param i (int) Atom index.
472  # @return Residue name as 3-character string.
473  # @exception AssertionError if index out of range or Structure not properly initialized
474  def residueName(self,i):
475  assert(i >= 0 and i < self.nAtoms())
476  assert(self._c_structure is not NULL)
478 
479  ## Get residue number for given atom
480  # @param i (int) Atom index.
481  # @return Residue number as 4-character string
482  # @exception AssertionError if index out of range or Structure not properly initialized
483  def residueNumber(self,i):
484  assert(i >= 0 and i < self.nAtoms())
485  assert(self._c_structure is not NULL)
487 
488  ## Get chain label for given atom.
489  # @param i (int) Atom index.
490  # @return Chain label as 1-character string
491  # @exception AssertionError if index out of range or Structure not properly initialized
492  def chainLabel(self,i):
493  assert(i >= 0 and i < self.nAtoms())
494  assert(self._c_structure is not NULL)
495  cdef char label[2]
496  label[0] = freesasa_structure_atom_chain(self._c_structure,i);
497  label[1] = '\0'
498  return label
499 
500  ## Get coordinates of given atom.
501  # @param i (int) Atom index.
502  # @return array of x, y, and z coordinates
503  # @exception AssertionError if index out of range or Structure not properly initialized
504  def coord(self, i):
505  assert(i >= 0 and i < self.nAtoms())
506  assert(self._c_structure is not NULL)
507  cdef const double *coord = freesasa_structure_coord_array(self._c_structure);
508  return [coord[3*i], coord[3*i+1], coord[3*i+2]]
509 
510  @staticmethod
511  def _get_structure_options(param):
512  options = 0
513 
514  # check validity of options
515  knownOptions = {'hetatm','hydrogen','join-models','separate-models',
516  'separate-chains','skip-unknown','halt-at-unknown'}
517  unknownOptions = []
518  for key in param:
519  if not key in knownOptions:
520  unknownOptions.append(key)
521  if len(unknownOptions) > 0:
522  raise AssertionError("Option(s): ",unknownOptions," unknown.");
523 
524  # calculate bitfield
525  if 'hetatm' in param and param['hetatm']:
526  options |= FREESASA_INCLUDE_HETATM
527  if 'hydrogen' in param and param['hydrogen']:
528  options |= FREESASA_INCLUDE_HYDROGEN
529  if 'join-models' in param and param['join-models']:
530  options |= FREESASA_JOIN_MODELS
531  if 'separate-models' in param and param['separate-models']:
532  options |= FREESASA_SEPARATE_MODELS
533  if 'separate-chains' in param and param['separate-chains']:
534  options |= FREESASA_SEPARATE_CHAINS
535  if 'skip-unknown' in param and param['skip-unknown']:
536  options |= FREESASA_SKIP_UNKNOWN
537  if 'halt-at-unknown' in param and param['halt-at-unknown']:
538  options |= FREESASA_HALT_AT_UNKNOWN
539  return options
540 
541  def _get_address(self, size_t ptr2ptr):
542  cdef freesasa_structure **p = <freesasa_structure**> ptr2ptr
543  p[0] = self._c_structure
544 
545  def _set_address(self, size_t ptr2ptr):
546  cdef freesasa_structure **p = <freesasa_structure**> ptr2ptr
547  self._c_structure = p[0]
548 
549  ## The destructor
550  def __dealloc__(self):
551  if self._c_structure is not NULL:
553 
554 ## Default options for structureArray.
555 # Defined separately for Doxygen's sake.
556 defaultStructureArrayOptions = {
557  'hetatm' : False,
558  'hydrogen' : False,
559  'separate-chains' : True,
560  'separate-models' : False
561 }
562 
563 ## Create array of structures from PDB file.
564 #
565 # Split PDB file into several structures by either by treating
566 # chains separately, by treating each MODEL as a separate
567 # structure, or both.
568 #
569 # @param fileName (str) The PDB file.
570 # @param options (dict) Specification for how to read the PDB-file
571 # (see def value for options).
572 # @param classifier: Classifier to assign atoms radii, default is used
573 # if none specified.
574 # @exception AssertionError if fileName is None
575 # @exception AssertionError if an option value is not recognized
576 # @exception AssertionError if neither of the options 'separate-chains'
577 # and 'separate-models' are specified.
578 # @exception IOError if can't open file
579 # @exception Exception: if there are problems parsing the input
580 # @return An array of Structures
581 def structureArray(fileName,
582  options = defaultStructureArrayOptions,
583  classifier = None):
584  assert fileName is not None
585  # we need to have at least one of these
586  assert(('separate-chains' in options and options['separate-chains'] is True)
587  or ('separate-models' in options and options['separate-models'] is True))
588  structure_options = Structure._get_structure_options(options)
589  cdef FILE *input
590  input = fopen(fileName,'r')
591  if input is NULL:
592  raise IOError("File '%s' could not be opened." % fileName)
593  cdef int n;
594  cdef freesasa_structure** sArray = freesasa_structure_array(input,&n,NULL,structure_options)
595  fclose(input)
596  if sArray is NULL:
597  raise Exception("Problems reading structures in '%s'." % fileName)
598  structures = []
599  for i in range(0,n):
600  structures.append(Structure())
601  structures[-1]._set_address(<size_t> &sArray[i])
602  if classifier is not None:
603  structures[-1].setRadiiWithClassifier(classifier)
604  free(sArray)
605  return structures
606 
607 
608 ## Calculate SASA of Structure
609 # @param structure Structure to be used
610 # @param parameters Parameters to use (if not specified defaults are used)
611 # @return A Result object
612 # @exception Exception: something went wrong in calculation (see C library error messages)
613 def calc(structure,parameters=None):
614  cdef const freesasa_parameters *p = NULL
615  cdef const freesasa_structure *s = NULL
616  if parameters is not None: parameters._get_address(<size_t>&p)
617  structure._get_address(<size_t>&s)
618  result = Result()
619  result._c_result = <freesasa_result*> freesasa_calc_structure(s,p)
620  if result._c_result is NULL:
621  raise Exception("Error calculating SASA.")
622  return result
623 
624 ## Calculate SASA for a set of coordinates and radii
625 # @param coord list of size 3*N with atomic coordinates (x1, y1, z1,
626 # x2, y2, z2, ..., x_N, y_N, z_N'.
627 # @param radii array of size N with atomic radii (r_1, r_2, ..., r_N)
628 # @param Parameters to use (if not specified, defaults are used)
629 # @exception AssertionError: mismatched array-sizes
630 # @exception Exception: Out of memory
631 # @exception Exception: something went wrong in calculation (see C library error messages)
632 def calcCoord(coord, radii, parameters=None):
633  assert(len(coord) == 3*len(radii))
634 
635  cdef const freesasa_parameters *p = NULL
636  cdef double *c = <double*> malloc(len(coord)*sizeof(double))
637  cdef double *r = <double*> malloc(len(radii)*sizeof(double))
638  if c is NULL or r is NULL:
639  raise Exception("Memory allocation error")
640 
641  for i in xrange(len(coord)):
642  c[i] = coord[i]
643  for i in xrange(len(radii)):
644  r[i] = radii[i]
645 
646  if parameters is not None: parameters._get_address(<size_t>&p)
647 
648  result = Result()
649  result._c_result = <freesasa_result*> freesasa_calc_coord(c, r, len(radii), p)
650 
651  if result._c_result is NULL:
652  raise Exception("Error calculating SASA.")
653 
654  free(c)
655  free(r)
656 
657  return result
658 
659 ## Break SASA result down into classes.
660 # @param result Result from sasa calculation.
661 # @param structure Structure used in calculation.
662 # @param classifier Classifier (if not specified default is used).
663 # @return Dictionary with names of classes as keys and their SASA values as values.
664 # @exception Exception: Problems with classification, see C library error messages
665 # (or Python exceptions if run with derived classifier).
666 def classifyResults(result,structure,classifier=None):
667  if classifier is None:
668  classifier = Classifier()
669  ret = dict()
670  for i in range(0,structure.nAtoms()):
671  name = classifier.classify(structure.residueName(i),structure.atomName(i))
672  if name not in ret:
673  ret[name] = 0
674  ret[name] += result.atomArea(i)
675  return ret
676 
677 ## Sum SASA result over a selection of atoms
678 # @param commands A list of commands with selections using Pymol
679 # syntax, e.g. `"s1, resn ala+arg"` or `"s2, chain A and resi 1-5"`
680 # (see @ref Selection).
681 # @param structure A Structure.
682 # @param result Result from sasa calculation on structure.
683 # @return Dictionary with names of selections (`"s1"`,`"s2"`,...) as
684 # keys, and the corresponding SASA values as values.
685 # @exception Exception: Parser failed (typically syntax error), see
686 # library error messages.
687 def selectArea(commands, structure, result):
688  cdef freesasa_structure *s
689  cdef freesasa_result *r
690  cdef freesasa_selection *selection
691  structure._get_address(<size_t> &s)
692  result._get_address(<size_t> &r)
693  value = dict()
694  for cmd in commands:
695  selection = freesasa_selection_new(cmd, s, r)
696  if selection == NULL:
697  raise Exception("Error parsing '%s'" % cmd)
698  value[freesasa_selection_name(selection)] = freesasa_selection_area(selection)
699  freesasa_selection_free(selection)
700  return value
701 
702 ## Set global verbosity
703 # @param verbosity Can have values freesasa.silent, freesasa.nowarnings or freesasa.normal
704 # @exception AssertionError if verbosity has illegal value
705 def setVerbosity(verbosity):
706  assert(verbosity in [silent, nowarnings, normal])
707  freesasa_set_verbosity(verbosity)
708 
709 
710 ## Get global verbosity
711 # @return Verbosity (freesasa.silent, freesasa.nowarnings or freesasa.normal)
713  return freesasa_get_verbosity()
714 
715 ## Create a freesasa structure from a Bio.PDB structure
716 #
717 # @remark Experimental, not thorougly tested yet
718 #
719 # @param bioPDBStructure a Bio.PDB structure
720 # @param classifier an optional classifier to specify atomic radii
721 # @param options Options supported are 'hetatm', 'skip-unknown' and 'halt-at-unknown'
722 # @return a freesasa.Structure
723 #
724 # @exception Exception if option 'halt-at-unknown' is selected and
725 # unknown atoms are encountered. Passes on exceptions from
726 # Structure.addAtom() and
727 # Structure.setRadiiWithClassifier().
728 def structureFromBioPDB(bioPDBStructure, classifier=None, options = Structure.defaultOptions):
729  structure = Structure()
730  if (classifier is None):
731  classifier = Classifier()
732  optbitfield = Structure._get_structure_options(options)
733 
734  atoms = bioPDBStructure.get_atoms()
735 
736  for a in atoms:
737  r = a.get_parent()
738  hetflag, resseq, icode = r.get_id()
739 
740  if (hetflag is not ' ' and not (optbitfield & FREESASA_INCLUDE_HETATM)):
741  continue
742 
743  c = r.get_parent()
744  v = a.get_vector()
745 
746  if (classifier.classify(r.get_resname(), a.get_fullname()) is 'Unknown'):
747  if (optbitfield & FREESASA_SKIP_UNKNOWN):
748  continue
749  if (optbitfield & FREESASA_HALT_AT_UNKNOWN):
750  raise Exception("Halting at unknown atom")
751 
752  structure.addAtom(a.get_fullname(), r.get_resname(), resseq, c.get_id(),
753  v[0], v[1], v[2])
754 
755  structure.setRadiiWithClassifier(classifier)
756  return structure
757 
758 ## Calc SASA from Bio.PDB structure
759 #
760 # Usage
761 #
762 # result, sasa_classes = calcBioPDB(structure, ...)
763 #
764 # @remark Experimental, not thorougly tested yet
765 #
766 # @param bioPDBStructure A Bio.PDB structure
767 # @param parameters A freesasa.Paremeters object
768 # @param classifier A freesasa.Classifier object
769 # @param options Options supported are 'hetatm', 'skip-unknown' and 'halt-at-unknown'
770 #
771 # @return A freesasa.Result object and a dictionary with classes
772 # defined by the classifier and associated areas
773 #
774 # @exception Exception if unknown atom is encountered and the option
775 # 'halt-at-unknown' is active. Passes on exceptions from
776 # calc(), classifyResults() and structureFromBioPDB().
777 def calcBioPDB(bioPDBStructure, parameters = Parameters(),
778  classifier = None, options = Structure.defaultOptions):
779  structure = structureFromBioPDB(bioPDBStructure, classifier, options)
780  result = calc(structure, parameters)
781  sasa_classes = classifyResults(result, structure, classifier)
782  return result, sasa_classes
783 
784 
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:728
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:224
Stores parameter values to be used by calculation.
Definition: freesasa.pyx:59
def calc
Calculate SASA of Structure.
Definition: freesasa.pyx:613
def __init__
Constructor.
Definition: freesasa.pyx:330
def probeRadius(self)
Get probe radius.
Definition: freesasa.pyx:115
def getVerbosity()
Get global verbosity.
Definition: freesasa.pyx:712
const char * freesasa_structure_atom_res_name(const freesasa_structure *structure, int i)
Get residue name.
def __dealloc__(self)
The destructor.
Definition: freesasa.pyx:550
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:687
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:583
def radius(self, residueName, atomName)
Radius of atom.
Definition: freesasa.pyx:285
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:301
def addAtom(self, atomName, residueName, residueNumber, chainLabel, x, y, z)
Add atom to structure.
Definition: freesasa.pyx:389
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:406
def radius(self, i)
Radius of atom.
Definition: freesasa.pyx:443
def nAtoms(self)
Number of atoms.
Definition: freesasa.pyx:435
def __cinit__(self)
The constructor.
Definition: freesasa.pyx:170
def calcBioPDB
Calc SASA from Bio.PDB structure.
Definition: freesasa.pyx:778
def atomName(self, i)
Get atom name.
Definition: freesasa.pyx:465
def setRadii(self, radiusArray)
Set atomic radii from an array.
Definition: freesasa.pyx:420
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:666
void freesasa_structure_free(freesasa_structure *structure)
Free structure.
def coord(self, i)
Get coordinates of given atom.
Definition: freesasa.pyx:504
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:483
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:455
def setVerbosity(verbosity)
Set global verbosity.
Definition: freesasa.pyx:705
freesasa_classifier * freesasa_classifier_from_file(FILE *file)
Generate a classifier from a config-file.
def __cinit__
Constructor.
Definition: freesasa.pyx:238
char freesasa_structure_atom_chain(const freesasa_structure *structure, int i)
Get chain label.
def __dealloc__(self)
The destructor.
Definition: freesasa.pyx:252
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:474
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:272
def totalArea(self)
Total SASA.
Definition: freesasa.pyx:190
def calcCoord
Calculate SASA for a set of coordinates and radii.
Definition: freesasa.pyx:632
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:492
def setNSlices(self, n)
Set the number of slices per atom in Lee & Richards algorithm.
Definition: freesasa.pyx:133
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:174