From 944ca45ab7fb8395a2d5c268884446aee2084677 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Vladim=C3=ADr=20Ulman?= <ulman@mpi-cbg.de> Date: Thu, 6 Oct 2016 14:53:35 +0200 Subject: [PATCH] Added initial (just contour yet) rendering of meshes into mask images. --- CMakeLists.txt | 13 ++++ cmath3d/TriangleMesh.cpp | 164 +++++++++++++++++++++++++++++++++++++++ cmath3d/TriangleMesh.h | 3 + src/main.cpp | 13 +++- src/params.h | 3 + 5 files changed, 193 insertions(+), 3 deletions(-) diff --git a/CMakeLists.txt b/CMakeLists.txt index 8c76f7c..ed99db4 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -55,6 +55,12 @@ else (WIN32) set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -fopenmp") endif (WIN32) +# +# due to Lapack+i3dalgo +# + +SET(CMAKE_EXE_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS} -Xlinker -defsym -Xlinker MAIN__=main") + # # user defined options # @@ -102,6 +108,13 @@ include_directories(${INC_GSL}) find_library(LIB_LAPACK lapack) set(LIBS ${LIBS} ${LIB_LAPACK}) +# i3dlibs libraries +find_path(INC_I3D "i3d/image3d.h") +include_directories(${INC_I3D}) +find_library(LIB_I3DCORE i3dcore) +find_library(LIB_I3DALGO i3dalgo) +set(LIBS ${LIBS} ${LIB_I3DCORE} ${LIB_I3DALGO}) + #-------------- # SOURCE FILES #-------------- diff --git a/cmath3d/TriangleMesh.cpp b/cmath3d/TriangleMesh.cpp index 179df5a..229f7ed 100644 --- a/cmath3d/TriangleMesh.cpp +++ b/cmath3d/TriangleMesh.cpp @@ -190,6 +190,170 @@ int ActiveMesh::ImportVTK(const char *filename) } +void ActiveMesh::RenderMask(i3d::Image3d<i3d::GRAY16>& mask) +{ + //time savers: resolution + const float xRes=mask.GetResolution().GetX(); + const float yRes=mask.GetResolution().GetY(); + const float zRes=mask.GetResolution().GetZ(); + + //time savers: offset + const float xOff=mask.GetOffset().x; + const float yOff=mask.GetOffset().y; + const float zOff=mask.GetOffset().z; + + //over all triangles + for (unsigned int i=0; i < ID.size()/3; ++i) + //for (unsigned int i=0; i < 15; ++i) + { + const Vector3F& v1=Pos[ID[3*i+0]]; + const Vector3F& v2=Pos[ID[3*i+1]]; + const Vector3F& v3=Pos[ID[3*i+2]]; + + //sweeping (and rendering) the triangle + for (float c=0.f; c <= 1.0f; c += 0.1f) + for (float b=0.f; b <= (1.0f-c); b += 0.1f) + { + float a=1.0f -b -c; + + Vector3F v=a*v1; + v+=b*v2; + v+=c*v3; + + //std::cout << "ID #" << i << ": v=(" << v.x << "," << v.y << "," << v.z << ")\n"; + + //nearest neighbor pixel coordinate + const int x=(int)roundf( (v.x-xOff) *xRes); + const int y=(int)roundf( (v.y-yOff) *yRes); + const int z=(int)roundf( (v.z-zOff) *zRes); + + if (mask.Include(x,y,z)) mask.SetVoxel(x,y,z,i); + } + } +} + + +void ActiveMesh::RenderMaskB(i3d::Image3d<i3d::GRAY16>& mask) +{ + //time savers: resolution + const float xRes=mask.GetResolution().GetX(); + const float yRes=mask.GetResolution().GetY(); + const float zRes=mask.GetResolution().GetZ(); + + //time savers: offset + const float xOff=mask.GetOffset().x; + const float yOff=mask.GetOffset().y; + const float zOff=mask.GetOffset().z; + + //over all triangles + //for (unsigned int i=0; i < ID.size()/3; ++i) + for (unsigned int i=0; i < 15; ++i) + { + const Vector3F& v1=Pos[ID[3*i+0]]; + const Vector3F& v2=Pos[ID[3*i+1]]; + const Vector3F& v3=Pos[ID[3*i+2]]; + + //Vector3F e12=v2-v1; //ID=0 + //Vector3F e13=v3-v1; //ID=1 + //Vector3F e23=v3-v2; //ID=2 + Vector3F edges[3]={v2-v1,v3-v1,v3-v2}; + + short order[3]={0,1,2}; + if (edges[1].LenQ() < edges[0].LenQ()) { order[0]=1; order[1]=0; } + if (edges[2].LenQ() < edges[order[1]].LenQ()) + { + order[2]=order[1]; + order[1]=2; + + if (edges[2].LenQ() < edges[order[0]].LenQ()) + { + order[1]=order[0]; + order[0]=2; + } + } + + //point within a triangle will of the form: vv + b*vb + c*vc + //vb is the shortest edge, vc is the longest edge + //vv is their common vertex + //vb,vc are oriented outward from this vertex + Vector3F vv,vb,vc; + if (order[0]==0 && order[2]==1) + { + //order[0] "points" at shortest edge + vv=v1; + vb=edges[0]; + vc=edges[1]; + } + else + if (order[0]==0 && order[2]==2) + { + vv=v2; + vb=-edges[0]; + vc=edges[2]; + } + else + if (order[0]==1 && order[2]==0) + { + vv=v1; + vb=edges[1]; + vc=edges[0]; + } + else + if (order[0]==1 && order[2]==2) + { + vv=v3; + vb=-edges[1]; + vc=-edges[2]; + } + else + if (order[0]==2 && order[2]==0) + { + vv=v2; + vb=edges[2]; + vc=-edges[0]; + } + else + if (order[0]==2 && order[2]==1) + { + vv=v3; + vb=-edges[2]; + vc=-edges[1]; + } + + //optimal increment for the short and long edge, respectively + float db=(vb.x*xRes > vb.y*yRes)? vb.x*xRes : vb.y*yRes; + db=(vb.z*zRes > db)? vb.z*zRes : db; + db=1.f/db; + + float dc=(vc.x*xRes > vc.y*yRes)? vc.x*xRes : vc.y*yRes; + dc=(vc.z*zRes > dc)? vc.z*zRes : dc; + dc=1.f/dc; + + std::cout << "\nID #" << i << ": v1=(" << v1.x << "," << v1.y << "," << v1.z << ")\n"; + std::cout << "ID #" << i << ": v2=(" << v2.x << "," << v2.y << "," << v2.z << ")\n"; + std::cout << "ID #" << i << ": v3=(" << v3.x << "," << v3.y << "," << v3.z << ")\n"; + + //sweeping (and rendering) the triangle + for (float c=0.f; c <= 1.0f; c += dc) + for (float b=0.f; b <= (1.0f-c); b += db) + { + Vector3F v=vv; + v+=b*vb; + v+=c*vc; + + std::cout << "ID #" << i << ": v=(" << v.x << "," << v.y << "," << v.z << ")\n"; + + //nearest neighbor pixel coordinate + const int x=(int)roundf( (v.x-xOff) *xRes); + const int y=(int)roundf( (v.y-yOff) *yRes); + const int z=(int)roundf( (v.z-zOff) *zRes); + + if (mask.Include(x,y,z)) mask.SetVoxel(x,y,z,i); + } + } +} + + void ActiveMesh::CenterMesh(const Vector3F& newCentre) { //calc geom. centre diff --git a/cmath3d/TriangleMesh.h b/cmath3d/TriangleMesh.h index 58eece6..0da81df 100644 --- a/cmath3d/TriangleMesh.h +++ b/cmath3d/TriangleMesh.h @@ -2,6 +2,7 @@ #define __TRIANGLEMESH__FOR_MT__ #include <vector> +#include <i3d/image3d.h> #include "../src/params.h" /** @@ -28,6 +29,8 @@ class ActiveMesh ///returns 0 on success, otherwise non-zero int ImportSTL(const char* filename); int ImportVTK(const char* filename); + void RenderMask(i3d::Image3d<i3d::GRAY16>& mask); + void RenderMaskB(i3d::Image3d<i3d::GRAY16>& mask); void CenterMesh(const Vector3F& newCentre); void ScaleMesh(const Vector3F& scale); diff --git a/src/main.cpp b/src/main.cpp index 46853cf..94d9981 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -1,6 +1,4 @@ -#include <list> -#include <string.h> -#include <stdlib.h> //for atoi() +#include <i3d/image3d.h> #include "params.h" #include "graphics.h" @@ -36,6 +34,15 @@ int main(void) mesh.CenterMesh(params.sceneCentre); //work with mesh + i3d::Image3d<i3d::GRAY16> mask; + mask.SetOffset(i3d::Offset(0.0,0.0,0.0)); + //mask.SetResolution(i3d::Resolution(20.0,20.0,20.0)); //for torus + //mask.MakeRoom(300,300,300); + mask.SetResolution(i3d::Resolution(2.0,2.0,2.0)); //for sample cell + mask.MakeRoom(600,600,600); + mask.GetVoxelData()=0; + mesh.RenderMask(mask); + mask.SaveImage("mesh.ics"); initializeGL(); loopGL(); diff --git a/src/params.h b/src/params.h index c42c48c..7caf761 100644 --- a/src/params.h +++ b/src/params.h @@ -71,6 +71,9 @@ class Vector3d { return( *this ); } + /// unary minus = negate + Vector3d<T> operator-() const { return Vector3d<T>(-x,-y,-z); } + Vector3d<T>& operator=(const T scalar) { this->x=scalar; -- GitLab