add code for gartenfeld scenario creation

src/main/java/org/matsim/policies/gartenfeld/CreateGartenfeldPopulation.java

@@ -0,0 +1,71 @@
+package org.matsim.policies.gartenfeld;
+
+import org.matsim.application.MATSimAppCommand;
+import org.matsim.application.prepare.population.MergePopulations;
+import org.matsim.prepare.population.CreateFixedPopulation;
+import org.matsim.prepare.population.InitLocationChoice;
+import org.matsim.prepare.population.RemoveUnavailableRoutes;
+import org.matsim.prepare.population.RunActivitySampling;
+import picocli.CommandLine;
+
+@CommandLine.Command(name = "create-gartenfeld-population", description = "Create the population for the Gartenfeld scenario.")
+public class CreateGartenfeldPopulation implements MATSimAppCommand {
+
+	private final static String SVN = "../shared-svn/projects/matsim-germany";
+
+	@CommandLine.Option(names = "--output", description = "Path to output population", defaultValue = "input/gartenfeld/gartenfeld-population-10pct.xml.gz")
+	private String output;
+
+	public static void main(String[] args) {
+		new CreateGartenfeldPopulation().execute(args);
+	}
+
+	@Override
+	public Integer call() throws Exception {
+
+		new CreateFixedPopulation().execute(
+			"--n", "7400",
+			"--sample", "0.1",
+			"--unemployed", "0.013",
+			"--age-dist", "0.149", "0.203",
+			"--facilities", "input/gartenfeld/DNG_residential.gpkg",
+			"--prefix", "dng",
+			"--output", output
+		);
+
+		new RunActivitySampling().execute(
+			"--seed", "1",
+			"--persons", "src/main/python/table-persons.csv",
+			"--activities", "src/main/python/table-activities.csv",
+			"--input", output,
+			"--output", output
+		);
+
+
+		new InitLocationChoice().execute(
+			"--input", output,
+			"--output", output,
+			"--k", "1",
+			"--facilities", "input/v6.3/berlin-v6.3-facilities.xml.gz",
+			"--network", "input/v6.3/berlin-v6.3-network.xml.gz",
+			"--shp", SVN + "/vg5000/vg5000_ebenen_0101/VG5000_GEM.shp",
+			"--commuter", SVN + "/regionalstatistik/commuter.csv",
+			"--commute-prob", "0.1",
+			"--sample", "0.1"
+		);
+
+		// Merge with calibrated plans into one
+		new MergePopulations().execute(
+			output, "input/v6.3/berlin-v6.3-10pct.plans.xml.gz",
+			"--output", output
+		);
+
+		new RemoveUnavailableRoutes().execute(
+			"--input", output,
+			"--network", "input/gartenfeld/gartenfeld-network.xml.gz",
+			"--output", output
+		);
+
+		return 0;
+	}
+}

src/main/java/org/matsim/policies/gartenfeld/RunGartenfeldScenario.java

@@ -0,0 +1,21 @@
+package org.matsim.policies.gartenfeld;
+
+import org.matsim.application.MATSimApplication;
+import org.matsim.run.OpenBerlinScenario;
+
+/**
+ * Run class for the Gartenfeld scenario.
+ */
+public final class RunGartenfeldScenario {
+
+	private RunGartenfeldScenario() {
+	}
+
+	public static void main(String[] args) {
+		MATSimApplication.runWithDefaults(OpenBerlinScenario.class, args,
+			"--config:plans.inputPlansFile", "../../input/gartenfeld/gartenfeld-population-10pct.xml.gz",
+			"--config:network.inputNetworkFile", "../../input/gartenfeld/gartenfeld-network.xml.gz"
+			);
+	}
+
+}

src/main/java/org/matsim/prepare/population/CreateBerlinPopulation.java

@@ -8,6 +8,7 @@
 import org.apache.logging.log4j.Logger;
 import org.geotools.api.feature.simple.SimpleFeature;
 import org.locationtech.jts.geom.Envelope;
+import org.locationtech.jts.geom.Geometry;
 import org.locationtech.jts.geom.MultiPolygon;
 import org.matsim.api.core.v01.Coord;
 import org.matsim.api.core.v01.Id;
@@ -85,7 +86,7 @@ public static Id<Person> generateId(Population population, String prefix, Splitt
 	/**
 	 * Samples a home coordinates from geometry and landuse.
 	 */
-	public static Coord sampleHomeCoordinate(MultiPolygon geometry, String crs, FacilityOptions facilities, SplittableRandom rnd) {
+	public static Coord sampleHomeCoordinate(Geometry geometry, String crs, FacilityOptions facilities, SplittableRandom rnd) {
 
 		Envelope bbox = geometry.getEnvelopeInternal();
 

src/main/java/org/matsim/prepare/population/CreateFixedPopulation.java

@@ -0,0 +1,151 @@
+package org.matsim.prepare.population;
+
+import org.apache.logging.log4j.LogManager;
+import org.apache.logging.log4j.Logger;
+import org.locationtech.jts.geom.Geometry;
+import org.matsim.api.core.v01.Coord;
+import org.matsim.api.core.v01.population.Person;
+import org.matsim.api.core.v01.population.Plan;
+import org.matsim.api.core.v01.population.Population;
+import org.matsim.api.core.v01.population.PopulationFactory;
+import org.matsim.application.MATSimAppCommand;
+import org.matsim.core.config.ConfigUtils;
+import org.matsim.core.population.PersonUtils;
+import org.matsim.core.population.PopulationUtils;
+import org.matsim.core.scenario.ProjectionUtils;
+import org.matsim.run.OpenBerlinScenario;
+import picocli.CommandLine;
+
+import java.nio.file.Path;
+import java.util.List;
+import java.util.Map;
+import java.util.SplittableRandom;
+import java.util.stream.IntStream;
+
+@CommandLine.Command(
+	name = "fixed-population",
+	description = "Create synthetic population for Gartenfeld."
+)
+public class CreateFixedPopulation implements MATSimAppCommand {
+
+	private static final Logger log = LogManager.getLogger(CreateFixedPopulation.class);
+
+	@CommandLine.Mixin
+	private FacilityOptions facilities = new FacilityOptions();
+
+	@CommandLine.Option(names = "--n", description = "Number of persons to generate", required = true)
+	private int n;
+
+	@CommandLine.Option(names = "--prefix", description = "Prefix for person ids", required = true)
+	private String prefix;
+
+	@CommandLine.Option(names = "--age-dist", description = "Age distribution for < 18 and 65+", arity = "2", required = true)
+	private List<Double> ageDist;
+
+	@CommandLine.Option(names = "--unemployed", description = "Unemployment rate", required = true)
+	private double unemployed;
+
+	@CommandLine.Option(names = "--gender-dist", description = "Proportion of women", defaultValue = "0.5")
+	private double genderDist;
+
+	@CommandLine.Option(names = "--output", description = "Path to output population", required = true)
+	private Path output;
+
+	@CommandLine.Option(names = "--sample", description = "Sample size to generate", defaultValue = "0.25")
+	private double sample;
+
+	private SplittableRandom rnd;
+	private Population population;
+
+	public static void main(String[] args) {
+		new CreateFixedPopulation().execute(args);
+	}
+
+	@Override
+	@SuppressWarnings("IllegalCatch")
+	public Integer call() throws Exception {
+
+		rnd = new SplittableRandom(0);
+		population = PopulationUtils.createPopulation(ConfigUtils.createConfig());
+
+		generatePersons();
+
+		log.info("Generated {} persons", population.getPersons().size());
+
+		PopulationUtils.sortPersons(population);
+
+		ProjectionUtils.putCRS(population, OpenBerlinScenario.CRS);
+		PopulationUtils.writePopulation(population, output.toString());
+
+		return 0;
+	}
+
+	private void generatePersons() {
+
+		double young = ageDist.get(0);
+		double old = ageDist.get(1);
+
+		// x women for 100 men
+		double quota = genderDist;
+
+		var sex = new EnumeratedAttributeDistribution<>(Map.of("f", quota, "m", 1 - quota));
+		var employment = new EnumeratedAttributeDistribution<>(Map.of(true, 1 - unemployed, false, unemployed));
+		var ageGroup = new EnumeratedAttributeDistribution<>(Map.of(
+			AgeGroup.YOUNG, young,
+			AgeGroup.MIDDLE, 1.0 - young - old,
+			AgeGroup.OLD, old
+		));
+
+		PopulationFactory f = population.getFactory();
+		Geometry geom = facilities.getGeometry().convexHull();
+
+		var youngDist = new UniformAttributeDistribution<>(IntStream.range(1, 18).boxed().toList());
+		var middleDist = new UniformAttributeDistribution<>(IntStream.range(18, 65).boxed().toList());
+		var oldDist = new UniformAttributeDistribution<>(IntStream.range(65, 100).boxed().toList());
+
+		for (int i = 0; i < n * sample; i++) {
+
+			Person person = f.createPerson(CreateBerlinPopulation.generateId(population, prefix, rnd));
+			PersonUtils.setSex(person, sex.sample());
+			PopulationUtils.putSubpopulation(person, "person");
+
+			AgeGroup group = ageGroup.sample();
+
+			if (group == AgeGroup.MIDDLE) {
+				PersonUtils.setAge(person, middleDist.sample());
+				PersonUtils.setEmployed(person, employment.sample());
+			} else if (group == AgeGroup.YOUNG) {
+				PersonUtils.setAge(person, youngDist.sample());
+				PersonUtils.setEmployed(person, false);
+			} else if (group == AgeGroup.OLD) {
+				PersonUtils.setAge(person, oldDist.sample());
+				PersonUtils.setEmployed(person, false);
+			}
+
+			Coord coord = CreateBerlinPopulation.sampleHomeCoordinate(geom, OpenBerlinScenario.CRS, facilities, rnd);
+
+			person.getAttributes().putAttribute(Attributes.HOME_X, coord.getX());
+			person.getAttributes().putAttribute(Attributes.HOME_Y, coord.getY());
+
+			// Currently hard-coded as berlin inhabitants
+			person.getAttributes().putAttribute(Attributes.GEM, 11000000);
+			person.getAttributes().putAttribute(Attributes.ARS, 110000000000L);
+			person.getAttributes().putAttribute(Attributes.RegioStaR7, 1);
+
+			Plan plan = f.createPlan();
+			plan.addActivity(f.createActivityFromCoord("home", coord));
+
+			person.addPlan(plan);
+			person.setSelectedPlan(plan);
+
+			population.addPerson(person);
+		}
+	}
+
+	private enum AgeGroup {
+		YOUNG,
+		MIDDLE,
+		OLD
+	}
+
+}

src/main/java/org/matsim/prepare/population/FacilityOptions.java

@@ -2,6 +2,7 @@
 
 import org.apache.logging.log4j.LogManager;
 import org.apache.logging.log4j.Logger;
+import org.locationtech.jts.geom.Geometry;
 import org.matsim.api.core.v01.Coord;
 import org.matsim.application.options.ShpOptions;
 import picocli.CommandLine;
@@ -45,13 +46,21 @@ public synchronized ShpOptions.Index getIndex(String queryCRS) {
 		ShpOptions shp = ShpOptions.ofLayer(facilityPath.toString(), null);
 
 		index = shp.createIndex(queryCRS, attr, ft -> Boolean.TRUE.equals(ft.getAttribute(attr))
-			|| Objects.equals(ft.getAttribute(attr), 1));
+			|| Objects.equals(ft.getAttribute(attr), 1) || Objects.equals(ft.getAttribute(attr), "1"));
 
 		log.info("Read {} features for {} facilities", index.size(), attr);
 
 		return index;
 	}
 
+	/**
+	 * Get the union geometry of the shape file.
+	 */
+	public Geometry getGeometry() {
+		ShpOptions shp = ShpOptions.ofLayer(facilityPath.toString(), null);
+		return shp.getGeometry();
+	}
+
 	/**
 	 * Tries to select a point that lies within one of the geometries of the index.
 	 * Will try at least {@link #iters} times, after which the last point is returned even if not within a geometry.

src/main/java/org/matsim/prepare/population/InitLocationChoice.java

@@ -70,6 +70,9 @@ public class InitLocationChoice implements MATSimAppCommand, PersonAlgorithm {
 	@CommandLine.Option(names = "--commuter", description = "Path to commuter.csv", required = true)
 	private Path commuterPath;
 
+	@CommandLine.Option(names = "--commute-prob", description = "Probability for commuting to a different zone", defaultValue = "1")
+	private double commuteProb;
+
 	@CommandLine.Option(names = "--facilities", description = "Path to facilities file", required = true)
 	private Path facilityPath;
 
@@ -312,8 +315,10 @@ private ActivityFacility sampleCommute(SplittableRandom rnd, double dist, Coord
 
 		ActivityFacility workPlace = null;
 
+		boolean commute = commuteProb >= 1 || rnd.nextDouble() < commuteProb;
+
 		// Only larger distances can be commuters to other zones
-		if (dist > 3000) {
+		if (dist > 3000 && commute) {
 			workPlace = commuter.selectTarget(rnd, ars, dist, MGC.coord2Point(refCoord), zone -> sampleZone(index, dist, refCoord, zone, rnd));
 		}
 

src/main/java/org/matsim/prepare/population/RemoveUnavailableRoutes.java

@@ -0,0 +1,74 @@
+package org.matsim.prepare.population;
+
+import org.matsim.api.core.v01.Id;
+import org.matsim.api.core.v01.network.Link;
+import org.matsim.api.core.v01.network.Network;
+import org.matsim.api.core.v01.population.*;
+import org.matsim.application.MATSimAppCommand;
+import org.matsim.core.network.NetworkUtils;
+import org.matsim.core.population.PopulationUtils;
+import org.matsim.core.population.algorithms.ParallelPersonAlgorithmUtils;
+import org.matsim.core.population.algorithms.PersonAlgorithm;
+import org.matsim.core.population.routes.NetworkRoute;
+import org.matsim.core.router.TripStructureUtils;
+import picocli.CommandLine;
+
+import java.util.Set;
+
+@CommandLine.Command(name = "remove-unavailable-routes", description = "Remove routes from the population that are not available in the network.")
+public class RemoveUnavailableRoutes implements MATSimAppCommand, PersonAlgorithm {
+
+	@CommandLine.Option(names = "--input", description = "Path to input population", required = true)
+	private String input;
+
+	@CommandLine.Option(names = "--network", description = "Path to network", required = true)
+	private String networkPath;
+
+	@CommandLine.Option(names = "--output", description = "Path to output population", required = true)
+	private String output;
+
+	private Network network;
+
+	@Override
+	public Integer call() throws Exception {
+
+		network = NetworkUtils.readNetwork(networkPath);
+
+		Population population = PopulationUtils.readPopulation(input);
+
+		ParallelPersonAlgorithmUtils.run(population, Runtime.getRuntime().availableProcessors(), this);
+
+		PopulationUtils.writePopulation(population, output);
+
+		return 0;
+	}
+
+	@Override
+	public void run(Person person) {
+
+		Set<Id<Link>> allLinks = network.getLinks().keySet();
+
+		for (Plan plan : person.getPlans()) {
+			for (Leg leg : TripStructureUtils.getLegs(plan.getPlanElements())) {
+				Route route = leg.getRoute();
+
+				if (route == null)
+					continue;
+
+				if (route instanceof NetworkRoute nr) {
+					if (!allLinks.containsAll(nr.getLinkIds()))
+						leg.setRoute(null);
+				}
+
+				if (!allLinks.contains(route.getStartLinkId()) || !allLinks.contains(route.getEndLinkId()))
+					leg.setRoute(null);
+			}
+
+			for (Activity act : TripStructureUtils.getActivities(plan.getPlanElements(), TripStructureUtils.StageActivityHandling.StagesAsNormalActivities)) {
+				if (act.getLinkId() != null && !allLinks.contains(act.getLinkId()))
+					act.setLinkId(null);
+
+			}
+		}
+	}
+}