Skip to content

Commit

Permalink
Merge pull request #331 from nicolas-f/swap_delaunay
Browse files Browse the repository at this point in the history 
Replace Poly2Tree by TinFour
  • Loading branch information
@nicolas-f
nicolas-f authored Apr 19, 2021
2 parents 8260f3f + d896d5e commit 6b6d8e9
Show file tree
Hide file tree
Showing 17 changed files with 670 additions and 664 deletions.
2 changes: 1 addition & 1 deletion h2gis-extension/pom.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
<parent>
<groupId>org.orbisgis</groupId>
<artifactId>noisemodelling-parent</artifactId>
<version>3.3.1</version>
<version>3.3.2-SNAPSHOT</version>
<relativePath>../pom.xml</relativePath>
</parent>
<description>Additional drivers and functions for H2GIS database</description>
Expand Down
2 changes: 1 addition & 1 deletion noisemodelling-emission/pom.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
<parent>
<groupId>org.orbisgis</groupId>
<artifactId>noisemodelling-parent</artifactId>
<version>3.3.1</version>
<version>3.3.2-SNAPSHOT</version>
<relativePath>../pom.xml</relativePath>
</parent>
<description>Translate light vehicle, heavy vehicle, public bus, tramway traffic into linear noise source. Then compute sound propagation and iso surface.</description>
Expand Down
2 changes: 1 addition & 1 deletion noisemodelling-jdbc/pom.xml
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,7 +10,7 @@
<parent>
<groupId>org.orbisgis</groupId>
<artifactId>noisemodelling-parent</artifactId>
<version>3.3.1</version>
<version>3.3.2-SNAPSHOT</version>
<relativePath>../pom.xml</relativePath>
</parent>
<description>Compute sound propagation rays.</description>
Expand Down
24 changes: 20 additions & 4 deletions noisemodelling-jdbc/src/main/java/org/noise_planet/noisemodelling/jdbc/JdbcNoiseMap.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -197,7 +197,16 @@ protected void fetchCellSoilAreas(Connection connection, Envelope fetchEnvelope,
}
}


void fetchCellBuildings(Connection connection, Envelope fetchEnvelope, MeshBuilder mesh) throws SQLException {
ArrayList<MeshBuilder.PolygonWithHeight> buildings = new ArrayList<>();
fetchCellBuildings(connection, fetchEnvelope, buildings);
for(MeshBuilder.PolygonWithHeight building : buildings) {
mesh.addGeometry(building);
}
}

void fetchCellBuildings(Connection connection, Envelope fetchEnvelope, List<MeshBuilder.PolygonWithHeight> buildings) throws SQLException {
Geometry envGeo = geometryFactory.toGeometry(fetchEnvelope);
boolean fetchAlpha = JDBCUtilities.hasField(connection, buildingsTableName, alphaFieldName);
String additionalQuery = "";
Expand Down Expand Up @@ -250,11 +259,14 @@ void fetchCellBuildings(Connection connection, Envelope fetchEnvelope, MeshBuild
alphaList.add(MeshBuilder.getWallAlpha(oldAlpha, freq));
}
}
MeshBuilder.PolygonWithHeight poly = mesh.addGeometry(intersectedGeometry,
heightField.isEmpty() ? Double.MAX_VALUE : rs.getDouble(heightField), alphaList);

MeshBuilder.PolygonWithHeight poly = new MeshBuilder.PolygonWithHeight(intersectedGeometry,
heightField.isEmpty() ? Double.MAX_VALUE : rs.getDouble(heightField),
alphaList);
if(columnIndex != 0) {
poly.setPrimaryKey(rs.getInt(columnIndex));
}
buildings.add(poly);
}
}
}
Expand Down Expand Up @@ -337,8 +349,12 @@ public void initialize(Connection connection, ProgressVisitor progression) throw
if(sourcesTableName.isEmpty()) {
throw new SQLException("A sound source table must be provided");
}
geometryFactory = new GeometryFactory(new PrecisionModel(),
SFSUtilities.getSRID(connection, TableLocation.parse(sourcesTableName)));
int srid = 0;
srid = SFSUtilities.getSRID(connection, TableLocation.parse(sourcesTableName));
if(srid == 0) {
srid = SFSUtilities.getSRID(connection, TableLocation.parse(buildingsTableName));
}
geometryFactory = new GeometryFactory(new PrecisionModel(), srid);

// Steps of execution
// Evaluation of the main bounding box (sourcesTableName+buildingsTableName)
Expand Down
125 changes: 53 additions & 72 deletions noisemodelling-jdbc/src/main/java/org/noise_planet/noisemodelling/jdbc/TriangleNoiseMap.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -11,6 +11,7 @@
import org.locationtech.jts.simplify.TopologyPreservingSimplifier;
import org.noise_planet.noisemodelling.pathfinder.*;
import org.noise_planet.noisemodelling.pathfinder.Triangle;
import org.noise_planet.noisemodelling.pathfinder.utils.Densifier3D;
import org.slf4j.Logger;
import org.slf4j.LoggerFactory;

Expand All @@ -33,15 +34,13 @@
*/
public class TriangleNoiseMap extends JdbcNoiseMap {
private static final int BATCH_MAX_SIZE = 100;
private final static double BUILDING_BUFFER = 0.5;
private Logger logger = LoggerFactory.getLogger(TriangleNoiseMap.class);
private double roadWidth = 2;
private double sourceDensification = 4;
private double maximumArea = 75;
private long nbreceivers = 0;
private double receiverHeight = 1.6;
private double buildingBuffer = 2;

private String exceptionDumpFolder = "";

/**
* @param buildingsTableName Buildings table
Expand All @@ -51,6 +50,20 @@ public TriangleNoiseMap(String buildingsTableName, String sourcesTableName) {
super(buildingsTableName, sourcesTableName);
}

/**
* @return When an exception occur, this folder with receiver the input data
*/
public String getExceptionDumpFolder() {
return exceptionDumpFolder;
}

/**
* @param exceptionDumpFolder When an exception occur, this folder with receiver the input data
*/
public void setExceptionDumpFolder(String exceptionDumpFolder) {
this.exceptionDumpFolder = exceptionDumpFolder;
}

/**
* @return Do not add receivers closer to specified distance
*/
Expand All @@ -67,15 +80,15 @@ public void setBuildingBuffer(double buildingBuffer) {
}

private void explodeAndAddPolygon(Geometry intersectedGeometry,
MeshBuilder delaunayTool)
LayerDelaunay delaunayTool)
throws LayerDelaunayError {
if (intersectedGeometry instanceof GeometryCollection) {
for (int j = 0; j < intersectedGeometry.getNumGeometries(); j++) {
Geometry subGeom = intersectedGeometry.getGeometryN(j);
explodeAndAddPolygon(subGeom, delaunayTool);
}
} else {
delaunayTool.addGeometry(intersectedGeometry);
} else if(intersectedGeometry instanceof Polygon && !intersectedGeometry.isEmpty()){
delaunayTool.addPolygon((Polygon)intersectedGeometry, 1);
}
}

Expand All @@ -90,9 +103,9 @@ private Geometry merge(LinkedList<Geometry> toUnite, double bufferSize) {
return bufferOp.getResultGeometry(bufferSize);
}

private void feedDelaunay(Collection<Geometry> buildings, MeshBuilder delaunayTool, Envelope boundingBoxFilter,
private void feedDelaunay(List<MeshBuilder.PolygonWithHeight> buildings, LayerDelaunay delaunayTool, Envelope boundingBoxFilter,
double srcDistance, LinkedList<LineString> delaunaySegments, double minRecDist,
double srcPtDist, double triangleSide, double buildingBuffer) throws LayerDelaunayError {
double buildingBuffer) throws LayerDelaunayError {
Envelope extendedEnvelope = new Envelope(boundingBoxFilter);
extendedEnvelope.expandBy(srcDistance * 2.);
Geometry linearRing = geometryFactory.toGeometry(boundingBoxFilter);
Expand All @@ -104,22 +117,18 @@ private void feedDelaunay(Collection<Geometry> buildings, MeshBuilder delaunayTo
Envelope fetchBox = new Envelope(boundingBoxFilter);
fetchBox.expandBy(buildingBuffer);
Geometry fetchGeometry = geometryFactory.toGeometry(fetchBox);
for(Geometry building : buildings) {
if(building.intersects(fetchGeometry)) {
toUnite.add(building);
for(MeshBuilder.PolygonWithHeight building : buildings) {
if(building.getGeometry().intersects(fetchGeometry)) {
toUnite.add(building.getGeometry());
}
}
// Reduce small artifacts to avoid, shortest geometry to be
// over-triangulated
LinkedList<Geometry> toUniteFinal = new LinkedList<>();
if (!toUnite.isEmpty()) {
Geometry bufferBuildings = merge(toUnite, buildingBuffer);
bufferBuildings = TopologyPreservingSimplifier.simplify(bufferBuildings,
minRecDist / 2);
// Remove small artifacts due to buildingsTableName buffer
if(triangleSide > 0) {
bufferBuildings = Densifier.densify(bufferBuildings, triangleSide);
}
//bufferBuildings = TopologyPreservingSimplifier.simplify(bufferBuildings,
// minRecDist / 2);
toUniteFinal.add(bufferBuildings); // Add buildingsTableName to triangulation
}
Geometry geom1 = geometryFactory.createPolygon();
Expand All @@ -134,17 +143,7 @@ private void feedDelaunay(Collection<Geometry> buildings, MeshBuilder delaunayTo
// Remove small artifacts due to multiple buffer crosses
bufferRoads = TopologyPreservingSimplifier.simplify(bufferRoads,
minRecDist / 2);
// Densify roads to set more receiver near roads.
if(srcPtDist > 0){
bufferRoads = Densifier.densify(bufferRoads, srcPtDist);
} else if (triangleSide > 0) {
bufferRoads = Densifier.densify(bufferRoads, triangleSide);
}
//Add points buffer to the final triangulation, this will densify sound level extraction near
//toUniteFinal.add(makeBufferSegmentsNearRoads(toUniteRoads,srcPtDist));
//roads, and helps to reduce over estimation due to inappropriate interpolation.
toUniteFinal.add(bufferRoads); // Merge roads with minRecDist m
// buffer
}
}
Geometry union = merge(toUniteFinal, 0.); // Merge roads and buildingsTableName
Expand All @@ -171,13 +170,12 @@ private void feedDelaunay(Collection<Geometry> buildings, MeshBuilder delaunayTo
* @param cellJ J cell index
* @param maxSrcDist Maximum propagation distance
* @param minRecDist Minimal distance receiver-source
* @param srcPtDist Densification distance of sources pts
* @param maximumArea Maximum area of triangles
* @throws LayerDelaunayError
*/
public void computeDelaunay(MeshBuilder cellMesh,
public void computeDelaunay(LayerDelaunay cellMesh,
Envelope mainEnvelope, int cellI, int cellJ, double maxSrcDist, Collection<Geometry> sources,
double minRecDist, double srcPtDist, double maximumArea, double buildingBuffer)
double minRecDist, double maximumArea, double buildingBuffer, List<MeshBuilder.PolygonWithHeight> buildings)
throws LayerDelaunayError {

Envelope cellEnvelope = getCellEnv(mainEnvelope, cellI, cellJ,
Expand All @@ -199,7 +197,7 @@ public void computeDelaunay(MeshBuilder cellMesh,
if (ptEnv.intersects(expandedCellEnvelop)) {
if (pt instanceof Point) {
// Add square in rendering
cellMesh.addGeometry(cellEnvelopeGeometry.intersection(pt.buffer(minRecDist, BufferParameters.CAP_SQUARE)));
cellMesh.addPolygon((Polygon)cellEnvelopeGeometry.intersection(pt.buffer(minRecDist, BufferParameters.CAP_SQUARE)), 1);
} else {
if (pt instanceof LineString) {
delaunaySegments.add((LineString) (pt));
Expand All @@ -215,26 +213,23 @@ public void computeDelaunay(MeshBuilder cellMesh,
}
}

// Compute equilateral triangle side from Area
double triangleSide = (2*Math.pow(maximumArea, 0.5)) / Math.pow(3, 0.25);
List<Geometry> buildings = new ArrayList<>(cellMesh.getPolygonWithHeight().size());
for(MeshBuilder.PolygonWithHeight poly : cellMesh.getPolygonWithHeight()) {
buildings.add(poly.getGeometry());
}
cellMesh.clearBuildings();
feedDelaunay(buildings, cellMesh, cellEnvelope, maxSrcDist, delaunaySegments,
minRecDist, srcPtDist, triangleSide, buildingBuffer);
minRecDist, buildingBuffer);

// Process delaunay
logger.info("Begin delaunay");
cellMesh.setComputeNeighbors(false);
cellMesh.setRetrieveNeighbors(false);
// Add cell envelope
if (maximumArea > 1) {
cellMesh.setMaximumArea(maximumArea);
Geometry densifiedEnvelope = Densifier.densify(new GeometryFactory().toGeometry(cellEnvelope), triangleSide);
cellMesh.finishPolygonFeeding(densifiedEnvelope);
cellMesh.setMaxArea(maximumArea);
double triangleSide = (2*Math.pow(maximumArea, 0.5)) / Math.pow(3, 0.25);
Polygon polygon = (Polygon)Densifier.densify(new GeometryFactory().toGeometry(cellEnvelope), triangleSide);
cellMesh.addLineString(polygon.getExteriorRing(), 0);
} else {
cellMesh.finishPolygonFeeding(cellEnvelope);
Polygon polygon = (Polygon) new GeometryFactory().toGeometry(cellEnvelope);
cellMesh.addLineString(polygon.getExteriorRing(), 0);
}
cellMesh.processDelaunay();
logger.info("End delaunay");
}

Expand All @@ -244,34 +239,33 @@ protected Envelope getComputationEnvelope(Connection connection) throws SQLExcep
}

public void generateReceivers(Connection connection, int cellI, int cellJ, String receiverTableName, String trianglesTableName, AtomicInteger receiverPK) throws SQLException, LayerDelaunayError, IOException {

int ij = cellI * gridDim + cellJ + 1;
if(verbose) {
logger.info("Begin processing of cell " + ij + " / " + gridDim * gridDim);
}
// Compute the first pass delaunay mesh
// The first pass doesn't take account of additional
// vertices of neighbor cells at the borders
// then, there are discontinuities in iso surfaces at each
// border of cell
MeshBuilder cellMesh = new MeshBuilder();
Envelope cellEnvelope = getCellEnv(mainEnvelope, cellI,
cellJ, getCellWidth(), getCellHeight());
// Fetch all source located in expandedCellEnvelop
PropagationProcessData data = new PropagationProcessData(null);
fetchCellSource(connection, cellEnvelope, data);

List<Geometry> sourceDelaunayGeometries = data.sourceGeometries;
fetchCellBuildings(connection, cellEnvelope, cellMesh);

MeshBuilder demMesh = new MeshBuilder();
FastObstructionTest freeFieldFinder = null;
if(!demTable.isEmpty()) {
fetchCellDem(connection, cellEnvelope, demMesh);
demMesh.finishPolygonFeeding(cellEnvelope);
freeFieldFinder = new FastObstructionTest(demMesh.getPolygonWithHeight(),
demMesh.getTriangles(), demMesh.getTriNeighbors(), demMesh.getVertices());
}

ArrayList<MeshBuilder.PolygonWithHeight> buildings = new ArrayList<>();
fetchCellBuildings(connection, cellEnvelope, buildings);

LayerTinfour cellMesh = new LayerTinfour();
cellMesh.setDumpFolder(exceptionDumpFolder);
try {
computeDelaunay(cellMesh, mainEnvelope, cellI,
cellJ,
maximumPropagationDistance, sourceDelaunayGeometries, roadWidth,
sourceDensification, maximumArea, buildingBuffer);
maximumPropagationDistance, sourceDelaunayGeometries, roadWidth, maximumArea, buildingBuffer, buildings);
} catch (LayerDelaunayError err) {
throw new SQLException(err.getLocalizedMessage(), err);
}
Expand All @@ -287,9 +281,6 @@ public void generateReceivers(Connection connection, int cellI, int cellJ, Strin
for(Coordinate vertex : cellMesh.getVertices()) {
Coordinate translatedVertex = new Coordinate(vertex);
double z = receiverHeight;
if(freeFieldFinder != null) {
z = freeFieldFinder.getHeightAtPosition(translatedVertex) + receiverHeight;
}
translatedVertex.setOrdinate(2, z);
vertices.add(translatedVertex);
}
Expand All @@ -312,13 +303,11 @@ public void generateReceivers(Connection connection, int cellI, int cellJ, Strin
}
int receiverPkOffset = receiverPK.get();
// Add vertices to receivers
PreparedStatement ps = connection.prepareStatement("INSERT INTO "+TableLocation.parse(receiverTableName)+" VALUES (?, ST_MAKEPOINT(?,?,?));");
PreparedStatement ps = connection.prepareStatement("INSERT INTO "+TableLocation.parse(receiverTableName)+" VALUES (?, ?);");
int batchSize = 0;
for(Coordinate v : vertices) {
ps.setInt(1, receiverPK.getAndAdd(1));
ps.setDouble(2, v.x);
ps.setDouble(3, v.y);
ps.setDouble(4, v.z);
ps.setObject(2, geometryFactory.createPoint(v));
ps.addBatch();
batchSize++;
if (batchSize >= BATCH_MAX_SIZE) {
Expand Down Expand Up @@ -361,14 +350,6 @@ public void setRoadWidth(double roadWidth) {
this.roadWidth = roadWidth;
}

public double getSourceDensification() {
return sourceDensification;
}

public void setSourceDensification(double sourceDensification) {
this.sourceDensification = sourceDensification;
}

public double getMaximumArea() {
return maximumArea;
}
Expand Down
29 changes: 29 additions & 0 deletions ...emodelling-jdbc/src/test/java/org/noise_planet/noisemodelling/jdbc/PointNoiseMapTest.java
View file
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,9 @@

import org.h2gis.api.EmptyProgressVisitor;
import org.h2gis.functions.factory.H2GISDBFactory;
import org.h2gis.functions.io.shp.SHPRead;
import org.h2gis.functions.io.shp.SHPWrite;
import org.h2gis.functions.io.shp.internal.SHPDriver;
import org.h2gis.utilities.SFSUtilities;
import org.junit.After;
import org.junit.Before;
Expand Down Expand Up @@ -162,6 +165,32 @@ public void testNoiseMapBuilding() throws Exception {
}
}

// @Test
// public void testNoiseMapBuilding2() throws Exception {
// try(Statement st = connection.createStatement()) {
// SHPRead.readShape(connection, LDENPointNoiseMapFactoryTest.class.getResource("roads_traff.shp").getFile(), "ROADS_GEOM");
// SHPRead.readShape(connection, LDENPointNoiseMapFactoryTest.class.getResource("buildings.shp").getFile(), " BUILDINGS");
// TriangleNoiseMap noisemap = new TriangleNoiseMap("BUILDINGS", "ROADS_GEOM");
// noisemap.setReceiverHasAbsoluteZCoordinates(false);
// noisemap.setSourceHasAbsoluteZCoordinates(false);
// noisemap.setHeightField("HEIGHT");
// noisemap.setMaximumArea(300);
// noisemap.setBuildingBuffer(0);
// noisemap.setMaximumPropagationDistance(800);
//
//
//
// noisemap.initialize(connection, new EmptyProgressVisitor());
// AtomicInteger pk = new AtomicInteger(0);
// for(int i=0; i < noisemap.getGridDim(); i++) {
// for(int j=0; j < noisemap.getGridDim(); j++) {
// noisemap.generateReceivers(connection, i, j, "NM_RECEIVERS", "TRIANGLES", pk);
// }
// }
// assertNotSame(0, pk.get());
// SHPWrite.exportTable(connection, "target/triangle.shp", "TRIANGLES");
// }
// }


}
Loading

0 comments on commit 6b6d8e9

Please sign in to comment.