Skip to content

Commit

Permalink
Merge pull request #2 from matsim-scenarios/commuter-analysis
Browse files Browse the repository at this point in the history
Commuter analysis
  • Loading branch information
simei94 authored Sep 4, 2024
2 parents 675f45f + 99c0bb0 commit c9bdfca
Show file tree
Hide file tree
Showing 5 changed files with 271 additions and 7 deletions.
52 changes: 52 additions & 0 deletions src/main/R/commuter-analysis.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,52 @@
library(matsim)
library(tidyverse)
library(readr)
library(sf)

#FILES
FILE_DIR = "../../shared-svn/projects/DiTriMo/data/commuters-by-town"
SIM <- paste0(FILE_DIR, "lausitz-v1.0-commuter.csv")
GEMEINDE <- paste0(FILE_DIR, "pgemeinden.csv")
COMMUTER <- paste0(FILE_DIR, "commuter.csv")
SHP <- paste0(FILE_DIR, "VG5000_GEM/VG5000_GEM.shp")
LAUSITZ.SHP <- paste0(FILE_DIR, "network-area/network-area.shp")

shp <- st_read(SHP)
lausitz.shp <- st_read(LAUSITZ.SHP)

sim <- read_csv(file = SIM)
gemeinden <- read_csv( file = GEMEINDE) %>%
mutate(code = str_remove(string = code, pattern = "P"))

commuter <- read_csv(file = COMMUTER) %>%
mutate(key = paste0(from, "-", to))

sim.1 <- sim %>%
filter(!is.na(from) & !is.na(to)) %>%
left_join(gemeinden, by = c("from" = "code")) %>%
left_join(gemeinden, by = c("from" = "code"), suffix = c("_from", "_to")) %>%
mutate(key = paste0(from, "-", to)) %>%
select(-c(from, to)) %>%
left_join(commuter, by = "key", suffix = c("_sim", "_real")) %>%
filter(!is.na(from)) %>%
# pivot_longer(cols = starts_with("n_"), names_to = "src", values_to = "n", names_prefix = "n_") %>%
arrange(desc(n_real))

breaks <- c(-Inf, 0.8, 1.2, Inf)
labels <- c("less", "exakt", "more")

sim.2 <- sim.1 %>%
filter(n_sim > 10) %>%
select(from, to, starts_with("n_"), starts_with("name_")) %>%
mutate(n_rel = n_sim / n_real,
quality = cut(n_rel, breaks = breaks, labels = labels))

ggplot(sim.2, aes(x = n_real, y = n_sim, col = quality)) +

geom_point() +

scale_x_log10() +

scale_y_log10() +

theme_bw()
34 changes: 28 additions & 6 deletions src/main/R/counts.R
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,10 @@ devtools::install_github("matsim-vsp/matsim-r",ref="counts")
library(matsim)
library(tidyverse)

COUNTS <- "Y:/matsim-lausitz/input/v1.0/lausitz-v1.0-counts-car-bast.xml.gz"
NETWORK <- "Y:/matsim-lausitz/input/v1.0/lausitz-v1.0-network-with-pt.xml.gz"
COUNTS <- "../../public-svn/matsim/scenarios/countries/de/lausitz/input/v1.0/lausitz-v1.0-counts-car-bast.xml.gz"
NETWORK <- "../../public-svn/matsim/scenarios/countries/de/lausitz/input/v1.0/lausitz-v1.0-network-with-pt.xml.gz"

linkstats <- readLinkStats(runId = "v1.0-uncalibrated", file = "Y:/matsim-lausitz/qa/output/lausitz-25pct.output_linkstats.csv.gz")
linkstats <- readLinkStats(runId = "v1.0-uncalibrated", file = "Y:/matsim-lausitz/qa/output/lausitz-25pct.output_linkstats.csv.gz", sampleSize = 1)

counts <- readCounts(COUNTS)
network <- loadNetwork(NETWORK)
Expand All @@ -17,7 +17,7 @@ join <- mergeCountsAndLinks(counts = counts, linkStats = list(linkstats), netwo

#### VIA-styled scatterplot ####

FILE_DIR = "C:/Users/ACER/Desktop/Uni/VSP/Lausitz-Plots/"
FILE_DIR = "../../shared-svn/projects/DiTriMo/data/commuters-by-town"

createCountScatterPlot(joinedFrame = join)
ggsave(filename = paste0(FILE_DIR, "Traffic_Count_Scatterplot_with_freight.jpg"))
Expand All @@ -43,7 +43,7 @@ rm(join.dtv.distribution)

#### Analysis of Estimation Quality ####

join.est.quality <- processDtvEstimationQuality(joinedFrame = join, aggr = T) %>%
join.est.quality <- processDtvEstimationQuality(joinedFrame = join, aggr = F) %>%
filter(!type %in% c("residential", "unclassified", NA))

ggplot(join.est.quality, aes(estimation, share, fill = type)) +
Expand All @@ -59,4 +59,26 @@ ggplot(join.est.quality, aes(estimation, share, fill = type)) +
theme(legend.position = "none", axis.text.x = element_text(angle = 90))

rm(join.est.quality)
ggsave(filename = paste0(FILE_DIR, "Estimation_quality_by_road_type_with_freight.jpg"))
ggsave(filename = paste0(FILE_DIR, "Estimation_quality_by_road_type_with_freight.jpg"))


#### network plot ####

library(tmap)
library(tmaptools)
library(OpenStreetMap)
library(sf)

link.geom <- join %>%
left_join(network$links, by = c("loc_id" = "id")) %>%
mutate(geom = sprintf("LINESTRING(%s %s, %s %s)", x.from, y.from, x.to, y.to)) %>%
st_as_sf(crs = 25832, wkt = "geom") %>%
transmute(loc_id, type.x, rel_vol = volume / count, geom)

tmap_mode("view")

tm_shape(shp = link.geom) +
tm_lines(col = "estimation", style = "cont", lwd = 3.5, palette = c("red", "green", "blue"))

tm_shape(shp = link.geom) +
tm_lines(col = "rel_vol", style = "cont", lwd = 5, palette = c("red", "yellow", "green"), breaks = c(0, 0.05, 0.8, 2))
4 changes: 3 additions & 1 deletion src/main/java/org/matsim/run/LausitzScenario.java
Original file line number Diff line number Diff line change
Expand Up @@ -34,7 +34,9 @@
import org.matsim.core.population.PopulationUtils;
import org.matsim.core.replanning.strategies.DefaultPlanStrategiesModule;
import org.matsim.core.scoring.functions.ScoringParametersForPerson;
import org.matsim.run.analysis.CommunityFilter;
import org.matsim.run.analysis.CommuterAnalysis;
import org.matsim.run.analysis.DistanceMatrix;
import org.matsim.run.prepare.PrepareNetwork;
import org.matsim.run.prepare.PreparePopulation;
import org.matsim.simwrapper.SimWrapperConfigGroup;
Expand All @@ -59,7 +61,7 @@
SplitActivityTypesDuration.class, CreateCountsFromBAStData.class, PreparePopulation.class, CleanPopulation.class, PrepareNetwork.class
})
@MATSimApplication.Analysis({
LinkStats.class, CheckPopulation.class, CommuterAnalysis.class,
LinkStats.class, CheckPopulation.class, CommuterAnalysis.class, CommunityFilter.class, DistanceMatrix.class
})
public class LausitzScenario extends MATSimApplication {

Expand Down
80 changes: 80 additions & 0 deletions src/main/java/org/matsim/run/analysis/CommunityFilter.java
Original file line number Diff line number Diff line change
@@ -0,0 +1,80 @@
package org.matsim.run.analysis;

import org.apache.commons.csv.CSVPrinter;
import org.geotools.api.feature.simple.SimpleFeature;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.Point;
import org.matsim.application.MATSimAppCommand;
import org.matsim.application.options.CsvOptions;
import org.matsim.core.utils.collections.Tuple;
import org.matsim.core.utils.gis.GeoFileReader;
import picocli.CommandLine;

import java.nio.file.Path;
import java.util.*;

@CommandLine.Command(name = "community-filter", description = "creates a csv with commuity keys within shape")
public class CommunityFilter implements MATSimAppCommand {

@CommandLine.Option(names = "--community-shp", description = "path to VG5000_GEM.shp", required = true)
private static Path communityShapePath;

@CommandLine.Option(names = "--dilution-area", description = "path to area-of-interest-shape file", required = true)
private static Path dilutionAreaShapePath;

@CommandLine.Option(names = "--output", description = "output csv filepath", required = true)
private static Path output;

@CommandLine.Mixin
CsvOptions csvOptions = new CsvOptions();
private final Map<String, Tuple<Double, Double>> filtered = new HashMap<>();

public static void main(String[] args) {
new CommunityFilter().execute(args);
}

@Override
public Integer call() throws Exception {

Collection<SimpleFeature> communities = readShapeFile(communityShapePath);
Collection<SimpleFeature> dilutionArea = readShapeFile(dilutionAreaShapePath);

List<Geometry> geometries = dilutionArea.stream().map(feature -> (Geometry) feature.getDefaultGeometry()).toList();

for(var community: communities){

String attribute = (String) community.getAttribute("ARS");

if(geometries.stream()
.anyMatch(geometry -> geometry.covers((Geometry) community.getDefaultGeometry()))){

Point centroid = ((Geometry) community.getDefaultGeometry()).getCentroid();
Tuple<Double, Double> coordinates = Tuple.of(centroid.getX(), centroid.getY());

filtered.put(attribute, coordinates);
}

}

try (CSVPrinter printer = csvOptions.createPrinter(output)) {
printer.print("ars");
printer.print("x");
printer.print("y");
printer.println();

for (Map.Entry<String, Tuple<Double, Double>> entry : filtered.entrySet()) {
printer.print(entry.getKey());
printer.print(entry.getValue().getFirst());
printer.print(entry.getValue().getSecond());
printer.println();
}
}

return 0;
}

private static Collection<SimpleFeature> readShapeFile(Path filepath){

return GeoFileReader.getAllFeatures(filepath.toString());
}
}
108 changes: 108 additions & 0 deletions src/main/java/org/matsim/run/analysis/DistanceMatrix.java
Original file line number Diff line number Diff line change
@@ -0,0 +1,108 @@
package org.matsim.run.analysis;

import org.apache.commons.csv.CSVPrinter;
import org.apache.logging.log4j.LogManager;
import org.apache.logging.log4j.Logger;
import org.geotools.api.feature.simple.SimpleFeature;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.Point;
import org.matsim.api.core.v01.Coord;
import org.matsim.application.MATSimAppCommand;
import org.matsim.application.options.CsvOptions;
import org.matsim.application.options.ShpOptions;
import org.matsim.core.utils.geometry.CoordUtils;
import org.matsim.core.utils.geometry.geotools.MGC;
import org.matsim.core.utils.gis.GeoFileReader;
import picocli.CommandLine;

import java.nio.file.Path;
import java.util.ArrayList;
import java.util.List;
import java.util.function.Predicate;

@CommandLine.Command(name = "distance-matrix", description = "creates a csv with commuity keys within shape")
public class DistanceMatrix implements MATSimAppCommand{

@CommandLine.Option(names = "--output", description = "output csv filepath", required = true)
private static Path output;

@CommandLine.Option(names = "--dilution-area", description = "shape to filter zones", required = false)
private static Path dilutionArea;

@CommandLine.Mixin
CsvOptions csvOptions = new CsvOptions();

@CommandLine.Mixin
ShpOptions shp = new ShpOptions();

private final List<String> distances = new ArrayList<>();
private static final Logger logger = LogManager.getLogger(DistanceMatrix.class);

public static void main(String[] args) {
new DistanceMatrix().execute(args);
}

@Override
public Integer call() throws Exception {

logger.info("Read features.");
List<SimpleFeature> communities = shp.readFeatures();

//to prevent RuntimeExceptions
ArrayList<SimpleFeature> copy = new ArrayList<>(communities);

Predicate<Geometry> filter = getFilter(dilutionArea);

logger.info("Calculate distance matrix.");
String delimiter = csvOptions.getFormat().getDelimiterString();
for(var community: communities){

if(!filter.test((Geometry) community.getDefaultGeometry()))
continue;
String nameFrom = (String) community.getAttribute("ARS");
Point centroid = ((Geometry) community.getDefaultGeometry()).getCentroid();
Coord from = MGC.point2Coord(centroid);
for(var target: copy){

String nameTo = (String) target.getAttribute("ARS");

Point centroid2 = ((Geometry) target.getDefaultGeometry()).getCentroid();
Coord to = MGC.point2Coord(centroid2);
double distance = CoordUtils.calcEuclideanDistance(from, to);

String distanceString = String.valueOf(distance).replace('.', ',');

distances.add(nameFrom + delimiter + nameTo + delimiter + distanceString);
}
}

logger.info("Print results to {}", output);
try (CSVPrinter printer = csvOptions.createPrinter(output)) {
printer.print("from");
printer.print("to");
printer.print("distance");
printer.println();

for (String entry : distances) {
for (String col : entry.split(csvOptions.getFormat().getDelimiterString()))
printer.print(col);
printer.println();
}
}

logger.info("Done!");
return 0;
}

private Predicate<Geometry> getFilter(Path path){

if(path == null)
return community -> true;

List<Geometry> geometries = GeoFileReader.getAllFeatures(path.toString()).stream()
.map(feature -> (Geometry) feature.getDefaultGeometry())
.toList();

return community -> geometries.stream().anyMatch(geometry -> geometry.covers(community));
}
}

0 comments on commit c9bdfca

Please sign in to comment.