From 3b0a861b93ae46a5737e85f8eff3f755989f06c5 Mon Sep 17 00:00:00 2001 From: GabijaBern Date: Wed, 6 Nov 2024 12:34:00 +0000 Subject: [PATCH 1/2] Option to crop interest zone and filter min dbh trees --- raycloudtools/rayextract/rayextract.cpp | 24 +- raylib/extraction/raytrees.cpp | 438 ++++++++++++++++++------ raylib/extraction/raytrees.h | 52 +-- raylib/raylaz.cpp | 10 + raylib/rayparse.cpp | 26 +- 5 files changed, 422 insertions(+), 128 deletions(-) diff --git a/raycloudtools/rayextract/rayextract.cpp b/raycloudtools/rayextract/rayextract.cpp index 5f0d25cd..d6c82888 100644 --- a/raycloudtools/rayextract/rayextract.cpp +++ b/raycloudtools/rayextract/rayextract.cpp @@ -55,12 +55,14 @@ void usage(int exit_code = 1) std::cout << " --crop_length 1.0 - (-p) crops small branches to this distance from end" << std::endl; std::cout << " --distance_limit 1 - (-d) maximum distance between neighbour points in a tree" << std::endl; std::cout << " --height_min 2 - (-h) minimum height counted as a tree" << std::endl; + std::cout << " --min_radius 0.01 - (-r) minimum tree radius to consider" << std::endl; std::cout << " --girth_height_ratio 0.12 - (-i) the amount up tree's height to estimate trunk girth" << std::endl; std::cout << " --global_taper 0.024 - (-a) force a taper value (diameter per length) for trees under global_taper_factor of max tree height. Use 0 to estimate global taper from the data" << std::endl; std::cout << " --global_taper_factor 0.3- (-o) 1 estimates same taper for whole scan, 0 is per-tree tapering. Like a soft cutoff at this amount of max tree height" << std::endl; std::cout << " --gravity_factor 0.3 - (-f) larger values preference vertical trees" << std::endl; std::cout << " --branch_segmentation- (-b) _segmented.ply is per branch segment" << std::endl; std::cout << " --grid_width 10 - (-w) crops results assuming cloud has been gridded with given width" << std::endl; + std::cout << " --grid_origin 0,0,0 - location of origin within grid cell that overlaps it. Defaults to a cell-centre origin (at grid_width/2 in each axis) matching raysplit grid. 0,0,0 is for a vertex origin." << std::endl; std::cout << " --use_rays - (-u) use rays to reduce trunk radius overestimation in noisy cloud data" << std::endl; std::cout << " (for internal constants -c -g -s see source file rayextract)" << std::endl; // These are the internal parameters that I don't expose as they are 'advanced' only, you shouldn't need to adjust them @@ -92,6 +94,9 @@ int rayExtract(int argc, char *argv[]) } ray::FileArgument cloud_file, mesh_file, trunks_file, trees_file, leaf_file; ray::TextArgument forest("forest"), trees("trees"), trunks("trunks"), terrain("terrain"), leaves("leaves"); + + ray::Vector3dArgument grid_origin; + ray::OptionalKeyValueArgument grid_origin_option("grid_origin", &grid_origin); ray::OptionalKeyValueArgument groundmesh_option("ground", 'g', &mesh_file); ray::OptionalKeyValueArgument trunks_option("trunks", 't', &trunks_file); ray::DoubleArgument gradient(0.001, 1000.0, 1.0), global_taper(0.0, 1.0), global_taper_factor(0.0, 1.0); @@ -102,7 +107,7 @@ int rayExtract(int argc, char *argv[]) ray::DoubleArgument max_diameter(0.01, 100.0), distance_limit(0.01, 10.0), height_min(0.01, 1000.0), min_diameter(0.01, 100.0), leaf_area(0.00001, 1.0, 0.002), leaf_droop(0.0, 10.0, 0.1), crop_length(0.01, 100.0);; ray::DoubleArgument girth_height_ratio(0.001, 0.5), length_to_radius(0.01, 10000.0), cylinder_length_to_width(0.1, 20.0), gap_ratio(0.01, 10.0), - span_ratio(0.01, 10.0); + span_ratio(0.01, 10.0), min_radius(0.01, 100.0); ray::DoubleArgument gravity_factor(0.0, 100.0), grid_width(1.0, 100000.0), grid_overlap(0.0, 0.9); ray::OptionalKeyValueArgument max_diameter_option("max_diameter", 'm', &max_diameter); @@ -112,6 +117,8 @@ int rayExtract(int argc, char *argv[]) ray::OptionalKeyValueArgument girth_height_ratio_option("girth_height_ratio", 'i', &girth_height_ratio); ray::OptionalKeyValueArgument cylinder_length_to_width_option("cylinder_length_to_width", 'c', &cylinder_length_to_width); + + ray::OptionalKeyValueArgument min_radius_option("min_radius", 'r', &min_radius); ray::OptionalKeyValueArgument gap_ratio_option("gap_ratio", 'g', &gap_ratio); ray::OptionalKeyValueArgument span_ratio_option("span_ratio", 's', &span_ratio); ray::OptionalKeyValueArgument gravity_factor_option("gravity_factor", 'f', &gravity_factor); @@ -137,7 +144,8 @@ int rayExtract(int argc, char *argv[]) argc, argv, { &trees, &cloud_file, &mesh_file }, { &max_diameter_option, &distance_limit_option, &height_min_option, &crop_length_option, &girth_height_ratio_option, &cylinder_length_to_width_option, &gap_ratio_option, &span_ratio_option, &gravity_factor_option, - &segment_branches, &grid_width_option, &global_taper_option, &global_taper_factor_option, &use_rays, &verbose }); + &segment_branches, &grid_width_option, &global_taper_option, &global_taper_factor_option, &use_rays, &verbose, + &grid_origin_option, &min_radius_option }); bool extract_leaves = ray::parseCommandLine(argc, argv, { &leaves, &cloud_file, &trees_file }, { &leaf_option, &leaf_area_option, &leaf_droop_option, &stalks }); @@ -226,8 +234,16 @@ int rayExtract(int argc, char *argv[]) if (global_taper_factor_option.isSet()) { params.global_taper_factor = global_taper_factor.value(); - } - params.use_rays = use_rays.isSet(); + } + if (grid_origin_option.isSet()) + { + params.grid_origin = grid_origin.value(); + } + if (min_radius_option.isSet()) + { + params.min_radius = min_radius.value(); + } + params.use_rays = use_rays.isSet(); params.segment_branches = segment_branches.isSet(); ray::Trees trees(cloud, offset, mesh, params, verbose.isSet()); diff --git a/raylib/extraction/raytrees.cpp b/raylib/extraction/raytrees.cpp index 99233cee..94b8485a 100644 --- a/raylib/extraction/raytrees.cpp +++ b/raylib/extraction/raytrees.cpp @@ -23,6 +23,7 @@ TreesParams::TreesParams() , segment_branches(false) , global_taper(0.012) , global_taper_factor(0.3) + , grid_origin(0, 0, 0) {} /// The main reconstruction algorithm @@ -75,19 +76,19 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons sections_[sec_].tree_height = tree_height; // Use a usr defined taper to control the height up the trunk to calculate the radius at - double girth_height = params_->girth_height_ratio * tree_height; // sections_[sec_].max_distance_to_end; + double girth_height = params_->girth_height_ratio * tree_height; // sections_[sec_].max_distance_to_end; double estimated_radius = 1e10; double best_dist = 0.0; Eigen::Vector3d best_tip; - + std::vector best_nodes; std::vector best_ends; - for (int j = 1; j<=3; j++) + for (int j = 1; j <= 3; j++) { - double max_dist = girth_height * (double)j / 2.0; // range from 0.5 to 1.5 times the specified height + double max_dist = girth_height * (double)j / 2.0; // range from 0.5 to 1.5 times the specified height nodes.clear(); sections_[sec_].ends.clear(); - extractNodesAndEndsFromRoots(nodes, base, children, max_dist * 2.0/3.0, max_dist); + extractNodesAndEndsFromRoots(nodes, base, children, max_dist * 2.0 / 3.0, max_dist); if (nodes.size() < 2) { continue; @@ -95,33 +96,36 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons sections_[sec_].tip = calculateTipFromVertices(nodes); if (verbose) { - for (auto &node: nodes) + for (auto &node : nodes) { - debug_cloud.addRay(Eigen::Vector3d(0,0,0), points_[node].pos, 0.0, ray::RGBA(j==1 ? 255 : 0, j==2 ? 255:0, j==3 ? 255:0, 255)); + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), points_[node].pos, 0.0, + ray::RGBA(j == 1 ? 255 : 0, j == 2 ? 255 : 0, j == 3 ? 255 : 0, 255)); } - } + } if (removeDistantPoints(nodes)) { sections_[sec_].tip = calculateTipFromVertices(nodes); } if (verbose) { - debug_cloud.addRay(Eigen::Vector3d(0,0,0), sections_[sec_].tip + Eigen::Vector3d(0,0,0.03), 0.0, ray::RGBA(255,255,0, 255)); + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), sections_[sec_].tip + Eigen::Vector3d(0, 0, 0.03), 0.0, + ray::RGBA(255, 255, 0, 255)); } // shift to cylinder's centre - Eigen::Vector3d up(0,0,1); + Eigen::Vector3d up(0, 0, 1); sections_[sec_].tip += vectorToCylinderCentre(nodes, up); // now find the segment radius double accuracy; double radius = estimateCylinderRadius(nodes, up, accuracy); - ray::RGBA col(j==1 ? 255 : 127, j==2 ? 255:127, j==3 ? 255:127, 255); + ray::RGBA col(j == 1 ? 255 : 127, j == 2 ? 255 : 127, j == 3 ? 255 : 127, 255); if (verbose) { - debug_cloud.addRay(Eigen::Vector3d(0,0,0), sections_[sec_].tip, 0.0, col); - for (double ang = 0; ang < 2.0*ray::kPi; ang += 0.1) + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), sections_[sec_].tip, 0.0, col); + for (double ang = 0; ang < 2.0 * ray::kPi; ang += 0.1) { - debug_cloud.addRay(Eigen::Vector3d(0,0,0), sections_[sec_].tip + radius * Eigen::Vector3d(std::sin(ang), std::cos(ang),0), 0.0, col); + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), + sections_[sec_].tip + radius * Eigen::Vector3d(std::sin(ang), std::cos(ang), 0), 0.0, col); } } if (radius < estimated_radius) @@ -136,21 +140,30 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons } if (best_dist == 0.0) { - std::cout << "warning: could not find any points on trunk " << sec_ << " at " << base.transpose() << " so removing the whole section" << std::endl; - sections_[sec_].tip = base + Eigen::Vector3d(0,0,0.01); + std::cout << "warning: could not find any points on trunk " << sec_ << " at " << base.transpose() + << " so removing the whole section" << std::endl; + sections_[sec_].tip = base + Eigen::Vector3d(0, 0, 0.01); + sections_[sec_].total_weight = 1e-10; + sections_[sec_].ends.clear(); // so this trunk is not ever used + continue; + } + + if (params_->min_radius && estimated_radius < params_->min_radius) + { + sections_[sec_].tip = base + Eigen::Vector3d(0, 0, 0.01); sections_[sec_].total_weight = 1e-10; - sections_[sec_].ends.clear(); // so this trunk is not ever used - continue; - } + sections_[sec_].ends.clear(); + continue; + } sections_[sec_].tip = best_tip; sections_[sec_].ends = best_ends; nodes = best_nodes; if (sections_[sec_].split_count < 2) { - double thickness = best_dist; + double thickness = best_dist; bool points_removed = false; - double gap = params_->gap_ratio * sections_[sec_].max_distance_to_end; // gap threshold for splitting - double span = params_->span_ratio * estimated_radius; // span threshold for splitting + double gap = params_->gap_ratio * sections_[sec_].max_distance_to_end; // gap threshold for splitting + double span = params_->span_ratio * estimated_radius; // span threshold for splitting std::vector> clusters = findPointClusters(base, points_removed, thickness, span, gap); if (clusters.size() > 1 || (points_removed && clusters.size() > 0)) // a bifurcation (or an alteration) @@ -163,21 +176,24 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons } if (verbose) { - for (auto &node: sections_[sec_].ends) + for (auto &node : sections_[sec_].ends) { - debug_cloud.addRay(Eigen::Vector3d(0,0,0), points_[node].pos + Eigen::Vector3d(0,0,0.02), 0.0, ray::RGBA(255, 0, 255, 255)); + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), points_[node].pos + Eigen::Vector3d(0, 0, 0.02), 0.0, + ray::RGBA(255, 0, 255, 255)); } - for (double ang = 0; ang < 2.0*ray::kPi; ang += 0.1) + for (double ang = 0; ang < 2.0 * ray::kPi; ang += 0.1) { - uint8_t shade = 255; // (uint8_t)(best_accuracy * 255.0); - debug_cloud.addRay(Eigen::Vector3d(0,0,0), best_tip + (estimated_radius + 0.01) * Eigen::Vector3d(std::sin(ang), std::cos(ang),0), 0.0, ray::RGBA(shade,shade,shade,255)); + uint8_t shade = 255; // (uint8_t)(best_accuracy * 255.0); + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), + best_tip + (estimated_radius + 0.01) * Eigen::Vector3d(std::sin(ang), std::cos(ang), 0), 0.0, + ray::RGBA(shade, shade, shade, 255)); } } nodes.clear(); sections_[sec_].ends.clear(); - extractNodesAndEndsFromRoots(nodes, base, children, 0.0, best_dist/2.0); // make it lower - estimateCylinderTaper(estimated_radius, best_accuracy, false); // update the expected taper + extractNodesAndEndsFromRoots(nodes, base, children, 0.0, best_dist / 2.0); // make it lower + estimateCylinderTaper(estimated_radius, best_accuracy, false); // update the expected taper } if (verbose) { @@ -196,11 +212,11 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons if (sections_[sec_].ends.size() > 0) { addChildSection(); - for (const auto &i: sections_[sec_].roots) + for (const auto &i : sections_[sec_].roots) { sections_[sec_].tip[2] = std::min(sections_[sec_].tip[2], points_[i].pos[2]); } - } + } else { std::cout << "weird, a trunk without end points! " << sec_ << std::endl; @@ -226,26 +242,27 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons if (!extract_from_ends) { Eigen::Vector3d base = getRootPosition(); - double span_rad = radius(sections_[sec_]); + double span_rad = radius(sections_[sec_]); double thickness = params_->cylinder_length_to_width * span_rad; extractNodesAndEndsFromRoots(nodes, base, children, 0.0, thickness); - + bool points_removed = false; - double gap = params_->gap_ratio * sections_[sec_].max_distance_to_end; // gap threshold for splitting - double span = params_->span_ratio * span_rad; // span thershold for splitting + double gap = params_->gap_ratio * sections_[sec_].max_distance_to_end; // gap threshold for splitting + double span = params_->span_ratio * span_rad; // span thershold for splitting std::vector> clusters = findPointClusters(base, points_removed, thickness, span, gap); if (clusters.size() > 1 || (points_removed && clusters.size() > 0)) // a bifurcation (or an alteration) { - extract_from_ends = true; // don't trust the found nodes as it is now two separate tree nodes - bool add_offshoots = true; // par != -1 && sections_[par].parent != -1; // when this is false it can lead to whole branches missing. + extract_from_ends = true; // don't trust the found nodes as it is now two separate tree nodes + bool add_offshoots = true; // par != -1 && sections_[par].parent != -1; // when this is false it can lead to + // whole branches missing. // if points have been removed then this only resets the current section's points // otherwise it creates new branch sections_ for each cluster and adds to the end of the sections_ list bifurcate(clusters, thickness, children, false, add_offshoots); } } - if (extract_from_ends) // we have split the ends, so we need to extract the set of nodes in a backwards manner + if (extract_from_ends) // we have split the ends, so we need to extract the set of nodes in a backwards manner { extractNodesFromEnds(nodes); } @@ -283,6 +300,11 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons removeOutOfBoundSections(cloud, min_bound, max_bound, offset); } + if (params_->min_radius) + { + removeSmallRadiusTrees(); + } + std::vector root_segs(cloud.ends.size(), -1); // now colour the ray cloud based on the segmentation segmentCloud(cloud, root_segs, section_ids); @@ -292,28 +314,29 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons removeOutOfBoundRays(cloud, min_bound, max_bound, root_segs); } - std::cout << "cloud's estimated mean taper ratio (diameter / length): " << 2.0 * forest_taper_ / forest_weight_ << std::endl; + std::cout << "cloud's estimated mean taper ratio (diameter / length): " << 2.0 * forest_taper_ / forest_weight_ + << std::endl; } -// If 1 tree plus some low lying foliage is in the node, then it won't be split if the foliage doesn't reach to the -// top of the section (so doesn't have an end point), but it can create a very large radius estimate. -// here we remove points that are too far from any end position. -// this allows us to use most nodes to get an accurate radius estimate, but still uses only the end nodes to decide whether -// to split +// If 1 tree plus some low lying foliage is in the node, then it won't be split if the foliage doesn't reach to the +// top of the section (so doesn't have an end point), but it can create a very large radius estimate. +// here we remove points that are too far from any end position. +// this allows us to use most nodes to get an accurate radius estimate, but still uses only the end nodes to decide +// whether to split bool Trees::removeDistantPoints(std::vector &nodes) { - Eigen::Vector3d up(0,0,1); - double max_rad = params_->gap_ratio * sections_[sec_].max_distance_to_end; // gap threshold for splitting + Eigen::Vector3d up(0, 0, 1); + double max_rad = params_->gap_ratio * sections_[sec_].max_distance_to_end; // gap threshold for splitting double max_rad_sqr = max_rad * max_rad; bool found = false; - for (int k = 0; k<(int)nodes.size(); k++) + for (int k = 0; k < (int)nodes.size(); k++) { if (nodes.size() <= 6) { break; } bool all_distant = true; - for (auto &end: sections_[sec_].ends) + for (auto &end : sections_[sec_].ends) { Eigen::Vector3d dif = points_[nodes[k]].pos - points_[end].pos; dif[2] = 0.0; @@ -345,18 +368,20 @@ double Trees::meanTaper(const BranchSection §ion) const mean_weight *= blend; weight *= 1.0 - blend; - double result = (mean_taper * mean_weight + taper*weight) / (mean_weight + weight); + double result = (mean_taper * mean_weight + taper * weight) / (mean_weight + weight); if (!(result == result)) { std::cout << "bad taper estimate" << std::endl; - std::cout << "root: " << root << " estimated taper: " << result << ", mt: " << mean_taper << ", mw: " << mean_weight << ", t: " << taper << ", w: " << weight << ", blend: " << blend << std::endl; + std::cout << "root: " << root << " estimated taper: " << result << ", mt: " << mean_taper << ", mw: " << mean_weight + << ", t: " << taper << ", w: " << weight << ", blend: " << blend << std::endl; } return result; } -double Trees::radius(const BranchSection §ion) const -{ - return sections_[section.root].tree_height * meanTaper(section) * section.radius_scale; // scaled down from root's estimated radius +double Trees::radius(const BranchSection §ion) const +{ + return sections_[section.root].tree_height * meanTaper(section) * + section.radius_scale; // scaled down from root's estimated radius } void Trees::calculatePointDistancesToEnd() @@ -447,7 +472,8 @@ Eigen::Vector3d Trees::getRootPosition() const // note that the nodes contain everything between min_dist and max_dist and also the adjacent nodes either side // i.e. a root node and an end node. This means the minimum number of nodes is 2 void Trees::extractNodesAndEndsFromRoots(std::vector &nodes, const Eigen::Vector3d &base, - const std::vector> &children, double min_dist, double max_dist) + const std::vector> &children, double min_dist, + double max_dist) { const int par = sections_[sec_].parent; std::vector iteration_nodes = sections_[sec_].roots; @@ -487,7 +513,7 @@ void Trees::extractNodesAndEndsFromRoots(std::vector &nodes, const Eigen::V // find clusters of points from the root points up the shortest paths, up to the cylinder length std::vector> Trees::findPointClusters(const Eigen::Vector3d &base, bool &points_removed, - double thickness, double span, double gap) + double thickness, double span, double gap) { const int par = sections_[sec_].parent; std::vector &all_ends = sections_[sec_].ends; @@ -548,7 +574,8 @@ std::vector> Trees::findPointClusters(const Eigen::Vector3d &ba // split into multiple branches and add as new branch sections to the end of the sections_ list // that is being iterated through. -void Trees::bifurcate(const std::vector> &clusters, double thickness, std::vector> &children, bool clip_tree, bool add_offshoots) +void Trees::bifurcate(const std::vector> &clusters, double thickness, + std::vector> &children, bool clip_tree, bool add_offshoots) { const int par = sections_[sec_].parent; // find the maximum distance (to tip) for each cluster @@ -570,14 +597,14 @@ void Trees::bifurcate(const std::vector> &clusters, double thic } double total_area = 0.0; - for (auto &dist: max_distances) + for (auto &dist : max_distances) { total_area += dist * dist; } // if clip_tree then we need to change the root nodes to be only those ones that - // each path goes to... - + // each path goes to... + // set the current branch section to the cluster with the // maximum distance to tip (i.e. the longest branch). sections_[sec_].ends = clusters[maxi]; @@ -640,9 +667,9 @@ void Trees::bifurcate(const std::vector> &clusters, double thic int child = end; while (node != -1) { - if (std::find(my_nodes.begin(), my_nodes.end(), node) != my_nodes.end()) + if (std::find(my_nodes.begin(), my_nodes.end(), node) != my_nodes.end()) { - break; // just hit my own tree + break; // just hit my own tree } if (std::find(all_nodes.begin(), all_nodes.end(), node) != all_nodes.end()) { @@ -652,7 +679,7 @@ void Trees::bifurcate(const std::vector> &clusters, double thic auto &list = children[node]; int find_count = 0; - for (int j = 0; j<(int)list.size(); j++) + for (int j = 0; j < (int)list.size(); j++) { if (list[j] == child) { @@ -665,8 +692,8 @@ void Trees::bifurcate(const std::vector> &clusters, double thic points_[child].parent = -1; // we need to tell all of the children what the new root is - std::vector child_list = {child}; - for (size_t j = 0; j child_list = { child }; + for (size_t j = 0; j < child_list.size(); j++) { points_[child_list[j]].root = child; child_list.insert(child_list.end(), children[child_list[j]].begin(), children[child_list[j]].end()); @@ -674,7 +701,8 @@ void Trees::bifurcate(const std::vector> &clusters, double thic break; } my_nodes.push_back(node); // fill in the nodes in this branch section - if (std::find(sections_[sec_].roots.begin(), sections_[sec_].roots.end(), node) != sections_[sec_].roots.end()) + if (std::find(sections_[sec_].roots.begin(), sections_[sec_].roots.end(), node) != + sections_[sec_].roots.end()) { my_roots.push_back(node); break; @@ -705,7 +733,8 @@ void Trees::bifurcate(const std::vector> &clusters, double thic } if (par > -1) { - sections_[sec_].radius_scale *= max_distances[maxi] / std::sqrt(total_area); // reduce branch radius due to the split + sections_[sec_].radius_scale *= + max_distances[maxi] / std::sqrt(total_area); // reduce branch radius due to the split } } @@ -832,7 +861,7 @@ double Trees::estimateCylinderRadius(const std::vector &nodes, const Eigen: { double rad = 0.0; // get the mean radius - double power = 0.25; // when 1 this is usual radius, but lower powers reduce outliers due to foliage + double power = 0.25; // when 1 this is usual radius, but lower powers reduce outliers due to foliage std::vector norms; norms.reserve(nodes.size()); if (params_->use_rays) @@ -840,10 +869,10 @@ double Trees::estimateCylinderRadius(const std::vector &nodes, const Eigen: for (auto &node : nodes) { Eigen::Vector3d offset = points_[node].pos - sections_[sec_].tip; - offset -= dir * offset.dot(dir); // flatten + offset -= dir * offset.dot(dir); // flatten Eigen::Vector3d ray = points_[node].start - points_[node].pos; - ray -= dir * ray.dot(dir); // flatten - offset += ray*std::max(0.0, std::min(-offset.dot(ray)/ray.dot(ray), 1.0)); // move to closest point + ray -= dir * ray.dot(dir); // flatten + offset += ray * std::max(0.0, std::min(-offset.dot(ray) / ray.dot(ray), 1.0)); // move to closest point norms.push_back(offset.norm()); } } @@ -852,18 +881,18 @@ double Trees::estimateCylinderRadius(const std::vector &nodes, const Eigen: for (auto &node : nodes) { Eigen::Vector3d offset = points_[node].pos - sections_[sec_].tip; - offset -= dir * offset.dot(dir); // flatten + offset -= dir * offset.dot(dir); // flatten norms.push_back(offset.norm()); } } - for (auto &norm: norms) + for (auto &norm : norms) { rad += std::pow(norm, power); } - const double eps = 1e-5; // prevent division by 0 + const double eps = 1e-5; // prevent division by 0 rad /= (double)nodes.size() + eps; - rad = std::pow(rad, 1.0/power); + rad = std::pow(rad, 1.0 / power); double e = 0.0; for (auto &norm : norms) { @@ -872,10 +901,12 @@ double Trees::estimateCylinderRadius(const std::vector &nodes, const Eigen: e /= (double)nodes.size() + eps; #define NEW_ACCURACY #if defined NEW_ACCURACY - const double sensor_noise = 0.02; // this prevents cylinders with a small number of very accurate points from totally dominating + const double sensor_noise = + 0.02; // this prevents cylinders with a small number of very accurate points from totally dominating accuracy = rad / (e + sensor_noise); -#else // this one goes down to 0, which could be bad if all cylinder points were low accuracy - accuracy = std::max(eps, rad - 2.0*e) / std::max(eps, rad); // so even distribution is accuracy 0, and cylinder shell is accuracy 1 +#else // this one goes down to 0, which could be bad if all cylinder points were low accuracy + accuracy = std::max(eps, rad - 2.0 * e) / + std::max(eps, rad); // so even distribution is accuracy 0, and cylinder shell is accuracy 1 accuracy *= std::min((double)nodes.size() / 3.0, 1.0); #endif accuracy = std::max(accuracy, eps); @@ -895,16 +926,17 @@ void Trees::estimateCylinderTaper(double radius, double accuracy, bool extract_f double junction_weight = 1.0; if (par > -1) { - junction_weight = extract_from_ends ? 0.25 : 1.0; // extracting from ends means we are less confident about the quality + junction_weight = + extract_from_ends ? 0.25 : 1.0; // extracting from ends means we are less confident about the quality } - double weight = l*l*l * accuracy * junction_weight; - // weight *= weight; // preference the strongest weight sections - double taper = (radius/L) * weight; + double weight = l * l * l * accuracy * junction_weight; + // weight *= weight; // preference the strongest weight sections + double taper = (radius / L) * weight; forest_taper_ += taper; forest_weight_ += weight; - forest_weight_squared_ += weight*weight; + forest_weight_squared_ += weight * weight; sections_[root].total_taper += taper; sections_[root].total_weight += weight; @@ -916,7 +948,7 @@ void Trees::estimateCylinderTaper(double radius, double accuracy, bool extract_f sections_[sec_].taper = taper; sections_[sec_].weight = weight; - sections_[sec_].len = l*l*l; + sections_[sec_].len = l * l * l; sections_[sec_].accuracy = accuracy; sections_[sec_].junction_weight = junction_weight; } @@ -931,12 +963,12 @@ void Trees::addChildSection() new_node.taper = sections_[sec_].taper; new_node.weight = sections_[sec_].weight; new_node.radius_scale = sections_[sec_].radius_scale; - + for (auto &root : new_node.roots) { new_node.max_distance_to_end = std::max(new_node.max_distance_to_end, points_[root].distance_to_end); } - if (new_node.max_distance_to_end > params_->crop_length) + if (new_node.max_distance_to_end > params_->crop_length) { sections_[sec_].children.push_back(static_cast(sections_.size())); new_node.tip = calculateTipFromVertices(new_node.roots); @@ -971,8 +1003,7 @@ void Trees::calculateSectionIds(std::vector §ion_ids, const std::vector break; } nodes.push_back(node); - if (std::find(sections_[sec_].roots.begin(), sections_[sec_].roots.end(), node) != - sections_[sec_].roots.end()) + if (std::find(sections_[sec_].roots.begin(), sections_[sec_].roots.end(), node) != sections_[sec_].roots.end()) { break; } @@ -1046,17 +1077,19 @@ void Trees::generateLocalSectionIds() std::cout << num << " trees saved" << std::endl; } -void Trees::removeOutOfBoundSections(const Cloud &cloud, Eigen::Vector3d &min_bound, Eigen::Vector3d &max_bound, const Eigen::Vector3d &offset) +void Trees::removeOutOfBoundSections(const Cloud &cloud, Eigen::Vector3d &min_bound, Eigen::Vector3d &max_bound, + const Eigen::Vector3d &offset) { const double width = params_->grid_width; cloud.calcBounds(&min_bound, &max_bound); - const Eigen::Vector3d mid = (min_bound + max_bound)/2.0 + offset; + const Eigen::Vector3d mid = (min_bound + max_bound) / 2.0 + offset + params_->grid_origin; const Eigen::Vector2d inds(std::round(mid[0] / width), std::round(mid[1] / width)); min_bound[0] = width * (inds[0] - 0.5) - offset[0]; min_bound[1] = width * (inds[1] - 0.5) - offset[1]; max_bound[0] = width * (inds[0] + 0.5) - offset[0]; max_bound[1] = width * (inds[1] + 0.5) - offset[1]; - std::cout << "min bound: " << (min_bound+offset).transpose() << ", max bound: " << (max_bound+offset).transpose() << std::endl; + std::cerr << "min bound: " << (min_bound + offset).transpose() << ", max bound: " << (max_bound + offset).transpose() + << " offset " << offset << std::endl; // disable trees out of bounds for (auto §ion : sections_) @@ -1076,7 +1109,8 @@ void Trees::removeOutOfBoundSections(const Cloud &cloud, Eigen::Vector3d &min_bo // colour the cloud by tree id, or by branch segment id void Trees::segmentCloud(Cloud &cloud, std::vector &root_segs, const std::vector §ion_ids) { - contiguous_section_ids_.resize(sections_.size(), -1); // these are different to the root section IDs as they exclude empty trees + contiguous_section_ids_.resize(sections_.size(), + -1); // these are different to the root section IDs as they exclude empty trees int num_trees = 0; bool first_non_root = false; for (size_t sec = 0; sec < sections_.size(); sec++) @@ -1106,21 +1140,22 @@ void Trees::segmentCloud(Cloud &cloud, std::vector &root_segs, const std::v const int root_id = points_[j].root; root_segs[i] = root_id == -1 ? -1 : section_ids[root_id]; - // This block places a cone of gradient 2 around the trunk direction, to remove the square of initial ground points + // This block places a cone of gradient 2 around the trunk direction, to remove the square of initial ground + // points if (seg != -1 && sections_[seg].parent == -1) { - Eigen::Vector3d dir(0,0,1e-10); + Eigen::Vector3d dir(0, 0, 1e-10); for (auto &child : sections_[seg].children) { dir += sections_[child].tip - sections_[seg].tip; } dir.normalize(); - const double grad = 2.0; // larger cuts out a steeper (narrower) cone - Eigen::Vector3d base = sections_[seg].tip - grad*dir*radius(sections_[seg]); + const double grad = 2.0; // larger cuts out a steeper (narrower) cone + Eigen::Vector3d base = sections_[seg].tip - grad * dir * radius(sections_[seg]); Eigen::Vector3d dif = cloud.ends[i] - base; double h = dif.dot(dir); - double w = (dif - dir*h).norm(); - if (grad*w > h) + double w = (dif - dir * h).norm(); + if (grad * w > h) { colour.red = colour.green = colour.blue = 0; continue; @@ -1187,13 +1222,15 @@ bool Trees::save(const std::string &filename, const Eigen::Vector3d &offset, boo return false; } ofs << std::setprecision(4) << std::fixed; - ofs << "# Tree file. Optional per-tree attributes (e.g. 'height,crown_radius, ') followed by 'x,y,z,radius' and any additional per-segment attributes:" << std::endl; - ofs << "x,y,z,radius,parent_id,section_id"; // simple format + ofs << "# Tree file. Optional per-tree attributes (e.g. 'height,crown_radius, ') followed by 'x,y,z,radius' and any " + "additional per-segment attributes:" + << std::endl; + ofs << "x,y,z,radius,parent_id,section_id"; // simple format if (verbose) { ofs << ",weight,len,accuracy,junction_weight"; } - ofs << std::endl; + ofs << std::endl; for (size_t sec = 0; sec < sections_.size(); sec++) { const auto §ion = sections_[sec]; @@ -1201,7 +1238,8 @@ bool Trees::save(const std::string &filename, const Eigen::Vector3d &offset, boo { continue; } - ofs << section.tip[0]+offset[0] << "," << section.tip[1]+offset[1] << "," << section.tip[2]+offset[2] << "," << radius(section) << ",-1," << sec; + ofs << section.tip[0] + offset[0] << "," << section.tip[1] + offset[1] << "," << section.tip[2] + offset[2] << "," + << radius(section) << ",-1," << sec; if (verbose) { ofs << "," << section.weight << "," << section.len << "," << section.accuracy << "," << section.junction_weight; @@ -1215,7 +1253,8 @@ bool Trees::save(const std::string &filename, const Eigen::Vector3d &offset, boo { std::cout << "bad format: " << node.root << " != " << root << std::endl; } - ofs << ", " << node.tip[0]+offset[0] << "," << node.tip[1]+offset[1] << "," << node.tip[2]+offset[2] << "," << radius(node) << "," << sections_[node.parent].id << "," << contiguous_section_ids_[children[c]]; + ofs << ", " << node.tip[0] + offset[0] << "," << node.tip[1] + offset[1] << "," << node.tip[2] + offset[2] << "," + << radius(node) << "," << sections_[node.parent].id << "," << contiguous_section_ids_[children[c]]; if (verbose) { ofs << "," << node.weight << "," << node.len << "," << node.accuracy << "," << node.junction_weight; @@ -1230,4 +1269,199 @@ bool Trees::save(const std::string &filename, const Eigen::Vector3d &offset, boo return true; } +void Trees::removeSmallRadiusTrees() +{ + for (sec_ = 0; sec_ < (int)sections_.size(); sec_++) + { + std::cerr << "------- processing " << sec_ << "th section" << std::endl; + double best_accuracy = -1.0; + std::vector nodes; // used + Eigen::Vector3d base = getRootPosition(); // used + + double best_dist = 0.0; + Eigen::Vector3d best_tip; + + std::vector best_nodes; + std::vector best_ends; + std::vector radius_estimates; + std::vector accuracy_estimates; + // set random section color + for (int j = 1; j <= 10; j++) + { + double initial_dbh_height = 1; + double max_dist = initial_dbh_height + 0.1 * (double)j; // range from 0.5 to 1.5 times the specified height + nodes.clear(); + sections_[sec_].ends.clear(); + extractNodesAndEndsFromRoots(nodes, base, children, max_dist - 0.5, max_dist); + if (nodes.size() < 2) + { + continue; + } + sections_[sec_].tip = calculateTipFromVertices(nodes); + if (verbose) + { + // choose random color for the section + + for (auto &node : nodes) + { + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), points_[node].pos, 0.0, section_color); + } + } + if (removeDistantPoints(nodes)) + { + sections_[sec_].tip = calculateTipFromVertices(nodes); + } + if (verbose) + { + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), sections_[sec_].tip + Eigen::Vector3d(0, 0, 0.03), 0.0, + ray::RGBA(0, 0, 0, 255)); + } + // shift to cylinder's centre + Eigen::Vector3d up(0, 0, 1); + sections_[sec_].tip += vectorToCylinderCentre(nodes, up); + // now find the segment radius + double accuracy = 0; + double circle_coverage = 0; + Eigen::Vector3d circle_centre; + double radius = estimateCylinderRadiusUsingRANSAC(nodes, up, accuracy, circle_coverage, circle_centre); + + std::cerr << sec_ << "Estimated radius " << radius << " and accuracy " << accuracy << " and coverage " + << circle_coverage << " j " << j << std::endl; + if (accuracy > 0.4 && circle_coverage > 0.5) + { + radius_estimates.push_back(radius); + accuracy_estimates.push_back(accuracy); + ray::RGBA col(uint8_t(circle_coverage * 10), uint8_t(10 * accuracy), uint8_t(sec_), 255); + if (verbose) + { + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), sections_[sec_].tip, 0.0, col); + for (double ang = 0; ang < 2.0 * ray::kPi; ang += 0.1) + { + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), + circle_centre + radius * Eigen::Vector3d(std::sin(ang), std::cos(ang), 0), 0.0, col); + } + } + } + } + + // take average of three smallest estimates + if (radius_estimates.size() > 4) + { + double radius_average = + std::accumulate(radius_estimates.begin(), radius_estimates.end(), 0.0) / radius_estimates.size(); + double accuracy_average = + std::accumulate(accuracy_estimates.begin(), accuracy_estimates.end(), 0.0) / accuracy_estimates.size(); + + std::cerr << sec_ << "radius_average " << radius_average << std::endl; + if (radius_average < 0.05 && accuracy_average > 0.4) + { + std::cerr << "Removing tree with average radius " << radius_average << " sec " << sec_ << std::endl; + sections_[sec_].children.clear(); + continue; + } + } + } +} + +double Trees::estimateCylinderRadiusUsingRANSAC(const std::vector &nodes, const Eigen::Vector3d &dir, + double &accuracy, double &circle_coverage) +{ + double best_rad = 0.0; + double best_accuracy = 0.0; + std::vector best_inliers; + int number_of_nodes = nodes.size(); + int max_iterations = 5; + + int sample_size = std::min(number_of_nodes / 2, 10000); + + for (int i = 0; i < max_iterations; i++) + { + // Step 1: Random Sampling + std::unordered_set sample_indices_set; + while (sample_indices_set.size() < sample_size) + { + sample_indices_set.insert(rand() % nodes.size()); + } + std::vector sample_indices(sample_indices_set.begin(), sample_indices_set.end()); + + // Step 2: Estimate the Radius from Sample + Eigen::Vector3d centroid = Eigen::Vector3d::Zero(); + for (auto &idx : sample_indices) + { + Eigen::Vector3d p = points_[nodes[idx]].pos - sections_[sec_].tip; + p -= dir * p.dot(dir); // Flatten + centroid += p; + } + centroid /= sample_size; + + double estimated_radius = 0.0; + for (auto &idx : sample_indices) + { + Eigen::Vector3d offset = points_[nodes[idx]].pos - sections_[sec_].tip; + offset -= dir * offset.dot(dir); // Flatten + estimated_radius += (offset - centroid).norm(); + } + estimated_radius /= sample_size; // Average radius from the sample + + double inlier_threshold = std::min(0.05, estimated_radius * 0.5); + + // Step 3: Identify Inliers Based on Estimated Radius + std::vector inliers; + for (auto &node : nodes) + { + Eigen::Vector3d offset = points_[node].pos - sections_[sec_].tip; + offset -= dir * offset.dot(dir); // Flatten + double radial_distance = (offset - centroid).norm(); + if (std::abs(radial_distance - estimated_radius) < inlier_threshold) + { + inliers.push_back(node); + } + } + + // Step 4: Recompute Radius Using Inliers + double total_radius = 0.0; + for (auto &node : inliers) + { + Eigen::Vector3d offset = points_[node].pos - sections_[sec_].tip; + offset -= dir * offset.dot(dir); // Flatten + total_radius += (offset - centroid).norm(); + } + double avg_radius = total_radius / inliers.size(); + + // Step 5: Update Circle Coverage, simulate points along the cirlce and see + // how many matching poitns + double coverage = 0.0; + double circle_coverage_th = 0.05; + for (int j = 0; j < 1000; j++) + { + double angle = 2.0 * M_PI * j / 1000; + Eigen::Vector3d p = centroid + avg_radius * Eigen::Vector3d(cos(angle), sin(angle), 0); + for (auto &node : inliers) + { + Eigen::Vector3d offset = points_[node].pos - sections_[sec_].tip; + offset -= dir * offset.dot(dir); // Flatten + if ((offset - p).norm() < circle_coverage_th) + { + coverage += 1.0; + break; + } + } + } + coverage /= 1000.0; + + // Step 6: Update Best Model if Necessary + if (best_accuracy < static_cast(inliers.size()) / nodes.size() && coverage > circle_coverage) + { + best_inliers = inliers; + best_rad = avg_radius; + best_accuracy = static_cast(inliers.size()) / nodes.size(); + circle_coverage = coverage; + } + } + + accuracy = best_accuracy; + + return best_rad; +} + } // namespace ray diff --git a/raylib/extraction/raytrees.h b/raylib/extraction/raytrees.h index 523101e4..e7f80859 100644 --- a/raylib/extraction/raytrees.h +++ b/raylib/extraction/raytrees.h @@ -18,20 +18,23 @@ namespace ray struct RAYLIB_EXPORT TreesParams { TreesParams(); - double max_diameter; // maximum tree diameter. Trees wider than this may be segmented into multiple trees - double crop_length; // distance to end of branch where it crops the branch, and doesn't generate further geometry - double distance_limit; // maximum distance between points that can count as connected - double height_min; // minimum height for a tree. Lower values are considered undergrowth and excluded + double max_diameter; // maximum tree diameter. Trees wider than this may be segmented into multiple trees + double crop_length; // distance to end of branch where it crops the branch, and doesn't generate further geometry + double distance_limit; // maximum distance between points that can count as connected + double height_min; // minimum height for a tree. Lower values are considered undergrowth and excluded double girth_height_ratio; // how far up tree to measure girth double cylinder_length_to_width; // the slenderness of the branch segment cylinders double gap_ratio; // max gap per branch length double span_ratio; // points that span a larger width determine that a branch has become two - double gravity_factor; // preferences branches that are less lateral, so penalises implausable horizontal branches - double grid_width; // used on a grid cell with overlap, to remove trees with a base in the overlap zone - bool segment_branches; // flag to output the ray cloud coloured by branch segment index rather than by tree index - double global_taper; // forced global taper, uses global_taper_factor to define how much it is applied - double global_taper_factor; // 0 estimates per-tree tapering, 1 uses per-scan tapering, 0.5 is mid-way on mid-weight trees - bool use_rays; // use the full rays in order to estimate a smaller radius when points are not all on the real branch + double gravity_factor; // preferences branches that are less lateral, so penalises implausable horizontal branches + double grid_width; // used on a grid cell with overlap, to remove trees with a base in the overlap zone + bool segment_branches; // flag to output the ray cloud coloured by branch segment index rather than by tree index + double global_taper; // forced global taper, uses global_taper_factor to define how much it is applied + double + global_taper_factor; // 0 estimates per-tree tapering, 1 uses per-scan tapering, 0.5 is mid-way on mid-weight trees + bool use_rays; // use the full rays in order to estimate a smaller radius when points are not all on the real branch + Eigen::Vector3d grid_origin; // if using a grid, then this is the left corner of the interest zone + double min_radius; // minimum radius for a branch }; struct BranchSection; // forwards declaration @@ -66,10 +69,11 @@ class RAYLIB_EXPORT Trees void extractNodesAndEndsFromRoots(std::vector &nodes, const Eigen::Vector3d &base, const std::vector> &children, double min_dist, double max_dist); /// find separate clusters of points within the branch section - std::vector> findPointClusters(const Eigen::Vector3d &base, bool &points_removed, - double thickness, double span, double gap); + std::vector> findPointClusters(const Eigen::Vector3d &base, bool &points_removed, double thickness, + double span, double gap); /// split the branch section to one branch for each cluster - void bifurcate(const std::vector> &clusters, double thickness, std::vector> &children, bool clip_tree, bool add_offshoots); + void bifurcate(const std::vector> &clusters, double thickness, + std::vector> &children, bool clip_tree, bool add_offshoots); /// find the points within the branch section from its end points void extractNodesFromEnds(std::vector &nodes); /// set the branch section tip position from the supplied list of Vertex IDs @@ -87,7 +91,8 @@ class RAYLIB_EXPORT Trees /// set ids that are locel (0-based) per tree void generateLocalSectionIds(); /// if using an overlapping grid, then remove trees with base outside the non-overlapping cell bounds - void removeOutOfBoundSections(const Cloud &cloud, Eigen::Vector3d &min_bound, Eigen::Vector3d &max_bound, const Eigen::Vector3d &offset); + void removeOutOfBoundSections(const Cloud &cloud, Eigen::Vector3d &min_bound, Eigen::Vector3d &max_bound, + const Eigen::Vector3d &offset); /// colour the cloud based on the section id for each point void segmentCloud(Cloud &cloud, std::vector &root_segs, const std::vector §ion_ids); /// remove points from the ray cloud if outside of the non-overlapping grid cell bounds @@ -99,16 +104,21 @@ class RAYLIB_EXPORT Trees double radius(const BranchSection §ion) const; /// remove elements of nodes that are too distant to the set of end points bool removeDistantPoints(std::vector &nodes); + /// remove trees with radius less than the param + void removeSmallRadiusTrees(); + /// Esimate cylinder radius using RANSAC + double estimateCylinderRadiusUsingRANSAC(const std::vector &nodes, const Eigen::Vector3d &dir, double &accuracy, + double &circle_coverage); // cached data that is used throughout the processing method int sec_; const TreesParams *params_; std::vector points_; - double forest_taper_{0}; - double forest_weight_{0}; - double forest_weight_squared_{0}; - std::vector contiguous_section_ids_; // converts section ids to contiguous (empty trees removed) sectino ids for output - + double forest_taper_{ 0 }; + double forest_weight_{ 0 }; + double forest_weight_squared_{ 0 }; + std::vector + contiguous_section_ids_; // converts section ids to contiguous (empty trees removed) sectino ids for output }; /// The structure for a single (cylindrical) branch section @@ -139,8 +149,8 @@ struct RAYLIB_EXPORT BranchSection int root; int id; // 0 based per tree double max_distance_to_end; - double total_taper; // a weighted taper estimate - double total_weight; // the weighting + double total_taper; // a weighted taper estimate + double total_weight; // the weighting double len, accuracy, junction_weight; double radius_scale; int split_count; diff --git a/raylib/raylaz.cpp b/raylib/raylaz.cpp index 46fa9e68..6370bf74 100644 --- a/raylib/raylaz.cpp +++ b/raylib/raylaz.cpp @@ -23,6 +23,16 @@ bool readLas(const std::string &file_name, { #if RAYLIB_WITH_LAS std::cout << "readLas: filename: " << file_name << std::endl; + std::ifstream file(file_name); + + if (file) + { + std::cout << "File exists!" << std::endl; + } + else + { + std::cout << "File does not exist." << std::endl; + } std::ifstream ifs; ifs.open(file_name.c_str(), std::ios::in | std::ios::binary); diff --git a/raylib/rayparse.cpp b/raylib/rayparse.cpp index 42d7d31e..c72a04f5 100644 --- a/raylib/rayparse.cpp +++ b/raylib/rayparse.cpp @@ -240,7 +240,7 @@ bool Vector3dArgument::parse(int argc, char *argv[], int &index, bool set_value) if (val < min_value_ || val > max_value_) { std::cout << "Please set argument " << index << " within the range: " << min_value_ << " to " << max_value_ - << std::endl; + << " provided " << val << std::endl; return false; } value_[i] = val; @@ -412,4 +412,28 @@ bool OptionalKeyValueArgument::parse(int argc, char *argv[], int &index, bool se return false; } +bool Optional3DVectorArgument::parse(int argc, char *argv[], int &index, bool set_value) +{ + if (index >= argc) + return false; + + std::string str(argv[index]); + if (str == ("--" + name_)) + { + if (set_value) + is_set_ = true; + + index++; + if (!value_->parse(argc, argv, index, set_value)) + { + index--; + return false; + } + + return true; + } + + return false; +} + } // namespace ray From 10fcc4a2e49e6d8ad202bdd7882b01deb4bec59c Mon Sep 17 00:00:00 2001 From: GabijaBern Date: Wed, 20 Nov 2024 10:00:25 +0000 Subject: [PATCH 2/2] Copying Tom's ideas --- raycloudtools/rayextract/rayextract.cpp | 77 +++++++++-------- raylib/extraction/raytrees.cpp | 108 +++++++++++++----------- raylib/extraction/raytrees.h | 9 +- raylib/rayparse.cpp | 28 +----- 4 files changed, 111 insertions(+), 111 deletions(-) diff --git a/raycloudtools/rayextract/rayextract.cpp b/raycloudtools/rayextract/rayextract.cpp index d6c82888..c8fe9bb9 100644 --- a/raycloudtools/rayextract/rayextract.cpp +++ b/raycloudtools/rayextract/rayextract.cpp @@ -5,10 +5,10 @@ // Author: Thomas Lowe #include "raylib/extraction/rayclusters.h" #include "raylib/extraction/rayforest.h" +#include "raylib/extraction/rayleaves.h" #include "raylib/extraction/rayterrain.h" #include "raylib/extraction/raytrees.h" #include "raylib/extraction/raytrunks.h" -#include "raylib/extraction/rayleaves.h" #include "raylib/raycloud.h" #include "raylib/rayforestgen.h" #include "raylib/rayforeststructure.h" @@ -25,7 +25,8 @@ static std::string extract_type; void usage(int exit_code = 1) { - const bool none = extract_type != "terrain" && extract_type != "trunks" && extract_type != "forest" && extract_type != "trees" && extract_type != "leaves"; + const bool none = extract_type != "terrain" && extract_type != "trunks" && extract_type != "forest" && + extract_type != "trees" && extract_type != "leaves"; // clang-format off std::cout << "Extract natural features into a text file or mesh file" << std::endl; std::cout << "usage:" << std::endl; @@ -54,21 +55,22 @@ void usage(int exit_code = 1) std::cout << " --max_diameter 0.9 - (-m) maximum trunk diameter in segmenting trees" << std::endl; std::cout << " --crop_length 1.0 - (-p) crops small branches to this distance from end" << std::endl; std::cout << " --distance_limit 1 - (-d) maximum distance between neighbour points in a tree" << std::endl; - std::cout << " --height_min 2 - (-h) minimum height counted as a tree" << std::endl; - std::cout << " --min_radius 0.01 - (-r) minimum tree radius to consider" << std::endl; + std::cout << " --height_min 2 - (-h) minimum height tree to reconstruct" << std::endl; + std::cout << " --radius_min 0 - (-r) minimum radius tree to reconstruct" << std::endl; std::cout << " --girth_height_ratio 0.12 - (-i) the amount up tree's height to estimate trunk girth" << std::endl; std::cout << " --global_taper 0.024 - (-a) force a taper value (diameter per length) for trees under global_taper_factor of max tree height. Use 0 to estimate global taper from the data" << std::endl; std::cout << " --global_taper_factor 0.3- (-o) 1 estimates same taper for whole scan, 0 is per-tree tapering. Like a soft cutoff at this amount of max tree height" << std::endl; std::cout << " --gravity_factor 0.3 - (-f) larger values preference vertical trees" << std::endl; std::cout << " --branch_segmentation- (-b) _segmented.ply is per branch segment" << std::endl; std::cout << " --grid_width 10 - (-w) crops results assuming cloud has been gridded with given width" << std::endl; - std::cout << " --grid_origin 0,0,0 - location of origin within grid cell that overlaps it. Defaults to a cell-centre origin (at grid_width/2 in each axis) matching raysplit grid. 0,0,0 is for a vertex origin." << std::endl; std::cout << " --use_rays - (-u) use rays to reduce trunk radius overestimation in noisy cloud data" << std::endl; - std::cout << " (for internal constants -c -g -s see source file rayextract)" << std::endl; + std::cout << " (for internal constants -c -g -s -d see source file rayextract)" << std::endl; // These are the internal parameters that I don't expose as they are 'advanced' only, you shouldn't need to adjust them // std::cout << " --cylinder_length_to_width 4- (-c) how slender the cylinders are" << std::endl; // std::cout << " --gap_ratio 0.016 - (-g) will split for lateral gaps at this multiple of branch length" << std::endl; // std::cout << " --span_ratio 4.5 - (-s) will split when branch width spans this multiple of radius" << std::endl; + // std::cout << " --grid_origin 0,0 - (-d) location of grid corner (any of them) when grid_width used, use 0,0 for grid with vertex at 0,0. + // Default is -grid_width/2,-grid_width/2 to match the grid in raysplit grid" << std::endl; } if (extract_type == "leaves" || none) { @@ -94,33 +96,33 @@ int rayExtract(int argc, char *argv[]) } ray::FileArgument cloud_file, mesh_file, trunks_file, trees_file, leaf_file; ray::TextArgument forest("forest"), trees("trees"), trunks("trunks"), terrain("terrain"), leaves("leaves"); - - ray::Vector3dArgument grid_origin; - ray::OptionalKeyValueArgument grid_origin_option("grid_origin", &grid_origin); ray::OptionalKeyValueArgument groundmesh_option("ground", 'g', &mesh_file); ray::OptionalKeyValueArgument trunks_option("trunks", 't', &trunks_file); ray::DoubleArgument gradient(0.001, 1000.0, 1.0), global_taper(0.0, 1.0), global_taper_factor(0.0, 1.0); ray::OptionalKeyValueArgument gradient_option("gradient", 'g', &gradient); - ray::OptionalFlagArgument exclude_rays("exclude_rays", 'e'), segment_branches("branch_segmentation", 'b'), stalks("stalks", 's'), use_rays("use_rays", 'u'); + ray::OptionalFlagArgument exclude_rays("exclude_rays", 'e'), segment_branches("branch_segmentation", 'b'), + stalks("stalks", 's'), use_rays("use_rays", 'u'); ray::DoubleArgument width(0.01, 10.0, 0.25), drop(0.001, 1.0), max_gradient(0.01, 5.0), min_gradient(0.01, 5.0); ray::DoubleArgument max_diameter(0.01, 100.0), distance_limit(0.01, 10.0), height_min(0.01, 1000.0), - min_diameter(0.01, 100.0), leaf_area(0.00001, 1.0, 0.002), leaf_droop(0.0, 10.0, 0.1), crop_length(0.01, 100.0);; - ray::DoubleArgument girth_height_ratio(0.001, 0.5), length_to_radius(0.01, 10000.0), cylinder_length_to_width(0.1, 20.0), gap_ratio(0.01, 10.0), - span_ratio(0.01, 10.0), min_radius(0.01, 100.0); - ray::DoubleArgument gravity_factor(0.0, 100.0), grid_width(1.0, 100000.0), - grid_overlap(0.0, 0.9); + radius_min(0.0, 1000.0), min_diameter(0.01, 100.0), leaf_area(0.00001, 1.0, 0.002), leaf_droop(0.0, 10.0, 0.1), + crop_length(0.01, 100.0); + ; + ray::DoubleArgument girth_height_ratio(0.001, 0.5), length_to_radius(0.01, 10000.0), + cylinder_length_to_width(0.1, 20.0), gap_ratio(0.01, 10.0), span_ratio(0.01, 10.0); + ray::DoubleArgument gravity_factor(0.0, 100.0), grid_width(1.0, 100000.0), grid_overlap(0.0, 0.9); + ray::Vector2dArgument grid_origin(-1e10, 1e10); ray::OptionalKeyValueArgument max_diameter_option("max_diameter", 'm', &max_diameter); ray::OptionalKeyValueArgument crop_length_option("crop_length", 'n', &crop_length); ray::OptionalKeyValueArgument distance_limit_option("distance_limit", 'd', &distance_limit); ray::OptionalKeyValueArgument height_min_option("height_min", 'h', &height_min); + ray::OptionalKeyValueArgument radius_min_option("radius_min", 'r', &radius_min); ray::OptionalKeyValueArgument girth_height_ratio_option("girth_height_ratio", 'i', &girth_height_ratio); ray::OptionalKeyValueArgument cylinder_length_to_width_option("cylinder_length_to_width", 'c', &cylinder_length_to_width); - - ray::OptionalKeyValueArgument min_radius_option("min_radius", 'r', &min_radius); ray::OptionalKeyValueArgument gap_ratio_option("gap_ratio", 'g', &gap_ratio); ray::OptionalKeyValueArgument span_ratio_option("span_ratio", 's', &span_ratio); + ray::OptionalKeyValueArgument grid_origin_option("grid_origin", 'd', &grid_origin); ray::OptionalKeyValueArgument gravity_factor_option("gravity_factor", 'f', &gravity_factor); ray::OptionalKeyValueArgument grid_width_option("grid_width", 'w', &grid_width); ray::OptionalKeyValueArgument global_taper_option("global_taper", 'a', &global_taper); @@ -142,11 +144,12 @@ int rayExtract(int argc, char *argv[]) { &groundmesh_option, &trunks_option, &width_option, &smooth_option, &drop_option, &verbose }); bool extract_trees = ray::parseCommandLine( argc, argv, { &trees, &cloud_file, &mesh_file }, - { &max_diameter_option, &distance_limit_option, &height_min_option, &crop_length_option, &girth_height_ratio_option, - &cylinder_length_to_width_option, &gap_ratio_option, &span_ratio_option, &gravity_factor_option, - &segment_branches, &grid_width_option, &global_taper_option, &global_taper_factor_option, &use_rays, &verbose, - &grid_origin_option, &min_radius_option }); - bool extract_leaves = ray::parseCommandLine(argc, argv, { &leaves, &cloud_file, &trees_file }, { &leaf_option, &leaf_area_option, &leaf_droop_option, &stalks }); + { &max_diameter_option, &distance_limit_option, &height_min_option, &radius_min_option, &crop_length_option, + &girth_height_ratio_option, &cylinder_length_to_width_option, &gap_ratio_option, &span_ratio_option, + &grid_origin_option, &gravity_factor_option, &segment_branches, &grid_width_option, &global_taper_option, + &global_taper_factor_option, &use_rays, &verbose }); + bool extract_leaves = ray::parseCommandLine(argc, argv, { &leaves, &cloud_file, &trees_file }, + { &leaf_option, &leaf_area_option, &leaf_droop_option, &stalks }); if (!extract_trunks && !extract_forest && !extract_terrain && !extract_trees && !extract_leaves) @@ -199,6 +202,10 @@ int rayExtract(int argc, char *argv[]) { params.height_min = height_min.value(); } + if (radius_min_option.isSet()) + { + params.radius_min = radius_min.value(); + } if (crop_length_option.isSet()) { params.crop_length = crop_length.value(); @@ -226,6 +233,12 @@ int rayExtract(int argc, char *argv[]) if (grid_width_option.isSet()) { params.grid_width = grid_width.value(); + params.grid_origin = + -Eigen::Vector2d(grid_width.value(), grid_width.value()) / 2.0; // the centred grid used by raysplit grid + } + if (grid_origin_option.isSet()) + { + params.grid_origin = grid_origin.value(); } if (global_taper_option.isSet()) { @@ -235,14 +248,6 @@ int rayExtract(int argc, char *argv[]) { params.global_taper_factor = global_taper_factor.value(); } - if (grid_origin_option.isSet()) - { - params.grid_origin = grid_origin.value(); - } - if (min_radius_option.isSet()) - { - params.min_radius = min_radius.value(); - } params.use_rays = use_rays.isSet(); params.segment_branches = segment_branches.isSet(); @@ -259,11 +264,15 @@ int rayExtract(int argc, char *argv[]) ray::ForestStructure forest; if (!forest.load(cloud_file.nameStub() + "_trees.txt")) { - usage(); + std::cerr << "Unable to load " + << cloud_file.nameStub() + + "_trees.txt to generate tree mesh file, this could mean that there were no trees output" + << std::endl; + exit(true); } ray::Mesh tree_mesh; forest.generateSmoothMesh(tree_mesh, -1, 1, 1, 1); - ray::writePlyMesh(cloud_file.nameStub() + "_trees_mesh.ply", tree_mesh, true); + ray::writePlyMesh(cloud_file.nameStub() + "_trees_mesh.ply", tree_mesh, true); } // extract the tree locations from a larger, aerial view of a forest else if (extract_forest) @@ -321,8 +330,8 @@ int rayExtract(int argc, char *argv[]) } else if (extract_leaves) { - ray::generateLeaves(cloud_file.nameStub(), trees_file.name(), leaf_file.name(), - leaf_area.value(), leaf_droop.value(), stalks.isSet()); + ray::generateLeaves(cloud_file.nameStub(), trees_file.name(), leaf_file.name(), leaf_area.value(), + leaf_droop.value(), stalks.isSet()); } else { diff --git a/raylib/extraction/raytrees.cpp b/raylib/extraction/raytrees.cpp index 94b8485a..934b15c2 100644 --- a/raylib/extraction/raytrees.cpp +++ b/raylib/extraction/raytrees.cpp @@ -7,6 +7,8 @@ #include #include "rayclusters.h" +#include + namespace ray { TreesParams::TreesParams() @@ -14,6 +16,7 @@ TreesParams::TreesParams() , crop_length(1.0) , distance_limit(1.0) , height_min(2.0) + , radius_min(0.0) , girth_height_ratio(0.12) , cylinder_length_to_width(4.0) , gap_ratio(0.016) @@ -23,7 +26,7 @@ TreesParams::TreesParams() , segment_branches(false) , global_taper(0.012) , global_taper_factor(0.3) - , grid_origin(0, 0, 0) + , grid_origin(0, 0) {} /// The main reconstruction algorithm @@ -148,13 +151,6 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons continue; } - if (params_->min_radius && estimated_radius < params_->min_radius) - { - sections_[sec_].tip = base + Eigen::Vector3d(0, 0, 0.01); - sections_[sec_].total_weight = 1e-10; - sections_[sec_].ends.clear(); - continue; - } sections_[sec_].tip = best_tip; sections_[sec_].ends = best_ends; nodes = best_nodes; @@ -300,9 +296,9 @@ Trees::Trees(Cloud &cloud, const Eigen::Vector3d &offset, const Mesh &mesh, cons removeOutOfBoundSections(cloud, min_bound, max_bound, offset); } - if (params_->min_radius) + if (params_->radius_min) { - removeSmallRadiusTrees(); + removeSmallRadiusTrees(verbose, offset, children); } std::vector root_segs(cloud.ends.size(), -1); @@ -1082,14 +1078,15 @@ void Trees::removeOutOfBoundSections(const Cloud &cloud, Eigen::Vector3d &min_bo { const double width = params_->grid_width; cloud.calcBounds(&min_bound, &max_bound); - const Eigen::Vector3d mid = (min_bound + max_bound) / 2.0 + offset + params_->grid_origin; - const Eigen::Vector2d inds(std::round(mid[0] / width), std::round(mid[1] / width)); - min_bound[0] = width * (inds[0] - 0.5) - offset[0]; - min_bound[1] = width * (inds[1] - 0.5) - offset[1]; - max_bound[0] = width * (inds[0] + 0.5) - offset[0]; - max_bound[1] = width * (inds[1] + 0.5) - offset[1]; - std::cerr << "min bound: " << (min_bound + offset).transpose() << ", max bound: " << (max_bound + offset).transpose() - << " offset " << offset << std::endl; + + const Eigen::Vector3d mid = (min_bound + max_bound) / 2.0 + offset; + Eigen::Vector2d mid_local(mid[0], mid[1]); + mid_local += params_->grid_origin; + const Eigen::Vector2d inds(std::floor(mid_local[0] / width), std::floor(mid_local[1] / width)); + min_bound[0] = width * (inds[0]) - offset[0]; + min_bound[1] = width * (inds[1]) - offset[1]; + max_bound[0] = width * (inds[0] + 1.0) - offset[0]; + max_bound[1] = width * (inds[1] + 1.0) - offset[1]; // disable trees out of bounds for (auto §ion : sections_) @@ -1269,14 +1266,23 @@ bool Trees::save(const std::string &filename, const Eigen::Vector3d &offset, boo return true; } -void Trees::removeSmallRadiusTrees() +void Trees::removeSmallRadiusTrees(bool verbose, const Eigen::Vector3d &offset, + const std::vector> &children) { + ray::Cloud debug_cloud; + double coverage_threshold = 0.7; + double accuracy_threshold = 0.5; for (sec_ = 0; sec_ < (int)sections_.size(); sec_++) { - std::cerr << "------- processing " << sec_ << "th section" << std::endl; + if (sections_[sec_].parent >= 0 || sections_[sec_].children.empty()) + { + continue; + } + if (verbose) + std::cerr << "------- Processing " << sec_ << "th section" << std::endl; double best_accuracy = -1.0; - std::vector nodes; // used - Eigen::Vector3d base = getRootPosition(); // used + std::vector nodes; + Eigen::Vector3d base = getRootPosition(); double best_dist = 0.0; Eigen::Vector3d best_tip; @@ -1284,12 +1290,12 @@ void Trees::removeSmallRadiusTrees() std::vector best_nodes; std::vector best_ends; std::vector radius_estimates; - std::vector accuracy_estimates; + std::vector coverage_estimates; // set random section color for (int j = 1; j <= 10; j++) { double initial_dbh_height = 1; - double max_dist = initial_dbh_height + 0.1 * (double)j; // range from 0.5 to 1.5 times the specified height + double max_dist = initial_dbh_height + 0.1 * (double)j; nodes.clear(); sections_[sec_].ends.clear(); extractNodesAndEndsFromRoots(nodes, base, children, max_dist - 0.5, max_dist); @@ -1301,7 +1307,7 @@ void Trees::removeSmallRadiusTrees() if (verbose) { // choose random color for the section - + ray::RGBA section_color(0, 0, 0, 255); for (auto &node : nodes) { debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), points_[node].pos, 0.0, section_color); @@ -1311,11 +1317,7 @@ void Trees::removeSmallRadiusTrees() { sections_[sec_].tip = calculateTipFromVertices(nodes); } - if (verbose) - { - debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), sections_[sec_].tip + Eigen::Vector3d(0, 0, 0.03), 0.0, - ray::RGBA(0, 0, 0, 255)); - } + // shift to cylinder's centre Eigen::Vector3d up(0, 0, 1); sections_[sec_].tip += vectorToCylinderCentre(nodes, up); @@ -1324,17 +1326,19 @@ void Trees::removeSmallRadiusTrees() double circle_coverage = 0; Eigen::Vector3d circle_centre; double radius = estimateCylinderRadiusUsingRANSAC(nodes, up, accuracy, circle_coverage, circle_centre); - - std::cerr << sec_ << "Estimated radius " << radius << " and accuracy " << accuracy << " and coverage " - << circle_coverage << " j " << j << std::endl; - if (accuracy > 0.4 && circle_coverage > 0.5) + if (verbose) + { + std::cerr << sec_ << ": estimated radius " << radius << " accuracy " << accuracy << " and coverage " + << circle_coverage << std::endl; + } + if (accuracy > accuracy_threshold && circle_coverage > coverage_threshold) { radius_estimates.push_back(radius); - accuracy_estimates.push_back(accuracy); + coverage_estimates.push_back(accuracy); ray::RGBA col(uint8_t(circle_coverage * 10), uint8_t(10 * accuracy), uint8_t(sec_), 255); if (verbose) { - debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), sections_[sec_].tip, 0.0, col); + debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), circle_centre, 0.0, col); for (double ang = 0; ang < 2.0 * ray::kPi; ang += 0.1) { debug_cloud.addRay(Eigen::Vector3d(0, 0, 0), @@ -1344,27 +1348,35 @@ void Trees::removeSmallRadiusTrees() } } - // take average of three smallest estimates - if (radius_estimates.size() > 4) + // take average of estimates + if (radius_estimates.size() > 3) { double radius_average = std::accumulate(radius_estimates.begin(), radius_estimates.end(), 0.0) / radius_estimates.size(); - double accuracy_average = - std::accumulate(accuracy_estimates.begin(), accuracy_estimates.end(), 0.0) / accuracy_estimates.size(); - - std::cerr << sec_ << "radius_average " << radius_average << std::endl; - if (radius_average < 0.05 && accuracy_average > 0.4) + double coverage_average = + std::accumulate(coverage_estimates.begin(), coverage_estimates.end(), 0.0) / coverage_estimates.size(); + if (verbose) + std::cerr << sec_ << "radius_average " << radius_average << " coverage_average " << coverage_average + << std::endl; + if (radius_average < params_->radius_min) { - std::cerr << "Removing tree with average radius " << radius_average << " sec " << sec_ << std::endl; + if (verbose) + std::cerr << "Removing tree with average radius " << radius_average << " sec " << sec_ << std::endl; sections_[sec_].children.clear(); continue; } } } + if (verbose) + { + debug_cloud.translate(offset); + debug_cloud.save("../data/debug.ply"); + } } double Trees::estimateCylinderRadiusUsingRANSAC(const std::vector &nodes, const Eigen::Vector3d &dir, - double &accuracy, double &circle_coverage) + double &accuracy, double &circle_coverage, + Eigen::Vector3d &circle_centre) { double best_rad = 0.0; double best_accuracy = 0.0; @@ -1403,7 +1415,7 @@ double Trees::estimateCylinderRadiusUsingRANSAC(const std::vector &nodes, c } estimated_radius /= sample_size; // Average radius from the sample - double inlier_threshold = std::min(0.05, estimated_radius * 0.5); + double inlier_threshold = std::min(0.03, estimated_radius * 0.6); // Step 3: Identify Inliers Based on Estimated Radius std::vector inliers; @@ -1431,7 +1443,7 @@ double Trees::estimateCylinderRadiusUsingRANSAC(const std::vector &nodes, c // Step 5: Update Circle Coverage, simulate points along the cirlce and see // how many matching poitns double coverage = 0.0; - double circle_coverage_th = 0.05; + double circle_coverage_th = std::min(0.03, estimated_radius * 0.2); for (int j = 0; j < 1000; j++) { double angle = 2.0 * M_PI * j / 1000; @@ -1456,6 +1468,8 @@ double Trees::estimateCylinderRadiusUsingRANSAC(const std::vector &nodes, c best_rad = avg_radius; best_accuracy = static_cast(inliers.size()) / nodes.size(); circle_coverage = coverage; + // save centroid in original space + circle_centre = sections_[sec_].tip + centroid; } } diff --git a/raylib/extraction/raytrees.h b/raylib/extraction/raytrees.h index e7f80859..2ab66db2 100644 --- a/raylib/extraction/raytrees.h +++ b/raylib/extraction/raytrees.h @@ -33,8 +33,8 @@ struct RAYLIB_EXPORT TreesParams double global_taper_factor; // 0 estimates per-tree tapering, 1 uses per-scan tapering, 0.5 is mid-way on mid-weight trees bool use_rays; // use the full rays in order to estimate a smaller radius when points are not all on the real branch - Eigen::Vector3d grid_origin; // if using a grid, then this is the left corner of the interest zone - double min_radius; // minimum radius for a branch + Eigen::Vector2d grid_origin; // if using a grid, then this is the left corner of the interest zone + double radius_min; // minimum radius for a branch }; struct BranchSection; // forwards declaration @@ -105,10 +105,11 @@ class RAYLIB_EXPORT Trees /// remove elements of nodes that are too distant to the set of end points bool removeDistantPoints(std::vector &nodes); /// remove trees with radius less than the param - void removeSmallRadiusTrees(); + void removeSmallRadiusTrees(bool verbose, const Eigen::Vector3d &offset, + const std::vector> &children); /// Esimate cylinder radius using RANSAC double estimateCylinderRadiusUsingRANSAC(const std::vector &nodes, const Eigen::Vector3d &dir, double &accuracy, - double &circle_coverage); + double &circle_coverage, Eigen::Vector3d &circle_centre); // cached data that is used throughout the processing method int sec_; diff --git a/raylib/rayparse.cpp b/raylib/rayparse.cpp index c72a04f5..77ecb97e 100644 --- a/raylib/rayparse.cpp +++ b/raylib/rayparse.cpp @@ -80,8 +80,8 @@ bool FileArgument::parse(int argc, char *argv[], int &index, bool set_value) std::string file = std::string(argv[index]); if (file.length() < 1) return false; - if (file[0] == '-') // no file should start with a dash. That is reserved for flag arguments - return false; + if (file[0] == '-') // no file should start with a dash. That is reserved for flag arguments + return false; if (check_extension_) { if (file.length() <= 2) @@ -412,28 +412,4 @@ bool OptionalKeyValueArgument::parse(int argc, char *argv[], int &index, bool se return false; } -bool Optional3DVectorArgument::parse(int argc, char *argv[], int &index, bool set_value) -{ - if (index >= argc) - return false; - - std::string str(argv[index]); - if (str == ("--" + name_)) - { - if (set_value) - is_set_ = true; - - index++; - if (!value_->parse(argc, argv, index, set_value)) - { - index--; - return false; - } - - return true; - } - - return false; -} - } // namespace ray