From ebb721c7e4788ddc420c63fc3292664ef83ac89e Mon Sep 17 00:00:00 2001 From: David Date: Fri, 8 Jan 2016 12:30:43 -0600 Subject: [PATCH] new repository for STIM lab --- 00_GT.obj | 29 +++++++++++++++++++++++++++++ 00_T.obj | 50 ++++++++++++++++++++++++++++++++++++++++++++++++++ 01_GT.obj | 9 +++++++++ 01_T.obj | 18 ++++++++++++++++++ 02_GT.obj | 21 +++++++++++++++++++++ 02_T.obj | 22 ++++++++++++++++++++++ 03_GT.obj | 40 ++++++++++++++++++++++++++++++++++++++++ 03_T.obj | 16 ++++++++++++++++ CMakeLists.txt | 59 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ DrawingFunctions.h | 919 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ ErrorMap_Fragment.glsl | 74 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ ErrorShader_Fragment.glsl | 72 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ FindANN.cmake | 42 ++++++++++++++++++++++++++++++++++++++++++ FindGLEW.cmake | 49 +++++++++++++++++++++++++++++++++++++++++++++++++ FindGLUT.cmake | 42 ++++++++++++++++++++++++++++++++++++++++++ GroundTruth.swc | 243 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ GroundTruthCopy.swc | 243 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ ImplicitCalculations.h | 436 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ NetworkComparison.cpp | 341 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ Node_SmoothError_Vertex.glsl | 11 +++++++++++ PerformanceData.h | 155 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ README | 9 +++++++++ SmoothShader_Fragment.glsl | 35 +++++++++++++++++++++++++++++++++++ SmoothShader_Vertex.glsl | 11 +++++++++++ anyoption.cpp | 1174 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ anyoption.h | 270 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ loop_proxy_GT.obj | 36 ++++++++++++++++++++++++++++++++++++ loop_proxy_T.obj | 41 +++++++++++++++++++++++++++++++++++++++++ molecule_GT.obj | 83 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ molecule_T.obj | 90 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/CHECK_OPENGL_ERROR.h | 15 +++++++++++++++ rts/objJedi.h | 1857 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rtsCamera.h | 191 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rtsFilename.h | 67 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rtsGraph.h | 142 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rtsLinearAlgebra.h | 15 +++++++++++++++ rts/rtsMatrix4x4.h | 196 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rtsPoint3d.h | 141 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rtsQuaternion.h | 144 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rtsSourceCode.h | 65 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rtsVector3d.h | 164 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rts_glShaderObject.h | 115 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rts_glShaderProgram.h | 295 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rts_glShaderUniform.h | 105 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rts_glTextureMap.h | 260 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rts/rts_glutRenderWindow.h | 155 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ rtsFiberNetwork.h | 2667 +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 47 files changed, 11234 insertions(+), 0 deletions(-) create mode 100644 00_GT.obj create mode 100644 00_T.obj create mode 100644 01_GT.obj create mode 100644 01_T.obj create mode 100644 02_GT.obj create mode 100644 02_T.obj create mode 100644 03_GT.obj create mode 100644 03_T.obj create mode 100644 CMakeLists.txt create mode 100644 DrawingFunctions.h create mode 100644 ErrorMap_Fragment.glsl create mode 100644 ErrorShader_Fragment.glsl create mode 100644 FindANN.cmake create mode 100644 FindGLEW.cmake create mode 100644 FindGLUT.cmake create mode 100644 GroundTruth.swc create mode 100644 GroundTruthCopy.swc create mode 100644 ImplicitCalculations.h create mode 100644 NetworkComparison.cpp create mode 100644 Node_SmoothError_Vertex.glsl create mode 100644 PerformanceData.h create mode 100644 README create mode 100644 SmoothShader_Fragment.glsl create mode 100644 SmoothShader_Vertex.glsl create mode 100644 anyoption.cpp create mode 100644 anyoption.h create mode 100644 loop_proxy_GT.obj create mode 100644 loop_proxy_T.obj create mode 100644 molecule_GT.obj create mode 100644 molecule_T.obj create mode 100644 rts/CHECK_OPENGL_ERROR.h create mode 100644 rts/objJedi.h create mode 100644 rts/rtsCamera.h create mode 100644 rts/rtsFilename.h create mode 100644 rts/rtsGraph.h create mode 100644 rts/rtsLinearAlgebra.h create mode 100644 rts/rtsMatrix4x4.h create mode 100644 rts/rtsPoint3d.h create mode 100644 rts/rtsQuaternion.h create mode 100644 rts/rtsSourceCode.h create mode 100644 rts/rtsVector3d.h create mode 100644 rts/rts_glShaderObject.h create mode 100644 rts/rts_glShaderProgram.h create mode 100644 rts/rts_glShaderUniform.h create mode 100644 rts/rts_glTextureMap.h create mode 100644 rts/rts_glutRenderWindow.h create mode 100644 rtsFiberNetwork.h diff --git a/00_GT.obj b/00_GT.obj new file mode 100644 index 0000000..0e66fe5 --- /dev/null +++ b/00_GT.obj @@ -0,0 +1,29 @@ +v 99 677 0 +v 229 510 0 +v 57 347 0 +v 344 724 0 +v 581 361 0 +v 341 199 0 +v 757 377 0 +v 981 152 0 +v 793 508 0 +v 941 600 0 +v 647 738 0 +v 88 563 0 +v 181 372 0 +v 316 596 0 +v 385 425 0 +v 481 165 0 +v 574 223 0 +v 818 229 0 +v 679 601 0 +l 1 12 2 +l 3 13 2 +l 2 14 4 +l 2 15 5 +l 6 16 17 5 +l 5 7 +l 8 18 7 +l 7 9 +l 9 10 +l 9 19 11 diff --git a/00_T.obj b/00_T.obj new file mode 100644 index 0000000..3eee399 --- /dev/null +++ b/00_T.obj @@ -0,0 +1,50 @@ +v 119 675 0 +v 343 730 0 +v 59 346 0 +v 577 361 0 +v 347 199 0 +v 756 375 0 +v 588 483 0 +v 789 514 0 +v 676 603 0 +v 652 744 0 +v 933 598 0 +v 105 611 0 +v 126 537 0 +v 227 519 0 +v 294 580 0 +v 319 677 0 +v 165 366 0 +v 192 438 0 +v 234 473 0 +v 294 463 0 +v 375 428 0 +v 441 404 0 +v 510 388 0 +v 419 188 0 +v 492 174 0 +v 551 205 0 +v 567 273 0 +v 572 334 0 +v 645 365 0 +v 695 410 0 +v 633 438 0 +v 776 424 0 +v 788 467 0 +v 742 552 0 +v 702 584 0 +v 663 659 0 +v 649 720 0 +v 839 539 0 +v 882 557 0 +v 804 606 0 +l 1 12 13 14 15 16 2 +l 3 17 18 19 20 21 22 23 4 +l 5 24 25 26 27 28 4 +l 4 29 6 +l 6 30 31 7 +l 6 32 33 8 +l 8 34 35 9 +l 9 36 37 10 +l 8 38 39 11 +l 9 40 11 diff --git a/01_GT.obj b/01_GT.obj new file mode 100644 index 0000000..6499cd4 --- /dev/null +++ b/01_GT.obj @@ -0,0 +1,9 @@ +v 291 412 0 +v 464 251 0 +v 635 405 0 +v 463 602 0 +v 635 405 0 +l 1 2 +l 2 3 +l 3 5 1 +l 3 4 diff --git a/01_T.obj b/01_T.obj new file mode 100644 index 0000000..10cc8a2 --- /dev/null +++ b/01_T.obj @@ -0,0 +1,18 @@ +v 292 411 0 +v 464 252 0 +v 531 312 0 +v 637 404 0 +v 537 405 0 +v 442 408 0 +v 368 409 0 +v 464 603 0 +v 531 312 0 +v 442 408 0 +l 1 2 +l 2 3 +l 3 4 +l 4 5 +l 3 9 5 +l 6 10 5 +l 1 7 +l 4 8 diff --git a/02_GT.obj b/02_GT.obj new file mode 100644 index 0000000..182e8d3 --- /dev/null +++ b/02_GT.obj @@ -0,0 +1,21 @@ +v 485 172 0 +v 681 360 0 +v 548 544 0 +v 300 356 0 +v 372 545 0 +v 481 378 0 +v 372 545 0 +v 481 378 0 +v 481 378 0 +v 481 378 0 +v 481 378 0 +l 1 2 +l 2 3 +l 1 4 +l 4 5 +l 5 7 3 +l 1 6 +l 6 8 3 +l 6 9 2 +l 6 10 5 +l 6 11 4 diff --git a/02_T.obj b/02_T.obj new file mode 100644 index 0000000..4b51745 --- /dev/null +++ b/02_T.obj @@ -0,0 +1,22 @@ +v 334 445 0 +v 370 544 0 +v 549 545 0 +v 681 360 0 +v 484 172 0 +v 469 477 0 +v 484 172 0 +v 469 477 0 +v 370 544 0 +v 549 545 0 +v 469 477 0 +v 484 172 0 +l 1 2 +l 2 3 +l 3 4 +l 4 5 +l 5 7 1 +l 6 8 1 +l 2 9 6 +l 3 10 6 +l 6 11 4 +l 5 12 6 diff --git a/03_GT.obj b/03_GT.obj new file mode 100644 index 0000000..3b968ba --- /dev/null +++ b/03_GT.obj @@ -0,0 +1,40 @@ +v 212 603 0 +v 365 427 0 +v 125 308 0 +v 457 253 0 +v 654 335 0 +v 578 630 0 +v 565 118 0 +v 813 145 0 +v 801 451 0 +v 949 54 0 +v 276 494 0 +v 208 348 0 +v 284 382 0 +v 404 352 0 +v 381 408 0 +v 414 444 0 +v 492 452 0 +v 554 451 0 +v 602 403 0 +v 662 429 0 +v 648 487 0 +v 617 570 0 +v 633 266 0 +v 614 219 0 +v 594 172 0 +v 705 322 0 +v 747 299 0 +v 796 243 0 +v 855 233 0 +v 825 359 0 +v 884 89 0 +l 1 11 2 +l 3 12 13 2 +l 4 14 15 2 +l 2 16 17 18 19 5 +l 5 20 21 22 6 +l 5 23 24 25 7 +l 5 26 27 28 8 +l 8 29 30 9 +l 8 31 10 diff --git a/03_T.obj b/03_T.obj new file mode 100644 index 0000000..981eaf3 --- /dev/null +++ b/03_T.obj @@ -0,0 +1,16 @@ +v 212 604 0 +v 123 310 0 +v 365 427 0 +v 576 631 0 +v 457 253 0 +v 563 116 0 +v 813 144 0 +v 800 450 0 +v 951 57 0 +v 800 450 0 +l 1 2 +l 3 4 +l 5 6 +l 6 7 +l 8 10 4 +l 8 9 diff --git a/CMakeLists.txt b/CMakeLists.txt new file mode 100644 index 0000000..ffa63fe --- /dev/null +++ b/CMakeLists.txt @@ -0,0 +1,59 @@ +#Specify the version being used aswell as the language +cmake_minimum_required(VERSION 2.8) +#Name your project here +project(netmets) + +#set the module directory +set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} "${CMAKE_SOURCE_DIR}") + +#find OpenGL +find_package(OpenGL REQUIRED) + +#find GLUT +set(GLUT_ROOT_PATH $ENV{GLUT_ROOT_PATH}) +find_package(GLUT REQUIRED) + +#find GLEW +find_package(GLEW REQUIRED) + +#find BOOST +find_package(Boost REQUIRED) + +#find the Approximate Nearest Neighbor Library +find_package(ANN REQUIRED) + +#set the include directories +include_directories( + ${CMAKE_CURRENT_BINARY_DIR} + ${OPENGL_INCLUDE_DIR} + ${GLEW_INCLUDE_PATH} + ${GLUT_INCLUDE_DIR} + ${ANN_INCLUDE_DIR} + ${Boost_INCLUDE_DIR} +) + +#Assign source files to the appropriate variables +file(GLOB SRC_CPP "*.cpp") +file(GLOB SRC_H "*.h") + +#set up copying data files +configure_file(00_GT.obj ${CMAKE_CURRENT_BINARY_DIR}/00_GT.obj @ONLY) +configure_file(00_T.obj ${CMAKE_CURRENT_BINARY_DIR}/00_T.obj @ONLY) +configure_file(01_GT.obj ${CMAKE_CURRENT_BINARY_DIR}/01_GT.obj @ONLY) +configure_file(01_T.obj ${CMAKE_CURRENT_BINARY_DIR}/01_T.obj @ONLY) +configure_file(02_GT.obj ${CMAKE_CURRENT_BINARY_DIR}/02_GT.obj @ONLY) +configure_file(02_T.obj ${CMAKE_CURRENT_BINARY_DIR}/02_T.obj @ONLY) +configure_file(03_GT.obj ${CMAKE_CURRENT_BINARY_DIR}/03_GT.obj @ONLY) +configure_file(03_T.obj ${CMAKE_CURRENT_BINARY_DIR}/03_T.obj @ONLY) +configure_file(SmoothShader_Vertex.glsl ${CMAKE_CURRENT_BINARY_DIR}/SmoothShader_Vertex.glsl @ONLY) +configure_file(SmoothShader_Fragment.glsl ${CMAKE_CURRENT_BINARY_DIR}/SmoothShader_Fragment.glsl @ONLY) +configure_file(ErrorMap_Fragment.glsl ${CMAKE_CURRENT_BINARY_DIR}/ErrorMap_Fragment.glsl @ONLY) + +#create an executable +add_executable(netmets ${SRC_CPP} ${SRC_H}) + +#set the link libraries +target_link_libraries(netmets ${ANN_LIBRARY} ${OPENGL_gl_LIBRARY} ${OPENGL_glu_LIBRARY} ${GLEW_LIBRARY} ${GLUT_glut_LIBRARY}) + + + diff --git a/DrawingFunctions.h b/DrawingFunctions.h new file mode 100644 index 0000000..eaba8f9 --- /dev/null +++ b/DrawingFunctions.h @@ -0,0 +1,919 @@ +#include "rtsFiberNetwork.h" +#include "rts/rts_glShaderProgram.h" +#include "GL/glut.h" +#include "rts/rts_glutRenderWindow.h" +#include + +extern void ComputeNetMets(); +extern rtsFiberNetwork* goldNetwork; +extern rtsFiberNetwork* testNetwork; +extern float sigmaG, sigmaC; +float network_span; +CoreGraphList coreGraph; +vector > sequenceColors; +int current_sequence = 0; + +//shader variables +rts_glShaderProgram Edge_ErrorShader; +rts_glShaderProgram Node_ErrorShader; +rts_glShaderProgram Smooth_Shader; + +//display lists +GLuint GT_FibersList=0; +GLuint T_FibersList=0; +GLuint GT_EndCaps=0; +GLuint T_EndCaps=0; +GLuint GT_NodesList=0; +GLuint T_NodesList=0; +GLuint T_PathList=0; +GLuint GT_PathList=0; + +//drawing variables +int tube_subdivisions = 20; +float node_radius_factor = 0.7; +float fiber_radius_factor = 0.5; +float cull_test_case_threshold = 1.0; + +#define DISPLAY_GT_NETWORK 1 +#define DISPLAY_T_NETWORK 2 +#define DISPLAY_GT_GRAPH 3 +#define DISPLAY_T_GRAPH 4 +#define DISPLAY_GT_SELECTED 5 +#define DISPLAY_T_SELECTED 6 + +//menu options +#define NETCOMP_EXIT 0 +#define DISPLAY_NETWORK 1 +#define DISPLAY_GRAPH 2 +#define DISPLAY_CONNECTED 3 +#define DISPLAY_SELECTED 4 +#define COLORMAP_ISOLUMINANT 5 +#define COLORMAP_BLACKBODY 6 +#define COLORMAP_BREWER 7 +#define COLORMAP_POLAR_CIELAB 8 +#define COLORMAP_RAINBOW 9 +#define CULL_TEST_CASE 10 +#define RECOMPUTE_METRIC 11 + + +//fibers to render in Graph Mode +//list T_DisplayEdges; +//list GT_DisplayEdges; + +int DisplayMode = DISPLAY_NETWORK; +float L0_pos[3]; +float L1_pos[3]; +//GLuint texColorMap=0; +rts_glTextureMap texColorMap; + + + + +void makeColormap(int ColorMapType = COLORMAP_BREWER) +{ + //if(texColorMap != 0) + // glDeleteTextures(1, &texColorMap); + + point3D* ctrlPts; + int num_points = 0; + if(ColorMapType == COLORMAP_ISOLUMINANT) + { + //allocate memory for the colormap + num_points = 2; + ctrlPts = new point3D[num_points]; + //memset(ctrlPts, 0, num_points*sizeof(point3D)); + + ctrlPts[0] = point3D(0.0, 1.0, 0.0); + ctrlPts[1] = point3D(1.0, 0.0, 0.0); + } + else if(ColorMapType == COLORMAP_RAINBOW) + { + //allocate memory for the colormap + num_points = 5; + ctrlPts = new point3D[num_points]; + //memset(ctrlPts, 0, num_points*sizeof(point3D)); + + //ctrlPts[0] = point3D(0.7, 0, 0.7); + ctrlPts[0] = point3D(0, 0, 1); + ctrlPts[1] = point3D(0, 0.7, 0.7); + ctrlPts[2] = point3D(0, 1, 0); + ctrlPts[3] = point3D(0.7, 0.7, 0); + ctrlPts[4] = point3D(1, 0, 0); + } + else if(ColorMapType == COLORMAP_BLACKBODY) + { + //allocate memory for the colormap + num_points = 4; + ctrlPts = new point3D[num_points]; + //memset(ctrlPts, 0, num_points*sizeof(point3D)); + + ctrlPts[0] = point3D(0.0, 0.0, 0.0); + ctrlPts[1] = point3D(1.0, 0.0, 0.0); + ctrlPts[2] = point3D(1.0, 1.0, 0.0); + ctrlPts[3] = point3D(1.0, 1.0, 1.0); + } + else if(ColorMapType == COLORMAP_BREWER) + { + //allocate memory for the colormap + num_points = 11; + ctrlPts = new point3D[num_points]; + //memset(ctrlPts, 0, num_points*sizeof(point3D)); + + ctrlPts[0] = point3D(0.192157, 0.211765, 0.584314); + ctrlPts[1] = point3D(0.270588, 0.458824, 0.705882); + ctrlPts[2] = point3D(0.454902, 0.678431, 0.819608); + ctrlPts[3] = point3D(0.670588, 0.85098, 0.913725); + ctrlPts[4] = point3D(0.878431, 0.952941, 0.972549); + ctrlPts[5] = point3D(1, 1, 0.74902); + ctrlPts[6] = point3D(0.996078, 0.878431, 0.564706); + ctrlPts[7] = point3D(0.992157, 0.682353, 0.380392); + ctrlPts[8] = point3D(0.956863, 0.427451, 0.262745); + ctrlPts[9] = point3D(0.843137, 0.188235, 0.152941); + ctrlPts[10] = point3D(0.647059, 0, 0.14902); + + } + else if(ColorMapType == COLORMAP_POLAR_CIELAB) + { + //allocate memory for the colormap + num_points = 33; + ctrlPts = new point3D[num_points]; + //memset(ctrlPts, 0, num_points*sizeof(point3D)); + + ctrlPts[0] = point3D(0.07514311, 0.468049805,1); + ctrlPts[1] = point3D(0.247872569, 0.498782363,1); + ctrlPts[2] = point3D(0.339526309, 0.528909511,1); + ctrlPts[3] = point3D(0.409505078, 0.558608486,1); + ctrlPts[4] = point3D(0.468487184, 0.588057293,1); + ctrlPts[5] = point3D(0.520796675, 0.617435078,1); + ctrlPts[6] = point3D(0.568724526, 0.646924167,1); + ctrlPts[7] = point3D(0.613686735, 0.676713218,1); + ctrlPts[8] = point3D(0.656658579, 0.707001303,1); + ctrlPts[9] = point3D(0.698372844, 0.738002964,1); + ctrlPts[10] = point3D(0.739424025, 0.769954435,1); + ctrlPts[11] = point3D(0.780330104, 0.803121429,1); + ctrlPts[12] = point3D(0.821573924, 0.837809045,1); + ctrlPts[13] = point3D(0.863634967, 0.874374691,1); + ctrlPts[14] = point3D(0.907017747, 0.913245283,1); + ctrlPts[15] = point3D(0.936129275, 0.938743558, 0.983038586); + ctrlPts[16] = point3D(0.943467973, 0.943498599, 0.943398095); + ctrlPts[17] = point3D(0.990146732, 0.928791426, 0.917447482); + ctrlPts[18] = point3D(1, 0.88332677, 0.861943246); + ctrlPts[19] = point3D(1, 0.833985467, 0.803839606); + ctrlPts[20] = point3D(1, 0.788626485, 0.750707739); + ctrlPts[21] = point3D(1, 0.746206642, 0.701389973); + ctrlPts[22] = point3D(1, 0.70590052, 0.654994046); + ctrlPts[23] = point3D(1, 0.667019783, 0.610806959); + ctrlPts[24] = point3D(1, 0.6289553, 0.568237474); + ctrlPts[25] = point3D(1, 0.591130233, 0.526775617); + ctrlPts[26] = point3D(1, 0.552955184, 0.485962266); + ctrlPts[27] = point3D(1, 0.513776083, 0.445364274); + ctrlPts[28] = point3D(1, 0.472800903, 0.404551679); + ctrlPts[29] = point3D(1, 0.428977855, 0.363073592); + ctrlPts[30] = point3D(1, 0.380759558, 0.320428137); + ctrlPts[31] = point3D(0.961891484, 0.313155629, 0.265499262); + ctrlPts[32] = point3D(0.916482116, 0.236630659, 0.209939162); + + } + + texColorMap.Init(ctrlPts, GL_TEXTURE_1D, num_points, 0, 0, GL_RGB, GL_RGB, GL_FLOAT); + + //glGenTextures(1, &texColorMap); + //glBindTexture(GL_TEXTURE_1D, texColorMap); + //glTexImage1D(GL_TEXTURE_1D, 0, GL_RGB, num_points, 0, GL_RGB, GL_FLOAT, ctrlPts); +} + + + +point3D HSLtoRGB(point3D HSL) +{ + float H = HSL.x; + float S = HSL.y; + float L = HSL.z; + + float temp2; + if(L < 0.5) + temp2 = L*(1.0+S); + else + temp2 = L+S - L*S; + + float temp1 = 2.0*L - temp2; + + point3D temp3(H+1.0/3.0, H, H-1.0/3.0); + if(temp3.x < 0.0) temp3.x = temp3.x+1.0; + if(temp3.y < 0.0) temp3.y = temp3.y+1.0; + if(temp3.z < 0.0) temp3.z = temp3.z+1.0; + + if(temp3.x > 1.0) temp3.x = temp3.x - 1.0; + if(temp3.y > 1.0) temp3.y = temp3.y - 1.0; + if(temp3.z > 1.0) temp3.z = temp3.z - 1.0; + + point3D result; + if(6.0*temp3.x < 1.0) result.x = temp1 +(temp2 - temp1)*6.0*temp3.x; + else if(2.0*temp3.x < 1.0) result.x = temp2; + else if(3.0*temp3.x < 2.0) result.x = temp1+(temp2-temp1)*((2.0/3.0) - temp3.x)*6.0; + else result.x = temp1; + + if(6.0*temp3.y < 1.0) result.y = temp1 +(temp2 - temp1)*6.0*temp3.y; + else if(2.0*temp3.y < 1.0) result.y = temp2; + else if(3.0*temp3.y < 2.0) result.y = temp1+(temp2-temp1)*((2.0/3.0) - temp3.y)*6.0; + else result.y = temp1; + + if(6.0*temp3.z < 1.0) result.z = temp1 +(temp2 - temp1)*6.0*temp3.z; + else if(2.0*temp3.z < 1.0) result.z = temp2; + else if(3.0*temp3.z < 2.0) result.z = temp1+(temp2-temp1)*((2.0/3.0) - temp3.z)*6.0; + else result.z = temp1; + + //result.a = 0.0; + return result; +} +void ColorFibers() +{ + //srand(time(NULL)); + sequenceColors.clear(); + //for each fiber + for(CoreGraphList::iterator i = coreGraph.begin(); i!=coreGraph.end(); i++) + { + float random_hue = (double)rand()/(double)RAND_MAX; + //cout<<"Random Hue: "< rgb = HSLtoRGB(point3D(random_hue, random_saturation, 0.5)); + //point3D rgb((double)rand()/(double)RAND_MAX, (double)rand()/(double)RAND_MAX, (double)rand()/(double)RAND_MAX); + sequenceColors.push_back(rgb); + } + + +} +void CenterCameraToSelected() +{ + if(coreGraph.size() == 0) + return; + + //center the fiber in both networks + point3D min_pt(9999, 9999, 9999); + point3D max_pt(-9999, -9999, -9999); + + //iterate through the first edge sequence + EdgeSequence::iterator i; + int node; + point3D test; + for(i=coreGraph[current_sequence].first.begin(); i!=coreGraph[current_sequence].first.end(); i++) + { + node = testNetwork->FiberList[*i].n0; + test = testNetwork->NodeList[node].p; + min_pt.x = min(test.x, min_pt.x); + min_pt.y = min(test.y, min_pt.y); + min_pt.z = min(test.z, min_pt.z); + max_pt.x = max(test.x, max_pt.x); + max_pt.y = max(test.y, max_pt.y); + max_pt.z = max(test.z, max_pt.z); + + node = testNetwork->FiberList[*i].n1; + test = testNetwork->NodeList[node].p; + min_pt.x = min(test.x, min_pt.x); + min_pt.y = min(test.y, min_pt.y); + min_pt.z = min(test.z, min_pt.z); + max_pt.x = max(test.x, max_pt.x); + max_pt.y = max(test.y, max_pt.y); + max_pt.z = max(test.z, max_pt.z); + } + point3D middle = min_pt+0.5*(max_pt - min_pt); + + rts_glut_camera.LookAt(middle); + +} +void IncrementSelectedFiber(int i) +{ + //get the currently selected fiber id + if(coreGraph.size() <= 0) + return; + + //get the number of fibers + int end_id = coreGraph.size(); + + current_sequence+=i; + if(current_sequence >= end_id) + current_sequence = 0; + if(current_sequence < 0) + current_sequence = coreGraph.size()-1; + + //print the selected edges + EdgeSequence::iterator EdgeI; + + for(EdgeI = coreGraph[current_sequence].first.begin(); EdgeI != coreGraph[current_sequence].first.end(); EdgeI++) + cout<<*EdgeI<<" "; + cout<<"--->"; + for(EdgeI = coreGraph[current_sequence].second.begin(); EdgeI != coreGraph[current_sequence].second.end(); EdgeI++) + cout<<*EdgeI<<" "; + cout< p; + //glColor3f(network->NodeList[n].error, 0.0, 0.0); + p = network->NodeList[n].p; + + glTranslatef(p.x, p.y, p.z); + //glutSolidSphere(node_radius*standard_deviation, 20, 20); + glTexCoord1f(network->NodeList[n].error); + if(network->NodeList[n].color < 0) + glColor3f(1.0, 0.0, 0.0); + else + glColor3f(1.0, 1.0, 1.0); + gluSphere(quadric,radius,32,32); + + glPopMatrix(); + +} +void DrawNodeSpheres(rtsFiberNetwork* network, float radius) +{ + + + unsigned int n; + for(n=0; n != network->FiberList.size(); n++) + { + if(!network->isCulled(n)) + { + DrawNodeSphere(network, network->FiberList[n].n0,radius); + DrawNodeSphere(network, network->FiberList[n].n1,radius); + } + } +} + +void FrenetFrame(vector3D &x, vector3D &y, vector3D &z) +{ + x = vector3D(0.0, 0.0, 1.0); + y = x.X(z); + x = z.X(y); + x.Normalize(); + y.Normalize(); + z.Normalize(); +} + +vector3D GetColor(float error) +{ + //This function converts an error value to a color + //The conversion is done by creating an HSV color from the error value and converting that HSV color to RGB + float H = (240.0/60.0)*(1.0 - error); + float S = 1.0; + float V = 1.0; + + int i = floor(H); + float f = H - i; + if(i%2 == 0) + f = 1-f; + float m = V*(1 - S); + float n = V*(1-S*f); + switch(i) + { + case 0: + return vector3D(V, n, m); + case 1: + return vector3D(n, V, m); + case 2: + return vector3D(m, V, n); + case 3: + return vector3D(m, n, V); + case 4: + return vector3D(n, m, V); + case 5: + return vector3D(V, m, n); + default: + return vector3D(0, 0, 0); + } + +} + +void DrawTube(point3D p0, vector3D d0, point3D p1, vector3D d1, float error0, float error1, float radius, int subdiv) +{ + + //draw the first circle + vector3D x0, y0, z0, x1, y1, z1; + z0 = d0; + FrenetFrame(x0, y0, z0); + + z1 = d1; + FrenetFrame(x1, y1, z1); + + float t_step = (2*3.14159)/subdiv; + + + float u, v; + + //get the RGB color + point3D circle0, circle1; + vector3D RGB0, RGB1; + vector3D normal; + + //RGB0 = GetColor(color0); + //RGB1 = GetColor(color1); + + glBegin(GL_TRIANGLE_STRIP); + for(int t=0; t<=subdiv; t++) + { + u = radius * cos(t*t_step); + v = radius * sin(t*t_step); + normal = u*x0 + v*y0; + circle0 = p0 + normal; + normal.Normalize(); + + glTexCoord1f(error0); + //glColor4f(error0, 0.0, 0.0, 1.0); + glNormal3f(normal.x, normal.y, normal.z); + glVertex3f(circle0.x, circle0.y, circle0.z); + + normal = u*x1 + v*y1; + circle1 = p1 + normal; + normal.Normalize(); + + glTexCoord1f(error1); + //glColor4f(error1, 0.0, 0.0, 1.0); + glNormal3f(normal.x, normal.y, normal.z); + glVertex3f(circle1.x, circle1.y, circle1.z); + + } + glEnd(); + CHECK_OPENGL_ERROR +} + +void ExtrudeFiber(rtsFiberNetwork* network, int fiber, float radius) +{ + vector3D x, y, z; + point3D p0, p1, p2, p3; + vector3D d1, d2; + float e1, e2; + + //get the first point + int node = network->FiberList[fiber].n0; + p1 = network->NodeList[node].p; + e1 = network->NodeList[node].error; + + //for each vertex in the fiber + int num_points = (int)network->FiberList[fiber].pointList.size(); + for(int v=0; vFiberList[fiber].pointList[v]; + e2 = network->FiberList[fiber].errorList[v]; + + if(vFiberList[fiber].pointList[v+1]; + else + { + node = network->FiberList[fiber].n1; + p3 = network->NodeList[node].p; + } + + d2 = p3-p1; + + //compute the fiber derivatives at p1 and p2 + if(v==0) //if this is the first fiber + d1 = p2 - p1; + else + { + d1 = p2 - p0; + } + + DrawTube(p1, d1, p2, d2, e1, e2, radius, tube_subdivisions); + + //shift + p0 = p1; + p1 = p2; + e1 = e2; + } + //make the last tube + + //if there were any points in the pointlist + if(num_points > 0) + { + p2 = p3; + node = network->FiberList[fiber].n1; + e2 = network->NodeList[node].error; + d1 = p2-p0; + d2 = p2-p1; + DrawTube(p1, d1, p2, d2, e1, e2, radius, tube_subdivisions); + } + //if there are only the two node points + else + { + node = network->FiberList[fiber].n1; + p2 = network->NodeList[node].p; + e2 = network->NodeList[node].error; + d1 = p2 - p1; + d2 = p2 - p1; + DrawTube(p1, d1, p2, d2, e1, e2, radius, tube_subdivisions); + } +} + +void DrawLineFiber(rtsFiberNetwork* network, int f) +{ + point3D p; + int node = network->FiberList[f].n0; + p = network->NodeList[node].p; + + glBegin(GL_LINE_STRIP); + glVertex3f(p.x, p.y, p.z); + for(int v=0; v!=network->FiberList[f].pointList.size(); v++) + { + p = network->FiberList[f].pointList[v]; + glVertex3f(p.x, p.y, p.z); + } + node = network->FiberList[f].n1; + p = network->NodeList[node].p; + glVertex3f(p.x, p.y, p.z); + glEnd(); + +} +void DrawLineNetwork(rtsFiberNetwork* network) +{ + + int num_fibers = network->FiberList.size(); + for(int f = 0; f < num_fibers; f++) + { + /*if(network->FiberList[f].mapped_to == -1) + glColor3f(0.0,0.0, 0.0); + else + glColor3f(1.0, 0.0, 0.0); + */ + DrawLineFiber(network, f); + CHECK_OPENGL_ERROR + } + + +} +void DrawGraphNodes(rtsFiberNetwork* network) +{ + //renders graph nodes, colored based on their node color + glMatrixMode(GL_MODELVIEW); + + unsigned int n; + for(n=0; n != network->NodeList.size(); n++) + { + glPushMatrix(); + + point3D p; + /*if(network->NodeList[n].color < 0) + glColor3f(1.0, 0.0, 0.0); + else + glColor3f(0.0, 1.0, 0.0); + */ + //glColor3f(network->NodeList[n].error, 0.0, 0.0); + p = network->NodeList[n].p; + + glTranslatef(p.x, p.y, p.z); + glutSolidSphere(node_radius_factor*sigmaC, 20, 20); + + glPopMatrix(); + } +} +void DrawFiberSequence(rtsFiberNetwork* network, EdgeSequence sequence, float fiber_radius, float node_radius) +{ + //glClear(GL_DEPTH_BUFFER_BIT); + for(EdgeSequence::iterator i = sequence.begin(); i != sequence.end(); i++) + { + ExtrudeFiber(network, *i, fiber_radius); + glPushAttrib(GL_CURRENT_BIT); + DrawNodeSphere(network, network->FiberList[*i].n0, node_radius); + DrawNodeSphere(network, network->FiberList[*i].n1, node_radius); + glPopAttrib(); + } + + +} +GLuint CreateFiberDisplayList(rtsFiberNetwork* network, float radius) +{ + GLuint result = glGenLists(1); + glNewList(result, GL_COMPILE); + + int num_fibers = network->FiberList.size(); + for(int f = 0; f < num_fibers; f++) + { + if(!network->isCulled(f)) + { + ExtrudeFiber(network, f, radius); + CHECK_OPENGL_ERROR + } + } + glEndList(); + return result; +} + +GLuint CreateNodeDisplayList(rtsFiberNetwork* network, float radius) +{ + GLuint result = glGenLists(1); + glNewList(result, GL_COMPILE); + DrawNodeSpheres(network, radius); + glEndList(); + return result; +} + + +void CreateFiberPathLists(float fiber_radius, float node_radius) +{ + GT_PathList = glGenLists(1); + glNewList(GT_PathList, GL_COMPILE); + + if(coreGraph.size() > 0) + { + for(unsigned int i=0; i rgb = sequenceColors[i]; + glColor3f(rgb.x, rgb.y, rgb.z); + DrawFiberSequence(goldNetwork, coreGraph[i].second, fiber_radius, node_radius); + } + } + glEndList(); + + T_PathList = glGenLists(1); + glNewList(T_PathList, GL_COMPILE); + + if(coreGraph.size() > 0) + { + for(unsigned int i=0; i rgb = sequenceColors[i]; + glColor3f(rgb.x, rgb.y, rgb.z); + DrawFiberSequence(testNetwork, coreGraph[i].first, fiber_radius, node_radius); + } + } + glEndList(); + +} +void CreateDisplayLists() +{ + if(GT_FibersList != 0) + glDeleteLists(GT_FibersList, 1); + if(T_FibersList != 0) + glDeleteLists(T_FibersList, 1); + if(GT_NodesList != 0) + glDeleteLists(GT_NodesList, 1); + if(T_NodesList != 0) + glDeleteLists(T_NodesList, 1); + if(GT_EndCaps != 0) + glDeleteLists(GT_EndCaps, 1); + if(T_EndCaps != 0) + glDeleteLists(T_EndCaps, 1); + + //create the display lists + GT_FibersList = CreateFiberDisplayList(goldNetwork, sigmaG*fiber_radius_factor); + T_FibersList = CreateFiberDisplayList(testNetwork, sigmaG*fiber_radius_factor); + GT_NodesList = CreateNodeDisplayList(goldNetwork, sigmaG*node_radius_factor); + T_NodesList = CreateNodeDisplayList(testNetwork, sigmaG*node_radius_factor); + GT_EndCaps = CreateNodeDisplayList(goldNetwork, sigmaG*fiber_radius_factor); + T_EndCaps = CreateNodeDisplayList(testNetwork, sigmaG*fiber_radius_factor); + + if(GT_PathList != 0) + glDeleteLists(GT_PathList,1); + if(T_PathList != 0) + glDeleteLists(T_PathList,1); + CreateFiberPathLists(fiber_radius_factor*sigmaG, node_radius_factor*sigmaG); +} +void RenderViewCamera() +{ + glMatrixMode(GL_PROJECTION); + glLoadIdentity(); + //compute the aspect ratio + float aspect_ratio = (float)glutGet(GLUT_WINDOW_WIDTH)/2.0/(float)glutGet(GLUT_WINDOW_HEIGHT); + gluPerspective(rts_glut_camera.getFOV(), aspect_ratio, network_span/10, network_span*10); + + //render the camera + glMatrixMode(GL_MODELVIEW); + glLoadIdentity(); + + point3D camera_position = rts_glut_camera.getPosition(); + vector3D camera_up = rts_glut_camera.getUp(); + point3D camera_lookat = rts_glut_camera.getLookAt(); + gluLookAt(camera_position.x, + camera_position.y, + camera_position.z, + camera_lookat.x, + camera_lookat.y, + camera_lookat.z, + camera_up.x, + camera_up.y, + camera_up.z); + + + //get the light positions (lights move with the camera) + vector3D up = rts_glut_camera.getUp(); + vector3D dir = rts_glut_camera.getDirection(); + vector3D side = dir.X(up); + L0_pos[0] = side.x; + L0_pos[1] = side.y; + L0_pos[2] = side.z; + + L1_pos[0] = -dir.x; + L1_pos[1] = -dir.y; + L1_pos[2] = -dir.z; + + //scale the viewport to the network + /*vector3D span = max_point - min_point; + float scale = 1.0/span.Length(); + //compute center point + point3D center = min_point + 0.5*span; + glScalef(scale, scale, scale); + glTranslatef(-center.x, -center.y, -center.z);*/ + + +} +void RecenterCamera() +{ + point3D min_point0 = goldNetwork->min_pos; + point3D min_point1 = testNetwork->min_pos; + + point3D max_point0 = goldNetwork->max_pos; + point3D max_point1 = testNetwork->max_pos; + + point3D min_point(min(min_point0.x, min_point1.x), min(min_point0.y, min_point1.y), min(min_point0.z, min_point1.z)); + point3D max_point(max(max_point0.x, max_point1.x), max(max_point0.y, max_point1.y), max(max_point0.z, max_point1.z)); + point3D center_point = min_point + 0.5*(max_point - min_point); + + + network_span = (max_point - point3D(0, 0, 0)).Length(); + rts_glut_camera.setPosition(0, 0, 3*network_span); + rts_glut_camera.LookAt(center_point, vector3D(0.0, 1.0, 0.0)); + +} +void DrawNetwork(int Display) +{ + //draw the network fibers + switch(Display) + { + case DISPLAY_GT_NETWORK: + Edge_ErrorShader.UpdateGlobalUniforms(); + Edge_ErrorShader.BeginProgram(); + glCallList(GT_FibersList); + Edge_ErrorShader.EndProgram(); + + Node_ErrorShader.UpdateGlobalUniforms(); + Node_ErrorShader.BeginProgram(); + glCallList(GT_EndCaps); + Node_ErrorShader.EndProgram(); + break; + case DISPLAY_T_NETWORK: + Edge_ErrorShader.UpdateGlobalUniforms(); + Edge_ErrorShader.BeginProgram(); + glCallList(T_FibersList); + Edge_ErrorShader.EndProgram(); + + Node_ErrorShader.UpdateGlobalUniforms(); + Node_ErrorShader.BeginProgram(); + glCallList(T_EndCaps); + Node_ErrorShader.EndProgram(); + break; + case DISPLAY_GT_GRAPH: + glColor3f(1.0, 1.0, 1.0); + Smooth_Shader.UpdateGlobalUniforms(); + Smooth_Shader.BeginProgram(); + glCallList(GT_PathList); + glColor3f(1.0, 1.0, 1.0); + glCallList(GT_FibersList); + glCallList(GT_NodesList); + Smooth_Shader.EndProgram(); + break; + case DISPLAY_T_GRAPH: + glColor3f(1.0, 1.0, 1.0); + Smooth_Shader.UpdateGlobalUniforms(); + Smooth_Shader.BeginProgram(); + glCallList(T_PathList); + glColor3f(1.0, 1.0, 1.0); + glCallList(T_FibersList); + glCallList(T_NodesList); + Smooth_Shader.EndProgram(); + break; + case DISPLAY_GT_SELECTED: + glColor3f(1.0, 1.0, 1.0); + Smooth_Shader.UpdateGlobalUniforms(); + Smooth_Shader.BeginProgram(); + glCallList(GT_FibersList); + glCallList(GT_NodesList); + + glClear(GL_DEPTH_BUFFER_BIT); + glColor3f(1.0, 0.0, 1.0); + DrawFiberSequence(goldNetwork, coreGraph[current_sequence].second, fiber_radius_factor*sigmaG, node_radius_factor*sigmaG); + Smooth_Shader.EndProgram(); + break; + case DISPLAY_T_SELECTED: + glColor3f(1.0, 1.0, 1.0); + Smooth_Shader.UpdateGlobalUniforms(); + Smooth_Shader.BeginProgram(); + glCallList(T_FibersList); + glCallList(T_NodesList); + + glClear(GL_DEPTH_BUFFER_BIT); + glColor3f(1.0, 0.0, 1.0); + DrawFiberSequence(testNetwork, coreGraph[current_sequence].first, fiber_radius_factor*sigmaG, node_radius_factor*sigmaG); + Smooth_Shader.EndProgram(); + break; + default: + break; + } + +} + + + + +void MyDisplayFunction() +{ + //glClearColor(0.95, 1.0, 0.95, 1.0); + glClearColor(1.0, 1.0, 1.0, 1.0); + glClear(GL_COLOR_BUFFER_BIT); + glClear(GL_DEPTH_BUFFER_BIT); + + //set left viewport + glMatrixMode(GL_PROJECTION); + glLoadIdentity(); + glViewport(0, 0, glutGet(GLUT_WINDOW_WIDTH)/2, glutGet(GLUT_WINDOW_HEIGHT)); + + RenderViewCamera(); + + + if(DisplayMode == DISPLAY_NETWORK) + DrawNetwork(DISPLAY_GT_NETWORK); + else if(DisplayMode == DISPLAY_GRAPH) + DrawNetwork(DISPLAY_GT_GRAPH); + else if(DisplayMode == DISPLAY_SELECTED) + DrawNetwork(DISPLAY_GT_SELECTED); + + //set right viewport + glMatrixMode(GL_PROJECTION); + glLoadIdentity(); + glViewport(glutGet(GLUT_WINDOW_WIDTH)/2, 0, glutGet(GLUT_WINDOW_WIDTH)/2, glutGet(GLUT_WINDOW_HEIGHT)); + + RenderViewCamera(); + //glClear(GL_COLOR_BUFFER_BIT); + //glClear(GL_DEPTH_BUFFER_BIT); + + if(DisplayMode == DISPLAY_NETWORK) + DrawNetwork(DISPLAY_T_NETWORK); + else if(DisplayMode == DISPLAY_GRAPH) + DrawNetwork(DISPLAY_T_GRAPH); + else if(DisplayMode == DISPLAY_SELECTED) + DrawNetwork(DISPLAY_T_SELECTED); + + + + //glutSwapBuffers(); +} + +void GlutMenuCallback(int option) +{ + if(option == NETCOMP_EXIT) + exit(1); + if(option >= DISPLAY_NETWORK && option <= DISPLAY_SELECTED) + DisplayMode = option; + if(option >= COLORMAP_ISOLUMINANT && option <= COLORMAP_RAINBOW) + makeColormap(option); + if(option == CULL_TEST_CASE) + { + //get the new threshold + cout<<"Enter new test case fiber threshold [0 1]: "; + float cull_value; + cin>>cull_value; + testNetwork->setCullValue(cull_value); + + //re-create the display lists + ComputeNetMets(); + CreateDisplayLists(); + + } + if(option == RECOMPUTE_METRIC) + { + cout<<"Please enter a sigma value: "; + cin>>sigmaG; + sigmaC = sigmaG; + + ComputeNetMets(); + CreateDisplayLists(); + } + + + + /* + switch(option) + { + case NETCOMP_EXIT: + exit(1); + break; + case NETCOMP_VIEW_GT: + DisplayMode = DISPLAY_GT_NETWORK; + break; + case NETCOMP_VIEW_T: + DisplayMode = DISPLAY_T_NETWORK; + break; + default: + break; + }*/ + +} \ No newline at end of file diff --git a/ErrorMap_Fragment.glsl b/ErrorMap_Fragment.glsl new file mode 100644 index 0000000..adb879f --- /dev/null +++ b/ErrorMap_Fragment.glsl @@ -0,0 +1,74 @@ +#extension GL_EXT_gpu_shader4 : enable + +varying vec3 fragment_normal; +//uniform vec3 L0_pos; +uniform vec3 L1_pos; +uniform sampler1D colorMap; + +vec4 HSLtoRGB(vec4 HSL) +{ + float H = HSL.x; + float S = HSL.y; + float L = HSL.z; + + float temp2; + if(L < 0.5) + temp2 = L*(1.0+S); + else + temp2 = L+S - L*S; + + float temp1 = 2.0*L - temp2; + + vec3 temp3 = vec3(H+1.0/3.0, H, H-1.0/3.0); + if(temp3.r < 0.0) temp3.r = temp3.r+1.0; + if(temp3.g < 0.0) temp3.g = temp3.g+1.0; + if(temp3.b < 0.0) temp3.b = temp3.b+1.0; + + if(temp3.r > 1.0) temp3.r = temp3.r - 1.0; + if(temp3.g > 1.0) temp3.g = temp3.g - 1.0; + if(temp3.b > 1.0) temp3.b = temp3.b - 1.0; + + vec4 result; + if(6.0*temp3.r < 1.0) result.r = temp1 +(temp2 - temp1)*6.0*temp3.r; + else if(2.0*temp3.r < 1.0) result.r = temp2; + else if(3.0*temp3.r < 2.0) result.r = temp1+(temp2-temp1)*((2.0/3.0) - temp3.r)*6.0; + else result.r = temp1; + + if(6.0*temp3.g < 1.0) result.g = temp1 +(temp2 - temp1)*6.0*temp3.g; + else if(2.0*temp3.g < 1.0) result.g = temp2; + else if(3.0*temp3.g < 2.0) result.g = temp1+(temp2-temp1)*((2.0/3.0) - temp3.g)*6.0; + else result.g = temp1; + + if(6.0*temp3.b < 1.0) result.b = temp1 +(temp2 - temp1)*6.0*temp3.b; + else if(2.0*temp3.b < 1.0) result.b = temp2; + else if(3.0*temp3.b < 2.0) result.b = temp1+(temp2-temp1)*((2.0/3.0) - temp3.b)*6.0; + else result.b = temp1; + + result.a = 0.0; + return result; +} + +float LightIntensity() +{ + //vec3 L0_pos = vec3(1.0, 0.0, 0.0); + //vec3 L1_pos = vec3(0.0, 0.0, 1.0); + //float L0 = max(dot(fragment_normal, L0_pos), 0.0); + float L1 = max(dot(fragment_normal, L1_pos), 0.0); + + //float total = L0 + L1; + //if(total > 1.0) + // total = 1.0; + + return L1; + //return 0.5; + +} + +void main(void) +{ + float error = gl_TexCoord[0].x; + + float light = LightIntensity(); + + gl_FragColor = texture1D(colorMap, error)*light; +} diff --git a/ErrorShader_Fragment.glsl b/ErrorShader_Fragment.glsl new file mode 100644 index 0000000..73d05f1 --- /dev/null +++ b/ErrorShader_Fragment.glsl @@ -0,0 +1,72 @@ +#extension GL_EXT_gpu_shader4 : enable + +varying vec3 fragment_normal; +uniform vec3 L0_pos; +uniform vec3 L1_pos; + +vec4 HSLtoRGB(vec4 HSL) +{ + float H = HSL.x; + float S = HSL.y; + float L = HSL.z; + + float temp2; + if(L < 0.5) + temp2 = L*(1.0+S); + else + temp2 = L+S - L*S; + + float temp1 = 2.0*L - temp2; + + vec3 temp3 = vec3(H+1.0/3.0, H, H-1.0/3.0); + if(temp3.r < 0.0) temp3.r = temp3.r+1.0; + if(temp3.g < 0.0) temp3.g = temp3.g+1.0; + if(temp3.b < 0.0) temp3.b = temp3.b+1.0; + + if(temp3.r > 1.0) temp3.r = temp3.r - 1.0; + if(temp3.g > 1.0) temp3.g = temp3.g - 1.0; + if(temp3.b > 1.0) temp3.b = temp3.b - 1.0; + + vec4 result; + if(6.0*temp3.r < 1.0) result.r = temp1 +(temp2 - temp1)*6.0*temp3.r; + else if(2.0*temp3.r < 1.0) result.r = temp2; + else if(3.0*temp3.r < 2.0) result.r = temp1+(temp2-temp1)*((2.0/3.0) - temp3.r)*6.0; + else result.r = temp1; + + if(6.0*temp3.g < 1.0) result.g = temp1 +(temp2 - temp1)*6.0*temp3.g; + else if(2.0*temp3.g < 1.0) result.g = temp2; + else if(3.0*temp3.g < 2.0) result.g = temp1+(temp2-temp1)*((2.0/3.0) - temp3.g)*6.0; + else result.g = temp1; + + if(6.0*temp3.b < 1.0) result.b = temp1 +(temp2 - temp1)*6.0*temp3.b; + else if(2.0*temp3.b < 1.0) result.b = temp2; + else if(3.0*temp3.b < 2.0) result.b = temp1+(temp2-temp1)*((2.0/3.0) - temp3.b)*6.0; + else result.b = temp1; + + result.a = 0.0; + return result; +} + +float LightIntensity() +{ + //vec3 L0_pos = vec3(1.0, 0.0, 0.0); + //vec3 L1_pos = vec3(0.0, 0.0, 1.0); + float L0 = 0.3*max(dot(fragment_normal, L0_pos), 0.0); + float L1 = 0.3*max(dot(fragment_normal, L1_pos), 0.0); + + return L0 + L1; + //return 0.5; + +} + +void main(void) +{ + float error = (1.0 - gl_TexCoord[0].x)*(240.0/360.0); + float light = LightIntensity(); + + //convert from HSL to RGB + vec4 HSL = vec4(error, 1.0, light, 1.0); + vec4 RGB = HSLtoRGB(HSL); + + gl_FragColor = vec4(RGB.r, RGB.g, RGB.b, 1.0); +} \ No newline at end of file diff --git a/FindANN.cmake b/FindANN.cmake new file mode 100644 index 0000000..c4b5649 --- /dev/null +++ b/FindANN.cmake @@ -0,0 +1,42 @@ +# - Try to find ANN +# Once done this will define +# +# ANN_FOUND - system has ANN +# ANN_INCLUDE_DIR - the ANN include directory +# ANN_LIBRARY - Link these to use ANN +# + +IF (ANN_INCLUDE_DIRS) + # Already in cache, be silent + SET(ANN_FIND_QUIETLY TRUE) +ENDIF (ANN_INCLUDE_DIRS) + +FIND_PATH( ANN_INCLUDE_DIR ANN/ANN.h + PATHS "/usr/include" "C:/libs/ANN/include") + +if( WIN32 ) + + set(ANN_LIBRARY $ENV{ANN_PATH}\\ANN.lib) + set(ANN_INCLUDE_DIR $ENV{ANN_PATH}) + + + # Store the library dir. May be used for linking to dll! + # GET_FILENAME_COMPONENT( ANN_LIBRARY_DIR ${ANN_LIBRARY} PATH ) + + find_package_handle_standard_args(ANN DEFAULT_MSG ANN_INCLUDE_DIR) + +else (WIN32) + + FIND_LIBRARY( ANN_LIBRARY + NAMES ann ANN + PATHS /lib /usr/lib /usr/lib64 /usr/local/lib ) + +endif( WIN32) + + +IF (ANN_INCLUDE_DIR AND ANN_LIBRARY) + SET(ANN_FOUND TRUE) +ELSE (ANN_INCLUDE_DIR AND ANN_LIBRARY) + SET( ANN_FOUND FALSE ) +ENDIF (ANN_INCLUDE_DIR AND ANN_LIBRARY) + diff --git a/FindGLEW.cmake b/FindGLEW.cmake new file mode 100644 index 0000000..a107814 --- /dev/null +++ b/FindGLEW.cmake @@ -0,0 +1,49 @@ +# +# Try to find GLEW library and include path. +# Once done this will define +# +# GLEW_FOUND +# GLEW_INCLUDE_PATH +# GLEW_LIBRARY +# + +IF (WIN32) + FIND_PATH( GLEW_INCLUDE_PATH GL/glew.h + $ENV{PROGRAMFILES}/GLEW/include + $ENV{GLEW_PATH}/include + DOC "The directory where GL/glew.h resides") + FIND_LIBRARY( GLEW_LIBRARY + NAMES glew GLEW glew32 glew32s + PATHS + $ENV{GLEW_PATH}/lib/Release/Win32 + DOC "The GLEW library") +ELSE (WIN32) + FIND_PATH( GLEW_INCLUDE_PATH GL/glew.h + /usr/include + /usr/local/include + /sw/include + /opt/local/include + DOC "The directory where GL/glew.h resides") + FIND_LIBRARY( GLEW_LIBRARY + NAMES GLEW glew + PATHS + /usr/lib64 + /usr/lib + /usr/local/lib64 + /usr/local/lib + /sw/lib + /opt/local/lib + DOC "The GLEW library") +ENDIF (WIN32) + +IF (GLEW_INCLUDE_PATH) + SET( GLEW_FOUND 1 CACHE STRING "Set to 1 if GLEW is found, 0 otherwise") +ELSE (GLEW_INCLUDE_PATH) + SET( GLEW_FOUND 0 CACHE STRING "Set to 1 if GLEW is found, 0 otherwise") +ENDIF (GLEW_INCLUDE_PATH) + +MARK_AS_ADVANCED( + GLEW_FOUND + GLEW_INCLUDE_PATH + GLEW_LIBRARY +) diff --git a/FindGLUT.cmake b/FindGLUT.cmake new file mode 100644 index 0000000..c4b5649 --- /dev/null +++ b/FindGLUT.cmake @@ -0,0 +1,42 @@ +# - Try to find ANN +# Once done this will define +# +# ANN_FOUND - system has ANN +# ANN_INCLUDE_DIR - the ANN include directory +# ANN_LIBRARY - Link these to use ANN +# + +IF (ANN_INCLUDE_DIRS) + # Already in cache, be silent + SET(ANN_FIND_QUIETLY TRUE) +ENDIF (ANN_INCLUDE_DIRS) + +FIND_PATH( ANN_INCLUDE_DIR ANN/ANN.h + PATHS "/usr/include" "C:/libs/ANN/include") + +if( WIN32 ) + + set(ANN_LIBRARY $ENV{ANN_PATH}\\ANN.lib) + set(ANN_INCLUDE_DIR $ENV{ANN_PATH}) + + + # Store the library dir. May be used for linking to dll! + # GET_FILENAME_COMPONENT( ANN_LIBRARY_DIR ${ANN_LIBRARY} PATH ) + + find_package_handle_standard_args(ANN DEFAULT_MSG ANN_INCLUDE_DIR) + +else (WIN32) + + FIND_LIBRARY( ANN_LIBRARY + NAMES ann ANN + PATHS /lib /usr/lib /usr/lib64 /usr/local/lib ) + +endif( WIN32) + + +IF (ANN_INCLUDE_DIR AND ANN_LIBRARY) + SET(ANN_FOUND TRUE) +ELSE (ANN_INCLUDE_DIR AND ANN_LIBRARY) + SET( ANN_FOUND FALSE ) +ENDIF (ANN_INCLUDE_DIR AND ANN_LIBRARY) + diff --git a/GroundTruth.swc b/GroundTruth.swc new file mode 100644 index 0000000..c2011b6 --- /dev/null +++ b/GroundTruth.swc @@ -0,0 +1,243 @@ +# Original neuron tracing performed in Neuromantic +# which was programmed by Darren Myatt (d.r.myatt@reading.ac.uk) +# See homepage at http://www.rdg.ac.uk/Neuromantic +# +# NEUROMANTIC V1.6.3 (1/15/2012 5:22:06 PM): Saved to GroundTruth.swc +1 1 1342.5 1069.5 0 6 -1 +2 1 1351.5 1063.5 0 6 1 +3 1 1361.5 1058.5 0 6 2 +4 1 1369.5 1052.5 0 4 3 +5 1 1378.5 1047.5 0 3.5 4 +6 1 1386.5 1041.5 0 3 5 +7 1 1394.5 1034.5 0 3 6 +8 1 1402.5 1027.5 0 3.5 7 +9 1 1412.5 1026.5 0 4 8 +10 1 1420.5 1018.5 0 4 9 +11 1 1428.5 1011.5 0 4 10 +12 1 1437.5 1005.5 0 4.5 11 +13 1 1447.5 1001.5 0 4.5 12 +14 1 1456.5 995.5 0 3 13 +15 1 1465.5 990.5 0 3 14 +16 1 1474.5 985.5 0 2.5 15 +17 1 1482.5 979.5 0 2.5 16 +18 1 1492.5 979.5 0 2.5 17 +19 1 1501.5 974.5 0 2 18 +20 1 1510.5 969.5 0 2 19 +21 1 1520.5 969.5 0 2 20 +22 1 1530.5 969.5 0 2.5 21 +23 1 1538.5 962.5 0 2.5 22 +24 1 1548.5 959.5 0 2.5 23 +25 1 1558.5 959.5 0 2 24 +26 1 1568.5 959.5 0 2 25 +27 1 1578.5 958.5 0 2 26 +28 1 1588.5 954.5 0 2.5 27 +29 1 1598.5 952.5 0 3.5 28 +30 1 1618.5 943.5 0 2 29 +31 1 1628.5 941.5 0 2 30 +32 1 1636.5 935.5 0 2 31 +33 1 1645.5 930.5 0 2 32 +34 1 1655.5 926.5 0 2 33 +35 1 1665.5 922.5 0 2 34 +36 1 1675.5 918.5 0 2 35 +37 1 1685.5 913.5 0 2 36 +38 1 1695.5 909.5 0 2 37 +39 1 1705.5 907.5 0 2 38 +40 1 1714.5 902.5 0 2 39 +41 1 1722.5 896.5 0 2 40 +42 1 1730.5 888.5 0 2 41 +43 1 1735.5 878.5 0 2 42 +44 1 1742.5 870.5 0 2 43 +45 1 1748.5 862.5 0 2 44 +46 1 1753.5 853.5 0 2 45 +47 1 1761.5 845.5 0 2 46 +48 1 1769.5 837.5 0 2.5 47 +49 1 1777.5 829.5 0 2.5 48 +50 1 1786.5 824.5 0 2.5 49 +51 1 1796.5 822.5 0 2.5 50 +52 1 1809.5 826.5 0 2 51 +53 1 1819.5 828.5 0 2 52 +54 1 1829.5 831.5 0 1.5 53 +55 1 1837.5 838.5 0 1.5 54 +56 1 1839.5 848.5 0 2.5 55 +57 1 1749.5 882.5 0 2 42 +58 1 1758.5 876.5 0 2 57 +59 1 1766.5 870.5 0 2 58 +60 1 1774.5 862.5 0 2.5 59 +61 1 1784.5 861.5 0 2.5 60 +62 1 1794.5 861.5 0 3 61 +63 1 1803.5 866.5 0 2 62 +64 1 1813.5 870.5 0 2 63 +65 1 1608.5 957.5 0 2 28 +66 1 1618.5 959.5 0 1.5 65 +67 1 1628.5 960.5 0 1.5 66 +68 1 1638.5 960.5 0 1.5 67 +69 1 1648.5 960.5 0 2 68 +70 1 1658.5 960.5 0 2.5 69 +71 1 1668.5 960.5 0 2.5 70 +72 1 1678.5 958.5 0 2 71 +73 1 1687.5 952.5 0 2 72 +74 1 1696.5 947.5 0 2 73 +75 1 1706.5 944.5 0 2 74 +76 1 1714.5 938.5 0 2.5 75 +77 1 1723.5 933.5 0 3 76 +78 1 1733.5 929.5 0 3 77 +79 1 1743.5 928.5 0 2.5 78 +80 1 1753.5 926.5 0 2.5 79 +81 1 1763.5 923.5 0 2.5 80 +82 1 1772.5 917.5 0 2.5 81 +83 1 1781.5 912.5 0 2 82 +84 1 1349.5 1054.5 0 4 1 +85 1 1358.5 1049.5 0 4 84 +86 1 1366.5 1041.5 0 3 85 +87 1 1371.5 1032.5 0 3 86 +88 1 1379.5 1024.5 0 3 87 +89 1 1387.5 1016.5 0 2.5 88 +90 1 1395.5 1008.5 0 2.5 89 +91 1 1403.5 1000.5 0 2.5 90 +92 1 1410.5 992.5 0 2.5 91 +93 1 1417.5 984.5 0 2 92 +94 1 1425.5 976.5 0 2.5 93 +95 1 1433.5 968.5 0 3 94 +96 1 1441.5 960.5 0 3 95 +97 1 1449.5 952.5 0 3 96 +98 1 1458.5 947.5 0 2.5 97 +99 1 1468.5 944.5 0 1.5 98 +100 1 1476.5 938.5 0 1.5 99 +101 1 1486.5 934.5 0 2 100 +102 1 1494.5 927.5 0 2 101 +103 1 1502.5 919.5 0 2 102 +104 1 1511.5 914.5 0 2 103 +105 1 1519.5 908.5 0 2 104 +106 1 1527.5 901.5 0 2 105 +107 1 1535.5 895.5 0 2 106 +108 1 1543.5 887.5 0 2.5 107 +109 1 1551.5 880.5 0 2.5 108 +110 1 1559.5 872.5 0 2.5 109 +111 1 1566.5 864.5 0 2 110 +112 1 1573.5 856.5 0 2 111 +113 1 1581.5 848.5 0 2.5 112 +114 1 1590.5 843.5 0 3.5 113 +115 1 1459.5 931.5 0 3 97 +116 1 1467.5 923.5 0 2 115 +117 1 1475.5 915.5 0 2 116 +118 1 1479.5 905.5 0 2 117 +119 1 1485.5 897.5 0 2.5 118 +120 1 1491.5 888.5 0 2.5 119 +121 1 1499.5 880.5 0 2.5 120 +122 1 1507.5 872.5 0 2.5 121 +123 1 1515.5 864.5 0 2.5 122 +124 1 1521.5 856.5 0 2 123 +125 1 1527.5 847.5 0 2.5 124 +126 1 1529.5 837.5 0 2.5 125 +127 1 1530.5 827.5 0 2.5 126 +128 1 1534.5 817.5 0 2.5 127 +129 1 1538.5 807.5 0 3 128 +130 1 1545.5 799.5 0 3 129 +131 1 1553.5 791.5 0 2.5 130 +132 1 1561.5 785.5 0 2.5 131 +133 1 1571.5 781.5 0 2 132 +134 1 1581.5 778.5 0 2 133 +135 1 1591.5 776.5 0 2 134 +136 1 1601.5 773.5 0 2 135 +137 1 1611.5 769.5 0 2 136 +138 1 1621.5 769.5 0 2.5 137 +139 1 1635.5 762.5 0 2 138 +140 1 1645.5 759.5 0 2 139 +141 1 1654.5 754.5 0 2.5 140 +142 1 1662.5 747.5 0 2 141 +143 1 1671.5 741.5 0 2 142 +144 1 1680.5 735.5 0 2 143 +145 1 1688.5 728.5 0 2 144 +146 1 1697.5 723.5 0 2 145 +147 1 1707.5 720.5 0 2 146 +148 1 1716.5 715.5 0 1.5 147 +149 1 1726.5 711.5 0 2 148 +150 1 1736.5 709.5 0 1.5 149 +151 1 1746.5 708.5 0 2 150 +152 1 1756.5 705.5 0 3 151 +153 1 1765.5 700.5 0 12 152 +154 1 1636.5 753.5 0 6 138 +155 1 1646.5 748.5 0 5.5 154 +156 1 1656.5 743.5 0 5.5 155 +157 1 1666.5 736.5 0 2 156 +158 1 1672.5 728.5 0 2 157 +159 1 1677.5 718.5 0 2 158 +160 1 1683.5 709.5 0 2 159 +161 1 1689.5 701.5 0 2 160 +162 1 1695.5 692.5 0 2 161 +163 1 1701.5 683.5 0 1.5 162 +164 1 1704.5 673.5 0 1.5 163 +165 1 1708.5 663.5 0 1.5 164 +166 1 1711.5 653.5 0 1.5 165 +167 1 1342.5 1038.5 0 3.5 1 +168 1 1339.5 1028.5 0 4 167 +169 1 1339.5 1018.5 0 4 168 +170 1 1339.5 1008.5 0 4 169 +171 1 1339.5 998.5 0 3 170 +172 1 1339.5 988.5 0 3 171 +173 1 1339.5 978.5 0 2.5 172 +174 1 1339.5 968.5 0 2.5 173 +175 1 1319.5 1084.5 0 4 1 +176 1 1310.5 1090.5 0 4 175 +177 1 1302.5 1098.5 0 4 176 +178 1 1294.5 1106.5 0 4 177 +179 1 1286.5 1114.5 0 4 178 +180 1 1278.5 1121.5 0 4 179 +181 1 1270.5 1128.5 0 4.5 180 +182 1 1262.5 1136.5 0 4 181 +183 1 1254.5 1144.5 0 4.5 182 +184 1 1246.5 1152.5 0 4.5 183 +185 1 1238.5 1160.5 0 4.5 184 +186 1 1231.5 1168.5 0 3.5 185 +187 1 1225.5 1176.5 0 3.5 186 +188 1 1221.5 1186.5 0 3.5 187 +189 1 1216.5 1196.5 0 3.5 188 +190 1 1210.5 1204.5 0 5 189 +191 1 1202.5 1212.5 0 5 190 +192 1 1197.5 1222.5 0 2.5 191 +193 1 1193.5 1232.5 0 2.5 192 +194 1 1188.5 1241.5 0 2.5 193 +195 1 1183.5 1251.5 0 2.5 194 +196 1 1177.5 1259.5 0 2.5 195 +197 1 1171.5 1268.5 0 2.5 196 +198 1 1165.5 1276.5 0 2.5 197 +199 1 1157.5 1284.5 0 2 198 +200 1 1149.5 1292.5 0 2 199 +201 1 1144.5 1302.5 0 2 200 +202 1 1133.5 1317.5 0 2 201 +203 1 1127.5 1326.5 0 2 202 +204 1 1122.5 1335.5 0 2.5 203 +205 1 1118.5 1345.5 0 2.5 204 +206 1 1115.5 1355.5 0 2 205 +207 1 1114.5 1365.5 0 2 206 +208 1 1114.5 1375.5 0 2 207 +209 1 1162.5 1298.5 0 2 197 +210 1 1158.5 1308.5 0 2 209 +211 1 1155.5 1318.5 0 2 210 +212 1 1153.5 1328.5 0 2 211 +213 1 1152.5 1338.5 0 2 212 +214 1 1151.5 1348.5 0 2 213 +215 1 1150.5 1358.5 0 2 214 +216 1 1149.5 1368.5 0 2 215 +217 1 1147.5 1378.5 0 2 216 +218 1 1144.5 1388.5 0 2 217 +219 1 1178.5 1232.5 0 2.5 191 +220 1 1172.5 1240.5 0 2.5 219 +221 1 1166.5 1248.5 0 2.5 220 +222 1 1159.5 1256.5 0 2.5 221 +223 1 1151.5 1264.5 0 2.5 222 +224 1 1142.5 1269.5 0 2.5 223 +225 1 1133.5 1274.5 0 2 224 +226 1 1125.5 1281.5 0 2 225 +227 1 1117.5 1289.5 0 2 226 +228 1 1109.5 1297.5 0 2 227 +229 1 1101.5 1303.5 0 2 228 +230 1 1093.5 1310.5 0 2 229 +231 1 1084.5 1315.5 0 2 230 +232 1 1076.5 1322.5 0 2 231 +233 1 1067.5 1327.5 0 2.5 232 +234 1 1059.5 1334.5 0 2.5 233 +235 1 1051.5 1342.5 0 2.5 234 +236 1 1043.5 1350.5 0 2 235 +237 1 1037.5 1359.5 0 2.5 236 +238 1 1031.5 1368.5 0 4.5 237 diff --git a/GroundTruthCopy.swc b/GroundTruthCopy.swc new file mode 100644 index 0000000..c2011b6 --- /dev/null +++ b/GroundTruthCopy.swc @@ -0,0 +1,243 @@ +# Original neuron tracing performed in Neuromantic +# which was programmed by Darren Myatt (d.r.myatt@reading.ac.uk) +# See homepage at http://www.rdg.ac.uk/Neuromantic +# +# NEUROMANTIC V1.6.3 (1/15/2012 5:22:06 PM): Saved to GroundTruth.swc +1 1 1342.5 1069.5 0 6 -1 +2 1 1351.5 1063.5 0 6 1 +3 1 1361.5 1058.5 0 6 2 +4 1 1369.5 1052.5 0 4 3 +5 1 1378.5 1047.5 0 3.5 4 +6 1 1386.5 1041.5 0 3 5 +7 1 1394.5 1034.5 0 3 6 +8 1 1402.5 1027.5 0 3.5 7 +9 1 1412.5 1026.5 0 4 8 +10 1 1420.5 1018.5 0 4 9 +11 1 1428.5 1011.5 0 4 10 +12 1 1437.5 1005.5 0 4.5 11 +13 1 1447.5 1001.5 0 4.5 12 +14 1 1456.5 995.5 0 3 13 +15 1 1465.5 990.5 0 3 14 +16 1 1474.5 985.5 0 2.5 15 +17 1 1482.5 979.5 0 2.5 16 +18 1 1492.5 979.5 0 2.5 17 +19 1 1501.5 974.5 0 2 18 +20 1 1510.5 969.5 0 2 19 +21 1 1520.5 969.5 0 2 20 +22 1 1530.5 969.5 0 2.5 21 +23 1 1538.5 962.5 0 2.5 22 +24 1 1548.5 959.5 0 2.5 23 +25 1 1558.5 959.5 0 2 24 +26 1 1568.5 959.5 0 2 25 +27 1 1578.5 958.5 0 2 26 +28 1 1588.5 954.5 0 2.5 27 +29 1 1598.5 952.5 0 3.5 28 +30 1 1618.5 943.5 0 2 29 +31 1 1628.5 941.5 0 2 30 +32 1 1636.5 935.5 0 2 31 +33 1 1645.5 930.5 0 2 32 +34 1 1655.5 926.5 0 2 33 +35 1 1665.5 922.5 0 2 34 +36 1 1675.5 918.5 0 2 35 +37 1 1685.5 913.5 0 2 36 +38 1 1695.5 909.5 0 2 37 +39 1 1705.5 907.5 0 2 38 +40 1 1714.5 902.5 0 2 39 +41 1 1722.5 896.5 0 2 40 +42 1 1730.5 888.5 0 2 41 +43 1 1735.5 878.5 0 2 42 +44 1 1742.5 870.5 0 2 43 +45 1 1748.5 862.5 0 2 44 +46 1 1753.5 853.5 0 2 45 +47 1 1761.5 845.5 0 2 46 +48 1 1769.5 837.5 0 2.5 47 +49 1 1777.5 829.5 0 2.5 48 +50 1 1786.5 824.5 0 2.5 49 +51 1 1796.5 822.5 0 2.5 50 +52 1 1809.5 826.5 0 2 51 +53 1 1819.5 828.5 0 2 52 +54 1 1829.5 831.5 0 1.5 53 +55 1 1837.5 838.5 0 1.5 54 +56 1 1839.5 848.5 0 2.5 55 +57 1 1749.5 882.5 0 2 42 +58 1 1758.5 876.5 0 2 57 +59 1 1766.5 870.5 0 2 58 +60 1 1774.5 862.5 0 2.5 59 +61 1 1784.5 861.5 0 2.5 60 +62 1 1794.5 861.5 0 3 61 +63 1 1803.5 866.5 0 2 62 +64 1 1813.5 870.5 0 2 63 +65 1 1608.5 957.5 0 2 28 +66 1 1618.5 959.5 0 1.5 65 +67 1 1628.5 960.5 0 1.5 66 +68 1 1638.5 960.5 0 1.5 67 +69 1 1648.5 960.5 0 2 68 +70 1 1658.5 960.5 0 2.5 69 +71 1 1668.5 960.5 0 2.5 70 +72 1 1678.5 958.5 0 2 71 +73 1 1687.5 952.5 0 2 72 +74 1 1696.5 947.5 0 2 73 +75 1 1706.5 944.5 0 2 74 +76 1 1714.5 938.5 0 2.5 75 +77 1 1723.5 933.5 0 3 76 +78 1 1733.5 929.5 0 3 77 +79 1 1743.5 928.5 0 2.5 78 +80 1 1753.5 926.5 0 2.5 79 +81 1 1763.5 923.5 0 2.5 80 +82 1 1772.5 917.5 0 2.5 81 +83 1 1781.5 912.5 0 2 82 +84 1 1349.5 1054.5 0 4 1 +85 1 1358.5 1049.5 0 4 84 +86 1 1366.5 1041.5 0 3 85 +87 1 1371.5 1032.5 0 3 86 +88 1 1379.5 1024.5 0 3 87 +89 1 1387.5 1016.5 0 2.5 88 +90 1 1395.5 1008.5 0 2.5 89 +91 1 1403.5 1000.5 0 2.5 90 +92 1 1410.5 992.5 0 2.5 91 +93 1 1417.5 984.5 0 2 92 +94 1 1425.5 976.5 0 2.5 93 +95 1 1433.5 968.5 0 3 94 +96 1 1441.5 960.5 0 3 95 +97 1 1449.5 952.5 0 3 96 +98 1 1458.5 947.5 0 2.5 97 +99 1 1468.5 944.5 0 1.5 98 +100 1 1476.5 938.5 0 1.5 99 +101 1 1486.5 934.5 0 2 100 +102 1 1494.5 927.5 0 2 101 +103 1 1502.5 919.5 0 2 102 +104 1 1511.5 914.5 0 2 103 +105 1 1519.5 908.5 0 2 104 +106 1 1527.5 901.5 0 2 105 +107 1 1535.5 895.5 0 2 106 +108 1 1543.5 887.5 0 2.5 107 +109 1 1551.5 880.5 0 2.5 108 +110 1 1559.5 872.5 0 2.5 109 +111 1 1566.5 864.5 0 2 110 +112 1 1573.5 856.5 0 2 111 +113 1 1581.5 848.5 0 2.5 112 +114 1 1590.5 843.5 0 3.5 113 +115 1 1459.5 931.5 0 3 97 +116 1 1467.5 923.5 0 2 115 +117 1 1475.5 915.5 0 2 116 +118 1 1479.5 905.5 0 2 117 +119 1 1485.5 897.5 0 2.5 118 +120 1 1491.5 888.5 0 2.5 119 +121 1 1499.5 880.5 0 2.5 120 +122 1 1507.5 872.5 0 2.5 121 +123 1 1515.5 864.5 0 2.5 122 +124 1 1521.5 856.5 0 2 123 +125 1 1527.5 847.5 0 2.5 124 +126 1 1529.5 837.5 0 2.5 125 +127 1 1530.5 827.5 0 2.5 126 +128 1 1534.5 817.5 0 2.5 127 +129 1 1538.5 807.5 0 3 128 +130 1 1545.5 799.5 0 3 129 +131 1 1553.5 791.5 0 2.5 130 +132 1 1561.5 785.5 0 2.5 131 +133 1 1571.5 781.5 0 2 132 +134 1 1581.5 778.5 0 2 133 +135 1 1591.5 776.5 0 2 134 +136 1 1601.5 773.5 0 2 135 +137 1 1611.5 769.5 0 2 136 +138 1 1621.5 769.5 0 2.5 137 +139 1 1635.5 762.5 0 2 138 +140 1 1645.5 759.5 0 2 139 +141 1 1654.5 754.5 0 2.5 140 +142 1 1662.5 747.5 0 2 141 +143 1 1671.5 741.5 0 2 142 +144 1 1680.5 735.5 0 2 143 +145 1 1688.5 728.5 0 2 144 +146 1 1697.5 723.5 0 2 145 +147 1 1707.5 720.5 0 2 146 +148 1 1716.5 715.5 0 1.5 147 +149 1 1726.5 711.5 0 2 148 +150 1 1736.5 709.5 0 1.5 149 +151 1 1746.5 708.5 0 2 150 +152 1 1756.5 705.5 0 3 151 +153 1 1765.5 700.5 0 12 152 +154 1 1636.5 753.5 0 6 138 +155 1 1646.5 748.5 0 5.5 154 +156 1 1656.5 743.5 0 5.5 155 +157 1 1666.5 736.5 0 2 156 +158 1 1672.5 728.5 0 2 157 +159 1 1677.5 718.5 0 2 158 +160 1 1683.5 709.5 0 2 159 +161 1 1689.5 701.5 0 2 160 +162 1 1695.5 692.5 0 2 161 +163 1 1701.5 683.5 0 1.5 162 +164 1 1704.5 673.5 0 1.5 163 +165 1 1708.5 663.5 0 1.5 164 +166 1 1711.5 653.5 0 1.5 165 +167 1 1342.5 1038.5 0 3.5 1 +168 1 1339.5 1028.5 0 4 167 +169 1 1339.5 1018.5 0 4 168 +170 1 1339.5 1008.5 0 4 169 +171 1 1339.5 998.5 0 3 170 +172 1 1339.5 988.5 0 3 171 +173 1 1339.5 978.5 0 2.5 172 +174 1 1339.5 968.5 0 2.5 173 +175 1 1319.5 1084.5 0 4 1 +176 1 1310.5 1090.5 0 4 175 +177 1 1302.5 1098.5 0 4 176 +178 1 1294.5 1106.5 0 4 177 +179 1 1286.5 1114.5 0 4 178 +180 1 1278.5 1121.5 0 4 179 +181 1 1270.5 1128.5 0 4.5 180 +182 1 1262.5 1136.5 0 4 181 +183 1 1254.5 1144.5 0 4.5 182 +184 1 1246.5 1152.5 0 4.5 183 +185 1 1238.5 1160.5 0 4.5 184 +186 1 1231.5 1168.5 0 3.5 185 +187 1 1225.5 1176.5 0 3.5 186 +188 1 1221.5 1186.5 0 3.5 187 +189 1 1216.5 1196.5 0 3.5 188 +190 1 1210.5 1204.5 0 5 189 +191 1 1202.5 1212.5 0 5 190 +192 1 1197.5 1222.5 0 2.5 191 +193 1 1193.5 1232.5 0 2.5 192 +194 1 1188.5 1241.5 0 2.5 193 +195 1 1183.5 1251.5 0 2.5 194 +196 1 1177.5 1259.5 0 2.5 195 +197 1 1171.5 1268.5 0 2.5 196 +198 1 1165.5 1276.5 0 2.5 197 +199 1 1157.5 1284.5 0 2 198 +200 1 1149.5 1292.5 0 2 199 +201 1 1144.5 1302.5 0 2 200 +202 1 1133.5 1317.5 0 2 201 +203 1 1127.5 1326.5 0 2 202 +204 1 1122.5 1335.5 0 2.5 203 +205 1 1118.5 1345.5 0 2.5 204 +206 1 1115.5 1355.5 0 2 205 +207 1 1114.5 1365.5 0 2 206 +208 1 1114.5 1375.5 0 2 207 +209 1 1162.5 1298.5 0 2 197 +210 1 1158.5 1308.5 0 2 209 +211 1 1155.5 1318.5 0 2 210 +212 1 1153.5 1328.5 0 2 211 +213 1 1152.5 1338.5 0 2 212 +214 1 1151.5 1348.5 0 2 213 +215 1 1150.5 1358.5 0 2 214 +216 1 1149.5 1368.5 0 2 215 +217 1 1147.5 1378.5 0 2 216 +218 1 1144.5 1388.5 0 2 217 +219 1 1178.5 1232.5 0 2.5 191 +220 1 1172.5 1240.5 0 2.5 219 +221 1 1166.5 1248.5 0 2.5 220 +222 1 1159.5 1256.5 0 2.5 221 +223 1 1151.5 1264.5 0 2.5 222 +224 1 1142.5 1269.5 0 2.5 223 +225 1 1133.5 1274.5 0 2 224 +226 1 1125.5 1281.5 0 2 225 +227 1 1117.5 1289.5 0 2 226 +228 1 1109.5 1297.5 0 2 227 +229 1 1101.5 1303.5 0 2 228 +230 1 1093.5 1310.5 0 2 229 +231 1 1084.5 1315.5 0 2 230 +232 1 1076.5 1322.5 0 2 231 +233 1 1067.5 1327.5 0 2.5 232 +234 1 1059.5 1334.5 0 2.5 233 +235 1 1051.5 1342.5 0 2.5 234 +236 1 1043.5 1350.5 0 2 235 +237 1 1037.5 1359.5 0 2.5 236 +238 1 1031.5 1368.5 0 4.5 237 diff --git a/ImplicitCalculations.h b/ImplicitCalculations.h new file mode 100644 index 0000000..2a5361c --- /dev/null +++ b/ImplicitCalculations.h @@ -0,0 +1,436 @@ +#include +#include +#include +//#include "DrawingFunctions.h" +#include "rtsFiberNetwork.h" +#include +#include +#include +#include +#include + +/*************************************************** +New Stuff +***************************************************/ +//networks for analysis +rtsFiberNetwork* goldNetwork; +rtsFiberNetwork* testNetwork; +rtsFiberNetwork* focusNetwork; + + +//Implicit network representations +rtsDTGrid3D* goldGrid; +rtsDTGrid3D* testGrid; +rtsDTGrid3D* focusGrid; + +//lexicographic lists, contains LexIDs for each sampled point as well as identifiers to the associated fiber +struct LexID +{ + unsigned long ID; + unsigned int fiber; + + bool operator<(LexID &rhs) + { + if(ID < rhs.ID) + return true; + return false; + } + bool operator==(LexID &rhs) + { + if(ID == rhs.ID) + return true; + return false; + } +}; +list goldLexList; +list testLexList; + +//standard deviation of the implicit envelope +float STD; + +//Size of each grid voxel. This value is defined in Network Space. +float VOXEL_SIZE; + +//Theoretical size of the grid. Used to compute the lexicographic identifiers. +vector3D GRID_DIM; + +//Transformation matrix that converts Network Space coordinates into Grid Space coordinates. +matrix4x4 toGRID; + +//Dilation parameter describes the radius of the envelope around the network (in voxels) +int DILATION; + +void ComputeGridProperties() +{ + /*Computes the following grid properties: + VOXEL_SIZE = the size (along one side) of each voxel "bin" + GRID_DIM = the theoretical dimensions of the grid used to represent the network + toGRID = transformation used to convert a point on the network into "grid space" + */ + + //pick a voxel size based on the standard daviation of the envelope that captures the shape of the Gaussian + //lets say that a voxel is 1 standard-deviation wide + VOXEL_SIZE = STD*1.0; + //dilate the DT Grid to 3 standard deviations + DILATION = (3.0*STD)/VOXEL_SIZE; + + + //find the minimum and maximum points in both networks to get the theoretical grid size + //find the bounding box for both networks + point3D net_min; + net_min.x = min(goldNetwork->min_pos.x, testNetwork->min_pos.x); + net_min.y = min(goldNetwork->min_pos.y, testNetwork->min_pos.y); + net_min.z = min(goldNetwork->min_pos.z, testNetwork->min_pos.z); + + point3D net_max; + net_max.x = max(goldNetwork->max_pos.x, testNetwork->max_pos.x); + net_max.y = max(goldNetwork->max_pos.y, testNetwork->max_pos.y); + net_max.z = max(goldNetwork->max_pos.z, testNetwork->max_pos.z); + + //compute the spatial size of the network + vector3D net_size = net_max - net_min; + + //compute the theoretical size of the grid + GRID_DIM.x = ceil(net_size.x / VOXEL_SIZE); + GRID_DIM.y = ceil(net_size.y / VOXEL_SIZE); + GRID_DIM.z = ceil(net_size.z / VOXEL_SIZE); + + //compute the transformation from Network Space to Grid Space + matrix4x4 scale; + scale. SetIdentity(); + scale.SetScale(1.0/VOXEL_SIZE, 1.0/VOXEL_SIZE, 1.0/VOXEL_SIZE); + + matrix4x4 trans; + trans.SetIdentity(); + trans.SetTranslation(-net_min.x, -net_min.y, -net_min.z); + + toGRID.SetIdentity(); + toGRID = scale*trans; + + cout<<"Voxel Size: "< p0, point3D p1, list* destList) +{ + /*Resamples a segment and stores the resulting points as + lexicographic values in LexicographicList. + */ + + //transform the Network Space coordinates to Grid Space + point3D grid_p0 = toGRID*p0; + point3D grid_p1 = toGRID*p1; + + //find the direction of travel + vector3D v = grid_p1 - grid_p0; + + //set the step size to the voxel size + int length = (int)v.Length(); + v.Normalize(); + + //start at p0, continue until p1 + point3D p = grid_p0; + + int l; + LexID lex; + + + for(l=0; l<=length; l++) + { + lex.ID = (unsigned long)p.x*GRID_DIM.y*GRID_DIM.z + (unsigned long)p.y*GRID_DIM.z + (unsigned long)p.z; + lex.fiber = f; + destList->push_back(lex); + + p = p + v; + } +} + +void RasterizeFiberLex(rtsFiberNetwork* sourceNet, int f, list* destList) +{ + /*resamples a fiber f and stores the resulting points as + lexicographic values in LexicographicList. + */ + + int num_points = sourceNet->FiberList[f].pointList.size(); + int p; + + point3D pos = sourceNet->getNodeCoord(f, 0); + + point3D coord; + for(p=0; pgetFiberPoint(f, p); + RasterizeSegmentLex(f, pos, coord, destList); + pos = coord; + } + RasterizeSegmentLex(f, pos, sourceNet->getNodeCoord(f, 1), destList); +} + +list RasterizeNetworkLex(rtsFiberNetwork* sourceNet, rtsDTGrid3D* destGrid) +{ + /*This function creates a list of points on the focus network + and stores them in a list as a lexicographic point: + lp = z*sx*sy + y*sx + x + */ + + //create a list that stores Lexicographic identifiers + list LexicographicList; + + int num_fibers = sourceNet->FiberList.size(); + int f; + for(f=0; f p; + list::iterator i; + unsigned long lexID; + for(i=LexicographicList.begin(); i!=LexicographicList.end(); i++) + { + lexID = (*i).ID; + p.x = (lexID)/(GRID_DIM.y*GRID_DIM.z); + p.y = (lexID - p.x*GRID_DIM.y*GRID_DIM.z)/GRID_DIM.z; + p.z = (lexID - p.x*GRID_DIM.y*GRID_DIM.z - p.y*GRID_DIM.z); + //cout<<"Pushing "; + //p.print(); + //focusField->set(p.x, p.y, p.z, 255); + destGrid->push(p.x, p.y, p.z, 255); + //p.print(); + + } + return LexicographicList; +} + +void EvaluateDistanceField(rtsDTGrid3D* initialField) +{ + + initialField->background = 255; + (*initialField) = 0; + FastSweeping3D(initialField, DILATION, VOXEL_SIZE); + +} + +float Gaussian(float x, float std) +{ + float prefix = 1.0/sqrt(2.0*3.14159*std*std); + float suffix = exp(-(x*x)/(2.0*std*std)); + return prefix*suffix; +} + +void EvaluateGaussianEnvelope(rtsDTGrid3D* distanceField) +{ + + rtsDTGrid3D::iterator i; + float G; + float zero_value = Gaussian(0.0, STD); + for(i=distanceField->begin(); i!=distanceField->after(); i++) + { + G = Gaussian(i.Value(), STD); + //normalize the gaussian + G /= zero_value; + //cout<background = 0.0; + +} + +rtsDTGrid3D CalculateDifference(rtsDTGrid3D* grid0, rtsDTGrid3D* grid1) +{ + rtsDTGrid3D result; + result = (*grid0) - (*grid1); + + + rts_itkVolume testVolume; + testVolume.InsertDTGrid(&result); + testVolume.SaveVOL("rasterized.vol"); + cout<<"Size: "< std_1) + test_val = 0.0; + if(test_val > std_1) + gold_val = 0.0; + */ + if(test_val > std_1 && gold_val > std_1) + return 0.0; + return gold_val - test_val; + +} + +rtsDTGrid3D CompareGrids(rtsDTGrid3D* goldGrid, rtsDTGrid3D* testGrid) +{ + rtsDTGrid3D result; + //create an iterator for each DT Grid + rtsDTGrid3D::iterator gold_i = goldGrid->begin(); + rtsDTGrid3D::iterator test_i = testGrid->begin(); + + //compute the value at one standard deviation + float std_1 = Gaussian(STD, STD); + cout<<"Value at 1 STD: "<after() && test_i != testGrid->after()) + { + //if the iterators are at the same coordinate + if(gold_i.Coord() == test_i.Coord()) + { + //insert result of the comparison into the new grid + result_value = ComparePoints(gold_i.Value(), test_i.Value(), std_1); + result.push(gold_i.X1(), gold_i.X2(), gold_i.X3(), result_value); + //increment both + gold_i++; test_i++; + } + //add the lowest (lexicographically) value to the background, insert the result, and increment + else if( (gold_i.Coord() < test_i.Coord()) ) + { + result_value = ComparePoints(gold_i.Value(), testGrid->background, std_1); + result.push(gold_i.X1(), gold_i.X2(), gold_i.X3(), result_value); + gold_i++; + } + else if( (test_i.Coord() < gold_i.Coord()) ) + { + result_value = ComparePoints(goldGrid->background, test_i.Value(), std_1); + result.push(test_i.X1(), test_i.X2(), test_i.X3(), result_value); + test_i++; + } + } + + cout<<"here"<after()) + { + result_value = ComparePoints(gold_i.Value(), testGrid->background, std_1); + result.push(gold_i.X1(), gold_i.X2(), gold_i.X3(), result_value); + gold_i++; + } + + while(test_i != testGrid->after()) + { + result_value = ComparePoints(goldGrid->background, test_i.Value(), std_1); + result.push(test_i.X1(), test_i.X2(), test_i.X3(), result_value); + test_i++; + } + + + return result; + + + rts_itkVolume testVolume; + testVolume.InsertDTGrid(&result); + testVolume.SaveVOL("rasterized.vol"); + cout<<"Size: "<* inputGrid, float &total_positive, float &total_negative) +{ + total_positive = total_negative = 0.0; + rtsDTGrid3D::iterator i; + float value; + for(i = inputGrid->begin(); i != inputGrid->after(); i++) + { + value = i.Value(); + if(value > 0.0) + total_positive += value; + if(value < 0.0) + total_negative += -value; + } +} + +float SumGrid(rtsDTGrid3D* inputGrid) +{ + rtsDTGrid3D::iterator i; + float result = 0.0; + for(i = inputGrid->begin(); i != inputGrid->after(); i++) + { + result += i.Value(); + } + return result; +} + +void StoreInFiber(rtsFiberNetwork* network, int fiber, float value) +{ + if(network->FiberList[fiber].FiberData == NULL) + { + network->FiberList[fiber].FiberData = (void*)(new float[2]); + ((float*)network->FiberList[fiber].FiberData)[0] = 0.0; + ((float*)network->FiberList[fiber].FiberData)[1] = 0.0; + } + + ((float*)network->FiberList[fiber].FiberData)[0] += value; + ((float*)network->FiberList[fiber].FiberData)[1] += 1.0; + +} +void MapToExplicit(int radius, rtsDTGrid3D* inGrid, list* inLexList, rtsFiberNetwork* inNetwork) +{ + //create a template iterator to move over the input grid + rtsMultiDirectionStencil3D i; + + //create an iterator to move through the lexicographic list + list::iterator L = inLexList->begin(); + + //set the stencil positions + int x, y, z; + for(x=-radius; x<=radius; x++) + for(y=-radius; y<=radius; y++) + for(z=-radius; z<=radius; z++) + { + i.addPosition(x, y, z); + } + + i.Initialize(inGrid); + unsigned long lex_id; + float positives; + int f; + int p; + for(i.StartPPP(); L != inLexList->end(); i.PPP()) + { + //convert the iterator position to a Lexicographic ID + lex_id = (unsigned long)i.X1()*GRID_DIM.y*GRID_DIM.z + (unsigned long)i.X2()*GRID_DIM.z + (unsigned long)i.X3(); + if(lex_id == (*L).ID) + { + //get the fiber ID associated with the point + f = (*L).fiber; + + //get the positive value + positives = 0.0; + for(p = 0; p 0) + positives += i.getValue(p); + } + //if(positives > 0) + // cout< +#include + +//AnyOption library for parsing command-line options: +//http://www.hackorama.com/anyoption/ +#include "anyoption.h" + +//#include "ImplicitCalculations.h" +#include "rts/objJedi.h" +#include "DrawingFunctions.h" +#include "rts/rtsGraph.h" +#include "rts/rtsFilename.h" + +#define PRINT_LOG 1 + + +using namespace std; + +bool displayGui = false; +bool displayVerbose = false; + +rtsFiberNetwork* goldNetwork; +GLint goldNodeDL; + +rtsFiberNetwork* testNetwork; + +//standard deviation of the gaussian envelope surrounding the networks +float sigmaG, sigmaC; + +int glutwindow_id; + +void GlutMenuCallback(int option); + +void LoadNetworks(string gold_filename, string test_filename) +{ + + goldNetwork = new rtsFiberNetwork(); + testNetwork = new rtsFiberNetwork(); + + cout<LoadSWC(gold_filename); + cout<<"Nodes: "<NodeList.size()<FiberList.size()<getTotalLength()<LoadSWC(test_filename); + cout<<"Nodes: "<NodeList.size()<FiberList.size()<getTotalLength()< min_point0 = goldNetwork->min_pos; + point3D min_point1 = testNetwork->min_pos; + + point3D max_point0 = goldNetwork->max_pos; + point3D max_point1 = testNetwork->max_pos; + + point3D min_point(min(min_point0.x, min_point1.x), min(min_point0.y, min_point1.y), min(min_point0.z, min_point1.z)); + point3D max_point(max(max_point0.x, max_point1.x), max(max_point0.y, max_point1.y), max(max_point0.z, max_point1.z)); + point3D center_point = min_point + 0.5*(max_point - min_point); + network_span = (max_point - point3D(0, 0, 0)).Length(); + rts_glut_camera.setPosition(0, 0, 3*network_span); + rts_glut_camera.LookAt(center_point, vector3D(0.0, 1.0, 0.0)); + + + //load the shaders + makeColormap(); + //create strings for all of the shader filenames + string SmoothShaderVertexFilename = ""; + SmoothShaderVertexFilename += "SmoothShader_Vertex.glsl"; + string SmoothShaderFragmentFilename = ""; + SmoothShaderFragmentFilename += "SmoothShader_Fragment.glsl"; + string ErrorMapFragmentFilename = ""; + ErrorMapFragmentFilename += "ErrorMap_Fragment.glsl"; + + //Edge_ErrorShader.Init(); + //ErrorShader.Clean(); + + Edge_ErrorShader.AttachShader(GL_VERTEX_SHADER, SmoothShaderVertexFilename.c_str()); + + //Edge_ErrorShader.AttachShader(GL_FRAGMENT_SHADER, "ErrorShader_Fragment.glsl"); + Edge_ErrorShader.AttachShader(GL_FRAGMENT_SHADER, ErrorMapFragmentFilename.c_str()); + Edge_ErrorShader.Compile(); + + //cout<<"Edge Error Shader Log-------------------------"<Resample(sigmaG); + testNetwork->SubdivideNetwork(sigmaG); + + goldNetwork->Resample(sigmaG); + goldNetwork->SubdivideNetwork(sigmaG); + + + //compare the two networks and get the core graph + float gFPR, gFNR, cFPR, cFNR; + coreGraph = goldNetwork->CompareNetworks(testNetwork, sigmaG, sigmaC, gFPR, gFNR, cFPR, cFNR); + + //output error metrics + if(displayVerbose) + { + cout<<"GEOMETRY++++++++++++"< +#include +int main(int argc, char* argv[]) +{ + //Create an AnyOption object + AnyOption* opt = new AnyOption(); + opt->addUsage( "" ); + opt->addUsage( "Usage: NetMets [OPTION]... [GROUND TRUTH FILE] [TEST CASE FILE]" ); + opt->addUsage( "" ); + opt->addUsage( " -h --help Prints this help " ); + opt->addUsage( " -s --sigma Input sigma value (default = 25)"); + opt->addUsage( " -g --gui Display NetMets GUI " ); + opt->addUsage( " -v --verbose Verbose output " ); + opt->addUsage( "" ); + + opt->setFlag( "help", 'h' ); /* a flag (takes no argument), supporting long and short form */ + opt->setFlag( "gui", 'g' ); + opt->setOption("sigma", 's'); + opt->setFlag( "verbose", 'v'); + + /* go through the command line and get the options */ + opt->processCommandArgs( argc, argv ); + + //display the help + if( opt->getFlag( "help" ) || opt->getFlag( 'h' ) ) + { + opt->printUsage(); + return 0; + } + + //set a flag to display the GUI if requested + if( opt->getFlag( "gui" ) || opt->getFlag( 'g' ) ) + displayGui = true; + if( opt->getFlag("verbose") || opt->getFlag('v')) + displayVerbose = true; + + //set the sigma value based on user input + sigmaG = 25; + sigmaC = 25; + + if( opt->getValue( 's' ) != NULL || opt->getValue( "sigma" ) != NULL ) + { + sigmaG = atof(opt->getValue("sigma")); + sigmaC = sigmaG; + + } + //cout << "size = " << atof(opt->getValue("sigma")) << endl ; + + //get the filename arguments + int nArgs = opt->getArgc(); + char** sArgs = (char**)malloc(nArgs); + for(int a=0; agetArgv(a); + + + string gold_filename; + string test_filename; + + int resolution = 512; + float epsilon = 0.1; + if(nArgs < 2) + { + gold_filename = "00_GT.obj"; test_filename = "00_T.obj"; sigmaG = sigmaC = 25.0; + //sigmaG = 25; sigmaC = 25.0; + + } + else + { + gold_filename = sArgs[0]; + test_filename = sArgs[1]; + } + + + //load the network files + testNetwork = new rtsFiberNetwork(); + testNetwork->LoadFile(test_filename); + + goldNetwork = new rtsFiberNetwork(); + goldNetwork->LoadFile(gold_filename); + + + + + + + //compute the metric + ComputeNetMets(); + + + + //if the user does not want the GUI displayed, just return + if(!displayGui) return 0; + + //initialize OpenGL + InitializeOpenGL(); + + //create display lists + CreateDisplayLists(); + + + //set up GLUT and start + rts_glutStart(MyDisplayFunction); + + return 0; + +} diff --git a/Node_SmoothError_Vertex.glsl b/Node_SmoothError_Vertex.glsl new file mode 100644 index 0000000..b93d8b1 --- /dev/null +++ b/Node_SmoothError_Vertex.glsl @@ -0,0 +1,11 @@ +varying vec3 fragment_normal; + +void main(void) +{ + //fragment_normal = normalize(gl_NormalMatrix * gl_Normal); + //fragment_normal = vec3(0.0, 0.0, 0.5); + fragment_normal = gl_Normal; + gl_Position = gl_ModelViewProjectionMatrix * gl_Vertex; + //gl_TexCoord[0] = gl_Color; + gl_FrontColor = gl_Color; +} diff --git a/PerformanceData.h b/PerformanceData.h new file mode 100644 index 0000000..3170ee5 --- /dev/null +++ b/PerformanceData.h @@ -0,0 +1,155 @@ +// add the following to a cpp file: +// PerformanceData PD; + +#pragma once + +enum PerformanceDataType +{ + PD_DISPLAY=0, + PD_SPS, + PD_UNUSED0, + + //my stuff + SUBDIVIDE_NETWORKS, + BRUTE_FORCE_DISTANCE, + LOG_N_DIST_BUILD0, + LOG_N_DIST_SEARCH0, + LOG_N_DIST_BUILD1, + LOG_N_DIST_SEARCH1, + LOG_N_DISTANCE, + MY_TOPOLOGY, + BOOST_TOPOLOGY, + + + //end my stuff + PERFORMANCE_DATA_TYPE_COUNT +}; + +static char PDTypeNames[][255] = { + "Display ", + "Simulation Total ", + " ----------------- ", + //my stuff + "Network Subdivision", + "Brute Force Dist.--", + " Build 0 ", + " Search 0 ", + " Build 1 ", + " Search 1 ", + "ANN Dist.----------", + "MY Topology--------", + "BOOST Topology-----", + //end my stuff + +}; + +//------------------------------------------------------------------------------- + +//#ifdef DISPLAY +// #define PERFORMANCE_DATA_MEASUE_ENABLE +//#endif + +//#ifdef PERFORMANCE_DATA_MEASUE_ENABLE + +//------------------------------------------------------------------------------- + +#include +#include +#include + +#include +#include + +//------------------------------------------------------------------------------- + +class PerformanceData +{ +public: + PerformanceData() { ClearAll(); QueryPerformanceFrequency(&cps); } + ~PerformanceData(){} + + void ClearAll() + { + for ( int i=0; i maxTime[type] ) maxTime[type] = t; + totalTime[type] -= times[type][ pos[type] ]; + times[type][ pos[type] ] = t; + totalTime[type] += t; + pos[type]++; + if ( pos[type] == 0 ) dataReady[type] = true; + } + + void PrintResult( ostream &os,int i=PERFORMANCE_DATA_TYPE_COUNT) + { + os.setf(ios::fixed); + if ((i=0)){ + double a = GetAvrgTime(i); + if ( a ) + os<< PDTypeNames[i]<<" : avrg="< MAX_LONG_PREFIX_LENGTH ){ + *( _prefix + MAX_LONG_PREFIX_LENGTH ) = '\0'; + } + + strcpy (long_opt_prefix, _prefix); +} + +void +AnyOption::setFileCommentChar( char _comment ) +{ + file_delimiter_char = _comment; +} + + +void +AnyOption::setFileDelimiterChar( char _delimiter ) +{ + file_comment_char = _delimiter ; +} + +bool +AnyOption::CommandSet() +{ + return( command_set ); +} + +bool +AnyOption::FileSet() +{ + return( file_set ); +} + +void +AnyOption::noPOSIX() +{ + posix_style = false; +} + +bool +AnyOption::POSIX() +{ + return posix_style; +} + + +void +AnyOption::setVerbose() +{ + verbose = true ; +} + +void +AnyOption::printVerbose() +{ + if( verbose ) + cout << endl ; +} +void +AnyOption::printVerbose( const char *msg ) +{ + if( verbose ) + cout << msg ; +} + +void +AnyOption::printVerbose( char *msg ) +{ + if( verbose ) + cout << msg ; +} + +void +AnyOption::printVerbose( char ch ) +{ + if( verbose ) + cout << ch ; +} + +bool +AnyOption::hasOptions() +{ + return hasoptions; +} + +void +AnyOption::autoUsagePrint(bool _autousage) +{ + autousage = _autousage; +} + +void +AnyOption::useCommandArgs( int _argc, char **_argv ) +{ + argc = _argc; + argv = _argv; + command_set = true; + appname = argv[0]; + if(argc > 1) hasoptions = true; +} + +void +AnyOption::useFiileName( const char *_filename ) +{ + filename = _filename; + file_set = true; +} + +/* + * set methods for options + */ + +void +AnyOption::setCommandOption( const char *opt ) +{ + addOption( opt , COMMAND_OPT ); + g_value_counter++; +} + +void +AnyOption::setCommandOption( char opt ) +{ + addOption( opt , COMMAND_OPT ); + g_value_counter++; +} + +void +AnyOption::setCommandOption( const char *opt , char optchar ) +{ + addOption( opt , COMMAND_OPT ); + addOption( optchar , COMMAND_OPT ); + g_value_counter++; +} + +void +AnyOption::setCommandFlag( const char *opt ) +{ + addOption( opt , COMMAND_FLAG ); + g_value_counter++; +} + +void +AnyOption::setCommandFlag( char opt ) +{ + addOption( opt , COMMAND_FLAG ); + g_value_counter++; +} + +void +AnyOption::setCommandFlag( const char *opt , char optchar ) +{ + addOption( opt , COMMAND_FLAG ); + addOption( optchar , COMMAND_FLAG ); + g_value_counter++; +} + +void +AnyOption::setFileOption( const char *opt ) +{ + addOption( opt , FILE_OPT ); + g_value_counter++; +} + +void +AnyOption::setFileOption( char opt ) +{ + addOption( opt , FILE_OPT ); + g_value_counter++; +} + +void +AnyOption::setFileOption( const char *opt , char optchar ) +{ + addOption( opt , FILE_OPT ); + addOption( optchar, FILE_OPT ); + g_value_counter++; +} + +void +AnyOption::setFileFlag( const char *opt ) +{ + addOption( opt , FILE_FLAG ); + g_value_counter++; +} + +void +AnyOption::setFileFlag( char opt ) +{ + addOption( opt , FILE_FLAG ); + g_value_counter++; +} + +void +AnyOption::setFileFlag( const char *opt , char optchar ) +{ + addOption( opt , FILE_FLAG ); + addOption( optchar , FILE_FLAG ); + g_value_counter++; +} + +void +AnyOption::setOption( const char *opt ) +{ + addOption( opt , COMMON_OPT ); + g_value_counter++; +} + +void +AnyOption::setOption( char opt ) +{ + addOption( opt , COMMON_OPT ); + g_value_counter++; +} + +void +AnyOption::setOption( const char *opt , char optchar ) +{ + addOption( opt , COMMON_OPT ); + addOption( optchar , COMMON_OPT ); + g_value_counter++; +} + +void +AnyOption::setFlag( const char *opt ) +{ + addOption( opt , COMMON_FLAG ); + g_value_counter++; +} + +void +AnyOption::setFlag( const char opt ) +{ + addOption( opt , COMMON_FLAG ); + g_value_counter++; +} + +void +AnyOption::setFlag( const char *opt , char optchar ) +{ + addOption( opt , COMMON_FLAG ); + addOption( optchar , COMMON_FLAG ); + g_value_counter++; +} + +void +AnyOption::addOption( const char *opt, int type ) +{ + if( option_counter >= max_options ){ + if( doubleOptStorage() == false ){ + addOptionError( opt ); + return; + } + } + options[ option_counter ] = opt ; + optiontype[ option_counter ] = type ; + optionindex[ option_counter ] = g_value_counter; + option_counter++; +} + +void +AnyOption::addOption( char opt, int type ) +{ + if( !POSIX() ){ + printVerbose("Ignoring the option character \""); + printVerbose( opt ); + printVerbose( "\" ( POSIX options are turned off )" ); + printVerbose(); + return; + } + + + if( optchar_counter >= max_char_options ){ + if( doubleCharStorage() == false ){ + addOptionError( opt ); + return; + } + } + optionchars[ optchar_counter ] = opt ; + optchartype[ optchar_counter ] = type ; + optcharindex[ optchar_counter ] = g_value_counter; + optchar_counter++; +} + +void +AnyOption::addOptionError( const char *opt ) +{ + cout << endl ; + cout << "OPTIONS ERROR : Failed allocating extra memory " << endl ; + cout << "While adding the option : \""<< opt << "\"" << endl; + cout << "Exiting." << endl ; + cout << endl ; + exit(0); +} + +void +AnyOption::addOptionError( char opt ) +{ + cout << endl ; + cout << "OPTIONS ERROR : Failed allocating extra memory " << endl ; + cout << "While adding the option: \""<< opt << "\"" << endl; + cout << "Exiting." << endl ; + cout << endl ; + exit(0); +} + +void +AnyOption::processOptions() +{ + if( ! valueStoreOK() ) + return; +} + +void +AnyOption::processCommandArgs(int max_args) +{ + max_legal_args = max_args; + processCommandArgs(); +} + +void +AnyOption::processCommandArgs( int _argc, char **_argv, int max_args ) +{ + max_legal_args = max_args; + processCommandArgs( _argc, _argv ); +} + +void +AnyOption::processCommandArgs( int _argc, char **_argv ) +{ + useCommandArgs( _argc, _argv ); + processCommandArgs(); +} + +void +AnyOption::processCommandArgs() +{ + if( ! ( valueStoreOK() && CommandSet() ) ) + return; + + if( max_legal_args == 0 ) + max_legal_args = argc; + new_argv = (int*) malloc( (max_legal_args+1) * sizeof(int) ); + for( int i = 1 ; i < argc ; i++ ){/* ignore first argv */ + if( argv[i][0] == long_opt_prefix[0] && + argv[i][1] == long_opt_prefix[1] ) { /* long GNU option */ + int match_at = parseGNU( argv[i]+2 ); /* skip -- */ + if( match_at >= 0 && i < argc-1 ) /* found match */ + setValue( options[match_at] , argv[++i] ); + }else if( argv[i][0] == opt_prefix_char ) { /* POSIX char */ + if( POSIX() ){ + char ch = parsePOSIX( argv[i]+1 );/* skip - */ + if( ch != '0' && i < argc-1 ) /* matching char */ + setValue( ch , argv[++i] ); + } else { /* treat it as GNU option with a - */ + int match_at = parseGNU( argv[i]+1 ); /* skip - */ + if( match_at >= 0 && i < argc-1 ) /* found match */ + setValue( options[match_at] , argv[++i] ); + } + }else { /* not option but an argument keep index */ + if( new_argc < max_legal_args ){ + new_argv[ new_argc ] = i ; + new_argc++; + }else{ /* ignore extra arguments */ + printVerbose( "Ignoring extra argument: " ); + printVerbose( argv[i] ); + printVerbose( ); + printAutoUsage(); + } + printVerbose( "Unknown command argument option : " ); + printVerbose( argv[i] ); + printVerbose( ); + printAutoUsage(); + } + } +} + +char +AnyOption::parsePOSIX( char* arg ) +{ + + for( unsigned int i = 0 ; i < strlen(arg) ; i++ ){ + char ch = arg[i] ; + if( matchChar(ch) ) { /* keep matching flags till an option */ + /*if last char argv[++i] is the value */ + if( i == strlen(arg)-1 ){ + return ch; + }else{/* else the rest of arg is the value */ + i++; /* skip any '=' and ' ' */ + while( arg[i] == whitespace + || arg[i] == equalsign ) + i++; + setValue( ch , arg+i ); + return '0'; + } + } + } + printVerbose( "Unknown command argument option : " ); + printVerbose( arg ); + printVerbose( ); + printAutoUsage(); + return '0'; +} + +int +AnyOption::parseGNU( char *arg ) +{ + int split_at = 0; + /* if has a '=' sign get value */ + for( unsigned int i = 0 ; i < strlen(arg) ; i++ ){ + if(arg[i] == equalsign ){ + split_at = i ; /* store index */ + i = strlen(arg); /* get out of loop */ + } + } + if( split_at > 0 ){ /* it is an option value pair */ + char* tmp = (char*) malloc( (split_at+1)*sizeof(char) ); + for( int i = 0 ; i < split_at ; i++ ) + tmp[i] = arg[i]; + tmp[split_at] = '\0'; + + if ( matchOpt( tmp ) >= 0 ){ + setValue( options[matchOpt(tmp)] , arg+split_at+1 ); + free (tmp); + }else{ + printVerbose( "Unknown command argument option : " ); + printVerbose( arg ); + printVerbose( ); + printAutoUsage(); + free (tmp); + return -1; + } + }else{ /* regular options with no '=' sign */ + return matchOpt(arg); + } + return -1; +} + + +int +AnyOption::matchOpt( char *opt ) +{ + for( int i = 0 ; i < option_counter ; i++ ){ + if( strcmp( options[i], opt ) == 0 ){ + if( optiontype[i] == COMMON_OPT || + optiontype[i] == COMMAND_OPT ) + { /* found option return index */ + return i; + }else if( optiontype[i] == COMMON_FLAG || + optiontype[i] == COMMAND_FLAG ) + { /* found flag, set it */ + setFlagOn( opt ); + return -1; + } + } + } + printVerbose( "Unknown command argument option : " ); + printVerbose( opt ) ; + printVerbose( ); + printAutoUsage(); + return -1; +} +bool +AnyOption::matchChar( char c ) +{ + for( int i = 0 ; i < optchar_counter ; i++ ){ + if( optionchars[i] == c ) { /* found match */ + if(optchartype[i] == COMMON_OPT || + optchartype[i] == COMMAND_OPT ) + { /* an option store and stop scanning */ + return true; + }else if( optchartype[i] == COMMON_FLAG || + optchartype[i] == COMMAND_FLAG ) { /* a flag store and keep scanning */ + setFlagOn( c ); + return false; + } + } + } + printVerbose( "Unknown command argument option : " ); + printVerbose( c ) ; + printVerbose( ); + printAutoUsage(); + return false; +} + +bool +AnyOption::valueStoreOK( ) +{ + int size= 0; + if( !set ){ + if( g_value_counter > 0 ){ + size = g_value_counter * sizeof(char*); + values = (char**)malloc( size ); + for( int i = 0 ; i < g_value_counter ; i++) + values[i] = NULL; + set = true; + } + } + return set; +} + +/* + * public get methods + */ +char* +AnyOption::getValue( const char *option ) +{ + if( !valueStoreOK() ) + return NULL; + + for( int i = 0 ; i < option_counter ; i++ ){ + if( strcmp( options[i], option ) == 0 ) + return values[ optionindex[i] ]; + } + return NULL; +} + +bool +AnyOption::getFlag( const char *option ) +{ + if( !valueStoreOK() ) + return false; + for( int i = 0 ; i < option_counter ; i++ ){ + if( strcmp( options[i], option ) == 0 ) + return findFlag( values[ optionindex[i] ] ); + } + return false; +} + +char* +AnyOption::getValue( char option ) +{ + if( !valueStoreOK() ) + return NULL; + for( int i = 0 ; i < optchar_counter ; i++ ){ + if( optionchars[i] == option ) + return values[ optcharindex[i] ]; + } + return NULL; +} + +bool +AnyOption::getFlag( char option ) +{ + if( !valueStoreOK() ) + return false; + for( int i = 0 ; i < optchar_counter ; i++ ){ + if( optionchars[i] == option ) + return findFlag( values[ optcharindex[i] ] ) ; + } + return false; +} + +bool +AnyOption::findFlag( char* val ) +{ + if( val == NULL ) + return false; + + if( strcmp( TRUE_FLAG , val ) == 0 ) + return true; + + return false; +} + +/* + * private set methods + */ +bool +AnyOption::setValue( const char *option , char *value ) +{ + if( !valueStoreOK() ) + return false; + for( int i = 0 ; i < option_counter ; i++ ){ + if( strcmp( options[i], option ) == 0 ){ + values[ optionindex[i] ] = (char*) malloc((strlen(value)+1)*sizeof(char)); + strcpy( values[ optionindex[i] ], value ); + return true; + } + } + return false; +} + +bool +AnyOption::setFlagOn( const char *option ) +{ + if( !valueStoreOK() ) + return false; + for( int i = 0 ; i < option_counter ; i++ ){ + if( strcmp( options[i], option ) == 0 ){ + values[ optionindex[i] ] = (char*) malloc((strlen(TRUE_FLAG)+1)*sizeof(char)); + strcpy( values[ optionindex[i] ] , TRUE_FLAG ); + return true; + } + } + return false; +} + +bool +AnyOption::setValue( char option , char *value ) +{ + if( !valueStoreOK() ) + return false; + for( int i = 0 ; i < optchar_counter ; i++ ){ + if( optionchars[i] == option ){ + values[ optcharindex[i] ] = (char*) malloc((strlen(value)+1)*sizeof(char)); + strcpy( values[ optcharindex[i] ], value ); + return true; + } + } + return false; +} + +bool +AnyOption::setFlagOn( char option ) +{ + if( !valueStoreOK() ) + return false; + for( int i = 0 ; i < optchar_counter ; i++ ){ + if( optionchars[i] == option ){ + values[ optcharindex[i] ] = (char*) malloc((strlen(TRUE_FLAG)+1)*sizeof(char)); + strcpy( values[ optcharindex[i] ] , TRUE_FLAG ); + return true; + } + } + return false; +} + + +int +AnyOption::getArgc( ) +{ + return new_argc; +} + +char* +AnyOption::getArgv( int index ) +{ + if( index < new_argc ){ + return ( argv[ new_argv[ index ] ] ); + } + return NULL; +} + +/* dotfile sub routines */ + +bool +AnyOption::processFile() +{ + if( ! (valueStoreOK() && FileSet()) ) + return false; + return ( consumeFile(readFile()) ); +} + +bool +AnyOption::processFile( const char *filename ) +{ + useFiileName(filename ); + return ( processFile() ); +} + +char* +AnyOption::readFile() +{ + return ( readFile(filename) ); +} + +/* + * read the file contents to a character buffer + */ + +char* +AnyOption::readFile( const char* fname ) +{ + int length; + char *buffer; + ifstream is; + is.open ( fname , ifstream::in ); + if( ! is.good() ){ + is.close(); + return NULL; + } + is.seekg (0, ios::end); + length = is.tellg(); + is.seekg (0, ios::beg); + buffer = (char*) malloc(length*sizeof(char)); + is.read (buffer,length); + is.close(); + return buffer; +} + +/* + * scans a char* buffer for lines that does not + * start with the specified comment character. + */ +bool +AnyOption::consumeFile( char *buffer ) +{ + + if( buffer == NULL ) + return false; + + char *cursor = buffer;/* preserve the ptr */ + char *pline = NULL ; + int linelength = 0; + bool newline = true; + for( unsigned int i = 0 ; i < strlen( buffer ) ; i++ ){ + if( *cursor == endofline ) { /* end of line */ + if( pline != NULL ) /* valid line */ + processLine( pline, linelength ); + pline = NULL; + newline = true; + }else if( newline ){ /* start of line */ + newline = false; + if( (*cursor != comment ) ){ /* not a comment */ + pline = cursor ; + linelength = 0 ; + } + } + cursor++; /* keep moving */ + linelength++; + } + free (buffer); + return true; +} + + +/* + * find a valid type value pair separated by a delimiter + * character and pass it to valuePairs() + * any line which is not valid will be considered a value + * and will get passed on to justValue() + * + * assuming delimiter is ':' the behaviour will be, + * + * width:10 - valid pair valuePairs( width, 10 ); + * width : 10 - valid pair valuepairs( width, 10 ); + * + * :::: - not valid + * width - not valid + * :10 - not valid + * width: - not valid + * :: - not valid + * : - not valid + * + */ + +void +AnyOption::processLine( char *theline, int length ) +{ + bool found = false; + char *pline = (char*) malloc( (length+1)*sizeof(char) ); + for( int i = 0 ; i < length ; i ++ ) + pline[i]= *(theline++); + pline[length] = nullterminate; + char *cursor = pline ; /* preserve the ptr */ + if( *cursor == delimiter || *(cursor+length-1) == delimiter ){ + justValue( pline );/* line with start/end delimiter */ + }else{ + for( int i = 1 ; i < length-1 && !found ; i++){/* delimiter */ + if( *cursor == delimiter ){ + *(cursor-1) = nullterminate; /* two strings */ + found = true; + valuePairs( pline , cursor+1 ); + } + cursor++; + } + cursor++; + if( !found ) /* not a pair */ + justValue( pline ); + } + free (pline); +} + +/* + * removes trailing and preceeding whitespaces from a string + */ +char* +AnyOption::chomp( char *str ) +{ + while( *str == whitespace ) + str++; + char *end = str+strlen(str)-1; + while( *end == whitespace ) + end--; + *(end+1) = nullterminate; + return str; +} + +void +AnyOption::valuePairs( char *type, char *value ) +{ + if ( strlen(chomp(type)) == 1 ){ /* this is a char option */ + for( int i = 0 ; i < optchar_counter ; i++ ){ + if( optionchars[i] == type[0] ){ /* match */ + if( optchartype[i] == COMMON_OPT || + optchartype[i] == FILE_OPT ) + { + setValue( type[0] , chomp(value) ); + return; + } + } + } + } + /* if no char options matched */ + for( int i = 0 ; i < option_counter ; i++ ){ + if( strcmp( options[i], type ) == 0 ){ /* match */ + if( optiontype[i] == COMMON_OPT || + optiontype[i] == FILE_OPT ) + { + setValue( type , chomp(value) ); + return; + } + } + } + printVerbose( "Unknown option in resourcefile : " ); + printVerbose( type ); + printVerbose( ); +} + +void +AnyOption::justValue( char *type ) +{ + + if ( strlen(chomp(type)) == 1 ){ /* this is a char option */ + for( int i = 0 ; i < optchar_counter ; i++ ){ + if( optionchars[i] == type[0] ){ /* match */ + if( optchartype[i] == COMMON_FLAG || + optchartype[i] == FILE_FLAG ) + { + setFlagOn( type[0] ); + return; + } + } + } + } + /* if no char options matched */ + for( int i = 0 ; i < option_counter ; i++ ){ + if( strcmp( options[i], type ) == 0 ){ /* match */ + if( optiontype[i] == COMMON_FLAG || + optiontype[i] == FILE_FLAG ) + { + setFlagOn( type ); + return; + } + } + } + printVerbose( "Unknown option in resourcefile : " ); + printVerbose( type ); + printVerbose( ); +} + +/* + * usage and help + */ + + +void +AnyOption::printAutoUsage() +{ + if( autousage ) printUsage(); +} + +void +AnyOption::printUsage() +{ + + if( once ) { + once = false ; + cout << endl ; + for( int i = 0 ; i < usage_lines ; i++ ) + cout << usage[i] << endl ; + cout << endl ; + } +} + + +void +AnyOption::addUsage( const char *line ) +{ + if( usage_lines >= max_usage_lines ){ + if( doubleUsageStorage() == false ){ + addUsageError( line ); + exit(1); + } + } + usage[ usage_lines ] = line ; + usage_lines++; +} + +void +AnyOption::addUsageError( const char *line ) +{ + cout << endl ; + cout << "OPTIONS ERROR : Failed allocating extra memory " << endl ; + cout << "While adding the usage/help : \""<< line << "\"" << endl; + cout << "Exiting." << endl ; + cout << endl ; + exit(0); + +} diff --git a/anyoption.h b/anyoption.h new file mode 100644 index 0000000..6eeaf46 --- /dev/null +++ b/anyoption.h @@ -0,0 +1,270 @@ +#ifndef _ANYOPTION_H +#define _ANYOPTION_H + +#include +#include +#include +#include + +#define COMMON_OPT 1 +#define COMMAND_OPT 2 +#define FILE_OPT 3 +#define COMMON_FLAG 4 +#define COMMAND_FLAG 5 +#define FILE_FLAG 6 + +#define COMMAND_OPTION_TYPE 1 +#define COMMAND_FLAG_TYPE 2 +#define FILE_OPTION_TYPE 3 +#define FILE_FLAG_TYPE 4 +#define UNKNOWN_TYPE 5 + +#define DEFAULT_MAXOPTS 10 +#define MAX_LONG_PREFIX_LENGTH 2 + +#define DEFAULT_MAXUSAGE 3 +#define DEFAULT_MAXHELP 10 + +#define TRUE_FLAG "true" + +using namespace std; + +class AnyOption +{ + +public: /* the public interface */ + AnyOption(); + AnyOption(int maxoptions ); + AnyOption(int maxoptions , int maxcharoptions); + ~AnyOption(); + + /* + * following set methods specifies the + * special characters and delimiters + * if not set traditional defaults will be used + */ + + void setCommandPrefixChar( char _prefix ); /* '-' in "-w" */ + void setCommandLongPrefix( char *_prefix ); /* '--' in "--width" */ + void setFileCommentChar( char _comment ); /* '#' in shellscripts */ + void setFileDelimiterChar( char _delimiter );/* ':' in "width : 100" */ + + /* + * provide the input for the options + * like argv[] for commndline and the + * option file name to use; + */ + + void useCommandArgs( int _argc, char **_argv ); + void useFiileName( const char *_filename ); + + /* + * turn off the POSIX style options + * this means anything starting with a '-' or "--" + * will be considered a valid option + * which alo means you cannot add a bunch of + * POIX options chars together like "-lr" for "-l -r" + * + */ + + void noPOSIX(); + + /* + * prints warning verbose if you set anything wrong + */ + void setVerbose(); + + + /* + * there are two types of options + * + * Option - has an associated value ( -w 100 ) + * Flag - no value, just a boolean flag ( -nogui ) + * + * the options can be either a string ( GNU style ) + * or a character ( traditional POSIX style ) + * or both ( --width, -w ) + * + * the options can be common to the commandline and + * the optionfile, or can belong only to either of + * commandline and optionfile + * + * following set methods, handle all the aboove + * cases of options. + */ + + /* options comman to command line and option file */ + void setOption( const char *opt_string ); + void setOption( char opt_char ); + void setOption( const char *opt_string , char opt_char ); + void setFlag( const char *opt_string ); + void setFlag( char opt_char ); + void setFlag( const char *opt_string , char opt_char ); + + /* options read from commandline only */ + void setCommandOption( const char *opt_string ); + void setCommandOption( char opt_char ); + void setCommandOption( const char *opt_string , char opt_char ); + void setCommandFlag( const char *opt_string ); + void setCommandFlag( char opt_char ); + void setCommandFlag( const char *opt_string , char opt_char ); + + /* options read from an option file only */ + void setFileOption( const char *opt_string ); + void setFileOption( char opt_char ); + void setFileOption( const char *opt_string , char opt_char ); + void setFileFlag( const char *opt_string ); + void setFileFlag( char opt_char ); + void setFileFlag( const char *opt_string , char opt_char ); + + /* + * process the options, registerd using + * useCommandArgs() and useFileName(); + */ + void processOptions(); + void processCommandArgs(); + void processCommandArgs( int max_args ); + bool processFile(); + + /* + * process the specified options + */ + void processCommandArgs( int _argc, char **_argv ); + void processCommandArgs( int _argc, char **_argv, int max_args ); + bool processFile( const char *_filename ); + + /* + * get the value of the options + * will return NULL if no value is set + */ + char *getValue( const char *_option ); + bool getFlag( const char *_option ); + char *getValue( char _optchar ); + bool getFlag( char _optchar ); + + /* + * Print Usage + */ + void printUsage(); + void printAutoUsage(); + void addUsage( const char *line ); + void printHelp(); + /* print auto usage printing for unknown options or flag */ + void autoUsagePrint(bool flag); + + /* + * get the argument count and arguments sans the options + */ + int getArgc(); + char* getArgv( int index ); + bool hasOptions(); + +private: /* the hidden data structure */ + int argc; /* commandline arg count */ + char **argv; /* commndline args */ + const char* filename; /* the option file */ + char* appname; /* the application name from argv[0] */ + + int *new_argv; /* arguments sans options (index to argv) */ + int new_argc; /* argument count sans the options */ + int max_legal_args; /* ignore extra arguments */ + + + /* option strings storage + indexing */ + int max_options; /* maximum number of options */ + const char **options; /* storage */ + int *optiontype; /* type - common, command, file */ + int *optionindex; /* index into value storage */ + int option_counter; /* counter for added options */ + + /* option chars storage + indexing */ + int max_char_options; /* maximum number options */ + char *optionchars; /* storage */ + int *optchartype; /* type - common, command, file */ + int *optcharindex; /* index into value storage */ + int optchar_counter; /* counter for added options */ + + /* values */ + char **values; /* common value storage */ + int g_value_counter; /* globally updated value index LAME! */ + + /* help and usage */ + const char **usage; /* usage */ + int max_usage_lines; /* max usage lines reseverd */ + int usage_lines; /* number of usage lines */ + + bool command_set; /* if argc/argv were provided */ + bool file_set; /* if a filename was provided */ + bool mem_allocated; /* if memory allocated in init() */ + bool posix_style; /* enables to turn off POSIX style options */ + bool verbose; /* silent|verbose */ + bool print_usage; /* usage verbose */ + bool print_help; /* help verbose */ + + char opt_prefix_char; /* '-' in "-w" */ + char long_opt_prefix[MAX_LONG_PREFIX_LENGTH + 1]; /* '--' in "--width" */ + char file_delimiter_char; /* ':' in width : 100 */ + char file_comment_char; /* '#' in "#this is a comment" */ + char equalsign; + char comment; + char delimiter; + char endofline; + char whitespace; + char nullterminate; + + bool set; //was static member + bool once; //was static member + + bool hasoptions; + bool autousage; + +private: /* the hidden utils */ + void init(); + void init(int maxopt, int maxcharopt ); + bool alloc(); + void cleanup(); + bool valueStoreOK(); + + /* grow storage arrays as required */ + bool doubleOptStorage(); + bool doubleCharStorage(); + bool doubleUsageStorage(); + + bool setValue( const char *option , char *value ); + bool setFlagOn( const char *option ); + bool setValue( char optchar , char *value); + bool setFlagOn( char optchar ); + + void addOption( const char* option , int type ); + void addOption( char optchar , int type ); + void addOptionError( const char *opt); + void addOptionError( char opt); + bool findFlag( char* value ); + void addUsageError( const char *line ); + bool CommandSet(); + bool FileSet(); + bool POSIX(); + + char parsePOSIX( char* arg ); + int parseGNU( char *arg ); + bool matchChar( char c ); + int matchOpt( char *opt ); + + /* dot file methods */ + char *readFile(); + char *readFile( const char* fname ); + bool consumeFile( char *buffer ); + void processLine( char *theline, int length ); + char *chomp( char *str ); + void valuePairs( char *type, char *value ); + void justValue( char *value ); + + void printVerbose( const char *msg ); + void printVerbose( char *msg ); + void printVerbose( char ch ); + void printVerbose( ); + + +}; + +#endif /* ! _ANYOPTION_H */ diff --git a/loop_proxy_GT.obj b/loop_proxy_GT.obj new file mode 100644 index 0000000..f13a7e6 --- /dev/null +++ b/loop_proxy_GT.obj @@ -0,0 +1,36 @@ +v 187 414 0 +v 405 293 0 +v 397 413 0 +v 483 496 0 +v 689 591 0 +v 801 705 0 +v 496 677 0 +v 187 293 0 +v 290 226 0 +v 288 462 0 +v 405 293 0 +v 468 232 0 +v 572 228 0 +v 602 303 0 +v 594 432 0 +v 397 413 0 +v 698 466 0 +v 810 415 0 +v 906 484 0 +v 902 608 0 +v 689 591 0 +v 493 561 0 +v 403 510 0 +v 307 597 0 +v 305 708 0 +v 410 758 0 +l 1 8 9 2 +l 1 10 3 +l 2 11 3 +l 2 12 13 14 15 4 +l 3 16 4 +l 4 5 +l 5 17 18 19 20 6 +l 5 21 6 +l 6 7 +l 7 22 23 24 25 26 7 diff --git a/loop_proxy_T.obj b/loop_proxy_T.obj new file mode 100644 index 0000000..626bd2c --- /dev/null +++ b/loop_proxy_T.obj @@ -0,0 +1,41 @@ +v 186 412 0 +v 185 294 0 +v 123 233 0 +v 343 257 0 +v 361 268 0 +v 404 293 0 +v 401 337 0 +v 396 413 0 +v 308 335 0 +v 443 459 0 +v 591 430 0 +v 514 510 0 +v 801 704 0 +v 497 675 0 +v 438 731 0 +v 291 227 0 +v 287 461 0 +v 467 232 0 +v 571 228 0 +v 600 302 0 +v 690 589 0 +v 699 468 0 +v 905 484 0 +v 900 609 0 +v 497 675 0 +v 424 533 0 +v 306 597 0 +v 304 705 0 +l 1 2 +l 2 3 +l 2 16 4 +l 5 6 +l 6 7 +l 7 8 +l 7 9 +l 8 17 1 +l 8 10 +l 6 18 19 20 11 +l 12 21 22 23 24 13 +l 14 25 13 +l 14 26 27 28 15 diff --git a/molecule_GT.obj b/molecule_GT.obj new file mode 100644 index 0000000..a99bf7d --- /dev/null +++ b/molecule_GT.obj @@ -0,0 +1,83 @@ +v 79 265 0 +v 138 305 0 +v 249 374 0 +v 197 417 0 +v 197 490 0 +v 278 545 0 +v 126 445 0 +v 193 598 0 +v 324 377 0 +v 358 317 0 +v 361 436 0 +v 468 512 0 +v 554 569 0 +v 608 718 0 +v 633 511 0 +v 742 512 0 +v 794 551 0 +v 869 493 0 +v 728 600 0 +v 872 608 0 +v 466 238 0 +v 551 173 0 +v 611 32 0 +v 633 237 0 +v 741 237 0 +v 796 198 0 +v 729 150 0 +v 870 140 0 +v 872 260 0 +v 195 264 0 +v 250 302 0 +v 143 381 0 +v 416 396 0 +v 468 433 0 +v 414 551 0 +v 361 516 0 +v 545 673 0 +v 635 433 0 +v 688 395 0 +v 741 431 0 +v 685 548 0 +v 409 355 0 +v 465 320 0 +v 358 238 0 +v 410 199 0 +v 549 74 0 +v 684 203 0 +v 632 319 0 +v 685 356 0 +v 743 316 0 +l 1 2 +l 2 30 31 3 +l 3 4 +l 4 32 2 +l 4 5 +l 5 6 +l 5 7 +l 5 8 +l 3 9 +l 9 10 +l 9 11 +l 11 33 34 12 +l 12 35 36 11 +l 12 13 +l 13 37 14 +l 13 15 +l 15 38 39 40 16 +l 16 41 15 +l 16 17 +l 17 18 +l 17 19 +l 17 20 +l 10 42 43 21 +l 10 44 45 21 +l 21 22 +l 22 46 23 +l 22 24 +l 24 47 25 +l 24 48 49 50 25 +l 25 26 +l 26 27 +l 26 28 +l 26 29 diff --git a/molecule_T.obj b/molecule_T.obj new file mode 100644 index 0000000..cf9f5ce --- /dev/null +++ b/molecule_T.obj @@ -0,0 +1,90 @@ +v 137 306 0 +v 75 304 0 +v 249 373 0 +v 197 417 0 +v 279 545 0 +v 126 444 0 +v 193 596 0 +v 324 377 0 +v 335 355 0 +v 359 316 0 +v 347 333 0 +v 466 239 0 +v 550 172 0 +v 611 32 0 +v 632 237 0 +v 633 319 0 +v 741 237 0 +v 688 354 0 +v 870 140 0 +v 729 151 0 +v 872 258 0 +v 362 435 0 +v 469 511 0 +v 634 510 0 +v 609 715 0 +v 547 627 0 +v 633 431 0 +v 688 393 0 +v 742 511 0 +v 794 550 0 +v 871 493 0 +v 871 606 0 +v 728 599 0 +v 249 302 0 +v 142 379 0 +v 197 453 0 +v 207 479 0 +v 224 507 0 +v 159 465 0 +v 180 492 0 +v 195 519 0 +v 382 221 0 +v 466 318 0 +v 412 353 0 +v 570 83 0 +v 686 202 0 +v 739 317 0 +v 443 414 0 +v 415 550 0 +v 362 516 0 +v 509 537 0 +v 547 553 0 +v 580 548 0 +v 567 690 0 +v 545 666 0 +v 741 431 0 +v 687 547 0 +v 688 393 0 +l 1 2 +l 1 34 3 +l 3 4 +l 4 35 1 +l 4 36 37 38 5 +l 6 39 40 41 7 +l 3 8 +l 8 9 +l 10 11 +l 10 42 12 +l 12 43 44 10 +l 12 13 +l 13 45 14 +l 13 15 +l 15 16 +l 15 46 17 +l 17 47 18 +l 17 19 +l 20 21 +l 8 22 +l 22 48 23 +l 23 49 50 22 +l 23 51 52 53 24 +l 25 54 55 26 +l 24 27 +l 28 56 29 +l 29 57 24 +l 29 30 +l 30 31 +l 30 32 +l 30 33 +l 28 58 18 diff --git a/rts/CHECK_OPENGL_ERROR.h b/rts/CHECK_OPENGL_ERROR.h new file mode 100644 index 0000000..92c2f7d --- /dev/null +++ b/rts/CHECK_OPENGL_ERROR.h @@ -0,0 +1,15 @@ +#ifndef RTS_OPENGL_ERROR +#define RTS_OPENGL_ERROR + +#include +#include +#include + +#define CHECK_OPENGL_ERROR \ +{ GLenum error; \ + while ( (error = glGetError()) != GL_NO_ERROR) { \ + printf( "OpenGL ERROR: %s\nCHECK POINT: %s (line %d)\n", gluErrorString(error), __FILE__, __LINE__ ); \ + } \ +} + +#endif \ No newline at end of file diff --git a/rts/objJedi.h b/rts/objJedi.h new file mode 100644 index 0000000..ecf5845 --- /dev/null +++ b/rts/objJedi.h @@ -0,0 +1,1857 @@ +#ifndef OBJJEDI_H +#define OBJJEDI_H + +/*OBJ reader and writer. + +One aspect of the writer is based on OpenGL output commands. You can send OpenGL commands +(replacing "gl" and "GL" with "obj" and "OBJ" respectively to render to an Wavefront OBJ file +*/ + +//#include "rtsMath.h" +#include "rtsLinearAlgebra.h" +#include +#include +#include +#include +#include +using namespace std; + +//masks for valid vertex components +#define OBJ_V_X 0x1 +#define OBJ_V_Y 0x2 +#define OBJ_V_Z 0x4 +#define OBJ_V_W 0x8 + +#define OBJ_VT_U 0x1 +#define OBJ_VT_V 0x2 +#define OBJ_VT_W 0x4 + +#define OBJ_VN_I 0x1 +#define OBJ_VN_J 0x2 +#define OBJ_VN_K 0x4 + +#define OBJ_V 0x1 +#define OBJ_VT 0x2 +#define OBJ_VN 0x4 + +//primitive types +typedef unsigned char primitive_type; + +#define OBJ_END 0x0 +#define OBJ_POINTS 0x1 +#define OBJ_LINES 0x2 +#define OBJ_FACE 0x3 + + +#define OBJ_INVALID ULONG_MAX +//define the point structures +struct vertex_position{float x,y,z,w; unsigned char mask;}; +struct vertex_texture{float u,v,w; unsigned char mask;}; +struct vertex_normal{float i,j,k; unsigned char mask;}; +//the point structure contains indices to the relative 3-vectors in the lists +struct vertex +{ + unsigned int v; + unsigned int vt; + unsigned int vn; +}; + +struct primitive +{ + vector p; //indices to the point + unsigned char mask; //mask describing which features (v, vt, vn) are used + primitive_type type; //the type of primitive (points, lines, face) +}; + +//axis-aligned bounding box +struct AABB +{ + vertex_position min; + vertex_position max; +}; + +//create variable types +typedef unsigned int OBJint; + +//define the OBJ data class +class rtsOBJ +{ +private: + + /*Current state variables. These variables are committed to the OBJ object + when a vertex is added. However, we are careful to only add them once (if they don't + change with every vertex. + */ + vertex_texture current_vt; + vertex_normal current_vn; + bool vt_changed; //true if a new vt or vn was inserted since the last vertex + bool vn_changed; + unsigned char current_primitive_mask; //defines what coordinates are being used by the current primitive + //global variable storing the current render mode + OBJint g_CurrentMode; + //output file stream + ofstream g_objFile; + //obj file object + //objData g_OBJ; + //number of vertices since the last BEGIN + unsigned int g_NumVertices; + /*Attribute mask. This indicates what attributes are stored for each vertex. + Only a single mask applies to each vertex between objBegin() and objEnd(). The + attribute mask is flipped the first time an attribute is set but it is fixed after + the first vertex is passed.*/ + unsigned char g_AttributeMask; + /*Attribute reset mask. This indicates whether or not an attribute has been + reset since the last vertex was rendered. This applies to OBJ_VT and OBJ_VN*/ + unsigned char g_AttributeResetMask; + //latest texture coordinate sent + unsigned int g_LatestVT; + //latest vertex normal sent + unsigned int g_LatestVN; + + //extents of the points in the OBJ file + AABB m_bounds; + + OBJint f_InsertPointVertex(); + OBJint f_InsertLineVertex(); + OBJint f_InsertLineLoopVertex(); + OBJint f_InsertLineStripVertex(); + OBJint f_InsertTriangleVertex(); + OBJint f_InsertTriangleStripVertex(); + OBJint f_InsertTriangleFanVertex(); + OBJint f_InsertQuadVertex(); + OBJint f_InsertQuadStripVertex(); + OBJint f_InsertPolygonVertex(); + OBJint f_InsertPreviousFaceVertex(unsigned int v_index); + OBJint f_InsertFirstFaceVertex(unsigned int v_index); + OBJint f_InsertNewFaceVertex(); + OBJint f_InsertNewLineVertex(); + OBJint f_CreateNewFace(); + OBJint f_CreateNewLine(); + OBJint f_TerminateLineLoop(); + OBJint f_LoadOBJ(const char* filename); + OBJint f_LoadSWC(const char* filename); + + //insert coordinate commands + //these are used to handle coordinates of different dimensions (and are called by the associated public function) + OBJint f_InsertVertexf(float x, float y, float z, float w, unsigned char mask); + OBJint f_InsertTexCoordf(float u, float v, float w, unsigned char mask); + OBJint f_InsertNormalf(float i, float j, float k, unsigned char mask); + //output functions + void f_OutputVertices(); + void f_OutputPoints(); + void f_OutputLines(); + void f_OutputFaces(); + void f_OutputVertexNormals(); + void f_OutputTextureCoordinates(); + void f_ClearAll(); + + //methods for reading from a file + OBJint f_ReadPosition(ifstream &infile); + OBJint f_ReadNormal(ifstream &infile); + OBJint f_ReadTexCoord(ifstream &infile); + OBJint f_ReadVertex(ifstream &infile, vertex &new_point, unsigned char &mask); + OBJint f_ReadPrimitive(ifstream &infile, primitive_type type); + OBJint f_AdjustExtents(vertex_position v); + +public: + /*These vectors store the vertex and primitive data from the obj file. + All vertices, texture coordinates, and normals are stored in m_v, m_vt, m_vn + respectively. The vectors for each primitive store an index into m_v, m_vt, + and m_vn identifying the associated coordinate. Basically, the data is stored + in a structure very similar to the OBJ file itself. + */ + vector v_list; + vector vt_list; + vector vn_list; + vector primitives; + + vector points; + vector lines; + vector faces; + +public: + + OBJint objOpen(const char* filename); + OBJint objBegin(OBJint obj_mode); + OBJint objEnd(); + OBJint objClose(); + OBJint objNormal3f(float i, float j, float k); + OBJint objNormal2f(float i, float j); + OBJint objNormal1f(float i); + OBJint objTexCoord3f(float u, float v, float w); + OBJint objTexCoord2f(float u, float v); + OBJint objTexCoord1f(float u); + OBJint objVertex1f(float x); + OBJint objVertex2f(float x, float y); + OBJint objVertex3f(float x, float y, float z); + OBJint objVertex4f(float x, float y, float z, float w); + OBJint LoadFile(const char* filename); + OBJint SaveFile(const char* filename); + + //direct insertion methods + void insertVertexPosition(float x, float y, float z, unsigned char mask = OBJ_V_X | OBJ_V_Y | OBJ_V_Z); + void insertLine(unsigned int num_points, unsigned int* pointlist, unsigned int* normallist, unsigned int* texturelist); + + //get methods + unsigned int getNumVertices(); + unsigned int getNumLines(); + unsigned int getNumFaces(); + unsigned int getNumPointLists(); + unsigned int getNumTexCoords(); + unsigned int getNumNormals(); + + //these functions return the coordinate index as well as the value + unsigned int getNumFaceVertices(unsigned int face); + unsigned int getFaceVertex(unsigned int face, unsigned int vertex); + unsigned int getFaceNormal(unsigned int face, unsigned int normal); + unsigned int getFaceTexCoord(unsigned int face, unsigned int texcoord); + unsigned int getNumLineVertices(unsigned int line); + unsigned int getLineVertex(unsigned int line, unsigned int vertex); + unsigned int getLineTexCoord(unsigned int line, unsigned int texcoord); + point3D getVertex3d(unsigned int index); + point3D getTexCoord3d(unsigned int index); + point3D getNormal3d(unsigned int index); + vertex_position getVertex(unsigned int index); + vertex_texture getTexCoord(unsigned int index); + vertex_normal getNormal(unsigned int index); + AABB getBoundingBox(){return m_bounds;} + void Scale(float scale_x, float scale_y, float scale_z); + void Translate(float trans_x, float trans_y, float trans_z); + + float GetDistance(float x, float y, float z); + + //return data about primitives + unsigned int getPrimitiveType(unsigned int primitive); + + //constructors + rtsOBJ(); + void CopyOBJ(const rtsOBJ& obj); + //assignment + rtsOBJ& operator=(const rtsOBJ& obj) + { + CopyOBJ(obj); + return *this; + } + //copy + rtsOBJ(const rtsOBJ& obj) + { + CopyOBJ(obj); + //return *this; + } + + + + //Iterator stuff + class iterator + { + friend class rtsOBJ; + private: + rtsOBJ* obj; + bool end_primitive; + bool end_object; + unsigned int primitive_index; + public: + unsigned int operator*(); + void operator++(); + bool operator==(rtsOBJ::iterator operand); + bool operator!=(rtsOBJ::iterator operand); + unsigned int size(){return obj->primitives[primitive_index].p.size();}; + void print(); + }; + + iterator begin(); + iterator end(); +}; + +//define error codes +#define OBJ_OK 0x0100 +#define OBJ_ERROR 0x0101 + +//define the different modes +#define OBJ_NONE 0x0 +#define OBJ_POINTS 0x1 +#define OBJ_LINES 0x2 +#define OBJ_LINE_STRIP 0x3 +#define OBJ_LINE_LOOP 0x4 +#define OBJ_TRIANGLES 0x5 +#define OBJ_TRIANGLE_STRIP 0x6 +#define OBJ_TRIANGLE_FAN 0x7 +#define OBJ_QUADS 0x8 +#define OBJ_QUAD_STRIP 0x9 +#define OBJ_POLYGON 0x10 + + + + + + +//initialize the OBJ file +OBJint objOpen(char* filename); +//close the obj file +OBJint objClose(); + +//start rendering in a certain mode +OBJint objBegin(OBJint mode); +//stop the current rendering sequence +OBJint objEnd(void); +//render a vertex to the file +OBJint objVertex1f(float x); +OBJint objVertex2f(float x, float y); +OBJint objVertex3f(float x, float y, float z); +OBJint objVertex4f(float x, float y, float z, float w); +//set a normal vector for a vertex +OBJint objNormal3f(float i, float j, float k); +OBJint objNormal2f(float i, float j); +OBJint objNormal1f(float i); +//set a texture coordinate for a vertex +OBJint objTexCoord3f(float u, float v, float w); +OBJint objTexCoord2f(float u, float v); +OBJint objTexCoord1f(float u); + +//global variable for global functions +extern rtsOBJ g_OBJ; + +/*BUG NOTES +The standard function calls for inserting vertices don't work anymore. I've fixed points +but everything beyond that has to be updated. +*/ + +#include "objJedi.h" +#include "rtsVector3d.h" +#include "rtsPoint3d.h" +#include +//variable for use in global functions +rtsOBJ g_OBJ; + +/********UTILITY METHODS*******************************/ +void rtsOBJ::Scale(float scale_x, float scale_y, float scale_z) +{ + vector::iterator i; + for(i = v_list.begin(); i!= v_list.end(); i++) + { + (*i).x *= scale_x; + (*i).y *= scale_y; + (*i).z *= scale_z; + } +} + +void rtsOBJ::Translate(float trans_x, float trans_y, float trans_z) +{ + vector::iterator i; + for(i = v_list.begin(); i!= v_list.end(); i++) + { + (*i).x += trans_x; + (*i).y += trans_y; + (*i).z += trans_z; + } + +} + +float rtsOBJ::GetDistance(float x, float y, float z) +{ + //gets the distance between the specified point and the nearest surface of the OBJ + //currently only works for lines + + //cout<<"Primitives: "< p0, p1, p2; + float min_dist = 255; + float dist; + unsigned int p, l; + vector3D v, w; + float c1, c2, b; + point3D Pb; + + //for each line + for(l=0; l OBJ_POLYGON) + return OBJ_ERROR; + + //otherwise, go ahead and set the mode + g_CurrentMode = mode; + //set the number of vertices to zero + g_NumVertices = 0; + + //reset the current state primitive state + current_primitive_mask = 0x0; + + return OBJ_OK; +} + +OBJint rtsOBJ::objEnd() +{ + OBJint error = OBJ_OK; + //check to make sure the number of rendered vertices is valid for the current mode + switch(g_CurrentMode) + { + case OBJ_NONE: + //can't quit if we haven't started + error = OBJ_ERROR; + break; + case OBJ_LINES: + //if less than two vertices or an odd number of vertices + if(g_NumVertices < 2 || g_NumVertices%2 != 0) + { + //if there wasn't a vertex at all + if(g_NumVertices == 0) + error = OBJ_ERROR; + //if there was at least one vertex + else + { + //pop the last line off the list + primitives.pop_back(); + lines.pop_back(); + error = OBJ_ERROR; + } + } + break; + case OBJ_LINE_STRIP: + //if less than two vertices + if(g_NumVertices < 2) + { + //if there wasn't a vertex at all + if(g_NumVertices == 0) + error = OBJ_ERROR; + //if there was at least one vertex + else + { + //pop the last line off the list + primitives.pop_back(); + lines.pop_back(); + error = OBJ_ERROR; + } + } + break; + case OBJ_LINE_LOOP: + //if less than three vertices + if(g_NumVertices < 3) + { + //pop the last line off the list + primitives.pop_back(); + lines.pop_back(); + error = OBJ_ERROR; + } + //connect the first and last points + else + { + error = f_TerminateLineLoop(); + } + break; + case OBJ_TRIANGLES: + //if less than three vertices or not a power of three + if(g_NumVertices < 3 || g_NumVertices%3 !=0) + { + primitives.pop_back(); + faces.pop_back(); + error = OBJ_ERROR; + } + break; + case OBJ_TRIANGLE_STRIP: + //if less than three vertices + if(g_NumVertices < 3) + { + primitives.pop_back(); + faces.pop_back(); + error = OBJ_ERROR; + } + break; + case OBJ_TRIANGLE_FAN: + //if less than three vertices + if(g_NumVertices < 3) + { + primitives.pop_back(); + faces.pop_back(); + error = OBJ_ERROR; + } + break; + case OBJ_QUADS: + if(g_NumVertices < 4 || g_NumVertices%4 != 0) + { + primitives.pop_back(); + faces.pop_back(); + error = OBJ_ERROR; + } + break; + case OBJ_QUAD_STRIP: + //has to be at least 4 vertices and an even number + if(g_NumVertices < 4 || g_NumVertices%2 != 0) + { + primitives.pop_back(); + faces.pop_back(); + error = OBJ_ERROR; + } + break; + case OBJ_POLYGON: + //has to be at least three vertices + if(g_NumVertices < 3) + { + primitives.pop_back(); + faces.pop_back(); + error = OBJ_ERROR; + } + break; + } + + //reset the attribute mask + g_AttributeMask = 0x0; + //just for closure, reset the attribute reset mask + g_AttributeResetMask = 0x0; + //stop rendering + g_CurrentMode = OBJ_NONE; + return error; +} + +OBJint rtsOBJ::f_InsertVertexf(float x, float y, float z, float w, unsigned char mask) +{ + //make sure we're rendering + if(g_CurrentMode == OBJ_NONE) + return OBJ_ERROR; + + //insert the vertex into the vertex vector + vertex_position v; + v.x = x; v.y = y; v.z=z; v.w=w; v.mask = mask; + v_list.push_back(v); + f_AdjustExtents(v); //set the bounding box + //insert texture coordinate and normal if specified for this primitive + if((current_primitive_mask & OBJ_VT) && (vt_changed)) + vt_list.push_back(current_vt); + if((current_primitive_mask & OBJ_VN) && (vn_changed)) + vn_list.push_back(current_vn); + + + //increment the number of vertices inserted + g_NumVertices++; + + //handle each case of the vertex creation individually + OBJint error = OBJ_OK; + switch(g_CurrentMode) + { + case OBJ_POINTS: + error = f_InsertPointVertex(); + break; + case OBJ_LINES: + error = f_InsertLineVertex(); + break; + case OBJ_LINE_LOOP: + error = f_InsertLineLoopVertex(); + break; + case OBJ_LINE_STRIP: + error = f_InsertLineStripVertex(); + break; + case OBJ_TRIANGLES: + error = f_InsertTriangleVertex(); + break; + case OBJ_TRIANGLE_STRIP: + error = f_InsertTriangleStripVertex(); + break; + case OBJ_TRIANGLE_FAN: + error = f_InsertTriangleFanVertex(); + break; + case OBJ_QUADS: + error = f_InsertQuadVertex(); + break; + case OBJ_QUAD_STRIP: + error = f_InsertQuadStripVertex(); + break; + case OBJ_POLYGON: + error = f_InsertPolygonVertex(); + break; + } + //set the reset mask to zero + g_AttributeResetMask = 0x0; + + //set the attribute mask to include vertex position + g_AttributeMask = g_AttributeMask | OBJ_V; + + //return the result of the insertion + return error; +} + +OBJint rtsOBJ::objVertex1f(float x) +{ + return f_InsertVertexf(x, 0.0, 0.0, 0.0, OBJ_V_X); +} + +OBJint rtsOBJ::objVertex2f(float x, float y) +{ + return f_InsertVertexf(x, y, 0.0, 0.0, OBJ_V_X | OBJ_V_Y); +} + +OBJint rtsOBJ::objVertex3f(float x, float y, float z) +{ + return f_InsertVertexf(x, y, z, 0.0, OBJ_V_X | OBJ_V_Y | OBJ_V_Z); + +} + +OBJint rtsOBJ::objVertex4f(float x, float y, float z, float w) +{ + return f_InsertVertexf(x, y, z, w, OBJ_V_X | OBJ_V_Y | OBJ_V_Z | OBJ_V_W); +} + +OBJint rtsOBJ::objNormal3f(float i, float j, float k) +{ + return f_InsertNormalf(i, j, k, OBJ_VN_I | OBJ_VN_J | OBJ_VN_K); +} +OBJint rtsOBJ::objNormal2f(float i, float j) +{ + return f_InsertNormalf(i, j, 0.0, OBJ_VN_I | OBJ_VN_J); +} + +OBJint rtsOBJ::objNormal1f(float i) +{ + return f_InsertNormalf(i, 0.0, 0.0, OBJ_VN_I); +} + +OBJint rtsOBJ::f_InsertNormalf(float i, float j, float k, unsigned char mask) +{ + /*DEPRECATED + //if the mode is not rendering faces, there is an error + if(g_CurrentMode < OBJ_TRIANGLES) + return OBJ_ERROR; + //if the normal attribute flag is not set, set it (as long as a vertex hasn't been written) + if(!(g_AttributeMask & OBJ_VN)) + { + //if a vertex has been rendered, then we can't change the attribute, so exit + if(g_NumVertices > 0) + return OBJ_ERROR; + else + //otherwise, just set the attribute flag + g_AttributeMask = g_AttributeMask | OBJ_VN; + } + + + //insert the new normal into the normal list for the file + vertex_normal new_vn; + new_vn.i = i; new_vn.j = j; new_vn.k = k; + new_vn.mask = mask; + vn_list.push_back(new_vn); + */ + current_vn.i = i; //set the current texture state to the given coordinates + current_vn.j = j; + current_vn.k = k; + current_vn.mask = mask; //set the mask as appropriate + vn_changed = true; //the texture coordinate state has changed + current_primitive_mask = current_primitive_mask | OBJ_VN; //the current primitive now uses texture coordinates (if it didn't previously) + + return OBJ_OK; +} + +OBJint rtsOBJ::objTexCoord3f(float u, float v, float w) +{ + return f_InsertTexCoordf(u, v, w, OBJ_VT_U | OBJ_VT_V | OBJ_VT_W); +} + +OBJint rtsOBJ::objTexCoord2f(float u, float v) +{ + return f_InsertTexCoordf(u, v, 0.0, OBJ_VT_U | OBJ_VT_V); +} + +OBJint rtsOBJ::objTexCoord1f(float u) +{ + return f_InsertTexCoordf(u, 0.0, 0.0, OBJ_VT_U); +} + +OBJint rtsOBJ::f_InsertTexCoordf(float u, float v, float w, unsigned char mask) +{ + /* + DEPRECATED + //if the normal attribute flag is not set, set it (as long as a vertex hasn't been written) + if(!(g_AttributeMask & OBJ_VT)) + { + //if a vertex has been rendered, then we can't change the attribute, so exit + if(g_NumVertices > 0) + return OBJ_ERROR; + else + //otherwise, just set the attribute flag + g_AttributeMask = g_AttributeMask | OBJ_VT; + } + + //insert the new texture coordinate into the list for the file + vertex_texture new_vt; + new_vt.u = u; new_vt.v = v; new_vt.w = w; + new_vt.mask = mask; + vt_list.push_back(new_vt); + */ + current_vt.u = u; //set the current texture state to the given coordinates + current_vt.v = v; + current_vt.w = w; + current_vt.mask = mask; //set the mask as appropriate + vt_changed = true; //the texture coordinate state has changed + current_primitive_mask = current_primitive_mask | OBJ_VT; //the current primitive now uses texture coordinates (if it didn't previously) + + return OBJ_OK; +} + + +void rtsOBJ::insertVertexPosition(float x, float y, float z, unsigned char mask) +{ + vertex_position v; + v.x = x; + v.y = y; + v.z = z; + v.mask = mask; + + v_list.push_back(v); +} + +void rtsOBJ::insertLine(unsigned int num_points, unsigned int* pointlist, unsigned int* normallist, unsigned int* texturelist) +{ + //create the new primitive + primitive line; + line.type = OBJ_LINES; + line.mask = 0; + + //set the mask based on the passed data + if(pointlist != NULL) + line.mask = line.mask | OBJ_V; + if(normallist != NULL) + line.mask = line.mask | OBJ_VN; + if(texturelist != NULL) + line.mask = line.mask | OBJ_VT; + + //insert the line + unsigned int v; + vertex new_vert; + for(v=0; v 1) + points.push_back(p_index); + + return OBJ_OK; +} + +OBJint rtsOBJ::f_InsertLineVertex() +{ + //if this is an odd vertex, create a new line + if(g_NumVertices%2 == 1) + { + f_CreateNewLine(); + } + + f_InsertNewLineVertex(); + + return OBJ_OK; +} +OBJint rtsOBJ::f_InsertLineLoopVertex() +{ + //technically, this is the same as inserting a line strip vertex + f_InsertLineStripVertex(); + return OBJ_OK; +} + +OBJint rtsOBJ::f_InsertLineStripVertex() +{ + if(g_NumVertices == 1) + { + f_CreateNewLine(); + } + + f_InsertNewLineVertex(); + + return OBJ_OK; +} + +OBJint rtsOBJ::f_InsertTriangleVertex() +{ + //if this is the first vertex in a triangle, create a new triangle + if(g_NumVertices%3 == 1) + { + f_CreateNewFace(); + } + + + f_InsertNewFaceVertex(); + + return OBJ_OK; +} +OBJint rtsOBJ::f_InsertTriangleStripVertex() +{ + //in the special case of the first three vertices, just create a triangle + if(g_NumVertices < 4) + f_InsertTriangleVertex(); + else + { + + //create a new face for the new triangle + f_CreateNewFace(); + + //make sure to get the ordering right: + // Even vertices provide the previous two vertices "in order": (v2, v3, v4) + // Odd vertices provide the previous two vertices in reverse order: (v1, v3, v2) + if(g_NumVertices % 2 == 1) + { + f_InsertPreviousFaceVertex(2); + f_InsertPreviousFaceVertex(1); + f_InsertNewFaceVertex(); + } + else + { + f_InsertPreviousFaceVertex(1); + f_InsertNewFaceVertex(); + f_InsertPreviousFaceVertex(2); + } + } + + return OBJ_OK; +} +OBJint rtsOBJ::f_InsertTriangleFanVertex() +{ + //in the special case of the first three vertices, just create a triangle + if(g_NumVertices <4) + f_InsertTriangleVertex(); + else + { + //create a new face for the new triangle + f_CreateNewFace(); + //add the previous vertices to the face + f_InsertFirstFaceVertex(0); + f_InsertPreviousFaceVertex(2); + //insert the current vertex + f_InsertNewFaceVertex(); + } + + return OBJ_OK; +} +OBJint rtsOBJ::f_InsertQuadVertex() +{ + //if this is the first vertex in a quad, create a new quad + if(g_NumVertices%4 == 1) + { + f_CreateNewFace(); + } + + f_InsertNewFaceVertex(); + + return OBJ_OK; +} +OBJint rtsOBJ::f_InsertQuadStripVertex() +{ + //in the case of one of the first four vertices, just create a quad + if(g_NumVertices < 5) + f_InsertQuadVertex(); + //if the vertex is odd (it will be the third vertex of a quad) + else if(g_NumVertices%2 == 1) + { + //create a new face for the new quad + f_CreateNewFace(); + //add the previous two vertices + f_InsertPreviousFaceVertex(2); + f_InsertPreviousFaceVertex(3); + //add the current vertex + f_InsertNewFaceVertex(); + } + else + { + //if this is the last vertex of the quad, just add it + f_InsertNewFaceVertex(); + + } + return OBJ_OK; +} +OBJint rtsOBJ::f_InsertPolygonVertex() +{ + //if this is the first vertex, create the quad + if(g_NumVertices == 1) + { + f_CreateNewFace(); + } + f_InsertNewFaceVertex(); + + return OBJ_OK; +} + +OBJint rtsOBJ::f_InsertPreviousFaceVertex(unsigned int v_index) +{ + /*Finds the vertex used in the previous face with the given index and + inserts it into the current face. This limits the redundancy in the file + for re-used vertices (in strips and fans). This also transfers texture + and normal information.*/ + + //find the index of the previous face + unsigned int prev_f_index = primitives.size() - 2; + //find the index of the current face + unsigned int f_index = prev_f_index +1; + //add the vertex information from the previous face to this face + primitives[f_index].p.push_back(primitives[prev_f_index].p[v_index]); + + return OBJ_OK; +} + +OBJint rtsOBJ::f_InsertFirstFaceVertex(unsigned int v_index) +{ + /*Finds the vertex used in the first face (since the last objBegin()) + with the given index and inserts it at the end of the current face. + This includes texture and normal information.*/ + + //The result depends on the type of face being rendered + //So far, this function only applies to triangle fans + if(g_CurrentMode != OBJ_TRIANGLE_FAN) + return OBJ_ERROR; + + //calculate the number of faces that have been rendered + unsigned int num_faces = g_NumVertices - 2; + //find the index of the first face + unsigned int first_f_index = primitives.size() - num_faces; + //find the index of the current face + unsigned int f_index = primitives.size() - 1; + //transfer the vertex information from the first face to this one + primitives[f_index].p.push_back(primitives[first_f_index].p[v_index]); + + + return OBJ_OK; +} + +OBJint rtsOBJ::f_InsertNewFaceVertex() +{ + /*This inserts information about the current vertex into the current face*/ + //find the new vertex index + vertex p; + p.v = v_list.size() -1; + p.vt = vt_list.size() - 1; + p.vn = vn_list.size() -1; + //find the face index + unsigned int f_index = primitives.size() -1; + //INSERT VERTEX AND ATTRIBUTE DATA + //just add the vertex to the face + primitives[f_index].p.push_back(p); + + return OBJ_OK; +} + +OBJint rtsOBJ::f_InsertNewLineVertex() +{ + /*This inserts information about the current vertex into the current line*/ + //find the new vertex index + vertex p; + p.v = v_list.size() -1; + p.vt = vt_list.size() - 1; + //find the line index + unsigned int l_index = primitives.size() -1; + + //ADD VERTEX AND ATTRIBUTE INFORMATION + //add the vertex to the line + primitives[l_index].p.push_back(p); + + return OBJ_OK; +} + +OBJint rtsOBJ::f_CreateNewFace() +{ + primitive new_f; + new_f.type = OBJ_FACE; + new_f.mask = g_AttributeMask; + faces.push_back(primitives.size()); + primitives.push_back(new_f); + + + return OBJ_OK; +} + +OBJint rtsOBJ::f_CreateNewLine() +{ + primitive new_l; + new_l.type = OBJ_LINES; + new_l.mask = g_AttributeMask; + lines.push_back(primitives.size()); + primitives.push_back(new_l); + + return OBJ_OK; +} + + +OBJint rtsOBJ::f_TerminateLineLoop() +{ + /*This function just terminates the line loop by setting the last vertex + to the first vertex.*/ + if(g_CurrentMode != OBJ_LINE_LOOP) + return OBJ_ERROR; + //find the index for the current line + unsigned int l_index = lines.size() -1; + //insert the first vertex as the last vertex + primitives[l_index].p.push_back(primitives[l_index].p[0]); + + return OBJ_OK; +} + +void rtsOBJ::f_OutputVertices() +{ + //get the number of vertices in the object + unsigned int v_num = v_list.size(); + for(unsigned int i=0; i>new_vertex.x; + new_vertex.mask = OBJ_V_X; + + if(infile.peek() == ' ') + { + infile>>new_vertex.y; + new_vertex.mask = new_vertex.mask | OBJ_V_Y; + } + if(infile.peek() == ' ') + { + infile>>new_vertex.z; + new_vertex.mask = new_vertex.mask | OBJ_V_Z; + } + if(infile.peek() == ' ') + { + infile>>new_vertex.w; + new_vertex.mask = new_vertex.mask | OBJ_V_W; + } + int c = infile.peek(); + //ignore the rest of the line + infile.ignore(1000, '\n'); + c = infile.peek(); + + //cout<<"vertex read: "<>new_normal.i; + new_normal.mask = OBJ_VN_I; + //get every other component + if(infile.peek() == ' ') + { + infile>>new_normal.j; + new_normal.mask = new_normal.mask | OBJ_VN_J; + } + if(infile.peek() == ' ') + { + infile>>new_normal.k; + new_normal.mask = new_normal.mask | OBJ_VN_K; + } + //ignore the rest of the line + infile.ignore(1000, '\n'); + //insert the normal + vn_list.push_back(new_normal); + + return OBJ_OK; +} + +OBJint rtsOBJ::f_ReadTexCoord(ifstream &infile) +{ + vertex_texture new_texcoord; + infile>>new_texcoord.u; + new_texcoord.mask = OBJ_VT_U; + //get every other component + if(infile.peek() == ' ') + { + infile>>new_texcoord.v; + new_texcoord.mask = new_texcoord.mask | OBJ_VT_V; + } + if(infile.peek() == ' ') + { + infile>>new_texcoord.w; + new_texcoord.mask = new_texcoord.mask | OBJ_VT_W; + } + + //ignore the rest of the line + infile.ignore(1000, '\n'); + //insert the texture coordinate + vt_list.push_back(new_texcoord); + + return OBJ_OK; + +} + +OBJint rtsOBJ::f_ReadVertex(ifstream &infile, vertex &new_point, unsigned char &mask) +{ + //store the line vertex + infile>>new_point.v; + new_point.v--; + mask = OBJ_V; + //new_face.v.push_back(new_v - 1); + if(infile.peek() == '/') + { + infile.get(); + //if there actually is a texcoord + if(infile.peek() != '/') + { + infile>>new_point.vt; //get the index + new_point.vt--; + mask = mask | OBJ_VT; //update the mask + //new_face.vt.push_back(new_vt - 1); + } + } + //check for a normal + if(infile.peek() == '/') + { + infile.get(); + infile>>new_point.vn; //get the index + new_point.vn--; + mask = mask | OBJ_VN; //update the mask + } + + return OBJ_OK; +} + +OBJint rtsOBJ::f_ReadPrimitive(ifstream &infile, primitive_type type) +{ + //create a new point list + primitive new_primitive; + new_primitive.type = type; + //until the end of the line + while(infile.peek() != '\n') + { + //read each point + if(infile.peek() == ' ') + infile.get(); + else + { + vertex new_point; + f_ReadVertex(infile, new_point, new_primitive.mask); + new_primitive.p.push_back(new_point); + } + } + //ignore the rest of the line + infile.ignore(1000, '\n'); + + //push the id of the primitive into the new list + //:DEBUG: + if(type == OBJ_POINTS) + points.push_back(primitives.size()); + else if(type == OBJ_LINES) + lines.push_back(primitives.size()); + else if(type == OBJ_FACE) + faces.push_back(primitives.size()); + + primitives.push_back(new_primitive); //push the new primitive + + return OBJ_OK; +} + + +OBJint rtsOBJ::f_AdjustExtents(vertex_position v) +{ + if(v.x < m_bounds.min.x) + m_bounds.min.x = v.x; + if(v.y < m_bounds.min.y) + m_bounds.min.y = v.y; + if(v.z < m_bounds.min.z) + m_bounds.min.z = v.z; + + if(v.x > m_bounds.max.x) m_bounds.max.x = v.x; + if(v.y > m_bounds.max.y) m_bounds.max.y = v.y; + if(v.z > m_bounds.max.z) m_bounds.max.z = v.z; + + return OBJ_OK; +} + +OBJint rtsOBJ::f_LoadOBJ(const char* filename) +{ + f_ClearAll(); + //open the file as a stream + ifstream infile; + infile.open(filename); + if(!infile) + exit(1); + + unsigned int vertices = 0; + + string token; + infile>>token; + while(!infile.eof()) + { + //if the token is some vertex property + if(token == "v") + f_ReadPosition(infile); + else if(token == "vn") + f_ReadNormal(infile); + else if(token == "vt") + f_ReadTexCoord(infile); + else if(token == "p") + f_ReadPrimitive(infile, OBJ_POINTS); + else if(token == "l") + f_ReadPrimitive(infile, OBJ_LINES); + else if(token == "f") + f_ReadPrimitive(infile, OBJ_FACE); + else + infile.ignore(9999, '\n'); + //vertices++; + + infile>>token; + } + +} + +OBJint rtsOBJ::f_LoadSWC(const char* filename) +{ + f_ClearAll(); + //open the file as a stream + ifstream infile; + infile.open(filename); + if(!infile.is_open()) + return OBJ_ERROR; + + vector swcVertices; + float token; + objBegin(OBJ_LINES); + while(!infile.eof()) + { + vertex_position v; + infile>>token; //get the id + infile>>token; //get the fiber type + infile>>v.x; //get the node position + infile>>v.y; + infile>>v.z; + infile>>token; //get the radius + infile>>token; //get the parent + + //insert the node into the swc vector + swcVertices.push_back(v); + //now draw the line from the parent to the current node + if(token != -1) + { + unsigned int i = (unsigned int)token - 1; + objVertex3f(swcVertices[i].x, swcVertices[i].y, swcVertices[i].z); + objVertex3f(v.x, v.y, v.z); + } + } + objEnd(); + + + return OBJ_OK; + +} +OBJint rtsOBJ::LoadFile(const char* filename) +{ + string strFilename = filename; + int length = strFilename.length(); + string extension = strFilename.substr(strFilename.length() - 3, 3); + if(!extension.compare(string("obj"))) + return f_LoadOBJ(filename); + else if(!extension.compare(string("swc"))) + return f_LoadSWC(filename); + else return f_LoadOBJ(filename); + +} + +OBJint rtsOBJ::SaveFile(const char* filename) +{ + //open the file as a stream + ofstream outfile; + outfile.open(filename); + if(!outfile.is_open()) + return OBJ_ERROR; + + //output vertex positions + vector::iterator v; + for(v=v_list.begin(); v!= v_list.end(); v++) + { + outfile<<"v"; + if((*v).mask & OBJ_V_X) + { + outfile<<' '<<(*v).x; + if((*v).mask & OBJ_V_Y) + { + outfile<<' '<<(*v).y; + if((*v).mask & OBJ_V_Z) + { + outfile<<' '<<(*v).z; + if((*v).mask & OBJ_V_W) + outfile<<' '<<(*v).w; + } + } + } + outfile<<'\n'; + } + + //output vertex texture coordinates + vector::iterator vt; + for(vt=vt_list.begin(); vt!= vt_list.end(); vt++) + { + outfile<<"vt"; + if((*vt).mask & OBJ_VT_U) + { + outfile<<' '<<(*vt).u; + if((*vt).mask & OBJ_VT_V) + { + outfile<<' '<<(*vt).v; + if((*vt).mask & OBJ_VT_W) + outfile<<' '<<(*vt).w; + } + } + outfile<<'\n'; + } + + //output vertex normal coordinates + vector::iterator vn; + for(vn=vn_list.begin(); vn!= vn_list.end(); vn++) + { + outfile<<"vn"; + if((*vn).mask & OBJ_VN_I) + { + outfile<<' '<<(*vn).i; + if((*vn).mask & OBJ_VN_J) + { + outfile<<' '<<(*vn).j; + if((*vn).mask & OBJ_VN_K) + outfile<<' '<<(*vn).k; + } + } + outfile<<'\n'; + } + + //output each primitive + vector::iterator p; + vector::iterator vert; + for(p=primitives.begin(); p!= primitives.end(); p++) + { + switch((*p).type) + { + case OBJ_POINTS: + outfile<<"p"; //output the points token + break; + case OBJ_LINES: + outfile<<"l"; //output the lines token + break; + case OBJ_FACE: + outfile<<"f"; //output the face token + break; + } + + //for each vertex in the list for the primitive + for(vert = (*p).p.begin(); vert != (*p).p.end(); vert++) + { + outfile<<' '<<(*vert).v + 1; + if((*p).mask & OBJ_VT) + outfile<<'/'<<(*vert).vt + 1; + if((*p).mask & OBJ_VN) + { + if(!((*p).mask & OBJ_VT)) + outfile<<'/'; + outfile<<'/'<<(*vert).vn + 1; + } + } + outfile<<'\n'; + + } + + + +} + +//get methods +unsigned int rtsOBJ::getNumVertices(){return v_list.size();} +unsigned int rtsOBJ::getNumLines(){return lines.size();} +unsigned int rtsOBJ::getNumFaces(){return faces.size();} +unsigned int rtsOBJ::getNumPointLists(){return points.size();} +unsigned int rtsOBJ::getNumTexCoords(){return vt_list.size();} +unsigned int rtsOBJ::getNumNormals(){return vn_list.size();} + +//these functions return the coordinate index as well as the value +unsigned int rtsOBJ::getNumFaceVertices(unsigned int face){return primitives[face].p.size();} +unsigned int rtsOBJ::getFaceVertex(unsigned int face, unsigned int vertex){return primitives[face].p[vertex].v;} +unsigned int rtsOBJ::getFaceNormal(unsigned int face, unsigned int normal){return primitives[face].p[normal].vn;} +unsigned int rtsOBJ::getFaceTexCoord(unsigned int face, unsigned int texcoord){return primitives[face].p[texcoord].vt;} +unsigned int rtsOBJ::getNumLineVertices(unsigned int line){return primitives[line].p.size();} +unsigned int rtsOBJ::getLineVertex(unsigned int line, unsigned int vertex){return primitives[line].p[vertex].v;} +unsigned int rtsOBJ::getLineTexCoord(unsigned int line, unsigned int texcoord){return primitives[line].p[texcoord].vt;} +point3D rtsOBJ::getVertex3d(unsigned int index) +{ + return point3D(v_list[index].x, v_list[index].y, v_list[index].z); +} +point3D rtsOBJ::getTexCoord3d(unsigned int index) +{ + return point3D(vt_list[index].u, vt_list[index].v, vt_list[index].w); +} +point3D rtsOBJ::getNormal3d(unsigned int index) +{ + return point3D(vn_list[index].i, vn_list[index].j, vn_list[index].k); +} +vertex_position rtsOBJ::getVertex(unsigned int index){return v_list[index];} +vertex_texture rtsOBJ::getTexCoord(unsigned int index){return vt_list[index];} +vertex_normal rtsOBJ::getNormal(unsigned int index){return vn_list[index];} + + +unsigned int rtsOBJ::getPrimitiveType(unsigned int primitive) +{ + /* + switch(primitives[i].type) + { + case OBJ_POINTS: + return OBJ_POINTS; + break; + case OBJ_LINES: + return OBJ_LINES; + break; + case OBJ_FACE: + f_RenderFace(i); + break; + }*/ + return 0; +} +/****************************************************/ +/*******Iterator Methods*****************************/ +/****************************************************/ + +rtsOBJ::iterator rtsOBJ::begin() +{ + //create an iterator that will be returned and assign it to this OBJ + iterator result; + result.obj = this; + result.end_object = false; + result.end_primitive = false; + + //if there are no primitives, return the end iterator + if(primitives.size() == 0) + return end(); + + //start at the beginning of the primitive array + result.primitive_index = 0; + + return result; +} + +rtsOBJ::iterator rtsOBJ::end() +{ + //create an end iterator to return + iterator result; + result.obj = this; + result.end_primitive = true; + result.primitive_index = result.obj->primitives.size(); + + return result; +} + +void rtsOBJ::iterator::operator++() +{ + primitive_index++; + if(primitive_index >= obj->primitives.size()) + (*this) = obj->end(); + +} + +bool rtsOBJ::iterator::operator ==(rtsOBJ::iterator operand) +{ + if(operand.primitive_index == primitive_index) + return true; + else return false; +} + +bool rtsOBJ::iterator::operator !=(rtsOBJ::iterator operand) +{ + if(operand.primitive_index != primitive_index) + return true; + else return false; +} + +unsigned int rtsOBJ::iterator::operator*() +{ + return primitive_index; +} + +void rtsOBJ::iterator::print() +{ + cout<<"This is a test"< d; //direction that the camera is pointing + point3D p; //position of the camera + vector3D up; //"up" direction + float focus; //focal length of the camera + float fov; + + //private function makes sure that the up vector is orthogonal to the direction vector and both are normalized + void stabalize() + { + vector3D side = up.X(d); + up = d.X(side); + up.Normalize(); + d.Normalize(); + } + +public: + void setPosition(point3D pos) + { + p = pos; + } + void setPosition(float x, float y, float z){setPosition(point3D(x, y, z));} + + void setFocalDistance(float distance){focus = distance;} + void setFOV(float field_of_view){fov = field_of_view;} + + void LookAt(point3D pos) + { + //find the new direction + d = pos - p; + + //find the distance from the look-at point to the current position + focus = d.Length(); + + //stabalize the camera + stabalize(); + } + void LookAt(float px, float py, float pz){LookAt(point3D(px, py, pz));} + void LookAt(point3D pos, vector3D new_up){up = new_up; LookAt(pos);} + void LookAt(float px, float py, float pz, float ux, float uy, float uz){LookAt(point3D(px, py, pz), vector3D(ux, uy, uz));} + void LookAtDolly(float lx, float ly, float lz) + { + //find the current focus point + point3D f = p + focus*d; + vector3D T = point3D(lx, ly, lz) - f; + p = p + T; + } + + void Dolly(vector3D direction) + { + p = p+direction; + } + void Dolly(float x, float y, float z){Dolly(vector3D(x, y, z));} + void Push(float delta) + { + if(delta > focus) + delta = focus; + + focus -= delta; + + Dolly(d*delta); + } + + void Pan(float theta_x, float theta_y, float theta_z) + { + //x rotation is around the up axis + rtsQuaternion qx; + qx.CreateRotation(theta_x, up.x, up.y, up.z); + + //y rotation is around the side axis + vector3D side = up.X(d); + rtsQuaternion qy; + qy.CreateRotation(theta_y, side.x, side.y, side.z); + + //z rotation is around the direction vector + rtsQuaternion qz; + qz.CreateRotation(theta_z, d.x, d.y, d.z); + + //combine the rotations in x, y, z order + rtsQuaternion final = qz*qy*qx; + + //get the rotation matrix + matrix4x4 rot_matrix = final.toMatrix(); + + //apply the rotation + d = rot_matrix*d; + up = rot_matrix*up; + + //stabalize the camera + stabalize(); + + } + void Pan(float theta_x){Pan(theta_x, 0, 0);} + void Tilt(float theta_y){Pan(0, theta_y, 0);} + void Twist(float theta_z){Pan(0, 0, theta_z);} + + void Zoom(float delta) + { + fov -= delta; + if(fov < 0.5) + fov = 0.5; + if(fov > 180) + fov = 180; + } + + void OrbitFocus(float theta_x, float theta_y) + { + //find the focal point + point3D focal_point = p + focus*d; + + //center the coordinate system on the focal point + point3D centered = p - (focal_point - point3D(0, 0, 0)); + + //create the x rotation (around the up vector) + rtsQuaternion qx; + qx.CreateRotation(theta_x, up.x, up.y, up.z); + centered = qx.toMatrix()*centered; + + //get a side vector for theta_y rotation + vector3D side = up.X((point3D(0, 0, 0) - centered).Normalize()); + + rtsQuaternion qy; + qy.CreateRotation(theta_y, side.x, side.y, side.z); + centered = qy.toMatrix()*centered; + + //perform the rotation on the centered camera position + //centered = final.toMatrix()*centered; + + //re-position the camera + p = centered + (focal_point - point3D(0, 0, 0)); + + //make sure we are looking at the focal point + LookAt(focal_point); + + //stabalize the camera + stabalize(); + + } + + void Slide(float u, float v) + { + vector3D V = up.Normalize(); + vector3D U = up.X(d).Normalize(); + + p = p + (V * v) + (U * u); + } + + //accessor methods + point3D getPosition(){return p;} + vector3D getUp(){return up;} + vector3D getDirection(){return d;} + point3D getLookAt(){return p + focus*d;} + float getFOV(){return fov;} + + //output the camera settings + const void print(ostream& output) + { + output<<"Position: "<(0, 0, 0); + d = vector3D(0, 0, 1); + up = vector3D(0, 1, 0); + focus = 1; + + } +}; + + + +#endif \ No newline at end of file diff --git a/rts/rtsFilename.h b/rts/rtsFilename.h new file mode 100644 index 0000000..3b4074a --- /dev/null +++ b/rts/rtsFilename.h @@ -0,0 +1,67 @@ +#include + +using namespace std; + +#ifndef RTS_FILENAME_H +#define RTS_FILENAME_H + +class rtsFilename +{ +private: + string filename; + +public: + string getFilename() + { + int pos = filename.find_last_of("/\\"); + string name = filename.substr(pos+1, filename.size()); + return name; + } + + rtsFilename& operator=(const string str) + { + filename = str; + return *this; + } + rtsFilename(const string str){filename = str;} + rtsFilename(){filename = "";} + + string getExtension() + { + int pos = filename.find_last_of("."); + string ext = filename.substr(pos+1, filename.size() - pos); + return ext; + } + + string getDirectory() + { + int pos = filename.find_last_of("/\\"); + string directory; + if(pos != -1) + directory = filename.substr(0, pos); + else + directory = ""; + return directory; + } + + string getPrefix() + { + int slash = filename.find_last_of("/\\"); + int dot = filename.find_last_of("."); + + string prefix; + prefix = filename.substr(slash+1, dot - slash - 1); + return prefix; + + } + + string getString() + { + return filename; + } + + + +}; + +#endif \ No newline at end of file diff --git a/rts/rtsGraph.h b/rts/rtsGraph.h new file mode 100644 index 0000000..b75a713 --- /dev/null +++ b/rts/rtsGraph.h @@ -0,0 +1,142 @@ +#ifndef RTS_GRAPH_H +#define RTS_GRAPH_H + +#include +#include + +using namespace std; + +template +class rtsGraph +{ + unsigned int current_id; +public: + class rtsGraphNode; + typedef typename list::iterator rtsGraphNodePtr; + + class rtsGraphNode + { + public: + list ConnectionList; + unsigned int ID; + T Data; + }; + + list NodeList; + + rtsGraphNodePtr addNode(T data) + { + //Adds a node to the graph. This node is unconnected by default. + + //create the new node + rtsGraphNode new_node; + //set the data to the data provided + new_node.Data = data; + //set the id + new_node.ID = current_id++; + + //add the node to the node list + NodeList.push_back(new_node); + rtsGraphNodePtr result = NodeList.end(); + result--; + return result; + } + bool cmpPtr(const rtsGraphNodePtr a, const rtsGraphNodePtr b) + { + return (*a).ID < (*b).ID; + } + + void addConnection(rtsGraphNodePtr n0, rtsGraphNodePtr n1) + { + //adds a connection between two nodes + //connect n0 to n1 + (*n0).ConnectionList.push_back(n1); + (*n1).ConnectionList.push_back(n0); + + } + + void removeNode(rtsGraphNodePtr n) + { + + typename list::iterator i; + typename list::iterator j; + typename list::iterator temp; + rtsGraphNodePtr m; + rtsGraphNodePtr k; + /*Remove all of the connections TO n*/ + + //for each node m connected to n + for(i = (*n).ConnectionList.begin(); i != (*n).ConnectionList.end(); i++) + { + m = (*i); + //for each node k connected to m + j = (*m).ConnectionList.begin(); + while(j != (*m).ConnectionList.end()) + { + k = (*j); + //if k is the same as n, remove it + if(k == n) + { + temp = j; + j++; + (*m).ConnectionList.erase(temp); + } + else + j++; + + } + } + + /*Add all connections to neighboring nodes*/ + + //for each node m connected to n + for(i = (*n).ConnectionList.begin(); i != (*n).ConnectionList.end(); i++) + { + m = (*i); + //for each node k connected to n + for(j = (*n).ConnectionList.begin(); j != (*n).ConnectionList.end(); j++) + { + k = (*j); + if(k != m) + (*m).ConnectionList.push_back(k); + + } + //(*m).ConnectionList.sort(&rtsGraph::cmpPtr); + //sort((*m).ConnectionList.begin(), (*m).ConnectionList.end(), rtsGraph::cmpPtr); + (*m).ConnectionList.unique(); + } + + + + + + + + } + + + void PrintGraph() + { + rtsGraphNodePtr i; + typename list::iterator c; + for(i = NodeList.begin(); i != NodeList.end(); i++) + { + cout<<(*i).Data<<": "; + for(c = (*i).ConnectionList.begin(); c != (*i).ConnectionList.end(); c++) + { + if(c != (*i).ConnectionList.begin()) + cout<<"--"; + cout<<(*(*c)).Data; + } + cout< +#include "rtsVector3d.h" +#include "rtsPoint3d.h" +#include "rtsMatrix4x4.h" \ No newline at end of file diff --git a/rts/rtsMatrix4x4.h b/rts/rtsMatrix4x4.h new file mode 100644 index 0000000..f0a8dd8 --- /dev/null +++ b/rts/rtsMatrix4x4.h @@ -0,0 +1,196 @@ +#include +#include +#include +#include + +#include "rtsVector3d.h" +#include "rtsPoint3d.h" + +using namespace std; +#ifndef RTSMATRIX4X4_H +#define RTSMATRIX4X4_H + +#define RTS_PI 3.14159 + +///This class represents a 4x4 matrix, which is often used to represent affine transformations on 3D points and vectors. + +/// +///This class is designed to work with point3d and vector3d. Data is stored internally in a column-major form which is compatible with OpenGL. +/// + +template +class matrix4x4 +{ + /*stored as: + | 0 4 8 12 | + | | + | 1 5 9 13 | + | | + | 2 6 10 14 | + | | + | 3 7 11 15 | + */ +public: + T m_matrix[16]; + + //constructors + matrix4x4(); /// operator*(vector3D rhs); /// operator*(point3D rhs); /// operator*(matrix4x4 rhs); /// +matrix4x4::matrix4x4() +{ + memset(m_matrix, 0, sizeof(T)*16); +} + +template +matrix4x4::matrix4x4(T c0r0, T c0r1, T c0r2, T c0r3, + T c1r0, T c1r1, T c1r2, T c1r3, + T c2r0, T c2r1, T c2r2, T c2r3, + T c3r0, T c3r1, T c3r2, T c3r3) +{ + + T new_matrix[16] = {c0r0,c0r1,c0r2,c0r3,c1r0,c1r1,c1r2,c1r3,c2r0,c2r1,c2r2,c2r3,c3r0,c3r1,c3r2,c3r3}; + memcpy(m_matrix, new_matrix, 16*sizeof(T)); +} + +template +vector3D matrix4x4::operator*(vector3D rhs) +{ + vector3D result; + result.x = m_matrix[0] * rhs.x + m_matrix[4] * rhs.y + m_matrix[8] * rhs.z; + result.y = m_matrix[1] * rhs.x + m_matrix[5] * rhs.y + m_matrix[9] * rhs.z; + result.z = m_matrix[2] * rhs.x + m_matrix[6] * rhs.y + m_matrix[10] * rhs.z; + return result; +} + +template +point3D matrix4x4::operator *(point3D rhs) +{ + point3D result; + result.x = m_matrix[0] * rhs.x + m_matrix[4] * rhs.y + m_matrix[8] * rhs.z + m_matrix[12]; + result.y = m_matrix[1] * rhs.x + m_matrix[5] * rhs.y + m_matrix[9] * rhs.z + m_matrix[13]; + result.z = m_matrix[2] * rhs.x + m_matrix[6] * rhs.y + m_matrix[10] * rhs.z + m_matrix[14]; + T w = m_matrix[3] + m_matrix[7] + m_matrix[11] + m_matrix[15]; + result.x/=w; + result.y/=w; + result.z/=w; + return result; +} + +template +matrix4x4 matrix4x4::operator *(matrix4x4 rhs) +{ + matrix4x4 result; + int i,j,r; + for(i = 0; i<4; i++) + for(j = 0; j<4; j++) + for(r=0; r<4; r++) + result(i, j) += (*this)(i, r) * rhs(r, j); + + return result; +} + + + + +template +void matrix4x4::SetIdentity() +{ + T new_matrix[16] = {1,0,0,0,0,1,0,0,0,0,1,0,0,0,0,1}; + memcpy(m_matrix, new_matrix, 16*sizeof(T)); +} + +template +void matrix4x4::SetTranslation(T x, T y, T z) +{ + T new_matrix[16] = {1,0,0,0,0,1,0,0,0,0,1,0,x,y,z,1}; + memcpy(m_matrix, new_matrix, 16*sizeof(T)); +} + +template +void matrix4x4::SetScale(T x, T y, T z) +{ + T new_matrix[16] = {x,0,0,0,0,y,0,0,0,0,z,0,0,0,0,1}; + memcpy(m_matrix, new_matrix, 16*sizeof(T)); +} + +template +void matrix4x4::SetRotation(T angle, T x, T y, T z) +{ + //create the axis of rotation + vector3D axis(x, y, z); + T length = axis.Length(); //make sure that it is normalized + if(length != 1) + axis.Normalize(); + + T c = cos(angle*RTS_PI/180.0); //compute the sine and cosine of angle + T s = sin(angle*RTS_PI/180.0); + + //create the matrix + m_matrix[0] = x*x*(1-c) + c; + m_matrix[1] = y*x*(1-c) + z*s; + m_matrix[2] = x*z*(1-c) - y*s; + m_matrix[4] = x*y*(1-c) - z*s; + m_matrix[5] = y*y*(1-c) + c; + m_matrix[6] = y*z*(1-c) + x*s; + m_matrix[8] = x*z*(1-c) + y*s; + m_matrix[9] = y*z*(1-c) - x*s; + m_matrix[10]= z*z*(1-c) + c; + m_matrix[15]= 1; +} + + +template +void matrix4x4::Print(int width, int precision) +{ + cout< for linear algebra operations. + +template class point3D +{ + //friend class vector3D; +public: + //data values + T x; + T y; + T z; + + point3D(); /// param); /// param); /// operator+(vector3D param); /// operator-(vector3D param); /// operator-(point3D param); /// operator point3D(); /// &rhs) + { + os<<"("< & operator=(const point3D & rhs) + { + if(this != &rhs) + { + x = rhs.x; + y = rhs.y; + z = rhs.z; + } + return *this; + }*/ + + void print(); +}; + +template +template +point3D::operator point3D() +{ + point3D cast_result(x, y, z); + return cast_result; +} + +template point3D::point3D() +{ + x=0; + y=0; + z=0; +} + +template point3D::point3D(T new_x, T new_y, T new_z) +{ + x=new_x; + y=new_y; + z=new_z; +} + +template +point3D point3D::operator -(vector3D param) +{ + point3D result; //create the result + + result.x = x-param.x; //perform the computation + result.y = y-param.y; + result.z = z-param.z; + + return result; +} + +template +vector3D point3D::operator -(point3D param) +{ + vector3D result; + result.x = x - param.x; + result.y = y - param.y; + result.z = z - param.z; + + return result; +} + +template +point3D point3D::operator+(vector3D param) +{ + point3D result; //create the result + + result.x = x+param.x; //perform the computation + result.y = y+param.y; + result.z = z+param.z; + + return result; +} + +template +bool point3D::operator ==(point3D param) +{ + if(x == param.x && y == param.y && z == param.z) + return true; + else + return false; +} + +template +bool point3D::operator !=(point3D param) +{ + if(x != param.x || y != param.y || z != param.z) + return true; + else + return false; +} + + + +template +void point3D::print() +{ + cout<<"("< +ostream& point3D::operator<<(ostream& os, point3D &rhs) +{ + os<<"("< +class rtsQuaternion +{ +public: + T w; + T x; + T y; + T z; + + void normalize(); + void CreateRotation(T theta, T axis_x, T axis_y, T axis_z); + rtsQuaternion operator*(rtsQuaternion &rhs); + matrix4x4 toMatrix(); + + + rtsQuaternion(); + rtsQuaternion(T w, T x, T y, T z); + +}; + +template +void rtsQuaternion::normalize() +{ + double length=sqrt(w*w + x*x + y*y + z*z); + w=w/length; + x=x/length; + y=y/length; + z=z/length; +} + +template +void rtsQuaternion::CreateRotation(T theta, T axis_x, T axis_y, T axis_z) +{ + //assign the given Euler rotation to this quaternion + w = cos(theta/2.0); + x = axis_x*sin(theta/2.0); + y = axis_y*sin(theta/2.0); + z = axis_z*sin(theta/2.0); + + +} + +template +rtsQuaternion rtsQuaternion::operator *(rtsQuaternion ¶m) +{ + float A, B, C, D, E, F, G, H; + + + A = (w + x)*(param.w + param.x); + B = (z - y)*(param.y - param.z); + C = (w - x)*(param.y + param.z); + D = (y + z)*(param.w - param.x); + E = (x + z)*(param.x + param.y); + F = (x - z)*(param.x - param.y); + G = (w + y)*(param.w - param.z); + H = (w - y)*(param.w + param.z); + + rtsQuaternion result; + result.w = B + (-E - F + G + H) /2; + result.x = A - (E + F + G + H)/2; + result.y = C + (E - F + G - H)/2; + result.z = D + (E - F - G + H)/2; + + return result; +} + +template +matrix4x4 rtsQuaternion::toMatrix() +{ + matrix4x4 result; + + + double wx, wy, wz, xx, yy, yz, xy, xz, zz, x2, y2, z2; + + + // calculate coefficients + x2 = x + x; y2 = y + y; + z2 = z + z; + xx = x * x2; xy = x * y2; xz = x * z2; + yy = y * y2; yz = y * z2; zz = z * z2; + wx = w * x2; wy = w * y2; wz = w * z2; + + result(0, 0) = 1.0 - (yy + zz); + result(0, 1) = xy - wz; + //m[0][0] = 1.0 - (yy + zz); m[1][0] = xy - wz; + + result(0, 2) = xz + wy; + result(0, 3) = 0.0; + //m[2][0] = xz + wy; m[3][0] = 0.0; + + result(1, 0) = xy + wz; + result(1, 1) = 1.0 - (xx + zz); + //m[0][1] = xy + wz; m[1][1] = 1.0 - (xx + zz); + + result(1, 2) = yz - wx; + result(1, 3) = 0.0; + //m[2][1] = yz - wx; m[3][1] = 0.0; + + result(2, 0) = xz - wy; + result(2, 1) = yz + wx; + //m[0][2] = xz - wy; m[1][2] = yz + wx; + + result(2, 2) = 1.0 - (xx + yy); + result(3, 2) = 0.0; + //m[2][2] = 1.0 - (xx + yy); m[3][2] = 0.0; + + result(3, 0) = 0.0; + result(3, 1) = 0.0; + result(3, 2) = 0.0; + result(3, 3) = 1.0; + //m[0][3] = 0; m[1][3] = 0; + //m[2][3] = 0; m[3][3] = 1; + /* + double* orientationmatrix=(double*)m; + char c; + + + double* result=new double[16]; + double* array=(double*)m; + for(int i=0; i<16; i++) + result[i]=array[i]; + */ + + return result; +} + +template +rtsQuaternion::rtsQuaternion() +{ + w=0.0; x=0.0; y=0.0; z=0.0; +} + +template +rtsQuaternion::rtsQuaternion(T c, T i, T j, T k) +{ + w=c; x=i; y=j; z=k; +} + +#endif diff --git a/rts/rtsSourceCode.h b/rts/rtsSourceCode.h new file mode 100644 index 0000000..9e53e45 --- /dev/null +++ b/rts/rtsSourceCode.h @@ -0,0 +1,65 @@ +#ifndef RTSSOURCECODE_H +#define RTSSOURCECODE_H + +#include +#include +#include +#include + +using namespace std; + +///This class defines generic source code that can be loaded from text files. It is primarily used by the rts_glShaderProgram class for GLSL programming. + +class rtsSourceCode +{ +public: + vector source; //the actual source code + void clear() /// +#include +using namespace std; + +#ifndef RTSVECTOR3D_H +#define RTSVECTOR3D_H + + +template class vector3D +{ + //friend class point3D; +public: + //data values + T x; + T y; + T z; + + vector3D(); /// operator+(vector3D param); /// operator-(vector3D param); /// param); /// operator*(T param); /// X(vector3D param); /// operator*(const T lhs, vector3D rhs){return rhs*lhs;} /// Times(vector3D param); /// Normalize(); /// operator vector3D(); /// &rhs) + { + os<<"("< +template +vector3D::operator vector3D() +{ + vector3D cast_result(x, y, z); + return cast_result; +} + +template vector3D::vector3D() +{ + x=0; + y=0; + z=0; +} + +template vector3D::vector3D(T new_x, T new_y, T new_z) +{ + x=new_x; + y=new_y; + z=new_z; +} + +template +vector3D vector3D::operator -(vector3D param) +{ + vector3D result; //create the result + + result.x = x-param.x; //perform the computation + result.y = y-param.y; + result.z = z-param.z; + + return result; +} + +template +vector3D vector3D::operator+(vector3D param) +{ + vector3D result; //create the result + + result.x = x+param.x; //perform the computation + result.y = y+param.y; + result.z = z+param.z; + + return result; +} + +template +vector3D vector3D::Times(vector3D param) +{ + vector3D result; //create the result + + result.x = x*param.x; //perform the computation + result.y = y*param.y; + result.z = z*param.z; + + return result; +} + +template +vector3D vector3D::operator*(T param) +{ + vector3D result; //create the result + + result.x = x*param; //perform the computation + result.y = y*param; + result.z = z*param; + + return result; +} + +template +T vector3D::operator *(vector3D param) +{ + //This function computes the dot product between two vectors + return x*param.x + y*param.y + z*param.z; +} + +template +vector3D vector3D::X(vector3D param) +{ + vector3D result; + result.x=y*param.z - z*param.y; + result.y=z*param.x - x*param.z; + result.z=x*param.y - y*param.x; + + return result; +} + +template +vector3D vector3D::Normalize() +{ + T length = Length(); + if(length ==0.0) + { + x = y = z = 0; + } + else + { + x/= length; + y /=length; + z /=length; + } + return *this; +} + +template +T vector3D::Length() +{ + return sqrt(x*x + y*y + z*z); +} + +template +void vector3D::print() +{ + cout<<"["< +//#include "windows.h" +#include +#include "rtsSourceCode.h" + +class rts_glShaderObject +{ +private: + void init() + { + id = 0; + compiled = false; + type = GL_FRAGMENT_SHADER; + } +public: + bool compiled; + GLenum type; + rtsSourceCode source; + GLuint id; + string log; + + rts_glShaderObject(GLenum type, const char* filename) + { + init(); //initialize the shader + SetType(type); //set the shader type + LoadSource(filename); //load the source code + } + rts_glShaderObject(GLenum type, rtsSourceCode sourceCode) + { + init(); //initialize the shader + SetType(type); //set the shader type + source = sourceCode; + } + rts_glShaderObject() + { + init(); + } + rts_glShaderObject(GLenum type) + { + init(); + SetType(type); + } + void LoadSource(const char* filename) + { + source = rtsSourceCode(filename); //get the shader source code + + } + void SetType(GLenum type) + { + if(id != 0) //if a shader currently exists, delete it + { + glDeleteShader(id); + id = 0; + } + type = type; + id = glCreateShader(type); //create a shader object + if(id == 0) //if a shader was not created, log an error + { + log = "Error getting shader ID from OpenGL"; + return; + } + } + void UploadSource() + { + //create the structure for the shader source code + GLsizei count = source.source.size(); + GLchar** code_string = new GLchar*[count]; + GLint* length = new GLint[count]; + for(int l = 0; l + +using namespace std; + +class rts_glShaderProgram +{ +private: + void get_uniforms() + { + GLint num_uniforms; + glGetProgramiv(id, GL_ACTIVE_UNIFORMS, &num_uniforms); //get the number of uniform variables + GLint max_name_length; + glGetProgramiv(id, GL_ACTIVE_UNIFORM_MAX_LENGTH, &max_name_length); //get the maximum uniform name length + GLchar* name_buffer = new GLchar[max_name_length]; //create a buffer to store the name + GLsizei length; //I'm not using these yet + GLint size; + GLenum type; //variable's data type + GLint location; //GPU location of the variable + for(int i=0; i shader_list; //list of opengl shaders + vector uniform_list; //list of active uniform variables + vector texture_list; //list of texture maps + + rts_glShaderProgram() + { + linked = false; + id = 0; + } + void AttachShader(rts_glShaderObject shader) + { + if(id == 0) + { + Init(); + } + if(shader.id == 0) //if the shader is invalid + { + log = "Shader is invalid"; + return; + } + + //attach the shader to the program + glAttachShader(id, shader.id); //attach the shader to the program in OpenGL + CHECK_OPENGL_ERROR + shader_list.push_back(shader); //push the shader onto our list for later access + } + //type = GL_FRAGMENT_SHADER or GL_VERTEX_SHADER + void AttachShader(GLenum type, const char* filename) + { + rts_glShaderObject shader(type, filename); + AttachShader(shader); + } + void AttachShader(GLenum type, rtsSourceCode source) + { + rts_glShaderObject shader(type, source); + AttachShader(shader); + } + void PrintLog() + { + cout<::iterator iter; + for(iter = shader_list.begin(); iter != shader_list.end(); iter++) + { + (*iter).Compile(); + //(*iter).PrintLog(); + } + } + void Link() + { + glLinkProgram(id); //link the current shader program + GLint link_status; //test to see if the link went alright + glGetProgramiv(id, GL_LINK_STATUS, &link_status); + if(link_status != GL_TRUE) + { + linked = false; + } + else + linked = true; + + GLsizei length; + GLchar buffer[1000]; + glGetProgramInfoLog(id, 1000, &length, buffer); + log = buffer; + + get_uniforms(); //create the list of active uniform variables + } + void BeginProgram() + { + CHECK_OPENGL_ERROR + if(id == 0) //if the program is invalid, return + { + log = "Invalid program, cannot use."; + return; + } + if(!linked) + { + cout<<"Shader Program used without being linked."< 0) + glActiveTexture(GL_TEXTURE0); + CHECK_OPENGL_ERROR + + //return to OpenGL default shading + glUseProgram(0); + CHECK_OPENGL_ERROR + } + void PrintUniforms() + { + cout<<"Shader Uniforms: "< +#include + +using namespace std; + +enum rtsUniformEnum {RTS_FLOAT, RTS_INT, RTS_BOOL, RTS_FLOAT_MATRIX}; + +///This class stores a single uniform variable for GLSL and is designed to be used by the rts_glShaderProgram class. +struct rts_glShaderUniform +{ +public: + string name; //the name of the variable + GLint location; //the location in the program + void* p_value; //pointer to the global data representing the value in main memory + GLenum type; //variable type (float, int, vec2, etc.) + //rtsUniformEnum rts_type; //type of variable in rts format + //unsigned int num; //the number of values required by the variable (1 for float, 2 for vec2, etc.) + string log; + + //void convert_type(GLenum gl_type); //converts the OpenGL data type to something useful for rts + void submit_to_gpu() + { + if(location < 0) + return; + if(p_value == NULL) + { + cout<<"Error in uniform address: "< +#include "rtsVector3d.h" +#include "CHECK_OPENGL_ERROR.h" +#include + +///This class stores an OpenGL texture map and is used by rts_glShaderProgram. +class rts_glTextureMap +{ +private: + void get_type() //guesses the texture type based on the size + { + if(size.y == 0) + texture_type = GL_TEXTURE_1D; + else if(size.z == 0) + texture_type = GL_TEXTURE_2D; + else + texture_type = GL_TEXTURE_3D; + } + void set_wrapping() //set the texture wrapping based on the dimensions + { + CHECK_OPENGL_ERROR + switch(texture_type) + { + case GL_TEXTURE_3D: + glTexParameteri(texture_type, GL_TEXTURE_WRAP_R_EXT, GL_REPEAT); + case GL_TEXTURE_2D: + glTexParameteri(texture_type, GL_TEXTURE_WRAP_T, GL_MIRRORED_REPEAT); + case GL_TEXTURE_1D: + glTexParameteri(texture_type, GL_TEXTURE_WRAP_S, GL_REPEAT); + break; + case GL_TEXTURE_RECTANGLE_ARB: + glTexParameteri(texture_type, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE); + glTexParameteri(texture_type, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE); + break; + + default: + break; + } + CHECK_OPENGL_ERROR + } + //void set_bits(GLvoid* bits); +public: + vector3D size; //vector representing the size of the texture + GLuint name; //texture name assigned by OpenGL + GLenum texture_type; //1D, 2D, 3D + GLint internal_format; //number of components (ex. 4 for RGBA) + GLenum pixel_format; //type of data (RGBA, LUMINANCE) + GLenum data_type; //data type of the bits (float, int, etc.) + + //constructor + rts_glTextureMap() + { + name = 0; + } + + void BeginTexture() + { + glEnable(texture_type); + glBindTexture(texture_type, name); + CHECK_OPENGL_ERROR + } + void EndTexture() + { + glDisable(texture_type); + CHECK_OPENGL_ERROR + } + + ///Creates an OpenGL texture map. This function requires basic information about the texture map as well as a pointer to the bit data describing the texture. + void Init(GLvoid *bits, + GLenum type = GL_TEXTURE_2D, + GLsizei width = 256, + GLsizei height = 256, + GLsizei depth = 0, + GLint internalformat = 1, + GLenum format = GL_LUMINANCE, + GLenum datatype = GL_UNSIGNED_BYTE, + GLint interpolation = GL_LINEAR) + { + CHECK_OPENGL_ERROR + if(name != 0) + glDeleteTextures(1, &name); + + + CHECK_OPENGL_ERROR + if(datatype == GL_FLOAT) + { + glPixelStorei(GL_PACK_ALIGNMENT, 4); + glPixelStorei(GL_UNPACK_ALIGNMENT, 4); //I honestly don't know what this does but it fixes problems + } + else if(datatype == GL_UNSIGNED_BYTE) + { + //glPixelStorei(GL_UNPACK_ALIGNMENT, 1); + //glPixelStorei(GL_PACK_ALIGNMENT, 1); + } + else if(datatype == GL_UNSIGNED_SHORT) + { + //glPixelStorei(GL_UNPACK_ALIGNMENT, 2); + //glPixelStorei(GL_PACK_ALIGNMENT, 2); + } + CHECK_OPENGL_ERROR + glGenTextures(1, &name); //get the texture name from OpenGL + //cout<<"OpenGL Name: "<(width, height, depth); //assign the texture size + //get_type(); //guess the type based on the size + texture_type = type; //set the type of texture + glEnable(texture_type); //enable the texture map + CHECK_OPENGL_ERROR + glBindTexture(texture_type, name); //bind the texture for editing + CHECK_OPENGL_ERROR + set_wrapping(); //set the texture wrapping parameters + CHECK_OPENGL_ERROR + glTexParameteri(texture_type, GL_TEXTURE_MAG_FILTER, interpolation); //set filtering + CHECK_OPENGL_ERROR + glTexParameteri(texture_type, GL_TEXTURE_MIN_FILTER, interpolation); + CHECK_OPENGL_ERROR + internal_format = internalformat; //set the number of components per pixel + pixel_format = format; //set the pixel format + data_type = datatype; //set the data type + SetBits(bits); //send the bits to the OpenGL driver + glTexEnvf(GL_TEXTURE_ENV, GL_TEXTURE_ENV_MODE, GL_REPLACE); //replace the specified vertex color + CHECK_OPENGL_ERROR + glDisable(texture_type); + + glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_S, GL_CLAMP_TO_EDGE); + glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_T, GL_CLAMP_TO_EDGE); + glTexParameteri(GL_TEXTURE_3D, GL_TEXTURE_WRAP_R, GL_CLAMP_TO_EDGE); + } + void Clean() + { + if(name != 0) + { + glDeleteTextures(1, &name); + CHECK_OPENGL_ERROR + name = 0; + } + } + void SetBits(GLvoid *bits) + { + glEnable(texture_type); //enable the texture map + CHECK_OPENGL_ERROR + glBindTexture(texture_type, name); + CHECK_OPENGL_ERROR + + switch(texture_type) + { + case GL_TEXTURE_3D: + glTexImage3D(texture_type, 0, internal_format, size.x, size.y, size.z, 0, pixel_format, data_type, bits); + CHECK_OPENGL_ERROR + break; + case GL_TEXTURE_2D: + case GL_TEXTURE_RECTANGLE_ARB: + glTexImage2D(texture_type, 0, internal_format, size.x, size.y, 0, pixel_format, data_type, bits); + CHECK_OPENGL_ERROR + break; + case GL_TEXTURE_1D: + glTexImage1D(texture_type, 0, internal_format, size.x, 0, pixel_format, data_type, bits); + CHECK_OPENGL_ERROR + break; + default: + //glTexImage2D(texture_type, 0, internal_format, size.x, size.y, 0, pixel_format, data_type, bits); + break; + } + CHECK_OPENGL_ERROR + } + void ResetBits(GLvoid *bits) + { + glEnable(texture_type); //enable the texture map + CHECK_OPENGL_ERROR + glBindTexture(texture_type, name); + CHECK_OPENGL_ERROR + + switch(texture_type) + { + case GL_TEXTURE_3D: + //glTexImage3D(texture_type, 0, internal_format, size.x, size.y, size.z, 0, pixel_format, data_type, bits); + break; + case GL_TEXTURE_2D: + case GL_TEXTURE_RECTANGLE_ARB: + glTexSubImage2D(texture_type, 0, 0, 0, size.x, size.y, pixel_format, data_type, bits); + CHECK_OPENGL_ERROR + break; + case GL_TEXTURE_1D: + //glTexImage1D(texture_type, 0, internal_format, size.x, 0, pixel_format, data_type, bits); + break; + default: + //glTexImage2D(texture_type, 0, internal_format, size.x, size.y, 0, pixel_format, data_type, bits); + break; + } + glDisable(texture_type); + CHECK_OPENGL_ERROR + } + void* GetBits(GLenum format, GLenum type) + { + //returns the texture data + + int components; + switch(format) + { + case GL_RED: + case GL_GREEN: + case GL_BLUE: + case GL_ALPHA: + case GL_LUMINANCE: + components = 1; + break; + case GL_LUMINANCE_ALPHA: + components = 2; + break; + case GL_RGB: + case GL_BGR: + components = 3; + break; + case GL_RGBA: + case GL_BGRA: + components = 4; + break; + } + + int type_size; + switch(type) + { + case GL_UNSIGNED_BYTE: + case GL_BYTE: + type_size = sizeof(char); + break; + case GL_UNSIGNED_SHORT: + case GL_SHORT: + type_size = sizeof(short); + break; + case GL_UNSIGNED_INT: + case GL_INT: + type_size = sizeof(int); + break; + case GL_FLOAT: + type_size = sizeof(float); + break; + } + + //allocate memory for the texture + void* result = malloc(components*type_size * size.x * size.y); + + BeginTexture(); + glGetTexImage(texture_type, 0, format, type, result); + + CHECK_OPENGL_ERROR + EndTexture(); + + + return result; + + } +}; + +#define RTS_UNKNOWN 0 + +#endif \ No newline at end of file diff --git a/rts/rts_glutRenderWindow.h b/rts/rts_glutRenderWindow.h new file mode 100644 index 0000000..091bc38 --- /dev/null +++ b/rts/rts_glutRenderWindow.h @@ -0,0 +1,155 @@ +#include +#include +#include "rtsQuaternion.h" +#include "rtsCamera.h" + +float d_angle = 0.05; + +rtsCamera rts_glut_camera; +unsigned int mouse_x; +unsigned int mouse_y; +int mouse_button; + + +//user display function pointer +void (*UserDisplay)(void); + + + +/*void KeyboardFunction(unsigned char key, int x, int y) +{ + if(key == 27) + exit(0); + +}*/ + +/*void SpecialKeys(int key, int x, int y) +{ + rtsQuaternion new_rotation; + + if(key == GLUT_KEY_UP) + rts_glut_camera.OrbitFocus(0, d_angle); + if(key == GLUT_KEY_DOWN) + rts_glut_camera.OrbitFocus(0, -d_angle); + if(key == GLUT_KEY_LEFT) + rts_glut_camera.OrbitFocus(d_angle, 0.0); + if(key == GLUT_KEY_RIGHT) + rts_glut_camera.OrbitFocus(-d_angle, 0.0); +}*/ + +void MouseDrag(int x, int y) +{ + int dx = x - mouse_x; + int dy = y - mouse_y; + + if(mouse_button == GLUT_LEFT_BUTTON) + rts_glut_camera.OrbitFocus(-dx*0.01, dy*0.01); + else if(mouse_button == GLUT_MIDDLE_BUTTON) + rts_glut_camera.Zoom(dy*rts_glut_camera.getFOV()/60.0); + + mouse_x = x; + mouse_y = y; + +} + +void MouseMove(int x, int y) +{ + mouse_x = x; + mouse_y = y; +} + +void MouseClick(int button, int state, int x, int y) +{ + if(state == GLUT_DOWN) + mouse_button = button; + + +} + +void IdleFunction() +{ + glutPostRedisplay(); +} + +void RenderCamera() +{ + //set viewport + glMatrixMode(GL_PROJECTION); + glLoadIdentity(); + glViewport(0, 0, glutGet(GLUT_WINDOW_WIDTH), glutGet(GLUT_WINDOW_HEIGHT)); + + //compute the aspect ratio + float aspect_ratio = (float)glutGet(GLUT_WINDOW_WIDTH)/(float)glutGet(GLUT_WINDOW_HEIGHT); + gluPerspective(rts_glut_camera.getFOV(), aspect_ratio, 0.01, 10.0); + + //render the camera + glMatrixMode(GL_MODELVIEW); + glLoadIdentity(); + + point3D camera_position = rts_glut_camera.getPosition(); + vector3D camera_up = rts_glut_camera.getUp(); + gluLookAt(camera_position.x, + camera_position.y, + camera_position.z, + 0.0, 0.0, 0.0, + camera_up.x, + camera_up.y, + camera_up.z); +} + +void rts_glutInitialize(const char* WindowName, int width = 1024, int height = 768) +{ + char *myargv [1]; + int myargc=1; + myargv [0]=strdup ("AppName"); + glutInit(&myargc, myargv); + + glutInitDisplayMode(GLUT_DOUBLE | GLUT_RGBA | GLUT_DEPTH); + glutInitWindowSize(width, height); + glutInitWindowPosition(200,0); + glutCreateWindow(WindowName); + + //initialize GLEW + GLenum err = glewInit(); + if(err != GLEW_OK) + fprintf(stderr, "Error: %s\n", glewGetErrorString(err)); + if (!glewIsSupported("GL_VERSION_2_0")) + { + printf("OpenGL 2.0 not supported\n"); + //exit(1); + } + glEnable(GL_DEPTH_TEST); + + glClearColor(1.0, 1.0, 1.0, 0.0); + + rts_glut_camera.setPosition(0.0, 0.0, -1.0); + rts_glut_camera.setFOV(60); + rts_glut_camera.LookAt(0.0, 0.0, 0.0, 0.0, 1.0, 0.0); +} + +void DisplayFunction() +{ + RenderCamera(); + UserDisplay(); + glutSwapBuffers(); +} + +void rts_glutStart(void (*display_func)(void)) +{ + //glutSpecialFunc(SpecialKeys); + //glutKeyboardFunc(KeyboardFunction); + glutDisplayFunc(DisplayFunction); + glutMotionFunc(MouseDrag); + glutMouseFunc(MouseClick); + glutPassiveMotionFunc(MouseMove); + UserDisplay = display_func; + glutIdleFunc(IdleFunction); + + + + //glutReshapeFunc(ReshapeFunction); + //glutMouseFunc(MouseFunction); + //glutMotionFunc(MotionFunction); + glutMainLoop(); + +} diff --git a/rtsFiberNetwork.h b/rtsFiberNetwork.h new file mode 100644 index 0000000..6f90cd0 --- /dev/null +++ b/rtsFiberNetwork.h @@ -0,0 +1,2667 @@ +#include +#include +#include +//#include "PerformanceData.h" +#include "rts/objJedi.h" +#include "rts/rtsPoint3d.h" +#include "rts/rtsFilename.h" +#include +//#include + +#include +#include +#include +using namespace boost; + + +//Performance +//PerformanceData PD; + +#ifndef RTS_FIBER_NETWORK_H +#define RTS_FIBER_NETWORK_H + +//definitions for topologyNode labels +#define RTS_TOPOLOGY_NODE_NOEXIST 0 +#define RTS_TOPOLOGY_NODE_INVALID 1 +#define RTS_TOPOLOGY_NODE_VALID 2 + +#define RTS_TOPOLOGY_EDGE_NOEXIST 0 +#define RTS_TOPOLOGY_EDGE_EXIST 1 + +//Properties for the topology graph +struct vertex_position_t +{ + typedef vertex_property_tag kind; +}; +typedef property, property > VertexProperties; +typedef list EdgeSequence; +typedef pair EdgeMapping; +typedef vector CoreGraphList; +typedef property > EdgeProperties; +typedef adjacency_list TopologyGraph; +typedef graph_traits::edge_descriptor TopologyEdge; +typedef graph_traits::vertex_descriptor TopologyVertex; + + +EdgeSequence global_EdgeSequence; +float global_Weight; +list global_NeighborWeights; +list global_EdgeDescriptorSequence; +vector global_Predecessors; +TopologyVertex global_Source; +pair BOOST_SmallestEdge(int v0, int v1, const TopologyGraph& G) +{ + pair edge_pair = edge(v0, v1, G); + //if the edge doesn't exist, return the false pair + if(!edge_pair.second) + return edge_pair; + + //cout<<"Smallest edge: "<::out_edge_iterator oi, oi_end; + TopologyEdge min_e = edge_pair.first; + float min_weight = get(edge_weight_t(), G, min_e); + + for(boost::tuples::tie(oi, oi_end) = out_edges(v0, G); oi!=oi_end; oi++) + { + if(target(*oi, G) == v1) + { + //cout<<"Edge weight for "<= 0 && v != global_Source) + { + TopologyVertex v0, v1; + v0 = v1 = v; + pair e; + global_Weight = 0.0; + while(v1 != global_Source) + { + v1 = global_Predecessors[v0]; + e = BOOST_SmallestEdge(v0, v1, G); + global_EdgeSequence.push_back(get(edge_color_t(), G, e.first).front()); + global_EdgeDescriptorSequence.push_back(e.first); + global_Weight += get(edge_weight_t(), G, e.first); + v0 = v1; + } + //cout<<"Weight test: "< e0, pair e1) +{ + if(e0.second < e1.second) + return true; + else return false; +} + + +struct Fiber +{ + vector > pointList; + vector errorList; + int n0; + int n1; + float length; + float error; + int mapped_to; //edge in another graph that maps to this one (-1 if a mapping doesn't exist) + void* FiberData; + + Fiber() + { + error = 1.0; + FiberData = NULL; + } +}; + +struct Node +{ + point3D p; + void* NodeData; + float error; + int color; + int incident; + + Node() + { + error = 1.0; + NodeData = NULL; + incident = 0; + } +}; + +struct geometryPoint +{ + point3D p; + unsigned int fiberIdx; + unsigned int pointIdx; + float dist; + int fiberID; +}; + +struct topologyEdge +{ + int label; + unsigned int n0; + unsigned int n1; + float error; +}; +struct topologyNode +{ + point3D p; + int label; + list connections; + unsigned int compatible; +}; + + + + +#define MAX_DIST 9999.0; + +class rtsFiberNetwork +{ +private: + bool fiber_started; + unsigned int num_points; + float cull_value; //used to cull fibers based on geometric error + + vector getNetPointSamples(float subdiv) + { + vector result; + + unsigned int f; + list > fiberPoints; + list >::iterator p; + for(f = 0; f* N0, vector* N1); + void KD_ComputeEnvelopeDistance(rtsFiberNetwork* network, vector* Samples, float sigma) + { + //build the point arrays + ANNpointArray dataPts0 = annAllocPts(Samples->size(), 3); + for(unsigned int i=0; isize(); i++) + { + dataPts0[i][0] = (*Samples)[i].p.x; + dataPts0[i][1] = (*Samples)[i].p.y; + dataPts0[i][2] = (*Samples)[i].p.z; + } + + //create ANN variables + ANNkd_tree* kdTree; + ANNpoint queryPt = annAllocPt(3); + ANNidxArray nearestIdx = new ANNidx[1]; + ANNdistArray nearestDist = new ANNdist[1]; + + //compare network 0 to network 1 + //PD.StartTimer(LOG_N_DIST_BUILD0); + kdTree = new ANNkd_tree(dataPts0, Samples->size(), 3); + //PD.EndTimer(LOG_N_DIST_BUILD0); + //PD.StartTimer(LOG_N_DIST_SEARCH0); + + //test each point in the network to the Samples list + unsigned int f, p; + int nodenum; + float gauss_dist; + for(f=0; fFiberList.size(); f++) + { + //clear the error list + network->FiberList[f].errorList.clear(); + + //compute the distance at the nodes + nodenum = network->FiberList[f].n0; + queryPt[0] = network->NodeList[nodenum].p.x; + queryPt[1] = network->NodeList[nodenum].p.y; + queryPt[2] = network->NodeList[nodenum].p.z; + kdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + gauss_dist = 1.0f - GaussianEnvelope(sqrtf((float)nearestDist[0]), sigma); + network->NodeList[nodenum].error = gauss_dist; + + nodenum = network->FiberList[f].n1; + queryPt[0] = network->NodeList[nodenum].p.x; + queryPt[1] = network->NodeList[nodenum].p.y; + queryPt[2] = network->NodeList[nodenum].p.z; + kdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + gauss_dist = 1.0f - GaussianEnvelope(sqrtf((float)nearestDist[0]), sigma); + network->NodeList[nodenum].error = gauss_dist; + + //compute the distance at each point along the fiber + for(p=0; pFiberList[f].pointList.size(); p++) + { + queryPt[0] = network->FiberList[f].pointList[p].x; + queryPt[1] = network->FiberList[f].pointList[p].y; + queryPt[2] = network->FiberList[f].pointList[p].z; + kdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + gauss_dist = 1.0f - GaussianEnvelope(sqrtf((float)nearestDist[0]), sigma); + network->FiberList[f].errorList.push_back(gauss_dist); + } + } + /*for(int i=0; isize(); i++) + { + queryPt[0] = (*N1)[i].p.x; + queryPt[1] = (*N1)[i].p.y; + queryPt[2] = (*N1)[i].p.z; + kdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + (*N1)[i].dist = sqrt(nearestDist[0]); + } */ + //delete kdTree; + //PD.EndTimer(LOG_N_DIST_SEARCH0); + + /* + //compare network 1 to network 0 + PD.StartTimer(LOG_N_DIST_BUILD1); + kdTree = new ANNkd_tree(dataPts1, N1->size(), 3); + PD.EndTimer(LOG_N_DIST_BUILD1); + PD.StartTimer(LOG_N_DIST_SEARCH1); + for(int i=0; isize(); i++) + { + queryPt[0] = (*N0)[i].p.x; + queryPt[1] = (*N0)[i].p.y; + queryPt[2] = (*N0)[i].p.z; + kdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + (*N0)[i].dist = sqrt(nearestDist[0]); + } + PD.EndTimer(LOG_N_DIST_SEARCH1); + //delete kdTree; + */ + annClose(); + //delete kdTree; + } + + void BD_ComputeL1Distance(vector* N0, vector* N1); + list > SubdivideSegment(point3D p0, point3D p1, float spacing) + { + + //find the direction of travel + vector3D v = p1 - p0; + + //set the step size to the voxel size + float length = v.Length(); + v.Normalize(); + + float l; + list > result; + point3D p; + for(l=0.0; l > SubdivideFiber(unsigned int f, float spacing) + { + list > result; + list > segment; + + point3D p0; + point3D p1; + + int node = FiberList[f].n0; + p0 = NodeList[node].p; + + for(unsigned int p=0; p initNodeList(rtsFiberNetwork* network); + void initTopologyGraph(vector* Nodes, vector* Edges, rtsFiberNetwork* Network); + void MapDeviationToNetwork(vector* source); + void topLabelNodes(vector* N0, vector* N1, float epsilon); + bool topMergeNode(vector* NodeList, vector* EdgeList, unsigned int node); + int topCollapse(vector* Nodelist, vector* EdgeList); + void topComputeMergePoints(vector* NodeList); + bool topDetectEdge(vector* NodeList, vector* EdgeList, unsigned int node0, unsigned int node1); + bool topDeleteEdge(vector* NodeList, vector* EdgeList, unsigned int node0, unsigned int node1); + unsigned int topGetNearestConnected(vector* NodeList, unsigned int node, bool must_be_compatible); + void ComputeBoundingVolume(); + float GeometryMetric(rtsFiberNetwork* network, float std) + { + //At this point each vertex in the network has an error in the range of [0 1] + //This function computes the average error of each fiber and the entire network based + //on the error at each vertex. + unsigned int f, p; + float fiber_length; + float fiber_metric; + float total_metric = 0.0; + float total_length = 0.0; + point3D p0, p1; + float e0, e1; + int node; + + //compute the metric for every fiber in the network + for(f=0; fFiberList.size(); f++) + { + fiber_metric = 0.0; + fiber_length = 0.0; + + //start at the first node + node = network->FiberList[f].n0; + p0 = network->NodeList[node].p; + e0 = network->NodeList[node].error; + + //iterate through every intermediary node + for(p=0; pFiberList[f].errorList.size(); p++) + { + p1 = network->FiberList[f].pointList[p]; + e1 = network->FiberList[f].errorList[p]; + //keep a running sum of both the fiber length and the average error + fiber_length += (p0 - p1).Length(); + fiber_metric += (e0+e1)/2 * (p0 - p1).Length(); + + //shift + p0 = p1; + e0 = e1; + + } + + //end at the last fiber node + node = network->FiberList[f].n1; + p1 = network->NodeList[node].p; + e1 = network->NodeList[node].error; + fiber_length += (p0 - p1).Length(); + fiber_metric += (e0+e1)/2 * (p0 - p1).Length(); + + //compute and store the average fiber error + network->FiberList[f].error = fiber_metric/fiber_length; + + //keep a running total of the network error + //do not include if the fiber is culled + if(network->FiberList[f].error <= network->cull_value) + { + total_length += fiber_length; + total_metric += fiber_metric; + } + } + + //compute the final error for the entire network + return total_metric / total_length; + } + + void MY_ComputeTopology(rtsFiberNetwork* testNetwork, float std); + void BOOST_MapEdgesToFibers(rtsFiberNetwork* network, TopologyGraph& G) + { + /* + //initialize the fiber mapping to -1 + for(int f=0; fFiberList.size(); f++) + network->FiberList[f].mapped_to = -1; + + //for each edge in the graph, set the mapping + graph_traits::edge_iterator ei, ei_end; + int map; + for(tie(ei, ei_end) = edges(G); ei != ei_end; ei++) + { + map = get(edge_color_t(), G, *ei); + network->FiberList[map].mapped_to = 1; + } + */ + + } + void BOOST_InvalidateValence2Vertices(TopologyGraph& G) + { + graph_traits::vertex_iterator vi, vi_end; + for(boost::tuples::tie(vi, vi_end) = vertices(G); vi!= vi_end; vi++) + { + if(degree(*vi, G) == 2) + { + put(vertex_color_t(), G, *vi, -1); + cout<<"invalidated"< > BOOST_FindNeighbors(TopologyGraph G, TopologyVertex node) + { + //Finds all colored vertices that can be reached from "node" + pair edge_map; + list > result; + + do{ + global_Predecessors.clear(); + global_Predecessors.resize(num_vertices(G)); + global_Source = node; + global_EdgeSequence.clear(); + global_EdgeDescriptorSequence.clear(); + global_Weight = 0.0; + boundary_bfs_visitor vis; + try{ + breadth_first_search(G, global_Source, visitor(vis)); + } + catch(TopologyVertex& vert_id) + { + edge_map.first = vert_id; + edge_map.second = global_EdgeSequence; + result.push_back(edge_map); + global_NeighborWeights.push_back(global_Weight); + //clear_vertex(vert_id, G); + remove_edge(global_EdgeDescriptorSequence.front(), G); + } + }while(!global_EdgeSequence.empty()); + return result; + } + + TopologyGraph BOOST_FindCoreGraph(TopologyGraph G) + { + TopologyGraph result; + //first insert all vertices + graph_traits::vertex_iterator vi, vi_end; + graph_traits::vertex_descriptor v; + for(boost::tuples::tie(vi, vi_end) = vertices(G); vi!=vi_end; vi++) + { + v = add_vertex(result); + put(vertex_color_t(), result, v, get(vertex_color_t(), G, *vi)); + put(vertex_position_t(), result, v, get(vertex_position_t(), G, *vi)); + } + + //for each vertex in the graph, find all of the neighbors + list > neighborhood; + list >::iterator ni; + list::iterator wi; + pair e; + for(boost::tuples::tie(vi, vi_end) = vertices(G); vi!=vi_end; vi++) + { + //only look for neighbors if the vertex has a color + if(get(vertex_color_t(), G, *vi) >= 0) + { + global_NeighborWeights.clear(); + neighborhood = BOOST_FindNeighbors(G, *vi); + for(ni = neighborhood.begin(), wi = global_NeighborWeights.begin(); ni != neighborhood.end(); ni++, wi++) + { + e = add_edge(*vi, (*ni).first, result); + put(edge_color_t(), result, e.first, (*ni).second); + //cout<<"Inserting weight: "<<*wi< e_Tc, e_GTc; + graph_traits::edge_iterator ei, ei_end, ei_temp; + graph_traits::vertex_descriptor v0, v1; + + //cout<<"Tc"<NodeList.size()<NodeList.size()); + + //for each fiber in the graph, create an edge + int n0, n1; + typedef std::pair Edge; + typedef std::pair::edge_descriptor, bool> EdgeIDType; + EdgeIDType edge_id; + Edge new_edge; + EdgeProperties ep; + EdgeSequence edge_sequence; + for(unsigned int f=0; fFiberList.size(); f++) + { + if(!network->isCulled(f)) + { + n0 = network->FiberList[f].n0; + n1 = network->FiberList[f].n1; + new_edge = Edge(n0, n1); + edge_id = add_edge(new_edge.first, new_edge.second, network->FiberList[f].error*network->FiberList[f].length, g); + + //add the starting edge color + edge_sequence.clear(); + edge_sequence.push_back(f); + put(edge_color_t(), g, edge_id.first, edge_sequence); + } + else + cout<<"culled"<::type PositionMap; + PositionMap positions = get(vertex_position_t(), g); + + for(unsigned int v=0; vNodeList.size(); v++) + positions[v] = network->NodeList[v].p; + + return g; + } + void BOOST_KD_ColorGraph(TopologyGraph& G, TopologyGraph& compareTo, float sigma) + { + //get the number of vertices in compareTo + int verts = num_vertices(compareTo); + //allocate enough space in an ANN array + ANNpointArray dataPts = annAllocPts(verts, 3); + + //get the vertex positions + typedef property_map::type PositionMap; + PositionMap positions = get(vertex_position_t(), compareTo); + //insert the positions into the ANN list + for(int v=0; v::type ColorMap; + ColorMap colors = get(vertex_color_t(), G); + //get the vertex compatibility map + //typedef property_map::type CompatibilityMap; + //CompatibilityMap compatibility = get(vertex_compatibility_t(), G); + //get the index property map + //typedef property_map::type IndexMap; + //IndexMap index = get(vertex_index, G); + + //query each vertex in G + typedef graph_traits::vertex_iterator vertex_iter; + std::pair vp; + point3D pos; + for (vp = vertices(G); vp.first != vp.second; ++vp.first) + { + pos = positions[*vp.first]; + queryPt[0] = pos.x; + queryPt[1] = pos.y; + queryPt[2] = pos.z; + //perform the 1-NN search + kdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + //if the distance is less than sigma, label as valid + if(sqrt(nearestDist[0]) < sigma) + { + colors[*vp.first] = nearestIdx[0]; + //compatibility[*vp.first] = nearestIdx[0]; + } + else + { + colors[*vp.first] = -1; + //compatibility[*vp.first] = -1; + } + } + } + void BOOST_KD_ColorGraphs(TopologyGraph& T, TopologyGraph& GT, float sigma) + { + /*Colors both graphs for the connectivity metric: + 1) Create a kd-tree using the vertices in GT + 2) Find a vertex in GT near each vertex in T + a) If a vertex in GT isn't found, the vertex in T is assigned a color of -1 + b) If a vertex in GT is found, the vertex in T is assigned the corresponding index + 3) Vertices in GT are assigned their own index if they are found (-1 if they aren't) + */ + + //initialize each vertex in T with a color of -1 + graph_traits::vertex_iterator vi, vi_end; + for (boost::tuples::tie(vi, vi_end) = vertices(T); vi != vi_end; vi++) + { + put(vertex_color_t(), T, *vi, -1); + } + + //initialize each vertex in GT with -1 + for (boost::tuples::tie(vi, vi_end) = vertices(GT); vi != vi_end; vi++) + { + put(vertex_color_t(), GT, *vi, -1); + } + //BOOST_PrintGraph(T); + + //CREATE THE KD TREE REPRESENTING GT + //get the number of vertices in GT + int verts = num_vertices(T); + //allocate enough space in an ANN array + ANNpointArray dataPts = annAllocPts(verts, 3); + + //get the vertex positions + typedef property_map::type PositionMap; + PositionMap positions = get(vertex_position_t(), T); + //insert the positions into the ANN list + for(int v=0; v::type ColorMap; + //ColorMap colors = get(vertex_color_t(), G); + + //query each vertex in T + //std::pair vp; + point3D pos; + int num_undetected = 0; + for (boost::tuples::tie(vi, vi_end) = vertices(GT); vi != vi_end; vi++) + { + //pos = positions[*vp.first]; + pos = get(vertex_position_t(), GT, *vi); + queryPt[0] = pos.x; + queryPt[1] = pos.y; + queryPt[2] = pos.z; + //perform the 1-NN search + kdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + //if the distance is less than sigma, label as valid + if(sqrt(nearestDist[0]) < sigma) + { + //colors[*vp.first] = nearestIdx[0]; + //color T + put(vertex_color_t(), T, nearestIdx[0], (int)*vi); + //color GT + put(vertex_color_t(), GT, *vi, (int)*vi); + //compatibility[*vp.first] = nearestIdx[0]; + } + else + { + //put(vertex_color_t(), T, *vi, -1); + num_undetected++; + //compatibility[*vp.first] = -1; + } + } + + //find undetected and falsely detected nodes + //cout<<"Number of Undetected Nodes: "<::vertex_iterator vi, vtmp, vi_end; + boost::tuples::tie(vi, vi_end) = vertices(G); + pair::edge_descriptor, bool> e_remove; + graph_traits::adjacency_iterator ai, ai_end; + graph_traits::vertex_descriptor v_remove; + point3D vp, ap; + vector3D v_dist; + float dist; + while(vi != vi_end) + { + //get the color of the vertex + vp = get(vertex_position_t(), G, *vi); + + //check every adjacent vertex + vertex_error = false; + for(boost::tuples::tie(ai, ai_end) = adjacent_vertices(*vi, G); ai != ai_end; ai++) + { + //get the position of the adjacent vertex + ap = get(vertex_position_t(), G, *ai); + //find the distance between the two points + v_dist = vp - ap; + dist = v_dist.Length(); + + if(dist <= sigma) + { + error_detected = true; + vertex_error = true; + v_remove = *vi; + e_remove = edge(*vi, *ai, G); + } + } + if(vertex_error) + { + //merge the vertices along edge "e_remove" + BOOST_MergeVertices(G, e_remove.first, v_remove); + + //refresh the vertex iterator + boost::tuples::tie(vtmp, vi_end) = vertices(G); + } + else + vi++; + } + + + } + void BOOST_SaveGTIndices(TopologyGraph& G) + { + //changes the ground truth colors to vertex indices + graph_traits::vertex_iterator vi, vi_end; + int index; + for(boost::tuples::tie(vi, vi_end) = vertices(G); vi != vi_end; vi++) + { + //only change the color if the vertex is valid + if(get(vertex_color_t(), G, *vi) >= 0) + { + index = get(vertex_index_t(), G, *vi); + put(vertex_color_t(), G, *vi, index); + } + } + } + void BOOST_ApplyColorsToNetwork(TopologyGraph& from, rtsFiberNetwork* to) + { + //get the index property + typedef property_map::type IndexMap; + IndexMap index = get(vertex_index, from); + //get the vertex color (valid/invalid) + typedef property_map::type ColorMap; + ColorMap colors = get(vertex_color_t(), from); + + //go through each vertex and apply the appropriate color to the associated fiber node + typedef graph_traits::vertex_iterator vertex_iter; + std::pair vp; + int idx; + for (vp = vertices(from); vp.first != vp.second; ++vp.first) + { + idx = index[*vp.first]; + to->NodeList[idx].color = colors[*vp.first]; + } + //cout << index[*vp.first] <<":"<::edge_descriptor e, graph_traits::vertex_descriptor v) + { + //cout<<"MERGE VERTEX TEST-------------------------------"<::vertex_descriptor v0, v1; + v0 = v; + v1 = target(e, G); + if(v1 == v) v1 = source(e, G); + + //we will merge v0 to v1 (although the order doesn't matter) + //well, okay, it does matter as far as the vertex position is concerned...maybe I'll fix that + + //get the current edge's colors + EdgeSequence removed_colors = get(edge_color_t(), G, e); + float removed_weight = get(edge_weight_t(), G, e); + + //first delete the current edge + remove_edge(e, G); + + //for each edge coming into v0, create a corresponding edge into v1 + pair::edge_descriptor, bool> new_edge; + graph_traits::in_edge_iterator ei, ei_end; + EdgeSequence new_colors; + float weight; + for(boost::tuples::tie(ei, ei_end) = in_edges(v0, G); ei != ei_end; ei++) + { + //create a new edge for each edge coming in to v0 + new_edge = add_edge(source(*ei, G), v1, G); + new_colors = get(edge_color_t(), G, *ei); + + //add the removed colors to the new edge + new_colors.insert(new_colors.end(), removed_colors.begin(), removed_colors.end()); + + weight = get(edge_weight_t(), G, *ei) + removed_weight; + put(edge_weight_t(), G, new_edge.first, weight); + put(edge_color_t(), G, new_edge.first, new_colors); + } + + //now remove v0 from the graph + clear_vertex(v0, G); + //remove_vertex(v0, G); + + //cout<<"END MERGE VERTEX TEST-------------------------------"<::vertex_iterator vi, vtmp, vi_end; + boost::tuples::tie(vi, vi_end) = vertices(G); + int color; + float min_weight; + float weight; + graph_traits::edge_descriptor e_remove; + graph_traits::in_edge_iterator ei, ei_end; + //graph_traits::vertex_descriptor v_remove; + bool vertex_removal; + + //merge neighbors that are both invalid + //for each vertex + for(boost::tuples::tie(vi, vi_end) = vertices(G); vi!=vi_end; vi++) + { + vertex_removal = false; + color = get(vertex_color_t(), G, *vi); + + //if the vertex is invalid + if(color < 0) + { + + //get the incident edge with the lowest weight + min_weight = 99999.0; + for(boost::tuples::tie(ei, ei_end) = in_edges(*vi, G); ei != ei_end; ei++) + { + if(get(vertex_color_t(), G, source(*ei, G)) == get(vertex_color_t(), G, target(*ei, G))) + { + vertex_removal = true; + weight = get(edge_weight_t(), G, *ei); + cout<<"weight: "<::vertex_iterator vi, vtmp, vi_end; + boost::tuples::tie(vi, vi_end) = vertices(G); + graph_traits::edge_descriptor max_edge; + graph_traits::in_edge_iterator ei, ei_end; + graph_traits::vertex_descriptor v_next, v_this; + graph_traits::adjacency_iterator ai, ai_end; + int v_degree; + + //for each vertex + for(boost::tuples::tie(vi, vi_end) = vertices(G); vi!=vi_end; vi++) + { + //if the vertex is invalid AND degree=1 + v_degree = degree(*vi, G); + if(get(vertex_color_t(), G, *vi) < 0 && v_degree == 1) + { + //delete the entire spur + v_this = *vi; + do{ + //get the adjacent vertex + boost::tuples::tie(ai, ai_end) = adjacent_vertices(v_this, G); + //cout<<*ai<::vertex_iterator vi, vtmp, vi_end; + boost::tuples::tie(vi, vi_end) = vertices(G); + int color; + pair::edge_descriptor, bool> e_remove; + graph_traits::adjacency_iterator ai, ai_end; + graph_traits::vertex_descriptor v_remove; + int num_incident; + bool merge_found; + while(vi != vi_end) + { + //get the color of the vertex + color = get(vertex_color_t(), G, *vi); + //if the vertex is valid + if(color >= 0) + { + //run through each adjacent vertex + merge_found = false; + for(boost::tuples::tie(ai, ai_end) = adjacent_vertices(*vi, G); ai != ai_end; ai++) + { + //if the two vertices are compatible + if(color == get(vertex_color_t(), G, *ai) && (*vi != *ai)) + { + v_remove = *vi; + e_remove = edge(*vi, *ai, G); + merge_found = true; + } + } + if(merge_found) + { + num_incident = degree(v_remove, G); + BOOST_MergeVertices(G, e_remove.first, v_remove); + if(num_incident == 1) merged_edges++; + boost::tuples::tie(vtmp, vi_end) = vertices(G); + } + else vi++; + } + else + vi++; + } + + return merged_edges; + } + + EdgeSequence BOOST_RemoveEdge(int v0, int v1, TopologyGraph& G) + { + //look for a direct edge + pair::edge_descriptor, bool> e0, e1; + e0 = BOOST_SmallestEdge(v0, v1, G); + + //if it exists, remove it and return the edge sequence + EdgeSequence sequence0, sequence1; + sequence0.clear(); + if(e0.second) + { + sequence0 = get(edge_color_t(), G, e0.first); + remove_edge(e0.first, G); + } + //otherwise look for an indirect edge + + else + { + graph_traits::adjacency_iterator ai, ai_end; + //for each vertex adjacent to v0 + for(boost::tuples::tie(ai, ai_end) = adjacent_vertices(v0, G); ai!=ai_end; ai++) + { + //if the adjacent vertex is invalid + if(get(vertex_color_t(), G, *ai) < 0) + { + //see if v1 is also connected + e1 = BOOST_SmallestEdge(v1, *ai, G); + //if it is + if(e1.second) + { + sequence0 = get(edge_color_t(), G, e1.first); + e0 = BOOST_SmallestEdge(v0, *ai, G); + sequence1 = get(edge_color_t(), G, e0.first); + sequence0.insert(sequence0.end(), sequence1.begin(), sequence1.end()); + remove_edge(e0.first, G); + remove_edge(e1.first, G); + return sequence0; + } + } + } + } + return sequence0; + + } + CoreGraphList BOOST_FindCore(TopologyGraph& GT, TopologyGraph& T) + { + CoreGraphList Gc; + + int FP = 0; + int FN = 0; + + //Find edges that exist in both GT and T + EdgeSequence T_Path, GT_Path; + EdgeMapping path_pair; + + //create iterators for the edges in T + graph_traits::edge_iterator ei, ei_end; + + //for each edge in T, remove the corresponding edge in GT (if it exists) and put the details in the core garph list + int v0, v1; + pair::edge_descriptor, bool> e; + pair::edge_descriptor, bool> e_T, e_GT; + + graph_traits::vertex_descriptor s, t; + list::edge_descriptor> to_remove; + list::edge_descriptor>::iterator remove_iter; + + /*//add all edges to a list + list> EdgeList; + for(tie(ei, ei_end) = edges(T); ei!=ei_end; ei++) + { + pair edge_data; + edge_data.first = *ei; + edge_data.second = get(edge_weight_t(), T, *ei); + EdgeList.push_back(edge_data); + } + EdgeList.sort(compare_edges);*/ + + //find direct connections + for(boost::tuples::tie(ei, ei_end) = edges(T); ei!=ei_end; ei++) + { + //get the ID for the corresponding vertices in GT + s = source(*ei, T); + t = target(*ei, T); + v0 = get(vertex_color_t(), T, s); + v1 = get(vertex_color_t(), T, t); + + //if both vertices are valid + if(v0 >= 0 && v1 >= 0) + { + GT_Path = BOOST_RemoveEdge(v0, v1, GT); + //if there was a corresponding edge in GT + if(GT_Path.size() > 0) + { + //get the edge path from T + T_Path = get(edge_color_t(), T, *ei); + //create the path pair and add it to the core graph list + path_pair.first = T_Path; + path_pair.second = GT_Path; + Gc.push_back(path_pair); + } + else + FP++; + //mark the edge for removal from T + to_remove.push_back(*ei); + + } + } + for(remove_iter = to_remove.begin(); remove_iter != to_remove.end(); remove_iter++) + remove_edge(*remove_iter, T); + + + //find indirect connections (these are connections that use one invalid point) + graph_traits::adjacency_iterator ai, ai_end; + bool match_found; + while(num_edges(T) > 0) + { + boost::tuples::tie(ei, ei_end) = edges(T); + + //s = valid point, t = invalid + s = source(*ei, T); + if(get(vertex_color_t(), T, s) < 0) + { + t = s; + s = target(*ei, T); + } + else + t = target(*ei, T); + //for each point adjacent to the invalid point + match_found = false; + v0 = get(vertex_color_t(), T, s); + for(boost::tuples::tie(ai, ai_end) = adjacent_vertices(t, T); ai!=ai_end; ai++) + { + //get the edge path from T + e = edge(t, *ai, T); + //make sure that the first and second edges are not the same + if(e.first != *ei) + { + v1 = get(vertex_color_t(), T, *ai); + GT_Path = BOOST_RemoveEdge(v0, v1, GT); + //if there was a corresponding edge in GT + if(GT_Path.size() > 0) + { + match_found = true; + T_Path = get(edge_color_t(), T, e.first); + EdgeSequence temp = get(edge_color_t(), T, *ei); + T_Path.insert(T_Path.end(), temp.begin(), temp.end()); + + //create the path pair and add it to the core graph list + path_pair.first = T_Path; + path_pair.second = GT_Path; + Gc.push_back(path_pair); + + //remove both edges and break + remove_edge(e.first, T); + remove_edge(*ei, T); + break; + } + } + } + if(!match_found) + { + FP++; + remove_edge(*ei, T); + } + } + + cout<<"False Positive Edges in Core: "<::type IndexMap; + IndexMap index = get(vertex_index, G); + + //get the position property + typedef property_map::type PositionMap; + PositionMap positions = get(vertex_position_t(), G); + //get the vertex color (valid/invalid) + typedef property_map::type ColorMap; + ColorMap colors = get(vertex_color_t(), G); + //get vertex compatibility + //typedef property_map::type CompatibilityMap; + //CompatibilityMap compatibility = get(vertex_compatibility_t(), G); + + std::cout << "vertices(g) = "<::vertex_iterator vertex_iter; + std::pair vp; + for (vp = vertices(G); vp.first != vp.second; ++vp.first) + cout << index[*vp.first] <<"("<::type WeightMap; + WeightMap weights = get(edge_weight_t(), G); + + std::cout << "edges(source, dest): weight to_fiber"<::edge_iterator ei, ei_end; + EdgeSequence edge_sequence; + for (boost::tuples::tie(ei, ei_end) = edges(G); ei != ei_end; ++ei) + { + std::cout << "(" << index[source(*ei, G)] << "," << index[target(*ei, G)] << "): "< p0 = NodeList[FiberList[f].n0].p; + point3D p1; + + double length = 0.0; + int num_points = FiberList[f].pointList.size(); + int p; + for(p=0; p FiberList; + vector NodeList; + point3D min_pos; + point3D max_pos; + + //constructors + rtsFiberNetwork() + { + fiber_started = false; + num_points = false; + cull_value = 1.0; + min_pos = point3D(9999, 9999, 9999); + max_pos = point3D(-9999, -9999, -9999); + } + + //get functions + point3D getNodeCoord(int node){return NodeList[node].p;} + point3D getNodeCoord(int fiber, bool node); + point3D getFiberPoint(unsigned int fiber, unsigned int point); + + //network culling + float getCullValue(){return cull_value;} + void setCullValue(float cull){cull_value = cull;} + bool isCulled(int f) + { + if(FiberList[f].error <= cull_value) + return false; + else return true; + } + + //statistics functions + double getTotalLength(); + double getFiberLength(int f) + { + return FiberList[f].length; + } + + + + //drawing functions + void StartFiber(float x, float y, float z) + { + fiber_started = true; + num_points = 1; + + //create a start node and an end node + Node n; + n.p = point3D(x, y, z); + NodeList.push_back(n); + NodeList.push_back(n); + + //create a fiber + Fiber f; + f.n0 = NodeList.size()-2; + f.n1 = NodeList.size()-1; + FiberList.push_back(f); + + } + void StartFiber(int node) + { + fiber_started = true; + num_points = 1; + + //the start node is specified + //specify the end node + Node n = NodeList[node]; + NodeList.push_back(n); + + //create a fiber + Fiber f; + f.n0 = node; + f.n1 = NodeList.size()-1; + FiberList.push_back(f); + } + + void ContinueFiber(float x, float y, float z) + { + if(!fiber_started) + { + StartFiber(x, y, z); + return; + } + num_points++; + + //store the last node coordinate in the fiber list + int f = FiberList.size() - 1; + int n = NodeList.size() - 1; + + if(num_points > 2) + FiberList[f].pointList.push_back(NodeList[n].p); + NodeList[n].p = point3D(x, y, z); + //NodeList[n].p.print(); + //cout<= NodeList.size() - 1) + return; + + int f = FiberList.size() - 1; + + //add the last point + FiberList[f].pointList.push_back(NodeList[FiberList[f].n1].p); + + //attach to the specified node + FiberList[f].n1 = node; + NodeList.pop_back(); + fiber_started = false; + } + void StartFiber(point3D p){StartFiber(p.x, p.y, p.z);} + void ContinueFiber(point3D p){ContinueFiber(p.x, p.y, p.z);} + + void Clear() + { + FiberList.clear(); + NodeList.clear(); + fiber_started = false; + } + + //loading functions + void MergeNodes(float sigma) + { + //create a KD tree for all nodes in the network + //get the number of vertices in GT + int verts = NodeList.size(); + //allocate enough space in an ANN array + ANNpointArray dataPts = annAllocPts(verts, 3); + + //insert the positions into the ANN list + for(int v=0; v NodeMap; + NodeMap.resize(NodeList.size(), -1); + + ANNpoint queryPt = annAllocPt(3); + ANNidxArray nearestIdx = new ANNidx[1]; + ANNdistArray nearestDist = new ANNdist[1]; + + for(unsigned int n=0; nannkSearch(queryPt, 1, nearestIdx, nearestDist); + + NodeMap[n] = nearestIdx[0]; + } + + //set the nodes for each fiber to those mapped + for(unsigned int f=0; f NodeMap; + NodeMap.resize(NodeList.size(), -1); + vector FiberNum; + FiberNum.resize(NodeList.size(), 0); + + vector newNodeList; + + //run through each fiber + int node; + for(unsigned int f=0; f branchIDList; + vector correspondingNodeList; + + vector::iterator iter; + int index; + + char c; + //first pass, get branch points and find bounds + while(!infile.eof()) + { + c = infile.peek(); + if((c <'0' || c > '9') && c != 32) + infile.ignore(9999, '\n'); + else + { + infile>>id; + infile>>type; + infile>>x; + infile>>y; + infile>>z; + infile>>r; + infile>>parent; + + + //root nodes and branch nodes + if(parent != id-1 && parent != -1) + branchIDList.push_back(parent); + if(x < min_pos.x) min_pos.x = x; + if(y < min_pos.y) min_pos.y = y; + if(z < min_pos.z) min_pos.z = z; + if(x > max_pos.x) max_pos.x = x; + if(y > max_pos.y) max_pos.y = y; + if(z > max_pos.z) max_pos.z = z; + } + }//end while + //sort the branch points + sort(branchIDList.begin(), branchIDList.end()); + //set the number of corresponding nodes + correspondingNodeList.resize(branchIDList.size()); + + //second pass, build the tree + infile.close(); + infile.clear(); + infile.open(filename.c_str()); + while(!infile.eof()) + { + c = infile.peek(); + if((c <'0' || c > '9') && c != 32) + infile.ignore(9999, '\n'); + else + { + infile>>id; + infile>>type; + infile>>x; + infile>>y; + infile>>z; + infile>>r; + infile>>parent; + + + //if we are starting a new fiber + if(parent == -1) + StartFiber(x, y, z); + else + { + //see if the parent is in the branch list + iter = find(branchIDList.begin(), branchIDList.end(), parent); + //if the parent point is a branch point + if(iter != branchIDList.end()) + { + index = iter - branchIDList.begin(); + StartBranch(correspondingNodeList[index]); + ContinueFiber(x, y, z); + } + else + ContinueFiber(x, y, z); + } + //if the current node is in the branch list + iter = find(branchIDList.begin(), branchIDList.end(), id); + if(iter != branchIDList.end()) + { + //get the id and store the corresponding node value + index = iter - branchIDList.begin(); + correspondingNodeList[index] = NodeList.size()-1; + } + } + }//end while + refreshFiberLengths(); + refreshIncidence(); + + } + + void LoadOBJ(string filename) + { + //first load the OBJ file + rtsOBJ objFile; + objFile.LoadFile(filename.c_str()); + + /*Create two lists. Each element represents a point in the OBJ file. We will + first go through each fiber (line) and find the vertex associated with the two ends of the fiber. + The validList value is "true" if the associated point in the OBJ file is a junction. It is "false" + if the point is just an intermediate point on a fiber. The pointList value stores the new value of the junction + in the NodeList of this FiberNetwork structure.*/ + + vector validList; + validList.assign(objFile.v_list.size(), false); + vector pointList; + pointList.assign(objFile.v_list.size(), 0); + + /*Run through each fiber: + 1) See if each fiber node has already been created by looking at validList (true = created) + 2) If the node has not been created, create it and set validList to true and pointList to the node index in this structure. + 3) Create an empty fiber with the appropriate node values assigned.*/ + unsigned int f; + unsigned int line_verts; + unsigned int v0, vn; + Node node_new; + for(f=0; f p){Translate(p.x, p.y, p.z);} + + void Oscillate(float frequency, float magnitude) + { + //impliment a sinusoidal oscillation along each fiber + vector3D side(0.0, 0.0, 1.0); + vector3D dir, normal; + point3D p0, p1; + float t; + + //for each fiber + for(unsigned int f=0; f newFiberList; + for(unsigned int f=0; f p0 = NodeList[n0].p; + point3D p1 = NodeList[n1].p; + + if(p0.x > px && p0.y > py && p0.z > pz && + p1.x > px && p1.y > py && p1.z > pz && + p0.x < px+sx && p0.y < py+sy && p0.z < pz+sz && + p1.x < px+sx && p1.y < py+sy && p1.z < pz+sz) + { + newFiberList.push_back(FiberList[f]); + } + } + FiberList.clear(); + FiberList = newFiberList; + RemoveExcessNodes(); + ComputeBoundingVolume(); + } + void ThresholdLength(float min_length, float max_length = 99999) + { + vector newFiberList; + refreshFiberLengths(); + for(unsigned int f=0; f min_length && FiberList[f].length < max_length) + { + newFiberList.push_back(FiberList[f]); + } + } + FiberList.clear(); + FiberList = newFiberList; + RemoveExcessNodes(); + ComputeBoundingVolume(); + } + void ThresholdSpines(float min_length) + { + vector newFiberList; + refreshIncidence(); + refreshFiberLengths(); + for(unsigned int f=0; f min_length || (NodeList[FiberList[f].n0].incident > 1 && NodeList[FiberList[f].n1].incident > 1)) + { + newFiberList.push_back(FiberList[f]); + } + } + FiberList.clear(); + FiberList = newFiberList; + RemoveExcessNodes(); + ComputeBoundingVolume(); + + } + //subdivision + void SubdivideNetwork(float spacing) + { + list > subdivided; + list >::iterator p; + for(unsigned int f=0; f p0, p1; + vector > newPointList; + for(unsigned int f=0; f= spacing ) + { + newPointList.push_back(p1); + } + } + FiberList[f].pointList = newPointList; + } + } + //network comparison + CoreGraphList CompareNetworks(rtsFiberNetwork* testNetwork, float sigmaG, float sigmaC, float &gFPR, float &gFNR, float &cFPR, float &cFNR) + { + //create point clouds that densely sample each network + vector netPointList0, netPointList1; + netPointList0 = getNetPointSamples(sigmaG); + netPointList1 = testNetwork->getNetPointSamples(sigmaG); + + //compute the L1 distance between vertices in one network to the point cloud representing the other network + KD_ComputeEnvelopeDistance(testNetwork, &netPointList0, sigmaG); + KD_ComputeEnvelopeDistance(this, &netPointList1, sigmaG); + + //compute the geometry metric using the distance values for each vertex + //float FPR, FNR; + gFNR = GeometryMetric(this, sigmaG); + gFPR = GeometryMetric(testNetwork, sigmaG); + + CoreGraphList core; + core = NEW_ComputeTopology(testNetwork, sigmaC); + + //Changes made by James Burck (Thanks James!)-------- + float TP = (float)core.size(); + + float TPandFP = (float)FiberList.size(); // formerly P, actaully TPandFN + float TPandFN = (float)testNetwork->FiberList.size(); // actually TPandFP + + cFNR = (TPandFN - TP) / TPandFN; + cFPR = (TPandFP - TP) / TPandFP; + //--------------------------------------------------- + + return core; + } + + +}; + +void rtsFiberNetwork::initTopologyGraph(vector* Nodes, vector* Edges, rtsFiberNetwork* network) +{ + /*This function constructs a graph based on the given network.*/ + Nodes->clear(); + Edges->clear(); + + topologyNode node; + topologyEdge edge; + //for each node in the fiber network, construct a topologyNode + for(unsigned int n=0; nNodeList.size(); n++) + { + node.compatible = 0; + node.label = RTS_TOPOLOGY_NODE_INVALID; + node.p = network->NodeList[n].p; + Nodes->push_back(node); + } + + //now fill in all the edges + for(unsigned int f=0; fFiberList.size(); f++) + { + edge.n0 = network->FiberList[f].n0; + edge.error = network->FiberList[f].error; + edge.n1 = network->FiberList[f].n1; + edge.label = RTS_TOPOLOGY_EDGE_EXIST; + + //attach the edge to each connected node in the node list + (*Nodes)[edge.n0].connections.push_back(Edges->size()); + (*Nodes)[edge.n1].connections.push_back(Edges->size()); + + //insert the edge into the list + Edges->push_back(edge); + } + +} + +void rtsFiberNetwork::BF_ComputeL1Distance(vector* N0, vector*N1) +{ + unsigned int i, j; + vector3D v; + float dist; + for(i=0; isize(); i++) + { + for(j=0; jsize(); j++) + { + v = (*N0)[i].p - (*N1)[j].p; + dist = v.Length(); + if(dist < (*N0)[i].dist) + (*N0)[i].dist = dist; + if(dist < (*N1)[j].dist) + (*N1)[j].dist = dist; + } + } +} + +void rtsFiberNetwork::BD_ComputeL1Distance(vector* N0, vector*N1) +{ + //build the point arrays + ANNpointArray dataPts0 = annAllocPts(N0->size(), 3); + for(unsigned int i=0; isize(); i++) + { + dataPts0[i][0] = (*N0)[i].p.x; + dataPts0[i][1] = (*N0)[i].p.y; + dataPts0[i][2] = (*N0)[i].p.z; + } + ANNpointArray dataPts1 = annAllocPts(N1->size(), 3); + for(unsigned int i=0; isize(); i++) + { + dataPts1[i][0] = (*N1)[i].p.x; + dataPts1[i][1] = (*N1)[i].p.y; + dataPts1[i][2] = (*N1)[i].p.z; + } + + //create ANN variables + ANNbd_tree* bdTree; + ANNpoint queryPt = annAllocPt(3); + ANNidxArray nearestIdx = new ANNidx[1]; + ANNdistArray nearestDist = new ANNdist[1]; + + //compare network 0 to network 1 + //bdTree = new ANNkd_tree(dataPts0, N0->size(), 3); + bdTree = new ANNbd_tree(dataPts0, N0->size(), 3); + for(unsigned int i=0; isize(); i++) + { + queryPt[0] = (*N1)[i].p.x; + queryPt[1] = (*N1)[i].p.y; + queryPt[2] = (*N1)[i].p.z; + bdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + (*N1)[i].dist = sqrtf((float)nearestDist[0]); + } + delete bdTree; + + //compare network 1 to network 0 + bdTree = new ANNbd_tree(dataPts1, N1->size(), 3); + for(unsigned int i=0; isize(); i++) + { + queryPt[0] = (*N0)[i].p.x; + queryPt[1] = (*N0)[i].p.y; + queryPt[2] = (*N0)[i].p.z; + bdTree->annkSearch(queryPt, 1, nearestIdx, nearestDist); + (*N0)[i].dist = sqrtf((float)nearestDist[0]); + } + delete bdTree; + + annClose(); +} + +void rtsFiberNetwork::MapDeviationToNetwork(vector* source) +{ + + +} + +void rtsFiberNetwork::topLabelNodes(vector* N0, vector* N1, float sigma) +{ + unsigned int i0, i1; + vector3D v; + float min_d; + unsigned int min_i; + for(i0=0; i0 < N0->size(); i0++) + { + v = (*N0)[i0].p - (*N1)[0].p; + min_d = v.Length(); + min_i = 0; + for(i1=0; i1 < N1->size(); i1++) + { + v = (*N0)[i0].p - (*N1)[i1].p; + if(v.Length() < min_d) + { + min_d = v.Length(); + min_i = i1; + } + } + //if the minimum distance from point i0 is less than sigma + if(min_d < sigma) + { + (*N0)[i0].label = RTS_TOPOLOGY_NODE_VALID; + (*N0)[i0].compatible = min_i; + } + } +} + +bool rtsFiberNetwork::topDetectEdge(vector* NodeList, vector* EdgeList, unsigned int node0, unsigned int node1) +{ + //This function determines if there is an edge linking node0 and node1 + list::iterator i; + for(i = (*NodeList)[node0].connections.begin(); i!=(*NodeList)[node0].connections.end(); i++) + { + if( ((*EdgeList)[*i].n0 == node0 && (*EdgeList)[*i].n1 == node1) || + ((*EdgeList)[*i].n0 == node1 && (*EdgeList)[*i].n1 == node0) ) + return true; + } + return false; +} +bool rtsFiberNetwork::topDeleteEdge(vector* NodeList, vector* EdgeList, unsigned int node0, unsigned int node1) +{ + //this function deletes the first edge found linking node0 and node1 + list::iterator i; + unsigned int edge_id; + for(i = (*NodeList)[node0].connections.begin(); i!=(*NodeList)[node0].connections.end(); i++) + { + if( (*EdgeList)[*i].n0 == node1 || (*EdgeList)[*i].n1 == node1 ) + { + //delete the edge + + edge_id = *i; + //remove the edge from node0 and node1 + (*NodeList)[node0].connections.remove(edge_id); + (*NodeList)[node1].connections.remove(edge_id); + //remove the edge + (*EdgeList)[edge_id].label = RTS_TOPOLOGY_EDGE_NOEXIST; + return true; + } + } + return false; + +} +bool rtsFiberNetwork::topMergeNode(vector* NodeList, vector* EdgeList, unsigned int node) +{ + /*this function merges a specific node with it's neighbor based on the following rules: + 1) If the node is invalid, remove adjacent edge with the highest error + 2) If the node is valid, merge with adjacent compatible node, removing the edge with the largest error + */ + + //if the node doesn't exist, just return. + if( (*NodeList)[node].label == RTS_TOPOLOGY_NODE_NOEXIST) return false; + //if this node isn't connected to anything, just remove the node + if( (*NodeList)[node].connections.size() == 0) + { + (*NodeList)[node].label = RTS_TOPOLOGY_NODE_NOEXIST; + return false; + } + + //FIND THE DESTINATION NODE + //create the destination node + unsigned int edge_to_remove; + + //if the node is invalid, find the edge with the highest error + if( (*NodeList)[node].label == RTS_TOPOLOGY_NODE_INVALID) + { + list::iterator i; + float highest_error = 0.0; + for(i=(*NodeList)[node].connections.begin(); i!=(*NodeList)[node].connections.end(); i++) + { + //if the current edge has a higher error, record it + if((*EdgeList)[(*i)].error >= highest_error) + { + highest_error = (*EdgeList)[(*i)].error; + edge_to_remove = (*i); + } + } + } + + //if the node is valid, find the compatible edge with the highest error + if( (*NodeList)[node].label == RTS_TOPOLOGY_NODE_VALID) + { + list::iterator i; + float highest_error = 0.0; + bool compatible_detected = false; + unsigned int node0, node1; + for(i=(*NodeList)[node].connections.begin(); i!=(*NodeList)[node].connections.end(); i++) + { + node0 = (*EdgeList)[(*i)].n0; + node1 = (*EdgeList)[(*i)].n1; + //find a compatible edge with the highest weight + if((*NodeList)[node0].label == RTS_TOPOLOGY_NODE_VALID && (*NodeList)[node1].label == RTS_TOPOLOGY_NODE_VALID && + (*NodeList)[node0].compatible == (*NodeList)[node1].compatible && (*EdgeList)[(*i)].error >= highest_error) + { + highest_error = (*EdgeList)[(*i)].error; + edge_to_remove = (*i); + compatible_detected = true; + } + } + //if a compatible node was not attached, just leave the node and return + if(!compatible_detected) return false; + } + + //PERFORM THE MERGE + + //find the node that we are merging to + unsigned int merge_to; + if((*EdgeList)[edge_to_remove].n0 == node) + merge_to = (*EdgeList)[edge_to_remove].n1; + else + merge_to = (*EdgeList)[edge_to_remove].n0; + + list::iterator i; + //remove the edge from 'node' + for(i = (*NodeList)[node].connections.begin(); i!=(*NodeList)[node].connections.end(); i++) + if((*i) == edge_to_remove) + { + (*NodeList)[node].connections.erase(i); + break; + } + //remove the edge from 'merge_to' + for(i = (*NodeList)[merge_to].connections.begin(); i!=(*NodeList)[merge_to].connections.end(); i++) + if((*i) == edge_to_remove) + { + (*NodeList)[merge_to].connections.erase(i); + break; + } + + //update all of the edges connected to 'node' + for(i = (*NodeList)[node].connections.begin(); i!=(*NodeList)[node].connections.end(); i++) + { + if((*EdgeList)[(*i)].n0 == node) + (*EdgeList)[(*i)].n0 = merge_to; + else + (*EdgeList)[(*i)].n1 = merge_to; + } + //add all edges in 'node' to the edges in 'merge_to' + for(i = (*NodeList)[node].connections.begin(); i!=(*NodeList)[node].connections.end(); i++) + { + (*NodeList)[merge_to].connections.push_back( (*i) ); + } + //sort the list and remove duplicates + //duplicates occur if two merged points were connected by multiple edges + (*NodeList)[merge_to].connections.sort(); + (*NodeList)[merge_to].connections.unique(); + + //get rid of 'node' + (*NodeList)[node].connections.clear(); + (*NodeList)[node].label = RTS_TOPOLOGY_NODE_NOEXIST; + + //remove the edge + (*EdgeList)[edge_to_remove].label = RTS_TOPOLOGY_EDGE_NOEXIST; + return true; +} +int rtsFiberNetwork::topCollapse(vector* NodeList, vector*EdgeList) +{ + unsigned int n; + unsigned int topology_changes = 0; + unsigned int num_connections; + bool node_merged = false; + for(n=0; nsize(); n++) + { + //if this node is the end of a barb, mark it as a topology change + num_connections = (*NodeList)[n].connections.size(); + node_merged = topMergeNode(NodeList, EdgeList, n); + if(num_connections == 1 && node_merged == true) + topology_changes++; + + } + + return topology_changes; +} +void rtsFiberNetwork::MY_ComputeTopology(rtsFiberNetwork* testNetwork, float sigma) +{ + //initialize the topology graphs + vector GT_nodes; + vector GT_edges; + initTopologyGraph(>_nodes, >_edges, this); + vector T_nodes; + vector T_edges; + initTopologyGraph(&T_nodes, &T_edges, testNetwork); + + //label the nodes in each list as VALID or INVALID + //this function also determines node compatibility in the Test array + topLabelNodes(>_nodes, &T_nodes, sigma); + topLabelNodes(&T_nodes, >_nodes, sigma); + + //copy the error to the fiber networks + for(unsigned int i=0; iNodeList[i].error = 0.0; + else + testNetwork->NodeList[i].error = 1.0; + } + for(unsigned int i=0; iNodeList[i].error = 0.0; + else + this->NodeList[i].error = 1.0; + } + + unsigned int FP_edges = topCollapse(&T_nodes, &T_edges); + unsigned int FN_edges = topCollapse(>_nodes, >_edges); + + //mark all the nodes in T as invalid again + //this will be used to compute topological errors + unsigned int n; + for(n=0; n validPoints; + for(n=0; n::iterator i; + unsigned int last; + unsigned int topErrors = 0; + for(i = validPoints.begin(); i != validPoints.end(); i++) + { + if(i != validPoints.begin()) + if(*i == last) + topErrors++; + last = *i; + } + + + + cout<<"False Positive Edges: "<NodeList.size()); + + //add all edges in G to result (based on compatibility (original index), NOT current index) + graph_traits::edge_iterator ei, ei_end; + graph_traits::vertex_descriptor v0, v1; + int idx0, idx1; + for(boost::tuples::tie(ei, ei_end) = edges(G); ei != ei_end; ei++) + { + v0 = source(*ei, G); + v1 = target(*ei, G); + idx0 = get(vertex_color_t(), G, v0); + idx1 = get(vertex_color_t(), G, v1); + add_edge(idx0, idx1, result); + } + return result; + +} +point3D rtsFiberNetwork::getNodeCoord(int fiber, bool node) +{ + //Return the node coordinate based on a fiber (edge) + //node value is false = source, true = dest + + if(node) + return getNodeCoord(FiberList[fiber].n1); + else + return getNodeCoord(FiberList[fiber].n0); + +} + +double rtsFiberNetwork::getTotalLength() +{ + //go through each fiber and measure the length + int num_fibers = FiberList.size(); + double length = 0.0; + int f; + for(f=0; f rtsFiberNetwork::getFiberPoint(unsigned int fiber, unsigned int point) +{ + if(point == 0) + return NodeList[FiberList[fiber].n0].p; + if(point == FiberList[fiber].pointList.size() + 1) + return NodeList[FiberList[fiber].n1].p; + + return FiberList[fiber].pointList[point-1]; + +} + +void rtsFiberNetwork::ComputeBoundingVolume() +{ + //find the bounding volume for the nodes + min_pos = NodeList[0].p; + max_pos = NodeList[0].p; + for(unsigned int n=0; n max_pos.x) + max_pos.x = NodeList[n].p.x; + if(NodeList[n].p.y > max_pos.y) + max_pos.y = NodeList[n].p.y; + if(NodeList[n].p.z > max_pos.z) + max_pos.z = NodeList[n].p.z; + } + + //combine with the bounding volume for the fibers + for(unsigned int f=0; f max_pos.x) + max_pos.x = FiberList[f].pointList[p].x; + if(FiberList[f].pointList[p].y > max_pos.y) + max_pos.y = FiberList[f].pointList[p].y; + if(FiberList[f].pointList[p].z > max_pos.z) + max_pos.z = FiberList[f].pointList[p].z; + } + } + +} +void rtsFiberNetwork::Translate(float x, float y, float z) +{ + //translates the network while maintaining all connectivity + + //create a translation vector + vector3D translate(x, y, z); + + //translate all nodes + int num_nodes = NodeList.size(); + int n; + for(n=0; n