diff options
-rw-r--r-- | source/blender/blenkernel/BKE_mesh_sample.hh | 59 | ||||
-rw-r--r-- | source/blender/blenkernel/intern/mesh_sample.cc | 172 | ||||
-rw-r--r-- | source/blender/editors/sculpt_paint/curves_sculpt_add.cc | 520 | ||||
-rw-r--r-- | source/blender/editors/sculpt_paint/curves_sculpt_brush.cc | 29 | ||||
-rw-r--r-- | source/blender/editors/sculpt_paint/curves_sculpt_intern.hh | 8 | ||||
-rw-r--r-- | source/blender/geometry/CMakeLists.txt | 2 | ||||
-rw-r--r-- | source/blender/geometry/GEO_add_curves_on_mesh.hh | 54 | ||||
-rw-r--r-- | source/blender/geometry/intern/add_curves_on_mesh.cc | 354 |
8 files changed, 701 insertions, 497 deletions
diff --git a/source/blender/blenkernel/BKE_mesh_sample.hh b/source/blender/blenkernel/BKE_mesh_sample.hh index 37c2ad7d3cb..f47201e89cc 100644 --- a/source/blender/blenkernel/BKE_mesh_sample.hh +++ b/source/blender/blenkernel/BKE_mesh_sample.hh @@ -6,12 +6,20 @@ * \ingroup bke */ +#include "BLI_function_ref.hh" #include "BLI_generic_virtual_array.hh" #include "BLI_math_vec_types.hh" +#include "DNA_meshdata_types.h" + #include "BKE_attribute.h" struct Mesh; +struct BVHTreeFromMesh; + +namespace blender { +struct RandomNumberGenerator; +} namespace blender::bke { struct ReadAttributeLookup; @@ -82,4 +90,55 @@ class MeshAttributeInterpolator { Span<float3> ensure_nearest_weights(); }; +/** + * Find randomly distributed points on the surface of a mesh within a 3D sphere. This does not + * sample an exact number of points because it comes with extra overhead to avoid bias that is only + * required in some cases. If an exact number of points is required, that has to be implemented at + * a higher level. + * + * \param approximate_density: Roughly the number of points per unit of area. + * \return The number of added points. + */ +int sample_surface_points_spherical(RandomNumberGenerator &rng, + const Mesh &mesh, + Span<int> looptri_indices_to_sample, + const float3 &sample_pos, + float sample_radius, + float approximate_density, + Vector<float3> &r_bary_coords, + Vector<int> &r_looptri_indices, + Vector<float3> &r_positions); + +/** + * Find randomly distributed points on the surface of a mesh within a circle that is projected on + * the mesh. This does not result in an exact number of points because that would come with extra + * overhead and is not always possible. If an exact number of points is required, that has to be + * implemented at a higher level. + * + * \param region_position_to_ray: Function that converts a 2D position into a 3D ray that is used + * to find positions on the mesh. + * \param mesh_bvhtree: BVH tree of the triangles in the mesh. Passed in so that it does not have + * to be retrieved again. + * \param tries_num: Number of 2d positions that are sampled. The maximum + * number of new samples. + * \return The number of added points. + */ +int sample_surface_points_projected( + RandomNumberGenerator &rng, + const Mesh &mesh, + BVHTreeFromMesh &mesh_bvhtree, + const float2 &sample_pos_re, + float sample_radius_re, + FunctionRef<void(const float2 &pos_re, float3 &r_start, float3 &r_end)> region_position_to_ray, + bool front_face_only, + int tries_num, + int max_points, + Vector<float3> &r_bary_coords, + Vector<int> &r_looptri_indices, + Vector<float3> &r_positions); + +float3 compute_bary_coord_in_triangle(const Mesh &mesh, + const MLoopTri &looptri, + const float3 &position); + } // namespace blender::bke::mesh_surface_sample diff --git a/source/blender/blenkernel/intern/mesh_sample.cc b/source/blender/blenkernel/intern/mesh_sample.cc index 7595c08a208..106c4c610ba 100644 --- a/source/blender/blenkernel/intern/mesh_sample.cc +++ b/source/blender/blenkernel/intern/mesh_sample.cc @@ -2,12 +2,15 @@ #include "BKE_attribute_access.hh" #include "BKE_attribute_math.hh" +#include "BKE_bvhutils.h" #include "BKE_mesh_runtime.h" #include "BKE_mesh_sample.hh" #include "DNA_mesh_types.h" #include "DNA_meshdata_types.h" +#include "BLI_rand.hh" + namespace blender::bke::mesh_surface_sample { template<typename T> @@ -259,4 +262,173 @@ void MeshAttributeInterpolator::sample_attribute(const ReadAttributeLookup &src_ } } +int sample_surface_points_spherical(RandomNumberGenerator &rng, + const Mesh &mesh, + const Span<int> looptri_indices_to_sample, + const float3 &sample_pos, + const float sample_radius, + const float approximate_density, + Vector<float3> &r_bary_coords, + Vector<int> &r_looptri_indices, + Vector<float3> &r_positions) +{ + const Span<MLoopTri> looptris{BKE_mesh_runtime_looptri_ensure(&mesh), + BKE_mesh_runtime_looptri_len(&mesh)}; + + const float sample_radius_sq = pow2f(sample_radius); + const float sample_plane_area = M_PI * sample_radius_sq; + /* Used for switching between two triangle sampling strategies. */ + const float area_threshold = sample_plane_area; + + const int old_num = r_bary_coords.size(); + + for (const int looptri_index : looptri_indices_to_sample) { + const MLoopTri &looptri = looptris[looptri_index]; + + const float3 &v0 = mesh.mvert[mesh.mloop[looptri.tri[0]].v].co; + const float3 &v1 = mesh.mvert[mesh.mloop[looptri.tri[1]].v].co; + const float3 &v2 = mesh.mvert[mesh.mloop[looptri.tri[2]].v].co; + + const float looptri_area = area_tri_v3(v0, v1, v2); + + if (looptri_area < area_threshold) { + /* The triangle is small compared to the sample radius. Sample by generating random + * barycentric coordinates. */ + const int amount = rng.round_probabilistic(approximate_density * looptri_area); + for ([[maybe_unused]] const int i : IndexRange(amount)) { + const float3 bary_coord = rng.get_barycentric_coordinates(); + const float3 point_pos = attribute_math::mix3(bary_coord, v0, v1, v2); + const float dist_to_sample_sq = math::distance_squared(point_pos, sample_pos); + if (dist_to_sample_sq > sample_radius_sq) { + continue; + } + + r_bary_coords.append(bary_coord); + r_looptri_indices.append(looptri_index); + r_positions.append(point_pos); + } + } + else { + /* The triangle is large compared to the sample radius. Sample by generating random points + * on the triangle plane within the sample radius. */ + float3 normal; + normal_tri_v3(normal, v0, v1, v2); + + float3 sample_pos_proj = sample_pos; + project_v3_plane(sample_pos_proj, normal, v0); + + const float proj_distance_sq = math::distance_squared(sample_pos_proj, sample_pos); + const float sample_radius_factor_sq = 1.0f - + std::min(1.0f, proj_distance_sq / sample_radius_sq); + const float radius_proj_sq = sample_radius_sq * sample_radius_factor_sq; + const float radius_proj = std::sqrt(radius_proj_sq); + const float circle_area = M_PI * radius_proj; + + const int amount = rng.round_probabilistic(approximate_density * circle_area); + + const float3 axis_1 = math::normalize(v1 - v0) * radius_proj; + const float3 axis_2 = math::normalize(math::cross(axis_1, math::cross(axis_1, v2 - v0))) * + radius_proj; + + for ([[maybe_unused]] const int i : IndexRange(amount)) { + const float r = std::sqrt(rng.get_float()); + const float angle = rng.get_float() * 2.0f * M_PI; + const float x = r * std::cos(angle); + const float y = r * std::sin(angle); + const float3 point_pos = sample_pos_proj + axis_1 * x + axis_2 * y; + if (!isect_point_tri_prism_v3(point_pos, v0, v1, v2)) { + /* Sampled point is not in the triangle. */ + continue; + } + + float3 bary_coord; + interp_weights_tri_v3(bary_coord, v0, v1, v2, point_pos); + + r_bary_coords.append(bary_coord); + r_looptri_indices.append(looptri_index); + r_positions.append(point_pos); + } + } + } + return r_bary_coords.size() - old_num; +} + +int sample_surface_points_projected( + RandomNumberGenerator &rng, + const Mesh &mesh, + BVHTreeFromMesh &mesh_bvhtree, + const float2 &sample_pos_re, + const float sample_radius_re, + const FunctionRef<void(const float2 &pos_re, float3 &r_start, float3 &r_end)> + region_position_to_ray, + const bool front_face_only, + const int tries_num, + const int max_points, + Vector<float3> &r_bary_coords, + Vector<int> &r_looptri_indices, + Vector<float3> &r_positions) +{ + const Span<MLoopTri> looptris{BKE_mesh_runtime_looptri_ensure(&mesh), + BKE_mesh_runtime_looptri_len(&mesh)}; + + int point_count = 0; + for ([[maybe_unused]] const int _ : IndexRange(tries_num)) { + if (point_count == max_points) { + break; + } + + const float r = sample_radius_re * std::sqrt(rng.get_float()); + const float angle = rng.get_float() * 2.0f * M_PI; + float3 ray_start, ray_end; + const float2 pos_re = sample_pos_re + r * float2(std::cos(angle), std::sin(angle)); + region_position_to_ray(pos_re, ray_start, ray_end); + const float3 ray_direction = math::normalize(ray_end - ray_start); + + BVHTreeRayHit ray_hit; + ray_hit.dist = FLT_MAX; + ray_hit.index = -1; + BLI_bvhtree_ray_cast(mesh_bvhtree.tree, + ray_start, + ray_direction, + 0.0f, + &ray_hit, + mesh_bvhtree.raycast_callback, + &mesh_bvhtree); + + if (ray_hit.index == -1) { + continue; + } + + if (front_face_only) { + const float3 normal = ray_hit.no; + if (math::dot(ray_direction, normal) >= 0.0f) { + continue; + } + } + + const int looptri_index = ray_hit.index; + const float3 pos = ray_hit.co; + + const float3 bary_coords = compute_bary_coord_in_triangle(mesh, looptris[looptri_index], pos); + + r_positions.append(pos); + r_bary_coords.append(bary_coords); + r_looptri_indices.append(looptri_index); + point_count++; + } + return point_count; +} + +float3 compute_bary_coord_in_triangle(const Mesh &mesh, + const MLoopTri &looptri, + const float3 &position) +{ + const float3 &v0 = mesh.mvert[mesh.mloop[looptri.tri[0]].v].co; + const float3 &v1 = mesh.mvert[mesh.mloop[looptri.tri[1]].v].co; + const float3 &v2 = mesh.mvert[mesh.mloop[looptri.tri[2]].v].co; + float3 bary_coords; + interp_weights_tri_v3(bary_coords, v0, v1, v2, position); + return bary_coords; +} + } // namespace blender::bke::mesh_surface_sample diff --git a/source/blender/editors/sculpt_paint/curves_sculpt_add.cc b/source/blender/editors/sculpt_paint/curves_sculpt_add.cc index f013fa05f4c..c475b90a0e7 100644 --- a/source/blender/editors/sculpt_paint/curves_sculpt_add.cc +++ b/source/blender/editors/sculpt_paint/curves_sculpt_add.cc @@ -22,9 +22,9 @@ #include "BKE_geometry_set.hh" #include "BKE_mesh.h" #include "BKE_mesh_runtime.h" +#include "BKE_mesh_sample.hh" #include "BKE_paint.h" #include "BKE_report.h" -#include "BKE_spline.hh" #include "DNA_brush_enums.h" #include "DNA_brush_types.h" @@ -38,10 +38,12 @@ #include "ED_screen.h" #include "ED_view3d.h" +#include "GEO_add_curves_on_mesh.hh" + #include "WM_api.h" /** - * The code below uses a prefix naming convention to indicate the coordinate space: + * The code below uses a suffix naming convention to indicate the coordinate space: * cu: Local space of the curves object that is being edited. * su: Local space of the surface object. * wo: World space. @@ -70,16 +72,6 @@ class AddOperation : public CurvesSculptStrokeOperation { void on_stroke_extended(const bContext &C, const StrokeExtension &stroke_extension) override; }; -static void initialize_straight_curve_positions(const float3 &p1, - const float3 &p2, - MutableSpan<float3> r_positions) -{ - const float step = 1.0f / (float)(r_positions.size() - 1); - for (const int i : r_positions.index_range()) { - r_positions[i] = math::interpolate(p1, p2, i * step); - } -} - /** * Utility class that actually executes the update when the stroke is updated. That's useful * because it avoids passing a very large number of parameters between functions. @@ -241,28 +233,28 @@ struct AddOperationExecutor { return; } - Array<NeighborsVector> neighbors_per_curve; if (use_interpolation_) { this->ensure_curve_roots_kdtree(); - neighbors_per_curve = this->find_curve_neighbors(added_points); - } - - /* Resize to add the new curves, building the offsets in the array owned by the curves. */ - const int tot_added_curves = added_points.bary_coords.size(); - curves_->resize(curves_->points_num(), curves_->curves_num() + tot_added_curves); - if (interpolate_point_count_) { - this->initialize_curve_offsets_with_interpolation(neighbors_per_curve); - } - else { - this->initialize_curve_offsets_without_interpolation(constant_points_per_curve_); } - /* Resize to add the correct point count calculated as part of building the offsets. */ - curves_->resize(curves_->offsets().last(), curves_->curves_num()); - - this->initialize_attributes(added_points, neighbors_per_curve); - - curves_->update_curve_types(); + geometry::AddCurvesOnMeshInputs add_inputs; + add_inputs.root_positions_cu = added_points.positions_cu; + add_inputs.bary_coords = added_points.bary_coords; + add_inputs.looptri_indices = added_points.looptri_indices; + add_inputs.interpolate_length = interpolate_length_; + add_inputs.interpolate_shape = interpolate_shape_; + add_inputs.interpolate_point_count = interpolate_point_count_; + add_inputs.fallback_curve_length = new_curve_length_; + add_inputs.fallback_point_count = constant_points_per_curve_; + add_inputs.surface = surface_; + add_inputs.surface_bvh = &surface_bvh_; + add_inputs.surface_looptris = surface_looptris_; + add_inputs.surface_uv_map = surface_uv_map_; + add_inputs.corner_normals_su = corner_normals_su_; + add_inputs.curves_to_surface_mat = curves_to_surface_mat_; + add_inputs.surface_to_curves_normal_mat = surface_to_curves_normal_mat_; + add_inputs.old_roots_kdtree = self_->curve_roots_kdtree_; + geometry::add_curves_on_mesh(*curves_, add_inputs); DEG_id_tag_update(&curves_id_->id, ID_RECALC_GEOMETRY); WM_main_add_notifier(NC_GEOM | ND_DATA, &curves_id_->id); @@ -312,7 +304,7 @@ struct AddOperationExecutor { const int looptri_index = ray_hit.index; const float3 brush_pos_su = ray_hit.co; - const float3 bary_coords = compute_bary_coord_in_triangle( + const float3 bary_coords = bke::mesh_surface_sample::compute_bary_coord_in_triangle( *surface_, surface_looptris_[looptri_index], brush_pos_su); const float3 brush_pos_cu = surface_to_curves_mat_ * brush_pos_su; @@ -339,60 +331,37 @@ struct AddOperationExecutor { const float4x4 &brush_transform) { const int old_amount = r_added_points.bary_coords.size(); - const int max_iterations = std::max(100'000, add_amount_ * 10); + const int max_iterations = 100; int current_iteration = 0; while (r_added_points.bary_coords.size() < old_amount + add_amount_) { if (current_iteration++ >= max_iterations) { break; } - - const float r = brush_radius_re_ * std::sqrt(rng.get_float()); - const float angle = rng.get_float() * 2.0f * M_PI; - const float2 pos_re = brush_pos_re_ + r * float2(std::cos(angle), std::sin(angle)); - - float3 ray_start_wo, ray_end_wo; - ED_view3d_win_to_segment_clipped( - ctx_.depsgraph, ctx_.region, ctx_.v3d, pos_re, ray_start_wo, ray_end_wo, true); - const float3 ray_start_cu = brush_transform * (world_to_curves_mat_ * ray_start_wo); - const float3 ray_end_cu = brush_transform * (world_to_curves_mat_ * ray_end_wo); - - const float3 ray_start_su = curves_to_surface_mat_ * ray_start_cu; - const float3 ray_end_su = curves_to_surface_mat_ * ray_end_cu; - const float3 ray_direction_su = math::normalize(ray_end_su - ray_start_su); - - BVHTreeRayHit ray_hit; - ray_hit.dist = FLT_MAX; - ray_hit.index = -1; - BLI_bvhtree_ray_cast(surface_bvh_.tree, - ray_start_su, - ray_direction_su, - 0.0f, - &ray_hit, - surface_bvh_.raycast_callback, - &surface_bvh_); - - if (ray_hit.index == -1) { - continue; - } - - if (use_front_face_) { - const float3 normal_su = ray_hit.no; - if (math::dot(ray_direction_su, normal_su) >= 0.0f) { - continue; - } + const int missing_amount = add_amount_ + old_amount - r_added_points.bary_coords.size(); + const int new_points = bke::mesh_surface_sample::sample_surface_points_projected( + rng, + *surface_, + surface_bvh_, + brush_pos_re_, + brush_radius_re_, + [&](const float2 &pos_re, float3 &r_start_su, float3 &r_end_su) { + float3 start_wo, end_wo; + ED_view3d_win_to_segment_clipped( + ctx_.depsgraph, ctx_.region, ctx_.v3d, pos_re, start_wo, end_wo, true); + const float3 start_cu = brush_transform * (world_to_curves_mat_ * start_wo); + const float3 end_cu = brush_transform * (world_to_curves_mat_ * end_wo); + r_start_su = curves_to_surface_mat_ * start_cu; + r_end_su = curves_to_surface_mat_ * end_cu; + }, + use_front_face_, + add_amount_, + missing_amount, + r_added_points.bary_coords, + r_added_points.looptri_indices, + r_added_points.positions_cu); + for (float3 &pos : r_added_points.positions_cu.as_mutable_span().take_back(new_points)) { + pos = surface_to_curves_mat_ * pos; } - - const int looptri_index = ray_hit.index; - const float3 pos_su = ray_hit.co; - - const float3 bary_coords = compute_bary_coord_in_triangle( - *surface_, surface_looptris_[looptri_index], pos_su); - - const float3 pos_cu = surface_to_curves_mat_ * pos_su; - - r_added_points.positions_cu.append(pos_cu); - r_added_points.bary_coords.append(bary_coords); - r_added_points.looptri_indices.append(looptri_index); } } @@ -505,9 +474,6 @@ struct AddOperationExecutor { const float brush_plane_area_su = M_PI * brush_radius_sq_su; const float approximate_density_su = add_amount_ / brush_plane_area_su; - /* Used for switching between two triangle sampling strategies. */ - const float area_threshold = brush_plane_area_su; - /* Usually one or two iterations should be enough. */ const int max_iterations = 5; int current_iteration = 0; @@ -517,78 +483,18 @@ struct AddOperationExecutor { if (current_iteration++ >= max_iterations) { break; } - - for (const int looptri_index : looptri_indices) { - const MLoopTri &looptri = surface_looptris_[looptri_index]; - - const float3 v0_su = surface_->mvert[surface_->mloop[looptri.tri[0]].v].co; - const float3 v1_su = surface_->mvert[surface_->mloop[looptri.tri[1]].v].co; - const float3 v2_su = surface_->mvert[surface_->mloop[looptri.tri[2]].v].co; - - const float looptri_area_su = area_tri_v3(v0_su, v1_su, v2_su); - - if (looptri_area_su < area_threshold) { - /* The triangle is small compared to the brush radius. Sample by generating random - * barycentric coordinates. */ - const int amount = rng.round_probabilistic(approximate_density_su * looptri_area_su); - for ([[maybe_unused]] const int i : IndexRange(amount)) { - const float3 bary_coord = rng.get_barycentric_coordinates(); - const float3 point_pos_su = attribute_math::mix3(bary_coord, v0_su, v1_su, v2_su); - const float distance_to_brush_sq_su = math::distance_squared(point_pos_su, - brush_pos_su); - if (distance_to_brush_sq_su > brush_radius_sq_su) { - continue; - } - - r_added_points.bary_coords.append(bary_coord); - r_added_points.looptri_indices.append(looptri_index); - r_added_points.positions_cu.append(surface_to_curves_mat_ * point_pos_su); - } - } - else { - /* The triangle is large compared to the brush radius. Sample by generating random points - * on the triangle plane within the brush radius. */ - float3 normal_su; - normal_tri_v3(normal_su, v0_su, v1_su, v2_su); - - float3 brush_pos_proj_su = brush_pos_su; - project_v3_plane(brush_pos_proj_su, normal_su, v0_su); - - const float proj_distance_sq_su = math::distance_squared(brush_pos_proj_su, - brush_pos_su); - const float brush_radius_factor_sq = 1.0f - - std::min(1.0f, - proj_distance_sq_su / brush_radius_sq_su); - const float radius_proj_sq_su = brush_radius_sq_su * brush_radius_factor_sq; - const float radius_proj_su = std::sqrt(radius_proj_sq_su); - const float circle_area_su = M_PI * radius_proj_su; - - const int amount = rng.round_probabilistic(approximate_density_su * circle_area_su); - - const float3 axis_1_su = math::normalize(v1_su - v0_su) * radius_proj_su; - const float3 axis_2_su = math::normalize(math::cross( - axis_1_su, math::cross(axis_1_su, v2_su - v0_su))) * - radius_proj_su; - - for ([[maybe_unused]] const int i : IndexRange(amount)) { - const float r = std::sqrt(rng.get_float()); - const float angle = rng.get_float() * 2.0f * M_PI; - const float x = r * std::cos(angle); - const float y = r * std::sin(angle); - const float3 point_pos_su = brush_pos_proj_su + axis_1_su * x + axis_2_su * y; - if (!isect_point_tri_prism_v3(point_pos_su, v0_su, v1_su, v2_su)) { - /* Sampled point is not in the triangle. */ - continue; - } - - float3 bary_coord; - interp_weights_tri_v3(bary_coord, v0_su, v1_su, v2_su, point_pos_su); - - r_added_points.bary_coords.append(bary_coord); - r_added_points.looptri_indices.append(looptri_index); - r_added_points.positions_cu.append(surface_to_curves_mat_ * point_pos_su); - } - } + const int new_points = bke::mesh_surface_sample::sample_surface_points_spherical( + rng, + *surface_, + looptri_indices, + brush_pos_su, + brush_radius_su, + approximate_density_su, + r_added_points.bary_coords, + r_added_points.looptri_indices, + r_added_points.positions_cu); + for (float3 &pos : r_added_points.positions_cu.as_mutable_span().take_back(new_points)) { + pos = surface_to_curves_mat_ * pos; } } @@ -613,312 +519,6 @@ struct AddOperationExecutor { BLI_kdtree_3d_balance(self_->curve_roots_kdtree_); } } - - void initialize_curve_offsets_with_interpolation(const Span<NeighborsVector> neighbors_per_curve) - { - MutableSpan<int> new_offsets = curves_->offsets_for_write().drop_front(tot_old_curves_); - - attribute_math::DefaultMixer<int> mixer{new_offsets}; - threading::parallel_for(neighbors_per_curve.index_range(), 1024, [&](IndexRange curves_range) { - for (const int i : curves_range) { - if (neighbors_per_curve[i].is_empty()) { - mixer.mix_in(i, constant_points_per_curve_, 1.0f); - } - else { - for (const NeighborInfo &neighbor : neighbors_per_curve[i]) { - const int neighbor_points_num = curves_->points_for_curve(neighbor.index).size(); - mixer.mix_in(i, neighbor_points_num, neighbor.weight); - } - } - } - }); - mixer.finalize(); - - bke::curves::accumulate_counts_to_offsets(new_offsets, tot_old_points_); - } - - void initialize_curve_offsets_without_interpolation(const int points_per_curve) - { - MutableSpan<int> new_offsets = curves_->offsets_for_write().drop_front(tot_old_curves_); - int offset = tot_old_points_; - for (const int i : new_offsets.index_range()) { - new_offsets[i] = offset; - offset += points_per_curve; - } - } - - void initialize_attributes(const AddedPoints &added_points, - const Span<NeighborsVector> neighbors_per_curve) - { - Array<float> new_lengths_cu(added_points.bary_coords.size()); - if (interpolate_length_) { - this->interpolate_lengths(neighbors_per_curve, new_lengths_cu); - } - else { - new_lengths_cu.fill(new_curve_length_); - } - - Array<float3> new_normals_su = this->compute_normals_for_added_curves_su(added_points); - if (!surface_uv_map_.is_empty()) { - this->initialize_surface_attachment(added_points); - } - - this->fill_new_selection(); - - if (interpolate_shape_) { - this->initialize_position_with_interpolation( - added_points, neighbors_per_curve, new_normals_su, new_lengths_cu); - } - else { - this->initialize_position_without_interpolation( - added_points, new_lengths_cu, new_normals_su); - } - } - - /** - * Select newly created points or curves in new curves if necessary. - */ - void fill_new_selection() - { - switch (curves_id_->selection_domain) { - case ATTR_DOMAIN_CURVE: { - const VArray<float> selection = curves_->selection_curve_float(); - if (selection.is_single() && selection.get_internal_single() >= 1.0f) { - return; - } - curves_->selection_curve_float_for_write().drop_front(tot_old_curves_).fill(1.0f); - break; - } - case ATTR_DOMAIN_POINT: { - const VArray<float> selection = curves_->selection_point_float(); - if (selection.is_single() && selection.get_internal_single() >= 1.0f) { - return; - } - curves_->selection_point_float_for_write().drop_front(tot_old_points_).fill(1.0f); - break; - } - default: - BLI_assert_unreachable(); - } - } - - Array<NeighborsVector> find_curve_neighbors(const AddedPoints &added_points) - { - const int tot_added_curves = added_points.bary_coords.size(); - Array<NeighborsVector> neighbors_per_curve(tot_added_curves); - threading::parallel_for(IndexRange(tot_added_curves), 128, [&](const IndexRange range) { - for (const int i : range) { - const float3 root_cu = added_points.positions_cu[i]; - std::array<KDTreeNearest_3d, max_neighbors> nearest_n; - const int found_neighbors = BLI_kdtree_3d_find_nearest_n( - self_->curve_roots_kdtree_, root_cu, nearest_n.data(), max_neighbors); - float tot_weight = 0.0f; - for (const int neighbor_i : IndexRange(found_neighbors)) { - KDTreeNearest_3d &nearest = nearest_n[neighbor_i]; - const float weight = 1.0f / std::max(nearest.dist, 0.00001f); - tot_weight += weight; - neighbors_per_curve[i].append({nearest.index, weight}); - } - /* Normalize weights. */ - for (NeighborInfo &neighbor : neighbors_per_curve[i]) { - neighbor.weight /= tot_weight; - } - } - }); - return neighbors_per_curve; - } - - void interpolate_lengths(const Span<NeighborsVector> neighbors_per_curve, - MutableSpan<float> r_lengths) - { - const Span<float3> positions_cu = curves_->positions(); - - threading::parallel_for(r_lengths.index_range(), 128, [&](const IndexRange range) { - for (const int added_curve_i : range) { - const Span<NeighborInfo> neighbors = neighbors_per_curve[added_curve_i]; - float length_sum = 0.0f; - for (const NeighborInfo &neighbor : neighbors) { - const IndexRange neighbor_points = curves_->points_for_curve(neighbor.index); - float neighbor_length = 0.0f; - for (const int segment_i : neighbor_points.drop_back(1)) { - const float3 &p1 = positions_cu[segment_i]; - const float3 &p2 = positions_cu[segment_i + 1]; - neighbor_length += math::distance(p1, p2); - } - length_sum += neighbor.weight * neighbor_length; - } - const float length = neighbors.is_empty() ? new_curve_length_ : length_sum; - r_lengths[added_curve_i] = length; - } - }); - } - - float3 compute_point_normal_su(const int looptri_index, const float3 &bary_coord) - { - const MLoopTri &looptri = surface_looptris_[looptri_index]; - const int l0 = looptri.tri[0]; - const int l1 = looptri.tri[1]; - const int l2 = looptri.tri[2]; - - const float3 &l0_normal_su = corner_normals_su_[l0]; - const float3 &l1_normal_su = corner_normals_su_[l1]; - const float3 &l2_normal_su = corner_normals_su_[l2]; - - const float3 normal_su = math::normalize( - attribute_math::mix3(bary_coord, l0_normal_su, l1_normal_su, l2_normal_su)); - return normal_su; - } - - Array<float3> compute_normals_for_added_curves_su(const AddedPoints &added_points) - { - Array<float3> normals_su(added_points.bary_coords.size()); - threading::parallel_for(normals_su.index_range(), 256, [&](const IndexRange range) { - for (const int i : range) { - const int looptri_index = added_points.looptri_indices[i]; - const float3 &bary_coord = added_points.bary_coords[i]; - normals_su[i] = compute_surface_point_normal( - surface_looptris_[looptri_index], bary_coord, corner_normals_su_); - } - }); - return normals_su; - } - - void initialize_surface_attachment(const AddedPoints &added_points) - { - MutableSpan<float2> surface_uv_coords = curves_->surface_uv_coords_for_write(); - threading::parallel_for( - added_points.bary_coords.index_range(), 1024, [&](const IndexRange range) { - for (const int i : range) { - const int curve_i = tot_old_curves_ + i; - const MLoopTri &looptri = surface_looptris_[added_points.looptri_indices[i]]; - const float2 &uv0 = surface_uv_map_[looptri.tri[0]]; - const float2 &uv1 = surface_uv_map_[looptri.tri[1]]; - const float2 &uv2 = surface_uv_map_[looptri.tri[2]]; - const float3 &bary_coords = added_points.bary_coords[i]; - const float2 uv = attribute_math::mix3(bary_coords, uv0, uv1, uv2); - surface_uv_coords[curve_i] = uv; - } - }); - } - - /** - * Initialize new curves so that they are just a straight line in the normal direction. - */ - void initialize_position_without_interpolation(const AddedPoints &added_points, - const Span<float> lengths_cu, - const MutableSpan<float3> normals_su) - { - MutableSpan<float3> positions_cu = curves_->positions_for_write(); - - threading::parallel_for( - added_points.bary_coords.index_range(), 256, [&](const IndexRange range) { - for (const int i : range) { - const IndexRange points = curves_->points_for_curve(tot_old_curves_ + i); - const float3 &root_cu = added_points.positions_cu[i]; - const float length = lengths_cu[i]; - const float3 &normal_su = normals_su[i]; - const float3 normal_cu = math::normalize(surface_to_curves_normal_mat_ * normal_su); - const float3 tip_cu = root_cu + length * normal_cu; - - initialize_straight_curve_positions(root_cu, tip_cu, positions_cu.slice(points)); - } - }); - } - - /** - * Use neighboring curves to determine the shape. - */ - void initialize_position_with_interpolation(const AddedPoints &added_points, - const Span<NeighborsVector> neighbors_per_curve, - const Span<float3> new_normals_su, - const Span<float> new_lengths_cu) - { - MutableSpan<float3> positions_cu = curves_->positions_for_write(); - - threading::parallel_for( - added_points.bary_coords.index_range(), 256, [&](const IndexRange range) { - for (const int i : range) { - const Span<NeighborInfo> neighbors = neighbors_per_curve[i]; - const IndexRange points = curves_->points_for_curve(tot_old_curves_ + i); - - const float length_cu = new_lengths_cu[i]; - const float3 &normal_su = new_normals_su[i]; - const float3 normal_cu = math::normalize(surface_to_curves_normal_mat_ * normal_su); - - const float3 &root_cu = added_points.positions_cu[i]; - - if (neighbors.is_empty()) { - /* If there are no neighbors, just make a straight line. */ - const float3 tip_cu = root_cu + length_cu * normal_cu; - initialize_straight_curve_positions(root_cu, tip_cu, positions_cu.slice(points)); - continue; - } - - positions_cu.slice(points).fill(root_cu); - - for (const NeighborInfo &neighbor : neighbors) { - const int neighbor_curve_i = neighbor.index; - const float3 &neighbor_first_pos_cu = - positions_cu[curves_->offsets()[neighbor_curve_i]]; - const float3 neighbor_first_pos_su = curves_to_surface_mat_ * neighbor_first_pos_cu; - - BVHTreeNearest nearest; - nearest.dist_sq = FLT_MAX; - BLI_bvhtree_find_nearest(surface_bvh_.tree, - neighbor_first_pos_su, - &nearest, - surface_bvh_.nearest_callback, - &surface_bvh_); - const int neighbor_looptri_index = nearest.index; - const MLoopTri &neighbor_looptri = surface_looptris_[neighbor_looptri_index]; - - const float3 neighbor_bary_coord = compute_bary_coord_in_triangle( - *surface_, neighbor_looptri, nearest.co); - - const float3 neighbor_normal_su = compute_surface_point_normal( - surface_looptris_[neighbor_looptri_index], - neighbor_bary_coord, - corner_normals_su_); - const float3 neighbor_normal_cu = math::normalize(surface_to_curves_normal_mat_ * - neighbor_normal_su); - - /* The rotation matrix used to transform relative coordinates of the neighbor curve - * to the new curve. */ - float normal_rotation_cu[3][3]; - rotation_between_vecs_to_mat3(normal_rotation_cu, neighbor_normal_cu, normal_cu); - - const IndexRange neighbor_points = curves_->points_for_curve(neighbor_curve_i); - const float3 &neighbor_root_cu = positions_cu[neighbor_points[0]]; - - /* Use a temporary #PolySpline, because that's the easiest way to resample an - * existing curve right now. Resampling is necessary if the length of the new curve - * does not match the length of the neighbors or the number of handle points is - * different. */ - PolySpline neighbor_spline; - neighbor_spline.resize(neighbor_points.size()); - neighbor_spline.positions().copy_from(positions_cu.slice(neighbor_points)); - neighbor_spline.mark_cache_invalid(); - - const float neighbor_length_cu = neighbor_spline.length(); - const float length_factor = std::min(1.0f, length_cu / neighbor_length_cu); - - const float resample_factor = (1.0f / (points.size() - 1.0f)) * length_factor; - for (const int j : IndexRange(points.size())) { - const Spline::LookupResult lookup = neighbor_spline.lookup_evaluated_factor( - j * resample_factor); - const float index_factor = lookup.evaluated_index + lookup.factor; - float3 p; - neighbor_spline.sample_with_index_factors<float3>( - neighbor_spline.positions(), {&index_factor, 1}, {&p, 1}); - const float3 relative_coord = p - neighbor_root_cu; - float3 rotated_relative_coord = relative_coord; - mul_m3_v3(normal_rotation_cu, rotated_relative_coord); - positions_cu[points[j]] += neighbor.weight * rotated_relative_coord; - } - } - } - }); - } }; void AddOperation::on_stroke_extended(const bContext &C, const StrokeExtension &stroke_extension) diff --git a/source/blender/editors/sculpt_paint/curves_sculpt_brush.cc b/source/blender/editors/sculpt_paint/curves_sculpt_brush.cc index 7e583773512..2dd0da9d21a 100644 --- a/source/blender/editors/sculpt_paint/curves_sculpt_brush.cc +++ b/source/blender/editors/sculpt_paint/curves_sculpt_brush.cc @@ -324,33 +324,4 @@ CurvesSculptCommonContext::CurvesSculptCommonContext(const bContext &C) this->rv3d = CTX_wm_region_view3d(&C); } -float3 compute_surface_point_normal(const MLoopTri &looptri, - const float3 &bary_coord, - const Span<float3> corner_normals) -{ - const int l0 = looptri.tri[0]; - const int l1 = looptri.tri[1]; - const int l2 = looptri.tri[2]; - - const float3 &l0_normal = corner_normals[l0]; - const float3 &l1_normal = corner_normals[l1]; - const float3 &l2_normal = corner_normals[l2]; - - const float3 normal = math::normalize( - attribute_math::mix3(bary_coord, l0_normal, l1_normal, l2_normal)); - return normal; -} - -float3 compute_bary_coord_in_triangle(const Mesh &mesh, - const MLoopTri &looptri, - const float3 &position) -{ - const float3 &v0 = mesh.mvert[mesh.mloop[looptri.tri[0]].v].co; - const float3 &v1 = mesh.mvert[mesh.mloop[looptri.tri[1]].v].co; - const float3 &v2 = mesh.mvert[mesh.mloop[looptri.tri[2]].v].co; - float3 bary_coords; - interp_weights_tri_v3(bary_coords, v0, v1, v2, position); - return bary_coords; -} - } // namespace blender::ed::sculpt_paint diff --git a/source/blender/editors/sculpt_paint/curves_sculpt_intern.hh b/source/blender/editors/sculpt_paint/curves_sculpt_intern.hh index 5c926b1a740..4aaf6671cdb 100644 --- a/source/blender/editors/sculpt_paint/curves_sculpt_intern.hh +++ b/source/blender/editors/sculpt_paint/curves_sculpt_intern.hh @@ -97,14 +97,6 @@ IndexMask retrieve_selected_curves(const Curves &curves_id, Vector<int64_t> &r_i void move_last_point_and_resample(MutableSpan<float3> positions, const float3 &new_last_position); -float3 compute_surface_point_normal(const MLoopTri &looptri, - const float3 &bary_coord, - const Span<float3> corner_normals); - -float3 compute_bary_coord_in_triangle(const Mesh &mesh, - const MLoopTri &looptri, - const float3 &position); - class CurvesSculptCommonContext { public: const Depsgraph *depsgraph = nullptr; diff --git a/source/blender/geometry/CMakeLists.txt b/source/blender/geometry/CMakeLists.txt index 1916c5f91f3..531487a45e2 100644 --- a/source/blender/geometry/CMakeLists.txt +++ b/source/blender/geometry/CMakeLists.txt @@ -15,6 +15,7 @@ set(INC ) set(SRC + intern/add_curves_on_mesh.cc intern/mesh_merge_by_distance.cc intern/mesh_primitive_cuboid.cc intern/mesh_to_curve_convert.cc @@ -24,6 +25,7 @@ set(SRC intern/reverse_uv_sampler.cc intern/uv_parametrizer.c + GEO_add_curves_on_mesh.hh GEO_mesh_merge_by_distance.hh GEO_mesh_primitive_cuboid.hh GEO_mesh_to_curve.hh diff --git a/source/blender/geometry/GEO_add_curves_on_mesh.hh b/source/blender/geometry/GEO_add_curves_on_mesh.hh new file mode 100644 index 00000000000..cf60a8e8ace --- /dev/null +++ b/source/blender/geometry/GEO_add_curves_on_mesh.hh @@ -0,0 +1,54 @@ +/* SPDX-License-Identifier: GPL-2.0-or-later */ + +#pragma once + +#include "BLI_float4x4.hh" +#include "BLI_kdtree.h" +#include "BLI_math_vector.hh" +#include "BLI_span.hh" + +#include "BKE_bvhutils.h" +#include "BKE_curves.hh" + +#include "DNA_mesh_types.h" +#include "DNA_meshdata_types.h" + +namespace blender::geometry { + +struct AddCurvesOnMeshInputs { + /** Information about the root points where new curves should be generated. */ + Span<float3> root_positions_cu; + Span<float3> bary_coords; + Span<int> looptri_indices; + + /** Determines shape of new curves. */ + bool interpolate_length = false; + bool interpolate_shape = false; + bool interpolate_point_count = false; + float fallback_curve_length = 0.0f; + int fallback_point_count = 0; + + /** Information about the surface that the new curves are attached to. */ + const Mesh *surface = nullptr; + BVHTreeFromMesh *surface_bvh = nullptr; + Span<MLoopTri> surface_looptris; + Span<float2> surface_uv_map; + Span<float3> corner_normals_su; + + /** Transformation matrices. */ + float4x4 curves_to_surface_mat; + float4x4 surface_to_curves_normal_mat; + + /** + * KD-Tree that contains the root points of existing curves. This is only necessary when + * interpolation is used. + */ + KDTree_3d *old_roots_kdtree = nullptr; +}; + +/** + * Generate new curves on a mesh surface with the given inputs. Existing curves stay intact. + */ +void add_curves_on_mesh(bke::CurvesGeometry &curves, const AddCurvesOnMeshInputs &inputs); + +} // namespace blender::geometry diff --git a/source/blender/geometry/intern/add_curves_on_mesh.cc b/source/blender/geometry/intern/add_curves_on_mesh.cc new file mode 100644 index 00000000000..34551bd474f --- /dev/null +++ b/source/blender/geometry/intern/add_curves_on_mesh.cc @@ -0,0 +1,354 @@ +/* SPDX-License-Identifier: GPL-2.0-or-later */ + +#include "BKE_mesh_sample.hh" +#include "BKE_spline.hh" + +#include "GEO_add_curves_on_mesh.hh" + +/** + * The code below uses a suffix naming convention to indicate the coordinate space: + * cu: Local space of the curves object that is being edited. + * su: Local space of the surface object. + */ + +namespace blender::geometry { + +using bke::CurvesGeometry; + +struct NeighborCurve { + /* Curve index of the neighbor. */ + int index; + /* The weights of all neighbors of a new curve add up to 1. */ + float weight; +}; + +static constexpr int max_neighbors = 5; +using NeighborCurves = Vector<NeighborCurve, max_neighbors>; + +static float3 compute_surface_point_normal(const MLoopTri &looptri, + const float3 &bary_coord, + const Span<float3> corner_normals) +{ + const int l0 = looptri.tri[0]; + const int l1 = looptri.tri[1]; + const int l2 = looptri.tri[2]; + + const float3 &l0_normal = corner_normals[l0]; + const float3 &l1_normal = corner_normals[l1]; + const float3 &l2_normal = corner_normals[l2]; + + const float3 normal = math::normalize( + attribute_math::mix3(bary_coord, l0_normal, l1_normal, l2_normal)); + return normal; +} + +static void initialize_straight_curve_positions(const float3 &p1, + const float3 &p2, + MutableSpan<float3> r_positions) +{ + const float step = 1.0f / (float)(r_positions.size() - 1); + for (const int i : r_positions.index_range()) { + r_positions[i] = math::interpolate(p1, p2, i * step); + } +} + +static Array<NeighborCurves> find_curve_neighbors(const Span<float3> root_positions, + const KDTree_3d &old_roots_kdtree) +{ + const int tot_added_curves = root_positions.size(); + Array<NeighborCurves> neighbors_per_curve(tot_added_curves); + threading::parallel_for(IndexRange(tot_added_curves), 128, [&](const IndexRange range) { + for (const int i : range) { + const float3 root = root_positions[i]; + std::array<KDTreeNearest_3d, max_neighbors> nearest_n; + const int found_neighbors = BLI_kdtree_3d_find_nearest_n( + &old_roots_kdtree, root, nearest_n.data(), max_neighbors); + float tot_weight = 0.0f; + for (const int neighbor_i : IndexRange(found_neighbors)) { + KDTreeNearest_3d &nearest = nearest_n[neighbor_i]; + const float weight = 1.0f / std::max(nearest.dist, 0.00001f); + tot_weight += weight; + neighbors_per_curve[i].append({nearest.index, weight}); + } + /* Normalize weights. */ + for (NeighborCurve &neighbor : neighbors_per_curve[i]) { + neighbor.weight /= tot_weight; + } + } + }); + return neighbors_per_curve; +} + +template<typename T, typename GetValueF> +void interpolate_from_neighbors(const Span<NeighborCurves> neighbors_per_curve, + const T &fallback, + const GetValueF &get_value_from_neighbor, + MutableSpan<T> r_interpolated_values) +{ + attribute_math::DefaultMixer<T> mixer{r_interpolated_values}; + threading::parallel_for(r_interpolated_values.index_range(), 512, [&](const IndexRange range) { + for (const int i : range) { + const NeighborCurves &neighbors = neighbors_per_curve[i]; + if (neighbors.is_empty()) { + mixer.mix_in(i, fallback, 1.0f); + } + else { + for (const NeighborCurve &neighbor : neighbors) { + const T neighbor_value = get_value_from_neighbor(neighbor.index); + mixer.mix_in(i, neighbor_value, neighbor.weight); + } + } + } + }); + mixer.finalize(); +} + +static void interpolate_position_without_interpolation( + CurvesGeometry &curves, + const int old_curves_num, + const Span<float3> root_positions_cu, + const Span<float> new_lengths_cu, + const Span<float3> new_normals_su, + const float4x4 &surface_to_curves_normal_mat) +{ + const int added_curves_num = root_positions_cu.size(); + MutableSpan<float3> positions_cu = curves.positions_for_write(); + threading::parallel_for(IndexRange(added_curves_num), 256, [&](const IndexRange range) { + for (const int i : range) { + const int curve_i = old_curves_num + i; + const IndexRange points = curves.points_for_curve(curve_i); + const float3 &root_cu = root_positions_cu[i]; + const float length = new_lengths_cu[i]; + const float3 &normal_su = new_normals_su[i]; + const float3 normal_cu = math::normalize(surface_to_curves_normal_mat * normal_su); + const float3 tip_cu = root_cu + length * normal_cu; + + initialize_straight_curve_positions(root_cu, tip_cu, positions_cu.slice(points)); + } + }); +} + +static void interpolate_position_with_interpolation(CurvesGeometry &curves, + const Span<float3> root_positions_cu, + const Span<NeighborCurves> neighbors_per_curve, + const int old_curves_num, + const Span<float> new_lengths_cu, + const Span<float3> new_normals_su, + const float4x4 &surface_to_curves_normal_mat, + const float4x4 &curves_to_surface_mat, + const BVHTreeFromMesh &surface_bvh, + const Span<MLoopTri> surface_looptris, + const Mesh &surface, + const Span<float3> corner_normals_su) +{ + MutableSpan<float3> positions_cu = curves.positions_for_write(); + const int added_curves_num = root_positions_cu.size(); + + threading::parallel_for(IndexRange(added_curves_num), 256, [&](const IndexRange range) { + for (const int i : range) { + const NeighborCurves &neighbors = neighbors_per_curve[i]; + const int curve_i = old_curves_num + i; + const IndexRange points = curves.points_for_curve(curve_i); + + const float length_cu = new_lengths_cu[i]; + const float3 &normal_su = new_normals_su[i]; + const float3 normal_cu = math::normalize(surface_to_curves_normal_mat * normal_su); + + const float3 &root_cu = root_positions_cu[i]; + + if (neighbors.is_empty()) { + /* If there are no neighbors, just make a straight line. */ + const float3 tip_cu = root_cu + length_cu * normal_cu; + initialize_straight_curve_positions(root_cu, tip_cu, positions_cu.slice(points)); + continue; + } + + positions_cu.slice(points).fill(root_cu); + + for (const NeighborCurve &neighbor : neighbors) { + const int neighbor_curve_i = neighbor.index; + const float3 &neighbor_first_pos_cu = positions_cu[curves.offsets()[neighbor_curve_i]]; + const float3 neighbor_first_pos_su = curves_to_surface_mat * neighbor_first_pos_cu; + + BVHTreeNearest nearest; + nearest.dist_sq = FLT_MAX; + BLI_bvhtree_find_nearest(surface_bvh.tree, + neighbor_first_pos_su, + &nearest, + surface_bvh.nearest_callback, + const_cast<BVHTreeFromMesh *>(&surface_bvh)); + const int neighbor_looptri_index = nearest.index; + const MLoopTri &neighbor_looptri = surface_looptris[neighbor_looptri_index]; + + const float3 neighbor_bary_coord = + bke::mesh_surface_sample::compute_bary_coord_in_triangle( + surface, neighbor_looptri, nearest.co); + + const float3 neighbor_normal_su = compute_surface_point_normal( + surface_looptris[neighbor_looptri_index], neighbor_bary_coord, corner_normals_su); + const float3 neighbor_normal_cu = math::normalize(surface_to_curves_normal_mat * + neighbor_normal_su); + + /* The rotation matrix used to transform relative coordinates of the neighbor curve + * to the new curve. */ + float normal_rotation_cu[3][3]; + rotation_between_vecs_to_mat3(normal_rotation_cu, neighbor_normal_cu, normal_cu); + + const IndexRange neighbor_points = curves.points_for_curve(neighbor_curve_i); + const float3 &neighbor_root_cu = positions_cu[neighbor_points[0]]; + + /* Use a temporary #PolySpline, because that's the easiest way to resample an + * existing curve right now. Resampling is necessary if the length of the new curve + * does not match the length of the neighbors or the number of handle points is + * different. */ + PolySpline neighbor_spline; + neighbor_spline.resize(neighbor_points.size()); + neighbor_spline.positions().copy_from(positions_cu.slice(neighbor_points)); + neighbor_spline.mark_cache_invalid(); + + const float neighbor_length_cu = neighbor_spline.length(); + const float length_factor = std::min(1.0f, length_cu / neighbor_length_cu); + + const float resample_factor = (1.0f / (points.size() - 1.0f)) * length_factor; + for (const int j : IndexRange(points.size())) { + const Spline::LookupResult lookup = neighbor_spline.lookup_evaluated_factor( + j * resample_factor); + const float index_factor = lookup.evaluated_index + lookup.factor; + float3 p; + neighbor_spline.sample_with_index_factors<float3>( + neighbor_spline.positions(), {&index_factor, 1}, {&p, 1}); + const float3 relative_coord = p - neighbor_root_cu; + float3 rotated_relative_coord = relative_coord; + mul_m3_v3(normal_rotation_cu, rotated_relative_coord); + positions_cu[points[j]] += neighbor.weight * rotated_relative_coord; + } + } + } + }); +} + +void add_curves_on_mesh(CurvesGeometry &curves, const AddCurvesOnMeshInputs &inputs) +{ + const bool use_interpolation = inputs.interpolate_length || inputs.interpolate_point_count || + inputs.interpolate_shape; + + Array<NeighborCurves> neighbors_per_curve; + if (use_interpolation) { + BLI_assert(inputs.old_roots_kdtree != nullptr); + neighbors_per_curve = find_curve_neighbors(inputs.root_positions_cu, *inputs.old_roots_kdtree); + } + + const int added_curves_num = inputs.root_positions_cu.size(); + const int old_points_num = curves.points_num(); + const int old_curves_num = curves.curves_num(); + const int new_curves_num = old_curves_num + added_curves_num; + + /* Grow number of curves first, so that the offsets array can be filled. */ + curves.resize(old_points_num, new_curves_num); + + /* Compute new curve offsets. */ + MutableSpan<int> curve_offsets = curves.offsets_for_write(); + MutableSpan<int> new_point_counts_per_curve = curve_offsets.take_back(added_curves_num); + if (inputs.interpolate_point_count) { + interpolate_from_neighbors<int>( + neighbors_per_curve, + inputs.fallback_point_count, + [&](const int curve_i) { return curves.points_for_curve(curve_i).size(); }, + new_point_counts_per_curve); + } + else { + new_point_counts_per_curve.fill(inputs.fallback_point_count); + } + for (const int i : IndexRange(added_curves_num)) { + curve_offsets[old_curves_num + i + 1] += curve_offsets[old_curves_num + i]; + } + + const int new_points_num = curves.offsets().last(); + curves.resize(new_points_num, new_curves_num); + MutableSpan<float3> positions_cu = curves.positions_for_write(); + + /* Determine length of new curves. */ + Array<float> new_lengths_cu(added_curves_num); + if (inputs.interpolate_length) { + interpolate_from_neighbors<float>( + neighbors_per_curve, + inputs.fallback_curve_length, + [&](const int curve_i) { + const IndexRange points = curves.points_for_curve(curve_i); + float length = 0.0f; + for (const int segment_i : points.drop_back(1)) { + const float3 &p1 = positions_cu[segment_i]; + const float3 &p2 = positions_cu[segment_i + 1]; + length += math::distance(p1, p2); + } + return length; + }, + new_lengths_cu); + } + else { + new_lengths_cu.fill(inputs.fallback_curve_length); + } + + /* Find surface normal at root points. */ + Array<float3> new_normals_su(added_curves_num); + threading::parallel_for(IndexRange(added_curves_num), 256, [&](const IndexRange range) { + for (const int i : range) { + const int looptri_index = inputs.looptri_indices[i]; + const float3 &bary_coord = inputs.bary_coords[i]; + new_normals_su[i] = compute_surface_point_normal( + inputs.surface_looptris[looptri_index], bary_coord, inputs.corner_normals_su); + } + }); + + /* Propagate attachment information. */ + if (!inputs.surface_uv_map.is_empty()) { + MutableSpan<float2> surface_uv_coords = curves.surface_uv_coords_for_write(); + bke::mesh_surface_sample::sample_corner_attribute( + *inputs.surface, + inputs.looptri_indices, + inputs.bary_coords, + GVArray::ForSpan(inputs.surface_uv_map), + IndexRange(added_curves_num), + surface_uv_coords.take_back(added_curves_num)); + } + + /* Update selection arrays when available. */ + const VArray<float> points_selection = curves.selection_point_float(); + if (points_selection.is_span()) { + MutableSpan<float> points_selection_span = curves.selection_point_float_for_write(); + points_selection_span.drop_front(old_points_num).fill(1.0f); + } + const VArray<float> curves_selection = curves.selection_curve_float(); + if (curves_selection.is_span()) { + MutableSpan<float> curves_selection_span = curves.selection_curve_float_for_write(); + curves_selection_span.drop_front(old_curves_num).fill(1.0f); + } + + /* Initialize position attribute. */ + if (inputs.interpolate_shape) { + interpolate_position_with_interpolation(curves, + inputs.root_positions_cu, + neighbors_per_curve, + old_curves_num, + new_lengths_cu, + new_normals_su, + inputs.surface_to_curves_normal_mat, + inputs.curves_to_surface_mat, + *inputs.surface_bvh, + inputs.surface_looptris, + *inputs.surface, + inputs.corner_normals_su); + } + else { + interpolate_position_without_interpolation(curves, + old_curves_num, + inputs.root_positions_cu, + new_lengths_cu, + new_normals_su, + inputs.surface_to_curves_normal_mat); + } + + curves.update_curve_types(); +} + +} // namespace blender::geometry |