NanoShaper
0.3.1
NanoShaper is a tool able to triangulate and inspect an arbitray triangulated surface or several types of molecular surfaces such as the Gaussian, Skin and the SES
|
00001 //--------------------------------------------------------- 00005 //--------------------------------------------------------- 00006 00007 #ifndef SkinSurface_h 00008 #define SkinSurface_h 00009 00010 #include "Surface.h" 00011 00012 //#define DEBUG_SKIN 00013 00014 #define AUX_GRID_DIM_SKIN 100 00015 #define MAX_MIXEDCELLS 400 00016 #define DEFAULT_S 0.45 // default shrink factor 00017 00018 // aggressive settings can be used in a 64 bit machine with 6 GB of memory 00019 00020 #define AUX_GRID_DIM_SKIN_2D 50 // aggressive setting is 150 00021 #define MAX_MIXEDCELLS_2D (400*AUX_GRID_DIM_SKIN_2D) // aggressive setting is (200*AUX_GRID_DIM_SKIN_2D) 00022 00024 class MixedCell{ 00025 public: 00026 double* quadric; 00027 int surface_type; 00028 vector<double*> mc_points; // these are only pointer, no explicit denstructor 00029 00030 virtual ~MixedCell() 00031 { 00032 deleteVector<double>(quadric); 00033 } 00034 }; 00035 00036 00037 // voronoi facet cell or del edge 00038 class Del1Cell : public MixedCell 00039 { 00040 public: 00041 int ids[2]; // the two atoms generating the del edge 00042 vector<double*> planes; // lateral planes 00043 00044 virtual ~Del1Cell() 00045 { 00046 for (unsigned int i=0;i<planes.size();i++) 00047 free(planes[i]); 00048 } 00049 00050 // #IMPROVE one could build pointers from here to the planes computed in Del2Cell or viceversa 00051 double upper[4]; // upper plane 00052 double lower[4]; // upper plane 00053 }; 00054 00055 // reduced voronoi cell 00056 class Del0Cell : public MixedCell 00057 { 00058 public: 00059 int id; // atom index 00060 vector<double*> planes; // planes pointer that points to upper/lower of Del1Cell 00061 00062 virtual ~Del0Cell() 00063 { 00064 // do nothing, planes are not managed by this class 00065 } 00066 }; 00067 00068 // del facet cell 00069 class Del2Cell : public MixedCell 00070 { 00071 public: 00072 int ids[3]; // atom indices 00073 double planes[3][4]; // lateral planes 00074 double upper[4]; 00075 double lower[4]; 00076 00077 virtual ~Del2Cell() 00078 { 00079 // do nothing, no pointers 00080 } 00081 }; 00082 00083 // reduced thethraedron 00084 class Del3Cell: public MixedCell 00085 { 00086 public: 00087 int ids[4]; // atom indices 00088 // planes of the reduced tethraedron 00089 double planes[4][4]; 00090 double vor[3]; // Voronoi center 00091 double reduced[4][3]; // reduced points in the same order of atom indices 00092 00093 virtual ~Del3Cell() 00094 { 00095 // do nothing, no pointers 00096 } 00097 }; 00098 00099 00100 #define GRID_MIXEDCELLMAP_2D(i,j,k,NA,NB) gridMixedCellMap2D[(k)+(MAX_MIXEDCELLS_2D-1)*((j)+(NB)*(i))] 00101 // fortran like 00102 //#define GRIDMIXEDCELLMAP(i,j,k,l,NX,NY,NZ) gridMixedCellMap[(i)+(NX)*((j)+(NY)*((k)+(NZ)*(l)))] 00103 // c like (NX is ignored) 00104 #define GRIDMIXEDCELLMAP(i,j,k,l,NX,NY,NZ) gridMixedCellMap[(l)+(MAX_MIXEDCELLS-1)*((k)+(NZ)*((j)+(NY)*(i)))] 00105 00106 #ifdef ENABLE_CGAL 00107 00108 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h> 00109 #include <CGAL/Regular_triangulation_3.h> 00110 #include <CGAL/Regular_triangulation_euclidean_traits_3.h> 00111 #include <CGAL/Regular_triangulation_filtered_traits_3.h> 00112 #include <CGAL/Triangulation_cell_base_with_info_3.h> 00113 #include <CGAL/Regular_triangulation_cell_base_3.h> 00114 #include <CGAL/Triangulation_vertex_base_with_info_3.h> 00115 #include <CGAL/Triangulation_data_structure_3.h> 00116 #include <CGAL/Polyhedron_3.h> 00117 #include <CGAL/convex_hull_3.h> 00118 #include <CGAL/Vector_3.h> 00119 #include <CGAL/Plane_3.h> 00120 #include <CGAL/Ray_3.h> 00121 #include <CGAL/Line_3.h> 00122 #include <CGAL/intersections.h> 00123 #include <CGAL/Linear_algebraCd.h> 00124 #include <cassert> 00125 #include <vector> 00126 #include <fstream> 00127 00129 #endif 00130 00131 #define DELAUNAY_POINT_CELL 0 00132 #define DELAUNAY_EDGE_CELL 1 00133 #define DELAUNAY_FACET_CELL 2 00134 #define DELAUNAY_TETRA_CELL 3 00135 00144 class SkinSurface: public Surface 00145 { 00146 00147 #ifdef ENABLE_CGAL 00148 private: 00149 00150 typedef CGAL::Exact_predicates_inexact_constructions_kernel K; 00151 typedef CGAL::Regular_triangulation_filtered_traits_3<K> Traits; 00152 00153 typedef CGAL::Triangulation_cell_base_with_info_3< Del3Cell*,Traits> Cell_with_info; 00154 typedef CGAL::Regular_triangulation_cell_base_3<Traits,Cell_with_info> RT_Cell_with_info; 00155 typedef CGAL::Triangulation_vertex_base_3<Traits> Vertex; 00156 typedef CGAL::Triangulation_data_structure_3<Vertex,RT_Cell_with_info> tds; 00157 typedef Traits::RT Weight; 00158 typedef Traits::Bare_point Point3; 00159 typedef Traits::Weighted_point Weighted_point; 00160 typedef CGAL::Polyhedron_3<K> Polyhedron; 00161 typedef CGAL::Regular_triangulation_3<Traits,tds> Rt; 00162 typedef Rt::Vertex_iterator Vertex_iterator; 00163 typedef Rt::Finite_vertices_iterator Finite_Vertex_Iterator; 00164 typedef Rt::Finite_cells_iterator Finite_Cells_Iterator; 00165 typedef Rt::Finite_edges_iterator Finite_Edges_Iterator; 00166 typedef Rt::Finite_facets_iterator Finite_Facets_Iterator; 00167 typedef Rt::Vertex_handle Vertex_handle; 00168 typedef tds::Cell_circulator Cell_circulator; 00169 typedef tds::Cell_handle Cell_handle; 00170 typedef Polyhedron::Halfedge_around_facet_circulator HF_circulator; 00171 typedef Rt::Facet Facet; 00172 typedef CGAL::Linear_algebraCd<double> LA; 00173 00174 #endif 00175 00176 private: 00178 int numMixedCells; 00184 int type[4]; 00186 int* gridMixedCellMap; 00188 int nx,ny,nz; 00189 double scale; 00190 double side; 00191 double*x,*y,*z; 00192 double s; 00193 double xmin,xmax,ymin,ymax,zmin,zmax; 00194 short*** ind; 00195 int* gridMixedCellMap2D; 00196 00198 int nx_2d,ny_2d,nz_2d; 00199 double scale_2d; 00200 double side_2d; 00201 double xmin_2d,xmax_2d,ymin_2d,ymax_2d,zmin_2d,zmax_2d; 00202 int** ind_2d; 00203 00205 MixedCell** mixedComplex; 00208 bool buildSkin(); 00210 bool buildAuxiliaryGrid(); 00211 #ifdef ENABLE_CGAL 00212 00213 bool buildSkinCGAL(); 00214 #endif 00215 00216 bool savePovRay; 00217 00218 public: 00220 SkinSurface(); 00222 SkinSurface(DelPhiShared* ds); 00223 00225 00226 virtual bool build(); 00228 virtual bool save(char* fileName); 00230 virtual bool load(char* fileName); 00232 virtual void printSummary(); 00234 virtual bool getProjection(double p[3],double* proj1,double* proj2, 00235 double* proj3,double* normal1,double* normal2,double* normal3); 00241 virtual void getRayIntersection(double p1[3],double p2[3],vector<pair<double,double*> >& intersections,int thdID); 00243 virtual void preProcessPanel(); 00245 00246 void setShrinking(double ss) 00247 { 00248 if (ss<=1.0 && ss>0.0) 00249 { 00250 s=ss; 00251 } 00252 else 00253 { 00254 cout << endl << WARN << "Cannot set " << ss << ". s parameter is in (0,1]. Setting "<< DEFAULT_S; 00255 s = DEFAULT_S; 00256 } 00257 } 00258 00259 double getShrinking() 00260 { 00261 return s; 00262 } 00263 00264 void setSavePovRay(bool ss) 00265 { 00266 savePovRay = ss; 00267 } 00268 00269 bool getSavePovRay() 00270 { 00271 return savePovRay; 00272 } 00273 00274 virtual ~SkinSurface(); 00275 00276 private: 00279 //bool rayQuadricIntersection(double*,double*,double**,double*,double*,int thdID,double* cache); 00280 bool rayQuadricIntersection(double*,double*,double*,double*,double*,int thdID,double* cache); 00282 bool isFeasible(MixedCell* mc,double* point); 00284 void projectToQuadric(double* y,double* Q,int type,double* proj,double* norm,double& dist); 00285 }; 00286 00287 #endif