segmentation / netmets

Commit 86a4e3734c36d9411c5d2fad2488e7e1a81d2d4a

Authored by David Mayerich
2 parents 5c83f00e 30641100

Merge branch 'JACK' into 'master' 

need some data to test on

See merge request !9
Showing 5 changed files with 130 additions and 102 deletions   Show diff stats
Inline Side-by-side
CMakeLists.txt
  Show/Hide comments   View file @86a4e37
... ... @@ -85,3 +85,5 @@ configure_file(data/04_Tb.obj ${CMAKE_CURRENT_BINARY_DIR}/04_Tb.obj @ONLY)
85 85 configure_file(data/1.swc ${CMAKE_CURRENT_BINARY_DIR}/1.swc @ONLY)
86 86 configure_file(data/2.swc ${CMAKE_CURRENT_BINARY_DIR}/2.swc @ONLY)
87 87 configure_file(data/00_GT.swc ${CMAKE_CURRENT_BINARY_DIR}/00_GT.swc @ONLY)
  88 +configure_file(data/01_GT.swc ${CMAKE_CURRENT_BINARY_DIR}/01_GT.swc @ONLY)
  89 +configure_file(data/01_T.swc ${CMAKE_CURRENT_BINARY_DIR}/01_T.swc @ONLY)
... ...
data/00_GT.swc
  Show/Hide comments   View file @86a4e37
  1 +# converted from 00_GT.obj
  2 +1 1 99 677 0 0.5 -1
  3 +2 1 88 563 0 0.5 1
  4 +3 1 229 510 0 0.5 2
  5 +4 1 181 372 0 0.5 3
  6 +5 1 57 347 0 0.5 4
  7 +6 1 316 596 0 0.5 3
  8 +7 1 344 724 0 0.5 6
  9 +8 1 385 425 0 0.5 3
  10 +9 1 581 361 0 0.5 8
  11 +10 1 574 223 0 0.5 9
  12 +11 1 481 165 0 0.5 10
  13 +12 1 341 199 0 0.5 11
  14 +13 1 757 377 0 0.5 9
  15 +14 1 818 229 0 0.5 13
  16 +15 1 981 152 0 0.5 14
  17 +16 1 793 508 0 0.5 13
  18 +17 1 941 600 0 0.5 16
  19 +18 1 679 601 0 0.5 16
  20 +19 1 647 738 0 0.5 18
... ...
data/01_GT.swc 0 → 100644
  Show/Hide comments   View file @86a4e37
  1 +1 1 0.0 0.0 0.0 0.5 -1
  2 +2 3 10 -12 11 0.5 1
  3 +3 3 42 -80 20 0.5 2
  4 +4 3 95 -55 4 0.5 2
  5 +5 3 4 9 0 0.5 1
  6 +6 3 46 55 -15 0.5 5
  7 +7 2 -20 11 30 0.5 1
  8 +8 4 -6 60 -4 0.5 7
  9 +9 4 -44 99 -8 0.5 8
... ...
data/01_T.swc 0 → 100644
  Show/Hide comments   View file @86a4e37
  1 +1 1 0.0 0.0 0.0 0.5 -1
  2 +2 3 10 -12 11 0.5 1
  3 +3 3 30 -40 10 0.5 2
  4 +4 3 4 9 0 0.5 1
  5 +5 3 46 55 -15 0.5 4
  6 +6 2 -20 11 30 0.5 1
  7 +7 4 -6 60 -4 0.5 6
  8 +8 4 -44 99 -8 0.5 7
  9 +9 4 -55 -10 8 0.5 7
  10 +10 4 -66 -10 8 0.5 9
  11 +11 3 40 -80 20 0.5 3
  12 +12 3 50 -90 30 0.5 11
  13 +13 1 20 20 20 0.5 1
... ...
... ... @@ -37,7 +37,6 @@ stim::gl_network<float> _T; //splitted T
37 37  
38 38 // indicator
39 39 unsigned ind = 0; //indicator of mapping
40   -unsigned swc_ind = 0; //indicator of rendering swc file as networks
41 40  
42 41 // relationships
43 42 std::vector<unsigned> _gt_t; // store indices of nearest edge points in _T for _GT
... ... @@ -47,7 +46,8 @@ std::vector&lt;unsigned&gt; _t_gt; // store indices of nearest edge points in _
47 46 float resample_rate = 0.5f; //sample rate for the network (fraction of sigma used as the maximum sample rate)
48 47 float camera_factor = 1.2f; //start point of the camera as a function of X and Y size
49 48 float orbit_factor = 0.01f; //degrees per pixel used to orbit the camera
50   -float zoom_factor = 10.0;
  49 +float zoom_factor = 10.0f;
  50 +float radius_factor = 0.5f;
51 51  
52 52 //mouse click
53 53 bool LButtonDown = false; // true when left button down
... ... @@ -74,7 +74,9 @@ GLuint dlist2;
74 74 //OpenGL objects
75 75 GLuint cmap_tex = 0; //texture name for the color map
76 76  
  77 +float delta;
77 78 float sigma = 3; //default sigma
  79 +float radius = 3; //equals to radius
78 80  
79 81 //sets an OpenGL viewport taking up the entire window
80 82 void glut_render_single_projection(){
... ... @@ -146,11 +148,12 @@ void glut_render(void) {
146 148 glTexEnvf(GL_TEXTURE_ENV, GL_TEXTURE_ENV_MODE, GL_REPLACE); //texture map will be used as the network color
147 149 glBindTexture(GL_TEXTURE_1D, cmap_tex); //bind the Brewer texture map
148 150  
149   - GT.glCylinder(sigma); //render the GT network
  151 + GT.glCylinder(sigma, radius); //render the GT network
150 152  
151 153 glut_render_right_projection(); //set up a projection for the right half of the window
152 154 glut_render_modelview(); //set up the modelview matrix using camera details
153   - T.glCylinder(sigma); //render the T network
  155 + T.glCylinder(sigma, radius); //render the T network
  156 + sigma = radius; // set sigma equal to radius
154 157 }
155 158 }
156 159  
... ... @@ -171,12 +174,12 @@ void glut_render(void) {
171 174 glTexEnvf(GL_TEXTURE_ENV, GL_TEXTURE_ENV_MODE, GL_REPLACE); //texture map will be used as the network color
172 175 glBindTexture(GL_TEXTURE_1D, cmap_tex); //bind the Brewer texture map
173 176  
174   - _GT.glCylinder(sigma); //render the GT network
  177 + _GT.glCylinder(sigma, radius); //render the GT network
175 178  
176 179 glut_render_right_projection(); //set up a projection for the right half of the window
177 180 glut_render_modelview(); //set up the modelview matrix using camera details
178   - _T.glCylinder(sigma); //render the T network
179   -
  181 + _T.glCylinder(sigma, radius); //render the T network
  182 + sigma = radius; //set sigma equal to radius
180 183 }
181 184 else {
182 185 glEnable(GL_DEPTH_TEST); //enable depth
... ... @@ -185,12 +188,13 @@ void glut_render(void) {
185 188 glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT); //clear the screen
186 189  
187 190 //_GT.glRandColorCenterlineGT(dlist1, _gt_t, colormap);
188   - _GT.glRandColorCylinder1(dlist1, _gt_t, colormap, sigma);
  191 + _GT.glRandColorCylinder1(dlist1, _gt_t, colormap, sigma, radius);
189 192  
190 193 glut_render_right_projection(); //set up a projection for the right half of the window
191 194 glut_render_modelview(); //set up the modelview matrix using camera details
192 195 //_T.glRandColorCenterlineT(dlist2, _t_gt, colormap);
193   - _T.glRandColorCylinder2(dlist2, _t_gt, colormap, sigma);
  196 + _T.glRandColorCylinder2(dlist2, _t_gt, colormap, sigma, radius);
  197 + sigma = radius; //set sigma equal to radius
194 198 }
195 199 }
196 200 }
... ... @@ -310,20 +314,10 @@ void glut_mouse(int button, int state, int x, int y){
310 314  
311 315 // define camera move based on mouse wheel move(actually we can combine this with glut_mouse)
312 316 void glut_wheel(int wheel, int direction, int x, int y) {
313   - /*float cam_move_fac; // camera move unit length
314   - stim::vec3<float> eye = cam.getPosition(); // get the camera position (eye point)
315   - stim::vec3<float> focus = cam.getLookAt(); // get the camera focal point
316   - cam_move_fac = fabs(focus[2] - eye[2]) / 50; // divided by 50
317   -
318   - if (direction > 0) // if it is button 3(up)
319   - des -= cam_move_fac;
320   - else // if it is button 4(down)
321   - des += cam_move_fac;*/
322   - float delta;
323 317  
324   - if (direction > 0) // if it is button 3(up)
  318 + if (direction > 0) // if it is button 3(up), move closer
325 319 delta = zoom_factor;
326   - else // if it is button 4(down)
  320 + else // if it is button 4(down), leave farther
327 321 delta = -zoom_factor;
328 322  
329 323 cam.Push(delta);
... ... @@ -332,19 +326,46 @@ void glut_wheel(int wheel, int direction, int x, int y) {
332 326  
333 327 // define keyboard inputs
334 328 void glut_keyboard(unsigned char key, int x, int y){
335   - if(key == 'm') // if keyboard 'm' is pressed
336   - {
337   - if(compareMode && !mappingMode){ // if it is in compare mode
338   - compareMode = false;
339   - mappingMode = true;
340   - }
341   - else{ // if it is in mapping mode
342   - compareMode = true;
343   - mappingMode = false;
344   - }
345   - }
346   - if (key == 27) { // if keyboard "ESC" is pressed
347   - exit(0); // exit
  329 +
  330 + // register different keyboard operation
  331 + switch (key) {
  332 +
  333 + // change render mode
  334 + case 'm': // if keyboard 'm' is pressed, then change render mode
  335 + if (compareMode && !mappingMode && ind) { // if current mode is comparing mode
  336 + compareMode = false;
  337 + mappingMode = true;
  338 + }
  339 + else if (!compareMode && mappingMode && ind) { // if current mode is mapping mode
  340 + compareMode = true;
  341 + mappingMode = false;
  342 + }
  343 + break;
  344 +
  345 + // zooming
  346 + case 'w': // if keyboard 'w' is pressed, then move closer
  347 + delta = zoom_factor;
  348 + cam.Push(delta);
  349 + break;
  350 + case 's': // if keyboard 's' is pressed, then leave farther
  351 + delta = -zoom_factor;
  352 + cam.Push(delta);
  353 + break;
  354 +
  355 + // resample and re-render the cylinder in different radius
  356 + case 'd': // if keyboard 'd' is pressed, then increase radius by radius_factor
  357 + radius += radius_factor;
  358 + break;
  359 + case 'a': // if keyboard 'a' is pressed, then decrease radius by radius_factor
  360 + radius -= radius_factor;
  361 + // get rid of the degenerated case when radius decrease below 0
  362 + if (radius < 0.001f)
  363 + radius = 0.2;
  364 + break;
  365 +
  366 + // close window and exit application
  367 + case 27: // if keyboard 'ESC' is pressed, then exit
  368 + exit(0);
348 369 }
349 370 glutPostRedisplay();
350 371 }
... ... @@ -396,13 +417,13 @@ void glut_initialize(){
396 417 glutMouseFunc(glut_mouse); //executed on a mouse click - sets starting mouse positions for rotations
397 418 glutMotionFunc(glut_motion); //executed when the mouse is moved while a button is pressed
398 419 if (ind == 1) { //only in mapping mode, keyboard will be used
399   - glutKeyboardFunc(glut_keyboard); //register keyboard callback
400 420 glutCreateMenu(glut_menu); //register menu option callback
401   - glutAddMenuEntry("Comparing Mode", 1);//register menu 1 as comparing mode
402   - glutAddMenuEntry("Mapping Mode", 2); //register menu 2 as mapping mode
403   - glutAddMenuEntry("Exit", 3); //register menu 3 as exiting
404   - glutAttachMenu(GLUT_RIGHT_BUTTON); //register right mouse to open menu option
405   - }
  421 + glutAddMenuEntry("Comparing Mode", 1); //register menu 1 as comparing mode
  422 + glutAddMenuEntry("Mapping Mode", 2); //register menu 2 as mapping mode
  423 + glutAddMenuEntry("Exit", 3); //register menu 3 as exiting
  424 + glutAttachMenu(GLUT_RIGHT_BUTTON); //register right mouse to open menu option
  425 + }
  426 + glutKeyboardFunc(glut_keyboard); //register keyboard callback
406 427 glutMouseWheelFunc(glut_wheel);
407 428  
408 429 texture_initialize(); //set up texture mapping (create texture maps, enable features)
... ... @@ -469,30 +490,6 @@ void map(float sigma, int device, float threshold){
469 490 std::cout << "FPR: " << FNR << std::endl;
470 491 }
471 492  
472   -void map_swc(float sigma, int device, float threshold) {
473   -
474   - // compare two networks
475   - _GT = GT.compare(T, sigma, device); //compare the ground truth to the test case - store errors in _GT
476   - _T = T.compare(_GT, sigma, device); //compare the test case to the ground truth - store errors in _T
477   -
478   - // mapping two networks and get their edge relation
479   - _GT.mapping(_T, _gt_t, device, threshold);
480   - _T.mapping(_GT, _t_gt, device, threshold);
481   -
482   - // generate random color set based on the number of edges in GT
483   - size_t num = _gt_t.size(); // also create random color for unmapping edge, but won't be used though
484   - colormap.resize(3 * num); // 3 portions compound RGB
485   - for (int i = 0; i < 3 * num; i++)
486   - colormap[i] = rand() / (float)RAND_MAX; // set to [0, 1]
487   -
488   - //calculate the metrics
489   - float FPR = _GT.average(0); //calculate the metrics
490   - float FNR = _T.average(0);
491   -
492   - std::cout << "FNR: " << FPR << std::endl; //print false alarms and misses
493   - std::cout << "FPR: " << FNR << std::endl;
494   -}
495   -
496 493 // writes features of the networks i.e average segment length, tortuosity, branching index, contraction, fractal dimension, number of end and branch points to a csv file
497 494 // Pranathi wrote this - saves network features to a CSV file
498 495 void features(std::string filename){
... ... @@ -563,10 +560,8 @@ int main(int argc, char* argv[])
563 560 if (args.nargs() >= 1) { // if at least one network file is specified
564 561 num_nets = 1; // set the number of networks to one
565 562 std::vector<std::string> tmp = stim::parser::split(args.arg(0), '.'); // split the filename at '.'
566   - if ("swc" == tmp[1]) { // loading swc file
  563 + if ("swc" == tmp[1]) // loading swc file
567 564 GT.load_swc(args.arg(0)); // load the specified file as the ground truth
568   - swc_ind = 1; // set the indicator of swc file to 1
569   - }
570 565 else if ("obj" == tmp[1]) // loading obj file
571 566 GT.load_obj(args.arg(0)); // load the specified file as the ground truth
572 567 else {
... ... @@ -575,43 +570,32 @@ int main(int argc, char* argv[])
575 570 }
576 571 }
577 572  
578   - if (args.nargs() == 2) { //if two files are specified, they will be displayed in neighboring viewports and compared
579   - if (1 == swc_ind) { //loading swc files
580   - int device = args["device"].as_int(); //get the device value from the user
581   - num_nets = 2; //set the number of networks to two
582   - //does it need to be resampled??
583   - sigma = args["sigma"].as_float(); //get the sigma value from the user
584   - T.load_swc(args.arg(1)); //load the second (test) network
585   - if (args["features"].is_set()) //if the user wants to save features
586   - features(args["features"].as_string());
587   - GT = GT.resample(resample_rate * sigma); //resample both networks based on the sigma value
588   - T = T.resample(resample_rate * sigma);
589   - if (args["mapping"].is_set()) {
590   - float threshold = args["mapping"].as_float();
591   - map_swc(sigma, device, threshold);
592   - //std::cout << "right now networks that are loaded from swc files do not need to be mapped with each other!" << std::endl;
593   - //exit(1);
594   - }
595   - else
596   - compare(sigma, device); //run the comparison algorithm
597   - }
  573 + if (args.nargs() == 2) { //if two files are specified, they will be displayed in neighboring viewports and compared
  574 + int device = args["device"].as_int(); //get the device value from the user
  575 + num_nets = 2; //set the number of networks to two
  576 + sigma = args["sigma"].as_float(); //get the sigma value from the user
  577 + radius = sigma;
  578 + std::vector<std::string> tmp = stim::parser::split(args.arg(1), '.'); // split the filename at '.'
  579 + if ("swc" == tmp[1]) //loading swc files
  580 + T.load_swc(args.arg(1)); //load the second (test) network
  581 + else if ("obj" == tmp[1]) //loading obj files
  582 + T.load_obj(args.arg(1));
598 583 else {
599   - int device = args["device"].as_int(); //get the device value from the user
600   - num_nets = 2; //set the number of networks to two
601   - sigma = args["sigma"].as_float(); //get the sigma value from the user
602   - T.load_obj(args.arg(1)); //load the second (test) network
603   - if (args["features"].is_set()) //if the user wants to save features
604   - features(args["features"].as_string());
605   - GT = GT.resample(resample_rate * sigma); //resample both networks based on the sigma value
606   - T = T.resample(resample_rate * sigma);
607   - if (args["mapping"].is_set()) {
608   - float threshold = args["mapping"].as_float();
609   - map(sigma, device, threshold);
610   - }
611   - else
612   - compare(sigma, device); //run the comparison algorithm
  584 + std::cout << "Invalid loading file" << std::endl;
  585 + exit(1);
  586 + }
  587 + if (args["features"].is_set()) //if the user wants to save features
  588 + features(args["features"].as_string());
  589 + //does it need to be resampled??
  590 + GT = GT.resample(resample_rate * sigma); //resample both networks based on the sigma value
  591 + T = T.resample(resample_rate * sigma);
  592 + if (args["mapping"].is_set()) {
  593 + float threshold = args["mapping"].as_float();
  594 + map(sigma, device, threshold);
  595 + }
  596 + else
  597 + compare(sigma, device); //run the comparison algorithm
613 598 }
614   - }
615 599  
616 600 //if a GUI is requested, display the network using OpenGL
617 601 if(args["gui"].is_set()){
... ...