Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Update radius #25

Merged
merged 6 commits into from
Oct 16, 2024
Merged
Show file tree
Hide file tree
Changes from 5 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -88,9 +88,9 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
# 0 - ajout de dimensions temporaires et de sortie
temporary_dimensions = [
"PT_VEG_DSM",
"PT_ON_BRIDGE",
"PT_ON_BUILDING",
"PT_ON_VEGET",
"PT_UNDER_BRIDGE",
"PT_CLOSED_BUILDING",
"PT_UNDER_VEGET",
"PT_ON_SOL",
"PT_ON_VIRT",
]
Expand Down Expand Up @@ -128,28 +128,28 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
condition_ref="Classification==2 && PT_VEG_DSM==0",
condition_out="PT_VEG_DSM=0",
)
# 1.3 Isolement en PT_ON_VEGET=1 des éléments sous la végétation (hors sol)
# 1.3 Isolement en PT_UNDER_VEGET=1 des éléments sous la végétation (hors sol)
pipeline = macro.add_radius_assign(
pipeline,
1,
False,
condition_src=macro.build_condition("Classification", [6, 9, 17, 67]),
condition_ref=macro.build_condition("Classification", [4, 5]),
condition_out="PT_ON_VEGET=1",
max2d_above=0,
max2d_below=-1,
condition_out="PT_UNDER_VEGET=1",
max2d_above=-1,
max2d_below=0,
)
pipeline = macro.add_radius_assign(
pipeline,
1,
False,
condition_src="PT_ON_VEGET==1 && ( "
condition_src="PT_UNDER_VEGET==1 && ( "
+ macro.build_condition("Classification", [6, 9, 17, 67])
+ " )",
condition_ref="PT_ON_VEGET==0 && ( "
condition_ref="PT_UNDER_VEGET==0 && ( "
+ macro.build_condition("Classification", [6, 9, 17, 67])
+ " )",
condition_out="PT_ON_VEGET=0",
condition_out="PT_UNDER_VEGET=0",
max2d_above=0.5,
max2d_below=0.5,
)
Expand Down Expand Up @@ -182,8 +182,8 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
condition_src="Classification==9",
condition_ref="Classification==2",
condition_out="PT_ON_SOL=1",
max2d_above=0,
max2d_below=-1,
max2d_above=-1,
max2d_below=0,
)
pipeline = macro.add_radius_assign(
pipeline,
Expand Down Expand Up @@ -232,7 +232,7 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
resolution=0.5,
output_dimension=dsm_dimension,
output_type="max",
where="(PT_ON_VEGET==0 && ("
where="(PT_UNDER_VEGET==0 && ("
+ macro.build_condition("Classification", [6, 17, 67])
+ f") || {dsm_dimension}==1)",
)
Expand All @@ -247,7 +247,7 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
resolution=0.5,
output_dimension=dsm_dimension,
output_type="min",
where="(PT_ON_VEGET==0 && PT_ON_SOL==0 && PT_ON_VIRT==0 && Classification==9)",
where="(PT_UNDER_VEGET==0 && PT_ON_SOL==0 && PT_ON_VIRT==0 && Classification==9)",
)
###################################################################################################################
# 4 - Gestion des points sol sous la veget,bâtis et ponts pour le MNS
Expand All @@ -271,14 +271,14 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
False,
condition_src="Classification==2 && PT_VEG_DSM==0",
condition_ref=macro.build_condition("Classification", [6, 67]),
condition_out="PT_ON_BUILDING=1",
condition_out="PT_CLOSED_BUILDING=1",
)
pipeline = macro.add_radius_assign(
pipeline,
1,
False,
condition_src=f"Classification==2 && {dsm_dimension}==0 && PT_ON_BUILDING==1 && {dtm_dimension}==1",
condition_ref="Classification==2 && PT_ON_BUILDING==0 && PT_VEG_DSM==0",
condition_src=f"Classification==2 && {dsm_dimension}==0 && PT_CLOSED_BUILDING==1 && {dtm_dimension}==1",
condition_ref="Classification==2 && PT_CLOSED_BUILDING==0 && PT_VEG_DSM==0",
condition_out=f"{dsm_dimension}=1",
)
###################################################################################################################
Expand All @@ -291,29 +291,29 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
False,
condition_src=macro.build_condition("Classification", [2, 3, 4, 5, 6, 9, 67]),
condition_ref="Classification==17",
condition_out="PT_ON_BRIDGE=1",
max2d_above=0, # ne pas prendre les points qui sont au dessus des points pont (condition_ref)
max2d_below=-1, # prendre tous les points qui sont en dessous des points pont (condition_ref)
condition_out="PT_UNDER_BRIDGE=1",
max2d_above=-1, # prendre les points (condition_src) qui on des points ponts au dessus d'eux (condition_ref)
max2d_below=0,
)
pipeline = macro.add_radius_assign(
pipeline,
1.25,
False,
condition_src="PT_ON_BRIDGE==1",
condition_ref="PT_ON_BRIDGE==0 && "
condition_src="PT_UNDER_BRIDGE==1",
condition_ref="PT_UNDER_BRIDGE==0 && "
+ macro.build_condition("Classification", [2, 3, 4, 5, 6, 9, 67]),
condition_out="PT_ON_BRIDGE=0",
max2d_above=0.5, # ne pas prendre les points qui sont au dessus des points pont (condition_ref)
max2d_below=0.5, # prendre tous les points qui sont en dessous des points pont (condition_ref)
condition_out="PT_UNDER_BRIDGE=0",
max2d_above=0.5,
max2d_below=0.5,
)
pipeline |= pdal.Filter.assign(value=[f"{dsm_dimension}=0 WHERE PT_ON_BRIDGE==1"])
pipeline |= pdal.Filter.assign(value=[f"{dsm_dimension}=0 WHERE PT_UNDER_BRIDGE==1"])
###################################################################################################################
# 6 - Ajout des point pour MNT (sol) qui servent au MNS également
###################################################################################################################

pipeline |= pdal.Filter.assign(
value=[
f"{dsm_dimension}=1 WHERE ({dtm_dimension}==1 && PT_VEG_DSM==0 && PT_ON_BRIDGE==0 && PT_ON_BUILDING==0 && PT_ON_VEGET==0)"
f"{dsm_dimension}=1 WHERE ({dtm_dimension}==1 && PT_VEG_DSM==0 && PT_UNDER_BRIDGE==0 && PT_CLOSED_BUILDING==0 && PT_UNDER_VEGET==0)"
]
)

Expand All @@ -333,15 +333,15 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
False,
condition_src="Classification==66",
condition_ref=macro.build_condition("Classification", [4, 5, 6, 17, 67]),
condition_out="PT_ON_VEGET=1",
condition_out="PT_UNDER_VEGET=1",
)
pipeline = macro.add_radius_assign(
pipeline,
0.5,
False,
condition_src="PT_ON_VEGET==1 && Classification==66",
condition_ref="PT_ON_VEGET==0 && Classification==66",
condition_out="PT_ON_VEGET=0",
condition_src="PT_UNDER_VEGET==1 && Classification==66",
condition_ref="PT_UNDER_VEGET==0 && Classification==66",
condition_out="PT_UNDER_VEGET=0",
)
# 7.3 Taguage pour les MNS des points virtuels eau seulement
pipeline = macro.add_radius_assign(
Expand All @@ -350,11 +350,11 @@ def define_marking_pipeline(input_las, output_las, dsm_dimension, dtm_dimension)
False,
condition_src="Classification==66",
condition_ref="Classification==17",
condition_out="PT_ON_BRIDGE=1",
condition_out="PT_UNDER_BRIDGE=1",
)
pipeline |= pdal.Filter.assign(
value=[
f"{dsm_dimension}=1 WHERE (Classification==66 && PT_ON_VEGET==0 && PT_ON_BRIDGE==0)"
f"{dsm_dimension}=1 WHERE (Classification==66 && PT_UNDER_VEGET==0 && PT_UNDER_BRIDGE==0)"
]
)

Expand Down
33 changes: 17 additions & 16 deletions src/filter_radius_assign/RadiusAssignFilter.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -68,33 +68,34 @@ void RadiusAssignFilter::ready(PointTableRef)
m_args->m_ptsToUpdate.clear();
}

void RadiusAssignFilter::doOneNoDomain(PointRef &point)
void RadiusAssignFilter::doOneNoDomain(PointRef &pointSrc)
{
// build3dIndex and build2dIndex are build once
PointIdList iNeighbors;
if (m_args->search3d) iNeighbors = refView->build3dIndex().radius(point, m_args->m_radius);
else iNeighbors = refView->build2dIndex().radius(point, m_args->m_radius);

if (m_args->search3d) iNeighbors = refView->build3dIndex().radius(pointSrc, m_args->m_radius);
else iNeighbors = refView->build2dIndex().radius(pointSrc, m_args->m_radius);
if (iNeighbors.size() == 0)
return;

if (!m_args->search3d)
if (!m_args->search3d && (m_args->m_max2d_below>=0 || m_args->m_max2d_above>=0))
{
double Zref = point.getFieldAs<double>(Dimension::Id::Z);
if (m_args->m_max2d_below>=0 || m_args->m_max2d_above>=0)
double Zsrc = pointSrc.getFieldAs<double>(Dimension::Id::Z);

bool take (false);
for (PointId ptId : iNeighbors)
{
bool take (true);
for (PointId ptId : iNeighbors)
{
double Zpt = refView->point(ptId).getFieldAs<double>(Dimension::Id::Z);
if (m_args->m_max2d_below>=0 && (Zref-Zpt)>m_args->m_max2d_below) {take=false; break;}
if (m_args->m_max2d_above>=0 && (Zpt-Zref)>m_args->m_max2d_above) {take=false; break;}
}
if (!take) return;
double Zref = refView->point(ptId).getFieldAs<double>(Dimension::Id::Z);
if (m_args->m_max2d_above>=0 && Zref>Zsrc && (Zref-Zsrc)>m_args->m_max2d_above) continue;
if (m_args->m_max2d_below>=0 && Zsrc>Zref && (Zsrc-Zref)>m_args->m_max2d_below) continue;
take = true;
break;
}

if (!take) return;
}

m_args->m_ptsToUpdate.push_back(point.pointId());
m_args->m_ptsToUpdate.push_back(pointSrc.pointId());
}

bool RadiusAssignFilter::doOne(PointRef& point)
Expand Down
Binary file added test/data/mnx/input/bat.laz
Binary file not shown.
Binary file added test/data/mnx/input/corse.laz
Binary file not shown.
Binary file added test/data/mnx/input/pont.laz
Binary file not shown.
Binary file added test/data/mnx/reference/bat.laz
Binary file not shown.
Binary file added test/data/mnx/reference/corse.laz
Binary file not shown.
Binary file modified test/data/mnx/reference/crop_1.laz
Binary file not shown.
Binary file modified test/data/mnx/reference/crop_2.laz
Binary file not shown.
Binary file modified test/data/mnx/reference/crop_3.laz
Binary file not shown.
Binary file added test/data/mnx/reference/pont.laz
Binary file not shown.
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,12 @@ def test_mark_points_to_use_for_digital_models_with_new_dimension_keep_dimension
assert all(
[
dim in output_dimensions
for dim in ["PT_VEG_DSM", "PT_ON_BRIDGE", "PT_ON_BUILDING", "PT_ON_VEGET"]
for dim in [
"PT_VEG_DSM",
"PT_UNDER_BRIDGE",
"PT_CLOSED_BUILDING",
"PT_UNDER_VEGET",
]
]
)

Expand Down Expand Up @@ -140,6 +145,9 @@ def test_parse_args():
"crop_1.laz",
# "crop_2.laz", ToDo : rebuild the reference for crop_2 which is false
"crop_3.laz",
"bat.laz",
"pont.laz",
"corse.laz",
],
)
def test_algo_mark_points_for_dm_with_reference(crop):
Expand Down
31 changes: 18 additions & 13 deletions test/test_radius_assign.py
Original file line number Diff line number Diff line change
Expand Up @@ -14,6 +14,7 @@
pt_ini = (pt_x, pt_y, pt_z, 1)

numeric_precision = 4
numeric_precision_z = 2
distance_radius = 1


Expand All @@ -28,6 +29,10 @@ def distance3d(pt1, pt2):
)


def distanceZ(pt1, pt2):
return round(pt1[2] - pt2[2], numeric_precision_z)


def run_filter(arrays_las, distance_radius, search_3d, limit_z_above=-1, limit_w_below=-1):

filter = "filters.radius_assign"
Expand Down Expand Up @@ -154,7 +159,7 @@ def test_radius_assign_2d_cylinder_below():

def func_test(pt):
distance_i = distance2d(pt_ini, pt)
distance_z = pt[2] - pt_ini[2]
distance_z = distanceZ(pt, pt_ini)
if distance_i < distance_radius and distance_z < limit_z_below:
return 1
return 0
Expand All @@ -177,7 +182,7 @@ def test_radius_assign_2d_cylinder_above():

def func_test(pt):
distance_i = distance2d(pt_ini, pt)
distance_z = pt_ini[2] - pt[2]
distance_z = distanceZ(pt_ini, pt)
if distance_i < distance_radius and distance_z < limit_z_above:
return 1
return 0
Expand All @@ -197,7 +202,7 @@ def test_radius_assign_2d_cylinder_above_below_null():

def func_test(pt):
distance_i = distance2d(pt_ini, pt)
distance_z = pt_ini[2] - pt[2]
distance_z = distanceZ(pt_ini, pt)
if distance_i < distance_radius and distance_z == 0:
return 1
return 0
Expand All @@ -217,7 +222,7 @@ def test_radius_assign_2d_cylinder_above_null_bellow_all():

def func_test(pt):
distance_i = distance2d(pt_ini, pt)
distance_z = pt[2] - pt_ini[2]
distance_z = distanceZ(pt, pt_ini)
if distance_i < distance_radius and distance_z <= 0:
return 1
return 0
Expand Down Expand Up @@ -261,13 +266,12 @@ def test_radius_assign_2d_cylinder_above_and_bellow(execution_number):

def func_test(pt):
distance_i = distance2d(pt_ini, pt)
distance_z = pt_ini[2] - pt[2]
if (
distance_i < distance_radius
and distance_z <= limit_z_above
and (-distance_z) <= limit_z_below
):
return 1
distance_z = distanceZ(pt, pt_ini) # src - ref
if distance_i < distance_radius:
if distance_z <= 0 and (-distance_z) <= limit_z_above: # src est sur ref
return 1
if distance_z >= 0 and distance_z <= limit_z_below: # src est sous ref
return 1
return 0

arrays_las, nb_points_take_2d = build_random_points_around_one_point(func_test, points)
Expand Down Expand Up @@ -298,9 +302,10 @@ def test_radius_assign_2d_cylinder(limit_z_above, limit_z_below):
def func_test(pt):
distance_i = distance2d(pt_ini, pt)
if distance_i < distance_radius:
if (limit_z_above >= 0) and ((pt_ini[2] - pt[2]) > limit_z_above):
distance_z = distanceZ(pt, pt_ini) # src - ref
if limit_z_above >= 0 and distance_z <= 0 and (-distance_z) > limit_z_above:
return 0
if (limit_z_below >= 0) and ((pt[2] - pt_ini[2]) > limit_z_below):
if limit_z_below >= 0 and distance_z >= 0 and distance_z > limit_z_below:
return 0
return 1
else:
Expand Down