Commit f4f7db17 authored by Sabrina Templeton's avatar Sabrina Templeton
Browse files

code cleanup in geometry shader

parent ae14fe01
......@@ -32,7 +32,7 @@ float random (vec2 uv) {
void
get_color( float u , out vec3 color ) {
u = random( vec2(v_Index,v_Index) ) * (u_nb_points+1);
//u = random( vec2(v_Index,v_Index) ) * (u_nb_points+1);
float umin = 0.;//u_umin;
float umax = u_nb_points +1;//u_umax;
......
// the version should be specified by a macro at compile time
// e.g. #version 410 or #version 330
//#version 410
#extension GL_ARB_gpu_shader_fp64 : enable
layout (points) in;
......@@ -17,26 +17,18 @@ const int MAX_VERTS = 30;
flat in int[] instance_ID;
out float v_Index;
//out vec3 v_Ind;
layout (triangle_strip, max_vertices = 40) out; // this was 40, upped it for debugging purposes
layout (triangle_strip, max_vertices = 40) out;
//uniform usamplerBuffer connectivity;
uniform samplerBuffer seeds;
uniform usamplerBuffer nn; //ask about why this has to be u I don't really understand the difference here but this is what made it work!
uniform usamplerBuffer nn;
// render full polygon
void render_poly(vec3[MAX_VERTS] poly, vec3 center, int space){
vec3 p0 = poly[0];
//vec3 p1;
for (int i = 0; i < space -2; i++){
gl_Position = u_ModelViewProjectionMatrix*vec4(p0, 1.0);
v_Index = float(instance_ID[0]);
v_Index = float(instance_ID[0]);
EmitVertex();
gl_Position = u_ModelViewProjectionMatrix*vec4(poly[i+1], 1.0);
......@@ -48,12 +40,15 @@ void render_poly(vec3[MAX_VERTS] poly, vec3 center, int space){
v_Index = float(instance_ID[0]);
EmitVertex();
EndPrimitive();
}
}
vec3 calc_intersect(vec3 p1, vec3 p2, vec3 q1, vec3 q2){ // p1 and p2 are seeds, q1 and q2 are borders
/*
* Originally from geogram (https://github.com/BrunoLevy/geogram)
* in src/lib/geogram/voronoi/generic_RVD_vertex.h
*/
//initialize variables
float d = 0;
float l1 = 0;
......@@ -94,17 +89,15 @@ vec3 calc_intersect(vec3 p1, vec3 p2, vec3 q1, vec3 q2){ // p1 and p2 are seeds,
}
void calc_poly(vec3 zi, vec3 zj, inout vec3[MAX_VERTS] curr_poly, inout int cspace){
//cspace, which is inout as it goes along with curr_poly and represents the actual length of the polygon information in curr_poly at any given moment (ideally)
//cspace is an argument, which is inout as it goes along with curr_poly and represents the actual length of the polygon information in curr_poly at any given moment (ideally)
int nspace = 0; // nspace is 'new space,' the variable which acts as a pointer to keep track of the next available spot to add a vertex to new_poly
//where zi is the pt and zj is the current nearest neighbor
//first step : start with just the initial polygon for now and perform one clip
vec3 new_poly[MAX_VERTS];//= vec3[4](vec3(0, 0, 0), vec3(0, 0, 0), vec3(0, 0, 0), vec3(0, 0, 0));
//where zi is the current site pt and zj is the current nearest neighbor
vec3 new_poly[MAX_VERTS];
for (int i = 0; i < cspace; i++) {
vec3 b0;
vec3 b1;// initialize b0 and b1
vec3 b1;
if (i == cspace-1){
b0 = curr_poly[cspace -1];
b1 = curr_poly[0];
......@@ -113,16 +106,16 @@ void calc_poly(vec3 zi, vec3 zj, inout vec3[MAX_VERTS] curr_poly, inout int cspa
b0 = curr_poly[i];
b1 = curr_poly[i+1];
}
int side1 = 0; // initialize side 1 and side 2
int side1 = 0;
int side2 = 0;
//calculate side 1
//calculate side 1: determine if the first border point is inside or outside
if (distance(b0, zi) < distance(b0, zj)) {
side1 = 1;
}
else {
side1 = -1;
}
//calculate side 2
//calculate side 2: determine if the second border point is inside or outside
if (distance(b1, zi) < distance(b1, zj)) {
side2 = 1;
}
......@@ -130,7 +123,6 @@ void calc_poly(vec3 zi, vec3 zj, inout vec3[MAX_VERTS] curr_poly, inout int cspa
side2 = -1;
}
vec3 intersect;
//again there's more to do here but I'm going to start by just calculating the intersection pt in the most hardcoded way possible
if (side1 != side2){
intersect = calc_intersect(zi, zj, b0, b1); //b0 and b1 are the current border points
if (side1 == 1){
......@@ -143,81 +135,48 @@ void calc_poly(vec3 zi, vec3 zj, inout vec3[MAX_VERTS] curr_poly, inout int cspa
nspace += 1;
}
}
else if (side1 == 1){ //meaning both sides are 1
new_poly[nspace] = b0; // append pt 1 essentially
else if (side1 == 1){ // meaning both sides are inside
new_poly[nspace] = b0; // we add in the first border point since we will add in the next point at the next step of the loop
nspace += 1;
}
//otherwise both sides are -1 and we take neither
// otherwise both sides are -1 and we take neither
}
// note to self: this is definitely outside of the loop ( and should remain to be)
curr_poly = new_poly; // does this need to be a deep copy? I don't think so.
curr_poly = new_poly;
cspace = nspace; // setting the current space to be the new space so that it is updated for the new polygon.
}
void
main() {
// lookup the seed points, taking in glPrimitiveIDIn
int site_idx = instance_ID[0];
//v_Index = instance_ID[0];
//v_Ind = vec3(instance_ID[0], instance_ID[0], instance_ID[0])
vec3 p0 = vec3(texelFetch( seeds , site_idx ).xy, 0);
//vec3 p1 = vec3(texelFetch( seeds , 0 ).xy, 0);
//if (site_idx == 0){ // this is for debugging purposes: REMOVE ONCE WE WANT TO EXPAND TO ACTUALLY COMPUTE FOR MORE THAN ONE SEED
// but also when we do that we will need to fix the nearest neighbors which right now are hardcoded in
//this is the initial poly
vec3 curr_poly1 = vec3(0, 0, 0);
vec3 curr_poly2 = vec3(1, 0, 0);
vec3 curr_poly3 = vec3(1, 1, 0);
vec3 curr_poly4 = vec3(0, 1, 0);
vec3 poly[MAX_VERTS];
// Lookup the seed points, taking in instance_ID which we output from the vertex shader. This should be equivalent to glPrimitiveID_In but is not
int site_idx = instance_ID[0];
vec3 zi = vec3(texelFetch( seeds , site_idx ).xy, 0);
// This is the initial polygon, we start by initializing it to the polygon which covers the total area of the diagram
vec3[MAX_VERTS] curr_poly;
// = vec3[MAX_VERTS](curr_poly1, curr_poly2, curr_poly3, curr_poly4);
curr_poly[0] = curr_poly1;
curr_poly[1] = curr_poly2;
curr_poly[2] = curr_poly3;
curr_poly[3] = curr_poly4;
int space = 4; //this is cspace
// this will obviously change once we are actually fetching from the nearest neighbor instead of the seed
float d = 0;
curr_poly[0] = vec3(0, 0, 0);
curr_poly[1] = vec3(1, 0, 0);
curr_poly[2] = vec3(1, 1, 0);
curr_poly[3] = vec3(0, 1, 0);
int space = 4;
for (int i = 1; i < u_nb_neighbors; i ++){
//works when i < 2 but not otherwise
//vec3 nn0 = vec3(texelFetch( seeds, i).xy, 0);
int nn0 = int(texelFetch( nn, i + (site_idx * u_nb_neighbors)).x); // i + (site_idx * 7)
//int nn0 = neighbor_hardcoded[i + (site_idx * 7)];//int(texelFetch( nn, i + (site_idx * 7)).x); // i + (site_idx * 7) // remember that this is just the index, we have to use this to lookup in seeds
vec3 nearest_neighbor = vec3(texelFetch( seeds, nn0).xy, 0);
float d1 = d;
d = distance(nearest_neighbor, p0);
vec3 zj= vec3(texelFetch( seeds, nn0).xy, 0); // nn0 is just the index so we need to look up the corresponding point in seeds
calc_poly(p0, nearest_neighbor, curr_poly, space);
calc_poly(zi, zj, curr_poly, space);
//calculate radius of security
// Below is the radius of security theorem
float radius = 0;
vec3 furthest_pt;
// loop over the vertices of the polygon to find the furthest one, the distance from zi to that point becomes the radius.
for (int j = 0; j < space; j++ ) {
vec3 curr_pt = curr_poly[j];
//render_square(curr_pt, .01);
if(distance(curr_pt, p0) > radius) {
radius = distance(curr_pt, p0);
furthest_pt = curr_pt;
if(distance(curr_pt, zi) > radius) {
radius = distance(curr_pt, zi);
}
}
if (distance(nearest_neighbor, p0) > (2.1 * radius)){
//render_square(vec3(0, 0, 0), .1);
// If the distance between the current nearest neighbor and the current site is greater than twice the radius, stop clipping
if (distance(zj, zi) > (2.1 * radius)){
break;
}
}
render_poly(curr_poly, p0, space);
render_poly(curr_poly, zi, space);
}
......@@ -35,7 +35,7 @@ UT_TEST_SUITE( voronoi_toy )
// set up the view
int width = 500, height = width;
graphics::vec3 eye = {0.5,0.5,2.0};
graphics::vec3 eye = {0.5,0.5,1.5};
graphics::vec3 up = {0,1,0};
graphics::vec3 center = {0.5,0.5,0.0};
......@@ -121,13 +121,18 @@ UT_TEST_CASE( test1 )
graphics::gl_index vertex_array;
GL_CALL( glGenVertexArrays( 1, &vertex_array ) );
GL_CALL( glBindVertexArray(vertex_array) );
// *************************************************************************
//generate random seeds desired times 3
int num_points = 1000;
int num_points = 1e6;
int size = num_points*3;
int DEBUG = 1;
int DEBUG = 0;
// *************************************************************************
std::vector<graphics::gl_float> seeds;
......@@ -168,76 +173,39 @@ else{
}
//cout << seeds;
//*/ //uncomment this to restore to random
// std::vector<graphics::gl_float> seeds = { // test case 1
// 0.25, 0.5 , 0.0,
// 0.25, 0.3 , 0.0,
// 0.15, 0.7 , 0.0,
// 0.55, 0.5 , 0.0,
// 0.01, 0.25 , 0.0,
// 0.8, 0.8 , 0.0,
// 0.9, 0.9 , 0.0,
// };
// std::vector<graphics::gl_float> seeds = { // test case 2
// 0.5, 0.5 , 0.0,
// 0.8, 0.8 , 0.0,
// 0.2, 0.8 , 0.0,
// 0.8, 0.2 , 0.0,
// 0.2, 0.2 , 0.0,
// 0.0, 0.0 , 0.0,
// 0.0, 0.0 , 0.0,
// };
std::vector<real_t> seeds_d( seeds.begin() , seeds.end() ); //convert the seeds to double
//compute the nearest neighbors here
// start the nearest neighbor timer around here,
clock_t TIME0 = clock();
// GEO::NearestNeighborSearch* nns = GEO::NearestNeighborSearch::create(3,"BNN"); // number of dimensions is really 2, but we have been using three
GEO::NearestNeighborSearch* nns = GEO::NearestNeighborSearch::create(3,"BNN"); // number of dimensions is really 2, but we have been using three
// index_t nb_points = 7;
// std::vector<real_t> x( nb_points*3 ); //3 even tho its 2d, equivalent to the seeds
// //std::vector<graphics::gl_float> seeds(nb_points*3); // the actual seeds
// for (index_t k=0;k<nb_points*3;k++)
// //float y = random_within(0.,1.);
// x[k] = random_within(0.,1.); //we will substitute seeds for this code chunk
// //x[k] = random_within(0.,1.);
nns->set_points( num_points , seeds_d.data() );
// nns->set_points( num_points , seeds_d.data() );
index_t nb_neighbors = 30; // should be 50 or 30
std::vector<index_t> neighbors(nb_neighbors);
std::vector<index_t> total_neighbors(nb_neighbors * num_points);
std::vector<double> neighbors_sq_dist(nb_neighbors);
// std::vector<index_t> neighbors(nb_neighbors);
// std::vector<index_t> total_neighbors(nb_neighbors * num_points);
// std::vector<double> neighbors_sq_dist(nb_neighbors);
index_t i = 0;
for (index_t k=0;k<num_points;k++)
{
nns->get_nearest_neighbors( nb_neighbors , k , neighbors.data() , neighbors_sq_dist.data() );
//UT_ASSERT_EQUALS( neighbors[0] , k );
//print_inline(neighbors);
// index_t i = 0;
// for (index_t k=0;k<num_points;k++)
// {
// nns->get_nearest_neighbors( nb_neighbors , k , neighbors.data() , neighbors_sq_dist.data() );
//here we need to resave the neighbors in a new full neighbors array before moving on
// //here we need to resave the neighbors in a new full neighbors array before moving on
for (index_t j = 0; j < nb_neighbors; j++){
// for (index_t j = 0; j < nb_neighbors; j++){
//total_neighbors[((k+1)*i)-1] = neighbors[i];
total_neighbors[i] = neighbors[j];
//std::cout << i << std::endl;
i++;
// //total_neighbors[((k+1)*i)-1] = neighbors[i];
// total_neighbors[i] = neighbors[j];
// //std::cout << i << std::endl;
// i++;
}
// }
}
// }
......@@ -249,28 +217,6 @@ else{
int nb_thread = ProcessCPU::maximum_concurrent_threads();
printf("--> nearest neighbor calculation time = %g seconds\n",real_t(TIME1-TIME0)/real_t(CLOCKS_PER_SEC)/nb_thread);
// std::vector<graphics::gl_float> seeds = {
// 0.0, 0.0 , 0.0,
// 1.0, 0.0 , 0.0,
// 0.0, 1.0 , 0.0,
// 0.0, 0.0 , 0.0,
// 1.0, 1.0 , 0.0,
// 0.0, -1.0, 0.0,
// };
// std::vector<GLuint> neighbors_x = {
// 0, 1, 2, 3, 4, 5, 6,
// 1, 0, 4, 3, 2, 5, 6,
// 2, 0, 1, 3, 4, 5, 6,
// 3, 0, 1, 2, 4, 5, 6,
// 4, 1, 0, 2, 3, 5, 6,
// 5, 6, 3, 2, 0, 1, 4,
// 6, 5, 3, 2, 0, 1, 4
// };
// convert nearest neighbor type
//std::vector<GLuint> neighbors_f( total_neighbors.begin() , total_neighbors.end() );
......@@ -299,13 +245,12 @@ else{
//get a colormaps and bind it to the colormap buffer
avro::Colormap colormap;
colormap.change_style("viridis");
colormap.change_style("jet");
index_t ncolor = 256*3;
GL_CALL( glBindBuffer( GL_TEXTURE_BUFFER , colormap_buffer) );
GL_CALL( glBufferData( GL_TEXTURE_BUFFER , sizeof(float) * ncolor , colormap.data() , GL_STATIC_DRAW) );
// generate a texture to hold the seeds
graphics::gl_index seeds_texture;
GL_CALL( glGenTextures( 1 , &seeds_texture) );
......@@ -313,9 +258,6 @@ else{
GL_CALL( glBindTexture( GL_TEXTURE_BUFFER , seeds_texture) );
GL_CALL( glTexBuffer( GL_TEXTURE_BUFFER , GL_RGB32F , points_buffer ) );
// generate a texture to hold the colormap
GLuint colormap_texture;
GL_CALL( glGenTextures( 1 , &colormap_texture) );
......@@ -331,10 +273,6 @@ else{
GL_CALL( glBindTexture( GL_TEXTURE_BUFFER, nn_texture) );
GL_CALL( glTexBuffer( GL_TEXTURE_BUFFER, GL_R32I, nn_buffer) );
// ensure we can capture the escape key being pressed
glfwSetInputMode(window, GLFW_STICKY_KEYS, GL_TRUE);
glfwPollEvents();
......@@ -373,30 +311,30 @@ else{
glUniform1i(nn_location, 2);
// draw
clock_t TIME2 = clock();
GL_CALL( glBindBuffer(GL_ARRAY_BUFFER, points_buffer ) ); // it doesn't really matter which buffer we bind here since it isn't used
// int points_remaining = num_points;
// int first_point_to_draw = 0;
// int call_size = 1e5;
// while (points_remaining > 0){
// //std::cout << "one iteration" << std::endl;
// std::cout << points_remaining << std::endl;
// // std::cout << first_point_to_draw << std::endl;
// if (points_remaining < call_size){
// // std::cout << "if statement" << std::endl;
// // std::cout << first_point_to_draw << std::endl;
// // std::cout << first_point_to_draw + points_remaining << std::endl;
// GL_CALL( glDrawArrays(GL_POINTS, first_point_to_draw , points_remaining ) );
// points_remaining = 0;
// }
// else{
// //std::cout << "test" << std::endl;
// GL_CALL( glDrawArrays(GL_POINTS, first_point_to_draw , call_size ) );
// points_remaining -= call_size;
// first_point_to_draw += call_size;
// }
/*
int points_remaining = num_points;
int first_point_to_draw = 0;
int call_size = 1e5;
while (points_remaining > 0){
//std::cout << "one iteration" << std::endl;
std::cout << points_remaining << std::endl;
// std::cout << first_point_to_draw << std::endl;
if (points_remaining < call_size){
// std::cout << "if statement" << std::endl;
// std::cout << first_point_to_draw << std::endl;
// std::cout << first_point_to_draw + points_remaining << std::endl;
GL_CALL( glDrawArrays(GL_POINTS, first_point_to_draw , points_remaining ) );
points_remaining = 0;
}
else{
//std::cout << "test" << std::endl;
GL_CALL( glDrawArrays(GL_POINTS, first_point_to_draw , call_size ) );
points_remaining -= call_size;
first_point_to_draw += call_size;
}
*/
GLuint query;
//GLuint64 elapsed_time;
......@@ -416,7 +354,7 @@ else{
//GL_CALL( glDrawArrays (GL_POINTS, 0, num_points) );
clock_t TIME3 = clock();
//std::cout << time;
FILE * pFile;
......@@ -427,24 +365,15 @@ else{
//fprintf(pFile, "--> nearest neighbor calculation time = %g seconds\n",real_t(TIME1-TIME0)/real_t(CLOCKS_PER_SEC));
//fprintf(pFile, "--> GPU run time = %f seconds\n", time/1e9 );
fprintf(pFile, "1e6 & %.2g & %.2f \\\\ \n", real_t(TIME1-TIME0)/real_t(CLOCKS_PER_SEC), time/1e9 );
//fprintf(pFile, "1e6 & %.2E & %.2E \\\\ \n", real_t((TIME1-TIME0)/real_t(CLOCKS_PER_SEC)\nb_thread, time/1e9 );
//std::cout << nb_thread;
std::cout << "should have written to a file.";
//std::cout << "should have written to a file.";
myfile.close();
//points_shader.use();
//points_shader.setUniform("u_ModelViewProjectionMatrix", mvp );
glfwSwapBuffers(window);
glfwSwapBuffers(window);
while (true) {
// swap buffers and wait for user input
//glfwSwapBuffers(window);
glfwPollEvents();
// determine if we should exit the render loop
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment