diff --git a/.gitignore b/.gitignore
index 0a7ae76b5..04dd911da 100644
--- a/.gitignore
+++ b/.gitignore
@@ -26,3 +26,4 @@ tripAttraction.csv
 tripProduction.csv
 
 src/main/java/META-INF/MANIFEST.MF
+.Rproj.user
diff --git a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimScenarioAssembler.java b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimScenarioAssembler.java
index 00825e99d..9bf3067c5 100644
--- a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimScenarioAssembler.java
+++ b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimScenarioAssembler.java
@@ -2,15 +2,28 @@
 
 import de.tum.bgu.msm.data.Day;
 import de.tum.bgu.msm.data.travelTimes.TravelTimes;
+import org.matsim.api.core.v01.Id;
 import org.matsim.api.core.v01.Scenario;
 import org.matsim.core.config.Config;
 
+import java.util.HashMap;
 import java.util.List;
 import java.util.Map;
 
 public interface MatsimScenarioAssembler {
 
     Scenario assembleScenario(Config initialMatsimConfig, int year, TravelTimes travelTimes);
+
+    // used default identifier to allow for using following method in only one class that implements MatsimScenarioAssembler.
+    default Scenario assembleScenarioWithPreviousPlans(
+            Scenario matsimScenario,
+            int year,
+            TravelTimes travelTimes,
+            HashMap<Id<org.matsim.api.core.v01.population.Person>, org.matsim.api.core.v01.population.Person> previousYearPopulation) {
+
+        throw new UnsupportedOperationException(
+                "This MatsimScenarioAssembler does not support previous-year plans.");
+    }
     Map<Day, Scenario> assembleMultiScenarios(Config initialMatsimConfig, int year, TravelTimes travelTimes);
 
 }
diff --git a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimSkimCreator.java b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimSkimCreator.java
index 0e12a397d..44b65d071 100644
--- a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimSkimCreator.java
+++ b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimSkimCreator.java
@@ -13,6 +13,7 @@
 import org.matsim.api.core.v01.Id;
 import org.matsim.api.core.v01.network.Network;
 import org.matsim.api.core.v01.network.Node;
+import org.matsim.api.core.v01.network.Link;
 import org.matsim.api.core.v01.population.Leg;
 import org.matsim.api.core.v01.population.PlanElement;
 import org.matsim.core.network.NetworkUtils;
@@ -24,6 +25,10 @@
 import org.matsim.facilities.Facility;
 import org.matsim.pt.transitSchedule.api.TransitStopFacility;
 
+import java.io.BufferedWriter;
+import java.io.IOException;
+import java.nio.file.Files;
+import java.nio.file.Path;
 import java.util.*;
 
 public class MatsimSkimCreator {
@@ -31,9 +36,25 @@ public class MatsimSkimCreator {
     private final static Logger logger = LogManager.getLogger(MatsimSkimCreator.class);
 
     private MatsimData matsimData;
+    private final Path csvOutputDirectory;
+    private final Integer csvOutputYear;
 
     public MatsimSkimCreator(MatsimData provider) {
+        this(provider, null, null);
+    }
+
+    /**
+     * Creates a skim creator that also writes every generated skim to CSV.
+     * The CSV files use the project standard FROM,TO,VALUE format.
+     */
+    public MatsimSkimCreator(MatsimData provider, Path csvOutputDirectory) {
+        this(provider, csvOutputDirectory, null);
+    }
+
+    public MatsimSkimCreator(MatsimData provider, Path csvOutputDirectory, Integer csvOutputYear) {
         this.matsimData = provider;
+        this.csvOutputDirectory = csvOutputDirectory;
+        this.csvOutputYear = csvOutputYear;
     }
 
     public IndexedDoubleMatrix2D createCarSkim(Collection<? extends de.tum.bgu.msm.data.Id> zones, int numberOfThreads, double peakHour_s) {
@@ -77,10 +98,29 @@ public IndexedDoubleMatrix2D createCarSkim(Collection<? extends de.tum.bgu.msm.d
             });
         }
         executor.execute();
+        logZeroCarSkimCells(skim);
         assignIntrazonals(5, Float.MAX_VALUE, 0.66f, skim);
+        writeSkimToCsv(skim, "car");
         return skim;
     }
 
+    private void logZeroCarSkimCells(IndexedDoubleMatrix2D skim) {
+        long zeroCells = 0;
+        long zeroOffDiagonalCells = 0;
+        for (int row = 0; row < skim.rows(); row++) {
+            int origin = skim.getIdForInternalRowIndex(row);
+            for (int column = 0; column < skim.columns(); column++) {
+                int destination = skim.getIdForInternalColumnIndex(column);
+                if (skim.getIndexed(origin, destination) == 0.) {
+                    zeroCells++;
+                    if (origin != destination) zeroOffDiagonalCells++;
+                }
+            }
+        }
+        logger.warn("Car skim contains " + zeroCells + " zero cells, including "
+                + zeroOffDiagonalCells + " zero off-diagonal cells, before intrazonal post-processing.");
+    }
+
     public IndexedDoubleMatrix2D createPtSkim(Collection<? extends de.tum.bgu.msm.data.Id> zones, int numberOfThreads,
                                               double peakHour_s) {
         final int partitionSize = (int) ((double) zones.size() / numberOfThreads) + 1;
@@ -160,6 +200,7 @@ public IndexedDoubleMatrix2D createPtSkim(Collection<? extends de.tum.bgu.msm.da
         }
         executor.execute();
         assignIntrazonals(5, Float.MAX_VALUE, 0.66f, skim);
+        writeSkimToCsv(skim, "pt");
         return skim;
     }
 
@@ -210,6 +251,7 @@ public IndexedDoubleMatrix2D createTeleportedSkim(Collection<? extends de.tum.bg
         }
         executor.execute();
         assignIntrazonals(5, Float.MAX_VALUE, 0.66f, skim);
+        writeSkimToCsv(skim, mode);
         return skim;
     }
 
@@ -259,21 +301,52 @@ public IndexedDoubleMatrix2D createFreeSpeedFactorSkim(Collection<? extends de.t
         }
         executor.execute();
         assignIntrazonals(5, Float.MAX_VALUE, 0.66f, skim);
+        writeSkimToCsv(skim, "freeSpeedFactor");
         return skim;
     }
 
+    /**
+     * Writes a skim in long CSV format with one row per origin-destination pair.
+     */
+    public void writeSkimToCsv(IndexedDoubleMatrix2D skim, String name) {
+        if (csvOutputDirectory == null) {
+            return;
+        }
+
+        String yearSuffix = csvOutputYear == null ? "" : "_" + csvOutputYear;
+        Path outputFile = csvOutputDirectory.resolve(name + "Skim" + yearSuffix + ".csv");
+        try {
+            Files.createDirectories(csvOutputDirectory);
+            try (BufferedWriter writer = Files.newBufferedWriter(outputFile)) {
+                writer.write("FROM,TO,VALUE");
+                writer.newLine();
+                for (int row = 0; row < skim.rows(); row++) {
+                    int origin = skim.getIdForInternalRowIndex(row);
+                    for (int column = 0; column < skim.columns(); column++) {
+                        int destination = skim.getIdForInternalColumnIndex(column);
+                        writer.write(origin + "," + destination + "," + skim.getIndexed(origin, destination));
+                        writer.newLine();
+                    }
+                }
+            }
+            logger.info("Wrote " + name + " skim to " + outputFile);
+        } catch (IOException e) {
+            throw new RuntimeException("Could not write skim CSV to " + outputFile, e);
+        }
+    }
+
     //TODO: copied from MITO car skim updater...maybe provide a utility function there
     private void assignIntrazonals(int numberOfNeighbours, float maximumMinutes, float proportionOfTime, IndexedDoubleMatrix2D skim) {
-        int nonIntrazonalCounter = 0;
-        for (int i = 1; i < skim.columns(); i++) {
+        for (int i = 0; i < skim.columns(); i++) {
             int i_id = skim.getIdForInternalColumnIndex(i);
             double[] minTimeValues = new double[numberOfNeighbours];
             for (int k = 0; k < numberOfNeighbours; k++) {
                 minTimeValues[k] = maximumMinutes;
             }
             //find the  n closest neighbors - the lower travel time values in the matrix column
-            for (int j = 1; j < skim.rows(); j++) {
+            for (int j = 0; j < skim.rows(); j++) {
                 int j_id = skim.getIdForInternalRowIndex(j);
+                if (i_id == j_id) continue;
                 int minimumPosition = 0;
                 while (minimumPosition < numberOfNeighbours) {
                     if (minTimeValues[minimumPosition] > skim.getIndexed(i_id, j_id) && skim.getIndexed(i_id, j_id) != 0) {
@@ -292,19 +365,12 @@ private void assignIntrazonals(int numberOfNeighbours, float maximumMinutes, flo
             }
             globalMinTime = globalMinTime / numberOfNeighbours * proportionOfTime;
 
-            //fill with the calculated value the cells with zero
-            for (int j = 1; j < skim.rows(); j++) {
-                int j_id = skim.getIdForInternalColumnIndex(j);
-                if (skim.getIndexed(i_id, j_id) == 0) {
-                    skim.setIndexed(i_id, j_id, globalMinTime);
-                    if (i != j) {
-                        nonIntrazonalCounter++;
-                    }
-                }
+            // Only the diagonal is intrazonal. Keep off-diagonal zero values
+            // visible because they indicate missing or unreachable routes.
+            if (skim.getIndexed(i_id, i_id) == 0) {
+                skim.setIndexed(i_id, i_id, globalMinTime);
             }
         }
         logger.info("Calculated intrazonal times and distances using the " + numberOfNeighbours + " nearest neighbours.");
-        logger.info("The calculation of intrazonals has also assigned values for cells with travel time equal to 0, that are not intrazonal: (" +
-                nonIntrazonalCounter + " cases).");
     }
-}
\ No newline at end of file
+}
diff --git a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimTransportModel.java b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimTransportModel.java
index 6a90af4fe..d3305275f 100644
--- a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimTransportModel.java
+++ b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimTransportModel.java
@@ -36,7 +36,6 @@
 import org.matsim.core.config.ConfigUtils;
 import org.matsim.core.controler.Controler;
 import org.matsim.core.controler.ControlerDefaults;
-import org.matsim.core.controler.OutputDirectoryHierarchy;
 import org.matsim.core.router.util.TravelDisutility;
 import org.matsim.core.router.util.TravelTime;
 import org.matsim.core.scenario.ScenarioUtils;
@@ -47,10 +46,7 @@
 
 import java.io.File;
 
-import java.util.Collections;
-import java.util.HashSet;
-import java.util.Objects;
-import java.util.Set;
+import java.util.*;
 
 import static de.tum.bgu.msm.properties.modules.TransportModelPropertiesModule.*;
 
@@ -71,6 +67,8 @@ public final class MatsimTransportModel implements TransportModel {
 
     private MatsimScenarioAssembler scenarioAssembler;
 
+    protected HashMap<Id<org.matsim.api.core.v01.population.Person>, org.matsim.api.core.v01.population.Person> previousYearPopulation = new HashMap<>();
+
     public MatsimTransportModel(DataContainer dataContainer, Config matsimConfig,
                                 Properties properties, MatsimScenarioAssembler scenarioAssembler,
                                 MatsimData matsimData) {
@@ -102,7 +100,7 @@ public void setup() {
             runTransportModel(properties.main.startYear);
         } else {
             String eventsFile = properties.main.baseDirectory + properties.transportModel.matsimInitialEventsFile;
-            replayFromEvents(eventsFile);
+            replayFromEvents(eventsFile, properties.main.startYear);
         }
         logger.warn("### About to leave MatsimTransportModel#setup()." );
     }
@@ -134,7 +132,7 @@ private void runTransportModel(int year) {
             //if using the SimpleCommuteModeChoiceScenarioAssembler, we need some initial travel times (this will use an unlodaded network)
             TravelTime myTravelTime = SiloMatsimUtils.getAnEmptyNetworkTravelTime();
             TravelDisutility myTravelDisutility = SiloMatsimUtils.getAnEmptyNetworkTravelDisutility();
-            updateTravelTimes(myTravelTime, myTravelDisutility);
+            updateTravelTimes(myTravelTime, myTravelDisutility, year);
         }
 		{
 			Set<Id<org.matsim.api.core.v01.population.Person>> toRemove = new HashSet<>( assembledScenario.getPopulation().getPersons().keySet() );
@@ -149,6 +147,7 @@ private void runTransportModel(int year) {
 			}
 		}
 
+        // no input plans file - we first add silo agents to matsim
         for (Household household: dataContainer.getHouseholdDataManager().getHouseholds()) {
             for (Person pp : household.getPersons().values()) {
                 final int noHHAUtos = (int) household.getVehicles().stream().filter( vv -> vv.getType().equals( de.tum.bgu.msm.data.vehicle.VehicleType.CAR ) ).count();
@@ -179,7 +178,28 @@ private void runTransportModel(int year) {
 
         // TODO remove config argument as it is duplicate (cf. above)
 		// yyyyyy the following replaces "assembledScenario"!!!
-        assembledScenario = scenarioAssembler.assembleScenario(initialMatsimConfig, year, travelTimes);
+//        assembledScenario = scenarioAssembler.assembleScenario(
+//                assembledScenario.getConfig(),
+//                year,
+//                travelTimes);
+
+//
+
+        if (year == properties.main.baseYear || properties.transportModel.matsimStartOption != MatsimStartOptionIdentifier.HOT) {
+            // this is where we generate plans for silo agents
+            assembledScenario = scenarioAssembler.assembleScenario(
+                    assembledScenario.getConfig(),
+                    year,
+                    travelTimes);
+        } else {
+            // this is what we call in subsequent years - tries to identify silo agents that already have previous plan
+            assembledScenario = scenarioAssembler.assembleScenarioWithPreviousPlans(
+                    assembledScenario,
+                    year,
+                    travelTimes,
+                    previousYearPopulation);
+        }
+
 
         finalizeConfig(assembledScenario.getConfig(), year, Integer.toString(year));
 
@@ -191,9 +211,21 @@ private void runTransportModel(int year) {
         // Get travel Times from MATSim
         logger.warn("###################################################");
         logger.warn("Using MATSim to compute travel times from zone to zone.");
+
         TravelTime travelTime = controler.getLinkTravelTimes();
         TravelDisutility travelDisutility = controler.getTravelDisutilityFactory().createTravelDisutility(travelTime);
-        updateTravelTimes(travelTime, travelDisutility);
+        updateTravelTimes(travelTime, travelDisutility, year);
+
+        //don't actually need to check if not base year, since we save this after the simulation of first year,
+        // which is correct to then save "previous" year's plans
+//        if (year != properties.main.baseYear && properties.transportModel.transportModelIdentifier == TransportModelIdentifier.MATSIM) {
+
+        // here is where we add our current year silo agents to a map (may run into memory issues, so design may change).
+        previousYearPopulation.clear();
+            for (org.matsim.api.core.v01.population.Person person : assembledScenario.getPopulation().getPersons().values()) {
+                previousYearPopulation.put(person.getId(), person);
+            }
+
     }
 
     private void finalizeConfig(Config config, int runId, String dir) {
@@ -213,16 +245,16 @@ private void finalizeConfig(Config config, int runId, String dir) {
     /**
      * @param eventsFile
      */
-    private void replayFromEvents(String eventsFile) {
+    private void replayFromEvents(String eventsFile, int year) {
         Scenario scenario = ScenarioUtils.loadScenario(initialMatsimConfig);
         TravelTime travelTime = TravelTimeUtils.createTravelTimesFromEvents(scenario.getNetwork(), scenario.getConfig(), eventsFile);
         TravelDisutility travelDisutility = ControlerDefaults.createDefaultTravelDisutilityFactory(scenario).createTravelDisutility(travelTime);
-        updateTravelTimes(travelTime, travelDisutility);
+        updateTravelTimes(travelTime, travelDisutility, year);
     }
 
-    private void updateTravelTimes(TravelTime travelTime, TravelDisutility disutility) {
+    private void updateTravelTimes(TravelTime travelTime, TravelDisutility disutility, int year) {
         matsimData.update(disutility, travelTime);
-        internalTravelTimes.update(matsimData);
+        internalTravelTimes.update(matsimData, year);
         final TravelTimes mainTravelTimes = dataContainer.getTravelTimes();
 
         if (mainTravelTimes != this.internalTravelTimes && mainTravelTimes instanceof SkimTravelTimes) {
@@ -234,4 +266,5 @@ private void updateTravelTimes(TravelTime travelTime, TravelDisutility disutilit
                     dataContainer.getGeoData().getZones().values());
     }
     }
+
 }
diff --git a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimTravelTimesAndCosts.java b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimTravelTimesAndCosts.java
index c2d00fccb..08215c709 100644
--- a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimTravelTimesAndCosts.java
+++ b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/MatsimTravelTimesAndCosts.java
@@ -38,6 +38,7 @@
 import org.matsim.vehicles.VehicleUtils;
 import org.matsim.vehicles.VehiclesFactory;
 
+import java.nio.file.Paths;
 import java.util.*;
 
 /**
@@ -59,6 +60,7 @@ public final class MatsimTravelTimesAndCosts implements TravelTimes {
     private Collection<Region> regions;
 
     private final Config config;
+    private Integer skimOutputYear;
 
     public MatsimTravelTimesAndCosts(Config config) {
         this.config = config;
@@ -71,7 +73,12 @@ public void initialize(GeoData geoData, MatsimData matsimData) {
     }
 
     public void update(MatsimData matsimData) {
+        update(matsimData, null);
+    }
+
+    public void update(MatsimData matsimData, Integer year) {
         this.matsimData = matsimData;
+        this.skimOutputYear = year;
         this.tripRouter = matsimData.createTripRouter();
         this.skimsByMode.clear();
         this.travelTimesFromRegion.clear();
@@ -282,7 +289,10 @@ public IndexedDoubleMatrix2D getPeakSkim(String mode) {
             logger.info("Calculating skim matrix for mode " + mode +
                     " using " + Properties.get().main.numberOfThreads + " threads.");
             IndexedDoubleMatrix2D skim;
-            final MatsimSkimCreator matsimSkimCreator = new MatsimSkimCreator(matsimData);
+            final String csvDirectory = Properties.get().main.baseDirectory + "scenOutput/"
+                    + Properties.get().main.scenarioName + "/matsim/skims";
+            final MatsimSkimCreator matsimSkimCreator = new MatsimSkimCreator(matsimData,
+                    Paths.get(csvDirectory), skimOutputYear);
             switch (mode) {
                 case TransportMode.car:
                     skim = matsimSkimCreator.createCarSkim(zones.values(), Properties.get().main.numberOfThreads,
@@ -318,6 +328,7 @@ public TravelTimes duplicate() {
         matsimTravelTimesAndCosts.zones = this.zones;
         matsimTravelTimesAndCosts.regions = this.regions;
         matsimTravelTimesAndCosts.matsimData = matsimData;
+        matsimTravelTimesAndCosts.skimOutputYear = skimOutputYear;
         matsimTravelTimesAndCosts.tripRouter = matsimData.createTripRouter();
         matsimTravelTimesAndCosts.skimsByMode.putAll(this.skimsByMode);
         matsimTravelTimesAndCosts.travelTimesFromRegion.putAll(travelTimesFromRegion);
diff --git a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/PreviousYearPlan.java b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/PreviousYearPlan.java
new file mode 100644
index 000000000..88c5a48ef
--- /dev/null
+++ b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/PreviousYearPlan.java
@@ -0,0 +1,21 @@
+package de.tum.bgu.msm.matsim;
+
+import org.matsim.api.core.v01.Coord;
+import org.matsim.api.core.v01.Id;
+import org.matsim.api.core.v01.population.Plan;
+import org.matsim.utils.objectattributes.attributable.Attributes;
+
+public class PreviousYearPlan {
+
+    Id agentId;
+    Plan previousYearPlan;
+    int jobId;
+    int dwellingId;
+
+    public PreviousYearPlan(Id agentId, Plan previousYearPlan, int jobId, int dwellingId) {
+        this.agentId = agentId;
+        this.previousYearPlan = previousYearPlan;
+        this.jobId =  jobId;
+        this.dwellingId = dwellingId;
+    }
+}
diff --git a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/SimpleCommuteModeChoiceMatsimScenarioAssembler.java b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/SimpleCommuteModeChoiceMatsimScenarioAssembler.java
index 5da1888cd..fd1804f8a 100644
--- a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/SimpleCommuteModeChoiceMatsimScenarioAssembler.java
+++ b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/SimpleCommuteModeChoiceMatsimScenarioAssembler.java
@@ -2,8 +2,7 @@
 
 import de.tum.bgu.msm.container.DataContainer;
 import de.tum.bgu.msm.data.Day;
-import de.tum.bgu.msm.data.dwelling.Dwelling;
-import de.tum.bgu.msm.data.dwelling.RealEstateDataManager;
+import de.tum.bgu.msm.data.dwelling.*;
 import de.tum.bgu.msm.data.household.Household;
 import de.tum.bgu.msm.data.job.Job;
 import de.tum.bgu.msm.data.job.JobDataManager;
@@ -12,36 +11,38 @@
 import de.tum.bgu.msm.data.travelTimes.TravelTimes;
 import de.tum.bgu.msm.data.vehicle.VehicleType;
 import de.tum.bgu.msm.models.modeChoice.CommuteModeChoice;
-import de.tum.bgu.msm.models.modeChoice.CommuteModeChoiceMapping;
 import de.tum.bgu.msm.properties.Properties;
 import de.tum.bgu.msm.utils.SiloUtil;
 import org.apache.logging.log4j.LogManager;
 import org.apache.logging.log4j.Logger;
 import org.locationtech.jts.geom.Coordinate;
 import org.matsim.api.core.v01.Coord;
+import org.matsim.api.core.v01.Id;
 import org.matsim.api.core.v01.Scenario;
-import org.matsim.api.core.v01.TransportMode;
-import org.matsim.api.core.v01.population.Activity;
-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.api.core.v01.population.*;
 import org.matsim.core.config.Config;
 import org.matsim.core.gbl.MatsimRandom;
+import org.matsim.core.router.TripStructureUtils;
 import org.matsim.core.scenario.ScenarioUtils;
 
-import java.util.Map;
-
-import java.util.Random;
+import java.nio.file.Paths;
+import java.util.*;
 
 public class SimpleCommuteModeChoiceMatsimScenarioAssembler implements MatsimScenarioAssembler {
 
-    private final static Logger logger = LogManager.getLogger(SimpleMatsimScenarioAssembler.class);
+    private final static Logger logger = LogManager.getLogger(SimpleCommuteModeChoiceMatsimScenarioAssembler.class);
 
     private final DataContainer dataContainer;
     private final Properties properties;
     private final HandlingOfRandomness handlingOfRandomObject;
     private final CommuteModeChoice commuteModeChoice;
     private final Random random;
+
+    static final String DWELLING_ID_ATTR = "dwellingId";
+    static final String JOB_ID_ATTR = "jobId";
+
+    private record AgentWithMode(org.matsim.api.core.v01.population.Person matsimPerson, String mode) {}
+
     // yyyy I found this using the regular silo random number sequence.  In consequence, it was using different random numbers every time it was
     // called, in consequence picking different agents from silo.  This is not what we want.  --  It also picked different agents every time it ran.
     // No idea why; by design, the silo rnd number sequence should be deterministic.  However, we also have a randomly occuring binarySearch error, so
@@ -73,82 +74,166 @@ public SimpleCommuteModeChoiceMatsimScenarioAssembler(DataContainer dataContaine
 
     }
 
-    @Override
-    public Scenario assembleScenario(Config matsimConfig, int year, TravelTimes travelTimes) {
-        logger.info("Starting creating (mode-respecting, home-work-home) MATSim scenario.");
+        @Override
+        public Scenario assembleScenarioWithPreviousPlans(
+                Scenario oldScenario,
+        int year,
+        TravelTimes travelTimes,
+        HashMap<Id<org.matsim.api.core.v01.population.Person>, org.matsim.api.core.v01.population.Person> previousYearPopulation) {
 
-        switch( handlingOfRandomObject ){
-            case fromSilo:
-                break;
-            case localInstanceFromMatsimWithAlwaysSameSeed:
-                random.setSeed( 4711 );
-                // (note that we WANT this with the same random seed for every year when matsim is called.  Could, however, be made dependent on the silo
-                // seed, so that with a change of the silo seed it also changes the random seed here.  kai' jun'23)
+            logger.info("Starting MATSim scenario using previous year's plans.");
+            initRandom();
 
-                // yyyy this will probably not work.  A person at a certain position will not be the same person some silo years later, since
-                // intermediate persons may have been deleted and/or added.  In matsim4urbansim, I therefore went through the remembered matsim
-                // population, and pulled the corresponding persons from the land use population.  And then added randomly drawn additional persons
-                // until the right amount was found (not necessarily the same as in the 1st iteration).  kai, jun'23
+            // fresh scenario, like in assembleScenario; matsimData's scenario stays untouched
+            Scenario scenario = ScenarioUtils.loadScenario(oldScenario.getConfig());
+            Population matsimPopulation = scenario.getPopulation();
+            PopulationFactory pf = matsimPopulation.getFactory();
 
-                break;
-            default:
-                throw new IllegalStateException( "Unexpected value: " + handlingOfRandomObject );
-        }
+            double populationScalingFactor = properties.transportModel.matsimScaleFactor;
+            SiloMatsimUtils.checkSiloPropertiesAndMatsimConfigConsistency(oldScenario.getConfig(), properties);
 
-        double populationScalingFactor = properties.transportModel.matsimScaleFactor;
-        SiloMatsimUtils.checkSiloPropertiesAndMatsimConfigConsistency(matsimConfig, properties);
+            JobDataManager jobDataManager = dataContainer.getJobDataManager();
+            RealEstateDataManager realEstateDataManager = dataContainer.getRealEstateDataManager();
+            Map<Id<org.matsim.api.core.v01.population.Person>, String> source = new HashMap<>();
 
-        Scenario scenario = ScenarioUtils.loadScenario(matsimConfig);
-        Population matsimPopulation = scenario.getPopulation();
+            // 1. expected pop size: all currently employed silo persons
+            int shouldPop = 0;
 
-        JobDataManager jobDataManager = dataContainer.getJobDataManager();
-        RealEstateDataManager realEstateDataManager = dataContainer.getRealEstateDataManager();
+            for (Person p : dataContainer.getHouseholdDataManager().getPersons()) {
+                if (isEmployedWithJob(p)) {
+                    shouldPop++;
+                }
+            }
 
-        for (Household household: dataContainer.getHouseholdDataManager().getHouseholds()) {
-//            if (SiloUtil.getRandomNumberAsDouble() > populationScalingFactor) {
-            if ( random.nextDouble() > populationScalingFactor ) {
-                // re random see above.  kai, jun'23
+            shouldPop *= populationScalingFactor;
 
-                continue;
-            }
+            // 2. carry over previous-year persons who still exist and still employed
+            List<org.matsim.api.core.v01.population.Person> previous = new ArrayList<>(previousYearPopulation.values());
+            Collections.shuffle(previous, random);   // so that, if we have too many, the ones dropped are random
 
-            for (Person person : household.getPersons().values()) {
-                if (person.getOccupation() != Occupation.EMPLOYED || person.getJobId() == -2) { // i.e. person does not work
+            int kept = 0;
+            int replanned = 0;
+            int dropped = 0;
+            for (org.matsim.api.core.v01.population.Person matsimPerson : previous) {
+                if (kept >= shouldPop) break;        // replaces the old "shouldPop < actualPop" trimming
+
+                Person siloPerson = dataContainer.getHouseholdDataManager()
+                        .getPersonFromId(Integer.parseInt(matsimPerson.getId().toString()));
+
+                boolean exists = siloPerson != null;
+                boolean employed = exists && isEmployedWithJob(siloPerson);
+
+                String mode = null;
+                boolean moved = false;
+                boolean modeChanged = false;
+
+                if (employed) {
+                    mode = determineMode(realEstateDataManager, siloPerson, travelTimes);
+
+                    TripStructureUtils.Trip firstTrip = TripStructureUtils.getTrips(matsimPerson.getSelectedPlan()).getFirst();
+                    String previousMode = TripStructureUtils.getRoutingModeIdentifier().identifyMainMode(firstTrip.getTripElements());
+
+                    moved = homeOrWorkChanged(
+                            (int) matsimPerson.getAttributes().getAttribute("dwellingId"),
+                            (int) matsimPerson.getAttributes().getAttribute("jobId"),
+                            siloPerson.getHousehold().getDwellingId(),
+                            siloPerson.getJobId());
+                    modeChanged = !Objects.equals(mode, previousMode);
+                }
+
+                Outcome outcome = classify(exists, employed, moved, modeChanged);
+                if (outcome == Outcome.DROP) {
+                    dropped++;
                     continue;
                 }
 
-                PopulationFactory pf = matsimPopulation.getFactory();
+                if (outcome == Outcome.REPLAN) {
+                    replanned++;
+                    matsimPerson.getPlans().clear();
+                    matsimPerson.setSelectedPlan(null);
+                    setHomeWorkAttributes(matsimPerson, siloPerson);
+                    createPlan(realEstateDataManager, jobDataManager, siloPerson, mode, pf, matsimPerson);
+                }
 
-                final int nHhCars = (int) household.getVehicles().stream().filter( vv -> vv.getType().equals( VehicleType.CAR ) ).count();
-                org.matsim.api.core.v01.population.Person matsimAlterEgo = SiloMatsimUtils.createMatsimAlterEgo(pf, person, nHhCars );
-                matsimPopulation.addPerson(matsimAlterEgo);
+                source.put(matsimPerson.getId(), outcome.label);
+                matsimPopulation.addPerson(matsimPerson);   // old object, with old or new plan
+                kept++;
+            }
 
-                Dwelling dwelling = realEstateDataManager.getDwelling(household.getDwellingId());
-                CommuteModeChoiceMapping commuteModeChoiceMapping = commuteModeChoice.assignCommuteModeChoice(dwelling, travelTimes, household);
+            // 3. top up with employed silo persons who weren't in last year's population
+            int added = 0;
+            if (kept < shouldPop) {
+                ArrayList<Person> candidates = new ArrayList<>();
+                for (Person p : dataContainer.getHouseholdDataManager().getPersons()) {
+                    if (!isEmployedWithJob(p)) continue;
+                    if (previousYearPopulation.containsKey(Id.createPersonId(Integer.toString(p.getId())))) continue;
+                    candidates.add(p);
+                }
 
-                String mode = commuteModeChoiceMapping.getMode(person).mode;
+                if (candidates.size() < shouldPop - kept) {
+                    logger.warn("not enough silo agents to top off matsim for year {}", year);
+                }
+
+                // this is probably a redundant duplicate check since we calc shouldPop earlier by filtering out unemployed.
+                int needed = Math.min(shouldPop - kept, candidates.size());
+                for (Person siloPerson : getRandomSample(candidates, needed)) {
+                    AgentWithMode agent = createAndAddAgent(realEstateDataManager, travelTimes, siloPerson, pf, matsimPopulation);
+                    source.put(agent.matsimPerson().getId(), "new");
+                    createPlan(realEstateDataManager, jobDataManager, siloPerson, agent.mode(), pf, agent.matsimPerson());
+                    added++;
+                }
+            }
 
-                if (mode.equals(TransportMode.car)) {
-                    Coord dwellingCoord = getOrRandomlyChooseDwellingCoord(dwelling);
+            int previousSize = previous.size();
+            logger.info("previous=" + previousSize + " shouldPop=" + shouldPop + " kept=" + kept
+                    + " replanned=" + replanned + " dropped=" + dropped + " added=" + added
+                    + " trimmed=" + (previousSize - kept - dropped) + " final=" + matsimPopulation.getPersons().size());
 
-                    Job job = jobDataManager.getJobFromId(person.getJobId());
-                    Coord jobCoord = getOrRandomlyChooseJobCoordinate(job);
+            java.nio.file.Path csvPath = java.nio.file.Paths.get(
+                    properties.main.baseDirectory + "scenOutput/" + properties.main.scenarioName
+                            + "/matsim/" + year + "/assembledPlans.csv");
 
-                    createHWHPlanAndAddToAlterEgo(pf, matsimAlterEgo, dwellingCoord, job, jobCoord, TransportMode.car);
-                } else {
-                    // TODO MATSim expects the plans to be there at least in VspPlansCleaner.notifyBeforeMobsim(VspPlansCleaner.java:54).
-                    // Therefore, the intended switch one line below may not work; maybe with other settings
-                    // if (!properties.transportModel.onlySimulateCarTrips) {
-                        Coord dwellingCoord = getOrRandomlyChooseDwellingCoord(dwelling);
+            writePlansToCsv(csvPath, source, matsimPopulation);
+            return scenario;
+        }
 
-                        Job job = jobDataManager.getJobFromId(person.getJobId());
-                        Coord jobCoord = getOrRandomlyChooseJobCoordinate(job);
+    @Override
+    public Scenario assembleScenario(Config matsimConfig, int year, TravelTimes travelTimes) {
+        logger.info("Starting creating (mode-respecting, home-work-home) MATSim scenario.");
 
-                        createHWHPlanAndAddToAlterEgo(pf, matsimAlterEgo, dwellingCoord, job, jobCoord, mode);
-                    // }
+        initRandom();
+        double populationScalingFactor = properties.transportModel.matsimScaleFactor;
+        SiloMatsimUtils.checkSiloPropertiesAndMatsimConfigConsistency(matsimConfig, properties);
+        Scenario scenario = ScenarioUtils.loadScenario(matsimConfig);
+        // at this pint - population is still 0!
+        Population matsimPopulation = scenario.getPopulation();
+
+        JobDataManager jobDataManager = dataContainer.getJobDataManager();
+        RealEstateDataManager realEstateDataManager = dataContainer.getRealEstateDataManager();
+        PopulationFactory pf = matsimPopulation.getFactory();
+        Map<Id<org.matsim.api.core.v01.population.Person>, String> source = new HashMap<>();
+
+        // first year: scale number of households based on populationScaleFactor
+        // for each person in this sampled num of houses, generate a plan with car or other mode
+        // decided by probability
+        for (Household household: dataContainer.getHouseholdDataManager().getHouseholds()) {
+            if ( random.nextDouble() > populationScalingFactor ) {
+                // re random see above.  kai, jun'23
+                continue;
+            }
+            for (Person person : household.getPersons().values()) {
+                if (!isEmployedWithJob(person)) {
+                    continue;
                 }
+                AgentWithMode agent = createAndAddAgent(realEstateDataManager, travelTimes, person, pf, matsimPopulation);
+                source.put(agent.matsimPerson().getId(), "new");
+                createPlan(realEstateDataManager, jobDataManager, person, agent.mode(), pf, agent.matsimPerson());
             }
         }
+
+        writePlansToCsv(Paths.get(properties.main.baseDirectory + "scenOutput/" + properties.main.scenarioName
+                + "/matsim/" + year + "/assembledPlans.csv"), source, matsimPopulation);
+
         logger.info("Finished creating MATSim scenario.");
         return scenario;
     }
@@ -221,4 +306,144 @@ private int defineWorkEndTime(Job job, int departureTime) {
         return departureTime + job.getWorkingTimeInSeconds().orElse(Math.max(0, (int) (8*3600 + random.nextGaussian() * 3600)));
         // re random see above.  kai, jun'23
     }
+
+    private <T> ArrayList<T> getRandomSample(ArrayList<T> source, int n) {
+        if (n > source.size()) throw new IllegalArgumentException("n > size");
+        List<T> copy = new ArrayList<>(source);
+        Collections.shuffle(copy, random);
+        return new ArrayList<>(copy.subList(0, n));
+    }
+
+    enum Outcome {
+        DROP("dropped"), KEEP("kept"), REPLAN("replanned");
+        final String label;
+        Outcome(String label) { this.label = label; }
+    }
+
+    static public Outcome classify(boolean stillExists, boolean employedWithJob,
+                            boolean homeOrWorkChanged, boolean modeChanged) {
+        if (!stillExists || !employedWithJob) return Outcome.DROP;
+        return (homeOrWorkChanged || modeChanged) ? Outcome.REPLAN : Outcome.KEEP;
+    }
+
+    static boolean homeOrWorkChanged(int oldDwellingId, int oldJobId, int newDwellingId, int newJobId) {
+        return oldDwellingId != newDwellingId || oldJobId != newJobId;
+    }
+
+
+    private void initRandom() {
+        switch( handlingOfRandomObject ){
+            case fromSilo:
+                break;
+            case localInstanceFromMatsimWithAlwaysSameSeed:
+                random.setSeed( 4711 );
+                // (note that we WANT this with the same random seed for every year when matsim is called.  Could, however, be made dependent on the silo
+                // seed, so that with a change of the silo seed it also changes the random seed here.  kai' jun'23)
+
+                // yyyy this will probably not work.  A person at a certain position will not be the same person some silo years later, since
+                // intermediate persons may have been deleted and/or added.  In matsim4urbansim, I therefore went through the remembered matsim
+                // population, and pulled the corresponding persons from the land use population.  And then added randomly drawn additional persons
+                // until the right amount was found (not necessarily the same as in the 1st iteration).  kai, jun'23
+
+                break;
+            default:
+                throw new IllegalStateException( "Unexpected value: " + handlingOfRandomObject );
+        }
+    }
+
+
+    private void createPlan(RealEstateDataManager realEstateDataManager,
+                            JobDataManager jobDataManager, Person siloPerson,
+                            String mode, PopulationFactory pf,
+                            org.matsim.api.core.v01.population.Person matsimPerson) {
+
+        Coord dwellingCoord = getOrRandomlyChooseDwellingCoord(realEstateDataManager.getDwelling(siloPerson.getHousehold().getDwellingId()));
+        Job job = jobDataManager.getJobFromId(siloPerson.getJobId());
+        Coord jobCoord = getOrRandomlyChooseJobCoordinate(job);
+        createHWHPlanAndAddToAlterEgo(pf, matsimPerson, dwellingCoord, job, jobCoord, mode);
+
+    }
+
+    private void writePlansToCsv(java.nio.file.Path csvPath, Map<Id<org.matsim.api.core.v01.population.Person>, String> source, Population matsimPopulation) {
+        try {
+            java.nio.file.Files.createDirectories(csvPath.getParent());
+            try (java.io.BufferedWriter w = java.nio.file.Files.newBufferedWriter(csvPath)) {
+                w.write("personId,source,nPlans,mode,homeX,homeY,workX,workY,homeId,workId,homeEndTime,workEndTime");
+                w.newLine();
+
+                for (org.matsim.api.core.v01.population.Person p : matsimPopulation.getPersons().values()) {
+                    Plan plan = p.getSelectedPlan();
+                    if (plan == null) {
+                        w.write(p.getId() + "," + source.getOrDefault(p.getId(), "unknown") + "," + p.getPlans().size() + ",,,,,,,,,");
+                        w.newLine();
+                        continue;
+                    }
+
+                    List<Activity> acts = new ArrayList<>();
+                    String legMode = "";
+                    for (PlanElement pe : plan.getPlanElements()) {
+                        if (pe instanceof Activity a) acts.add(a);
+                        else if (pe instanceof Leg l && legMode.isEmpty()) legMode = l.getMode();
+                    }
+                    Activity home = acts.size() > 0 ? acts.get(0) : null;
+                    Activity work = acts.size() > 1 ? acts.get(1) : null;
+
+                    String dwellingId = String.valueOf(p.getAttributes().getAttribute("dwellingId"));
+                    String jobId      = String.valueOf(p.getAttributes().getAttribute("jobId"));
+
+                    w.write(String.join(",",
+                            p.getId().toString(),
+                            source.getOrDefault(p.getId(), "unknown"),
+                            Integer.toString(p.getPlans().size()),
+                            legMode,
+                            home != null ? Double.toString(home.getCoord().getX()) : "",
+                            home != null ? Double.toString(home.getCoord().getY()) : "",
+                            work != null ? Double.toString(work.getCoord().getX()) : "",
+                            work != null ? Double.toString(work.getCoord().getY()) : "",
+                            dwellingId,
+                            jobId,
+                            home != null ? Double.toString(home.getEndTime().orElse(-1)) : "",
+                            work != null ? Double.toString(work.getEndTime().orElse(-1)) : ""));
+                    w.newLine();
+                }
+            }
+            logger.info("Wrote assembled plans to " + csvPath);
+        } catch (java.io.IOException e) {
+            logger.error("Could not write plans csv", e);
+        }
+    }
+
+    private AgentWithMode createAndAddAgent(RealEstateDataManager realEstateDataManager, TravelTimes travelTimes, Person siloPerson, PopulationFactory pf, Population matsimPopulation) {
+        Household household = siloPerson.getHousehold();
+        int nHhCars = (int) household.getVehicles().stream()
+                .filter(v -> v.getType().equals(VehicleType.CAR)).count();
+        org.matsim.api.core.v01.population.Person matsimPerson =
+                SiloMatsimUtils.createMatsimAlterEgo(pf, siloPerson, nHhCars);
+        setHomeWorkAttributes(matsimPerson, siloPerson);
+        matsimPopulation.addPerson(matsimPerson);
+
+        String mode = determineMode(realEstateDataManager, siloPerson, travelTimes);
+
+        return new AgentWithMode(matsimPerson, mode);
+    }
+
+    static boolean isEmployedWithJob(Occupation occupation, int jobId) {
+        return occupation == Occupation.EMPLOYED && jobId != -2;
+    }
+
+    private boolean isEmployedWithJob(Person p) {          // keeps call sites unchanged
+        return isEmployedWithJob(p.getOccupation(), p.getJobId());
+    }
+
+    private String determineMode(RealEstateDataManager realEstateDataManager, Person siloPerson, TravelTimes travelTimes) {
+        return commuteModeChoice.assignCommuteModeChoice(
+                realEstateDataManager.getDwelling(siloPerson.getHousehold().getDwellingId()),
+                travelTimes, siloPerson.getHousehold()).getMode(siloPerson).mode;
+    }
+
+    private void setHomeWorkAttributes(org.matsim.api.core.v01.population.Person matsimPerson, Person siloPerson) {
+        matsimPerson.getAttributes().putAttribute(DWELLING_ID_ATTR, siloPerson.getHousehold().getDwellingId());
+        matsimPerson.getAttributes().putAttribute(JOB_ID_ATTR, siloPerson.getJobId());
+    }
+
 }
diff --git a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/StandaloneSkimCreator.java b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/StandaloneSkimCreator.java
index 2cca3ed69..49274ab5c 100644
--- a/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/StandaloneSkimCreator.java
+++ b/extensions/matsim2silo/src/main/java/de/tum/bgu/msm/matsim/StandaloneSkimCreator.java
@@ -23,6 +23,8 @@
 import org.matsim.pt.transitSchedule.api.TransitSchedule;
 import org.matsim.pt.transitSchedule.api.TransitScheduleReader;
 
+import java.io.File;
+import java.nio.file.Path;
 import java.util.Collection;
 import java.util.List;
 import java.util.Map;
@@ -78,7 +80,11 @@ public static void createSkims(int threads, String shapePath, String networkPath
         FreespeedTravelTimeAndDisutility freespeed = new FreespeedTravelTimeAndDisutility(config.scoring());
         data.update(freespeed, freespeed);
 
-        final MatsimSkimCreator matsimSkimCreator = new MatsimSkimCreator(data);
+        final File outputFile = new File(outputPath).getAbsoluteFile();
+        final Path csvOutputDirectory = outputFile.getParentFile() == null
+                ? new File(".").toPath()
+                : outputFile.getParentFile().toPath();
+        final MatsimSkimCreator matsimSkimCreator = new MatsimSkimCreator(data, csvOutputDirectory);
         final SkimTravelTimes skimTravelTimes = new SkimTravelTimes();
 
         logger.info("Creating car skim...");
diff --git a/extensions/matsim2silo/src/test/java/de/tum/bgu/msm/matsim/AssembledPlansYearOverYearTest.java b/extensions/matsim2silo/src/test/java/de/tum/bgu/msm/matsim/AssembledPlansYearOverYearTest.java
new file mode 100644
index 000000000..5d5c7b05a
--- /dev/null
+++ b/extensions/matsim2silo/src/test/java/de/tum/bgu/msm/matsim/AssembledPlansYearOverYearTest.java
@@ -0,0 +1,123 @@
+package de.tum.bgu.msm.matsim;
+
+import org.junit.jupiter.api.Test;
+import org.junit.jupiter.api.Tag;
+
+import java.io.IOException;
+import java.nio.file.Files;
+import java.nio.file.Path;
+import java.util.*;
+import java.util.stream.Stream;
+
+import static org.junit.jupiter.api.Assertions.*;
+
+@Tag("integration")
+public class AssembledPlansYearOverYearTest {
+    enum AutofreiScenario {
+        BASE,
+        SMALL,
+        LARGE
+    }
+
+    enum PtScenario {
+        //        NONE,
+        NES,
+        RINGSX
+    }
+
+    //    static String scenario = "base";
+    static AutofreiScenario autofreiScenario = AutofreiScenario.LARGE;
+    static PtScenario ptScenario = PtScenario.NES;
+
+    static String prefix = "2026-10-23a__";
+    static String scenario = prefix + "autofrei-" + autofreiScenario + "__pt-" + ptScenario;
+
+    record Row(String id, String source, String mode, String homeX, String homeY,
+               String workX, String workY, String homeId, String workId) {}
+
+    // point this at a finished run (or run the model in @BeforeAll)
+    // AT THE MOMENT JUST SETTING AS ABSOLUTE PATH: PLEASE CHANGE!
+    static final Path MATSIM_DIR = Path.of("/home/brendan/git/silo/useCases/fabiland/scenario/scenOutput/" + scenario + "/matsim");
+
+    /**
+     * year -> rows, sorted by year
+     */
+    static SortedMap<Integer, Map<String, Row>> readAllYears() throws IOException {
+//        assertTrue(Files.isDirectory(MATSIM_DIR),
+//                "matsim dir not found: " + MATSIM_DIR.toAbsolutePath());
+        SortedMap<Integer, Map<String, Row>> result = new TreeMap<>();
+        try (Stream<Path> dirs = Files.list(MATSIM_DIR)) {
+            for (Path dir : dirs.toList()) {
+                Path csv = dir.resolve("assembledPlans.csv");
+                if (!Files.isRegularFile(csv)) continue;
+                try {
+                    result.put(Integer.parseInt(dir.getFileName().toString()), readCsv(csv));
+                } catch (NumberFormatException ignored) { /* not a year folder */ }
+            }
+        }
+        return result;
+    }
+
+    static Map<String, Row> readCsv(Path csv) throws IOException {
+        List<String> lines = Files.readAllLines(csv);
+        Map<String, Row> rows = new LinkedHashMap<>();
+        for (String line : lines.subList(1, lines.size())) {           // skip header
+            String[] c = line.split(",", -1);                          // keep empty trailing fields
+            // header: personId,source,nPlans,mode,homeX,homeY,workX,workY,homeEndTime,workEndTime
+
+            Row row = new Row(c[0], c[1], c[3], c[4], c[5], c[6], c[7], c[8], c[9]);
+            assertNull(rows.put(row.id(), row), "duplicate person id " + row.id() + " in " + csv);
+        }
+        return rows;
+    }
+
+    @Test
+    public void everyYearIsConsistentWithThePreviousOne() throws IOException {
+        SortedMap<Integer, Map<String, Row>> years = readAllYears();
+        assertTrue(years.containsKey(0), "year 0 has no assembledPlans.csv");
+        assertEquals(11, years.size());
+
+        Integer prevYear = null;
+        for (var entry : years.entrySet()) {
+            int year = entry.getKey();
+            Map<String, Row> curr = entry.getValue();
+
+            if (prevYear == null) {
+                curr.values().forEach(r -> assertEquals("new", r.source()));
+            }
+
+            // no plan-less rows, every row has a known source
+            for (Row r : curr.values()) {
+                assertTrue(Set.of("kept", "replanned", "new").contains(r.source()),
+                        year + ": unexpected source " + r.source() + " for " + r.id());
+                assertFalse(r.mode().isEmpty(), year + ": no plan for " + r.id());
+            }
+
+            if (prevYear != null) {
+                Map<String, Row> prev = years.get(prevYear);
+                for (Row r : curr.values()) {
+                    switch (r.source()) {
+                        case "kept", "replanned" ->
+                                assertTrue(prev.containsKey(r.id()), year + ": " + r.id() + " not in " + prevYear);
+                        case "new" ->
+                                assertFalse(prev.containsKey(r.id()), year + ": " + r.id() + " is new but existed");
+                    }
+                    if (r.source().equals("kept")) {
+                        Row p = prev.get(r.id());
+                        // workID stays the same, but Coords change? these coord assertions always fail.
+                        // I assume that's because of getOrRandomlyChooseDwellingCoord & getOrRandomlyChooseJobCoord
+//                        assertEquals(p.homeX(), r.homeX(), year + ": home moved for kept " + r.id());
+//                        assertEquals(p.homeY(), r.homeY());
+//                        assertEquals(p.workX(), r.workX(), year + ": work moved for kept " + r.id());
+//                        assertEquals(p.workY(), r.workY());
+                        assertEquals(p.homeId(), r.homeId(), year + ": dwelling changed for kept " + r.id());
+                        assertEquals(p.workId(), r.workId(), year + ": job changed for kept " + r.id());
+                    }
+                }
+                long carried = curr.values().stream().filter(r -> !r.source().equals("new")).count();
+                assertTrue(carried <= prev.size());
+            }
+            prevYear = year;
+        }
+    }
+}
diff --git a/extensions/matsim2silo/src/test/java/de/tum/bgu/msm/matsim/SimpleCommuteModeChoiceMatsimScenarioAssemblerTest.java b/extensions/matsim2silo/src/test/java/de/tum/bgu/msm/matsim/SimpleCommuteModeChoiceMatsimScenarioAssemblerTest.java
new file mode 100644
index 000000000..522f4dff6
--- /dev/null
+++ b/extensions/matsim2silo/src/test/java/de/tum/bgu/msm/matsim/SimpleCommuteModeChoiceMatsimScenarioAssemblerTest.java
@@ -0,0 +1,45 @@
+package de.tum.bgu.msm.matsim;
+
+import de.tum.bgu.msm.data.person.Occupation;
+import org.junit.jupiter.api.Test;
+
+import static de.tum.bgu.msm.matsim.SimpleCommuteModeChoiceMatsimScenarioAssembler.*;
+import static org.junit.jupiter.api.Assertions.*;
+
+class SimpleCommuteModeChoiceMatsimScenarioAssemblerTest {
+
+    @Test
+    void onlyEmployedWithRealJobCounts() {
+        assertTrue(isEmployedWithJob(Occupation.EMPLOYED, 5));
+        assertFalse(isEmployedWithJob(Occupation.EMPLOYED, -2));
+        assertFalse(isEmployedWithJob(Occupation.UNEMPLOYED, 5));
+        assertFalse(isEmployedWithJob(Occupation.UNEMPLOYED, -2));
+    }
+
+    @Test
+    void homeOrWorkChange() {
+        assertFalse(homeOrWorkChanged(1, 2, 1, 2));
+        assertTrue(homeOrWorkChanged(1, 2, 9, 2));
+        assertTrue(homeOrWorkChanged(1, 2, 1, 9));
+    }
+
+    @Test
+    void classifyCoversAllCases() {
+        assertEquals(Outcome.DROP,   classify(false, true,  false, false));
+        assertEquals(Outcome.DROP,   classify(true,  false, false, false));
+        assertEquals(Outcome.KEEP,   classify(true,  true,  false, false));
+        assertEquals(Outcome.REPLAN, classify(true,  true,  true,  false));
+        assertEquals(Outcome.REPLAN, classify(true,  true,  false, true));
+    }
+
+
+//    @Test
+//    void randomSampleIsDistinctSubsetAndDeterministic() {
+//        List<Integer> src = IntStream.range(0, 50).boxed().toList();
+//        List<Integer> a = getRandomSample(src, 10, new Random(1));
+//        assertEquals(10, new HashSet<>(a).size());
+//        assertTrue(src.containsAll(a));
+//        assertEquals(a, getRandomSample(src, 10, new Random(1)));
+//        assertThrows(IllegalArgumentException.class, () -> getRandomSample(src, 51, new Random(1)));
+//    }
+}
\ No newline at end of file
diff --git a/pom.xml b/pom.xml
index df3e1422a..3c4ffd9a9 100644
--- a/pom.xml
+++ b/pom.xml
@@ -102,9 +102,9 @@
 		<module>useCases</module>
 		<module>siloCore</module>
 		<module>synthetic-population</module>
-        <module>analysis</module>
-        <module>extensions</module>
-    </modules>
+		<module>analysis</module>
+		<module>extensions</module>
+	</modules>
 
 	<build>
 		<defaultGoal>compiler:compile</defaultGoal>
@@ -154,7 +154,7 @@
 
 					<!--execute each test in separate JVM (slow, but otherwise it does not-->
 					<!--work, and no problem as long as there are only integration tests):-->
-<!--					<forkMode>always</forkMode>-->
+					<!--					<forkMode>always</forkMode>-->
 					<forkCount>1</forkCount>
 					<reuseForks>false</reuseForks>
 
@@ -180,11 +180,17 @@
 			<version>2.24.0</version> <!-- Replace with the correct version -->
 		</dependency>
 
+		<dependency>
+			<groupId>tech.tablesaw</groupId>
+			<artifactId>tablesaw-core</artifactId>
+			<version>0.43.1</version>
+		</dependency>
+
 	</dependencies>
 
 	<properties>
 		<geotools.version>31.3</geotools.version>
-<!--		<matsimVersion>2025.0-PR3643</matsimVersion>-->
+		<!--		<matsimVersion>2025.0-PR3643</matsimVersion>-->
 		<matsimVersion>2025.0</matsimVersion>
 		<jackson.version>2.9.5</jackson.version>
 		<junit.version>5.10.3</junit.version>
diff --git a/siloCore/src/main/java/de/tum/bgu/msm/data/job/JobDataManagerImpl.java b/siloCore/src/main/java/de/tum/bgu/msm/data/job/JobDataManagerImpl.java
index 0993cc784..29c286ad4 100644
--- a/siloCore/src/main/java/de/tum/bgu/msm/data/job/JobDataManagerImpl.java
+++ b/siloCore/src/main/java/de/tum/bgu/msm/data/job/JobDataManagerImpl.java
@@ -366,6 +366,8 @@ public Job findVacantJob(Zone homeZone, Collection<Region> regions) {
                 }
             }
             if (regionSampler.getCumulatedProbability() == 0) {
+                // Discard zero-weight entries before adding regions with fallback weights.
+                regionSampler = new Sampler<>(regions.size(), Region.class, SiloUtil.getRandomObject());
                 // could not find job in reasonable distance. Person will have to commute far and is likely to relocate in the future
                 for (Region reg : regions) {
                     if (getNumberOfVacantJobsByRegion(reg.getId()) > 0) {
diff --git a/siloCore/src/main/java/de/tum/bgu/msm/models/demography/marriage/MarriageModelImpl.java b/siloCore/src/main/java/de/tum/bgu/msm/models/demography/marriage/MarriageModelImpl.java
index c7012cf79..3e94af44f 100644
--- a/siloCore/src/main/java/de/tum/bgu/msm/models/demography/marriage/MarriageModelImpl.java
+++ b/siloCore/src/main/java/de/tum/bgu/msm/models/demography/marriage/MarriageModelImpl.java
@@ -316,7 +316,7 @@ private Household chooseRelocationTarget(Person partner1, Person partner2, House
         return moveTo;
     }
 
-    private boolean moveTogether(Person person1, Person person2, Household moveTo) {
+    protected boolean moveTogether(Person person1, Person person2, Household moveTo) {
 
         movePerson(person1, moveTo);
         movePerson(person2, moveTo);
diff --git a/siloCore/src/main/java/de/tum/bgu/msm/models/relocation/moves/MovesModelImpl.java b/siloCore/src/main/java/de/tum/bgu/msm/models/relocation/moves/MovesModelImpl.java
index 6f2e09f6a..59ea97e6a 100644
--- a/siloCore/src/main/java/de/tum/bgu/msm/models/relocation/moves/MovesModelImpl.java
+++ b/siloCore/src/main/java/de/tum/bgu/msm/models/relocation/moves/MovesModelImpl.java
@@ -15,7 +15,6 @@
 import de.tum.bgu.msm.io.output.YearByYearCsvModelTracker;
 import de.tum.bgu.msm.models.AbstractModel;
 import de.tum.bgu.msm.properties.Properties;
-import de.tum.bgu.msm.properties.modules.TransportModelPropertiesModule;
 import de.tum.bgu.msm.util.concurrent.ConcurrentExecutor;
 import de.tum.bgu.msm.utils.SampleException;
 import de.tum.bgu.msm.utils.Sampler;
@@ -25,6 +24,7 @@
 import org.apache.logging.log4j.Logger;
 import org.locationtech.jts.geom.Coordinate;
 
+import java.io.PrintWriter;
 import java.nio.file.Path;
 import java.nio.file.Paths;
 import java.util.*;
@@ -32,6 +32,7 @@
 import java.util.concurrent.ConcurrentHashMap;
 import java.util.concurrent.ConcurrentLinkedQueue;
 import java.util.concurrent.CyclicBarrier;
+import java.util.stream.Collectors;
 
 /**
  * @author Nico
@@ -60,6 +61,16 @@ public class MovesModelImpl extends AbstractModel implements MovesModel {
     private final Map<Integer, Double > sumOfSatisfactionsByZone = new LinkedHashMap<>();
     private YearByYearCsvModelTracker relocationTracker;
 
+    PrintWriter pwd;
+    boolean headerWritten = false;
+    public final Map<Integer, Double> sumAvgSatisfactionByZone = new LinkedHashMap<>();
+    public final Map<Integer, Double> sumCurrentUtilByZone = new LinkedHashMap<>();
+    public final Map<Integer, Double> sumPropByZone = new LinkedHashMap<>();
+    public final Map<Integer, Integer> cntByZone = new LinkedHashMap<>();
+    public final Map<Integer, Integer> cntDefaultSatisfaction = new LinkedHashMap<>();
+
+
+
 
     public MovesModelImpl(DataContainer dataContainer, Properties properties, MovesStrategy movesStrategy,
                           HousingStrategy housingStrategy, Random random) {
@@ -82,6 +93,9 @@ public void setup() {
         String header = new StringJoiner(",").add("hh").add("oldDdd").add("newDd").add("oldX").add("oldY").add("newX").add("newY").add("oldZone").add("newZone").add("autos").add("licenses").add("workers").toString();
         Path basePath = Paths.get(properties.main.baseDirectory).resolve("scenOutput").resolve(properties.main.scenarioName).resolve("siloResults/relocation");
         relocationTracker = new YearByYearCsvModelTracker(basePath, "relocation", header);
+        pwd = SiloUtil.openFileForSequentialWriting( Paths.get(properties.main.baseDirectory).resolve("scenOutput").resolve(properties.main.scenarioName).resolve("avgSatisfactionByZone.csv").toString(), false);
+//        pwd.println("zone,avgSatisfaction1,avgSatisfaction2,avgCurrentUtil,avgProp,defaultProp,year");
+
     }
 
     @Override
@@ -111,16 +125,73 @@ public List<MoveEvent> getEventsForCurrentYear(int year) {
     @Override
     public void endYear(int year) {
         UtilityUtils.endYear();
+
+        assert housingStrategy instanceof SimpleCommuteModeChoiceHousingStrategyImpl;
+
+        SimpleCommuteModeChoiceHousingStrategyImpl housingStrategy1 = (SimpleCommuteModeChoiceHousingStrategyImpl) housingStrategy;
+
+
+
+        for (int zoneNumber = 1; zoneNumber <= sumOfSatisfactionsByZone.size(); zoneNumber++) {
+            Map<String, Object> row = new LinkedHashMap<>();
+            row.put("zone", zoneNumber);
+            row.put("avgSatisfaction1", sumOfSatisfactionsByZone.getOrDefault(zoneNumber,0.)/householdsByZone.getOrDefault(zoneNumber,0));
+            row.put("avgSatisfaction2", sumAvgSatisfactionByZone.getOrDefault(zoneNumber,0.) / cntByZone.getOrDefault(zoneNumber,0));
+            row.put("avgCurrentUtil", sumCurrentUtilByZone.getOrDefault(zoneNumber,0.) / cntByZone.getOrDefault(zoneNumber,0));
+            row.put("avgProp", sumPropByZone.getOrDefault(zoneNumber,0.) / cntByZone.getOrDefault(zoneNumber,0));
+            row.put("defaultProp", ((double) cntDefaultSatisfaction.getOrDefault(zoneNumber,0)) / cntByZone.getOrDefault(zoneNumber,0));
+            row.put("ddQualityUtility",  housingStrategy1.ddQualityUtilityMap.getOrDefault(zoneNumber,0.) / housingStrategy1.cntMap.getOrDefault(zoneNumber,0));
+            row.put("ddSizeUtility", housingStrategy1.ddSizeUtilityMap.getOrDefault(zoneNumber,0.) / housingStrategy1.cntMap.getOrDefault(zoneNumber,0));
+            row.put("ddAutoAccessibilityUtility", housingStrategy1.ddAutoAccessibilityUtilityMap.getOrDefault(zoneNumber,0.) / housingStrategy1.cntMap.getOrDefault(zoneNumber,0));
+            row.put("transitAccessibilityUtility", housingStrategy1.transitAccessibilityUtilityMap.getOrDefault(zoneNumber,0.) / housingStrategy1.cntMap.getOrDefault(zoneNumber,0));
+            row.put("ddPriceUtility", housingStrategy1.ddPriceUtilityMap.getOrDefault(zoneNumber,0.) / housingStrategy1.cntMap.getOrDefault(zoneNumber,0));
+            row.put("workDistanceUtility", housingStrategy1.workDistanceUtilityMap.getOrDefault(zoneNumber,0.) / housingStrategy1.cntMap.getOrDefault(zoneNumber,0));
+            row.put("commuteModeCarCnt", ((double) housingStrategy1.commuteModeCarCntMap.getOrDefault(zoneNumber,0)) / housingStrategy1.commuteModeAllCntMap.getOrDefault(zoneNumber,0));
+            row.put("commuteModePtCnt", ((double) housingStrategy1.commuteModePtCntMap.getOrDefault(zoneNumber, 0)) / housingStrategy1.commuteModeAllCntMap.getOrDefault(zoneNumber, 0));
+            row.put("commuteModeCarUtilSum",  housingStrategy1.commuteModeCarUtilSumMap.getOrDefault(zoneNumber,0.) / housingStrategy1.commuteModeCarCntMap.getOrDefault(zoneNumber,0));
+            row.put("commuteModePtUtilSum", housingStrategy1.commuteModePtUtilSumMap.getOrDefault(zoneNumber,0.) / housingStrategy1.commuteModePtCntMap.getOrDefault(zoneNumber,0));
+
+            row.put("year", year);
+
+
+            if (!headerWritten) {
+                pwd.println(String.join(",", row.keySet()));
+                headerWritten = true;
+            }
+
+            pwd.println(row.values().stream()
+                    .map(String::valueOf)
+                    .collect(Collectors.joining(",")));
+        }
+
+
+        sumAvgSatisfactionByZone.clear();
+        sumCurrentUtilByZone.clear();
+        sumPropByZone.clear();
+        cntByZone.clear();
+        cntDefaultSatisfaction.clear();
+
+        housingStrategy1.ddQualityUtilityMap.clear();
+        housingStrategy1.ddSizeUtilityMap.clear();
+        housingStrategy1.ddAutoAccessibilityUtilityMap.clear();
+        housingStrategy1.transitAccessibilityUtilityMap.clear();
+        housingStrategy1.ddPriceUtilityMap.clear();
+        housingStrategy1.workDistanceUtilityMap.clear();
+        housingStrategy1.cntMap.clear();
+
+        housingStrategy1.commuteModeAllCntMap.clear();
+        housingStrategy1.commuteModePtCntMap.clear();
+        housingStrategy1.commuteModeCarCntMap.clear();
+        housingStrategy1.commuteModeCarUtilSumMap.clear();
+        housingStrategy1.commuteModePtUtilSumMap.clear();
+
     }
 
     @Override
     public void endSimulation() {
         relocationTracker.end();
-//        try {
-//            fileWriter.close();
-//        } catch (IOException e) {
-//            e.printStackTrace();
-//        }
+
+        pwd.close();
     }
 
     /**
@@ -225,7 +296,7 @@ public int searchForNewDwelling(Household household) {
         }
 
         // Step 2: select vacant dwelling in selected region
-        List<Dwelling> vacantDwellings
+            List<Dwelling> vacantDwellings
                 = new ArrayList<>(dataContainer.getRealEstateDataManager().getListOfVacantDwellingsInRegion(selectedRegion.getId()));
         if (vacantDwellings.isEmpty()) {
             return -1;
@@ -276,13 +347,28 @@ public int searchForNewDwelling(Household household) {
     private boolean moveOrNot(Household household) {
         HouseholdType hhType = household.getHouseholdType();
         Dwelling dd = dataContainer.getRealEstateDataManager().getDwelling(household.getDwellingId());
+        int zone = dd.getZoneId();
         if (!housingStrategy.isHouseholdEligibleToLiveHere(household, dd)) {
             return true;
         }
         final double currentUtil = satisfactionByHousehold.get(household.getId());
+
+        if (!averageHousingSatisfaction.containsKey(hhType)) {
+            cntDefaultSatisfaction.merge(dd.getZoneId(), 1, Integer::sum);
+        }
         final double avgSatisfaction = averageHousingSatisfaction.getOrDefault(hhType, currentUtil);
 
+
+
         final double prop = movesStrategy.getMovingProbability(avgSatisfaction, currentUtil);
+
+
+        cntByZone.merge(dd.getZoneId(), 1, Integer::sum);
+
+        sumAvgSatisfactionByZone.merge(dd.getZoneId(), avgSatisfaction, Double::sum);
+        sumCurrentUtilByZone.merge(dd.getZoneId(), currentUtil, Double::sum);
+        sumPropByZone.merge(dd.getZoneId(), prop, Double::sum);
+
         return this.random.nextDouble() <= prop;
     }
 
@@ -468,7 +554,10 @@ public Map<Integer, Integer> getHouseholdsByZone() {
         return householdsByZone;
     }
 
+
     public Map<Integer, Double> getSumOfSatisfactionsByZone() {
         return sumOfSatisfactionsByZone;
     }
 }
+
+
diff --git a/siloCore/src/main/java/de/tum/bgu/msm/models/relocation/moves/SimpleCommuteModeChoiceHousingStrategyImpl.java b/siloCore/src/main/java/de/tum/bgu/msm/models/relocation/moves/SimpleCommuteModeChoiceHousingStrategyImpl.java
index bf1ccbe54..794afb218 100644
--- a/siloCore/src/main/java/de/tum/bgu/msm/models/relocation/moves/SimpleCommuteModeChoiceHousingStrategyImpl.java
+++ b/siloCore/src/main/java/de/tum/bgu/msm/models/relocation/moves/SimpleCommuteModeChoiceHousingStrategyImpl.java
@@ -33,6 +33,21 @@ public class SimpleCommuteModeChoiceHousingStrategyImpl implements HousingStrate
     // used the hardcoded "Dwelling" type.  So we can as well specify it.  kai, jun'26
 
     private final static Logger logger = LogManager.getLogger(SimpleCommuteModeChoiceHousingStrategyImpl.class);
+    public Map<Integer, Double> ddQualityUtilityMap = new HashMap<>();
+    public Map<Integer, Double> workDistanceUtilityMap = new HashMap<>();
+    public Map<Integer, Double> ddPriceUtilityMap = new HashMap<>();
+    public Map<Integer, Double> transitAccessibilityUtilityMap = new HashMap<>();
+    public Map<Integer, Double> ddAutoAccessibilityUtilityMap = new HashMap<>();
+    public Map<Integer, Double> ddSizeUtilityMap = new HashMap<>();
+    public Map<Integer, Integer> cntMap = new HashMap<>();
+
+
+    public Map<Integer, Integer> commuteModeCarCntMap = new HashMap<>();
+    public Map<Integer, Double> commuteModeCarUtilSumMap = new HashMap<>();
+    public Map<Integer, Integer> commuteModePtCntMap = new HashMap<>();
+    public Map<Integer, Double> commuteModePtUtilSumMap = new HashMap<>();
+    public Map<Integer, Integer> commuteModeAllCntMap = new HashMap<>();
+
 
     private enum Normalizer {
         /**
@@ -145,10 +160,36 @@ public double calculateHousingUtility(Household hh, Dwelling dwelling) {
         for (Person pp : hh.getPersons().values()) {
             if (pp.getOccupation() == Occupation.EMPLOYED && pp.getJobId() != -2) {
 
-                workDistanceUtility *= commuteModeChoiceMapping.getMode(pp).utility;
+
+                CommuteModeChoiceMapping.CommuteMode mode = commuteModeChoiceMapping.getMode(pp);
+                workDistanceUtility *= mode.utility;
+
+
+                commuteModeAllCntMap.merge(dwelling.getZoneId(), 1, Integer::sum);
+                if (mode.mode.equals("car")) {
+                    commuteModeCarCntMap.merge(dwelling.getZoneId(), 1, Integer::sum);
+                    commuteModeCarUtilSumMap.merge(dwelling.getZoneId(), mode.utility, Double::sum);
+                } else if (mode.mode.equals("pt")) {
+                    commuteModePtCntMap.merge(dwelling.getZoneId(), 1, Integer::sum);
+                    commuteModePtUtilSumMap.merge(dwelling.getZoneId(), mode.utility, Double::sum);
+                } else {
+                    throw new RuntimeException();
+                }
 
             }
         }
+
+        cntMap.merge(dwelling.getZoneId(), 1, Integer::sum);
+        ddQualityUtilityMap.merge(dwelling.getZoneId(), ddQualityUtility, Double::sum);
+        ddSizeUtilityMap.merge(dwelling.getZoneId(), ddSizeUtility, Double::sum);
+        ddAutoAccessibilityUtilityMap.merge(dwelling.getZoneId(), ddAutoAccessibilityUtility, Double::sum);
+        transitAccessibilityUtilityMap.merge(dwelling.getZoneId(), transitAccessibilityUtility, Double::sum);
+        ddPriceUtilityMap.merge(dwelling.getZoneId(), ddPriceUtility, Double::sum);
+
+
+        workDistanceUtilityMap.merge(dwelling.getZoneId(), workDistanceUtility, Double::sum);
+
+
         return dwellingUtilityStrategy.calculateSelectDwellingUtility(ht, ddSizeUtility, ddPriceUtility,
                 ddQualityUtility, ddAutoAccessibilityUtility,
                 transitAccessibilityUtility, workDistanceUtility);
diff --git a/siloCore/src/main/java/de/tum/bgu/msm/properties/modules/TransportModelPropertiesModule.java b/siloCore/src/main/java/de/tum/bgu/msm/properties/modules/TransportModelPropertiesModule.java
index 71d1ac1e9..11408df37 100644
--- a/siloCore/src/main/java/de/tum/bgu/msm/properties/modules/TransportModelPropertiesModule.java
+++ b/siloCore/src/main/java/de/tum/bgu/msm/properties/modules/TransportModelPropertiesModule.java
@@ -24,6 +24,11 @@ public class TransportModelPropertiesModule {
      */
     public final TransportModelIdentifier transportModelIdentifier;
 
+    /**
+     * Identifier for start type in MATSim {@link MatsimStartOptionIdentifier}: COLD, HOT.
+     */
+    public final MatsimStartOptionIdentifier matsimStartOption;
+
     /**
      * Identifier for the transport model {@link TransportModelIdentifier}: MITO_MATSIM, MATSIM, NONE.
      */
@@ -53,6 +58,10 @@ public enum TransportModelIdentifier {
         MITO_MATSIM, MATSIM, NONE;
     }
 
+    public enum MatsimStartOptionIdentifier {
+        COLD, HOT;
+    }
+
     /**
      * Identifier for which data structure to use for travel times.
      */
@@ -76,6 +85,9 @@ public TransportModelPropertiesModule(ResourceBundle bundle) {
         PropertiesUtil.newPropertySubmodule("Transport model identifier (MITO_MATSIM, MATSIM, NONE, or empty)");
         transportModelIdentifier = TransportModelIdentifier.valueOf(PropertiesUtil.getStringProperty(bundle, "transport.model", "NONE").toUpperCase());
 
+        PropertiesUtil.newPropertySubmodule("Transport model identifier (MITO_MATSIM, MATSIM, NONE, or empty)");
+        matsimStartOption = MatsimStartOptionIdentifier.valueOf(PropertiesUtil.getStringProperty(bundle, "matsim.start", "COLD").toUpperCase());
+
         PropertiesUtil.newPropertySubmodule("Travel time data structure identifier (MATSIM or SKIM)");
         travelTimeImplIdentifier = TravelTimeImplIdentifier.valueOf(PropertiesUtil.getStringProperty(bundle, "travel.time", "SKIM").toUpperCase());
 
diff --git a/siloCore/src/main/java/de/tum/bgu/msm/simulator/Simulator.java b/siloCore/src/main/java/de/tum/bgu/msm/simulator/Simulator.java
index f69b8e01e..3d430a9bc 100644
--- a/siloCore/src/main/java/de/tum/bgu/msm/simulator/Simulator.java
+++ b/siloCore/src/main/java/de/tum/bgu/msm/simulator/Simulator.java
@@ -8,11 +8,14 @@
 import de.tum.bgu.msm.models.AbstractModel;
 import de.tum.bgu.msm.models.EventModel;
 import de.tum.bgu.msm.models.ModelUpdateListener;
+import de.tum.bgu.msm.models.relocation.moves.MovesModelImpl;
 import de.tum.bgu.msm.utils.SiloUtil;
 import de.tum.bgu.msm.utils.TimeTracker;
 import org.apache.logging.log4j.LogManager;
 import org.apache.logging.log4j.Logger;
 
+import java.io.PrintWriter;
+import java.nio.file.Path;
 import java.util.*;
 
 /**
@@ -31,7 +34,6 @@ public final class Simulator {
 
     private final List<MicroEvent> events = new ArrayList<>();
     private final TimeTracker timeTracker;
-
     private final Set<ResultsMonitor> resultsMonitors = new HashSet<>() ;
 
     public Simulator(TimeTracker timeTracker) {
@@ -144,5 +146,7 @@ public void endSimulation() {
             resultsMonitor.endSimulation();
         }
 
+
+
     }
 }
diff --git a/useCases/fabiland/scenario/25r_ae.properties b/useCases/fabiland/scenario/25r_ae.properties
index 932e6c1cd..e23561747 100644
--- a/useCases/fabiland/scenario/25r_ae.properties
+++ b/useCases/fabiland/scenario/25r_ae.properties
@@ -6,6 +6,7 @@ zonal.data.file = input/zoneSystem_25-reg.csv
 employment.types = IND
 transport.model = MATSIM
 travel.time = matsim
+matsim.start = hot
 # transport.model.years = 5,10
 transport.model.years = 1,2,3,4,5,6,7,8,9,10
 matsim.simulate.car.trips.only = false
diff --git a/useCases/fabiland/scenario/config_cap5_1-l_nes_smc.xml b/useCases/fabiland/scenario/config_cap5_1-l_nes_smc.xml
new file mode 100644
index 000000000..ff555f73e
--- /dev/null
+++ b/useCases/fabiland/scenario/config_cap5_1-l_nes_smc.xml
@@ -0,0 +1,65 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<!DOCTYPE config SYSTEM "http://www.matsim.org/files/dtd/config_v2.dtd">
+<config>
+	<module name="global" >
+		<param name="numberOfThreads" value="1" />
+	</module>
+	<module name="controler" >
+		<param name="lastIteration" value="5" />
+		<param name="overwriteFiles" value="deleteDirectoryIfExists" />
+	</module>
+	<module name="network" >
+		<param name="inputNetworkFile" value="matsimInput/nw_cap30_1-l_nes.xml" />
+	</module>
+	<module name="transit">
+		<param name="transitScheduleFile" value="matsimInput/ts_1-l_nes.xml" />
+		<param name="useTransit" value="true" />
+	</module>
+	<module name="qsim" >
+		<param name="trafficDynamics" value="withHoles" />
+		<!-- Note that this setting is independent from the "global" threads setting -->
+		<param name="numberOfThreads" value="1" />
+	</module>
+	<module name="parallelEventHandling" >
+		<param name="numberOfThreads" value="1" />
+	</module>
+	<module name="planCalcScore" >
+		<parameterset type="activityParams" >
+			<param name="activityType" value="home" />
+			<param name="typicalDuration" value="12:00:00" />
+		</parameterset>
+		<parameterset type="activityParams" >
+			<param name="activityType" value="work" />
+			<param name="typicalDuration" value="08:00:00" />
+		</parameterset>
+	</module>
+	<module name="strategy" >
+		<param name="fractionOfIterationsToDisableInnovation" value="0.8" />
+		<parameterset type="strategysettings" >
+			<param name="strategyName" value="ChangeExpBeta" />
+			<param name="weight" value="0.8" />
+		</parameterset>
+		<parameterset type="strategysettings" >
+			<param name="strategyName" value="ReRoute" />
+			<param name="weight" value="0.1" />
+		</parameterset>
+		<parameterset type="strategysettings" >
+			<param name="strategyName" value="SubtourModeChoice" />
+			<param name="weight" value="0.1" />
+		</parameterset>
+	</module>
+	<module name="subtourModeChoice" >
+		<param name="considerCarAvailability" value="true" />
+		<param name="modes" value="car,pt" />
+	</module>
+	<module name="transit" >
+		<param name="usingTransitInMobsim" value="false" />
+	</module>
+	<module name="transitRouter" >
+		<param name="extensionRadius" value="500.0" />
+	</module>
+	<module name="vspExperimental" >
+		<!-- Options: | ignore | info | warn | abort | .  When violating VSP defaults, this results in nothing, logfile infos, logfile warnings, or aborts.  Members of VSP should use `abort' or talk to kai. -->
+		<param name="vspDefaultsCheckingLevel" value="warn" />
+	</module>
+</config>
\ No newline at end of file
diff --git a/useCases/fabiland/scenario/matsimInput/_nw_jr_tmp.xml b/useCases/fabiland/scenario/matsimInput/_nw_jr_tmp.xml
new file mode 100644
index 000000000..a0316726e
--- /dev/null
+++ b/useCases/fabiland/scenario/matsimInput/_nw_jr_tmp.xml
@@ -0,0 +1,245 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<!DOCTYPE network SYSTEM "http://www.matsim.org/files/dtd/network_v2.dtd">
+<network>
+
+<!-- ====================================================================== -->
+
+	<nodes>
+		<node id="1" x="-10000.0" y="10000.0" >
+		</node>
+		<node id="10" x="10000.0" y="5000.0" >
+		</node>
+		<node id="11" x="-10000.0" y="0.0" >
+		</node>
+		<node id="12" x="-5000.0" y="0.0" >
+		</node>
+		<node id="13" x="0.0" y="0.0" >
+		</node>
+		<node id="14" x="5000.0" y="0.0" >
+		</node>
+		<node id="15" x="10000.0" y="0.0" >
+		</node>
+		<node id="16" x="-10000.0" y="-5000.0" >
+		</node>
+		<node id="17" x="-5000.0" y="-5000.0" >
+		</node>
+		<node id="18" x="0.0" y="-5000.0" >
+		</node>
+		<node id="19" x="5000.0" y="-5000.0" >
+		</node>
+		<node id="2" x="-5000.0" y="10000.0" >
+		</node>
+		<node id="20" x="10000.0" y="-5000.0" >
+		</node>
+		<node id="21" x="-10000.0" y="-10000.0" >
+		</node>
+		<node id="22" x="-5000.0" y="-10000.0" >
+		</node>
+		<node id="23" x="0.0" y="-10000.0" >
+		</node>
+		<node id="24" x="5000.0" y="-10000.0" >
+		</node>
+		<node id="25" x="10000.0" y="-10000.0" >
+		</node>
+		<node id="3" x="0.0" y="10000.0" >
+		</node>
+		<node id="4" x="5000.0" y="10000.0" >
+		</node>
+		<node id="5" x="10000.0" y="10000.0" >
+		</node>
+		<node id="6" x="-10000.0" y="5000.0" >
+		</node>
+		<node id="7" x="-5000.0" y="5000.0" >
+		</node>
+		<node id="8" x="0.0" y="5000.0" >
+		</node>
+		<node id="9" x="5000.0" y="5000.0" >
+		</node>
+		<node id="pt-13" x="0.0" y="0.0" >
+		</node>
+		<node id="pt-5" x="10000.0" y="10000.0" >
+		</node>
+		<node id="pt-9" x="5000.0" y="5000.0" >
+		</node>
+	</nodes>
+
+<!-- ====================================================================== -->
+
+	<links capperiod="01:00:00" effectivecellsize="7.5" effectivelanewidth="3.75">
+		<link id="0102" from="1" to="2" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0106" from="1" to="6" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0201" from="2" to="1" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0203" from="2" to="3" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0207" from="2" to="7" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0302" from="3" to="2" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0304" from="3" to="4" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0308" from="3" to="8" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0403" from="4" to="3" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0405" from="4" to="5" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0409" from="4" to="9" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0504" from="5" to="4" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0510" from="5" to="10" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0601" from="6" to="1" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0607" from="6" to="7" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0611" from="6" to="11" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0702" from="7" to="2" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0706" from="7" to="6" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0708" from="7" to="8" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0712" from="7" to="12" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0803" from="8" to="3" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0807" from="8" to="7" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0809" from="8" to="9" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0813" from="8" to="13" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0904" from="9" to="4" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0908" from="9" to="8" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0910" from="9" to="10" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0914" from="9" to="14" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1005" from="10" to="5" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1009" from="10" to="9" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1015" from="10" to="15" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1106" from="11" to="6" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1112" from="11" to="12" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1116" from="11" to="16" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1207" from="12" to="7" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1211" from="12" to="11" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1213" from="12" to="13" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1217" from="12" to="17" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1308" from="13" to="8" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1312" from="13" to="12" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1314" from="13" to="14" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1318" from="13" to="18" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1409" from="14" to="9" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1413" from="14" to="13" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1415" from="14" to="15" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1419" from="14" to="19" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1510" from="15" to="10" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1514" from="15" to="14" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1520" from="15" to="20" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1611" from="16" to="11" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1617" from="16" to="17" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1621" from="16" to="21" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1712" from="17" to="12" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1716" from="17" to="16" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1718" from="17" to="18" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1722" from="17" to="22" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1813" from="18" to="13" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1817" from="18" to="17" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1819" from="18" to="19" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1823" from="18" to="23" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1914" from="19" to="14" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1918" from="19" to="18" length="5000.0" freespeed="0.1388888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1920" from="19" to="20" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1924" from="19" to="24" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2015" from="20" to="15" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2019" from="20" to="19" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2025" from="20" to="25" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2116" from="21" to="16" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2122" from="21" to="22" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2217" from="22" to="17" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2221" from="22" to="21" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2223" from="22" to="23" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2318" from="23" to="18" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2322" from="23" to="22" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2324" from="23" to="24" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2419" from="24" to="19" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2423" from="24" to="23" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2425" from="24" to="25" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2520" from="25" to="20" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2524" from="25" to="24" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="pt-0" from="pt-5" to="pt-5" length="50.0" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-1" from="pt-5" to="pt-9" length="7071.067811865475" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-2" from="pt-9" to="pt-13" length="7071.067811865475" freespeed="0.08333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-3" from="pt-13" to="pt-13" length="50.0" freespeed="0.08333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-4" from="pt-13" to="pt-9" length="7071.067811865475" freespeed="0.08333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-5" from="pt-9" to="pt-5" length="7071.067811865475" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+	</links>
+
+<!-- ====================================================================== -->
+
+</network>
diff --git a/useCases/fabiland/scenario/matsimInput/autofrei.xml b/useCases/fabiland/scenario/matsimInput/autofrei.xml
new file mode 100644
index 000000000..8602d28f9
--- /dev/null
+++ b/useCases/fabiland/scenario/matsimInput/autofrei.xml
@@ -0,0 +1,245 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<!DOCTYPE network SYSTEM "http://www.matsim.org/files/dtd/network_v2.dtd">
+<network>
+
+<!-- ====================================================================== -->
+
+	<nodes>
+		<node id="1" x="-10000.0" y="10000.0" >
+		</node>
+		<node id="10" x="10000.0" y="5000.0" >
+		</node>
+		<node id="11" x="-10000.0" y="0.0" >
+		</node>
+		<node id="12" x="-5000.0" y="0.0" >
+		</node>
+		<node id="13" x="0.0" y="0.0" >
+		</node>
+		<node id="14" x="5000.0" y="0.0" >
+		</node>
+		<node id="15" x="10000.0" y="0.0" >
+		</node>
+		<node id="16" x="-10000.0" y="-5000.0" >
+		</node>
+		<node id="17" x="-5000.0" y="-5000.0" >
+		</node>
+		<node id="18" x="0.0" y="-5000.0" >
+		</node>
+		<node id="19" x="5000.0" y="-5000.0" >
+		</node>
+		<node id="2" x="-5000.0" y="10000.0" >
+		</node>
+		<node id="20" x="10000.0" y="-5000.0" >
+		</node>
+		<node id="21" x="-10000.0" y="-10000.0" >
+		</node>
+		<node id="22" x="-5000.0" y="-10000.0" >
+		</node>
+		<node id="23" x="0.0" y="-10000.0" >
+		</node>
+		<node id="24" x="5000.0" y="-10000.0" >
+		</node>
+		<node id="25" x="10000.0" y="-10000.0" >
+		</node>
+		<node id="3" x="0.0" y="10000.0" >
+		</node>
+		<node id="4" x="5000.0" y="10000.0" >
+		</node>
+		<node id="5" x="10000.0" y="10000.0" >
+		</node>
+		<node id="6" x="-10000.0" y="5000.0" >
+		</node>
+		<node id="7" x="-5000.0" y="5000.0" >
+		</node>
+		<node id="8" x="0.0" y="5000.0" >
+		</node>
+		<node id="9" x="5000.0" y="5000.0" >
+		</node>
+		<node id="pt-13" x="0.0" y="0.0" >
+		</node>
+		<node id="pt-5" x="10000.0" y="10000.0" >
+		</node>
+		<node id="pt-9" x="5000.0" y="5000.0" >
+		</node>
+	</nodes>
+
+<!-- ====================================================================== -->
+
+	<links capperiod="01:00:00" effectivecellsize="7.5" effectivelanewidth="3.75">
+		<link id="0102" from="1" to="2" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0106" from="1" to="6" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0201" from="2" to="1" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0203" from="2" to="3" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0207" from="2" to="7" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0302" from="3" to="2" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0304" from="3" to="4" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0308" from="3" to="8" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0403" from="4" to="3" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0405" from="4" to="5" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0409" from="4" to="9" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0504" from="5" to="4" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0510" from="5" to="10" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0601" from="6" to="1" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0607" from="6" to="7" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0611" from="6" to="11" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0702" from="7" to="2" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0706" from="7" to="6" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0708" from="7" to="8" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0712" from="7" to="12" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0803" from="8" to="3" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0807" from="8" to="7" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0809" from="8" to="9" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0813" from="8" to="13" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0904" from="9" to="4" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0908" from="9" to="8" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0910" from="9" to="10" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="0914" from="9" to="14" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1005" from="10" to="5" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1009" from="10" to="9" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1015" from="10" to="15" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1106" from="11" to="6" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1112" from="11" to="12" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1116" from="11" to="16" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1207" from="12" to="7" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1211" from="12" to="11" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1213" from="12" to="13" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1217" from="12" to="17" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1308" from="13" to="8" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1312" from="13" to="12" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1314" from="13" to="14" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1318" from="13" to="18" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1409" from="14" to="9" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1413" from="14" to="13" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1415" from="14" to="15" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1419" from="14" to="19" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1510" from="15" to="10" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1514" from="15" to="14" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1520" from="15" to="20" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1611" from="16" to="11" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1617" from="16" to="17" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1621" from="16" to="21" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1712" from="17" to="12" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1716" from="17" to="16" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1718" from="17" to="18" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1722" from="17" to="22" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1813" from="18" to="13" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1817" from="18" to="17" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1819" from="18" to="19" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1823" from="18" to="23" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1914" from="19" to="14" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1918" from="19" to="18" length="5000.0" freespeed="0.1" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1920" from="19" to="20" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="1924" from="19" to="24" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2015" from="20" to="15" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2019" from="20" to="19" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2025" from="20" to="25" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2116" from="21" to="16" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2122" from="21" to="22" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2217" from="22" to="17" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2221" from="22" to="21" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2223" from="22" to="23" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2318" from="23" to="18" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2322" from="23" to="22" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2324" from="23" to="24" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2419" from="24" to="19" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2423" from="24" to="23" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2425" from="24" to="25" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2520" from="25" to="20" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="2524" from="25" to="24" length="5000.0" freespeed="13.88888888888889" capacity="30.0" permlanes="1.0" oneway="1" modes="car" >
+		</link>
+		<link id="pt-0" from="pt-5" to="pt-5" length="50.0" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-1" from="pt-5" to="pt-9" length="7071.067811865475" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-2" from="pt-9" to="pt-13" length="7071.067811865475" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-3" from="pt-13" to="pt-13" length="50.0" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-4" from="pt-13" to="pt-9" length="7071.067811865475" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+		<link id="pt-5" from="pt-9" to="pt-5" length="7071.067811865475" freespeed="8.333333333333334" capacity="500.0" permlanes="1.0" oneway="1" modes="pt" >
+		</link>
+	</links>
+
+<!-- ====================================================================== -->
+
+</network>
diff --git a/useCases/fabiland/scenario/test2.properties b/useCases/fabiland/scenario/test2.properties
new file mode 100644
index 000000000..32ad3e754
--- /dev/null
+++ b/useCases/fabiland/scenario/test2.properties
@@ -0,0 +1,14 @@
+scenario.name = test
+base.year = 0
+end.year = 2
+random.seed = 1
+zones.shapefile = input/zonesShapefile/fabiland.shp
+zonal.data.file = input/zoneSystem_1-reg.csv
+employment.types = IND
+transport.model = MATSIM
+travel.time = matsim
+#transport.model.years=1
+#transport.model.years = 5,10
+transport.model.years = 1,2
+matsim.simulate.car.trips.only = false
+matsim.include.access.egress = true
diff --git a/useCases/fabiland/scripts/AutoFreiTests.ipynb b/useCases/fabiland/scripts/AutoFreiTests.ipynb
new file mode 100644
index 000000000..00ac23f90
--- /dev/null
+++ b/useCases/fabiland/scripts/AutoFreiTests.ipynb
@@ -0,0 +1,311 @@
+{
+ "cells": [
+  {
+   "cell_type": "raw",
+   "id": "d4f71dbc45fd0850",
+   "metadata": {},
+   "source": []
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 1,
+   "id": "fc795b2fb2a17b82",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "import pandas as pd"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 58,
+   "id": "74720f06",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "base = pd.read_csv(\"/home/brendan/git/silo/useCases/fabiland/scenario/scenOutput/base/avgSatisfactionByZone_base.csv\")\n",
+    "policy = pd.read_csv(\"/home/brendan/git/silo/useCases/fabiland/scenario/scenOutput/base/avgSatisfactionByZone_policy.csv\")\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 50,
+   "id": "dc3521fc",
+   "metadata": {},
+   "outputs": [
+    {
+     "data": {
+      "text/html": [
+       "<div>\n",
+       "<style scoped>\n",
+       "    .dataframe tbody tr th:only-of-type {\n",
+       "        vertical-align: middle;\n",
+       "    }\n",
+       "\n",
+       "    .dataframe tbody tr th {\n",
+       "        vertical-align: top;\n",
+       "    }\n",
+       "\n",
+       "    .dataframe thead th {\n",
+       "        text-align: right;\n",
+       "    }\n",
+       "</style>\n",
+       "<table border=\"1\" class=\"dataframe\">\n",
+       "  <thead>\n",
+       "    <tr style=\"text-align: right;\">\n",
+       "      <th></th>\n",
+       "      <th>zone</th>\n",
+       "      <th>avgSatisfaction</th>\n",
+       "      <th>avgProp</th>\n",
+       "      <th>year</th>\n",
+       "    </tr>\n",
+       "  </thead>\n",
+       "  <tbody>\n",
+       "    <tr>\n",
+       "      <th>0</th>\n",
+       "      <td>1</td>\n",
+       "      <td>0.12280260233844255</td>\n",
+       "      <td>0.08719006839701185</td>\n",
+       "      <td>1</td>\n",
+       "    </tr>\n",
+       "    <tr>\n",
+       "      <th>1</th>\n",
+       "      <td>2</td>\n",
+       "      <td>0.08515294864097894</td>\n",
+       "      <td>0.14104434132920057</td>\n",
+       "      <td>1</td>\n",
+       "    </tr>\n",
+       "    <tr>\n",
+       "      <th>2</th>\n",
+       "      <td>3</td>\n",
+       "      <td>0.2650501217684014</td>\n",
+       "      <td>0.0913587355063226</td>\n",
+       "      <td>1</td>\n",
+       "    </tr>\n",
+       "    <tr>\n",
+       "      <th>3</th>\n",
+       "      <td>4</td>\n",
+       "      <td>0.2448101344060386</td>\n",
+       "      <td>0.07201396410230788</td>\n",
+       "      <td>1</td>\n",
+       "    </tr>\n",
+       "    <tr>\n",
+       "      <th>4</th>\n",
+       "      <td>5</td>\n",
+       "      <td>0.09848246455445277</td>\n",
+       "      <td>0.09321934310869487</td>\n",
+       "      <td>1</td>\n",
+       "    </tr>\n",
+       "  </tbody>\n",
+       "</table>\n",
+       "</div>"
+      ],
+      "text/plain": [
+       "  zone      avgSatisfaction              avgProp year\n",
+       "0    1  0.12280260233844255  0.08719006839701185    1\n",
+       "1    2  0.08515294864097894  0.14104434132920057    1\n",
+       "2    3   0.2650501217684014   0.0913587355063226    1\n",
+       "3    4   0.2448101344060386  0.07201396410230788    1\n",
+       "4    5  0.09848246455445277  0.09321934310869487    1"
+      ]
+     },
+     "execution_count": 50,
+     "metadata": {},
+     "output_type": "execute_result"
+    }
+   ],
+   "source": [
+    "policy.head()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 78,
+   "id": "9f56517f",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "zone 13 base Satisfaction:  0.19873045559881253\n",
+      "zone 13 base Prop to move:  0.08308358488306133\n"
+     ]
+    },
+    {
+     "name": "stderr",
+     "output_type": "stream",
+     "text": [
+      "/tmp/ipykernel_88934/2214987494.py:2: SettingWithCopyWarning: \n",
+      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
+      "Try using .loc[row_indexer,col_indexer] = value instead\n",
+      "\n",
+      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
+      "  base_13_S[\"avgSatisfaction\"] = pd.to_numeric(base_13_S[\"avgSatisfaction\"], errors=\"coerce\")\n",
+      "/tmp/ipykernel_88934/2214987494.py:6: SettingWithCopyWarning: \n",
+      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
+      "Try using .loc[row_indexer,col_indexer] = value instead\n",
+      "\n",
+      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
+      "  base_13_P[\"avgProp\"] = pd.to_numeric(base_13_P[\"avgProp\"], errors=\"coerce\")\n"
+     ]
+    }
+   ],
+   "source": [
+    "base_13_S = base[base[\"zone\"] == \"13\"]\n",
+    "base_13_S[\"avgSatisfaction\"] = pd.to_numeric(base_13_S[\"avgSatisfaction\"], errors=\"coerce\")\n",
+    "# base_13_S[\"avgSatisfaction\"].mean()\n",
+    "\n",
+    "base_13_P = base[base[\"zone\"] == \"13\"]\n",
+    "base_13_P[\"avgProp\"] = pd.to_numeric(base_13_P[\"avgProp\"], errors=\"coerce\")\n",
+    "# base_13_P[\"avgProp\"].mean()\n",
+    "\n",
+    "print(\"zone 13 base Satisfaction: \", base_13_S[\"avgSatisfaction\"].mean())\n",
+    "print(\"zone 13 base Prop to move: \", base_13_S[\"avgProp\"].mean())"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 74,
+   "id": "e64f5f99",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "base[\"avgSatisfaction\"] = pd.to_numeric(base[\"avgSatisfaction\"], errors=\"coerce\")\n",
+    "base_S_allzones = base[\"avgSatisfaction\"]\n",
+    "# base_S_allzones.mean()\n",
+    "\n",
+    "\n",
+    "base[\"avgProp\"] = pd.to_numeric(base[\"avgProp\"], errors=\"coerce\")\n",
+    "base_P_allzones = base[\"avgProp\"]\n",
+    "# base_P_allzones.mean()\n",
+    "\n",
+    "# base_allzones[base_allzones[\"avgSatisfaction\"]].mean()"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 75,
+   "id": "34c7f93f",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "0.9172196139377837\n",
+      "0.8135163248642484\n"
+     ]
+    }
+   ],
+   "source": [
+    "print(base_13_S[\"avgSatisfaction\"].mean() / base_S_allzones.mean())\n",
+    "\n",
+    "print(base_13_P[\"avgProp\"].mean() / base_P_allzones.mean())"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 81,
+   "id": "d15ccb0c",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "zone 13 policy  Satisfaction:  0.1570603020779263\n",
+      "zone 13 policy Prop to move:  0.08525122821083374\n"
+     ]
+    },
+    {
+     "name": "stderr",
+     "output_type": "stream",
+     "text": [
+      "/tmp/ipykernel_88934/2745029612.py:2: SettingWithCopyWarning: \n",
+      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
+      "Try using .loc[row_indexer,col_indexer] = value instead\n",
+      "\n",
+      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
+      "  policy_13_S[\"avgSatisfaction\"] = pd.to_numeric(policy_13_S[\"avgSatisfaction\"], errors=\"coerce\")\n",
+      "/tmp/ipykernel_88934/2745029612.py:6: SettingWithCopyWarning: \n",
+      "A value is trying to be set on a copy of a slice from a DataFrame.\n",
+      "Try using .loc[row_indexer,col_indexer] = value instead\n",
+      "\n",
+      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
+      "  policy_13_P[\"avgProp\"] = pd.to_numeric(policy_13_P[\"avgProp\"], errors=\"coerce\")\n"
+     ]
+    }
+   ],
+   "source": [
+    "policy_13_S = policy[policy[\"zone\"] == \"13\"]\n",
+    "policy_13_S[\"avgSatisfaction\"] = pd.to_numeric(policy_13_S[\"avgSatisfaction\"], errors=\"coerce\")\n",
+    "# policy_13_S[\"avgSatisfaction\"].mean()\n",
+    "\n",
+    "policy_13_P = policy[policy[\"zone\"] == \"13\"]\n",
+    "policy_13_P[\"avgProp\"] = pd.to_numeric(policy_13_P[\"avgProp\"], errors=\"coerce\")\n",
+    "# policy_13_P[\"avgProp\"].mean()\n",
+    "\n",
+    "print(\"zone 13 policy  Satisfaction: \", policy_13_S[\"avgSatisfaction\"].mean())\n",
+    "print(\"zone 13 policy Prop to move: \", policy_13_S[\"avgProp\"].mean())"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 71,
+   "id": "cd92afd6",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "policy[\"avgSatisfaction\"] = pd.to_numeric(policy[\"avgSatisfaction\"], errors=\"coerce\")\n",
+    "policy_S_allzones = policy[\"avgSatisfaction\"]\n",
+    "# policy_S_allzones.mean()\n",
+    "\n",
+    "policy[\"avgProp\"] = pd.to_numeric(policy[\"avgProp\"], errors=\"coerce\")\n",
+    "policy_P_allzones = policy[\"avgProp\"]\n",
+    "# policy_P_allzones.mean()\n"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": 72,
+   "id": "3f73c7ba",
+   "metadata": {},
+   "outputs": [
+    {
+     "name": "stdout",
+     "output_type": "stream",
+     "text": [
+      "0.7620574214414428\n",
+      "0.892060397948594\n"
+     ]
+    }
+   ],
+   "source": [
+    "print(policy_13_S[\"avgSatisfaction\"].mean() / policy_S_allzones.mean())\n",
+    "print(policy_13_P[\"avgProp\"].mean() / policy_P_allzones.mean())"
+   ]
+  }
+ ],
+ "metadata": {
+  "kernelspec": {
+   "display_name": "3.10.13",
+   "language": "python",
+   "name": "python3"
+  },
+  "language_info": {
+   "codemirror_mode": {
+    "name": "ipython",
+    "version": 3
+   },
+   "file_extension": ".py",
+   "mimetype": "text/x-python",
+   "name": "python",
+   "nbconvert_exporter": "python",
+   "pygments_lexer": "ipython3",
+   "version": "3.10.13"
+  }
+ },
+ "nbformat": 4,
+ "nbformat_minor": 5
+}
diff --git a/useCases/fabiland/scripts/AutofreiAnalysis.Rmd b/useCases/fabiland/scripts/AutofreiAnalysis.Rmd
new file mode 100644
index 000000000..0cb14456f
--- /dev/null
+++ b/useCases/fabiland/scripts/AutofreiAnalysis.Rmd
@@ -0,0 +1,666 @@
+---
+title: "xxx"
+output: html_document
+---
+
+
+```{r}
+library(tidyverse)
+# library(sf)
+# library(tmap)
+```
+
+# FABILAND
+
+## Setup
+networks for ggplots
+```{r}
+net_total <-  list(annotate("segment", x = -10000, xend =10000,
+           y = c(-10000,-5000, 0, 5000,10000),
+           yend = c(-10000,-5000, 0, 5000,10000)),
+  annotate("segment", y = -10000, yend =10000,
+           x = c(-10000,-5000, 0, 5000,10000),
+           xend = c(-10000,-5000, 0, 5000,10000)) )
+
+
+net_autofrei_large <-  list(annotate("segment", x = -5000, xend =5000,
+           y = c(-5000, 0, 5000),
+           yend = c(-5000, 0, 5000), col = "red"),
+  annotate("segment", y = -5000, yend =5000,
+           x = c(-5000, 0, 5000),
+           xend = c(-5000, 0, 5000),col = "red") )
+
+net_autofrei_small <-  list(annotate("segment", x = -5000, xend =5000,
+           y = 0, col = "red"),
+  annotate("segment", y = -5000, yend =5000, x = 0,col = "red") )
+
+
+pt_color <- "#D8BFD8"
+
+
+net_pt_nes <- list(annotate("segment", x = 0, xend = 10000, y = 0, yend = 10000, col = pt_color,lwd = 2))
+
+
+net_pt_ringsx <- list(annotate("segment", x = -10000, xend = 10000, y = -10000, yend = 10000, col = pt_color, lwd = 2),
+                      annotate("segment", x = -10000, xend = 10000, y = 10000, yend = -10000, col = pt_color, lwd = 2),
+                      annotate("segment", x = c(-10000,10000), y = -10000, yend = 10000, col = pt_color, lwd = 2),
+                      annotate("segment", y = c(-10000,10000), x = -10000, xend = 10000, col = pt_color, lwd = 2),
+                      annotate("segment", x = c(-5000,5000), y = -5000, yend = 5000, col = pt_color, lwd = 2),
+                      annotate("segment", y = c(-5000,5000), x = -5000, xend = 5000, col = pt_color, lwd = 2))
+
+
+
+
+```
+
+
+Define Directory for Scenario
+```{r}
+
+setwd("/home/brendan/svn/")
+root_dir <- "shared-svn/projects/fabilut/fabiland/"
+
+# 1:
+# Carfree: Small
+# PT: NES
+
+# net_autofrei <- net_autofrei_small
+# net_pt <- net_pt_nes
+
+base_dir <- paste0(root_dir,"2026-09-15b__autofrei-BASE__pt-NES/")
+policy_dir <- paste0(root_dir,"2026-09-15b__autofrei-SMALL__pt-NES/")
+
+# 2:
+# Carfree: Large
+# PT: NES
+
+# net_autofrei <- net_autofrei_large
+# net_pt <- net_pt_nes
+# 
+# base_dir <- paste0(root_dir,"2026-09-15b__autofrei-BASE__pt-NES/")
+# policy_dir <- paste0(root_dir,"2026-09-15b__autofrei-LARGE__pt-NES/")
+
+
+# 3:
+# Carfree: Small
+# PT: RINGSX
+# 
+# net_autofrei <- net_autofrei_small
+# net_pt <- net_pt_ringsx
+# 
+# base_dir <- paste0(root_dir,"2026-09-15b__autofrei-BASE__pt-RINGSX/")
+# policy_dir <- paste0(root_dir,"2026-09-15b__autofrei-SMALL__pt-RINGSX/")
+
+# 4:
+# Carfree: Large
+# PT: RINGSX
+
+# net_autofrei <- net_autofrei_large
+# net_pt <- net_pt_ringsx
+# 
+# base_dir <- paste0(root_dir,"2026-09-15b__autofrei-BASE__pt-RINGSX/")
+# policy_dir <- paste0(root_dir,"2026-09-15b__autofrei-LARGE__pt-RINGSX/")
+
+
+
+
+###
+year <- 10
+base_pp <- read_csv(paste0(base_dir,"microData/pp_",year,".csv"))
+base_hh <- read_csv(paste0(base_dir,"microData/hh_",year,".csv"))
+base_dd <- read_csv(paste0(base_dir,"microData/dd_",year,".csv"))
+base_jj <- read_csv(paste0(base_dir,"microData/jj_",year,".csv"))
+
+
+policy_pp <- read_csv(paste0(policy_dir,"microData/pp_",year,".csv"))
+policy_hh <- read_csv(paste0(policy_dir,"microData/hh_",year,".csv"))
+policy_dd <- read_csv(paste0(policy_dir,"microData/dd_",year,".csv"))
+policy_jj <- read_csv(paste0(policy_dir,"microData/jj_",year,".csv"))
+
+
+
+base_merged <- base_pp |> left_join(base_hh, by = c("hhid" = "id")) |>  left_join(base_dd, by = c("dwelling"="id")) |> left_join(base_jj |> select(id = personId, jobX = coordX, jobY = coordY)) |> mutate(coordXrounded = round(coordX/5000) *5000, coordYrounded = round(coordY/5000) *5000) 
+
+
+
+policy_merged <- policy_pp |> left_join(policy_hh, by = c("hhid" = "id")) |> left_join(policy_dd, by = c("dwelling"="id")) |> 
+  left_join(policy_jj |> select(id = personId, jobX = coordX, jobY = coordY))  |> mutate(coordXrounded = round(coordX/5000) *5000, coordYrounded = round(coordY/5000) *5000) 
+
+
+# get agents who left zone 13 from base to policy
+policy_from13 <- policy_merged |>
+  filter(!(coordXrounded == 0 & coordYrounded == 0)) |>
+  filter(relationShip != "child") |>
+  inner_join(
+    base_merged |>
+      filter(coordXrounded == 0, coordYrounded == 0) |>
+      filter(relationShip != "child"),
+    by = "id",
+    suffix = c("_policy", "_base")
+  )
+
+policy_to13 <- policy_merged |>
+  filter(coordXrounded == 0, coordYrounded == 0) |>
+  filter(relationShip != "child") |>
+  inner_join(
+    base_merged |>
+      filter(!(coordXrounded == 0 & coordYrounded == 0)) |>
+      filter(relationShip != "child"),
+    by = "id",
+    suffix = c("_policy", "_base")
+  )
+
+mean(policy_from13$income_policy)
+mean(policy_to13$income_policy)
+
+employment_from13 <- table(policy_from13$occupation_policy) |> 
+  setNames(c("employed", "unemployed"))
+rate <- employment_from13["employed"] / (employment_from13["employed"] + employment_from13["unemployed"])
+rate
+
+employment_to13 <- table(policy_to13$occupation_policy) |> 
+  setNames(c("employed", "unemployed"))
+rate <- employment_to13["employed"] / (employment_to13["employed"] + employment_to13["unemployed"])
+rate
+
+# print("Average Income of Adults")
+# print(paste0("base: ", base_merged |> filter(coordXrounded == 0 & coordYrounded == 0) |> filter(relationShip!="child") |> summarise(income = mean(income))))
+# print(paste0("policy: ", policy_merged |> filter(coordXrounded == 0 & coordYrounded == 0) |> filter(relationShip!="child") |> summarise(income = mean(income))))
+# 
+# zonal_counts <- base_merged |>
+#   select(id, coordXrounded, coordYrounded) |> 
+#   mutate(scenario = "base") |> 
+#   bind_rows(policy_merged |> select(id,coordXrounded,coordYrounded) |> mutate(scenario = "policy")) |>
+#   group_by(coordXrounded,coordYrounded,scenario) |> count()
+# 
+# zonal_mean_incomes  <- base_merged |>
+#   select(id, coordXrounded, coordYrounded, income) |> 
+#   mutate(scenario = "base") |> 
+#   bind_rows(policy_merged |> select(id,coordXrounded,coordYrounded, income) |> mutate(scenario = "policy")) |>
+#   group_by(coordXrounded,coordYrounded,scenario) |>
+#   summarise(mean_income = mean(income, na.rm = TRUE), .groups = "drop")
+# 
+# 
+# xxx <-  zonal_counts |>
+#   pivot_wider(names_from = scenario, values_from = n) |>
+#   mutate(diff = policy - base)
+# 
+# xxx |> filter(coordXrounded == 0 & coordYrounded == 0)
+# 
+# yyy <- zonal_mean_incomes |> 
+#   pivot_wider(names_from = scenario, values_from = mean_income) |>
+#   mutate(diff = policy - base) 
+# 
+# income_pop_diff_table <- xxx |> left_join(yyy, by = c("coordXrounded" = "coordXrounded", "coordYrounded" = "coordYrounded"))
+# 
+# 
+# # diif plot
+#   # ggplot(xxx) +
+#   # net_pt +
+#   # net_total +
+#   # net_autofrei +
+#   # geom_point(aes(coordXrounded, coordYrounded, fill = diff, size = base),shape = 21) +
+#   # scale_fill_gradient2(
+#   #   low = "red",
+#   #   mid = "white",
+#   #   high = "blue",
+#   #   midpoint = 0,
+#   #   limits = c(-100, 100),
+#   #   oob = scales::squish
+#   # ) +
+#   # scale_size_continuous(range = c(1, 10)) +
+#   # theme_minimal(base_size = 14) +
+#   # labs(title = "Autofrei Policy vs. Base Case (at year 10)", fill = "pop diff", size = "Base Pop", x="", y="")
+#   
+#   ggplot(income_pop_diff_table) +
+#   net_pt +
+#   net_total +
+#   net_autofrei +
+#   geom_point(aes(coordXrounded, coordYrounded, fill = base.y, size = base.x),shape = 21) +
+#   scale_fill_gradient2(
+#     low = "red",
+#     mid = "white",
+#     high = "blue",
+#     midpoint = 13500,
+#     limits = c(10000, 17000),
+#     oob = scales::squish
+#   ) +
+#   scale_size_continuous(range = c(5, 15)) +
+#   theme_minimal(base_size = 14) +
+#   labs(title = "Base Case (at year 10)", fill = "mean income", size = "Base Pop", x="", y="")
+
+
+
+```
+
+Read Data
+```{r}
+
+```
+
+Merge Dataset
+```{r}
+
+```
+
+
+
+
+```{r}
+
+```
+
+
+
+
+## Plot Home Coords
+```{r}
+
+
+# sanity check
+# base_merged |> full_join(policy_merged, by = "id", suffix = c(".base",".policy")) |> 
+#   mutate(age_same = age.base == age.policy, gender_same = gender.base == gender.policy) |> filter(!age_same)
+
+
+# People moving away from carfree zone
+base_merged |> full_join(policy_merged, by = "id", suffix = c(".base",".policy")) |> 
+  filter(age.base == age.policy & gender.base == gender.policy) |> 
+  mutate(coordX.base = round(coordX.base/5000) *5000, coordY.base = round(coordY.base/5000) *5000) |>
+  mutate(coordX.policy = round(coordX.policy/5000) *5000, coordY.policy = round(coordY.policy/5000) *5000) |>
+  mutate(jobX.base = round(jobX.base/5000) *5000, jobY.base = round(jobY.base/5000) *5000) |> 
+  mutate(jobX.policy = round(jobX.policy/5000) *5000, jobY.policy = round(jobY.policy/5000) *5000) |> 
+ filter(coordX.base == 0 & coordY.base == 0) |> 
+  mutate(category = case_when(
+  is.na(coordX.policy) ~ "disappeared",
+  coordX.policy == 0 & coordY.policy == 0 ~ "stayed",
+  .default = "movedFromCBA"
+)) |> 
+  mutate(job_in_cba = jobX.base == 0 & jobY.base == 0) |>
+  group_by(category) |> 
+  summarise(n = n(),  income.base = mean(income.base), drivers_license = mean(driversLicense.base == TRUE), job_in_cba =mean(job_in_cba == TRUE, na.rm = TRUE),has_car = mean(autos.base > 0))
+
+# People moving to carfree zone
+base_merged |> full_join(policy_merged, by = "id", suffix = c(".base",".policy")) |> 
+  filter(age.base == age.policy & gender.base == gender.policy) |> 
+  mutate(coordX.base = round(coordX.base/5000) *5000, coordY.base = round(coordY.base/5000) *5000) |> 
+  mutate(coordX.policy = round(coordX.policy/5000) *5000, coordY.policy = round(coordY.policy/5000) *5000) |>
+  mutate(jobX.base = round(jobX.base/5000) *5000, jobY.base = round(jobY.base/5000) *5000) |> 
+  mutate(jobX.policy = round(jobX.policy/5000) *5000, jobY.policy = round(jobY.policy/5000) *5000) |> 
+ filter(coordX.policy == 0 & coordY.policy == 0) |> 
+  mutate(category = case_when(
+  is.na(coordX.base) ~ "appeared",
+  coordX.base == 0 & coordY.base == 0 ~ "stayed",
+  .default = "moved2CBA"
+)) |>  mutate(job_in_cba = jobX.base == 0 & jobY.base == 0) |> 
+  group_by(category) |> 
+  summarise(n = n(),  income.base = mean(income.base), drivers_license = mean(driversLicense.base == TRUE), job_in_cba =mean(job_in_cba == TRUE, na.rm = TRUE), has_car = mean(autos.base > 0))
+
+
+```
+
+
+## Skims
+
+```{r}
+# base_dir <- "/Users/jakob/git/silo/useCases/fabiland/scenario/scenOutput/2026-09-10a-base/"
+# policy_dir <- "/Users/jakob/git/silo/useCases/fabiland/scenario/scenOutput/2026-09-10a-policy/"
+
+# base_dir <- "/Users/jakob/git/silo/useCases/fabiland/scenario/scenOutput/2026-09-11a-base/"
+# policy_dir <- "/Users/jakob/git/silo/useCases/fabiland/scenario/scenOutput/2026-09-11b-policy/"
+
+year <- 10
+base_car_skim <- read_csv(paste0(base_dir,"matsim/skims/carSkim_",year,".csv"))
+policy_car_skim <- read_csv(paste0(policy_dir, "matsim/skims/carSkim_",year,".csv"))
+base_pt_skim <- read_csv(paste0(base_dir, "matsim/skims/ptSkim_",year,".csv"))
+policy_pt_skim <- read_csv(paste0(policy_dir, "matsim/skims/ptSkim_",year,".csv"))
+
+
+fabiland_shp <- read_sf("/Users/jakob/git/silo/useCases/fabiland/scenario/input/zonesShapefile/fabiland.shp")
+
+from_node = 13
+
+base_car_skim |> filter(FROM  == from_node) |> left_join(fabiland_shp, by = c("TO" = "Id")) |> st_as_sf() |> tm_shape() + tm_polygons(fill = "VALUE") + tm_title("Base - Car")
+
+policy_car_skim |> filter(FROM  == from_node) |> left_join(fabiland_shp, by = c("TO" = "Id")) |> st_as_sf() |> tm_shape() + tm_polygons(fill = "VALUE") + tm_title("Policy - Car")
+
+
+base_pt_skim |> filter(FROM  == from_node) |> left_join(fabiland_shp, by = c("TO" = "Id")) |>  
+st_as_sf() |> tm_shape() + tm_polygons(fill = "VALUE") + tm_title("Base - PT")
+
+policy_pt_skim |> filter(FROM  == from_node) |> left_join(fabiland_shp, by = c("TO" = "Id")) |> st_as_sf() |> tm_shape() + tm_polygons(fill = "VALUE") + tm_title("Policy - PT")
+
+```
+
+
+## Zonal System
+## Analytics
+```{r}
+# base_dir <- "/Users/jakob/git/silo/useCases/fabiland/scenario/scenOutput/2026-09-14e-base/"
+# policy_dir <- "/Users/jakob/git/silo/useCases/fabiland/scenario/scenOutput/2026-09-14e-policy/"
+
+base_analytics <- read_csv(paste0(base_dir,"avgSatisfactionByZone.csv")) |> left_join(fabiland_shp, by = c("zone" = "Id")) |> st_as_sf() 
+
+policy_analytics <- read_csv(paste0(policy_dir,"avgSatisfactionByZone.csv")) |> left_join(fabiland_shp, by = c("zone" = "Id")) |> st_as_sf()
+```
+
+
+```{r}
+year1 <- 0
+
+
+base_analytics |> filter(year == year1)  |> tm_shape() + tm_polygons(fill = "avgProp") + tm_title("Probability of Moving, Base, year 0")
+
+policy_analytics |> filter(year == year1)  |> tm_shape() + tm_polygons(fill = "avgProp") + tm_title("Probability of Moving, Policy, year 0")
+
+# avgSatisfaction1
+ggplot(base_analytics) + geom_line( aes(year, avgSatisfaction1,colour = as.factor(zone))) + geom_line(data = base_analytics |> filter(zone == 13), mapping = aes(year, avgSatisfaction1), col = "red", lwd = 2) + lims(y = c(0.0,0.6))+ labs(title = "avgSatisfaction1 + Base")
+
+
+ggplot(policy_analytics) + geom_line( aes(year, avgSatisfaction1,colour = as.factor(zone))) + geom_line(data = policy_analytics |> filter(zone == 13), mapping = aes(year, avgSatisfaction1), col = "red", lwd = 2)+ lims(y = c(0.0,0.6)) +labs(title = "avgSatisfaction1 + Policy")
+
+
+# avgSatisfaction2
+
+ggplot(base_analytics) + geom_line( aes(year, avgSatisfaction2,colour = as.factor(zone))) + geom_line(data = base_analytics |> filter(zone == 13), mapping = aes(year, avgSatisfaction2), col = "red", lwd = 2) + lims(y = c(0.0,0.6)) + labs(title = "avgSatisfaction2 + Base")
+
+
+ggplot(policy_analytics) + geom_line( aes(year, avgSatisfaction2,colour = as.factor(zone))) + geom_line(data = policy_analytics |> filter(zone == 13), mapping = aes(year, avgSatisfaction2), col = "red", lwd = 2)+ lims(y = c(0.0,0.6))+labs(title = "avgSatisfaction2 + Policy")
+
+
+# avgCurrentUtil
+
+ggplot(base_analytics) + geom_line( aes(year, avgCurrentUtil,colour = as.factor(zone))) + geom_line(data = base_analytics |> filter(zone == 13), mapping = aes(year, avgCurrentUtil), col = "red", lwd = 2) + lims(y = c(0.0,0.6)) + labs(title = "avgCurrentUtil + Base")
+
+
+ggplot(policy_analytics) + geom_line( aes(year, avgCurrentUtil,colour = as.factor(zone))) + geom_line(data = policy_analytics |> filter(zone == 13), mapping = aes(year, avgCurrentUtil), col = "red", lwd = 2)+ lims(y = c(0.0,0.6))+ labs(title = "avgCurrentUtil + Policy")
+
+# PROP
+
+ggplot(base_analytics) + geom_line( aes(year, avgProp,colour = as.factor(zone))) + geom_line(data = base_analytics |> filter(zone == 13), mapping = aes(year, avgProp), col = "red", lwd = 2)+ lims(y = c(0.0,0.4))+ labs(title = "avgProp + Base")
+
+
+ggplot(policy_analytics) + geom_line( aes(year, avgProp,colour = as.factor(zone))) + geom_line(data = policy_analytics |> filter(zone == 13), mapping = aes(year, avgProp), col = "red", lwd = 2)+ lims(y = c(0.0,0.4))+ labs(title = "avgProp + Policy")
+
+
+base_analytics |> mutate(scenario = "base") |> bind_rows(policy_analytics |> mutate(scenario = "policy")) |> st_drop_geometry() |>  filter(zone ==13) |>  ggplot() + geom_col(aes(year,avgCurrentUtil, fill = scenario), position = "dodge") + labs(title = "avgCurrentUtil + zone13")
+
+base_analytics |> mutate(scenario = "base") |> bind_rows(policy_analytics |> mutate(scenario = "policy")) |> st_drop_geometry() |>  filter(zone ==13) |>  ggplot() + geom_col(aes(year,avgSatisfaction1, fill = scenario), position = "dodge")+ labs(title = "avgSatisfaction1 + zone13")
+
+base_analytics |> mutate(scenario = "base") |> bind_rows(policy_analytics |> mutate(scenario = "policy")) |> st_drop_geometry() |>  filter(zone ==13) |>  ggplot() + geom_col(aes(year,avgProp, fill = scenario), position = "dodge")+ labs(title = "avgProp + zone13")
+
+
+### 
+
+base_analytics |> mutate(scenario = "base") |> bind_rows(policy_analytics |> mutate(scenario = "policy")) |> st_drop_geometry() |>  group_by(zone, scenario) |> summarise(avgCurrentUtil = mean(avgCurrentUtil)) |>  ggplot() + geom_col(aes(zone,avgCurrentUtil, fill = scenario), position = "dodge")+ labs(title = "avgCurrentUtil + average over years")+ geom_vline(xintercept = 13)
+
+base_analytics |> mutate(scenario = "base") |> bind_rows(policy_analytics |> mutate(scenario = "policy")) |> st_drop_geometry() |>  group_by(zone, scenario) |> summarise(avgSatisfaction1 = mean(avgSatisfaction1)) |>  ggplot() + geom_col(aes(zone,avgSatisfaction1, fill = scenario), position = "dodge")+ labs(title = "avgSatisfaction1 + average over years")+ geom_vline(xintercept = 13)
+
+base_analytics |> mutate(scenario = "base") |> bind_rows(policy_analytics |> mutate(scenario = "policy")) |> st_drop_geometry() |>  group_by(zone, scenario) |> summarise(avgProp = mean(avgProp)) |>  ggplot() + geom_col(aes(zone,avgProp, fill = scenario), position = "dodge")+ labs(title = "avgProp + average over years")+ geom_vline(xintercept = 13)
+```
+
+
+```{r}
+
+variables <- base_analytics |> st_drop_geometry() |> select(-c(zone,year)) |> colnames()
+# variables <- c("avgSatisfaction1","avgSatisfaction2","avgCurrentUtil","avgProp","defaultProp","ddQualityUtility","ddSizeUtility","ddAutoAccessibilityUtility","transitAccessibilityUtility","ddPriceUtility","workDistanceUtility") # Add your variable names here
+
+analytics <- base_analytics |>
+  mutate(scenario = "base") |>
+  bind_rows(policy_analytics |> mutate(scenario = "policy")) |>
+  st_drop_geometry()
+
+for (variable in variables) {
+  p <- analytics |>
+    group_by(zone, scenario) |>
+    summarise(value = mean(.data[[variable]]), .groups = "drop") |>
+    ggplot(aes(zone, value, fill = scenario)) +
+    geom_col(position = "dodge") +
+    labs(
+      title = paste(variable, "— average over years"),
+      y = variable
+    ) +
+    geom_vline(xintercept = 13)
+
+  print(p)
+}
+
+
+
+
+
+
+```
+## Plot Jobs
+```{r}
+base_merged |> ggplot() + geom_point(aes(jobX,jobY), alpha = 0.1) + theme_minimal(base_size = 14) + labs(title = "Base")
+policy_merged |> ggplot() + geom_point(aes(jobX,jobY), alpha = 0.1) + theme_minimal(base_size = 14) + labs(title = "Policy")
+
+
+zonal_counts_jobs <- base_merged |>
+  select(id, jobX, jobY) |> 
+ 
+  mutate(scenario = "base") |> 
+  bind_rows(policy_merged |> select(id,jobX, jobY) |> mutate(scenario = "policy")) |>
+  filter(!is.na(jobX)) |>
+  mutate(coordX = round(jobX/5000) *5000, coordY = round(jobY/5000) *5000) |> 
+  group_by(coordX, coordY,scenario) |> count()
+
+
+zonal_counts_jobs |> 
+  filter(scenario == "base") |> 
+  ggplot() + 
+  geom_point(aes(coordX,coordY, size = n)) + 
+  theme_minimal(base_size = 14) +
+  labs(title = "Base Counts")
+  
+zonal_counts_jobs |> 
+  filter(scenario == "policy") |> 
+  ggplot() + 
+  geom_point(aes(coordX,coordY, size = n)) + 
+  theme_minimal(base_size = 14) +
+  labs(title = "Policy Counts")
+
+
+zonal_counts_jobs |>
+  pivot_wider(names_from = scenario, values_from = n) |>
+  mutate(diff = policy - base) |>
+  ggplot() +
+  net_pt +
+  net_total +
+  net_autofrei +
+  geom_point(aes(coordX, coordY, fill = diff), size = 10,shape = 21) +
+  scale_fill_gradient2(
+    low = "red",
+    mid = "white",
+    high = "blue",
+    midpoint = 0
+  ) +
+  theme_minimal(base_size = 14) +
+  labs(title = "Policy - Base", fill = "Job Diff") 
+
+
+
+# jobs for to13 agents:
+policy_to13_jobs <- policy_to13 |> select(jobX_policy, jobY_policy) |> 
+    filter(!is.na(jobX_policy)) |>
+  mutate(coordX_policy = round(jobX_policy/5000) *5000, coordY_policy = round(jobY_policy/5000) *5000) |> 
+  group_by(coordX_policy, coordY_policy) |> count()
+
+policy_to13_jobs
+
+policy_from13_jobs <- policy_from13 |> select(jobX_policy, jobY_policy) |> 
+    filter(!is.na(jobX_policy)) |>
+  mutate(coordX_policy = round(jobX_policy/5000) *5000, coordY_policy = round(jobY_policy/5000) *5000) |> 
+  group_by(coordX_policy, coordY_policy) |> count()
+
+policy_from13_jobs
+
+policy_from13_jobs |> 
+  ggplot() + 
+  geom_point(aes(coordX_policy, coordY_policy, size = n)) +
+    scale_size_continuous(
+    breaks = c(10, 20, 30, 40, 50)
+  ) +
+  theme_minimal(base_size = 14) +
+  labs(title = "Job counts - from zone 13")
+
+policy_to13_jobs |> 
+  ggplot() + 
+  geom_point(aes(coordX_policy, coordY_policy, size = n)) + 
+  scale_size_continuous(
+    breaks = c(10, 20, 30, 40, 50)
+  ) +
+  theme_minimal(base_size = 14) +
+  labs(title = "Job counts - to zone 13")
+
+
+
+
+
+```
+
+
+## Deprecated Homes Plots
+```{r}
+
+base_merged |> ggplot() + geom_point(aes(coordX,coordY), alpha = 0.1) + theme_minimal(base_size = 14) + labs(title = "Base")
+  
+
+policy_merged |> ggplot() + geom_point(aes(coordX,coordY), alpha = 0.1) + theme_minimal(base_size = 14) + labs(title = "Policy")
+
+
+
+
+base_merged |>
+  select(id, coordX, coordY) |> 
+  left_join(policy_merged |> select(id, coordX, coordY), by = "id") |>
+  filter(!is.na(coordX.x) & !is.na(coordX.y)) |>
+  mutate(diff = coordX.x !=coordX.y) |> 
+  filter(diff) |> 
+  ggplot() + 
+  geom_point(aes(coordX.x,coordY.x), alpha = 0.1) +
+  theme_minimal(base_size = 14) +
+  labs(title = "Diff (Base)") 
+
+base_merged |>
+  select(id, coordX, coordY) |> 
+  left_join(policy_merged |> select(id, coordX, coordY), by = "id") |>
+  filter(!is.na(coordX.x) & !is.na(coordX.y)) |>
+  mutate(diff = coordX.x !=coordX.y) |> 
+  filter(diff) |> 
+  ggplot() + 
+  geom_point(aes(coordX.y,coordY.y), alpha = 0.1) +
+  theme_minimal(base_size = 14) +
+  labs(title = "Diff (Policy)") 
+zonal_counts |> 
+  filter(scenario == "base") |> 
+  ggplot() + 
+  geom_point(aes(coordXrounded,coordYrounded, col = n), size = 10) + 
+  theme_minimal(base_size = 14) +
+  labs(title = "Base Counts")
+  
+zonal_counts |> 
+  filter(scenario == "policy") |> 
+  ggplot() + 
+  geom_point(aes(coordXrounded,coordYrounded, col = n), size = 10) + 
+  theme_minimal(base_size = 14) +
+  labs(title = "Policy Counts")
+
+
+
+# net_total <-  list(annotate("segment", x = -10000, xend =10000,
+#            y = c(-10000,-5000, 0, 5000,10000),
+#            yend = c(-10000,-5000, 0, 5000,10000)),
+#   annotate("segment", y = -10000, yend =10000,
+#            x = c(-10000,-5000, 0, 5000,10000),
+#            xend = c(-10000,-5000, 0, 5000,10000)) )
+# 
+# 
+# net_autofrei <-  list(annotate("segment", x = -5000, xend =5000,
+#            y = c(-5000, 0, 5000),
+#            yend = c(-5000, 0, 5000), col = "red"),
+#   annotate("segment", y = -5000, yend =5000,
+#            x = c(-5000, 0, 5000),
+#            xend = c(-5000, 0, 5000),col = "red") )
+
+
+```
+
+# (UNUSED) BERLIN
+```{r}
+
+zones_sf <- read_sf("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/input/zonesShapefile/zones_Berlin-Brandenburg_size_3715_popCap_5000_crs_31468.shp")
+
+od_car_base <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08b-base/matsim/skims/carSkim_2012.csv")
+
+od_pt_base <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08b-base/matsim/skims/ptSkim_2012.csv")
+
+od_car_policy <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08-policy/matsim/skims/carSkim_2012.csv")
+
+od_pt_policy <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08-policy/matsim/skims/ptSkim_2012.csv")
+
+car_tt_by_origin_base <- od_car_base |> group_by(FROM) |> summarise(tt_mean = mean(VALUE)) |> select(id = FROM, tt_mean)
+pt_tt_by_origin_base <- od_pt_base |> group_by(FROM) |> summarise(tt_mean = mean(VALUE)) |> select(id = FROM, tt_mean)
+
+
+tmap_mode("view")
+od_car_base |> filter(TO == 1713) |> left_join(zones_sf, by = c("FROM" = "id")) |> st_as_sf() |> tm_shape() + tm_polygons(fill = "VALUE")
+
+
+# CAR
+tmap_mode("view")
+zones_sf |>
+  left_join(
+    od_car_base |> filter(TO == 1713),
+    by = c("id" = "FROM")
+  ) |>
+  tm_shape() +
+  tm_polygons(
+    fill = "VALUE",
+    fill.scale = tm_scale_intervals(
+      breaks = c(0, 50,100,500, 1000, 5000, 10000, Inf)
+    )
+  )
+
+
+# PT
+
+
+zones_sf |>
+  left_join(
+    od_pt_base |> filter(TO == 1713),
+    by = c("id" = "FROM")
+  ) |>
+  tm_shape() +
+  tm_polygons(
+    fill = "VALUE",
+    fill.scale = tm_scale_intervals(
+      breaks = c(0, 50,100,500, 1000, 5000, 10000, Inf)
+    )
+  )
+
+
+pt_tt_by_origin_base |> left_join(zones_sf) |> st_as_sf() |> tm_shape() + tm_polygons(fill = "tt_mean")
+
+```
+
+```{r}
+#base
+base_pp <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08-base/microData/pp_2050.csv")
+base_hh <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08-base/microData/hh_2050.csv")
+base_dd <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08-base/microData/dd_2050.csv")
+
+
+policy_pp <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08-policy/microData/pp_2050.csv")
+policy_hh <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08-policy/microData/hh_2050.csv")
+policy_dd <- read_csv("/Users/jakob/git/silo/useCases/berlinBrandenburg/scenario/scenOutput/2026-09-08-policy/microData/dd_2050.csv")
+
+tmap_mode("plot")
+base_pp |> left_join(base_hh, by = c("hhid" = "id")) |> left_join(base_dd, by = c("hhid"="id")) |> st_as_sf(coords = c("coordX","coordY"), crs = 25832) |> tm_shape() + tm_dots() + tm_basemap("Esri.WorldGrayCanvas")
+
+
+
+policy_pp |> left_join(base_hh, by = c("hhid" = "id")) |> left_join(base_dd, by = c("hhid"="id")) |> st_as_sf(coords = c("coordX","coordY"), crs = 25832) |> tm_shape() + tm_dots() + tm_basemap("Esri.WorldGrayCanvas")
+```
+
+
diff --git a/useCases/fabiland/scripts/simplifiedSiloAnalysis.Rmd b/useCases/fabiland/scripts/simplifiedSiloAnalysis.Rmd
new file mode 100644
index 000000000..c68fabf0d
--- /dev/null
+++ b/useCases/fabiland/scripts/simplifiedSiloAnalysis.Rmd
@@ -0,0 +1,60 @@
+---
+title: "simplifiedSiloAnalysis"
+output: html_document
+---
+
+```{r}
+# install.packages("tidyverse")
+library(tidyverse)
+
+```
+
+```{r}
+dir <- "/home/brendan/git/silo/useCases/fabiland/scenario/scenOutput/base/microData/"
+
+```
+
+Persons
+```{r}
+pp_start <- read_csv(paste0(dir,"pp_0.csv"))
+pp_end <- read_csv(paste0(dir,"pp_10.csv"))
+
+pp_comb <- full_join(pp_start,pp_end, by = "id")
+
+# how many people were born? = 498
+pp_comb %>% filter(is.na(age.x))
+
+# how many people died = 596
+pp_comb %>% filter(is.na(age.y))
+
+# is the age diff always the same? age diff is always 9.
+
+pp_comb_normal <- pp_comb %>%  filter(!is.na(age.x) & !is.na(age.y)) 
+pp_comb_normal %>% mutate(age_diff = age.y - age.x)
+
+
+#  how many people changed their workplace? = Half changed their workplace, 
+pp_comb_normal %>% mutate(job_changed = case_when(
+  workplace.x == -1 & workplace.y == -1 ~ "stay_unemployed",
+  workplace.x == -1 & workplace.y != -1 ~ "got_job",
+  workplace.x != -1 & workplace.y == -1 ~ "retired",
+  workplace.x != workplace.y ~ "changed_job",
+  .default = "stay_job"
+)) %>% group_by(job_changed) %>% count()
+
+# how many people got married or divorced (HH changed)?
+pp_comb_normal %>% mutate(relationship_changed = case_when(
+  relationShip.x != "married" & relationShip.y == "married" ~ "got married",
+  relationShip.x == "married" & relationShip.y != "married" ~ "got divorced",
+  .default = "no relationship change"
+))  %>% group_by(relationship_changed) %>% count()
+
+
+# how many people moved? 
+pp_comb_normal %>% mutate(household_changed = case_when(
+  hhid.x != hhid.y ~ "moved household",
+  .default = "stayed put"
+)) %>% group_by(household_changed) %>% count()
+
+```
+
diff --git a/useCases/fabiland/src/main/java/inputGeneration/PTScheduleCreatorRingsAndDiagonals.java b/useCases/fabiland/src/main/java/inputGeneration/PTScheduleCreatorRingsAndDiagonals.java
new file mode 100644
index 000000000..7a42eb983
--- /dev/null
+++ b/useCases/fabiland/src/main/java/inputGeneration/PTScheduleCreatorRingsAndDiagonals.java
@@ -0,0 +1,124 @@
+package inputGeneration;
+
+import org.apache.logging.log4j.LogManager;
+import org.apache.logging.log4j.Logger;
+import org.matsim.api.core.v01.Id;
+import org.matsim.api.core.v01.Scenario;
+import org.matsim.api.core.v01.TransportMode;
+import org.matsim.api.core.v01.network.Node;
+import org.matsim.core.config.ConfigUtils;
+import org.matsim.core.network.io.MatsimNetworkReader;
+import org.matsim.core.network.io.NetworkWriter;
+import org.matsim.core.scenario.ScenarioUtils;
+import org.matsim.pt.transitSchedule.api.*;
+import org.matsim.pt.utils.CreatePseudoNetwork;
+import org.matsim.pt.utils.CreateVehiclesForSchedule;
+import org.matsim.pt.utils.TransitScheduleValidator;
+import org.matsim.vehicles.MatsimVehicleWriter;
+
+import java.nio.file.Path;
+import java.util.ArrayList;
+import java.util.List;
+
+/** Generates both rings and diagonals, including a combined road/PT network and vehicles. */
+public class PTScheduleCreatorRingsAndDiagonals {
+    private static final Logger LOG = LogManager.getLogger(PTScheduleCreatorRingsAndDiagonals.class);
+    private static final Path DEFAULT_NETWORK = Path.of(
+            "useCases/fabiland/scenario/matsimInput/nw_cap30.xml");
+    private static final double PT_HEADWAY = 600;
+    private static final double PT_OPERATION_START_TIME = 6 * 3600;
+    private static final double PT_OPERATION_END_TIME = 24 * 3600;
+    private static final double RING_SEGMENT_TIME = 360;
+    private static final double DIAGONAL_SEGMENT_TIME = 500;
+
+    /** Run from the repository root, or supply the road network path as the only argument. */
+    public static void main(String[] args) {
+        if (args.length > 1) {
+            throw new IllegalArgumentException("Usage: PTScheduleCreatorRingsAndDiagonals [road-network.xml]");
+        }
+        Path inputNetwork = (args.length == 0 ? DEFAULT_NETWORK : Path.of(args[0])).toAbsolutePath();
+        Path outputDirectory = inputNetwork.getParent();
+        Scenario scenario = createScenario(inputNetwork);
+
+        new NetworkWriter(scenario.getNetwork()).write(
+                outputDirectory.resolve("nw_cap30_rings_x.xml").toString());
+        // Vehicle generation must precede schedule writing to persist departure vehicle IDs.
+        new TransitScheduleWriter(scenario.getTransitSchedule()).writeFile(
+                outputDirectory.resolve("ts_rings_x.xml").toString());
+        new MatsimVehicleWriter(scenario.getTransitVehicles()).writeFile(
+                outputDirectory.resolve("tv_rings_x.xml").toString());
+        LOG.info("Wrote rings and diagonals network, schedule and vehicles to {}", outputDirectory);
+    }
+
+    public static Scenario createScenario(Path inputNetwork) {
+        Scenario scenario = ScenarioUtils.createScenario(ConfigUtils.createConfig());
+        new MatsimNetworkReader(scenario.getNetwork()).readFile(inputNetwork.toString());
+        if (scenario.getNetwork().getNodes().keySet().stream().anyMatch(id -> id.toString().startsWith("pt-"))
+                || scenario.getNetwork().getLinks().keySet().stream().anyMatch(id -> id.toString().startsWith("pt-"))) {
+            throw new IllegalArgumentException("Input must be a road network without the reserved pt- prefix");
+        }
+        TransitSchedule schedule = scenario.getTransitSchedule();
+        TransitScheduleFactory factory = schedule.getFactory();
+        for (int number = 1; number <= 25; number++) {
+            Node node = scenario.getNetwork().getNodes().get(Id.createNodeId(number));
+            if (node == null) {
+                throw new IllegalArgumentException("Road network is missing grid node " + number);
+            }
+            TransitStopFacility stop = factory.createTransitStopFacility(
+                    Id.create(number, TransitStopFacility.class), node.getCoord(), false);
+            stop.setStopAreaId(Id.create(number, TransitStopArea.class));
+            schedule.addStopFacility(stop);
+        }
+
+        addLine(schedule, "outer-ring", RING_SEGMENT_TIME,
+                21, 22, 23, 24, 25, 20, 15, 10, 5, 4, 3, 2, 1, 6, 11, 16, 21);
+        addLine(schedule, "inner-ring", RING_SEGMENT_TIME,
+                7, 8, 9, 14, 19, 18, 17, 12, 7);
+        addLine(schedule, "x-ne", DIAGONAL_SEGMENT_TIME, 5, 9, 13, 17, 21);
+        addLine(schedule, "x-nw", DIAGONAL_SEGMENT_TIME, 1, 7, 13, 19, 25);
+
+        new CreatePseudoNetwork(schedule, scenario.getNetwork(), "pt-").createNetwork();
+        new CreateVehiclesForSchedule(schedule, scenario.getTransitVehicles()).run(TransportMode.pt);
+
+        TransitScheduleValidator.ValidationResult validation =
+                TransitScheduleValidator.validateAll(schedule, scenario.getNetwork());
+        validation.getWarnings().forEach(LOG::warn);
+        if (!validation.getErrors().isEmpty()) {
+            throw new IllegalStateException("Invalid transit schedule: " + validation.getErrors());
+        }
+        return scenario;
+    }
+
+    private static void addLine(TransitSchedule schedule, String lineId, double segmentTime, int... nodes) {
+        TransitScheduleFactory factory = schedule.getFactory();
+        TransitLine line = factory.createTransitLine(Id.create(lineId, TransitLine.class));
+        for (int direction = 0; direction < 2; direction++) {
+            List<TransitRouteStop> stops = new ArrayList<>();
+            for (int index = 0; index < nodes.length; index++) {
+                int node = nodes[direction == 0 ? index : nodes.length - 1 - index];
+                TransitStopFacility facility = schedule.getFacilities().get(Id.create(node, TransitStopFacility.class));
+                if (index == 0 && nodes[0] == nodes[nodes.length - 1]) {
+                    // A ring revisits its origin. Sharing its departure platform with another
+                    // line makes Raptor discard transfers there as an earlier stop on the ring.
+                    TransitStopFacility departurePlatform = factory.createTransitStopFacility(
+                            Id.create(lineId + "-" + (direction + 1) + "-start", TransitStopFacility.class),
+                            facility.getCoord(), false);
+                    departurePlatform.setStopAreaId(facility.getStopAreaId());
+                    schedule.addStopFacility(departurePlatform);
+                    facility = departurePlatform;
+                }
+                double offset = index * segmentTime;
+                TransitRouteStop stop = factory.createTransitRouteStop(facility, offset, offset);
+                stop.setAwaitDepartureTime(false);
+                stops.add(stop);
+            }
+            TransitRoute route = factory.createTransitRoute(
+                    Id.create(lineId + "-" + (direction + 1), TransitRoute.class), null, stops, TransportMode.pt);
+            for (double time = PT_OPERATION_START_TIME; time < PT_OPERATION_END_TIME; time += PT_HEADWAY) {
+                route.addDeparture(factory.createDeparture(Id.create(String.valueOf(time), Departure.class), time));
+            }
+            line.addRoute(route);
+        }
+        schedule.addTransitLine(line);
+    }
+}
diff --git a/useCases/fabiland/src/main/java/run/ModelBuilderFabilandSimplified.java b/useCases/fabiland/src/main/java/run/ModelBuilderFabilandSimplified.java
index 5a29bc9d8..4d7bf15ab 100644
--- a/useCases/fabiland/src/main/java/run/ModelBuilderFabilandSimplified.java
+++ b/useCases/fabiland/src/main/java/run/ModelBuilderFabilandSimplified.java
@@ -2,11 +2,10 @@
 
 import de.tum.bgu.msm.container.DataContainer;
 import de.tum.bgu.msm.container.ModelContainer;
-import de.tum.bgu.msm.data.Region;
-import de.tum.bgu.msm.data.dwelling.Dwelling;
 import de.tum.bgu.msm.data.dwelling.DwellingFactory;
 import de.tum.bgu.msm.data.household.Household;
 import de.tum.bgu.msm.data.household.HouseholdFactory;
+import de.tum.bgu.msm.data.person.Person;
 import de.tum.bgu.msm.data.person.PersonFactory;
 import de.tum.bgu.msm.events.impls.household.MoveEvent;
 import de.tum.bgu.msm.matsim.*;
@@ -38,24 +37,18 @@
 import de.tum.bgu.msm.models.jobmography.JobMarketUpdate;
 import de.tum.bgu.msm.models.jobmography.JobMarketUpdateImpl;
 import de.tum.bgu.msm.models.realEstate.construction.*;
-import de.tum.bgu.msm.models.realEstate.demolition.DefaultDemolitionStrategy;
-import de.tum.bgu.msm.models.realEstate.demolition.DemolitionModel;
-import de.tum.bgu.msm.models.realEstate.demolition.DemolitionModelImpl;
 import de.tum.bgu.msm.models.realEstate.pricing.DefaultPricingStrategy;
 import de.tum.bgu.msm.models.realEstate.pricing.PricingModel;
 import de.tum.bgu.msm.models.realEstate.pricing.PricingModelImpl;
-import de.tum.bgu.msm.models.realEstate.renovation.DefaultRenovationStrategy;
-import de.tum.bgu.msm.models.realEstate.renovation.RenovationModel;
-import de.tum.bgu.msm.models.realEstate.renovation.RenovationModelImpl;
 import de.tum.bgu.msm.models.relocation.migration.InOutMigration;
 import de.tum.bgu.msm.models.relocation.migration.InOutMigrationImpl;
 import de.tum.bgu.msm.models.relocation.moves.*;
 import de.tum.bgu.msm.models.transportModel.TransportModel;
 import de.tum.bgu.msm.properties.Properties;
 import de.tum.bgu.msm.utils.SiloUtil;
-import models.FabilandConstructionLocationStrategy;
 import org.matsim.api.core.v01.Scenario;
 import org.matsim.core.config.Config;
+import org.matsim.core.scenario.MutableScenario;
 import org.matsim.core.scenario.ScenarioUtils;
 
 import java.util.Collection;
@@ -68,10 +61,24 @@
 public class ModelBuilderFabilandSimplified{
 
     public static ModelContainer getModelContainer(DataContainer dataContainer, Properties properties, Config config) {
+        // As a base assumption, let us say that people work at "facilities", which we can imagine as Berlin-type blocks.  These
+        // would have something like 20k or even more sqm.  We could now assume that everybody just stays where they are, but
+        // dwellings grow or shring as necessary.  E.g. if people get children, the dwelling grows to a plausible size; and if
+        // people die, the dwelling vanisches and leaves sqm for others.  We could even imagine children moving out but staying in
+        // the same facility.  Also "male" partners would not move in.
+
+        // My intuition is that the first order error is not so bad, e.g.:
+        // * For every "male" partner not moving in, there will be some other "male" partner staying in the facility with having its "female" partner elsewhere.
+        // * For each family not moving out when getting children there will be some other person dying and making space.
+        // * Etc.
+
+        // I think that for such an approach we can simplify even more than what I have done below.
 
         final PersonFactory ppFactory = dataContainer.getHouseholdDataManager().getPersonFactory();
         final HouseholdFactory hhFactory = dataContainer.getHouseholdDataManager().getHouseholdFactory();
         final DwellingFactory ddFactory = dataContainer.getRealEstateDataManager().getDwellingFactory();
+        // (in most if not all cases, the dataContainer is handed over anyways.  --> add separate constructor w/o those factories;
+        // set old constructor to deprecated (but do not make effort to remove).
 
         final BirthModel birthModel = new BirthModelImpl(dataContainer, ppFactory, properties, new DefaultBirthStrategy(), SiloUtil.provideNewRandom());
 
@@ -109,18 +116,23 @@ public static ModelContainer getModelContainer(DataContainer dataContainer, Prop
         };
 
         final CreateCarOwnershipModel carOwnershipModel = new FabilandCarOwnership();
-        // yy (for VSP purposes, car ownership could also be done by matsim)
+        // yy (for VSP purposes, car ownership could  be done by matsim)
 
         final DivorceModel divorceModel = new DivorceModelImpl(
                 dataContainer, movesModel, carOwnershipModel, hhFactory,
                 properties, new DefaultDivorceStrategy(), SiloUtil.provideNewRandom());
 
         final DriversLicenseModel driversLicenseModel = new DriversLicenseModelImpl(dataContainer, properties, new DefaultDriversLicenseStrategy(), SiloUtil.provideNewRandom());
-        // yy for VSP purposes, this might not be needed
+        // yy (for VSP purposes, this might not be needed)
 
         final EducationModel educationModel = new EducationModelImpl(dataContainer, properties, SiloUtil.provideNewRandom());
 
-        final EmploymentModel employmentModel = new EmploymentModelImpl(dataContainer, properties, SiloUtil.provideNewRandom());
+        EmploymentModel employmentModel = new EmploymentModelImpl(dataContainer, properties, SiloUtil.provideNewRandom()) {
+            @Override
+            public boolean lookForJob(int perId) {
+                return false;   // stub this one
+            }
+        };
 
         final LeaveParentHhModel leaveParentsModel = new LeaveParentHhModelImpl(dataContainer, movesModel,
                 carOwnershipModel, hhFactory, properties, new DefaultLeaveParentalHouseholdStrategy(), SiloUtil.provideNewRandom());
@@ -133,9 +145,17 @@ public static ModelContainer getModelContainer(DataContainer dataContainer, Prop
 
         final InOutMigration inOutMigration = new InOutMigrationImpl(dataContainer, employmentModel, movesModel,
                 carOwnershipModel, driversLicenseModel, properties, SiloUtil.provideNewRandom());
+        // (do we need this at VSP?)
+        // (if we need it, the car ownership model can be null.  Presumably, also the drivers licence model.)
 
         final MarriageModel marriageModel = new MarriageModelImpl(dataContainer, movesModel, inOutMigration,
-                carOwnershipModel, hhFactory, properties, new DefaultMarriageStrategy(), SiloUtil.provideNewRandom());
+                carOwnershipModel, hhFactory, properties, new DefaultMarriageStrategy(), SiloUtil.provideNewRandom()) {
+            @Override
+            protected boolean moveTogether(Person person1, Person person2, Household moveTo) {
+                return true;
+            };
+        };
+        // (do we need this at VSP?  We could also have women have children.)
 
         TransportModel transportModel;
         MatsimScenarioAssembler scenarioAssembler;
@@ -143,8 +163,13 @@ public static ModelContainer getModelContainer(DataContainer dataContainer, Prop
         MatsimData matsimData = null;
         if (config != null) {
             final Scenario scenario = ScenarioUtils.loadScenario(config);
-            matsimData = new MatsimData(config, properties, ZoneConnectorMethod.WEIGHTED_BY_POPULATION, dataContainer, scenario.getNetwork(), scenario.getTransitSchedule());
+            //  seems to be breaking
+//            matsimData = new MatsimData(config, properties, ZoneConnectorMethod.WEIGHTED_BY_POPULATION, dataContainer, scenario.getNetwork(), scenario.getTransitSchedule());
             // (only the constructor is deprecated)
+
+            // old version
+            matsimData = new MatsimData(properties, ZoneConnectorMethod.WEIGHTED_BY_POPULATION, dataContainer, (MutableScenario) scenario );
+
         }
         switch (properties.transportModel.transportModelIdentifier) {
             case MATSIM:
@@ -163,16 +188,16 @@ public static ModelContainer getModelContainer(DataContainer dataContainer, Prop
         final ModelContainer modelContainer = new ModelContainer(
                 birthModel, birthdayModel,
                 deathModel, marriageModel,
-                divorceModel, driversLicenseModel,
-                educationModel, employmentModel,
-                leaveParentsModel, jobMarketUpdateModel,
-                null, null, pricing, null,
-                constructionOverwrite, inOutMigration, movesModel, transportModel);
+                divorceModel, null,
+                null, employmentModel,
+                null, null,
+                null, null, null, null,
+                null, inOutMigration, movesModel, transportModel);
 
         return modelContainer;
     }
 
-    private static class FabilandCarOwnership implements CreateCarOwnershipModel {
+    protected static class FabilandCarOwnership implements CreateCarOwnershipModel {
 
         private final Random random;
 
diff --git a/useCases/fabiland/src/main/java/run/RunFabiland.java b/useCases/fabiland/src/main/java/run/RunFabiland.java
index ebe039970..5d5c0707e 100644
--- a/useCases/fabiland/src/main/java/run/RunFabiland.java
+++ b/useCases/fabiland/src/main/java/run/RunFabiland.java
@@ -6,7 +6,6 @@
 import de.tum.bgu.msm.io.output.DefaultResultsMonitor;
 import de.tum.bgu.msm.io.output.HouseholdSatisfactionMonitor;
 import de.tum.bgu.msm.io.output.MultiFileResultsMonitor;
-import de.tum.bgu.msm.io.output.ResultsMonitor;
 import de.tum.bgu.msm.properties.Properties;
 import de.tum.bgu.msm.utils.SiloUtil;
 import org.apache.logging.log4j.LogManager;
@@ -25,43 +24,49 @@ public class RunFabiland {
     public static void main(String[] args) {
         // see regression test
 
-        Properties properties = SiloUtil.siloInitialization(args[0]);
+        // args: SILO config, MATSim config
+        // e.g., "useCases/fabiland/scenario/1r_ae.properties useCases/fabiland/scenario/config_cap30_1-l_nes_smc.xml"
+        // or, to match regression test ...
+        // "useCases/fabiland/scenario/test.properties useCases/fabiland/scenario/config_cap30_1-l_nes_smc.xml --config:controller.lastIteration 1 "
+
+        Properties siloConfig = SiloUtil.siloInitialization(args[0]);
 
         String[] matsimArgs = Arrays.copyOfRange( args, 1, args.length );
 
-        Config config = null;
+        Config matsimConfig = null;
 //        if (args.length > 1 && args[1] != null) {
-            config = ConfigUtils.loadConfig(matsimArgs);
+        matsimConfig = ConfigUtils.loadConfig(matsimArgs);
 //        }
         logger.info("Started SILO Fabiland sandbox model");
 
         // The following is obviously just a dirty quickfix until access/egress is default in MATSim
-        if (properties.transportModel.includeAccessEgress) {
+        if (siloConfig.transportModel.includeAccessEgress) {
 ////            config.plansCalcRoute().setInsertingAccessEgressWalk(true); // in matsim-12
-            config.routing().setAccessEgressType(RoutingConfigGroup.AccessEgressType.accessEgressModeToLink); // in matsim-13-w37
+            matsimConfig.routing().setAccessEgressType(RoutingConfigGroup.AccessEgressType.accessEgressModeToLink); // in matsim-13-w37
         }
 //		config.routing().setAccessEgressType( RoutingConfigGroup.AccessEgressType.none );
-		// yyyyyy Silo uses a re-implementation of a lot of matsim infrastructure, and that is outside injection.  The more advanced access/egress types are not implemented there.
-		// kai, apr'26
+        // yyyyyy Silo uses a re-implementation of a lot of matsim infrastructure, and that is outside injection.  The more advanced access/egress types are not implemented there.
+        // kai, apr'26
 
-		config.controller().setOverwriteFileSetting( OutputDirectoryHierarchy.OverwriteFileSetting.overwriteExistingFiles );
-		// Somehow, some version matsim is starting again for the accessibility computation, and that wipes the directory after the main run.
-		// --> did not help
+        matsimConfig.controller().setOverwriteFileSetting( OutputDirectoryHierarchy.OverwriteFileSetting.overwriteExistingFiles );
+        // Somehow, some version matsim is starting again for the accessibility computation, and that wipes the directory after the main run.
+        // --> did not help
 
         logger.warn("Constructing data container ...");
-        DataContainer dataContainer = DataBuilderFabiland.buildDataContainer(properties, config);
-        DataBuilderFabiland.readInput(properties, dataContainer);
+        DataContainer dataContainer = DataBuilderFabiland.buildDataContainer(siloConfig, matsimConfig);
+        DataBuilderFabiland.readInput(siloConfig, dataContainer);
         logger.warn("... done with constructing data container.");
 
         logger.warn("Constructing model container ...");
-        ModelContainer modelContainer = ModelBuilderFabiland.getModelContainer(dataContainer, properties, config);
+        ModelContainer modelContainer = ModelBuilderFabiland.getModelContainer(dataContainer, siloConfig, matsimConfig);
+//        ModelContainer modelContainer = ModelBuilderFabilandSimplified.getModelContainer(dataContainer, siloConfig, matsimConfig);
         logger.warn("... done with constructing model container.");
 
         logger.warn("Constructing silo model ...");
-		SiloModel model = new SiloModel(properties, dataContainer, modelContainer);
-        model.addResultMonitor( new DefaultResultsMonitor(dataContainer, properties) );
-        model.addResultMonitor( new MultiFileResultsMonitor(dataContainer, properties) );
-        model.addResultMonitor( new HouseholdSatisfactionMonitor(dataContainer, properties, modelContainer) );
+        SiloModel model = new SiloModel(siloConfig, dataContainer, modelContainer);
+        model.addResultMonitor( new DefaultResultsMonitor(dataContainer, siloConfig) );
+        model.addResultMonitor( new MultiFileResultsMonitor(dataContainer, siloConfig) );
+        model.addResultMonitor( new HouseholdSatisfactionMonitor(dataContainer, siloConfig, modelContainer) );
         logger.warn("... done with constructing silo model.");
 
         logger.warn("Running silo model ...");
diff --git a/useCases/fabiland/src/main/java/run/RunFabilandAutofrei.java b/useCases/fabiland/src/main/java/run/RunFabilandAutofrei.java
new file mode 100644
index 000000000..30f829ae0
--- /dev/null
+++ b/useCases/fabiland/src/main/java/run/RunFabilandAutofrei.java
@@ -0,0 +1,168 @@
+package run;
+
+import de.tum.bgu.msm.SiloModel;
+import de.tum.bgu.msm.container.DataContainer;
+import de.tum.bgu.msm.container.ModelContainer;
+import de.tum.bgu.msm.io.output.DefaultResultsMonitor;
+import de.tum.bgu.msm.io.output.HouseholdSatisfactionMonitor;
+import de.tum.bgu.msm.io.output.MultiFileResultsMonitor;
+import de.tum.bgu.msm.properties.Properties;
+import de.tum.bgu.msm.utils.SiloUtil;
+import org.apache.logging.log4j.LogManager;
+import org.apache.logging.log4j.Logger;
+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.network.NetworkWriter;
+import org.matsim.api.core.v01.network.Node;
+import org.matsim.core.config.Config;
+import org.matsim.core.config.ConfigUtils;
+import org.matsim.core.config.groups.RoutingConfigGroup;
+import org.matsim.core.controler.OutputDirectoryHierarchy;
+import org.matsim.core.network.NetworkUtils;
+import org.matsim.core.network.io.MatsimNetworkReader;
+
+import java.io.IOException;
+import java.nio.file.Files;
+import java.nio.file.Path;
+import java.util.*;
+import java.util.stream.Collectors;
+
+public class RunFabilandAutofrei {
+
+    enum AutofreiScenario {
+        BASE,
+        SMALL,
+        LARGE
+    }
+
+    enum PtScenario {
+//        NONE,
+        NES,
+        RINGSX
+    }
+
+
+    //    static String scenario = "base";
+    static AutofreiScenario autofreiScenario = AutofreiScenario.LARGE;
+    static PtScenario ptScenario = PtScenario.NES; // nes / ringsX
+
+    private final static Logger logger = LogManager.getLogger(RunFabilandAutofrei.class);
+
+    public static void main(String[] args) throws IOException {
+
+
+
+        // Create a temporary SILO config, wherein we overwrite the scenario name:
+        String prefix = "2026-10-23a__";
+        String scenario = prefix + "autofrei-" + autofreiScenario + "__pt-" + ptScenario;
+        Path baseSiloPropertiesFile = Path.of("useCases/fabiland/scenario/25r_ae.properties");
+        Path tempSiloPropertiesFile = Files.createTempFile(baseSiloPropertiesFile.toAbsolutePath().getParent(), "tmp-", ".properties");
+        Files.writeString(tempSiloPropertiesFile,
+                Files.readString(baseSiloPropertiesFile) + "\nscenario.name = " + scenario + "\n");
+        tempSiloPropertiesFile.toFile().deleteOnExit();
+
+
+        // load Silo properties
+        Properties siloConfig = SiloUtil.siloInitialization(tempSiloPropertiesFile.toString());
+
+        // load MATSim config & modify certain options
+        Config matsimConfig = ConfigUtils.loadConfig("useCases/fabiland/scenario/config_cap5_1-l_nes_smc.xml");
+
+        logger.info("Started SILO Fabiland sandbox model");
+
+        // The following is obviously just a dirty quickfix until access/egress is default in MATSim
+        if (siloConfig.transportModel.includeAccessEgress) {
+            matsimConfig.routing().setAccessEgressType(RoutingConfigGroup.AccessEgressType.accessEgressModeToLink); // in matsim-13-w37
+        }
+//		config.routing().setAccessEgressType( RoutingConfigGroup.AccessEgressType.none );
+        // yyyyyy Silo uses a re-implementation of a lot of matsim infrastructure, and that is outside injection.  The more advanced access/egress types are not implemented there.
+        // kai, apr'26
+
+        matsimConfig.controller().setOverwriteFileSetting( OutputDirectoryHierarchy.OverwriteFileSetting.overwriteExistingFiles );
+//        matsimConfig.controller().setLastIteration(50);
+
+
+//        if (ptScenario.equals(PtScenario.NONE)) {
+//            matsimConfig.network().setInputFile("matsimInput/nw_cap30.xml");
+//            matsimConfig.transit().setTransitScheduleFile(null);
+//        } else
+        if (ptScenario.equals(PtScenario.NES)) {
+            // do nothing
+        } else if (ptScenario.equals(PtScenario.RINGSX)) {
+            matsimConfig.network().setInputFile("matsimInput/nw_cap30_rings_x.xml");
+            matsimConfig.transit().setTransitScheduleFile("matsimInput/ts_rings_x.xml");
+//            matsimConfig.transit().setVehiclesFile("matsimInput/tv_rings_x.xml");
+        } else {
+            throw new RuntimeException();
+        }
+        
+
+        if (autofreiScenario.equals(AutofreiScenario.BASE)) {
+            // do nothing
+        } else if (autofreiScenario.equals(AutofreiScenario.SMALL)) {
+            Network network = NetworkUtils.createNetwork();
+            new MatsimNetworkReader(network).readFile("useCases/fabiland/scenario/" + matsimConfig.network().getInputFile());
+
+            // Carfree
+            Node node13 = network.getNodes().get(Id.createNodeId(13));
+            Set<Link> carfreeLinks = new HashSet<>(node13.getInLinks().values());
+            carfreeLinks.addAll(node13.getOutLinks().values());
+
+            for (Link carfreeLink : carfreeLinks) {
+                carfreeLink.setFreespeed(carfreeLink.getFreespeed() / 100);
+            }
+
+
+            new NetworkWriter(network).write("useCases/fabiland/scenario/matsimInput/_nw_jr_tmp.xml");
+            matsimConfig.network().setInputFile("matsimInput/_nw_jr_tmp.xml");
+        } else if (autofreiScenario.equals(AutofreiScenario.LARGE)) {
+            Network network = NetworkUtils.createNetwork();
+            new MatsimNetworkReader(network).readFile("useCases/fabiland/scenario/" + matsimConfig.network().getInputFile());
+
+            // all links must have a centroid within +-5000m
+            Set<Link> carfreeLinks = network.getLinks().values().stream().filter(link -> Math.abs(link.getCoord().getX()) <= 5000 & Math.abs(link.getCoord().getY()) <= 5000).collect(Collectors.toSet());
+            for (Link carfreeLink : carfreeLinks) {
+                carfreeLink.setFreespeed(carfreeLink.getFreespeed() / 100);
+            }
+
+            new NetworkWriter(network).write("useCases/fabiland/scenario/matsimInput/_nw_jr_tmp.xml");
+            matsimConfig.network().setInputFile("matsimInput/_nw_jr_tmp.xml");
+
+        } else {
+            throw new RuntimeException();
+        }
+
+//        RoutingConfigGroup.TeleportedModeParams walkParams = matsimConfig.routing().getTeleportedModeParams().get(TransportMode.walk);
+//        walkParams.setTeleportedModeSpeed(walkParams.getTeleportedModeSpeed() * 10);
+//
+//        RoutingConfigGroup.TeleportedModeParams nNwalkParams = matsimConfig.routing().getTeleportedModeParams().get(TransportMode.non_network_walk);
+//        nNwalkParams.setTeleportedModeSpeed(nNwalkParams.getTeleportedModeSpeed() * 10);
+        // Somehow, some version matsim is starting again for the accessibility computation, and that wipes the directory after the main run.
+        // --> did not help
+
+        logger.warn("Constructing data container ...");
+        DataContainer dataContainer = DataBuilderFabiland.buildDataContainer(siloConfig, matsimConfig);
+        DataBuilderFabiland.readInput(siloConfig, dataContainer);
+        logger.warn("... done with constructing data container.");
+
+
+        logger.warn("Constructing model container ...");
+        ModelContainer modelContainer = ModelBuilderFabiland.getModelContainer(dataContainer, siloConfig, matsimConfig);
+
+//        ModelContainer modelContainer = ModelBuilderFabilandSimplified.getModelContainer(dataContainer, siloConfig, matsimConfig);
+        logger.warn("... done with constructing model container.");
+
+        logger.warn("Constructing silo model ...");
+        SiloModel model = new SiloModel(siloConfig, dataContainer, modelContainer);
+
+        model.addResultMonitor( new DefaultResultsMonitor(dataContainer, siloConfig) );
+        model.addResultMonitor( new MultiFileResultsMonitor(dataContainer, siloConfig) );
+        model.addResultMonitor( new HouseholdSatisfactionMonitor(dataContainer, siloConfig, modelContainer) );
+        logger.warn("... done with constructing silo model.");
+
+        logger.warn("Running silo model ...");
+        model.runModel();
+        logger.warn("Finished SILO.");
+    }
+}
diff --git a/useCases/fabiland/src/test/java/run/RunFabilandModelCombinationTests.java b/useCases/fabiland/src/test/java/run/RunFabilandModelCombinationTests.java
new file mode 100644
index 000000000..1fe189b8b
--- /dev/null
+++ b/useCases/fabiland/src/test/java/run/RunFabilandModelCombinationTests.java
@@ -0,0 +1,789 @@
+package run;
+
+import de.tum.bgu.msm.SiloModel;
+import de.tum.bgu.msm.container.DataContainer;
+import de.tum.bgu.msm.container.ModelContainer;
+import de.tum.bgu.msm.data.dwelling.DwellingFactory;
+import de.tum.bgu.msm.data.household.Household;
+import de.tum.bgu.msm.data.household.HouseholdFactory;
+import de.tum.bgu.msm.data.person.Person;
+import de.tum.bgu.msm.data.person.PersonFactory;
+import de.tum.bgu.msm.events.impls.household.MoveEvent;
+import de.tum.bgu.msm.io.output.DefaultResultsMonitor;
+import de.tum.bgu.msm.io.output.HouseholdSatisfactionMonitor;
+import de.tum.bgu.msm.io.output.MultiFileResultsMonitor;
+import de.tum.bgu.msm.matsim.*;
+import de.tum.bgu.msm.models.autoOwnership.CreateCarOwnershipModel;
+import de.tum.bgu.msm.models.demography.birth.BirthModel;
+import de.tum.bgu.msm.models.demography.birth.BirthModelImpl;
+import de.tum.bgu.msm.models.demography.birth.DefaultBirthStrategy;
+import de.tum.bgu.msm.models.demography.birthday.BirthdayModel;
+import de.tum.bgu.msm.models.demography.birthday.BirthdayModelImpl;
+import de.tum.bgu.msm.models.demography.death.DeathModel;
+import de.tum.bgu.msm.models.demography.death.DeathModelImpl;
+import de.tum.bgu.msm.models.demography.death.DefaultDeathStrategy;
+import de.tum.bgu.msm.models.demography.divorce.DefaultDivorceStrategy;
+import de.tum.bgu.msm.models.demography.divorce.DivorceModel;
+import de.tum.bgu.msm.models.demography.divorce.DivorceModelImpl;
+import de.tum.bgu.msm.models.demography.driversLicense.DefaultDriversLicenseStrategy;
+import de.tum.bgu.msm.models.demography.driversLicense.DriversLicenseModel;
+import de.tum.bgu.msm.models.demography.driversLicense.DriversLicenseModelImpl;
+import de.tum.bgu.msm.models.demography.education.EducationModel;
+import de.tum.bgu.msm.models.demography.education.EducationModelImpl;
+import de.tum.bgu.msm.models.demography.employment.EmploymentModel;
+import de.tum.bgu.msm.models.demography.employment.EmploymentModelImpl;
+import de.tum.bgu.msm.models.demography.leaveParentalHousehold.DefaultLeaveParentalHouseholdStrategy;
+import de.tum.bgu.msm.models.demography.leaveParentalHousehold.LeaveParentHhModel;
+import de.tum.bgu.msm.models.demography.leaveParentalHousehold.LeaveParentHhModelImpl;
+import de.tum.bgu.msm.models.demography.marriage.DefaultMarriageStrategy;
+import de.tum.bgu.msm.models.demography.marriage.MarriageModel;
+import de.tum.bgu.msm.models.demography.marriage.MarriageModelImpl;
+import de.tum.bgu.msm.models.jobmography.JobMarketUpdate;
+import de.tum.bgu.msm.models.jobmography.JobMarketUpdateImpl;
+import de.tum.bgu.msm.models.modeChoice.CommuteModeChoice;
+import de.tum.bgu.msm.models.realEstate.construction.*;
+import de.tum.bgu.msm.models.realEstate.demolition.DefaultDemolitionStrategy;
+import de.tum.bgu.msm.models.realEstate.demolition.DemolitionModel;
+import de.tum.bgu.msm.models.realEstate.demolition.DemolitionModelImpl;
+import de.tum.bgu.msm.models.realEstate.pricing.DefaultPricingStrategy;
+import de.tum.bgu.msm.models.realEstate.pricing.PricingModel;
+import de.tum.bgu.msm.models.realEstate.pricing.PricingModelImpl;
+import de.tum.bgu.msm.models.realEstate.renovation.DefaultRenovationStrategy;
+import de.tum.bgu.msm.models.realEstate.renovation.RenovationModel;
+import de.tum.bgu.msm.models.realEstate.renovation.RenovationModelImpl;
+import de.tum.bgu.msm.models.relocation.migration.InOutMigration;
+import de.tum.bgu.msm.models.relocation.migration.InOutMigrationImpl;
+import de.tum.bgu.msm.models.relocation.moves.*;
+import de.tum.bgu.msm.models.transportModel.TransportModel;
+import de.tum.bgu.msm.properties.Properties;
+import de.tum.bgu.msm.utils.SiloUtil;
+import models.FabilandConstructionLocationStrategy;
+import org.apache.logging.log4j.LogManager;
+import org.apache.logging.log4j.Logger;
+import org.checkerframework.checker.nullness.qual.NonNull;
+import org.junit.jupiter.api.AfterEach;
+import org.junit.jupiter.api.Assertions;
+import org.junit.jupiter.api.BeforeEach;
+import org.junit.jupiter.api.Test;
+import org.junit.jupiter.api.extension.RegisterExtension;
+import org.matsim.api.core.v01.Scenario;
+import org.matsim.core.config.Config;
+import org.matsim.core.config.ConfigUtils;
+import org.matsim.core.config.groups.RoutingConfigGroup;
+import org.matsim.core.controler.OutputDirectoryHierarchy;
+import org.matsim.core.scenario.MutableScenario;
+import org.matsim.core.scenario.ScenarioUtils;
+import org.matsim.testcases.MatsimTestUtils;
+import tech.tablesaw.api.DoubleColumn;
+import tech.tablesaw.api.Row;
+import tech.tablesaw.api.StringColumn;
+import tech.tablesaw.api.Table;
+
+import java.util.Arrays;
+import java.util.Collection;
+import java.util.Collections;
+import java.util.Random;
+
+//@ExtendWith(MatsimTestUtils.class)
+public class RunFabilandModelCombinationTests {
+
+    //region Define global variables
+    private static final Logger log = LogManager.getLogger(RunFabilandModelCombinationTests.class);
+    private static MovesModel movesModel;
+    private static BirthModel birthModel;
+    private static PersonFactory ppFactory;
+    private static HouseholdFactory hhFactory;
+    private static DwellingFactory ddFactory;
+    private static BirthdayModel birthdayModel;
+    private static DeathModel deathModel;
+    private static EducationModel educationModel;
+    private static CreateCarOwnershipModel carOwnershipModel;
+    private static DivorceModel divorceModel;
+    private static DriversLicenseModel driversLicenseModel;
+    private static LeaveParentHhModel leaveParentsModel;
+    private static JobMarketUpdate jobMarketUpdateModel;
+    private static EmploymentModel employmentModel;
+    private static MarriageModel marriageModel;
+    private static PricingModel pricingModel;
+    private static TransportModel transportModel;
+    private static ConstructionOverwrite constructionOverwrite;
+    private static InOutMigration inOutMigration;
+    private static ConstructionModel constructionModel;
+    private static DemolitionModel demolitionModel;
+    private static RenovationModel renovationModel;
+
+    private static DataContainer dataContainer = null;
+    private static Properties siloConfig = null;
+    private static Config matsimConfig = null;
+    //endregion
+
+    @RegisterExtension
+	public MatsimTestUtils utils = new MatsimTestUtils();
+    private String popfiles;
+    private double expectedAgeDiff;
+
+    @AfterEach
+    void tearDown() {
+        Properties.reset();
+    }
+
+    @BeforeEach
+    public void setup() {
+
+
+        // SILO configuration
+        siloConfig = SiloUtil.siloInitialization("./scenario/test2.properties");
+        popfiles = siloConfig.main.baseDirectory + "scenOutput/" + siloConfig.main.scenarioName + "/microData/";
+
+        expectedAgeDiff = siloConfig.main.endYear - siloConfig.main.startYear - 1;
+
+        // MATSim configuration
+        matsimConfig = ConfigUtils.loadConfig( "./scenario/config_cap30_1-l_nes_smc.xml");
+        matsimConfig.controller().setLastIteration(1);
+
+        // The following is obviously just a dirty quickfix until access/egress is default in MATSim
+        if (siloConfig.transportModel.includeAccessEgress) {
+            matsimConfig.routing().setAccessEgressType(RoutingConfigGroup.AccessEgressType.accessEgressModeToLink); // in matsim-13-w37
+        }
+
+        matsimConfig.controller().setOverwriteFileSetting(OutputDirectoryHierarchy.OverwriteFileSetting.overwriteExistingFiles);
+
+        // data containers
+        dataContainer = DataBuilderFabiland.buildDataContainer(siloConfig, matsimConfig);
+        DataBuilderFabiland.readInput(siloConfig, dataContainer);
+
+        //sets up all default models
+        setupModels(dataContainer, siloConfig, matsimConfig);
+    }
+
+    /**
+     * BASE TEST: This test runs FABILAND with all models activated. This reflects RunFabiland or RunFabilandTest
+     * we test whether ages change, people are born, people die, people change households
+     */
+	@Test
+	public void testAllModels() {
+
+    ModelContainer modelContainer = new ModelContainer(
+            birthModel, birthdayModel,
+            deathModel, marriageModel,
+            divorceModel, driversLicenseModel,
+            educationModel, employmentModel,
+            leaveParentsModel, jobMarketUpdateModel,
+            constructionModel, demolitionModel, pricingModel, renovationModel,
+            constructionOverwrite, inOutMigration, movesModel, transportModel);
+
+
+//        ModelContainer modelContainer = ModelBuilderFabiland.getModelContainer(dataContainer, siloConfig, matsimConfig);
+        SiloModel model = new SiloModel(siloConfig, dataContainer, modelContainer);
+
+        model.runModel();
+
+        // PREPARE PERSONS TABLE FOR TESTS
+        // join first and last year
+        Table personsJoined = readAndJoinPp();
+
+        String ageStartCol = "age";
+        String ageEndCol = "T2.age"; // <-- verify this against the logged column names above
+
+        // people present in first and last years
+        Table personsAliveContinuously = personsJoined.dropWhere(
+                personsJoined.column(ageStartCol).isMissing()
+                        .or(personsJoined.column(ageEndCol).isMissing())
+        );
+
+        // column indicating whether people changed households or stayed put
+        StringColumn hhChanged = getHhChangedColumn(personsAliveContinuously);
+        personsAliveContinuously.addColumns(hhChanged);
+
+        // TEST NUMBER BORN > 0
+        // people born = present in pp_10 but not pp_0 -> age.x is missing
+        long numBorn = personsJoined.where(personsJoined.column(ageStartCol).isMissing()).rowCount();
+        Assertions.assertTrue(numBorn > 0);
+        // TEST NUMBER DIED > 0
+        // people who died = present in pp_0 but not pp_10 -> age.y is missing
+        long numDied = personsJoined.where(personsJoined.column(ageEndCol).isMissing()).rowCount();
+        Assertions.assertTrue(numDied > 0);
+
+
+        // TEST AGE DIFFERENCE = 9 FOR ALL CONTINUOUSLY ALIVE AGENTS
+        double[] uniqueDiffs = findAgeDifferences(personsAliveContinuously, ageStartCol, ageEndCol);
+
+        Assertions.assertEquals(1, uniqueDiffs.length,
+                "Age difference should be identical for every person, but found: " + Arrays.toString(uniqueDiffs));
+
+        Assertions.assertEquals(expectedAgeDiff, uniqueDiffs[0], 1e-6,
+                "Expected age difference of 9 years between pp_0 and pp_10");
+
+        // TEST THAT PEOPLE MOVED
+        int movedCnt =  hhChanged.isEqualTo("moved household").size();
+        Assertions.assertTrue(movedCnt > 0);
+    }
+
+
+    /**
+     * Fabiland without construction, demolition, renovation, models.
+     * we test whether ages change, people are born, people die, people change households (same as previous)
+     */
+    @Test
+    public void testModelsLarge() {
+
+        ModelContainer modelContainer = new ModelContainer(
+                birthModel, birthdayModel,
+                deathModel, marriageModel,
+                divorceModel, driversLicenseModel,
+                educationModel, employmentModel,
+                leaveParentsModel, jobMarketUpdateModel,
+                null, null, pricingModel, null,
+                null, inOutMigration, movesModel, transportModel);
+
+
+        SiloModel model = new SiloModel(siloConfig, dataContainer, modelContainer);
+
+        model.runModel();
+
+
+        // PREPARE PERSONS TABLE FOR TESTS
+        // join first and last year
+        Table personsJoined = readAndJoinPp();
+
+        String ageStartCol = "age";
+        String ageEndCol = "T2.age"; // <-- verify this against the logged column names above
+
+        // people present in first and last years
+        Table personsAliveContinuously = personsJoined.dropWhere(
+                personsJoined.column(ageStartCol).isMissing()
+                        .or(personsJoined.column(ageEndCol).isMissing())
+        );
+
+        // column indicating whether people changed households or stayed put
+        StringColumn hhChanged = getHhChangedColumn(personsAliveContinuously);
+        personsAliveContinuously.addColumns(hhChanged);
+
+        // TEST NUMBER BORN > 0
+        // people born = present in pp_10 but not pp_0 -> age.x is missing
+        long numBorn = personsJoined.where(personsJoined.column(ageStartCol).isMissing()).rowCount();
+        Assertions.assertTrue(numBorn > 0);
+        // TEST NUMBER DIED > 0
+        // people who died = present in pp_0 but not pp_10 -> age.y is missing
+        long numDied = personsJoined.where(personsJoined.column(ageEndCol).isMissing()).rowCount();
+        Assertions.assertTrue(numDied > 0);
+
+
+        // TEST AGE DIFFERENCE = 9 FOR ALL CONTINUOUSLY ALIVE AGENTS
+        double[] uniqueDiffs = findAgeDifferences(personsAliveContinuously, ageStartCol, ageEndCol);
+
+        Assertions.assertEquals(1, uniqueDiffs.length,
+                "Age difference should be identical for every person, but found: " + Arrays.toString(uniqueDiffs));
+
+        Assertions.assertEquals(expectedAgeDiff, uniqueDiffs[0], 1e-6,
+                "Expected age difference of 9 years between pp_0 and pp_10");
+
+        // TEST THAT PEOPLE MOVED
+        int movedCnt =  hhChanged.isEqualTo("moved household").size();
+        Assertions.assertTrue(movedCnt > 0);
+    }
+
+
+    /**
+     * Additionally to models removed in previous test, we also take out
+     * driversLicense, education, employment, and jobMarket models
+     * we test whether ages change, people are born, people die, people change households (same as previous)
+     */
+    @Test
+    public void testModelsMedium() {
+
+        ModelContainer modelContainer = new ModelContainer(
+                birthModel, birthdayModel,
+                deathModel, marriageModel,
+                divorceModel, null,
+                null, null,
+                leaveParentsModel, null,
+                null, null, pricingModel, null,
+                null, inOutMigration, movesModel, transportModel);
+
+
+//        ModelContainer modelContainer = ModelBuilderFabiland.getModelContainer(dataContainer, siloConfig, matsimConfig);
+        SiloModel model = new SiloModel(siloConfig, dataContainer, modelContainer);
+
+        model.runModel();
+
+
+        // PREPARE PERSONS TABLE FOR TESTS
+        // join first and last year
+        Table personsJoined = readAndJoinPp();
+
+        String ageStartCol = "age";
+        String ageEndCol = "T2.age"; // <-- verify this against the logged column names above
+
+        // people present in first and last years
+        Table personsAliveContinuously = personsJoined.dropWhere(
+                personsJoined.column(ageStartCol).isMissing()
+                        .or(personsJoined.column(ageEndCol).isMissing())
+        );
+
+        // column indicating whether people changed households or stayed put
+        StringColumn hhChanged = getHhChangedColumn(personsAliveContinuously);
+        personsAliveContinuously.addColumns(hhChanged);
+
+        // TEST NUMBER BORN > 0
+        // people born = present in pp_10 but not pp_0 -> age.x is missing
+        long numBorn = personsJoined.where(personsJoined.column(ageStartCol).isMissing()).rowCount();
+        Assertions.assertTrue(numBorn > 0);
+        // TEST NUMBER DIED > 0
+        // people who died = present in pp_0 but not pp_10 -> age.y is missing
+        long numDied = personsJoined.where(personsJoined.column(ageEndCol).isMissing()).rowCount();
+        Assertions.assertTrue(numDied > 0);
+
+
+        // TEST AGE DIFFERENCE = 9 FOR ALL CONTINUOUSLY ALIVE AGENTS
+        double[] uniqueDiffs = findAgeDifferences(personsAliveContinuously, ageStartCol, ageEndCol);
+
+        Assertions.assertEquals(1, uniqueDiffs.length,
+                "Age difference should be identical for every person, but found: " + Arrays.toString(uniqueDiffs));
+
+        Assertions.assertEquals(expectedAgeDiff, uniqueDiffs[0], 1e-6,
+                "Expected age difference of 9 years between pp_0 and pp_10");
+
+        // TEST THAT PEOPLE MOVED
+        int movedCnt =  hhChanged.isEqualTo("moved household").size();
+        Assertions.assertTrue(movedCnt > 0);
+    }
+
+    /**
+     * Additionally to models removed in previous test, we also take out
+     * marriage, divorce, & leaveParentsHH
+     * we test whether ages change, people are born, people die, people change households (same as previous)
+     * PEOPLE DO NOT MOVE, I guess since others models are turned off
+     */
+    @Test
+    public void testModelsSmall() {
+
+        ModelContainer modelContainer = new ModelContainer(
+                birthModel, birthdayModel,
+                deathModel, null,
+                null, null,
+                null, null,
+                null, null,
+                null, null, pricingModel, null,
+                null, inOutMigration, movesModel, transportModel);
+
+
+        SiloModel model = new SiloModel(siloConfig, dataContainer, modelContainer);
+        model.runModel();
+
+
+        // PREPARE PERSONS TABLE FOR TESTS
+        // join first and last year
+        Table personsJoined = readAndJoinPp();
+
+        String ageStartCol = "age";
+        String ageEndCol = "T2.age"; // <-- verify this against the logged column names above
+
+        // people present in first and last years
+        Table personsAliveContinuously = personsJoined.dropWhere(
+                personsJoined.column(ageStartCol).isMissing()
+                        .or(personsJoined.column(ageEndCol).isMissing())
+        );
+
+        // column indicating whether people changed households or stayed put
+        StringColumn hhChanged = getHhChangedColumn(personsAliveContinuously);
+        personsAliveContinuously.addColumns(hhChanged);
+
+        // TEST NUMBER BORN > 0
+        // people born = present in pp_10 but not pp_0 -> age.x is missing
+        long numBorn = personsJoined.where(personsJoined.column(ageStartCol).isMissing()).rowCount();
+        Assertions.assertTrue(numBorn > 0);
+        // TEST NUMBER DIED > 0
+        // people who died = present in pp_0 but not pp_10 -> age.y is missing
+        long numDied = personsJoined.where(personsJoined.column(ageEndCol).isMissing()).rowCount();
+        Assertions.assertTrue(numDied > 0);
+
+
+        // TEST AGE DIFFERENCE = 9 FOR ALL CONTINUOUSLY ALIVE AGENTS
+        double[] uniqueDiffs = findAgeDifferences(personsAliveContinuously, ageStartCol, ageEndCol);
+
+        Assertions.assertEquals(1, uniqueDiffs.length,
+                "Age difference should be identical for every person, but found: " + Arrays.toString(uniqueDiffs));
+
+        Assertions.assertEquals(expectedAgeDiff, uniqueDiffs[0], 1e-6,
+                "Expected age difference of 9 years between pp_0 and pp_10");
+
+        // TEST THAT PEOPLE MOVED
+        int movedCnt =  hhChanged.isEqualTo("moved household").size();
+        Assertions.assertEquals(0, movedCnt);
+    }
+
+    /**
+     * Additionally to models removed in previous test, we also take out
+     * birth, death, pricing, migration
+     * ALl we are left with is birthday and move and transport
+     * we test whether ages change, people are NOT born, people DO NOT die, people DO NOT change households
+     */
+    @Test
+    public void testModelsOnlyBirthdays() {
+
+        ModelContainer modelContainer = new ModelContainer(
+                null, birthdayModel,
+                null, null,
+                null, null,
+                null, null,
+                null, null,
+                null, null, null, null,
+                null, null, null, transportModel);
+
+
+        SiloModel model = new SiloModel(siloConfig, dataContainer, modelContainer);
+
+        model.runModel();
+
+
+        // PREPARE PERSONS TABLE FOR TESTS
+        // join first and last year
+        Table personsJoined = readAndJoinPp();
+
+        String ageStartCol = "age";
+        String ageEndCol = "T2.age"; // <-- verify this against the logged column names above
+
+        // people present in first and last years
+        Table personsAliveContinuously = personsJoined.dropWhere(
+                personsJoined.column(ageStartCol).isMissing()
+                        .or(personsJoined.column(ageEndCol).isMissing())
+        );
+
+        // column indicating whether people changed households or stayed put
+        StringColumn hhChanged = getHhChangedColumn(personsAliveContinuously);
+        personsAliveContinuously.addColumns(hhChanged);
+
+        // TEST NUMBER BORN = 0
+        // people born = present in pp_10 but not pp_0 -> age.x is missing
+        long numBorn = personsJoined.where(personsJoined.column(ageStartCol).isMissing()).rowCount();
+        Assertions.assertEquals( 0, numBorn);
+        // TEST NUMBER DIED = 0
+        // people who died = present in pp_0 but not pp_10 -> age.y is missing
+        long numDied = personsJoined.where(personsJoined.column(ageEndCol).isMissing()).rowCount();
+        Assertions.assertEquals( 0, numDied);
+
+
+        // TEST AGE DIFFERENCE = 9 FOR ALL CONTINUOUSLY ALIVE AGENTS
+        double[] uniqueDiffs = findAgeDifferences(personsAliveContinuously, ageStartCol, ageEndCol);
+
+        Assertions.assertEquals(1, uniqueDiffs.length,
+                "Age difference should be identical for every person, but found: " + Arrays.toString(uniqueDiffs));
+
+        Assertions.assertEquals(expectedAgeDiff, uniqueDiffs[0], 1e-6,
+                "Expected age difference of 9 years between pp_0 and pp_10");
+
+        // TEST THAT PEOPLE MOVED
+        int movedCnt =  hhChanged.isEqualTo("moved household").size();
+        Assertions.assertEquals(0, movedCnt);
+    }
+
+
+    /**
+     * Tests custom set of models: birth, birthday, death, marriage*, drivers license, education, employment, transport
+     * moves model is technically present, but is overridden to do nothing
+     * marriage model is overriden such that newly wed couple doesn't move in together
+     * we test whether ages change, people are born, people die, people DO NOT change households
+     */
+
+    @Test
+    public void testModelsKai() {
+        marriageModel = new MarriageModelImpl(dataContainer, movesModel, inOutMigration,
+                carOwnershipModel, hhFactory, siloConfig, new DefaultMarriageStrategy(), SiloUtil.provideNewRandom()) {
+            @Override
+            protected boolean moveTogether(Person person1, Person person2, Household moveTo) {
+                return true;
+            }
+        };
+
+        // breaks if we take moves out completely
+        // taken from ModelBuilderFabilandSimplified.
+        movesModel = new MovesModel(){
+            @Override public int searchForNewDwelling( Household household ){
+                return -1; // means no dwelling was found
+            }
+            @Override public void moveHousehold( Household hh, int idOldDD, int idNewDD ){
+                // do nothing
+            }
+            @Override public Collection<MoveEvent> getEventsForCurrentYear(int year ){
+                return Collections.emptyList();
+            }
+            @Override public boolean handleEvent( MoveEvent event ){
+                // this should not happen.  Maybe test?
+                // If this is needed, one cold use MovesModelImpl as a delegate and then go from there.
+                return false;
+            }
+            @Override public void setup(){
+                // do nothing
+            }
+            @Override public void prepareYear( int year ){
+                // do nothing
+            }
+            @Override public void endYear( int year ){
+                // do nothing
+            }
+            @Override public void endSimulation(){
+                // do nothing
+            }
+        };
+
+        ModelContainer modelContainer = new ModelContainer(
+                birthModel, birthdayModel,
+                deathModel, marriageModel,
+                null, driversLicenseModel,
+                educationModel, employmentModel,
+                null, null,
+                null, null, null, null,
+                null, null, movesModel, transportModel);
+
+        SiloModel model = new SiloModel(siloConfig, dataContainer, modelContainer);
+        model.runModel();
+
+
+        // PREPARE PERSONS TABLE FOR TESTS
+        // join first and last year
+        Table personsJoined = readAndJoinPp();
+
+        String ageStartCol = "age";
+        String ageEndCol = "T2.age"; // <-- verify this against the logged column names above
+
+        // people present in first and last years
+        Table personsAliveContinuously = personsJoined.dropWhere(
+                personsJoined.column(ageStartCol).isMissing()
+                        .or(personsJoined.column(ageEndCol).isMissing())
+        );
+
+        // column indicating whether people changed households or stayed put
+        StringColumn hhChanged = getHhChangedColumn(personsAliveContinuously);
+        personsAliveContinuously.addColumns(hhChanged);
+
+        // TEST NUMBER BORN > 0
+        // people born = present in pp_10 but not pp_0 -> age.x is missing
+        long numBorn = personsJoined.where(personsJoined.column(ageStartCol).isMissing()).rowCount();
+        Assertions.assertTrue(numBorn > 0);
+        // TEST NUMBER DIED > 0
+        // people who died = present in pp_0 but not pp_10 -> age.y is missing
+        long numDied = personsJoined.where(personsJoined.column(ageEndCol).isMissing()).rowCount();
+        Assertions.assertTrue(numDied > 0);
+
+
+        // TEST AGE DIFFERENCE = 9 FOR ALL CONTINUOUSLY ALIVE AGENTS
+        double[] uniqueDiffs = findAgeDifferences(personsAliveContinuously, ageStartCol, ageEndCol);
+
+        Assertions.assertEquals(1, uniqueDiffs.length,
+                "Age difference should be identical for every person, but found: " + Arrays.toString(uniqueDiffs));
+
+        Assertions.assertEquals(expectedAgeDiff, uniqueDiffs[0], 1e-6,
+                "Expected age difference of 9 years between pp_0 and pp_10");
+
+        // TEST THAT PEOPLE MOVED = 0
+        int movedCnt =  hhChanged.isEqualTo("moved household").size();
+        Assertions.assertEquals(0, movedCnt);
+    }
+
+
+    private static @NonNull StringColumn getHhChangedColumn(Table personsAliveContinuously) {
+        StringColumn hhChanged = StringColumn.create("household_changed", personsAliveContinuously.rowCount());
+
+        for (Row row : personsAliveContinuously) {
+            Integer hhx = row.getInt("hhid");
+            Integer hhy = row.getInt("T2.hhid");
+
+            String category = !hhx.equals(hhy) ? "moved household" : "stayed put";
+            hhChanged.set(row.getRowNumber(), category);
+        }
+        return hhChanged;
+    }
+
+    private static double[] findAgeDifferences(Table personsAliveContinuously, String ageStartCol, String ageEndCol) {
+        DoubleColumn ageStart = personsAliveContinuously.numberColumn(ageStartCol).asDoubleColumn();
+        DoubleColumn ageEnd = personsAliveContinuously.numberColumn(ageEndCol).asDoubleColumn();
+        DoubleColumn ageDiff = ageEnd.subtract(ageStart).setName("age_diff");
+        personsAliveContinuously.addColumns(ageDiff);
+
+        double[] uniqueDiffs = ageDiff.unique().asDoubleArray();
+        return uniqueDiffs;
+    }
+
+    private Table readAndJoinPp() {
+        int lastYear = siloConfig.main.endYear - siloConfig.main.startYear;
+        Table ppStart = Table.read().csv(popfiles + "pp_0.csv");
+        Table ppEnd = Table.read().csv(popfiles + "pp_" + lastYear + ".csv");
+
+        // full outer join on "id" — Tablesaw will rename colliding columns from the
+        // right-hand table (check joined.columnNames() to confirm the actual suffix,
+        // it's typically something like "age_2" but can vary by version)
+        Table joined = ppStart.joinOn("id").fullOuter(ppEnd, true, false, "id");
+        return joined;
+    }
+
+    public void setupModels(DataContainer dataContainer, Properties properties, Config config) {
+
+        try {
+        ppFactory = dataContainer.getHouseholdDataManager().getPersonFactory();
+        hhFactory = dataContainer.getHouseholdDataManager().getHouseholdFactory();
+        ddFactory = dataContainer.getRealEstateDataManager().getDwellingFactory();
+
+        birthModel = new BirthModelImpl(dataContainer, ppFactory, properties, new DefaultBirthStrategy(), SiloUtil.provideNewRandom());
+
+        birthdayModel = new BirthdayModelImpl(dataContainer, properties, SiloUtil.provideNewRandom());
+
+        deathModel = new DeathModelImpl(dataContainer, properties, new DefaultDeathStrategy(), SiloUtil.provideNewRandom());
+
+        {
+
+            final DwellingUtilityStrategy dwellingUtilityStrategy = new DwellingUtilityStrategyImpl();
+            // (This is something like
+            // [alpha * sizeUtility + beta * autoAccessibilityUtility + gamma * transitAccessibilityUtility + (1.0 - alpha - beta - gamma) * qualityUtility]^delta   * priceUtl^eps * workDistanceUtl^{1-delta-eps}
+            // That is, some kind of Cobbs Douglas function with contributions
+            // * workDistanceUtil
+            // * priceUtl
+            // * a weighted sum of sizeUtl, autoAccUtl, transitAccUtl, qualityUtl.
+            // All of the params depend on income and hh size, and are given for corresponding categories.
+            // )
+
+            final DwellingProbabilityStrategy dwellingProbabilityStrategy = new DefaultDwellingProbabilityStrategy();
+            // (is just exp(beta*util); should be called "weight" instead of "probability" since it is not normalized.  kai, apr'26)
+
+            final RegionUtilityStrategy regionUtilityStrategy = new RegionUtilityStrategyImpl();
+            // (This is something like (1 - alpha) * price + alpha * accessibility, with alpha depending on the income category.)
+
+            final RegionProbabilityStrategy regionProbabilityStrategy = new RegionProbabilityStrategyImpl();
+            //( same as DefaultDwellingProbabilityStrategy, see above. kai, apr'26)
+
+            final CommuteModeChoice commuteModeChoice1 = new SimpleMatsimCommuteModeChoice( dataContainer, properties, SiloUtil.provideNewRandom() );
+            // (there is a comment in SimpleCommuteModeChoiceHousingStrategyImpl that the constructor should actually work w/o providing the CommuteChoiceModel.  It then provides
+            // CommuteModeChoice internally, as SimpleCommuteModeChoice.  That model looks similar; presumably some of the lookups (e.g. travel time) used to have different arguments.
+            // kai, apr'26)
+            // (I think that the existing implementations do the following: compute the logit probas for car and pt if car is an option; then go from HH member with largest car
+            // proba down and select car with logit proba as long as another car is available in the HH.)
+            // (yy this implies, as long as nothing else comes in, that moves consider new residencies under the assumption that they will not change the number of vehicles)
+
+            final HousingStrategy housingStrategy = new SimpleCommuteModeChoiceHousingStrategyImpl( dataContainer,
+                    properties,
+                    dataContainer.getTravelTimes(),
+                    dwellingUtilityStrategy,
+                    dwellingProbabilityStrategy,
+                    regionUtilityStrategy,
+                    regionProbabilityStrategy,
+                    commuteModeChoice1
+            );
+            /// (If I see this right, this is computing the necessary inputs to {@link DwellingUtilityStrategy} using all the other dependencies, and then computing the dwelling
+            /// probabilities using the {@link DwellingUtilityStrategy}.
+
+            final MovesStrategy movesStrategy = new DefaultMovesStrategy();
+            // (This says that the moving proba is  1 - 1/(1+0.03 * Math.exp(10*(householdSatisfaction - currentDwellingUtility)))
+
+            movesModel = new MovesModelImpl( dataContainer, properties, movesStrategy, housingStrategy, SiloUtil.provideNewRandom() );
+        }
+        // (Overall, I think that life events are not explicitly triggering a move.  They will, however, shift the respective utilities, and so a move becomes more probable.)
+
+        // (For our own issue, which (only) is to "age" the population, this seems like a lot of overhead.  On the other hand, it is not immediately clear how else this should be
+        // resolved; location choice needs to be something like logit based on location utility.  What might end up being a bit of a problem is (commute) mode choice ... we would
+        // need to say that maybe the commute mode choice to make a residence decision is based on what is defined above, but the final mode choice may be different once people
+        // optimize into their new environment. kai, apr'26)
+
+        // (Also, the MatsimScenarioAssembler takes the commute mode from the Silo mode choice model (although I am not sure that it will give the same results are for housing
+        // choice, since the random numbers are different--????).  Evidently, we can override this in MATSim.  Clearly, this will be two models pull into different directions.
+        // Unfortunately, I fail to see how we calibrade mode choice on the MATSim side if we rely on upstream mode choice.  Maybe possible in principle, but will not be sensitive
+        // to (the details of) many transport policies.  kai, apr'26)
+
+        // (As we have known for some time now, this is two competing modelling paradigms: mode (and time) choice as part of the DTA vs mode (and time?) choice as part of the
+        // upstream model.  Rolf seems to have been a proponent of the former.  However, from a MATSim side and MATSim research side, this is not something we can continue since we
+        // have experience with mode choice in MATSim but not in upstream models, and this is also a path that we do not want to take (would need to be done by someone else).
+        // We should consider to re-code the commute mode choice for our purposes, i.e. always allow for car and pt during housing search, using income-dependent utilities and
+        // fixed costs for cars (and pt), etc., and then sort out the actual choice later.  However, this may end up not so different from the current approach, so why not leave
+        // the current approach in place, but change as follows: Only if a silo person has moved residency or job, then we use the mode choice from silo as an initial suggestion;
+        // otherwise we keep the plan (including mode) from previous iterations. -- Such an approach would make a lot of sense anyways; it was already encoded as "hot start" in
+        // both matsim-urbansim couplings (by KN and by TN); and it would be far enough downstream of SILO, i.e. only in the adapter class. --  A bit disappointing that 20 year
+        // later we are not any further here. :-(   kai, apr'26)
+
+        carOwnershipModel = new FabilandCarOwnership();
+
+        divorceModel = new DivorceModelImpl( dataContainer, movesModel, carOwnershipModel, hhFactory, properties, new DefaultDivorceStrategy(), SiloUtil.provideNewRandom());
+
+        driversLicenseModel = new DriversLicenseModelImpl(dataContainer, properties, new DefaultDriversLicenseStrategy(), SiloUtil.provideNewRandom());
+
+        educationModel = new EducationModelImpl(dataContainer, properties, SiloUtil.provideNewRandom());
+
+        employmentModel = new EmploymentModelImpl(dataContainer, properties, SiloUtil.provideNewRandom());
+
+        leaveParentsModel = new LeaveParentHhModelImpl(dataContainer, movesModel,
+                carOwnershipModel, hhFactory, properties, new DefaultLeaveParentalHouseholdStrategy(), SiloUtil.provideNewRandom());
+
+        jobMarketUpdateModel = new JobMarketUpdateImpl(dataContainer, properties, SiloUtil.provideNewRandom());
+
+        constructionModel = new ConstructionModelImpl(dataContainer, ddFactory,
+                properties, new FabilandConstructionLocationStrategy(), new DefaultConstructionDemandStrategy(), SiloUtil.provideNewRandom());
+
+        pricingModel = new PricingModelImpl(dataContainer, properties, new DefaultPricingStrategy(), SiloUtil.provideNewRandom());
+
+        renovationModel = new RenovationModelImpl(dataContainer, properties, new DefaultRenovationStrategy(), SiloUtil.provideNewRandom());
+
+        constructionOverwrite = new ConstructionOverwriteImpl(dataContainer, ddFactory, properties, SiloUtil.provideNewRandom());
+
+        inOutMigration = new InOutMigrationImpl(dataContainer, employmentModel, movesModel,
+                carOwnershipModel, driversLicenseModel, properties, SiloUtil.provideNewRandom());
+
+        demolitionModel = new DemolitionModelImpl(dataContainer, movesModel,
+                inOutMigration, properties, new DefaultDemolitionStrategy(), SiloUtil.provideNewRandom());
+
+        marriageModel = new MarriageModelImpl(dataContainer, movesModel, inOutMigration,
+                carOwnershipModel, hhFactory, properties, new DefaultMarriageStrategy(), SiloUtil.provideNewRandom());
+
+        MatsimScenarioAssembler scenarioAssembler;
+
+        MatsimData matsimData = null;
+        if (config != null) {
+            final Scenario scenario = ScenarioUtils.loadScenario(config);
+//			matsimData = new MatsimData(config, properties, ZoneConnectorManagerImpl.ZoneConnectorMethod.WEIGHTED_BY_POPULATION, dataContainer, scenario.getNetwork(), scenario.getTransitSchedule());
+            matsimData = new MatsimData(properties, ZoneConnectorManager.ZoneConnectorMethod.WEIGHTED_BY_POPULATION, dataContainer, (MutableScenario) scenario );
+            // yyyyyy the scenario is now generated in the MatsimData class, and in consequence the above constructor could/should be removed again
+        }
+        switch (properties.transportModel.transportModelIdentifier) {
+            case MATSIM:
+//                SimpleCommuteModeChoice commuteModeChoice = new SimpleCommuteModeChoice(dataContainer, properties, SiloUtil.provideNewRandom());
+                SimpleMatsimCommuteModeChoice commuteModeChoice = new SimpleMatsimCommuteModeChoice(dataContainer, properties, SiloUtil.provideNewRandom());
+                // (yyyy is also instantiated above.  why not re-use?)
+                // --> fails the regression test.  Which is no wonder, since it changes the random number sequence.
+
+                scenarioAssembler = new SimpleCommuteModeChoiceMatsimScenarioAssembler(dataContainer, properties, commuteModeChoice, SimpleCommuteModeChoiceMatsimScenarioAssembler.HandlingOfRandomness.localInstanceFromMatsimWithAlwaysSameSeed);
+                // yyyyyy the above needs to be re-coded from a MATSim perspective ... i.e. (1) select the same agents, and (2) take the mode "suggestion" from silo only if home or
+                // job location has changed, otherwise keep existing plan with existing mode.  kai, apr'26
+
+                transportModel = new MatsimTransportModel(dataContainer, config, properties, scenarioAssembler, matsimData);
+                break;
+            case NONE:
+            default:
+                transportModel = null;
+        }
+        } catch (Exception ee) {
+            log.fatal("there was an exception: \n" + ee);
+
+            ee.printStackTrace();
+
+            // if one catches an exception, then one needs to explicitly fail the test:
+            Assertions.fail();
+        }
+    }
+
+    private static class FabilandCarOwnership implements CreateCarOwnershipModel {
+
+        private final Random random;
+
+        FabilandCarOwnership() {
+            this.random = SiloUtil.provideNewRandom();
+        }
+
+        @Override
+        public void run() {
+
+        }
+
+        @Override
+        public void simulateCarOwnership(Household hh) {
+            hh.setAutos(random.nextInt(3)+1);
+        }
+    }
+}
+