add fare zone based pt pricing

contribs/vsp/src/main/java/playground/vsp/pt/fare/DistanceBasedPtFareParams.java

@@ -22,6 +22,7 @@ public class DistanceBasedPtFareParams extends ReflectiveConfigGroup {
     public static final String LONG_DISTANCE_TRIP_THRESHOLD = "longDistanceTripThreshold";
     public static final String LONG_DISTANCE_TRIP_SLOPE = "longDistanceTripSlope";
     public static final String LONG_DISTANCE_TRIP_INTERCEPT = "longDistanceTripIntercept";
+    public static final String FARE_ZONE_SHP = "fareZoneShp";
 
     @PositiveOrZero
     private double minFare = 2.0;
@@ -35,6 +36,7 @@ public class DistanceBasedPtFareParams extends ReflectiveConfigGroup {
     private double longDistanceTripIntercept = 30.0;
     @PositiveOrZero
     private double longDistanceTripSlope = 0.00025;
+	private String fareZoneShp;
 
     public DistanceBasedPtFareParams() {
         super(SET_NAME);
@@ -52,6 +54,7 @@ public Map<String, String> getComments() {
         map.put(LONG_DISTANCE_TRIP_THRESHOLD, "Threshold of the long trips in meters. Below this value, " +
                 "the trips are considered as normal trips. Above this value, the trips are considered as " +
                 "inter-city trips");
+		map.put(FARE_ZONE_SHP, "Shp file with fare zone(s). This parameter is only used for PtFareCalculationModel 'fareZoneBased'.");
         return map;
     }
 
@@ -114,4 +117,14 @@ public double getLongDistanceTripThreshold() {
     public void setLongDistanceTripThreshold(double longDistanceTripThreshold) {
         this.longDistanceTripThreshold = longDistanceTripThreshold;
     }
+
+	@StringGetter(FARE_ZONE_SHP)
+	public String getFareZoneShp() {
+		return fareZoneShp;
+	}
+
+	@StringSetter(FARE_ZONE_SHP)
+	public void setFareZoneShp(String fareZoneShp) {
+		this.fareZoneShp = fareZoneShp;
+	}
 }

contribs/vsp/src/main/java/playground/vsp/pt/fare/FareZoneBasedPtFareHandler.java

@@ -0,0 +1,163 @@
+package playground.vsp.pt.fare;
+
+import com.google.inject.Inject;
+import org.apache.commons.math.stat.regression.SimpleRegression;
+import org.geotools.api.feature.simple.SimpleFeature;
+import org.locationtech.jts.geom.Geometry;
+import org.matsim.api.core.v01.Coord;
+import org.matsim.api.core.v01.Id;
+import org.matsim.api.core.v01.events.ActivityStartEvent;
+import org.matsim.api.core.v01.events.PersonMoneyEvent;
+import org.matsim.api.core.v01.events.handler.ActivityStartEventHandler;
+import org.matsim.api.core.v01.population.Person;
+import org.matsim.application.options.ShpOptions;
+import org.matsim.core.api.experimental.events.EventsManager;
+import org.matsim.core.router.StageActivityTypeIdentifier;
+import org.matsim.core.utils.geometry.CoordUtils;
+import org.matsim.core.utils.geometry.geotools.MGC;
+import org.matsim.pt.PtConstants;
+
+import java.util.AbstractMap;
+import java.util.HashMap;
+import java.util.List;
+import java.util.Map;
+
+public class FareZoneBasedPtFareHandler implements ActivityStartEventHandler {
+	@Inject
+	private EventsManager events;
+
+	public static final String FARE = "fare";
+	public static final String PT_FARE_ZONE_BASED = "fare zone based pt fare";
+	public static final String PT_GERMANWIDE_FARE_BASED = "german-wide fare based pt fare";
+
+	private final ShpOptions shp;
+
+	private final Map<Id<Person>, Coord> personDepartureCoordMap = new HashMap<>();
+	private final Map<Id<Person>, Coord> personArrivalCoordMap = new HashMap<>();
+
+	public FareZoneBasedPtFareHandler(DistanceBasedPtFareParams params) {
+		this.shp = new ShpOptions(params.getFareZoneShp(), null, null);
+	}
+
+	@Override
+	public void handleEvent(ActivityStartEvent event) {
+
+		if (event.getActType().equals(PtConstants.TRANSIT_ACTIVITY_TYPE)) {
+			personDepartureCoordMap.computeIfAbsent(event.getPersonId(), c -> event.getCoord()); // The departure place is fixed to the place of first pt interaction an agent has in the whole leg
+			personArrivalCoordMap.put(event.getPersonId(), event.getCoord()); // The arrival stop will keep updating until the agent start a real activity (i.e. finish the leg)
+		}
+
+		if (!StageActivityTypeIdentifier.isStageActivity(event.getActType())) {
+			Id<Person> personId = event.getPersonId();
+			if (personDepartureCoordMap.containsKey(personId)) {
+				double distance = CoordUtils.calcEuclideanDistance
+					(personDepartureCoordMap.get(personId), personArrivalCoordMap.get(personId));
+
+				SimpleFeature departureZone = determineFareZone(personDepartureCoordMap.get(personId), shp.readFeatures());
+				SimpleFeature arrivalZone = determineFareZone(personArrivalCoordMap.get(personId), shp.readFeatures());
+
+				Map.Entry<String, Double> fareEntry = computeFare(distance, departureZone, arrivalZone);
+				// charge fare to the person
+				events.processEvent(
+					new PersonMoneyEvent(event.getTime(), event.getPersonId(), -fareEntry.getValue(),
+						PtFareConfigGroup.PT_FARE, fareEntry.getKey(), event.getPersonId().toString()));
+
+				personDepartureCoordMap.remove(personId);
+				personArrivalCoordMap.remove(personId);
+			}
+		}
+	}
+
+	public static Map.Entry<String, Double> computeFare(double distance, SimpleFeature departureZone, SimpleFeature arrivalZone) {
+
+		if (departureZone != null && arrivalZone != null) {
+//			if both zones are not null -> departure and arrival point are inside of one of the tarifzonen
+			if (departureZone.getID().equals(arrivalZone.getID())) {
+				return new AbstractMap.SimpleEntry<>(PT_FARE_ZONE_BASED ,(double) departureZone.getAttribute(FARE));
+			}
+		}
+//		in every other case return german wide fare / Deutschlandtarif
+		return getGermanWideFare(distance);
+	}
+
+	private static Map.Entry<String, Double> getGermanWideFare(double distance) {
+
+		SimpleRegression regression = new SimpleRegression();
+
+//		in Deutschlandtarif, the linear function for the prices above 100km seem to have a different steepness
+//		hence the following difference in data points
+//		prices taken from https://deutschlandtarifverbund.de/wp-content/uploads/2024/07/20231201_TBDT_J_10_Preisliste_V07.pdf
+		if (distance / 1000 <= 100.) {
+			regression.addData(1, 1.70);
+			regression.addData(2,1.90);
+			regression.addData(3,2.00);
+			regression.addData(4,2.10);
+			regression.addData(5,2.20);
+			regression.addData(6,3.20);
+			regression.addData(7,3.70);
+			regression.addData(8,3.80);
+			regression.addData(9,3.90);
+			regression.addData(10,4.10);
+			regression.addData(11,5.00);
+			regression.addData(12,5.40);
+			regression.addData(13,5.60);
+			regression.addData(14,5.80);
+			regression.addData(15,5.90);
+			regression.addData(16,6.40);
+			regression.addData(17,6.50);
+			regression.addData(18,6.60);
+			regression.addData(19,6.70);
+			regression.addData(20,6.90);
+			regression.addData(30,9.90);
+			regression.addData(40,13.70);
+			regression.addData(50,16.30);
+			regression.addData(60,18.10);
+			regression.addData(70,20.10);
+			regression.addData(80,23.20);
+			regression.addData(90,26.20);
+			regression.addData(100,28.10);
+		} else {
+			regression.addData(100,28.10);
+			regression.addData(200,47.20);
+			regression.addData(300,59.70);
+			regression.addData(400,71.70);
+			regression.addData(500,83.00);
+			regression.addData(600,94.60);
+			regression.addData(700,106.30);
+			regression.addData(800,118.20);
+			regression.addData(900,130.10);
+			regression.addData(1000,141.00);
+			regression.addData(1100,148.60);
+			regression.addData(1200,158.10);
+			regression.addData(1300,169.20);
+			regression.addData(1400,179.80);
+			regression.addData(1500,190.10);
+			regression.addData(1600,201.50);
+			regression.addData(1700,212.80);
+			regression.addData(1800,223.30);
+			regression.addData(1900,233.90);
+			regression.addData(2000,244.00);
+		}
+		return new AbstractMap.SimpleEntry<>(PT_GERMANWIDE_FARE_BASED, regression.getSlope() * distance / 1000 + regression.getIntercept());
+	}
+
+	static SimpleFeature determineFareZone(Coord coord, List<SimpleFeature> features) {
+		SimpleFeature zone = null;
+
+		for (SimpleFeature ft : features) {
+			Geometry geom = (Geometry) ft.getDefaultGeometry();
+
+			if (MGC.coord2Point(coord).within(geom)) {
+				zone = ft;
+				break;
+			}
+		}
+		return zone;
+	}
+
+	@Override
+	public void reset(int iteration) {
+		personArrivalCoordMap.clear();
+		personDepartureCoordMap.clear();
+	}
+}

contribs/vsp/src/main/java/playground/vsp/pt/fare/PtFareConfigGroup.java

@@ -13,7 +13,7 @@ public class PtFareConfigGroup extends ReflectiveConfigGroup {
     public static final String APPLY_UPPER_BOUND = "applyUpperBound";
     public static final String UPPER_BOUND_FACTOR = "upperBoundFactor";
 
-    public enum PtFareCalculationModels {distanceBased} // More to come (e.g. zone based, hybrid...)
+    public enum PtFareCalculationModels {distanceBased, fareZoneBased} // More to come (e.g. zone based, hybrid...)
 
     private static final String PT_FARE_CALCULATION_CMT = "PT fare calculation scheme. Current implementation: distanceBased (more to come...)";
     public static final String UPPER_BOUND_FACTOR_CMT = "When upper bound is applied, upperBound  = upperBoundFactor * max Fare of the day. " +

contribs/vsp/src/main/java/playground/vsp/pt/fare/PtFareModule.java

@@ -11,12 +11,15 @@ public void install() {
         getConfig().scoring().getModes().get(TransportMode.pt).setDailyMonetaryConstant(0);
         getConfig().scoring().getModes().get(TransportMode.pt).setMarginalUtilityOfDistance(0);
         PtFareConfigGroup ptFareConfigGroup = ConfigUtils.addOrGetModule(this.getConfig(), PtFareConfigGroup.class);
+		DistanceBasedPtFareParams distanceBasedPtFareParams = ConfigUtils.addOrGetModule(this.getConfig(), DistanceBasedPtFareParams.class);
+
         if (ptFareConfigGroup.getPtFareCalculation() == PtFareConfigGroup.PtFareCalculationModels.distanceBased) {
-            DistanceBasedPtFareParams distanceBasedPtFareParams = ConfigUtils.addOrGetModule(this.getConfig(), DistanceBasedPtFareParams.class);
-            addEventHandlerBinding().toInstance(new DistanceBasedPtFareHandler(distanceBasedPtFareParams));
-        } else {
+			addEventHandlerBinding().toInstance(new DistanceBasedPtFareHandler(distanceBasedPtFareParams));
+		} else if (ptFareConfigGroup.getPtFareCalculation() == PtFareConfigGroup.PtFareCalculationModels.fareZoneBased) {
+			addEventHandlerBinding().toInstance(new FareZoneBasedPtFareHandler(distanceBasedPtFareParams));
+		} else {
             throw new RuntimeException("Please choose from the following fare Calculation method: [" +
-                    PtFareConfigGroup.PtFareCalculationModels.distanceBased + "]");
+                    PtFareConfigGroup.PtFareCalculationModels.distanceBased + ", " + PtFareConfigGroup.PtFareCalculationModels.fareZoneBased + "]");
         }
 
         if (ptFareConfigGroup.getApplyUpperBound()) {

contribs/vsp/src/test/java/playground/vsp/pt/fare/FareZoneBasedPtFareHandlerTest.java

@@ -0,0 +1,124 @@
+package playground.vsp.pt.fare;
+
+import org.junit.jupiter.api.Assertions;
+import org.junit.jupiter.api.Test;
+import org.junit.jupiter.api.extension.RegisterExtension;
+import org.matsim.analysis.personMoney.PersonMoneyEventsAnalysisModule;
+import org.matsim.api.core.v01.Coord;
+import org.matsim.api.core.v01.Id;
+import org.matsim.api.core.v01.TransportMode;
+import org.matsim.api.core.v01.population.*;
+import org.matsim.core.config.Config;
+import org.matsim.core.config.ConfigUtils;
+import org.matsim.core.config.groups.ScoringConfigGroup;
+import org.matsim.core.controler.AbstractModule;
+import org.matsim.core.controler.Controler;
+import org.matsim.core.controler.OutputDirectoryHierarchy;
+import org.matsim.core.scenario.MutableScenario;
+import org.matsim.core.scenario.ScenarioUtils;
+import org.matsim.core.utils.io.IOUtils;
+import org.matsim.examples.ExamplesUtils;
+import org.matsim.testcases.MatsimTestUtils;
+
+import java.io.BufferedReader;
+import java.io.IOException;
+import java.net.URL;
+import java.nio.file.Files;
+import java.nio.file.Path;
+import java.util.ArrayList;
+import java.util.List;
+
+import static org.matsim.application.ApplicationUtils.globFile;
+
+public class FareZoneBasedPtFareHandlerTest {
+
+	@RegisterExtension
+	private MatsimTestUtils utils = new MatsimTestUtils();
+
+	@Test
+	void testFareZoneBasedPtFareHandler() {
+
+		URL context = ExamplesUtils.getTestScenarioURL("kelheim");
+		Config config = ConfigUtils.loadConfig(IOUtils.extendUrl(context, "config.xml"));
+		config.controller().setOverwriteFileSetting(OutputDirectoryHierarchy.OverwriteFileSetting.deleteDirectoryIfExists);
+		config.controller().setOutputDirectory(utils.getOutputDirectory());
+		config.controller().setLastIteration(0);
+
+		PtFareConfigGroup ptFareConfigGroup = ConfigUtils.addOrGetModule(config, PtFareConfigGroup.class);
+		ptFareConfigGroup.setPtFareCalculationModel(PtFareConfigGroup.PtFareCalculationModels.fareZoneBased);
+
+		DistanceBasedPtFareParams fareParams = ConfigUtils.addOrGetModule(config, DistanceBasedPtFareParams.class);
+		fareParams.setFareZoneShp(IOUtils.extendUrl(context, "ptTestArea/pt-area.shp").toString());
+
+
+		ScoringConfigGroup scoring = ConfigUtils.addOrGetModule(config, ScoringConfigGroup.class);
+
+		ScoringConfigGroup.ActivityParams homeParams = new ScoringConfigGroup.ActivityParams("home");
+		ScoringConfigGroup.ActivityParams workParams = new ScoringConfigGroup.ActivityParams("work");
+		homeParams.setTypicalDuration(8 * 3600.);
+		workParams.setTypicalDuration(8 * 3600.);
+		scoring.addActivityParams(homeParams);
+		scoring.addActivityParams(workParams);
+
+
+		MutableScenario scenario = (MutableScenario) ScenarioUtils.loadScenario(config);
+
+		Population population = ScenarioUtils.createScenario(ConfigUtils.createConfig()).getPopulation();
+		PopulationFactory fac = population.getFactory();
+
+		Person person = fac.createPerson(Id.createPersonId("fareTestPerson"));
+		Plan plan = fac.createPlan();
+
+		Activity home = fac.createActivityFromCoord("home", new Coord(710300.624,5422165.737));
+//		bus to Saal (Donau) work location departs at 09:14
+		home.setEndTime(9 * 3600.);
+		Activity work = fac.createActivityFromCoord("work", new Coord(714940.65,5420707.78));
+//		rb17 to regensburg 2nd home location departs at 13:59
+		work.setEndTime(13 * 3600. + 45 * 60);
+		Activity home2 = fac.createActivityFromCoord("home", new Coord(726634.40,5433508.07));
+
+		Leg leg = fac.createLeg(TransportMode.pt);
+
+		plan.addActivity(home);
+		plan.addLeg(leg);
+		plan.addActivity(work);
+		plan.addLeg(leg);
+		plan.addActivity(home2);
+
+		person.addPlan(plan);
+		population.addPerson(person);
+		scenario.setPopulation(population);
+
+		Controler controler = new Controler(scenario);
+		controler.addOverridingModule(new AbstractModule() {
+			@Override
+			public void install() {
+				install(new PtFareModule());
+				install(new PersonMoneyEventsAnalysisModule());
+			}
+		});
+		controler.run();
+
+		Assertions.assertTrue(Files.exists(Path.of(utils.getOutputDirectory())));
+
+//		read personMoneyEvents.tsv and check if both fare entries do have the correct fare type
+		String filePath = globFile(Path.of(utils.getOutputDirectory()), "*output_personMoneyEvents.tsv*").toString();
+		String line;
+		List<String[]> events = new ArrayList<>();
+
+		try (BufferedReader br = IOUtils.getBufferedReader(filePath)) {
+//			skip header
+			br.readLine();
+
+			while ((line = br.readLine()) != null) {
+				events.add(line.split(";"));
+			}
+		} catch (IOException e) {
+			e.printStackTrace();
+		}
+
+		Assertions.assertEquals(2, events.size());
+		Assertions.assertEquals(FareZoneBasedPtFareHandler.PT_FARE_ZONE_BASED, events.get(0)[4]);
+		Assertions.assertEquals(FareZoneBasedPtFareHandler.PT_GERMANWIDE_FARE_BASED, events.get(1)[4]);
+	}
+}

examples/scenarios/kelheim/ptTestArea/pt-area.cpg

@@ -0,0 +1 @@
+UTF-8

examples/scenarios/kelheim/ptTestArea/pt-area.prj

@@ -0,0 +1 @@
+PROJCS["ETRS_1989_UTM_Zone_32N",GEOGCS["GCS_ETRS_1989",DATUM["D_ETRS_1989",SPHEROID["GRS_1980",6378137.0,298.257222101]],PRIMEM["Greenwich",0.0],UNIT["Degree",0.0174532925199433]],PROJECTION["Transverse_Mercator"],PARAMETER["False_Easting",500000.0],PARAMETER["False_Northing",0.0],PARAMETER["Central_Meridian",9.0],PARAMETER["Scale_Factor",0.9996],PARAMETER["Latitude_Of_Origin",0.0],UNIT["Meter",1.0]]