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 if (quadric!=NULL) 00033 deleteVector<double>(quadric); 00034 } 00035 }; 00036 00037 00038 // voronoi facet cell or del edge 00039 class Del1Cell : public MixedCell 00040 { 00041 public: 00042 int ids[2]; // the two atoms generating the del edge 00043 vector<double*> planes; // lateral planes 00044 00045 virtual ~Del1Cell() 00046 { 00047 for (unsigned int i=0;i<planes.size();i++) 00048 free(planes[i]); 00049 } 00050 00051 // #IMPROVE one could build pointers from here to the planes computed in Del2Cell or viceversa 00052 double upper[4]; // upper plane 00053 double lower[4]; // upper plane 00054 }; 00055 00056 // reduced voronoi cell 00057 class Del0Cell : public MixedCell 00058 { 00059 public: 00060 int id; // atom index 00061 vector<double*> planes; // planes pointer that points to upper/lower of Del1Cell 00062 00063 virtual ~Del0Cell() 00064 { 00065 // do nothing, planes are not managed by this class 00066 } 00067 }; 00068 00069 // del facet cell 00070 class Del2Cell : public MixedCell 00071 { 00072 public: 00073 int ids[3]; // atom indices 00074 double planes[3][4]; // lateral planes 00075 double upper[4]; 00076 double lower[4]; 00077 00078 virtual ~Del2Cell() 00079 { 00080 // do nothing, no pointers 00081 } 00082 }; 00083 00084 // reduced thethraedron 00085 class Del3Cell: public MixedCell 00086 { 00087 public: 00088 int ids[4]; // atom indices 00089 // planes of the reduced tethraedron 00090 double planes[4][4]; 00091 double vor[3]; // Voronoi center 00092 double reduced[4][3]; // reduced points in the same order of atom indices 00093 00094 virtual ~Del3Cell() 00095 { 00096 // do nothing, no pointers 00097 } 00098 }; 00099 00100 00101 #define GRID_MIXEDCELLMAP_2D(i,j,k,NA,NB) gridMixedCellMap2D[(k)+(MAX_MIXEDCELLS_2D-1)*((j)+(NB)*(i))] 00102 // fortran like 00103 //#define GRIDMIXEDCELLMAP(i,j,k,l,NX,NY,NZ) gridMixedCellMap[(i)+(NX)*((j)+(NY)*((k)+(NZ)*(l)))] 00104 // c like (NX is ignored) 00105 #define GRIDMIXEDCELLMAP(i,j,k,l,NX,NY,NZ) gridMixedCellMap[(l)+(MAX_MIXEDCELLS-1)*((k)+(NZ)*((j)+(NY)*(i)))] 00106 00107 #ifdef ENABLE_CGAL 00108 00109 #include <CGAL/Exact_predicates_inexact_constructions_kernel.h> 00110 #include <CGAL/Regular_triangulation_3.h> 00111 #include <CGAL/Regular_triangulation_euclidean_traits_3.h> 00112 #include <CGAL/Regular_triangulation_filtered_traits_3.h> 00113 #include <CGAL/Triangulation_cell_base_with_info_3.h> 00114 #include <CGAL/Regular_triangulation_cell_base_3.h> 00115 #include <CGAL/Triangulation_vertex_base_with_info_3.h> 00116 #include <CGAL/Triangulation_data_structure_3.h> 00117 #include <cassert> 00118 #include <vector> 00119 #include <fstream> 00120 00121 // use to translate in Pov-Ray 00122 #include <CGAL/Polyhedron_3.h> 00123 #include <CGAL/convex_hull_3.h> 00124 00126 #endif 00127 00128 #define DELAUNAY_POINT_CELL 0 00129 #define DELAUNAY_EDGE_CELL 1 00130 #define DELAUNAY_FACET_CELL 2 00131 #define DELAUNAY_TETRA_CELL 3 00132 00141 class SkinSurface: public Surface 00142 { 00143 00144 #ifdef ENABLE_CGAL 00145 private: 00146 00147 typedef CGAL::Exact_predicates_inexact_constructions_kernel K; 00148 typedef CGAL::Regular_triangulation_filtered_traits_3<K> Traits; 00149 typedef CGAL::Triangulation_cell_base_with_info_3< Del3Cell*,Traits> Cell_with_info; 00150 typedef CGAL::Regular_triangulation_cell_base_3<Traits,Cell_with_info> RT_Cell_with_info; 00151 typedef CGAL::Triangulation_vertex_base_3<Traits> Vertex; 00152 typedef CGAL::Triangulation_data_structure_3<Vertex,RT_Cell_with_info> tds; 00153 typedef Traits::RT Weight; 00154 typedef Traits::Bare_point Point3; 00155 typedef Traits::Weighted_point Weighted_point; 00156 typedef CGAL::Regular_triangulation_3<Traits,tds> Rt; 00157 typedef Rt::Vertex_iterator Vertex_iterator; 00158 typedef Rt::Finite_vertices_iterator Finite_Vertex_Iterator; 00159 typedef Rt::Finite_cells_iterator Finite_Cells_Iterator; 00160 typedef Rt::Finite_edges_iterator Finite_Edges_Iterator; 00161 typedef Rt::Finite_facets_iterator Finite_Facets_Iterator; 00162 typedef Rt::Vertex_handle Vertex_handle; 00163 typedef tds::Cell_circulator Cell_circulator; 00164 typedef tds::Cell_handle Cell_handle; 00165 typedef Rt::Facet Facet; 00166 00167 // used to translate in Pov-Ray 00168 typedef CGAL::Polyhedron_3<K> Polyhedron; 00169 typedef Polyhedron::Halfedge_around_facet_circulator HF_circulator; 00170 00171 //a functor computing the plane containing a triangular facet 00172 struct Plane_from_facet { 00173 Polyhedron::Plane_3 operator()(Polyhedron::Facet& f) { 00174 Polyhedron::Halfedge_handle h = f.halfedge(); 00175 return Polyhedron::Plane_3( h->vertex()->point(), 00176 h->next()->vertex()->point(), 00177 h->opposite()->vertex()->point()); 00178 } 00179 }; 00180 00181 00182 #endif 00183 00184 private: 00186 int numMixedCells; 00192 int type[4]; 00194 int* gridMixedCellMap; 00196 int nx,ny,nz; 00197 double scale; 00198 double side; 00199 double*x,*y,*z; 00200 double s; 00201 double xmin,xmax,ymin,ymax,zmin,zmax; 00202 short*** ind; 00203 int* gridMixedCellMap2D; 00204 00206 int nx_2d,ny_2d,nz_2d; 00207 double scale_2d; 00208 double side_2d; 00209 double xmin_2d,xmax_2d,ymin_2d,ymax_2d,zmin_2d,zmax_2d; 00210 int** ind_2d; 00211 00213 MixedCell** mixedComplex; 00216 bool buildSkin(); 00218 bool buildAuxiliaryGrid(); 00219 #ifdef ENABLE_CGAL 00220 00221 bool buildSkinCGAL(); 00222 #endif 00223 00224 bool savePovRay; 00225 00226 public: 00228 SkinSurface(); 00230 SkinSurface(DelPhiShared* ds); 00231 00233 00234 virtual bool build(); 00236 virtual bool save(char* fileName); 00238 virtual bool load(char* fileName); 00240 virtual void printSummary(); 00242 virtual bool getProjection(double p[3],double* proj1,double* proj2, 00243 double* proj3,double* normal1,double* normal2,double* normal3); 00249 virtual void getRayIntersection(double p1[3],double p2[3],vector<pair<double,double*> >& intersections,int thdID); 00251 virtual void preProcessPanel(); 00253 00254 void setShrinking(double ss) 00255 { 00256 if (ss<=1.0 && ss>0.0) 00257 { 00258 s=ss; 00259 } 00260 else 00261 { 00262 cout << endl << WARN << "Cannot set " << ss << ". s parameter is in (0,1]. Setting "<< DEFAULT_S; 00263 s = DEFAULT_S; 00264 } 00265 } 00266 00267 double getShrinking() 00268 { 00269 return s; 00270 } 00271 00272 void setSavePovRay(bool ss) 00273 { 00274 savePovRay = ss; 00275 } 00276 00277 bool getSavePovRay() 00278 { 00279 return savePovRay; 00280 } 00281 00282 virtual ~SkinSurface(); 00283 00284 private: 00287 //bool rayQuadricIntersection(double*,double*,double**,double*,double*,int thdID,double* cache); 00288 bool rayQuadricIntersection(double*,double*,double*,double*,double*,int thdID,double* cache); 00290 bool isFeasible(MixedCell* mc,double* point); 00292 void projectToQuadric(double* y,double* Q,int type,double* proj,double* norm,double& dist); 00294 void saveSkinPatch(ofstream& of,MixedCell* mc,int index,vector<Weighted_point>& la); 00295 }; 00296 00297 #endif