diff --git a/contribs/application/src/main/java/org/matsim/application/analysis/traffic/CountComparisonAnalysis.java b/contribs/application/src/main/java/org/matsim/application/analysis/traffic/CountComparisonAnalysis.java index adc9a86f435..e0c02a41831 100644 --- a/contribs/application/src/main/java/org/matsim/application/analysis/traffic/CountComparisonAnalysis.java +++ b/contribs/application/src/main/java/org/matsim/application/analysis/traffic/CountComparisonAnalysis.java @@ -76,6 +76,16 @@ private static int[] sum(int[] a, int[] b) { return counts; } + /** + * Calculate the geh value for simulated and reference count + */ + private static double geh(double simulated, double observed) { + final double diff = simulated - observed; + final double sum = simulated + observed; + + return Math.sqrt(2 * diff * diff / sum); + } + @Override public Integer call() throws Exception { @@ -114,14 +124,16 @@ private Table writeOutput(Counts counts, Network network, VolumesAnalyzer StringColumn.create("road_type"), IntColumn.create("hour"), DoubleColumn.create("observed_traffic_volume"), - DoubleColumn.create("simulated_traffic_volume") + DoubleColumn.create("simulated_traffic_volume"), + DoubleColumn.create("geh") ); Table dailyTrafficVolume = Table.create(StringColumn.create("link_id"), StringColumn.create("name"), StringColumn.create("road_type"), DoubleColumn.create("observed_traffic_volume"), - DoubleColumn.create("simulated_traffic_volume") + DoubleColumn.create("simulated_traffic_volume"), + DoubleColumn.create("geh") ); for (Map.Entry, MeasurementLocation> entry : counts.getMeasureLocations().entrySet()) { @@ -170,6 +182,7 @@ private Table writeOutput(Counts counts, Network network, VolumesAnalyzer row.setInt("hour", hour); row.setDouble("observed_traffic_volume", observedTrafficVolumeAtHour); row.setDouble("simulated_traffic_volume", simulatedTrafficVolumeAtHour); + row.setDouble("geh", geh(simulatedTrafficVolumeAtHour, observedTrafficVolumeAtHour)); } } else { // Get the daily values @@ -183,6 +196,7 @@ private Table writeOutput(Counts counts, Network network, VolumesAnalyzer row.setString("road_type", type); row.setDouble("observed_traffic_volume", observedTrafficVolumeByDay); row.setDouble("simulated_traffic_volume", simulatedTrafficVolumeByDay); + row.setDouble("geh", geh(simulatedTrafficVolumeByDay, observedTrafficVolumeByDay)); } DoubleColumn relError = dailyTrafficVolume.doubleColumn("simulated_traffic_volume") diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/Predictor.java b/contribs/application/src/main/java/org/matsim/application/prepare/Predictor.java index 460d1a41abe..473048a99ed 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/Predictor.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/Predictor.java @@ -11,6 +11,7 @@ public interface Predictor { /** * Predict value from given features. + * @return predicted value, maybe NaN if no prediction is possible. */ double predict(Object2DoubleMap features, Object2ObjectMap categories); diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/SampleNetwork.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/SampleNetwork.java index 82a7a278e32..242667bec49 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/SampleNetwork.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/SampleNetwork.java @@ -144,6 +144,10 @@ public Integer call() throws Exception { List list = e.getValue(); + // for now only consider traffic lights + if (!e.getKey().equals("traffic_light")) + continue; + log.info("Sampling {} out of {} intersections for type {}", sample, list.size(), e.getKey()); for (int i = 0; i < sample && !list.isEmpty(); i++) { @@ -193,7 +197,7 @@ public Integer call() throws Exception { RandomizedTravelTime tt = new RandomizedTravelTime(rnd); - LeastCostPathCalculator router = createRandomizedRouter(network, tt); + LeastCostPathCalculator router = createRandomizedRouter(cityNetwork, tt); sampleCityRoutes(cityNetwork, router, tt, rnd); @@ -223,9 +227,16 @@ private void sampleCityRoutes(Network network, LeastCostPathCalculator router, R Link to = NetworkUtils.getNearestLink(network, dest); + // Links could be on the very edge so that nodes are outside the network + if (to == null || !network.getNodes().containsKey(link.getFromNode().getId()) || + !network.getNodes().containsKey(to.getToNode().getId())) { + i--; + continue; + } + LeastCostPathCalculator.Path path = router.calcLeastCostPath(link.getFromNode(), to.getToNode(), 0, null, null); - if (path.nodes.size() < 2) { + if (path == null || path.nodes.size() < 2) { i--; continue; } diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ApplyNetworkParams.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ApplyNetworkParams.java index deef01dfe5e..048e2fa81b7 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ApplyNetworkParams.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ApplyNetworkParams.java @@ -4,7 +4,6 @@ import com.fasterxml.jackson.annotation.JsonInclude; import com.fasterxml.jackson.annotation.PropertyAccessor; import com.fasterxml.jackson.databind.ObjectMapper; -import it.unimi.dsi.fastutil.objects.Object2DoubleMap; import org.apache.logging.log4j.LogManager; import org.apache.logging.log4j.Logger; import org.matsim.api.core.v01.Id; @@ -14,12 +13,16 @@ import org.matsim.application.MATSimAppCommand; import org.matsim.application.options.InputOptions; import org.matsim.application.options.OutputOptions; +import org.matsim.application.prepare.Predictor; import org.matsim.application.prepare.network.params.NetworkParamsOpt.Feature; import org.matsim.core.network.NetworkUtils; import org.matsim.core.utils.io.IOUtils; import picocli.CommandLine; import java.io.BufferedReader; +import java.math.BigDecimal; +import java.math.RoundingMode; +import java.util.List; import java.util.Map; import java.util.Set; @@ -40,6 +43,7 @@ public class ApplyNetworkParams implements MATSimAppCommand { @CommandLine.Mixin private final OutputOptions output = OutputOptions.ofCommand(ApplyNetworkParams.class); + @CommandLine.Parameters(arity = "1..*", description = "Type of parameters to apply. Available: ${COMPLETION-CANDIDATES}") private Set params; @@ -52,6 +56,19 @@ public class ApplyNetworkParams implements MATSimAppCommand { @CommandLine.Option(names = "--factor-bounds", split = ",", description = "Speed factor limits (lower,upper bound)", defaultValue = NetworkParamsOpt.DEFAULT_FACTOR_BOUNDS) private double[] speedFactorBounds; + @CommandLine.Option(names = "--capacity-bounds", split = ",", defaultValue = "0.4,0.6,0.8", + description = "Minimum relative capacity against theoretical max (traffic light,right before left, priority)") + private List capacityBounds; + + @CommandLine.Option(names = "--road-types", split = ",", description = "Road types to apply changes to") + private Set roadTypes; + + @CommandLine.Option(names = "--junction-types", split = ",", description = "Junction types to apply changes to") + private Set junctionTypes; + + @CommandLine.Option(names = "--decrease-only", description = "Only set values if the are lower than the current value", defaultValue = "false") + private boolean decrease; + private NetworkModel model; private NetworkParams paramsOpt; @@ -62,14 +79,14 @@ public static void main(String[] args) { } /** - * Theoretical capacity. + * Theoretical capacity assuming fixed car length and headway. This should usually not be exceeded. */ - private static double capacityEstimate(double v) { + public static double capacityEstimate(double v) { // headway double tT = 1.2; - // car length + // car length + buffer double lL = 7.0; double Qc = v / (v * tT + lL); @@ -94,13 +111,19 @@ public Integer call() throws Exception { } } - Map, Feature> features = NetworkParamsOpt.readFeatures(input.getPath("features.csv"), network.getLinks().size()); + Map, Feature> features = NetworkParamsOpt.readFeatures(input.getPath("features.csv"), network.getLinks()); for (Link link : network.getLinks().values()) { Feature ft = features.get(link.getId()); + if (roadTypes != null && !roadTypes.isEmpty() && !roadTypes.contains(ft.highwayType())) + continue; + + if (junctionTypes != null && !junctionTypes.isEmpty() && !junctionTypes.contains(ft.junctionType())) + continue; + try { - applyChanges(link, ft.junctionType(), ft.features()); + applyChanges(link, ft); } catch (IllegalArgumentException e) { warn++; log.warn("Error processing link {}", link.getId(), e); @@ -117,67 +140,126 @@ public Integer call() throws Exception { /** * Apply speed and capacity models and apply changes. */ - private void applyChanges(Link link, String junctionType, Object2DoubleMap features) { - - String type = NetworkUtils.getHighwayType(link); + private void applyChanges(Link link, Feature ft) { boolean modified = false; if (params.contains(NetworkAttribute.capacity)) { + modified = applyCapacity(link, ft); + } + + if (params.contains(NetworkAttribute.freespeed)) { + modified |= applyFreeSpeed(link, ft); + } - FeatureRegressor capacity = model.capacity(junctionType); + if (modified) + warn++; + } - double perLane = capacity.predict(features); + private boolean applyCapacity(Link link, Feature ft) { - double cap = capacityEstimate(features.getDouble("speed")); + Predictor capacity = model.capacity(ft.junctionType(), ft.highwayType()); + // No operation performed if not supported + if (capacity == null) { + return false; + } - // Minimum thresholds - double threshold = switch (junctionType) { - // traffic light can reduce capacity at least to 50% (with equal green split) - case "traffic_light" -> 0.4; - case "right_before_left" -> 0.6; - // Motorways are kept at their max theoretical capacity - case "priority" -> type.startsWith("motorway") ? 1 : 0.8; - default -> 0; - }; + double perLane = capacity.predict(ft.features(), ft.categories()); + if (Double.isNaN(perLane)) { + return true; + } - if (perLane < cap * threshold) { - log.warn("Increasing capacity per lane on {} ({}, {}) from {} to {}", link.getId(), type, junctionType, perLane, cap * threshold); - perLane = cap * threshold; - modified = true; - } + double cap = capacityEstimate(ft.features().getDouble("speed")); - link.setCapacity(link.getNumberOfLanes() * perLane); + if (capacityBounds.isEmpty()) + capacityBounds.add(0.0); + + // Fill up to 3 elements if not provided + if (capacityBounds.size() < 3) { + capacityBounds.add(capacityBounds.get(0)); + capacityBounds.add(capacityBounds.get(1)); } + // Minimum thresholds + double threshold = switch (ft.junctionType()) { + case "traffic_light" -> capacityBounds.get(0); + case "right_before_left" -> capacityBounds.get(1); + case "priority" -> capacityBounds.get(2); + default -> 0; + }; - if (params.contains(NetworkAttribute.freespeed)) { + boolean modified = false; - double speedFactor = 1.0; - FeatureRegressor speedModel = model.speedFactor(junctionType); + if (perLane < cap * threshold) { + log.warn("Increasing capacity per lane on {} ({}, {}) from {} to {}", + link.getId(), ft.highwayType(), ft.junctionType(), perLane, cap * threshold); + perLane = cap * threshold; + modified = true; + } - speedFactor = paramsOpt != null ? - speedModel.predict(features, paramsOpt.getParams(junctionType)) : - speedModel.predict(features); + if (perLane > cap) { + log.warn("Reducing capacity per lane on {} ({}, {}) from {} to {}", + link.getId(), ft.highwayType(), ft.junctionType(), perLane, cap); + perLane = cap; + modified = true; + } - if (speedFactor > speedFactorBounds[1]) { - log.warn("Reducing speed factor on {} from {} to {}", link.getId(), speedFactor, speedFactorBounds[1]); - speedFactor = speedFactorBounds[1]; - modified = true; - } + if (ft.features().getOrDefault("num_lanes", link.getNumberOfLanes()) != link.getNumberOfLanes()) + log.warn("Number of lanes for link {} does not match the feature file", link.getId()); - // Threshold for very low speed factors - if (speedFactor < speedFactorBounds[0]) { - log.warn("Increasing speed factor on {} from {} to {}", link, speedFactor, speedFactorBounds[0]); - speedFactor = speedFactorBounds[0]; - modified = true; - } + int totalCap = BigDecimal.valueOf(link.getNumberOfLanes() * perLane).setScale(0, RoundingMode.HALF_UP).intValue(); + + if (decrease && totalCap > link.getCapacity()) + return false; + + link.setCapacity(totalCap); + + return modified; + } + + private boolean applyFreeSpeed(Link link, Feature ft) { - link.setFreespeed((double) link.getAttributes().getAttribute("allowed_speed") * speedFactor); - link.getAttributes().putAttribute("speed_factor", speedFactor); + Predictor speedModel = model.speedFactor(ft.junctionType(), ft.highwayType()); + + // No operation performed if not supported + if (speedModel == null) { + return false; } - if (modified) - warn++; + double speedFactor = paramsOpt != null ? + speedModel.predict(ft.features(), ft.categories(), paramsOpt.getParams(ft.junctionType())) : + speedModel.predict(ft.features(), ft.categories()); + + if (Double.isNaN(speedFactor)) { + return false; + } + + boolean modified = false; + + if (speedFactor > speedFactorBounds[1]) { + log.warn("Reducing speed factor on {} from {} to {}", link.getId(), speedFactor, speedFactorBounds[1]); + speedFactor = speedFactorBounds[1]; + modified = true; + } + + // Threshold for very low speed factors + if (speedFactor < speedFactorBounds[0]) { + log.warn("Increasing speed factor on {} from {} to {}", link, speedFactor, speedFactorBounds[0]); + speedFactor = speedFactorBounds[0]; + modified = true; + } + + double freeSpeed = (double) link.getAttributes().getAttribute("allowed_speed") * speedFactor; + + freeSpeed = BigDecimal.valueOf(freeSpeed).setScale(3, RoundingMode.HALF_EVEN).doubleValue(); + + if (decrease && freeSpeed > link.getFreespeed()) + return false; + + link.setFreespeed(freeSpeed); + link.getAttributes().putAttribute("speed_factor", freeSpeed); + + return modified; } + } diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/EvalFreespeedParams.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/EvalFreespeedParams.java index 014d2615549..cd2ba4e016a 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/EvalFreespeedParams.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/EvalFreespeedParams.java @@ -18,6 +18,7 @@ import org.matsim.application.analysis.traffic.traveltime.SampleValidationRoutes; import org.matsim.application.options.InputOptions; import org.matsim.application.options.OutputOptions; +import org.matsim.application.prepare.Predictor; import org.matsim.contrib.osm.networkReader.LinkProperties; import org.matsim.core.network.NetworkUtils; import picocli.CommandLine; @@ -104,7 +105,7 @@ static Result applyAndEvaluateParams( continue; } - FeatureRegressor speedModel = model.speedFactor(ft.junctionType()); + Predictor speedModel = model.speedFactor(ft.junctionType(), ft.highwayType()); if (speedModel == null) { link.setFreespeed(allowedSpeed); @@ -115,15 +116,15 @@ static Result applyAndEvaluateParams( if (request.hasParams()) { double[] p = request.getParams(ft.junctionType()); - speedFactor = speedModel.predict(ft.features(), p); + speedFactor = speedModel.predict(ft.features(), ft.categories(), p); } else - speedFactor = speedModel.predict(ft.features()); + speedFactor = speedModel.predict(ft.features(), ft.categories()); // apply lower and upper bound speedFactor = Math.max(speedFactorBounds[0], speedFactor); speedFactor = Math.min(speedFactorBounds[1], speedFactor); - attributes.put(link.getId(), speedModel.getData(ft.features())); + attributes.put(link.getId(), speedModel.getData(ft.features(), ft.categories())); link.setFreespeed(allowedSpeed * speedFactor); link.getAttributes().putAttribute("speed_factor", speedFactor); @@ -151,7 +152,7 @@ public Integer call() throws Exception { mapper.setSerializationInclusion(JsonInclude.Include.NON_DEFAULT); validationSet = readValidation(validationFiles, refHours); - features = readFeatures(input.getPath("features.csv"), network.getLinks().size()); + features = readFeatures(input.getPath("features.csv"), network.getLinks()); CSVPrinter csv; Path out = output.getPath("eval.csv"); diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/FeatureRegressor.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/FeatureRegressor.java index 051135e89b6..cade7eb1693 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/FeatureRegressor.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/FeatureRegressor.java @@ -1,11 +1,15 @@ package org.matsim.application.prepare.network.params; import it.unimi.dsi.fastutil.objects.Object2DoubleMap; +import it.unimi.dsi.fastutil.objects.Object2ObjectMap; +import org.matsim.application.prepare.Predictor; /** * Predictor interface for regression. + * @deprecated Use {@link Predictor} instead. */ -public interface FeatureRegressor { +@Deprecated +public interface FeatureRegressor extends Predictor { /** @@ -28,4 +32,22 @@ default double[] getData(Object2DoubleMap ft) { throw new UnsupportedOperationException("Not implemented"); } + default double predict(Object2DoubleMap features, Object2ObjectMap categories) { + return predict(features); + } + + /** + * Predict values with adjusted model params. + */ + default double predict(Object2DoubleMap features, Object2ObjectMap categories, double[] params) { + return predict(features, params); + } + + /** + * Return data that is used for internal prediction function (normalization already applied). + */ + default double[] getData(Object2DoubleMap features, Object2ObjectMap categories) { + return getData(features); + } + } diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/FreespeedOptServer.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/FreespeedOptServer.java index 387824839f4..fef19850305 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/FreespeedOptServer.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/FreespeedOptServer.java @@ -102,7 +102,7 @@ public Integer call() throws Exception { } validationSet = NetworkParamsOpt.readValidation(validationFiles, refHours); - features = NetworkParamsOpt.readFeatures(input.getPath("features.csv"), network.getLinks().size()); + features = NetworkParamsOpt.readFeatures(input.getPath("features.csv"), network.getLinks()); log.info("Initial score:"); applyAndEvaluateParams(null, "init"); diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/NetworkModel.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/NetworkModel.java index f98bb9bd41e..f88c4183c58 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/NetworkModel.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/NetworkModel.java @@ -1,22 +1,24 @@ package org.matsim.application.prepare.network.params; +import org.matsim.application.prepare.Predictor; + /** * A model for estimating network parameters. */ public interface NetworkModel { /** - * Flow Capacity (per lane) + * Flow Capacity (per lane). */ - default FeatureRegressor capacity(String junctionType) { - return null; + default Predictor capacity(String junctionType, String highwayType) { + throw new UnsupportedOperationException("Capacity model not implemented for class: " + getClass().getName()); } /** - * Speed factor (relative to free flow speed). + * Speed factor (relative to allowed speed). */ - default FeatureRegressor speedFactor(String junctionType) { - return null; + default Predictor speedFactor(String junctionType, String highwayType) { + throw new UnsupportedOperationException("Speed factor model not implemented for class: " + getClass().getName()); } } diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/NetworkParamsOpt.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/NetworkParamsOpt.java index 2c69e7c6498..e9185fc7a0f 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/NetworkParamsOpt.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/NetworkParamsOpt.java @@ -4,6 +4,8 @@ import it.unimi.dsi.fastutil.ints.Int2ObjectMap; import it.unimi.dsi.fastutil.objects.Object2DoubleMap; import it.unimi.dsi.fastutil.objects.Object2DoubleOpenHashMap; +import it.unimi.dsi.fastutil.objects.Object2ObjectMap; +import it.unimi.dsi.fastutil.objects.Object2ObjectOpenHashMap; import org.apache.commons.csv.CSVFormat; import org.apache.commons.csv.CSVParser; import org.apache.commons.csv.CSVPrinter; @@ -15,6 +17,7 @@ import org.matsim.api.core.v01.network.Network; import org.matsim.api.core.v01.network.Node; import org.matsim.application.analysis.traffic.traveltime.SampleValidationRoutes; +import org.matsim.core.network.NetworkUtils; import org.matsim.core.router.DijkstraFactory; import org.matsim.core.router.costcalculators.OnlyTimeDependentTravelDisutility; import org.matsim.core.router.util.LeastCostPathCalculator; @@ -55,9 +58,14 @@ static NetworkModel load(Class modelClazz) { /** * Read network edge features from csv. */ - static Map, Feature> readFeatures(String input, int expectedLinks) throws IOException { + static Map, Feature> readFeatures(String input, Map, ? extends Link> links) throws IOException { - Map, Feature> features = new IdMap<>(Link.class, expectedLinks); + Map, Feature> features = new IdMap<>(Link.class, links.size()); + + // Create features from link attributes + for (Link link : links.values()) { + features.put(link.getId(), createDefaultFeature(link)); + } try (CSVParser reader = new CSVParser(IOUtils.getBufferedReader(input), CSVFormat.DEFAULT.builder().setHeader().setSkipHeaderRecord(true).build())) { @@ -67,27 +75,65 @@ static Map, Feature> readFeatures(String input, int expectedLinks) thro for (CSVRecord row : reader) { Id id = Id.createLinkId(row.get("linkId")); - - Object2DoubleOpenHashMap ft = new Object2DoubleOpenHashMap<>(); - ft.defaultReturnValue(Double.NaN); + Link link = links.get(id); + Feature ft = features.computeIfAbsent(id, (k) -> createDefaultFeature(link)); for (String column : header) { String v = row.get(column); try { - ft.put(column, Double.parseDouble(v)); + ft.features.put(column, Double.parseDouble(v)); } catch (NumberFormatException e) { // every not equal to True will be false - ft.put(column, Boolean.parseBoolean(v) ? 1 : 0); + ft.features.put(column, Boolean.parseBoolean(v) ? 1 : 0); + ft.categories.put(column, v); } } - features.put(id, new Feature(row.get("junction_type"), ft)); + String highwayType = header.contains(NetworkUtils.TYPE) ? row.get(NetworkUtils.TYPE) : + (link != null ? NetworkUtils.getHighwayType(link) : null); + + features.put(id, new Feature(row.get("junction_type").intern(), highwayType, ft.features, ft.categories)); } } return features; } + /** + * Create default feature based on link attributes. + */ + private static Feature createDefaultFeature(Link link) { + Object2DoubleOpenHashMap ft = new Object2DoubleOpenHashMap<>(); + ft.defaultReturnValue(Double.NaN); + Object2ObjectMap categories = new Object2ObjectOpenHashMap<>(); + + // Link might not be present in the network + if (link == null) + return new Feature("", "", ft, categories); + + String highwayType = NetworkUtils.getHighwayType(link); + categories.put("highway_type", highwayType); + ft.put("speed", NetworkUtils.getAllowedSpeed(link)); + ft.put("num_lanes", link.getNumberOfLanes()); + ft.put("length", link.getLength()); + ft.put("capacity", link.getCapacity()); + ft.put("freespeed", link.getFreespeed()); + + for (Map.Entry e : link.getAttributes().getAsMap().entrySet()) { + String key = e.getKey(); + Object value = e.getValue(); + if (value instanceof Number) { + ft.put(key, ((Number) value).doubleValue()); + } else if (value instanceof Boolean) { + ft.put(key, (Boolean) value ? 1 : 0); + } else { + categories.put(key, value.toString()); + } + } + + return new Feature("", highwayType, ft, categories); + } + /** * Read validation files and calc target speed. */ @@ -170,7 +216,7 @@ static Result evaluate(Network network, Object2DoubleMap features) { + record Feature(String junctionType, String highwayType, Object2DoubleMap features, Object2ObjectMap categories) { } record Result(double rmse, double mae, Map> data) { diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSMotorwayCapacity.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSMotorwayCapacity.java new file mode 100644 index 00000000000..4a8bb92e219 --- /dev/null +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSMotorwayCapacity.java @@ -0,0 +1,44 @@ +package org.matsim.application.prepare.network.params.hbs; + +import it.unimi.dsi.fastutil.objects.Object2DoubleMap; +import it.unimi.dsi.fastutil.objects.Object2ObjectMap; +import org.matsim.application.prepare.Predictor; + +/** + * Capacity for motorways. + */ +public class HBSMotorwayCapacity implements Predictor { + @Override + public double predict(Object2DoubleMap features, Object2ObjectMap categories) { + + // speed in km/h + int speed = (int) Math.round(features.getDouble("speed") * 3.6); + int lanes = (int) features.getOrDefault("num_lanes", 1); + + // Capacity for 1 lane motorways is not defined in HBS + double capacity = 2000; + if (lanes == 2) { + if (speed >= 130) + capacity = 3700.0; + else if (speed >= 120) + capacity = 3800; + else + capacity = 3750; + } else if (lanes == 3) { + if (speed >= 130) + capacity = 5300; + else if (speed >= 120) + capacity = 5400; + else + capacity = 5350; + } else if (lanes >= 4) { + if (speed >= 130) + capacity = 7300; + else + capacity = 7400; + } + + // Return capacity per lane + return capacity / lanes; + } +} diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSNetworkParams.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSNetworkParams.java new file mode 100644 index 00000000000..faff8379a3c --- /dev/null +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSNetworkParams.java @@ -0,0 +1,31 @@ +package org.matsim.application.prepare.network.params.hbs; + +import org.matsim.application.prepare.Predictor; +import org.matsim.application.prepare.network.params.NetworkModel; + +/** + * Capacity params calculated according to + * "Handbuch für die Bemessung von Straßenverkehrsanlagen“ (HBS) + */ +public class HBSNetworkParams implements NetworkModel { + + private static final Predictor MOTORWAY = new HBSMotorwayCapacity(); + private static final Predictor ROAD = new HBSRoadCapacity(); + + @Override + public Predictor capacity(String junctionType, String highwayType) { + + // Traffic lights are not considered + if (junctionType.equals("traffic_light")) { + return null; + } + + if (highwayType.startsWith("motorway")) { + return MOTORWAY; + } + + // All lower category roads + return ROAD; + } + +} diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSRoadCapacity.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSRoadCapacity.java new file mode 100644 index 00000000000..3e07362db8f --- /dev/null +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/hbs/HBSRoadCapacity.java @@ -0,0 +1,144 @@ +package org.matsim.application.prepare.network.params.hbs; + +import it.unimi.dsi.fastutil.objects.Object2DoubleMap; +import it.unimi.dsi.fastutil.objects.Object2ObjectMap; +import org.matsim.application.prepare.Predictor; + +/** + * Capacity for general roads, that are not motorways or residential roads. + */ +public class HBSRoadCapacity implements Predictor { + /** + * Capacity on "Landstraße", often osm secondary. These numbers are taken from the HBS (see contribs/application/src/main/python/capacity/hbs.py) + */ + private static double capacityLandStr(int lanes, int curvature) { + + + if (lanes == 1) { + if (curvature == 1) return 1369.532383465354; + + if (curvature == 2) return 1117.1498589355958; + + if (curvature == 3) return 1048.5840399296935; + + if (curvature == 4) return 956.0314100959505; + } + + if (lanes == 2) return 1956.6719999999998; + + + // Own assumption of increasing capacity with more lanes + // This is not covered by the HBS and is a very rare case + return (1956.6719999999998 * 1.3) / lanes; + } + + /** + * Bundesstraße with at least 70km/h, often osm primary or trunk + */ + private static double capacityBundesStr(int lanes) { + + if (lanes == 1) + return 2033.868926820213; + + if (lanes == 2) + return 3902.4390243902435 / 2; + + return (3902.4390243902435 * 1.3) / lanes; + } + + /** + * Capacity on a side road merging into a main road at a junction. + * + * @param qP the vehicle volume of the main road + */ + private static double capacityMerging(double qP) { + + // See HBS page 5-20, eq. S5-12 table S5-5 + + // mean Folgezeitlücken of the different combinations + double tf = 3.5; + + // mean Grenzzeitlücke + double tg = 6.36; + + + return Math.exp((-qP / 3600) * (tg - tf / 2)) * 3600 / tf; + } + + private static int curvatureCategory(double length, double curvature) { + + // for too short segment, curvature is not relevant + if (length < 50 || curvature == 0) + return 1; + + double sumChanges = curvature * (length / 1000); + + // Scale length of segment to at least 300m, because the manual recommends at least a certain segment length + double ku = sumChanges / Math.max(0.3, length / 1000); + + if (ku > 150) + return 4; + if (ku > 100) + return 3; + if (ku > 50) + return 2; + return 1; + + } + + /** + * Capacity of a side road of type merging into a main road. + * + * @param roadType type of the target road + * @param mainType type of the higher priority road + */ + private static double capacityMerging(String roadType, String mainType) { + + if (mainType.equals("trunk") || mainType.equals("primary") || roadType.contains("residential")) + // ~600 veh/h + return capacityMerging(600); + + // ~800 veh/h + return capacityMerging(400); + } + + @Override + public double predict(Object2DoubleMap features, Object2ObjectMap categories) { + + // Speed in km/h + int speed = (int) Math.round(features.getDouble("speed") * 3.6); + int lanes = (int) features.getOrDefault("num_lanes", 1); + int curvature = curvatureCategory(features.getDouble("length"), features.getOrDefault("curvature", 0)); + String type = categories.get("highway_type"); + + // Primary and trunk roads are often Bundesstraßen, + // but only if they have a speed limit of 70 or more, this calculation is valid + // From OSM alone it is not possible to distinguish anbaufreie Hauptverkehrsstraßen and Landstraße clearly + if (speed >= 90 || ((type.contains("primary") || type.contains("trunk")) && speed >= 70)) { + return capacityBundesStr(lanes); + } else if (speed >= 70) { + return capacityLandStr(lanes, curvature); + } + + String merging = categories.get("is_merging_into"); + + // Only merging with a single lane road is considered + if (!merging.isEmpty() && lanes == 1) { + return capacityMerging(type, merging); + } + + // Capacity for city roads + if (speed >= 40 || lanes >= 2 || features.getDouble("is_secondary_or_higher") == 1) { + return switch (lanes) { + case 1 -> 1139.0625; + case 2 -> 2263.438914027149 / 2; + // Own assumption, rare edge-case + default -> 2263.438914027149 * 1.3 / lanes; + }; + } + + // Remaining are residential which are assumed to have high urbanisation + return 800.0 / lanes; + } + +} diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ref/DecisionTreeParams.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ref/DecisionTreeParams.java index 4bb84fa6577..7dab8d49135 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ref/DecisionTreeParams.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ref/DecisionTreeParams.java @@ -18,7 +18,7 @@ public final class DecisionTreeParams implements NetworkModel { private static final FeatureRegressor INSTANCE = new Model(); @Override - public FeatureRegressor speedFactor(String junctionType) { + public FeatureRegressor speedFactor(String junctionType, String highwayType) { return INSTANCE; } diff --git a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ref/GermanyNetworkParams.java b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ref/GermanyNetworkParams.java index ffe7250b2f4..80541a82a7e 100644 --- a/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ref/GermanyNetworkParams.java +++ b/contribs/application/src/main/java/org/matsim/application/prepare/network/params/ref/GermanyNetworkParams.java @@ -9,7 +9,7 @@ */ public final class GermanyNetworkParams implements NetworkModel { @Override - public FeatureRegressor capacity(String junctionType) { + public FeatureRegressor capacity(String junctionType, String highwayType) { return switch (junctionType) { case "traffic_light" -> GermanyNetworkParams_capacity_traffic_light.INSTANCE; case "right_before_left" -> GermanyNetworkParams_capacity_right_before_left.INSTANCE; @@ -19,7 +19,7 @@ public FeatureRegressor capacity(String junctionType) { } @Override - public FeatureRegressor speedFactor(String junctionType) { + public FeatureRegressor speedFactor(String junctionType, String highwayType) { return switch (junctionType) { case "traffic_light" -> GermanyNetworkParams_speedRelative_traffic_light.INSTANCE; case "right_before_left" -> GermanyNetworkParams_speedRelative_right_before_left.INSTANCE; diff --git a/contribs/application/src/main/python/capacity/hbs.py b/contribs/application/src/main/python/capacity/hbs.py new file mode 100644 index 00000000000..aad41f1a45e --- /dev/null +++ b/contribs/application/src/main/python/capacity/hbs.py @@ -0,0 +1,234 @@ +import math +import matplotlib.pyplot as plt +import plotly.graph_objects as go + + +def capacity_estimate(v): + # headway + tT = 1.2 + # car length + lL = 7.0 + Qc = v / (v * tT + lL) + return 3600 * Qc + + +def merge(qp=400): + tf = 3.5 + tg = 6.36 + return math.exp((-qp / 3600) * (tg - tf / 2)) * 3600 / tf + + +class Street: + def __init__(self, osm_type, lanes, speed, curvature): + self.osm_type = osm_type + self.lanes = lanes + self.speed = speed + self.curvature = curvature + + +def calc_capacity_stadtstrasse(street): + # Source: HSB S3 + k = 0 + a = 54 + b = 0.850 + + if street.lanes == 1: + f = 1.0 + # Assume middle-high "Erschließungsintensität" + if street.speed == 30: + k = 45 + a = 38 + b = 0.715 + elif street.speed == 50: + k = 45 + a = 54 + b = 0.850 + elif street.speed == 70: + k = 40 + a = 89 + b = 0.846 + + return (-1.0 * b * math.pow(k, 3 / 2) * math.sqrt(4 * a * f + b * b * k) + 2 * a * f * k + b * b * k * k) / ( + 2 * f * f) + elif street.lanes == 2: + f = 0.7 + if street.speed == 50: + k = 45 + a = -0.009 + b = 55.58 + elif street.speed == 70: + f = 0.5 + k = 40 + a = -0.008 + b = 80 + + return (b * k) / (f - a * k) + + +def calc_capacity_landstrasse(street): + # Source: HSB table L3-4 + if street.lanes == 1: + k = 20 + + a = 0 + b = 0 + + # Table L3-4 + if street.curvature == 1: + a = 98.73 + b = 0.8175 + elif street.curvature == 2: + a = 83.88 + b = 0.8384 + elif street.curvature == 3: + a = 74.41 + b = 0.6788 + elif street.curvature == 4: + a = 68.02 + b = 0.6539 + else: + raise ValueError(f"Unknown curvature {street.curvature}") + + return 0.5 * (-b * math.pow(k, 3 / 2) * math.sqrt(4 * a + b * b * k) + 2 * a * k + b * b * k * k) + if street.lanes == 2: + k = 48 + a = 55.5 + b = -0.614 + + # Divide by 2 because the formula is for both directions + return k * (2 * a + b * k) / 2 + + +def calc_capacity_autobahn(street): + # Source: HSB A3 + anzahl_fahrstreifen = 1 + if street.lanes <= 2: + if street.speed >= 130: + return 3700 / anzahl_fahrstreifen + elif street.speed == 120: + return 3800 / anzahl_fahrstreifen + elif street.speed <= 120: + return 3750 / anzahl_fahrstreifen + elif street.lanes == 3: + if street.speed >= 130: + return 5300 / anzahl_fahrstreifen + elif street.speed == 120: + return 5400 / anzahl_fahrstreifen + elif street.speed <= 120: + return 5350 / anzahl_fahrstreifen + else: + if street.speed >= 130: + return 7300 / anzahl_fahrstreifen + elif street.speed == 120: + return 7400 / anzahl_fahrstreifen + elif street.speed <= 120: + return 7400 / anzahl_fahrstreifen + + +osm_types = { + 'residential': 'Stadtstraßen', + 'unclassified': 'Landstraßen', + 'motorway': 'Autobahnen' +} + +speeds = [30, 50, 70, 80, 90, 100, 120, 130] +curvatures = [1, 2, 3, 4] + +capacities = {} +streets = [] + +for osm_type, street_type in osm_types.items(): + capacities[osm_type] = {} + for lanes in range(1, 5): + capacities[osm_type][lanes] = {} + for curvature in curvatures: + capacities[osm_type][lanes][curvature] = {} + for speed in speeds: + if (osm_type == 'residential' and speed > 70) or (osm_type == 'residential' and speed < 30) or ( + osm_type == 'residential' and lanes > 2): + continue + if (osm_type == 'unclassified' and speed > 100) or (osm_type == 'unclassified' and speed <= 50) or ( + osm_type == 'unclassified' and lanes > 2): + continue + if (osm_type == 'motorway' and lanes == 1) or (osm_type == 'motorway' and speed < 80): + continue + street = Street(osm_type, lanes, speed, curvature) + streets.append(street) + capacities[osm_type][lanes][curvature][speed] = 0 + +for street in streets: + capacity = 0 + if street.osm_type == 'residential': + capacity = calc_capacity_stadtstrasse(street) + elif street.osm_type == 'motorway': + capacity = calc_capacity_autobahn(street) + else: + capacity = calc_capacity_landstrasse(street) + capacities[street.osm_type][street.lanes][street.curvature][street.speed] = capacity + +for street_type, lanes_capacity in capacities.items(): + for lanes, speeds_capacity in lanes_capacity.items(): + print(f"{osm_types[street_type]}, Fahrstreifen: {lanes}:") + for curvature, capacity in speeds_capacity.items(): + for speed, capacity in capacity.items(): + print(f"Geschwindigkeit {speed} km/h: Kurvigkeit {curvature} Verkehrskapazität {capacity} ") + +plt.figure(figsize=(14, 10)) +for street_type, lanes_capacity in capacities.items(): + for lanes, speeds_capacity in lanes_capacity.items(): + for curvature, capacity in speeds_capacity.items(): + plt.plot(list(capacity.keys()), list(capacity.values()), + label=f"{osm_types[street_type]}, Fahrstreifen: {lanes} Kurvigkeit: {curvature}") + +# plt.figure(figsize=(14, 10)) +# for street_type, lanes_capacity in capacities.items(): +# for lanes, speeds_capacity in lanes_capacity.items(): +# for speed, curvatures_capacity in speeds_capacity.items(): +# for curvature, capacity in curvatures_capacity.items(): +# plt.plot(speed, capacity, +# label=f"{osm_types[street_type]}, Fahrstreifen: {lanes}, Kurvigkeit: {curvature}") + +speeds_compare = [] +capacity_compare = [] +for v in range(30, 130, 10): + speeds_compare.append(v) + capacity_compare.append(capacity_estimate(v)) +plt.plot(speeds_compare, capacity_compare, label="Referenz") + +plt.xlabel('Geschwindigkeit (km/h)') +plt.ylabel('Verkehrskapazität') +plt.title('Verkehrskapazität nach Straßentyp, Fahrstreifen und Geschwindigkeit') +plt.legend(loc='upper center', bbox_to_anchor=(0.5, -0.2), fancybox=True, ncol=4) +plt.subplots_adjust(bottom=0.3) +plt.show() + +# fig = go.Figure() +# +# for street_type, lanes_capacity in capacities.items(): +# for lanes, speeds_capacity in lanes_capacity.items(): +# fig.add_trace(go.Scatter(x=list(speeds_capacity.keys()), y=list(speeds_capacity.values()), +# mode='lines', +# name=f"{osm_types[street_type]}, Fahrstreifen: {lanes}")) +# +# speeds_compare = [] +# capacity_compare = [] +# for v in range(50, 140, 10): +# speeds_compare.append(v) +# capacity_compare.append(capacity_estimate(v)) +# fig.add_trace(go.Scatter(x=speeds_compare, y=capacity_compare, mode='lines', showlegend=True, name="REFERENZ")) +# +# for qp in (0, 200, 400, 600): +# kontenpunkt_y = [] +# kontenpunkt_y.append(merge(qp)) +# kontenpunkt_y.append(merge(qp)) +# fig.add_trace(go.Scatter(x=(30, 50), y=kontenpunkt_y, mode='lines', showlegend=True, name="Knotenpunkt %d" % qp)) +# +# fig.update_layout( +# title='Verkehrskapazität nach Straßentyp, Fahrstreifen und Geschwindigkeit', +# xaxis_title='Geschwindigkeit (km/h)', +# yaxis_title='Verkehrskapazität', +# legend=dict(x=0.5, y=-0.2, orientation='h'), +# margin=dict(b=100) +# ) +# +# fig.show() diff --git a/contribs/application/src/test/java/org/matsim/application/prepare/network/params/ApplyNetworkParamsTest.java b/contribs/application/src/test/java/org/matsim/application/prepare/network/params/ApplyNetworkParamsTest.java index b25737d6c40..e9f723d6db8 100644 --- a/contribs/application/src/test/java/org/matsim/application/prepare/network/params/ApplyNetworkParamsTest.java +++ b/contribs/application/src/test/java/org/matsim/application/prepare/network/params/ApplyNetworkParamsTest.java @@ -40,6 +40,5 @@ void apply() throws Exception { ); assertThat(output.resolve("network-opt.xml")).exists(); - } } diff --git a/contribs/application/src/test/java/org/matsim/application/prepare/network/params/HBSNetworkParamsTest.java b/contribs/application/src/test/java/org/matsim/application/prepare/network/params/HBSNetworkParamsTest.java new file mode 100644 index 00000000000..ab6172db208 --- /dev/null +++ b/contribs/application/src/test/java/org/matsim/application/prepare/network/params/HBSNetworkParamsTest.java @@ -0,0 +1,44 @@ +package org.matsim.application.prepare.network.params; + +import org.junit.jupiter.api.Test; +import org.junit.jupiter.api.extension.RegisterExtension; +import org.matsim.contrib.sumo.SumoNetworkConverter; +import org.matsim.testcases.MatsimTestUtils; + +import java.nio.file.Path; +import java.util.List; + +import static org.assertj.core.api.Assertions.assertThat; + +class HBSNetworkParamsTest { + + @RegisterExtension + MatsimTestUtils utils = new MatsimTestUtils(); + + @Test + void apply() throws Exception { + + Path networkPath = Path.of(utils.getPackageInputDirectory()).resolve("osm.net.xml"); + + Path output = Path.of(utils.getOutputDirectory()); + + SumoNetworkConverter converter = SumoNetworkConverter.newInstance(List.of(networkPath), + output.resolve("network.xml"), + "EPSG:4326", "EPSG:4326"); + + converter.call(); + + assertThat(output.resolve("network.xml")).exists(); + assertThat(output.resolve("network-ft.csv")).exists(); + + new ApplyNetworkParams().execute( + "capacity", + "--network", output.resolve("network.xml").toString(), + "--input-features", output.resolve("network-ft.csv").toString(), + "--output", output.resolve("network-hbs.xml").toString(), + "--model", "org.matsim.application.prepare.network.params.hbs.HBSNetworkParams" + ); + + assertThat(output.resolve("network-hbs.xml")).exists(); + } +} diff --git a/contribs/sumo/src/main/java/org/matsim/contrib/sumo/SumoNetworkFeatureExtractor.java b/contribs/sumo/src/main/java/org/matsim/contrib/sumo/SumoNetworkFeatureExtractor.java index 7bf1154bf2d..ee8c9a45b87 100644 --- a/contribs/sumo/src/main/java/org/matsim/contrib/sumo/SumoNetworkFeatureExtractor.java +++ b/contribs/sumo/src/main/java/org/matsim/contrib/sumo/SumoNetworkFeatureExtractor.java @@ -1,5 +1,7 @@ package org.matsim.contrib.sumo; +import it.unimi.dsi.fastutil.objects.Object2IntMap; +import it.unimi.dsi.fastutil.objects.Object2IntOpenHashMap; import org.apache.commons.csv.CSVPrinter; import org.matsim.contrib.osm.networkReader.LinkProperties; @@ -28,7 +30,7 @@ class SumoNetworkFeatureExtractor { incomingEdges = new HashMap<>(); for (SumoNetworkHandler.Edge edge : this.handler.edges.values()) { - incomingEdges.computeIfAbsent(edge.to, (k) -> new ArrayList<>()) + incomingEdges.computeIfAbsent(edge.to, k -> new ArrayList<>()) .add(edge); } } @@ -77,6 +79,60 @@ private static Set directionSet(SumoNetworkHandler.Connection c) { return set; } + /** + * Calculate the curvature of an edge. One gon is 1/400 of a full circle. + * The formula is: KU = (Sum of the curvature of the subsegments) / (Length of the edge) + * + * @return curvature in gon/km + */ + static double calcCurvature(SumoNetworkHandler.Edge edge) { + double totalGon = 0; + List coordinates = edge.shape; + + for (int i = 2; i < coordinates.size(); i++) { + double[] pointA = coordinates.get(i - 2); + double[] pointB = coordinates.get(i - 1); + double[] pointC = coordinates.get(i); + + double[] vectorAB = {pointB[0] - pointA[0], pointB[1] - pointA[1]}; + double[] vectorBC = {pointC[0] - pointB[0], pointC[1] - pointB[1]}; + + double dotProduct = calcDotProduct(vectorAB, vectorBC); + double magnitudeAB = calcMagnitude(vectorAB); + double magnitudeBC = calcMagnitude(vectorBC); + + double cosine = dotProduct / (magnitudeAB * magnitudeBC); + double angleRadians = Math.acos(cosine); + double angleDegrees = Math.toDegrees(angleRadians); + + totalGon += Math.abs((angleDegrees / 360) * 400); + } + + return totalGon / (edge.getLength() / 1000); + } + + /** + * Calculate the dot product of two vectors. + * + * @param vec1 vector 1 + * @param vec2 vector 2 + * @return dot product + */ + private static double calcDotProduct(double[] vec1, double[] vec2) { + return vec1[0] * vec2[0] + vec1[1] * vec2[1]; + } + + /** + * Calculate the magnitude of a vector. + * + * @param vec vector + * @return magnitude + */ + private static double calcMagnitude(double[] vec) { + return Math.sqrt(vec[0] * vec[0] + vec[1] * vec[1]); + } + + /** * Get priority. Higher is more important. */ @@ -86,9 +142,12 @@ private int getPrio(SumoNetworkHandler.Edge edge) { public List getHeader() { return List.of("linkId", "highway_type", "speed", "length", "num_lanes", "change_num_lanes", "change_speed", "num_to_links", "num_conns", - "num_response", "num_foes", "dir_multiple_s", "dir_l", "dir_r", "dir_s", "dir_exclusive", + "num_response", "num_foes", "dir_multiple_s", "dir_l", "dir_r", "dir_s", "dir_exclusive", "curvature", "junction_type", "junction_inc_lanes", "priority_higher", "priority_equal", "priority_lower", - "is_secondary_or_higher", "is_primary_or_higher", "is_motorway", "is_link"); + "is_secondary_or_higher", "is_primary_or_higher", "is_motorway", + "is_link", "has_merging_link", "is_merging_into", + "num_left", "num_right", "num_straight" + ); } public void print(CSVPrinter out) { @@ -106,7 +165,7 @@ public void print(CSVPrinter out, String linkId, SumoNetworkHandler.Edge edge) t String highwayType = getHighwayType(edge.type); SumoNetworkHandler.Junction junction = handler.junctions.get(edge.to); - List connections = handler.connections.computeIfAbsent(edge.id, (k) -> new ArrayList<>()); + List connections = handler.connections.computeIfAbsent(edge.id, k -> new ArrayList<>()); Set toEdges = connections.stream() .filter(c -> !c.dir.equals("t")) @@ -129,12 +188,18 @@ public void print(CSVPrinter out, String linkId, SumoNetworkHandler.Edge edge) t Set dirs = new HashSet<>(); boolean multipleDirS = false; boolean exclusiveDirs = true; + + // Number of connections per direction + Object2IntMap numConnections = new Object2IntOpenHashMap<>(); + for (SumoNetworkHandler.Connection c : connections) { Set d = directionSet(c); if (dirs.contains('s') && d.contains('s')) multipleDirS = true; + d.forEach(dir -> numConnections.mergeInt(dir, 1, Integer::sum)); + Set intersection = new HashSet<>(dirs); // use the copy constructor intersection.retainAll(d); if (!intersection.isEmpty()) @@ -152,12 +217,27 @@ public void print(CSVPrinter out, String linkId, SumoNetworkHandler.Edge edge) t .mapToInt(e -> e.lanes.size()) .sum(); - boolean geq_secondary = switch (highwayType) { + boolean merging = incomingEdges.get(junction.id).stream() + .anyMatch(e -> e.type.contains("link")); + + OptionalInt highestPrio = incomingEdges.get(junction.id).stream() + .mapToInt(this::getPrio).max(); + + // Find category of the highest merging lane + String mergingHighest = ""; + if (highestPrio.isPresent() && highestPrio.getAsInt() > getPrio(edge)) { + Optional> m = osm.entrySet().stream().filter(e -> e.getValue().getHierarchyLevel() == -highestPrio.getAsInt()) + .findFirst(); + if (m.isPresent()) + mergingHighest = m.get().getKey(); + } + + boolean geqSecondary = switch (highwayType) { case "secondary", "primary", "trunk", "motorway" -> true; default -> false; }; - boolean geq_primary = switch (highwayType) { + boolean geqPrimary = switch (highwayType) { case "primary", "trunk", "motorway" -> true; default -> false; }; @@ -170,7 +250,7 @@ public void print(CSVPrinter out, String linkId, SumoNetworkHandler.Edge edge) t out.print(Math.max(-3, Math.min(3, maxLanes - edge.lanes.size()))); out.print(maxSpeed - edge.lanes.get(0).speed); out.print(toEdges.size()); - out.print(Math.min(6, handler.connections.get(edge.id).size())); + out.print(Math.min(6, connections.size())); out.print(Math.min(12, aggr.response().cardinality())); out.print(Math.min(12, aggr.foes().cardinality())); out.print(bool(multipleDirS)); @@ -178,15 +258,21 @@ public void print(CSVPrinter out, String linkId, SumoNetworkHandler.Edge edge) t out.print(bool(dirs.contains('r'))); out.print(bool(dirs.contains('s'))); out.print(bool(exclusiveDirs)); + out.print(calcCurvature(edge)); out.print(junction.type); out.print(Math.min(12, incomingLanes)); out.print(bool("higher".equals(prio))); out.print(bool("equal".equals(prio))); out.print(bool("lower".equals(prio))); - out.print(bool(geq_secondary)); - out.print(bool(geq_primary)); + out.print(bool(geqSecondary)); + out.print(bool(geqPrimary)); out.print(bool("motorway".equals(highwayType))); out.print(bool(highwayType.contains("link"))); + out.print(bool(merging)); + out.print(mergingHighest); + out.print(numConnections.getInt('l')); + out.print(numConnections.getInt('r')); + out.print(numConnections.getInt('s')); out.println(); } diff --git a/contribs/sumo/src/test/java/org/matsim/contrib/sumo/SumoNetworkFeatureExtractorTest.java b/contribs/sumo/src/test/java/org/matsim/contrib/sumo/SumoNetworkFeatureExtractorTest.java new file mode 100644 index 00000000000..e8b41e359da --- /dev/null +++ b/contribs/sumo/src/test/java/org/matsim/contrib/sumo/SumoNetworkFeatureExtractorTest.java @@ -0,0 +1,81 @@ +package org.matsim.contrib.sumo; + +import org.junit.jupiter.api.Assertions; +import org.junit.jupiter.api.Test; + +class SumoNetworkFeatureExtractorTest { + + @Test + void twoCoordsBackAndForth() { + String[] coords = {"0,0", "0,100", "0,0"}; + int lentgh = 200; + + SumoNetworkHandler.Edge edge = new SumoNetworkHandler.Edge("1", "2", "3", "highway", 0, "name", coords); + edge.lanes.add(new SumoNetworkHandler.Lane("1", 0, lentgh, 50 / 3.6, null, null)); + + double ku = SumoNetworkFeatureExtractor.calcCurvature(edge); + + // Length: 0.2 km + // Gon: 200 + // Curvature: 200 / 0.2 = 1000 + + Assertions.assertEquals(1000, ku, 0.000001); + + System.out.println(ku); + } + + + @Test + void nintyDegreeCorner() { + + String[] coords = {"0,0", "0,100", "100,100"}; + int lentgh = 200; + + SumoNetworkHandler.Edge edge = new SumoNetworkHandler.Edge("1", "2", "3", "highway", 0, "name", coords); + edge.lanes.add(new SumoNetworkHandler.Lane("1", 0, lentgh, 50 / 3.6, null, null)); + + double ku = SumoNetworkFeatureExtractor.calcCurvature(edge); + + // Length: 0.2 km + // Gon: 100 + // Curvature: 100 / 0.2 = 500 + + Assertions.assertEquals(500, ku, 0.000001); + } + + @Test + void twoCorners() { + + String[] coords = {"0,0", "0,100", "100,100", "100,200"}; + int lentgh = 300; + + SumoNetworkHandler.Edge edge = new SumoNetworkHandler.Edge("1", "2", "3", "highway", 0, "name", coords); + edge.lanes.add(new SumoNetworkHandler.Lane("1", 0, lentgh, 50 / 3.6, null, null)); + + double ku = SumoNetworkFeatureExtractor.calcCurvature(edge); + + // Length: 0.3 km + // Gon: 100 + 100 = 200 + // Curvature: 200 / 0.3 = 666.6666666666666 + + Assertions.assertEquals((200 / 0.3), ku, 0.000001); + } + + @Test + void rectangle() { + + String[] coords = {"0,0", "0,100", "100,100", "100,0", "0,0"}; + int lentgh = 400; + + SumoNetworkHandler.Edge edge = new SumoNetworkHandler.Edge("1", "2", "3", "highway", 0, "name", coords); + edge.lanes.add(new SumoNetworkHandler.Lane("1", 0, lentgh, 50 / 3.6, null, null)); + + double ku = SumoNetworkFeatureExtractor.calcCurvature(edge); + + // Length: 0.4 km + // Gon: 100 + 100 + 100 = 300 + // Curvature: 300 / 0.4 = 666.6666666666666 + + Assertions.assertEquals((300 / 0.4), ku, 0.000001); + } +}