diff --git a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/path/Scene.java b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/path/Scene.java index c8a7a833f..26ffe64d8 100644 --- a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/path/Scene.java +++ b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/path/Scene.java @@ -38,6 +38,7 @@ public class Scene { public static final double DEFAULT_RECEIVER_DIST = 1.0; public static final double DEFAULT_GS = 0.0; public static final double DEFAULT_G = 0.0; + public static final double DEFAULT_G_BUILDING = 0.0; public static final String YAW_DATABASE_FIELD = "YAW"; public static final String PITCH_DATABASE_FIELD = "PITCH"; public static final String ROLL_DATABASE_FIELD = "ROLL"; diff --git a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/CutPoint.java b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/CutPoint.java index 526cfec94..4efb2408a 100644 --- a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/CutPoint.java +++ b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/CutPoint.java @@ -18,24 +18,24 @@ public class CutPoint implements Comparable { /** {@link Coordinate} of the cut point. */ - Coordinate coordinate; + Coordinate coordinate = new Coordinate(); /** Intersection type. */ ProfileBuilder.IntersectionType type; /** Identifier of the cut element. */ - int id; + int id = -1; /** Identifier of the building containing the point. -1 if no building. */ - int buildingId; + int buildingId = -1; /** Identifier of the wall containing the point. -1 if no wall. */ - int wallId; + int wallId = -1; /** Height of the building containing the point. NaN of no building. */ - double height; + double height = Double.NaN; /** Topographic height of the point. */ double zGround = Double.NaN; /** Ground effect coefficient. 0 if there is no coefficient. */ - double groundCoef; + double groundCoef = Double.NaN; /** Wall alpha. NaN if there is no coefficient. */ List wallAlpha = Collections.emptyList(); - boolean corner; + boolean corner = false; //todo with horizontal plane diffraction rework: remove, replace with intersection type-> DIFFRACTION_POINT /** * Constructor using a {@link Coordinate}. @@ -44,14 +44,9 @@ public class CutPoint implements Comparable { * @param id Identifier of the cut element. */ public CutPoint(Coordinate coord, ProfileBuilder.IntersectionType type, int id, boolean corner) { - this.coordinate = new Coordinate(coord.x, coord.y, coord.z); + this.coordinate = new Coordinate(coord); this.type = type; this.id = id; - this.buildingId = -1; - this.wallId = -1; - this.groundCoef = 0; - this.wallAlpha = new ArrayList<>(); - this.height = 0; this.corner = corner; } public CutPoint(Coordinate coord, ProfileBuilder.IntersectionType type, int id) { @@ -59,7 +54,6 @@ public CutPoint(Coordinate coord, ProfileBuilder.IntersectionType type, int id) } public CutPoint() { - coordinate = new Coordinate(); } /** diff --git a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/CutProfile.java b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/CutProfile.java index 52f243239..82d33d738 100644 --- a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/CutProfile.java +++ b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/CutProfile.java @@ -110,7 +110,7 @@ public void addTopoCutPt(Coordinate coord, int id) { /** * In order to reduce the number of reallocation, reserve the provided points size - * @param numberOfPointsToBePushed + * @param numberOfPointsToBePushed Number of items to preallocate */ public void reservePoints(int numberOfPointsToBePushed) { pts.ensureCapacity(pts.size() + numberOfPointsToBePushed); @@ -118,12 +118,15 @@ public void reservePoints(int numberOfPointsToBePushed) { /** * Add a ground effect cutting point. - * @param coord Coordinate of the cutting point. + * @param coordinate Coordinate of the cutting point. * @param id Id of the cut topography. */ - public void addGroundCutPt(Coordinate coord, int id) { - pts.add(new CutPoint(coord, ProfileBuilder.IntersectionType.GROUND_EFFECT, id)); + public CutPoint addGroundCutPt(Coordinate coordinate, int id, double groundCoefficient) { + CutPoint pt = new CutPoint(coordinate, ProfileBuilder.IntersectionType.GROUND_EFFECT, id); + pt.setGroundCoef(groundCoefficient); + pts.add(pt); hasGroundEffectInter = true; + return pt; } /** diff --git a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/ProfileBuilder.java b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/ProfileBuilder.java index e023f8a09..e4bba6b48 100644 --- a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/ProfileBuilder.java +++ b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/ProfileBuilder.java @@ -22,6 +22,8 @@ import org.locationtech.jts.geom.Polygon; import org.locationtech.jts.index.ItemVisitor; import org.locationtech.jts.index.strtree.STRtree; +import org.locationtech.jts.math.Vector2D; +import org.locationtech.jts.math.Vector3D; import org.locationtech.jts.operation.distance.DistanceOp; import org.locationtech.jts.triangulate.quadedge.Vertex; import org.noise_planet.noisemodelling.pathfinder.delaunay.LayerDelaunay; @@ -29,6 +31,7 @@ import org.noise_planet.noisemodelling.pathfinder.delaunay.LayerTinfour; import org.noise_planet.noisemodelling.pathfinder.delaunay.Triangle; import org.noise_planet.noisemodelling.pathfinder.path.PointPath; +import org.noise_planet.noisemodelling.pathfinder.path.Scene; import org.noise_planet.noisemodelling.pathfinder.utils.IntegerTuple; import org.noise_planet.noisemodelling.pathfinder.utils.geometry.JTSUtility; import org.slf4j.Logger; @@ -66,6 +69,8 @@ */ public class ProfileBuilder { public static final double epsilon = 1e-7; + public static final double CENTIMETER = 0.01; + public static final double LEFT_SIDE = Math.PI / 2; /** Class {@link java.util.logging.Logger}. */ private static final Logger LOGGER = LoggerFactory.getLogger(ProfileBuilder.class); /** Default RTree node capacity. */ @@ -98,7 +103,7 @@ public class ProfileBuilder { /** Building RTree. */ private STRtree wallTree = new STRtree(TREE_NODE_CAPACITY); /** RTree with Buildings's walls linestrings, walls linestring, GroundEffect linestrings - * The object is an instance of Wall */ + * The object is an integer. It's an index of the array {@link #processedWalls} */ private STRtree rtree; private STRtree groundEffectsRtree = new STRtree(TREE_NODE_CAPACITY); @@ -474,7 +479,7 @@ public ProfileBuilder addWall(LineString geom, double height, List alpha } for(int i=0; i splitSegment(Coordinate c0, Coordinate c1, doubl /** * Retrieve the cutting profile following the line build from the given coordinates. - * @param c0 Starting point. - * @param c1 Ending point. - * @param gS Default ground effect value + * @param sourceCoordinate Starting point. + * @param receiverCoordinate Ending point. + * @param gS Default source absorption ground effect value if no ground absorption value is found * @return Cutting profile. */ - public CutProfile getProfile(Coordinate c0, Coordinate c1, double gS) { + public CutProfile getProfile(Coordinate sourceCoordinate, Coordinate receiverCoordinate, double gS) { CutProfile profile = new CutProfile(); - //Fetch topography evolution between c0 and c1 - if(topoTree != null) { - addTopoCutPts(c0, c1, profile); + // Add sourceCoordinate + CutPoint sourcePoint = profile.addSource(sourceCoordinate); + int groundAbsorptionIndex = getIntersectingGroundAbsorption(FACTORY.createPoint(sourceCoordinate)); + if(groundAbsorptionIndex >= 0) { + sourcePoint.setGroundCoef(groundAbsorptions.get(groundAbsorptionIndex).getCoefficient()); + } else { + sourcePoint.setGroundCoef(gS); } - // Split line into segments for structures based on RTree in order to limit the number of queries - // (for large area of the line segment envelope) - LineSegment fullLine = new LineSegment(c0, c1); - List lines = splitSegment(c0, c1, maxLineLength); + //Fetch topography evolution between sourceCoordinate and receiverCoordinate + if(topoTree != null) { + addTopoCutPts(sourceCoordinate, receiverCoordinate, profile); + } //Add Buildings/Walls and Ground effect transition points if(rtree != null) { - addGroundBuildingCutPts(lines, fullLine, profile); + LineSegment fullLine = new LineSegment(sourceCoordinate, receiverCoordinate); + addGroundBuildingCutPts(fullLine, profile); + } + + // Add receiver point + CutPoint receiverPoint = profile.addReceiver(receiverCoordinate); + + //Sort all the cut point from sourceCoordinate to receiverCoordinate positions + profile.sort(sourceCoordinate); + + // Propagate ground coefficient for unknown coefficients + double currentCoefficient = sourcePoint.groundCoef; + for (CutPoint cutPoint : profile.pts) { + if(Double.isNaN(cutPoint.groundCoef)) { + cutPoint.setGroundCoef(currentCoefficient); + } else if (cutPoint.getType().equals(GROUND_EFFECT)) { + currentCoefficient = cutPoint.getGroundCoef(); + } } - //Sort all the cut point from c0 to c1 positions - profile.sort(c0); return profile; } @@ -1042,131 +1035,89 @@ public CutProfile getProfile(Coordinate c0, Coordinate c1, double gS) { * @param query The geometry object to check for intersection * @return The ground absorption object or null if nothing is found here */ - public GroundAbsorption getIntersectingGroundAbsorption(Geometry query) { + public int getIntersectingGroundAbsorption(Geometry query) { if(groundEffectsRtree != null) { var res = groundEffectsRtree.query(query.getEnvelopeInternal()); for (Object groundEffectAreaIndex : res) { if(groundEffectAreaIndex instanceof Integer) { GroundAbsorption groundAbsorption = groundAbsorptions.get((Integer) groundEffectAreaIndex); if(groundAbsorption.geom.intersects(query)) { - return groundAbsorption; + return (Integer) groundEffectAreaIndex; } } } } - return null; - } - /** - * In the cut profile, when reaching buildings sides, - * add a point to reach the bottom of the wall before/after the building intersection - * @param profile - * @param c0 - * @param c1 - * @param gS Default ground coefficient - */ - private void addBuildingBaseCutPts(CutProfile profile, Coordinate c0, Coordinate c1, double gS) { - ArrayList pts = new ArrayList<>(profile.pts.size()); - int buildId = -1; - CutPoint lastBuild = null; - for(int i=0; i lines, LineSegment fullLine, CutProfile profile) { - List indexes = new ArrayList<>(); + private void addGroundBuildingCutPts(LineSegment fullLine, CutProfile profile) { + Vector2D directionBefore = Vector2D.create(fullLine.p1, fullLine.p0).normalize().multiply(CENTIMETER); + Vector2D directionAfter = Vector2D.create(fullLine.p0, fullLine.p1).normalize().multiply(CENTIMETER); + // Collect all objects where envelope intersects all sub-segments of fullLine + Set indexes = new HashSet<>(); + + // Segmented fullLine, this is the query for rTree indexes + // Split line into segments for structures based on RTree in order to limit the number of queries + // (for large area of the line segment envelope) + List lines = splitSegment(fullLine.p0, fullLine.p1, maxLineLength); for (LineSegment line : lines) { indexes.addAll(rtree.query(new Envelope(line.p0, line.p1))); } - indexes = indexes.stream().distinct().collect(Collectors.toList()); - Map processedGround = new HashMap<>(); for (int i : indexes) { Wall facetLine = processedWalls.get(i); Coordinate intersection = fullLine.intersection(facetLine.ls); if (intersection != null) { intersection = new Coordinate(intersection); if(!isNaN(facetLine.p0.z) && !isNaN(facetLine.p1.z)) { - if(facetLine.p0.z == facetLine.p1.z) { + // same z in the line, so useless to compute interpolation between points + if(Double.compare(facetLine.p0.z, facetLine.p1.z) == 0) { intersection.z = facetLine.p0.z; + } else { + intersection.z = Vertex.interpolateZ(intersection, facetLine.p0, facetLine.p1); } - else { - intersection.z = facetLine.p0.z + ((intersection.x - facetLine.p0.x) / (facetLine.p1.x - facetLine.p0.x) * (facetLine.p1.z - facetLine.p0.z)); - } - } - else if(topoTree == null) { - intersection.z = NaN; - } - else { - intersection.z = getZGround(intersection); } if(facetLine.type == IntersectionType.BUILDING) { - profile.addBuildingCutPt(intersection, facetLine.originId, i, facetLine.p0.equals(intersection)||facetLine.p1.equals(intersection)); + CutPoint pt = profile.addBuildingCutPt(intersection, facetLine.originId, i,false); + pt.setGroundCoef(Scene.DEFAULT_G_BUILDING); + pt.setWallAlpha(buildings.get(facetLine.getOriginId()).alphas); + // add a point at the bottom of the building on the exterior side of the building + Vector2D facetVector = Vector2D.create(facetLine.p0, facetLine.p1); + // exterior polygon segments are CW, so the exterior of the polygon is on the left side of the vector + // it works also with polygon holes as interiors are CCW + Vector2D exteriorVector = facetVector.rotate(LEFT_SIDE).normalize().multiply(CENTIMETER); + Coordinate exteriorPoint = exteriorVector.add(Vector2D.create(intersection)).toCoordinate(); + CutPoint exteriorPointCutPoint = profile.addBuildingCutPt(exteriorPoint, facetLine.originId, i,false); + if(topoTree == null) { + exteriorPointCutPoint.coordinate.setZ(0.0); + } else { + exteriorPointCutPoint.coordinate.setZ(getZGround(exteriorPointCutPoint)); + pt.zGround = exteriorPointCutPoint.coordinate.z; + exteriorPointCutPoint.zGround = exteriorPointCutPoint.coordinate.z; + } } else if(facetLine.type == IntersectionType.WALL) { - profile.addWallCutPt(intersection, facetLine.originId, facetLine.p0.equals(intersection)||facetLine.p1.equals(intersection), facetLine.alphas); + profile.addWallCutPt(intersection, facetLine.originId, false, facetLine.alphas); } else if(facetLine.type == GROUND_EFFECT) { - if(!intersection.equals(facetLine.p0) && !intersection.equals(facetLine.p1)) { - // we hit the border of a ground effect - // we need to add a new point with the new value of the ground effect - // we will query for the point that lie after the intersection with the ground effect border - // in order to have the new value of the ground effect, if there is nothing at this location - // we fall back to the default value of ground effect - - // TODO + // we hit the border of a ground effect + // we need to add a new point with the new value of the ground effect + // we will query for the point that lie after the intersection with the ground effect border + // in order to have the new value of the ground effect, if there is nothing at this location + // we fall back to the default value of ground effect + // if this is another ground effect it will be processed in another loop (two intersections on the same coordinate) + // retrieve the ground coefficient after the intersection in the direction of the profile + // this method will solve the question if we enter a new ground absorption or we will leave one + Point afterIntersectionPoint = FACTORY.createPoint(Vector2D.create(intersection).add(directionAfter).toCoordinate()); + GroundAbsorption groundAbsorption = groundAbsorptions.get(facetLine.getOriginId()); + if(groundAbsorption.geom.intersects(afterIntersectionPoint)) { + // we enter a new ground effect + profile.addGroundCutPt(intersection, facetLine.getOriginId(), groundAbsorption.getCoefficient()); + } else { + // no new ground effect, we fall back to default G + profile.addGroundCutPt(intersection, facetLine.getOriginId(), Scene.DEFAULT_G); } } } diff --git a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/Wall.java b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/Wall.java index af8a610ca..289364900 100644 --- a/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/Wall.java +++ b/noisemodelling-pathfinder/src/main/java/org/noise_planet/noisemodelling/pathfinder/profilebuilder/Wall.java @@ -23,8 +23,6 @@ public class Wall implements ProfileBuilder.Obstacle { List alphas; /** Wall height, if -1, use z coordinate. */ double height; - boolean hasP0Neighbour = false; - boolean hasP1Neighbour = false; public Coordinate p0; public Coordinate p1; LineSegment ls; @@ -76,24 +74,6 @@ public Wall(Coordinate p0, Coordinate p1, int originId, ProfileBuilder.Intersect this.alphas = new ArrayList<>(); } - /** - * Constructor using start/end point and id. - * @param p0 Start point of the segment. - * @param p1 End point of the segment. - * @param originId Id or index of the source building or topographic triangle. - */ - public Wall(Coordinate p0, Coordinate p1, int originId, ProfileBuilder.IntersectionType type, boolean hasP0Neighbour, boolean hasP1Neighbour) { - this.line = FACTORY.createLineString(new Coordinate[]{p0, p1}); - this.p0 = p0; - this.p1 = p1; - this.ls = new LineSegment(p0, p1); - this.originId = originId; - this.type = type; - this.alphas = new ArrayList<>(); - this.hasP0Neighbour = hasP0Neighbour; - this.hasP1Neighbour = hasP1Neighbour; - } - /** * @return Index of this wall in the ProfileBuild list */ @@ -169,14 +149,6 @@ public ProfileBuilder.IntersectionType getType() { return type; } - /*public boolean hasP0Neighbour() { - return hasP0Neighbour; - } - - public boolean hasP1Neighbour() { - return hasP1Neighbour; - }*/ - public ProfileBuilder.Obstacle getObstacle() { return obstacle; }