#[compute]
#version 450

// Workgroup size
layout(local_size_x = 8, local_size_y = 8, local_size_z = 1) in;

// Storage buffer
layout(set = 0, binding = 0, std430) buffer DataBuffer {
    float sample_points[];
} voxels;

layout(set = 0, binding = 1) uniform Params {
    int world_size;
    float threshold;
} params;

layout(set = 0, binding = 2, std430) buffer SurfaceBuffer {
    float surface_points[]; 
} surface;

uint index3(uint x, uint y, uint z) {
    return x + y * params.world_size + z * params.world_size * params.world_size;
}

void store_surface_point(uint idx, vec3 pos) {
    uint base = idx * 3u;
    surface.surface_points[base + 0u] = pos.x;
    surface.surface_points[base + 1u] = pos.y;
    surface.surface_points[base + 2u] = pos.z;
}

void main() {
    uvec3 id = gl_GlobalInvocationID;
    if (id.x >= uint(params.world_size) ||
        id.y >= uint(params.world_size) ||
        id.z >= uint(params.world_size))
        return;

    const ivec3 SURFACE_AXIS[8] = ivec3[](
	      ivec3(1,0,0),
	      ivec3(0,1,0),
	      ivec3(0,0,1),
	      ivec3(1,0,1),
	      ivec3(0,1,1),
	      ivec3(1,1,0),
	      ivec3(1,1,1),
	      ivec3(0,0,0)
      );

    vec3 p = vec3(id) - params.world_size / 2 + 1;
    float d = length(p) - params.threshold;
    uint idx = index3(id.x, id.y, id.z);
    voxels.sample_points[idx] = d;

    int previous_sign = 0;
    vec3 previous_sample_point_coords = vec3(id.x, id.y, id.z);
    float previous_sample_point_distance = 0.0;
    vec3 intersection_points_sum = vec3(0.0, 0.0, 0.0);
    uint intersection_points_count = 0u;
    for (int sample_point_to_check_index = 0; sample_point_to_check_index < 8; sample_point_to_check_index++){
      ivec3 sample_point_to_check = ivec3(id) + SURFACE_AXIS[sample_point_to_check_index];
      
      // bounds check
      if (sample_point_to_check.x < 0 || sample_point_to_check.y < 0 || sample_point_to_check.z < 0) continue;
      if (uint(sample_point_to_check.x) >= params.world_size || uint(sample_point_to_check.y) >= params.world_size || uint(sample_point_to_check.z) >= params.world_size) continue;
      
      uint buffer_index = index3(uint(sample_point_to_check.x), uint(sample_point_to_check.y), uint(sample_point_to_check.z));
      float sample_point_distance_from_sdf = voxels.sample_points[buffer_index];
      
      int current_sign = (sample_point_distance_from_sdf >= 0.0 ? 1 : -1);
      
      if (previous_sign != 0 && current_sign != previous_sign) {
        float t = previous_sample_point_distance / (previous_sample_point_distance - sample_point_distance_from_sdf);
        vec3 intersect = mix(vec3(previous_sample_point_coords), vec3(sample_point_to_check), t);
        intersection_points_sum += intersect;
        intersection_points_count++;
      }

      previous_sign = current_sign;
      previous_sample_point_distance = sample_point_distance_from_sdf;
      previous_sample_point_coords = sample_point_to_check;
    } 

    if (intersection_points_count > 0u) {
        vec3 avg = intersection_points_sum / float(intersection_points_count);
        store_surface_point(idx, avg);
    } else {
        store_surface_point(idx, vec3(0.0)); // no surface
    }



}