Commit 887575e6 authored by karius's avatar karius

first working implementation of chebychev nodes

parent d83ae270
This diff is collapsed.
......@@ -57,41 +57,31 @@ public:
__host__ static void read_matrix_SO3_sampling(const char * file_name,float4 *& h_rotations, uint *& h_num_rotations);
__host__ static void shift_and_scale_envelope_score(float * d_score, float * d_normalized_score,float *d_min, float *d_max, size_t num_scores);
__host__ static void shift_and_scale_chamfer_score(float * d_score, float * d_normalized_score,float *d_min, float *d_max, size_t num_scores);
TransformationGrid(const std::pair<float3,float3> & bounding_box_outer,const std::pair<float3,float3> & bounding_box_inner,
float * const& translation_offset,const float4 & angle_box, const float & alpha, const float & tau_in);
//Testgrid constructor
TransformationGrid(const std::pair<float3,float3> & bounding_box_outer, // @suppress("Class members should be properly initialized")
const std::pair<float3,float3> & bounding_box_inner,float * const& translation_offset,float tau_in, uint n_testrots);
//empty grid with only one transformation (identity) for test purposes
TransformationGrid::TransformationGrid(const std::pair<float3,float3> & bounding_box_outer, // @suppress("Class members should be properly initialized")
const std::pair<float3,float3> & bounding_box_inner,float * const& translation_offset,float tau_in,const char * rotation_samples_file);
TransformationGrid(float * const& translation_offset);
TransformationGrid(const std::pair<float3,float3> bounding_box, const float tau, float * const& translation_offset);
virtual ~TransformationGrid();
__host__ void print_transformations(void);
bool * h_TEST;
bool * d_TEST;
uint n_testrot;
uint * h_num_translations;
uint * d_num_translations;
uint * h_num_rotations;
uint * d_num_rotations;
int3 * d_translation_index_dim;
int3 * h_translation_index_dim;
uint * h_translation_refinement;
uint * d_translation_refinement;
uint * h_translation_refinement_dim;
uint * d_translation_refinement_dim;
float3 * d_translation_total_shift;
float3 * h_translation_total_shift;
float3 * d_translation_step_shift;
float3 * h_translation_step_shift;
float4 * h_translation_offset;
float4 * d_translation_offset;
float4 * h_rotations;
float4 * d_rotations;
thrust::host_vector<float> th_scores;
thrust::device_vector<float> td_scores;
float * h_chebychev_nodes;
float * d_chebychev_nodes;
float * h_chebychev_weights;
float * d_chebychev_weights;
float * h_scores;
float * d_scores;
uint * h_num_transformations;
uint * d_num_transformations;
float tau;
__device__ int4 d_transformation_index_from_linear(const uint& linear_index);
__host__ int4 h_transformation_index_from_linear(const uint& linear_index);
__device__ void d_transformation_to_memory(const uint& i, float4 * const& rotation, float4 * const& translation, float4 * const& translation_offset);
......@@ -99,12 +89,13 @@ public:
__host__ void write_to_csv(const char * file_name);
__host__ int write_to_png(const char * filename);
__host__ void normalize_chamfer_score(void);
__host__ void refine_translation_grid(uint & i, float * h_translation);
__host__ void calculate_current_chebychev_nodes(void);
// static __host__ void evaluate(const char * file_name);
__host__ static thrust::device_vector<thrust::tuple<float4,float3,float>> n_best_from_file(const char * file_name, uint n);
__host__ static float3 transform(const float4& t);
void print(void);
private:
std::pair<float4,float4> bounding_box;
float alpha_actual;
std::pair<float3,float3> bounding_box;
};
#endif /* TRANSFORMATIONGRID_H_ */
No preview for this file type
......@@ -165,7 +165,7 @@ void envelope_score(Density binarized_map, float * d_workspaces, float * d_value
//transformation to shared
if(threadIdx.x == 0){
grid.d_transformation_to_memory(b,&rotation,&translation,grid.d_translation_offset);
// printf("Block %i handling, iteration %i transformation %f %f %f %f - %f %f %f\n",blockIdx.x,b,rotation.x,rotation.y,rotation.z,rotation.w,translation.x,translation.y,translation.z);
printf("Block %i handling, iteration %i transformation %f %f %f %f - %f %f %f\n",blockIdx.x,b,rotation.x,rotation.y,rotation.z,rotation.w,translation.x,translation.y,translation.z);
}
__syncthreads();
//imprint particle ids, iterate over particles
......@@ -533,94 +533,96 @@ void kai_score(float4 * d_coords, int * d_segmentations, uint * d_labels, float
// printf("\n\n");
//}
//
//BOOST_AUTO_TEST_CASE(envelope_test_simpler)
// {
// printf(" ============================= Testing Envelope score ==========================\n");
// printf(" ============================= ++++++++++++++++++++++ ==========================\n");
// //memory management
// nvmlInit_v2();
// //checkout free memory
// nvmlMemory_t memory_info;
// nvmlDevice_t device_handle;
// nvmlDeviceGetHandleByIndex_v2(0,&device_handle);
// nvmlDeviceGetMemoryInfo(device_handle, &memory_info);
//
// float above_value = -1.0;
// float below_value = 0.0;
// float coord_dim_outer[3] = {70,70,70};
// float coord_dim[3] = {50,50,50};
// float distance_threshold = 0.5;
// uint pixel_dim_outer[3] = {70,70,70};
// uint pixel_dim[3] = {50,50,50};
// uint3 cube_pixel_dim= {30,30,30};
// Density density;
// Density::from_bounds(density, &coord_dim_outer[0],&pixel_dim_outer[0]);
//// Density density = Density::from_mrc("test/tauA_Brf1_TBP_negstain_30A.mrc");
// density.init_device();
// density.to_device();
// density.mid_cube_with_surface(cube_pixel_dim,1,1.0,1.5);
//
// Density binarized_map;
// density.cut_and_binarize(binarized_map,1.0,above_value,below_value,{20,20,20});
//
// Density to_fit_density;
// Density::from_bounds(to_fit_density, &coord_dim[0],&pixel_dim[0]);
// to_fit_density.mid_cube_with_surface(cube_pixel_dim,1,1.0,1.5);
//
// //every pixel a particle, for test reasons
// thrust::device_vector<uint> td_on_surface_to_fit;
// size_t on_surface_num_to_fit;
// to_fit_density.indices_in_range(td_on_surface_to_fit,true,1.25,2,on_surface_num_to_fit);
//
// Particles pixel_coords = to_fit_density.particles_from_density(td_on_surface_to_fit);
// pixel_coords.origin_to_center();
//
// float memory_usage_factor = 0.5;
// size_t density_memory_usage = pixel_coords.particle_count()*sizeof(uint) + binarized_map.h_pixel_vol()*sizeof(float);
// unsigned long int max = 100;
// uint work_spaces = std::min(max,((uint)(((float)memory_info.free)*memory_usage_factor))/density_memory_usage);
//
// //set up TransformationGrid
// float tau_in = 1.0;
// std::pair<float3,float3> bounding_box_outer = density.get_bounding_box();
// std::pair<float3,float3> bounding_box_inner = pixel_coords.get_rotationally_safe_bounding_box();
// TransformationGrid transformation_grid(bounding_box_outer,bounding_box_inner,binarized_map.h_mid_point,tau_in,"rot_samples/N23_M5880_IcoC7.dat");
// CudaCheckError();
// printf("Transformation num: %u\n",transformation_grid.h_num_transformations[0]);
// printf("Translation num: %u\n",transformation_grid.h_num_translations[0]);
// printf("Translational shift: %f %f %f\n",transformation_grid.h_translation_step_shift->x,transformation_grid.h_translation_step_shift->y,transformation_grid.h_translation_step_shift->z);
// printf("Rotation num: %u\n",transformation_grid.h_num_rotations[0]);
//
//// TransformationGrid transformation_grid(binarized_map.h_mid_point);
// const size_t envelope_score_block_dim = 256;
// size_t envelope_score_grid_dim = work_spaces;
// size_t envelope_workspace_replaced_indices_volume = work_spaces*pixel_coords.particle_count()*sizeof(uint);
// size_t envelope_workspace_replaced_values_buffer_volume = work_spaces*pixel_coords.particle_count()*sizeof(float);
// //allocate workspace resources
//
// float * d_replaced_value_buffer;
// cudaMalloc((void **)&d_replaced_value_buffer,envelope_workspace_replaced_values_buffer_volume);
// uint * d_replaced_linear_indices;
// cudaMalloc((void **)&d_replaced_linear_indices,envelope_workspace_replaced_indices_volume);
// CudaCheckError();
// float * d_workspaces;
// cudaMalloc((void **)&d_workspaces,work_spaces*binarized_map.h_pixel_vol()*sizeof(float));
// printf("Allocating %u workspaces for %u pixels each occupying %u bytes of memory ...\n",work_spaces,binarized_map.h_pixel_vol(),envelope_workspace_replaced_indices_volume + envelope_workspace_replaced_values_buffer_volume);
// CudaCheckError();
// //(Density binarized_map, float * d_workspaces, uint * d_linear_offsets, float4 * d_coords,size_t num_particles , TransformationGrid grid)
// envelope_score<envelope_score_block_dim><<<envelope_score_grid_dim,envelope_score_block_dim>>>
// (binarized_map, d_workspaces,d_replaced_value_buffer, d_replaced_linear_indices,pixel_coords.d_data,pixel_coords.particle_count(),transformation_grid);
//// transformation_grid.write_to_csv("test/env_score.csv");
// float * d_normalized_scores;
// cudaMalloc((void **)&d_normalized_scores,transformation_grid.h_num_transformations[0]*sizeof(*d_normalized_scores));
// float * d_min;
// cudaMalloc((void **)&d_min,sizeof(*d_min));
// float * d_max;
// cudaMalloc((void **)&d_max,sizeof(*d_max));
// transformation_grid.shift_and_scale_envelope_score(transformation_grid.d_scores,d_normalized_scores,d_min,d_max,transformation_grid.h_num_transformations[0]);
// transformation_grid.write_to_csv("test/env_score_5880.csv");
// printf("Finished testing envelope score ... \n\n");
//}
BOOST_AUTO_TEST_CASE(envelope_test_simpler)
{
printf(" ============================= Testing Envelope score ==========================\n");
printf(" ============================= ++++++++++++++++++++++ ==========================\n");
//memory management
nvmlInit_v2();
//checkout free memory
nvmlMemory_t memory_info;
nvmlDevice_t device_handle;
nvmlDeviceGetHandleByIndex_v2(0,&device_handle);
nvmlDeviceGetMemoryInfo(device_handle, &memory_info);
float above_value = -1.0;
float below_value = 0.0;
float coord_dim_outer[3] = {70,70,70};
float coord_dim[3] = {50,50,50};
float distance_threshold = 0.5;
uint pixel_dim_outer[3] = {70,70,70};
uint pixel_dim[3] = {50,50,50};
uint3 cube_pixel_dim= {30,30,30};
Density density;
Density::from_bounds(density, &coord_dim_outer[0],&pixel_dim_outer[0]);
// Density density = Density::from_mrc("test/tauA_Brf1_TBP_negstain_30A.mrc");
density.init_device();
density.to_device();
density.mid_cube_with_surface(cube_pixel_dim,1,1.0,1.5);
Density binarized_map;
density.cut_and_binarize(binarized_map,1.0,above_value,below_value,{20,20,20});
Density to_fit_density;
Density::from_bounds(to_fit_density, &coord_dim[0],&pixel_dim[0]);
to_fit_density.mid_cube_with_surface(cube_pixel_dim,1,1.0,1.5);
//every pixel a particle, for test reasons
thrust::device_vector<uint> td_on_surface_to_fit;
size_t on_surface_num_to_fit;
to_fit_density.indices_in_range(td_on_surface_to_fit,true,1.25,2,on_surface_num_to_fit);
Particles pixel_coords = to_fit_density.particles_from_density(td_on_surface_to_fit);
pixel_coords.origin_to_center();
float memory_usage_factor = 0.5;
size_t density_memory_usage = pixel_coords.particle_count()*sizeof(uint) + binarized_map.h_pixel_vol()*sizeof(float);
unsigned long int max = 100;
uint work_spaces = std::min(max,((uint)(((float)memory_info.free)*memory_usage_factor))/density_memory_usage);
//set up TransformationGrid
float tau_in = 1.0;
std::pair<float3,float3> bounding_box_outer = density.get_bounding_box();
std::pair<float3,float3> bounding_box_inner = pixel_coords.get_rotationally_safe_bounding_box();
std::pair<float3,float3> effective_translational_shift;
effective_translational_shift.first = bounding_box_outer.first - bounding_box_inner.first;
effective_translational_shift.second = bounding_box_outer.second - bounding_box_inner.second;
TransformationGrid transformation_grid(effective_translational_shift,tau_in,binarized_map.h_mid_point);
CudaCheckError();
printf("Transformation num: %u\n",transformation_grid.h_num_transformations[0]);
printf("Translation num: %u\n",transformation_grid.h_num_translations[0]);
printf("Rotation num: %u\n",transformation_grid.h_num_rotations[0]);
// TransformationGrid transformation_grid(binarized_map.h_mid_point);
const size_t envelope_score_block_dim = 256;
size_t envelope_score_grid_dim = work_spaces;
size_t envelope_workspace_replaced_indices_volume = work_spaces*pixel_coords.particle_count()*sizeof(uint);
size_t envelope_workspace_replaced_values_buffer_volume = work_spaces*pixel_coords.particle_count()*sizeof(float);
//allocate workspace resources
float * d_replaced_value_buffer;
cudaMalloc((void **)&d_replaced_value_buffer,envelope_workspace_replaced_values_buffer_volume);
uint * d_replaced_linear_indices;
cudaMalloc((void **)&d_replaced_linear_indices,envelope_workspace_replaced_indices_volume);
CudaCheckError();
float * d_workspaces;
cudaMalloc((void **)&d_workspaces,work_spaces*binarized_map.h_pixel_vol()*sizeof(float));
printf("Allocating %u workspaces for %u pixels each occupying %u bytes of memory ...\n",work_spaces,binarized_map.h_pixel_vol(),envelope_workspace_replaced_indices_volume + envelope_workspace_replaced_values_buffer_volume);
CudaCheckError();
//(Density binarized_map, float * d_workspaces, uint * d_linear_offsets, float4 * d_coords,size_t num_particles , TransformationGrid grid)
envelope_score<envelope_score_block_dim><<<envelope_score_grid_dim,envelope_score_block_dim>>>
(binarized_map, d_workspaces,d_replaced_value_buffer, d_replaced_linear_indices,pixel_coords.d_data,pixel_coords.particle_count(),transformation_grid);
// transformation_grid.write_to_csv("test/env_score.csv");
float * d_normalized_scores;
cudaMalloc((void **)&d_normalized_scores,transformation_grid.h_num_transformations[0]*sizeof(*d_normalized_scores));
float * d_min;
cudaMalloc((void **)&d_min,sizeof(*d_min));
float * d_max;
cudaMalloc((void **)&d_max,sizeof(*d_max));
transformation_grid.shift_and_scale_envelope_score(transformation_grid.d_scores,d_normalized_scores,d_min,d_max,transformation_grid.h_num_transformations[0]);
transformation_grid.write_to_csv("test/env_score_5880.csv");
printf("Finished testing envelope score ... \n\n");
}
//
//BOOST_AUTO_TEST_CASE(chamfer_test_simpler)
//{
......@@ -928,39 +930,39 @@ void kai_score(float4 * d_coords, int * d_segmentations, uint * d_labels, float
// }
//}
BOOST_AUTO_TEST_CASE(chamfer_test)
{
float coord_dim[3] = {150,150,150};
uint pixel_dim[3] = {100,100,100};
Density density;
printf("Here\n");
Density::from_bounds(density,&coord_dim[0],&pixel_dim[0]);
printf("Here2\n");
Particles particles = Particles::from_pdb("test/tauA.pdb",true);
printf("Here3\n");
Density::from_particles_on_grid(density,particles,5.0);
density.to_mrc("test/tauA.mrc");
// thrust::device_vector<uint> td_on_surface;
// Density difference_map = sub_density.difference_map_density_threshold(threshold-deviation,threshold + deviation,td_on_surface);
// Particles pixel_coords = structure_density.particles_from_density(td_on_surface);
// pixel_coords.center_on_density(sub_density.h_coord_dim);
//// TransformationGrid::TransformationGrid(const std::pair<float3,float3> & bounding_box_outer,
//// const std::pair<float3,float3> & bounding_box_inner,const float4 & angle_box, const float & alpha, const float & tau_in)
// TransformationGrid transformation_grid(sub_density.get_bounding_box(),particles.get_bounding_box(),{0,0,0,0},20,2.5);
// printf("Transformation num: %u\n",transformation_grid.h_num_transformations[0]);
// printf("Translation num: %u\n",transformation_grid.h_num_translations[0]);
// printf("Rotation num: %u\n",transformation_grid.h_num_rotations[0]);
//// void chamfer_score(float4 * d_coords, float * d_difference_map, uint * d_pixel_dim,
//// float pixel_size, size_t num_particles,TransformationGrid grid)
// size_t chamfer_kernel_block_dim = 128;
// size_t chamfer_kernel_grid_dim = 32;
// size_t chamfer_kernel_shared_mem = sizeof(float4) + sizeof(float3);
// chamfer_score<<<chamfer_kernel_grid_dim,chamfer_kernel_block_dim,chamfer_kernel_shared_mem>>>
// (pixel_coords.d_data,difference_map.d_data,difference_map.d_pixel_dim,difference_map.pixel_size,
// pixel_coords.particle_count(),transformation_grid);
// printf("");
}
//BOOST_AUTO_TEST_CASE(chamfer_test)
//{
// float coord_dim[3] = {150,150,150};
// uint pixel_dim[3] = {100,100,100};
// Density density;
// printf("Here\n");
// Density::from_bounds(density,&coord_dim[0],&pixel_dim[0]);
// printf("Here2\n");
// Particles particles = Particles::from_pdb("test/tauA.pdb",true);
// printf("Here3\n");
//
// Density::from_particles_on_grid(density,particles,5.0);
// density.to_mrc("test/tauA.mrc");
//// thrust::device_vector<uint> td_on_surface;
//// Density difference_map = sub_density.difference_map_density_threshold(threshold-deviation,threshold + deviation,td_on_surface);
//// Particles pixel_coords = structure_density.particles_from_density(td_on_surface);
//// pixel_coords.center_on_density(sub_density.h_coord_dim);
////// TransformationGrid::TransformationGrid(const std::pair<float3,float3> & bounding_box_outer,
////// const std::pair<float3,float3> & bounding_box_inner,const float4 & angle_box, const float & alpha, const float & tau_in)
//// TransformationGrid transformation_grid(sub_density.get_bounding_box(),particles.get_bounding_box(),{0,0,0,0},20,2.5);
//// printf("Transformation num: %u\n",transformation_grid.h_num_transformations[0]);
//// printf("Translation num: %u\n",transformation_grid.h_num_translations[0]);
//// printf("Rotation num: %u\n",transformation_grid.h_num_rotations[0]);
////// void chamfer_score(float4 * d_coords, float * d_difference_map, uint * d_pixel_dim,
////// float pixel_size, size_t num_particles,TransformationGrid grid)
//// size_t chamfer_kernel_block_dim = 128;
//// size_t chamfer_kernel_grid_dim = 32;
//// size_t chamfer_kernel_shared_mem = sizeof(float4) + sizeof(float3);
//// chamfer_score<<<chamfer_kernel_grid_dim,chamfer_kernel_block_dim,chamfer_kernel_shared_mem>>>
//// (pixel_coords.d_data,difference_map.d_data,difference_map.d_pixel_dim,difference_map.pixel_size,
//// pixel_coords.particle_count(),transformation_grid);
//// printf("");
//}
//BOOST_AUTO_TEST_CASE(difference_map_test)
//{
......
Markdown is supported
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