BALL  1.4.1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
reducedSurface.h
Go to the documentation of this file.
1 // -*- Mode: C++; tab-width: 2; -*-
2 // vi: set ts=2:
3 //
4 
5 #ifndef BALL_STRUCTURE_REDUCEDSURFACE_H
6 #define BALL_STRUCTURE_REDUCEDSURFACE_H
7 
8 #ifndef BALL_MATHC_COMMON_H
9 # include <BALL/MATHS/common.h>
10 #endif
11 
12 #ifndef BALL_MATHS_SIMPLEBOX3_H
13 # include <BALL/MATHS/simpleBox3.h>
14 #endif
15 
16 #ifndef BALL_MATHS_CIRCLE3_H
17 # include <BALL/MATHS/circle3.h>
18 #endif
19 
20 #ifndef BALL_MATHS_SPHERE_H
21 # include <BALL/MATHS/sphere3.h>
22 #endif
23 
24 #ifndef BALL_MATHS_VECTOR3_H
25 # include <BALL/MATHS/vector3.h>
26 #endif
27 
28 #ifndef BALL_DATATYPE_HASHSET_H
29 # include <BALL/DATATYPE/hashMap.h>
30 #endif
31 
32 #ifndef BALL_DATATYPE_HASHSET_H
33 # include <BALL/DATATYPE/hashSet.h>
34 #endif
35 
36 #ifndef BALL_COMMON_EXCEPTION_H
37 # include <BALL/COMMON/exception.h>
38 #endif
39 
40 #ifndef BALL_STRUCTURE_RSEDGE_H
41 # include <BALL/STRUCTURE/RSEdge.h>
42 #endif
43 
44 #ifndef BALL_STRUCTURE_RSFACE_H
45 # include <BALL/STRUCTURE/RSFace.h>
46 #endif
47 
48 #ifndef BALL_STRUCTURE_RSVERTEX_H
49 # include <BALL/STRUCTURE/RSVertex.h>
50 #endif
51 
52 #include <set>
53 #include <list>
54 #include <deque>
55 #include <vector>
56 
57 namespace BALL
58 {
60  {
62  : a(a1), b(a2)
63  {
64  if (a > b) std::swap(a, b);
65  }
66 
67  bool operator==(const SortedPosition2& pos) const
68  {
69  return (a == pos.a) && (b == pos.b);
70  }
71 
74  };
75 
77  {
79  : a(a1), b(a2), c(a3)
80  {
81  if (a > b) std::swap(a, b);
82  if (a > c) std::swap(a, c);
83  if (b > c) std::swap(b, c);
84  }
85 
86  bool operator==(const SortedPosition3& pos) const
87  {
88  return (a == pos.a) && (b == pos.b) && (c == pos.c);
89  }
90 
94  };
95 }
96 
97 #ifdef BALL_HAS_BOOST_UNORDERED_MAP
98 namespace boost
99 {
100 #endif
101 
102 #ifdef BALL_EXTEND_HASH_IN_STD_NS
103 namespace std
104 {
105 #endif // BALL_EXTEND_HASH_IN_STD_NS
106 
107 #ifdef BALL_HAS_TR1_UNORDERED_MAP
108  namespace tr1
109  {
110 #endif // BALL_HAS_TR1_UNORDERED_MAP
111  template<>
112  struct hash<BALL::SortedPosition2> : public std::unary_function<BALL::SortedPosition2, size_t>
113  {
114  inline size_t operator()(const BALL::SortedPosition2& p) const
115  {
116  return 5323 * p.a + 1847 * p.b;
117  }
118  };
119 
120  template<>
121  struct hash<BALL::SortedPosition3> : public std::unary_function<BALL::SortedPosition3, size_t>
122  {
123  inline size_t operator()(const BALL::SortedPosition3& p) const
124  {
125  return 5323 * p.a + 1847 * p.b + 2347 * p.c;
126  }
127  };
128 #ifdef BALL_HAS_TR1_UNORDERED_MAP
129  }
130 #endif // BALL_HAS_TR1_UNORDERED_MAP
131 
132 #ifdef BALL_EXTEND_HASH_IN_STD_NS
133 }
134 #endif // BALL_EXTEND_HASH_IN_STD_NS
135 
136 #ifdef BALL_HAS_BOOST_UNORDERED_MAP
137 }
138 #endif
139 
140 namespace BALL
141 {
142  class RSComputer;
143  class SolventExcludedSurface;
144  class SESComputer;
145  class SESSingularityCleaner;
146  class TriangulatedSES;
147  class SolventAccessibleSurface;
148  class TriangulatedSAS;
149  class SESTriangulator;
150 
155  {
156  public:
157 
170  friend class RSComputer;
172  friend class SESComputer;
173  friend class SESSingularityCleaner;
175  friend class TriangulatedSES;
176  friend class TriangulatedSAS;
177  friend class SESTriangulator;
178 
180 
181 
184 
189  ReducedSurface();
190 
195  ReducedSurface(const ReducedSurface& reduced_surface, bool = true);
196 
200  ReducedSurface(const std::vector< TSphere3<double> >& spheres,
201  const double& probe_radius);
202 
205  virtual ~ReducedSurface();
206 
208 
211 
215  void operator = (const ReducedSurface& reduced_surface);
216 
220  void set(const ReducedSurface& reduced_surface);
221 
224  void clear();
225 
228  void clean();
229 
231 
234 
238  Size numberOfAtoms() const;
239 
243  Size numberOfVertices() const;
244 
248  Size numberOfEdges() const;
249 
253  Size numberOfFaces() const;
254 
258  double getProbeRadius() const;
259 
264  TSphere3<double> getSphere(Position i) const
265  throw(Exception::IndexOverflow);
266 
271  RSVertex* getVertex(Position i) const
272  throw(Exception::IndexOverflow);
273 
278  RSEdge* getEdge(Position i) const
279  throw(Exception::IndexOverflow);
280 
285  RSFace* getFace(Position i) const
286  throw(Exception::IndexOverflow);
287 
291  void insert(RSVertex* rsvertex);
292 
296  void insert(RSEdge* rsedge);
297 
301  void insert(RSFace* rsface);
302 
306  double getMaximalRadius() const;
307 
311  TSimpleBox3<double> getBoundingBox() const;
312 
317  void deleteSimilarFaces(RSFace* face1, RSFace* face2);
318 
327  bool getAngle(RSFace* face1, RSFace* face2,
328  RSVertex* vertex1, RSVertex* vertex2,
329  TAngle<double>& angle, bool check = false) const;
330 
333  void compute()
334  throw(Exception::GeneralException,
335  Exception::DivisionByZero,
336  Exception::IndexOverflow);
337 
339 
340  private:
341 
342  /*_ Test whether a ReducedSurface object can be copied.
343  */
344  bool canBeCopied(const ReducedSurface& reduced_surface);
345 
346  /*_ Copy a ReducedSurface object.
347  */
348  void copy(const ReducedSurface& reduced_surface);
349 
350  /*_
351  */
352  void correctEdges(RSFace* face1, RSFace* face2,
353  RSEdge* edge1, RSEdge* edge2);
354 
355  /*_
356  */
357  void joinVertices(RSFace* face1, RSFace* face2,
358  RSVertex* vertex1, RSVertex* vertex2);
359 
360  /*_
361  */
362  void findSimilarVertices(RSFace* face1, RSFace* face2,
363  std::vector<RSVertex*>& rsvertex1,
364  std::vector<RSVertex*>& rsvertex2);
365 
366  /*_
367  */
368  void findSimilarEdges(RSFace* face1, RSFace* face2,
369  std::vector<RSEdge*>& rsedge1,
370  std::vector<RSEdge*>& rsedge2);
371 
372  protected:
373 
374  /*_ the number of atoms of the reduced surface
375  */
376  Size number_of_atoms_;
377 
378  /*_ the atoms of the molecule
379  */
380 
381  std::vector< TSphere3<double> > atom_;
382 
383  /*_ probe radius
384  */
385  double probe_radius_;
386 
387  /*_ the number of vertices of the reduced surface
388  */
389  Size number_of_vertices_;
390 
391  /*_ the vertices of the reduced surface
392  */
393  std::vector< RSVertex* > vertices_;
394 
395  /*_ the number of edges of the reduced surface
396  */
397  Size number_of_edges_;
398 
399  /*_ the edges of the reduced surface
400  */
401  std::vector< RSEdge* > edges_;
402 
403  /*_ the number of faces of the reduced surface
404  */
405  Size number_of_faces_;
406 
407  /*_ the faces of the reduced surface
408  */
409  std::vector< RSFace* > faces_;
410 
411  /*_ maximal radius of all atoms
412  */
413  double r_max_;
414 
415  /*_ bounding SimpleBox of the atom centers of the molecule
416  */
417  TSimpleBox3<double> bounding_box_;
418  };
419 
423 
427  BALL_EXPORT std::ostream& operator << (std::ostream& s, const ReducedSurface& rs);
428 
430 
435  {
436  public:
437 
438  BALL_CREATE(RSComputer)
439 
440 
443 
450  {
451  STATUS_OK = 0,
453  STATUS_NOT_TESTED
454  };
455 
462  {
463  STATUS_ON_SURFACE = 0,
465  STATUS_UNKNOWN
466  };
468 
470  {
471  ProbeStatus status[2];
473  };
474 
478 
483  RSComputer();
484 
488 
491  virtual ~RSComputer();
492 
494 
497 
500  void run()
501  throw(Exception::GeneralException,
502  Exception::DivisionByZero,
503  Exception::IndexOverflow);
504 
506 
507 
508  private:
509 
510  /*_ @name Computing reduced surface
511  */
513 
514  /*_
515  */
516  void preProcessing();
517 
518  /*_ Compute a RSComponent.
519  */
520  void getRSComponent()
521  throw(Exception::GeneralException,
522  Exception::DivisionByZero,
523  Exception::IndexOverflow);
524 
525  /*_ Treat all edges of a face.
526  @param face the RSFace to be treated
527  */
528  bool treatFace(RSFace* face)
529  throw(Exception::GeneralException,
530  Exception::DivisionByZero,
531  Exception::IndexOverflow);
532 
533  /*_ Roll over an edge that belongs to only one face and find the other one.
534  @param edge the RSEdge to be treated
535  */
536  bool treatEdge(RSEdge* edge)
537  throw(Exception::GeneralException,
538  Exception::DivisionByZero,
539  Exception::IndexOverflow);
540 
541  /*_ Treat an ambiguous situation.
542  All vertices on an ambiguous atom are deleted with all its edges and
543  faces. The radius of the atom is decreased by 10 EPSILON.
544  @param atom the index of the atom
545  */
546  void correct(Index atom);
547 
548  /*_ Check all new created vertices for extensions
549  */
550  void extendComponent()
551  throw(Exception::GeneralException,
552  Exception::DivisionByZero,
553  Exception::IndexOverflow);
554 
555  /*_ Find a third atom rolling over two vertices starting on a face.
556  From all atoms which can be touched by the probe sphere when it
557  touches the given two vertices we choose the one with smallest
558  rotation angle.
559  If the rotation angle equals zero, the probe sphere can touch four
560  atoms and an exception is thrown.
561  If no atom can be found an exception is thrown.
562  @param vertex1 the first vertex
563  @param vertex2 the second vertex
564  @param face the starting face
565  @param probe the new probe sphere
566  @param phi the rotation angle
567  @return Index index of the found atom
568  */
569  Index thirdAtom(RSVertex* vertex1, RSVertex* vertex2,
570  RSFace* face, TSphere3<double>& probe, TAngle<double>& phi)
571  throw(Exception::GeneralException,
572  Exception::DivisionByZero,
573  Exception::IndexOverflow);
574 
576  /*_ @name Finding a start position
577  */
579 
580  /*_ Find a start position
581  @param vertex a pointer to the found vertex, if only a vertex
582  can be found
583  @param edge a pointer to the found edge, if only an edge can be
584  found
585  @param face a pointer to the found face, if a face can be found
586  @return Position 0, if no start position is found,
587  1, if a single vertex is found,
588  2, if an edge is found,
589  3, if a face is found
590  */
591  Position getStartPosition()
592  throw(Exception::DivisionByZero);
593 
595  /*_ @name Finding a first face
596  */
598 
599  /*_ Try to find a starting face
600  @return RSFace* a pointer to the found face, if a face can be found,
601  NULL otherwise
602  */
603  RSFace* findFirstFace()
604  throw(Exception::DivisionByZero);
605 
606  /*_ Try to find a starting face in a given direction
607  @param direction search in x-direction, if direction is 0,
608  search in y-direction, if direction is 1,
609  search in z-direction, if direction is 2
610  @param extrem search in min direction, if extrem is 0,
611  search in max direction, if extrem is 1
612  @return RSFace* a pointer to the found face, if a face can be found,
613  NULL otherwise
614  */
615  RSFace* findFace(Position direction, Position extrem)
616  throw(Exception::DivisionByZero);
617 
619  /*_ @name Finding a first edge
620  */
622 
623  /*_ Try to find a starting edge
624  @return RSEdge* a pointer to the found edge, if a face can be found,
625  NULL otherwise
626  */
627  RSEdge* findFirstEdge();
628 
629  /*_ Try to find a starting edge in a given direction
630  @param direction search in x-direction, if direction is 0,
631  search in y-direction, if direction is 1,
632  search in z-direction, if direction is 2
633  @param extrem search in min direction, if extrem is 0,
634  search in max direction, if extrem is 1
635  @return RSEdge* a pointer to the found edge, if a face can be found,
636  NULL otherwise
637  */
638  RSEdge* findEdge(Position direction, Position extrem);
639 
641  /*_ @name Finding a first vertex
642  */
644 
645  /*_ Try to find a single atom
646  @return RSVertex* a pointer to the found vertex, if a vertex can be
647  found, NULL otherwise
648  */
649  RSVertex* findFirstVertex();
650 
651  /*_ Find a single atom in a given direction
652  @param direction search in x-direction, if direction is 0,
653  search in y-direction, if direction is 1,
654  search in z-direction, if direction is 2
655  @param extrem search in min direction, if extrem is 0,
656  search in max direction, if extrem is 1
657  @return Index the index of the found atom
658  */
659  Index findFirstAtom(Position direction, Position extrem);
660 
661  /*_ Find a second atom close enougth to the first atom in a given direction
662  @param atom1 the index of the first atom
663  @param direction search in x-direction, if direction is 0,
664  search in y-direction, if direction is 1,
665  search in z-direction, if direction is 2
666  @param extrem search in min direction, if extrem is 0,
667  search in max direction, if extrem is 1
668  @return Index the index of the found atom
669  */
670  Index findSecondAtom(Index atom, Position direction, Position extrem);
671 
672  /*_ Find a second atom close enougth to the first two atoms
673  @param atom1 the index of the first atom
674  @param atom2 the index of the second atom
675  @param atom_list a HashSet of the indices of all candidate atoms
676  @return ::std::list< ::std::pair< Index,TSphere3<double> > >
677  a list of all candidates with their probe spheres
678  */
679  void findThirdAtom(Index atom1, Index atom2, const std::deque<Index>& third,
680  std::deque< std::pair< Index,TSphere3<double> > >& atoms);
681 
683  /*_ @name Some utilities
684  */
686 
687  /*_ Find all atoms close enougth to two given atoms.
688  The indices of all atoms which can be touched by the probe sphere when
689  it touches the given atoms are computed.
690  @param atom1 the index of the first given atom
691  @param atom2 the index of the second given atom
692  @param output_list list of all atoms close enougth to the given atoms
693  */
694  const std::deque<Index>& neighboursOfTwoAtoms(const SortedPosition2& pos);
695 
696  /*_ Find all atoms close enougth to three given atoms.
697  The indices of all atoms which can be touched by the probe sphere when
698  it touches the given atoms are computed.
699  @param atom1 the index of the first given atom
700  @param atom2 the index of the second given atom
701  @param atom3 the index of the third given atom
702  @param output_list list of all atoms close enougth to the given atoms
703  */
704  void neighboursOfThreeAtoms(Index atom1, Index atom2, Index atom3,
705  std::deque<Index>& output_list);
706 
707  /*_ Get the extrem coordinate of a circle in a given direction
708  @param circle the circle
709  @param direction search in x-direction, if direction is 0,
710  search in y-direction, if direction is 1,
711  search in z-direction, if direction is 2
712  @param extrem search in min direction, if extrem is 0,
713  search in max direction, if extrem is 1
714  @return double the extrem coordinate
715  */
716  double getCircleExtremum(const TCircle3<double>& circle,
717  Position direction, Position extrem);
718 
720  /*_ @name Creating / updating edges / faces
721  */
723 
724  /*_ Create a free edge from two vertices if it is a free edge
725  @param vertex1 a pointer to the first vertex
726  @param vertex2 a pointer to the second vertex
727  @return RSEdge* a pointer to the created free edge, if there is one,
728  NULL otherwise
729  */
730  RSEdge* createFreeEdge(RSVertex* vertex1, RSVertex* vertex2);
731 
732  /*_ Get the circle described by the center of the probe sphere and the two
733  contact circles with the atoms when the probe sphere rolls over two
734  atoms
735  @param atom1 the index of the first atom
736  @param atom2 the index of the second atom
737  @param circle1 the circle described by the center of the probe sphere
738  @param circle2 the contact circle with atom1
739  @param circle3 the contact circle with atom2
740  @return bool true, if the probe sphere can touch both atoms,
741  false, otherwise
742  */
743  bool getCircles(Index atom1, Index atom2, TCircle3<double>& circle1,
744  TCircle3<double>& circle2, TCircle3<double>& circle3);
745 
746  /*_ Get the normal vector of the face described by three atoms and a probe
747  @param atom1 the index of the first atom
748  @param atom2 the index of the second atom
749  @param atom3 the index of the third atom
750  @param probe the probe sphere lying on atom1, atom2, atom3
751  @return TVector3<double> the normal vector
752  */
753  TVector3<double> getFaceNormal(const TSphere3<double>& atom1, const TSphere3<double>& atom2,
754  const TSphere3<double>& atom3, const TSphere3<double>& probe);
755 
756  /*_ Update a face and it's edges
757  @param v1 the first vertex of the face
758  @param v2 the second vertex of the face
759  @param v3 the third vertex of the face
760  @param e1 the first edge
761  @param e2 the second edge
762  @param e3 the third edge
763  @param f the face
764  @param probe the probe sphere of the face
765  */
766  void updateFaceAndEdges(RSVertex* v1, RSVertex* v2, RSVertex* v3,
767  RSEdge* e1, RSEdge* e2, RSEdge* e3,
768  RSFace* f, const TSphere3<double>& probe);
769 
770  /*_ Test, whether a face exists or not
771  @param face a pointer to the face to be tested
772  @return RSFace* a pointer to the face, if it exists, otherwise NULL
773  */
774  RSFace* faceExists(RSFace* face, const std::list< RSVertex* >& vertices);
775 
777  /*_ @name Finding a probe sphere
778  */
780 
781  /*_ Get the centers of the probe sphere when it lies on three atoms
782  @param pos the three atom's indices
783  @param c1 the first center
784  @param c2 the second center
785  @return bool true, if the probe sphere can touch the three atoms,
786  false, otherwise
787  */
788  bool centerOfProbe(const SortedPosition3& pos, TVector3<double>& c1, TVector3<double>& c2);
789 
790  /*_ Check,weather a probe sphere is inside an atom
791  @param probe the probe sphere to be tested
792  @return bool true, if the probe sphere is not intersecting any atom
793  false, otherwise
794  */
795  bool checkProbe(const TSphere3<double>& probe, const SortedPosition3& pos);
796 
797  /*_
798  */
799  void correctProbePosition(Position atom);
800 
801  /*_
802  */
803  void correctProbePosition(const SortedPosition3& pos);
804 
805  /*_
806  */
807  void insert(RSVertex* vertex);
808 
809  /*_
810  */
811  void insert(RSEdge* edge);
812 
813  /*_
814  */
815  void insert(RSFace* face);
816 
818 
819  protected:
820 
821 
822  /*_ a pointer to the reduced surface to compute
823  */
825 
826  /*_ for each atom a list of its neighbours
827  */
828  std::vector< std::deque<Index> > neighbours_;
829 
830  /*_ for each atom a status
831  */
832  std::vector< AtomStatus > atom_status_;
833 
834  /*_ for each pair of atoms a list of its neighbours
835  */
836  HashMap< SortedPosition2, std::deque<Index> > neighbours_of_two_;
837 
838  /*_ for each triple of atoms its probe positions
839  */
840  HashMap< SortedPosition3, ProbePosition* > probe_positions_;
841 
842  /*_ all new created vertices which are not yet checked for extensions
843  */
844  HashSet<RSVertex*> new_vertices_;
845 
846  /*_ all new created faces which are not completely treated yet
847  */
848  HashSet<RSFace*> new_faces_;
849 
850  /*_ for each atom a list of the rsvertices of the atom
851  */
852  std::vector< std::list<RSVertex*> > vertices_;
853  };
854 } // namespace BALL
855 
856 #endif // BALL_STRUCTURE_REDUCEDSURFACE_H