From e4667cee90ab54837ae9450fcd788539c0a4aa0f Mon Sep 17 00:00:00 2001 From: David Mayerich Date: Wed, 17 Apr 2019 15:45:31 -0500 Subject: [PATCH] first commit --- CMakeLists.txt | 69 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ FindGLEW.cmake | 126 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ FindGLFW.cmake | 68 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ FindSTIM.cmake | 27 +++++++++++++++++++++++++++ main.cpp | 490 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ tensorfield.cpp | 0 tensorfield.h | 122 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 7 files changed, 902 insertions(+), 0 deletions(-) create mode 100644 CMakeLists.txt create mode 100644 FindGLEW.cmake create mode 100644 FindGLFW.cmake create mode 100644 FindSTIM.cmake create mode 100644 main.cpp create mode 100644 tensorfield.cpp create mode 100644 tensorfield.h diff --git a/CMakeLists.txt b/CMakeLists.txt new file mode 100644 index 0000000..b631b17 --- /dev/null +++ b/CMakeLists.txt @@ -0,0 +1,69 @@ +#Specify the version being used aswell as the language +cmake_minimum_required(VERSION 3.12) + +#Name your project here +project(tensorview) + +#set the module directory +set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}") + +#default to release mode +if(NOT CMAKE_BUILD_TYPE) + set(CMAKE_BUILD_TYPE Release) +endif(NOT CMAKE_BUILD_TYPE) + +#build the executable in the binary directory on MS Visual Studio +if ( MSVC ) + SET( CMAKE_RUNTIME_OUTPUT_DIRECTORY_DEBUG "${OUTPUT_DIRECTORY}") + SET( CMAKE_RUNTIME_OUTPUT_DIRECTORY_RELEASE "${OUTPUT_DIRECTORY}") + SET( LIBRARY_OUTPUT_DIRECTORY_DEBUG "${OUTPUT_DIRECTORY}") + SET( LIBRARY_OUTPUT_DIRECTORY_RELEASE "${OUTPUT_DIRECTORY}") + add_definitions(-D_CRT_SECURE_NO_WARNINGS) + add_definitions(-D_SCL_SECURE_NO_WARNINGS) +endif ( MSVC ) +#MAYBE REMOVE----------------- +#set C++11 flags if using GCC +if( CMAKE_COMPILER_IS_GNUCC ) + SET( CMAKE_CXX_FLAGS "-std=c++11") + SET( CUDA_NVCC_FLAGS "-std=c++11") +endif( CMAKE_COMPILER_IS_GNUCC ) +#----------------------------- + +#find packages----------------------------------- + +#find the pthreads package +find_package(Threads) + +#find the X11 package +find_package(X11) + +find_package(GLFW REQUIRED) +if(WIN32) + find_package(GLEW REQUIRED) + include_directories(${GLEW_INCLUDE_DIR}) +endif(WIN32) + +#find and set up OpenGL +find_package(OpenGL REQUIRED) + +find_package(STIM REQUIRED) + +add_executable(tensorview + main.cpp + tensorfield.h +) + +#include include directories +include_directories(${GLFW_INCLUDE_DIRS} + ${OpenGL_INCLUDE_DIRS} + ${STIM_INCLUDE_DIRS} +) + +target_link_libraries(tensorview + ${GLFW_LIBRARIES} + ${OPENGL_gl_LIBRARY} + ${OPENGL_glu_LIBRARY} +) +if(WIN32) + target_link_libraries(tensorview ${GLEW_GLEW_LIBRARY}) +endif(WIN32) \ No newline at end of file diff --git a/FindGLEW.cmake b/FindGLEW.cmake new file mode 100644 index 0000000..6ae9a27 --- /dev/null +++ b/FindGLEW.cmake @@ -0,0 +1,126 @@ +# Copyright (c) 2012-2016 DreamWorks Animation LLC +# +# All rights reserved. This software is distributed under the +# Mozilla Public License 2.0 ( http://www.mozilla.org/MPL/2.0/ ) +# +# Redistributions of source code must retain the above copyright +# and license notice and the following restrictions and disclaimer. +# +# * Neither the name of DreamWorks Animation nor the names of +# its contributors may be used to endorse or promote products derived +# from this software without specific prior written permission. +# +# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS +# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT +# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR +# A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT +# OWNER OR CONTRIBUTORS BE LIABLE FOR ANY INDIRECT, INCIDENTAL, +# SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT +# LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, +# DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY +# THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT +# (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE +# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. +# IN NO EVENT SHALL THE COPYRIGHT HOLDERS' AND CONTRIBUTORS' AGGREGATE +# LIABILITY FOR ALL CLAIMS REGARDLESS OF THEIR BASIS EXCEED US$250.00. +# + +#-*-cmake-*- +# - Find GLEW +# +# Author : Nicholas Yue yue.nicholas@gmail.com +# +# This auxiliary CMake file helps in find the GLEW headers and libraries +# +# GLEW_FOUND set if Glew is found. +# GLEW_INCLUDE_DIR GLEW's include directory +# GLEW_GLEW_LIBRARY GLEW libraries +# GLEW_glewmx_LIBRARY GLEWmx libraries (Mulitple Rendering Context) + +FIND_PACKAGE ( PackageHandleStandardArgs ) + +FIND_PATH( GLEW_LOCATION include/GL/glew.h + "$ENV{GLEW_ROOT}" + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + +FIND_PACKAGE_HANDLE_STANDARD_ARGS ( GLEW + REQUIRED_VARS GLEW_LOCATION + ) + +IF ( GLEW_LOCATION ) + + SET( GLEW_INCLUDE_DIR "${GLEW_LOCATION}/include" CACHE STRING "GLEW include path") + + SET ( ORIGINAL_CMAKE_FIND_LIBRARY_SUFFIXES ${CMAKE_FIND_LIBRARY_SUFFIXES}) + IF (GLEW_USE_STATIC_LIBS) + IF (APPLE) + SET(CMAKE_FIND_LIBRARY_SUFFIXES ".a") + FIND_LIBRARY ( GLEW_LIBRARY_PATH GLEW PATHS ${GLEW_LOCATION}/lib + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + FIND_LIBRARY ( GLEWmx_LIBRARY_PATH GLEWmx PATHS ${GLEW_LOCATION}/lib + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + # MESSAGE ( "APPLE STATIC" ) + # MESSAGE ( "GLEW_LIBRARY_PATH = " ${GLEW_LIBRARY_PATH} ) + ELSEIF (WIN32) + # Link library + SET(CMAKE_FIND_LIBRARY_SUFFIXES ".lib") + FIND_LIBRARY ( GLEW_LIBRARY_PATH GLEW32S PATHS ${GLEW_LOCATION}/lib ) + FIND_LIBRARY ( GLEWmx_LIBRARY_PATH GLEW32MXS PATHS ${GLEW_LOCATION}/lib ) + ELSE (APPLE) + SET(CMAKE_FIND_LIBRARY_SUFFIXES ".a") + FIND_LIBRARY ( GLEW_LIBRARY_PATH GLEW PATHS ${GLEW_LOCATION}/lib + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + FIND_LIBRARY ( GLEWmx_LIBRARY_PATH GLEWmx PATHS ${GLEW_LOCATION}/lib + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + # MESSAGE ( "LINUX STATIC" ) + # MESSAGE ( "GLEW_LIBRARY_PATH = " ${GLEW_LIBRARY_PATH} ) + ENDIF (APPLE) + ELSE () + IF (APPLE) + SET(CMAKE_FIND_LIBRARY_SUFFIXES ".dylib") + FIND_LIBRARY ( GLEW_LIBRARY_PATH GLEW PATHS ${GLEW_LOCATION}/lib ) + FIND_LIBRARY ( GLEWmx_LIBRARY_PATH GLEWmx PATHS ${GLEW_LOCATION}/lib ) + ELSEIF (WIN32) + # Link library + SET(CMAKE_FIND_LIBRARY_SUFFIXES ".lib") + FIND_LIBRARY ( GLEW_LIBRARY_PATH GLEW32 PATHS ${GLEW_LOCATION}/lib ) + FIND_LIBRARY ( GLEWmx_LIBRARY_PATH GLEW32mx PATHS ${GLEW_LOCATION}/lib ) + # Load library + SET(CMAKE_FIND_LIBRARY_SUFFIXES ".dll") + FIND_LIBRARY ( GLEW_DLL_PATH GLEW32 PATHS ${GLEW_LOCATION}/bin + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + FIND_LIBRARY ( GLEWmx_DLL_PATH GLEW32mx PATHS ${GLEW_LOCATION}/bin + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + ELSE (APPLE) + # Unices + FIND_LIBRARY ( GLEW_LIBRARY_PATH GLEW PATHS ${GLEW_LOCATION}/lib + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + FIND_LIBRARY ( GLEWmx_LIBRARY_PATH GLEWmx PATHS ${GLEW_LOCATION}/lib + NO_DEFAULT_PATH + NO_SYSTEM_ENVIRONMENT_PATH + ) + ENDIF (APPLE) + ENDIF () + # MUST reset + SET(CMAKE_FIND_LIBRARY_SUFFIXES ${ORIGINAL_CMAKE_FIND_LIBRARY_SUFFIXES}) + + SET( GLEW_GLEW_LIBRARY ${GLEW_LIBRARY_PATH} CACHE STRING "GLEW library") + SET( GLEW_GLEWmx_LIBRARY ${GLEWmx_LIBRARY_PATH} CACHE STRING "GLEWmx library") + +ENDIF () diff --git a/FindGLFW.cmake b/FindGLFW.cmake new file mode 100644 index 0000000..d84a2e7 --- /dev/null +++ b/FindGLFW.cmake @@ -0,0 +1,68 @@ +# Find GLFW 3 +# +# GLFW_LIBRARIES +# GLFW_INCLUDE_DIRS. +# GLFW_FOUND +set(GLFW_ROOT $ENV{GLFW_ROOT}) +IF(NOT UNIX) + IF(NOT GLFW_ROOT) + MESSAGE("ERROR: GLFW_ROOT must be set!") + ENDIF(NOT GLFW_ROOT) + + FIND_PATH(GLFW_INCLUDE_DIRS DOC "Path to GLFW include directory." + NAMES GLFW/glfw3.h + PATHS ${GLFW_ROOT}/include) + IF(MSVC15) + FIND_LIBRARY(GLFW_LIBRARIES DOC "Absolute path to GLFW library." + NAMES glfw3.lib + PATHS ${GLFW_ROOT}/lib-vc2015) + ELSEIF(MSVC13) + FIND_LIBRARY(GLFW_LIBRARIES DOC "Absolute path to GLFW library." + NAMES glfw3.lib + PATHS ${GLFW_ROOT}/lib-vc2013) + ELSEIF(MSVC12) + FIND_LIBRARY(GLFW_LIBRARIES DOC "Absolute path to GLFW library." + NAMES glfw3.lib + PATHS ${GLFW_ROOT}/lib-vc2012) + ELSEIF(MSVC10) + FIND_LIBRARY(GLFW_LIBRARIES DOC "Absolute path to GLFW library." + NAMES glfw3.lib + PATHS ${GLFW_ROOT}/lib-vc2010) + ELSEIF(MINGW) + IF(CMAKE_CL_64) + FIND_LIBRARY(GLFW_LIBRARIES DOC "Absolute path to GLFW library." + NAMES glfw3.dll + PATHS ${GLFW_ROOT}/lib-mingw-w64) + ELSE(CMAKE_CL_64) + FIND_LIBRARY(GLFW_LIBRARIES DOC "Absolute path to GLFW library." + NAMES glfw3.dll + PATHS ${GLFW_ROOT}/lib-mingw) + ENDIF(CMAKE_CL_64) + ELSE(MINGW) + # Default to latest version of VC libs + FIND_LIBRARY(GLFW_LIBRARIES DOC "Absolute path to GLFW library." + NAMES glfw3.lib + PATHS ${GLFW_ROOT}/lib-vc2015) + ENDIF(MSVC15) +ELSE(NOT UNIX) + FIND_PATH(GLFW_INCLUDE_DIRS DOC "Path to GLFW include directory." + NAMES GLFW/glfw3.h + PATHS + /usr/include + /usr/local/include + /usr/target/include + /sw/include + /opt/local/include) + + FIND_LIBRARY(GLFW_LIBRARIES DOC "Absolute path to GLFW library." + NAMES glfw3.dll + PATHS + /usr/local/lib + /usr/lib + /lib) +ENDIF(NOT UNIX) + +include(FindPackageHandleStandardArgs) +find_package_handle_standard_args(GLFW DEFAULT_MSG GLFW_LIBRARIES GLFW_INCLUDE_DIRS) + +mark_as_advanced(GLFW_INCLUDE_DIRS GLFW_LIBRARIES) diff --git a/FindSTIM.cmake b/FindSTIM.cmake new file mode 100644 index 0000000..c05f9f6 --- /dev/null +++ b/FindSTIM.cmake @@ -0,0 +1,27 @@ +# finds the STIM library (downloads it if it isn't present) +# set STIMLIB_PATH to the directory containing the stim subdirectory (the stim repository) + +include(FindPackageHandleStandardArgs) + +set(STIM_ROOT $ENV{STIM_ROOT}) + +IF(NOT UNIX) + IF(NOT STIM_ROOT) + MESSAGE("ERROR: STIM_ROOT environment variable must be set!") + ENDIF(NOT STIM_ROOT) + + FIND_PATH(STIM_INCLUDE_DIRS DOC "Path to GLFW include directory." + NAMES stim/image/image.h + PATHS ${STIM_ROOT}) +ENDIF(NOT UNIX) + +find_package_handle_standard_args(STIM DEFAULT_MSG STIM_INCLUDE_DIRS) + +if(STIM_FOUND) + set(STIM_INCLUDE_DIRS ${STIM_INCLUDE_DIRS}) +elseif(STIM_FOUND) + message("STIM library not found. Set the STIM_ROOT environment variable to the STIM location.") + message("STIMLIB can be found here: https://git.stim.ee.uh.edu/codebase/stimlib") +endif(STIM_FOUND) + +find_package_handle_standard_args(STIM DEFAULT_MSG STIM_INCLUDE_DIRS) diff --git a/main.cpp b/main.cpp new file mode 100644 index 0000000..9bed744 --- /dev/null +++ b/main.cpp @@ -0,0 +1,490 @@ +#include +#include + +#include +#include +const double PI = 3.14159265358979323846; + +#include "tensorfield.h" +#include +//#include + +// Global variable to store the GLFW window +GLFWwindow* window; +tensorfield T; +size_t z_slice = 0; +bool cout_frame = false; + +float diffuse_intensity = 0.7f; +float ambient_intensity = 0.3f; + +int glyph_type = 1; +float gamma = 3; +int glyph_resolution = 10; +float glyph_scale = 0.5; + +//glyph vertices +bool glyph_calculated = false; + +struct float3 { + float x; + float y; + float z; +}; + +struct float2 { + float u; + float v; +}; + +std::vector vertices; +std::vector normals; +std::vector texcoords; +std::vector sin_theta; +std::vector cos_theta; +std::vector sin_phi; +std::vector cos_phi; + +struct eigendecomposition { + + stim::vec3 v0; + stim::vec3 v2; + stim::vec3 lambda; +}; + +inline float sgn(float x) { + if (x < 0) + return -1; + else if (x > 0) + return 1; + return 0; +} + +stim::vec3 triangle_norm(stim::vec3 p0, stim::vec3 p1, stim::vec3 p2) { + stim::vec3 a = p1 - p0; + stim::vec3 b = p2 - p0; + stim::vec3 n = a.cross(b); + return n; +} + + +// returns the fractional anisotropy and stores the spherical, linear, and planar components in cs, cl, and cp +float get_anisotropy(stim::vec3 lambda, float& cs, float& cl, float& cp) { + + //calculate the denominator for each specific anisotropy + float denom = lambda[0] + lambda[1] + lambda[2]; + + cl = (lambda[2] - lambda[1]) / denom; + cp = 2 * (lambda[1] - lambda[0]) / denom; + cs = 3 * lambda[0] / denom; + + return 0; +} + +//calculate the eigenvectors and eigenvalues for the input pixel (x, y, z) +// eigenvalues are embedded in the length of the eigenvector +void get_pixel_eigen(eigendecomposition* e, size_t x, size_t y, size_t z) { + + e->v0[0] = T(x, y, z, 0, 0); + e->v0[1] = T(x, y, z, 1, 0); + e->v0[2] = T(x, y, z, 2, 0); + + //temporarily store v1 in order to get the eigenvalue + // (the eigenvector is redundant because it's the cross product of v0 and v1) + stim::vec3 v1; + v1[0] = T(x, y, z, 0, 1); + v1[1] = T(x, y, z, 1, 1); + v1[2] = T(x, y, z, 2, 1); + + e->v2[0] = T(x, y, z, 0, 2); + e->v2[1] = T(x, y, z, 1, 2); + e->v2[2] = T(x, y, z, 2, 2); + + e->lambda[0] = e->v0.len(); + e->lambda[1] = v1.len(); + e->lambda[2] = e->v2.len(); + e->v0 = e->v0 / e->lambda[0]; + e->v2 = e->v2 / e->lambda[2]; +} + +/// Create a rotation matrix to orient a glyph along the tensor direction. Orientations are based on the input vectors +/// v0, v1, and v2 +void get_glyph_rotation_matrix(float* R, stim::vec3 v0, stim::vec3 v2) { + + stim::matrix_sq M; + M(0, 0) = v0[0]; + M(1, 0) = v0[1]; + M(2, 0) = v0[2]; + M(3, 0) = 0; + + M(0, 2) = v2[0]; + M(1, 2) = v2[1]; + M(2, 2) = v2[2]; + M(3, 2) = 0; + + stim::vec3 a(v0[0], v0[1], v0[2]); + stim::vec3 b(v2[0], v2[1], v2[2]); + stim::vec3 c = a.cross(b); + + M(0, 1) = c[0]; + M(1, 1) = c[1]; + M(2, 1) = c[2]; + M(3, 1) = 0; + + M(0, 3) = 0; + M(1, 3) = 0; + M(2, 3) = 0; + M(3, 3) = 1; + + memcpy(R, M.M, 16 * sizeof(float)); + return; +} + + +// Perform the necessary updates when the user reshapes the window +void window_reshape() { + int width, height; + + //get the context size + glfwGetFramebufferSize(window, &width, &height); + + + //create an OpenGL viewport + glViewport(0, 0, width, height); + + //set the default orthographic view (assuming that the window and image aspect ratio are identical) + float left = 0.0f; + float right = (float)T.shape[2]; + float bottom = 0.0f; + float top = (float)T.shape[1]; + + //create an orthographic projection + glMatrixMode(GL_PROJECTION); + glLoadIdentity(); + gluOrtho2D(left, right, bottom, top); + + glMatrixMode(GL_MODELVIEW); + glLoadIdentity(); +} + +void init() { + glEnable(GL_DEPTH_TEST); + glEnable(GL_CULL_FACE); + glFrontFace(GL_CCW); +} + +void generate_glyph_points(int resolution = 5) { + int sectorCount = 2 * resolution; + int stackCount = resolution; + + int num_vertices = (sectorCount + 1) * (stackCount + 1); + sin_theta.resize(num_vertices); + sin_phi.resize(num_vertices); + cos_theta.resize(num_vertices); + cos_phi.resize(num_vertices); + + float theta, phi; // vertex texCoord + + float sectorStep = 2 * (float)PI / sectorCount; + float stackStep = (float)PI / stackCount; + + size_t i = 0; + //calculate the vertex positions + for (int phi_i = 0; phi_i <= stackCount; ++phi_i) + { + + // add (sectorCount+1) vertices per stack + // the first and last vertices have same position and normal, but different tex coords + for (int theta_i = 0; theta_i <= sectorCount; ++theta_i) + { + + // calculate the spherical coordinates of the vertex + theta = (float)theta_i / sectorCount * 2 * (float)PI; + phi = (float)phi_i / stackCount * (float)PI; + + //pre-compute the sine and cosine values used to create the superquadric + sin_theta[i] = sinf(theta); + sin_phi[i] = sinf(phi); + cos_theta[i] = cosf(theta); + cos_phi[i] = cosf(phi); + i++; + } + } +} + +inline float weird_exp(float x, float e) { + return sgn(x) * powf(fabs(x), e); +} + +stim::vec3 qx(size_t i, float alpha = 1, float beta = 1) { + + float sin_phi_beta = weird_exp(sin_phi[i], beta); + float sin_theta_alpha = weird_exp(sin_theta[i], alpha); + float cos_theta_alpha = weird_exp(cos_theta[i], alpha); + float cos_phi_beta = weird_exp(cos_phi[i], beta); + + stim::vec3 p; + p[0] = cos_phi_beta; + p[1] = -sin_theta_alpha * sin_phi_beta; + p[2] = cos_theta_alpha * sin_phi_beta; + return p; +} + +void render_triangle(stim::vec3 p0, stim::vec3 p1, stim::vec3 p2) { + stim::vec3 n = triangle_norm(p0, p1, p2); + + glNormal3f(n[0], n[1], n[2]); + glVertex3f(p0[0], p0[1], p0[2]); + glVertex3f(p1[0], p1[1], p1[2]); + glVertex3f(p2[0], p2[1], p2[2]); +} + +//render a glyph (0 = ellipsoid, 1 = superquadric) +void render_glyph(eigendecomposition* e, int glyph_type = 0, int resolution = 5) { + glMatrixMode(GL_MODELVIEW); + glPushMatrix(); + + stim::vec3 norm_lambda = e->lambda.norm(); + glScalef(glyph_scale, glyph_scale, glyph_scale); + glScalef(norm_lambda[2], norm_lambda[1], norm_lambda[0]); + + if (glyph_type == 0) { + GLUquadric* q = gluNewQuadric(); + gluSphere(q, 1, 2 * resolution, resolution); + } + else if (glyph_type == 1) { + float radius = 0.5; + int sectorCount = 2 * resolution; + int stackCount = resolution; + + if (!glyph_calculated) { + generate_glyph_points(resolution); + glyph_calculated = true; + } + + //get the anisotropy values + float fa, cs, cl, cp, alpha, beta; + fa = get_anisotropy(e->lambda, cs, cl, cp); + + if (cl >= cp) { + alpha = powf(1 - cp, gamma); + beta = powf(1 - cl, gamma); + } + else { + alpha = powf(1 - cl, gamma); + beta = powf(1 - cp, gamma); + } + + //draw the sphere + int k1, k2; + float3 p[3]; + float2 s[3]; + float3 n; + + //draw the glyph + glBegin(GL_TRIANGLES); + for (int i = 0; i < stackCount; ++i) + { + k1 = i * (sectorCount + 1); // beginning of current stack + k2 = k1 + sectorCount + 1; // beginning of next stack + + for (int j = 0; j < sectorCount; ++j, ++k1, ++k2) + { + // 2 triangles per sector excluding first and last stacks + // k1 => k2 => k1+1 + if (i != 0) + { + + stim::vec3 p0 = qx(k1, alpha, beta); + stim::vec3 p1 = qx(k2, alpha, beta); + stim::vec3 p2 = qx(k1 + 1, alpha, beta); + + render_triangle(p0, p1, p2); + } + + // k1+1 => k2 => k2+1 + if (i != (stackCount - 1)) + { + stim::vec3 p0 = qx(k1 + 1, alpha, beta); + stim::vec3 p1 = qx(k2, alpha, beta); + stim::vec3 p2 = qx(k2 + 1, alpha, beta); + + render_triangle(p0, p1, p2); + } + } + } + glEnd(); + } + glPopMatrix(); +} + +void lighting() { + + GLfloat light_position[] = { 1, 1, -1.0, 0.0 }; + GLfloat light_ambient[] = { ambient_intensity, ambient_intensity, ambient_intensity, 1.0 }; + GLfloat light_diffuse[] = { diffuse_intensity, diffuse_intensity, diffuse_intensity, 1.0 }; + glShadeModel(GL_SMOOTH); + + glLightfv(GL_LIGHT0, GL_POSITION, light_position); + glLightfv(GL_LIGHT0, GL_DIFFUSE, light_diffuse); + glLightfv(GL_LIGHT0, GL_AMBIENT, light_ambient); + + glEnable(GL_LIGHTING); + glEnable(GL_LIGHT0); + + glColorMaterial(GL_FRONT, GL_AMBIENT_AND_DIFFUSE); + glEnable(GL_COLOR_MATERIAL); + glEnable(GL_NORMALIZE); + +} +void render() { + size_t zi = z_slice; + + float dc = 1.0f / T.shape[1]; + float c[3]; + eigendecomposition e; + float R[16]; + + + glMatrixMode(GL_MODELVIEW_MATRIX); + + lighting(); + + for (size_t yi = 0; yi < T.shape[1]; yi++) { + for (size_t xi = 0; xi < T.shape[2]; xi++) { + get_pixel_eigen(&e, xi, yi, zi); + + c[0] = abs(e.v0[0]); + c[1] = abs(e.v0[1]); + c[2] = abs(e.v0[2]); + + glPushMatrix(); + + + glTranslatef((float)xi + 0.5f, (float)yi + 0.5f, 0); + + stim::vec3 v0(e.v0[0], e.v0[1], e.v0[2]); + stim::vec3 v2(e.v2[0], e.v2[1], e.v2[2]); + get_glyph_rotation_matrix(R, v0,v2); + glMultMatrixf((GLfloat*)R); + + glColor3f(c[0], c[1], c[2]); + render_glyph(&e, glyph_type, glyph_resolution); + + + glPopMatrix(); + } + } + cout_frame = false; + +} + +void key_callback(GLFWwindow* window, int key, int scancode, int action, int mods) +{ + if (key == GLFW_KEY_RIGHT && action != GLFW_RELEASE) + z_slice++; + else if (key == GLFW_KEY_LEFT && action != GLFW_RELEASE) + z_slice--; + if (z_slice >= T.shape[0]) z_slice = 0; + if (z_slice < 0) z_slice = T.shape[0] - 1; +} + +void display_rotation_matrix(eigendecomposition* e) { + float R[16]; + stim::vec3 v0(e->v0[0], e->v0[1], e->v0[2]); + stim::vec3 v2(e->v2[0], e->v2[1], e->v2[2]); + get_glyph_rotation_matrix(R, v0, v2); + for (int r = 0; r < 4; r++) { + for (int c = 0; c < 4; c++) { + std::cout << R[r * 4 + c] << " "; + } + std::cout << std::endl; + } +} + + +void display_eigen(eigendecomposition* e) { + std::cout << "v0 = (" << e->v0[0] << ", " << e->v0[1] << ", " << e->v0[2] << ")" << std::endl; +} + +void mouse_button_callback(GLFWwindow* window, int button, int action, int mods) { + //if the user clicks inside the window, display information about the tensor field + if (button == GLFW_MOUSE_BUTTON_LEFT && action == GLFW_PRESS) { + double xpos, ypos; + glfwGetCursorPos(window, &xpos, &ypos); //get the position of the mouse pointer + int width, height; + glfwGetFramebufferSize(window, &width, &height); + int pixel_x = (int)(xpos / width * T.shape[2]); + int pixel_y = T.shape[1] - (int)(ypos / height * (int)T.shape[1]) - 1; + std::cout << "----------------------------" << std::endl; + std::cout << "x: " << pixel_x << " y: " << pixel_y << std::endl; + std::cout << "----------------------------" << std::endl; + eigendecomposition e; + get_pixel_eigen(&e, pixel_x, pixel_y, z_slice); + display_eigen(&e); + std::cout << "===========" << std::endl; + display_rotation_matrix(&e); + } +} + +void reshape_callback(GLFWwindow* window, int width, int height) { + window_reshape(); + render(); +} + + +int main(int argc, char** argv) { + std::string filename(argv[1]); + if (T.load_tira(filename) != TIRA_SUCCESS) { + return -1; + } + + /* Initialize the library */ + if (!glfwInit()) + return -1; + + + /* Create a windowed mode window and its OpenGL context */ + window = glfwCreateWindow(700, 700, "Hello World", NULL, NULL); + if (!window) + { + glfwTerminate(); + return -1; + } + + /* Make the window's context current */ + glfwMakeContextCurrent(window); + + glfwSetKeyCallback(window, key_callback); + glfwSetFramebufferSizeCallback(window, reshape_callback); + glfwSetMouseButtonCallback(window, mouse_button_callback); + + GLenum err = glewInit(); + //deal with a GLEW initialization failure + if (GLEW_OK != err) + std::cout << "GLEW Error: " << glewGetErrorString(err) << std::endl; + + //set up the window viewport + window_reshape(); + + //initialize the OpenGL render details + init(); + + /* Loop until the user closes the window */ + while (!glfwWindowShouldClose(window)) + { + /* Render here */ + glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT); + + render(); + + /* Swap front and back buffers */ + glfwSwapBuffers(window); + + /* Poll for and process events */ + glfwPollEvents(); + } + + glfwTerminate(); + return 0; +} \ No newline at end of file diff --git a/tensorfield.cpp b/tensorfield.cpp new file mode 100644 index 0000000..e69de29 --- /dev/null +++ b/tensorfield.cpp diff --git a/tensorfield.h b/tensorfield.h new file mode 100644 index 0000000..959986a --- /dev/null +++ b/tensorfield.h @@ -0,0 +1,122 @@ +#ifndef TIRA_TENSORFIELD +#define TIRA_TENSORFIELD + +#include +#include +#include +#include + +enum tira_floattype {unused0, unused1, float16, unused3, float32, unused4, unused5, unused6, unused7, float64}; +enum tira_error {TIRA_SUCCESS, TIRA_UNDEFINED, TIRA_FILE_OPEN_ERROR, TIRA_UNKNOWN_FILE_TYPE}; +enum tira_filetype {unused, spectrum, tensor2}; +struct tira_header { + std::valarray shape; //stores the resolution of the tensor field along each dimension + tira_filetype filetype; //stores the file type + tira_floattype datatype; //stores the data type of the tensor field + float version; //file version +}; + +tira_error tira_load_header(std::ifstream &infile, tira_header &header) { + + //see if the file represents a TIRA structure + unsigned char TIRA[5]; + TIRA[4] = 0; + infile.read((char*)TIRA, 4); + if (strcmp((char*)TIRA, "TIRA") != 0) + return TIRA_UNKNOWN_FILE_TYPE; + + //get the file type + infile.read((char*)&(header.filetype), sizeof(tira_filetype)); + + //get the version + infile.read((char*)&(header.version), sizeof(float)); + + //get the number of dimensions in the field + size_t ndim; + infile.read((char*)&ndim, sizeof(size_t)); + header.shape.resize(ndim); + + //read each dimension + for (int d = 0; d < ndim; d++) + infile.read((char*)&(header.shape[d]), sizeof(size_t)); + + //read the data type + infile.read((char*)&(header.datatype), sizeof(tira_floattype)); + + return TIRA_SUCCESS; +} + +template +class tensorfield { +public: + + std::valarray shape; //number of samples in each dimension + std::valarray cum_offsets; //cumulative factors for fast array accesses + tira_filetype filetype; //stores the file type + tira_floattype datatype; //stores the data type of the tensor field + T* ptr; //pointer to the raw tensor data + + void init() { + size_t ndims = shape.size(); //get the number of dimensions + cum_offsets.resize(ndims); //allocate space for the cumulative offsets + size_t cum = 1; + for (int i = 0; i < ndims; i++) { + cum_offsets[ndims - 1 - i] = cum; + cum *= shape[ndims - 1 - i]; + } + } + + size_t bytes() { + size_t product = 1; + for (size_t i = 0; i < shape.size(); i++) + product *= shape[i]; + return product * (size_t)datatype; + } + + // Function turns a series of array indices into a 1D index to the corresponding element + inline size_t idx(std::initializer_list coords) { + size_t index = 0; + size_t i = 0; + for (auto c : coords) { + index += c * cum_offsets[i]; + i++; + } + return index; + } + + tira_error load_tira(std::string filename) { + + //open the file for binary reading + std::ifstream infile(filename, std::ios::binary); + if (!infile) return TIRA_FILE_OPEN_ERROR; + + tira_header header; + tira_load_header(infile, header); + shape = header.shape; + filetype = header.filetype; + datatype = header.datatype; + + //allocate space in memory for the tensor field + size_t size_bytes = bytes(); + ptr = (T*)malloc(size_bytes); + + //copy data from the file to memory + infile.read((char*)ptr, bytes()); + + //close the file + infile.close(); + + //initialize all internal variables + init(); + + return TIRA_SUCCESS; + } + + T& operator()(size_t x, size_t y, size_t z, size_t r, size_t c) { + size_t i = idx({ z, y, x, r, c }); + return ptr[i]; + } + +}; + +#endif \ No newline at end of file -- libgit2 0.21.4