From 36bc4997a3b9ed28079e7039337a482d599f583c Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 10 Mar 2026 16:20:57 +0100
Subject: [PATCH 01/29] Add note and rule functions to printer.py

---
 printer.py | 36 ++++++++++++++++++++++++++++++++++++
 1 file changed, 36 insertions(+)

diff --git a/printer.py b/printer.py
index 1517c9c..9cb772f 100644
--- a/printer.py
+++ b/printer.py
@@ -217,6 +217,37 @@ def _log(self, text: str) -> None:
                         f.write(line + "\n")
         return None
 
+    def note(self, text: str, prefix: str = "Note: ", hard_wrap_chars: str = None) -> None:
+        """
+        Handles a note message. Printed in yellow (same as warning but semantically different).
+
+        :param text: Text to be printed
+        :param prefix: Prefix (default: "Note: ")
+        :param hard_wrap_chars: Chars to be added at the end of the text if it exceeds the console width,
+        None if text should not be truncated
+        :return: None
+        """
+        text = self.handle_hard_wrap_chars(text, prefix, hard_wrap_chars)
+        if len(prefix) > 0:
+            self.console.print(f"[yellow]{escape(prefix)}[/yellow]{escape(text)}")
+        else:
+            self.console.print(f"[yellow]{escape(text)}[/yellow]")
+        self._log(f"{prefix}{text}")
+        return None
+
+    def rule(self, title: str = "", style: str = "white") -> None:
+        """
+        Prints a horizontal rule (line) to the console, optionally with a centered title.
+        Delegates to rich.console.Console.rule().
+
+        :param title: Optional title to display in the center of the rule
+        :param style: Rich style string for the rule (default: "white")
+        :return: None
+        """
+        self.console.rule(title, style=style)
+        self._log(f"--- {title} ---" if title else "-" * self.console.width)
+        return None
+
     def separator(self) -> None:
         """
         Prints a separator line to the console. The line is made up of dashes and
@@ -280,3 +311,8 @@ def pprint_zoi_var(var, zoi, index_positions: list = None, decimals: int = 2):
     # Iterate over all lists and print the values
     for i in range(len(value_list)):
         print(f"    {key_list[i]:>{key_spacer}} : {lower_list[i]:>{lower_spacer}} : {value_list[i]:>{value_spacer}} : {upper_list[i]:>{upper_spacer}} : {fixed_list[i]:>{fixed_spacer}} : {stale_list[i]:>{stale_spacer}} : {domain_list[i]:>{domain_spacer}}")
+
+
+def timestamp() -> str:
+    """Returns a timestamp string in the format 'YYYYMMDD_HHMMSS'."""
+    return datetime.datetime.now().strftime('%Y%m%d_%H%M%S')

From 87e9bc76ca2f2c978bc4a7ec38aa8ba8c219dc1e Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 10 Mar 2026 16:37:38 +0100
Subject: [PATCH 02/29] Move Caller.py from LEGO-Pyomo to InOutModule

---
 Caller.py | 128 ++++++++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 128 insertions(+)
 create mode 100644 Caller.py

diff --git a/Caller.py b/Caller.py
new file mode 100644
index 0000000..21d6b5d
--- /dev/null
+++ b/Caller.py
@@ -0,0 +1,128 @@
+import argparse
+import datetime
+import os
+import subprocess
+import time
+
+from rich.highlighter import ReprHighlighter
+from rich.live import Live
+from rich.text import Text
+
+_highlighter = ReprHighlighter()
+
+from printer import Printer
+
+printer = Printer.getInstance()
+
+parser = argparse.ArgumentParser(description='Calls the exact lines from the given file, can be called multiple times.')
+
+parser.add_argument('jobs', type=str, help='Path to the text-file containing the commands to be called.')
+parser.add_argument("--spawn", type=int, help='Number of jobs to spawn (if this is specified, it will call itself multiple times)', nargs='?', default=0)
+args = parser.parse_args()
+printer.information(f"Using jobs from '{args.jobs}'")
+
+if args.spawn >= 1:
+    printer.information(f"Spawning {args.spawn} parallel jobs")
+    for i in range(args.spawn):
+        subprocess.Popen([
+            "cmd", "/c", "start", f"Caller {i}: {args.jobs}", "cmd", "/k",
+            f"set POST_ACTIVATE_COMMAND=python {os.path.abspath(__file__)} {args.jobs} && call Conda-Activation-Scripts/activate_environment_windows.bat"
+        ])
+    printer.information(f"Spawned {args.spawn} parallel jobs, exiting... ")
+    exit(0)
+
+barrier_waited = {}  # tracks total seconds waited per barrier line index; resets when barrier passes or a job is picked up
+
+
+def _all_previous_done(jobs_file, lines, barrier_index):
+    return all(
+        os.path.exists(f"{jobs_file}.finished{j}") or os.path.exists(f"{jobs_file}.error{j}")
+        for j in range(barrier_index) if lines[j].strip() not in ("---", "")
+    )
+
+
+def _any_previous_unclaimed(jobs_file, lines, barrier_index):
+    return any(
+        not os.path.exists(f"{jobs_file}.started{j}")
+        and not os.path.exists(f"{jobs_file}.finished{j}")
+        and not os.path.exists(f"{jobs_file}.error{j}")
+        for j in range(barrier_index) if lines[j].strip() not in ("---", "")
+    )
+
+
+while True:
+    with open(args.jobs, 'r') as f:
+        lines = f.readlines()
+
+    found_one = False
+    restart = False
+    for i, line in enumerate(lines):
+        started_job_flag = f"{args.jobs}.started{i}"
+        finished_job_flag = f"{args.jobs}.finished{i}"
+        error_job_flag = f"{args.jobs}.error{i}"
+        if line.strip() == "---":
+            dot_i = 0
+            found_unclaimed = False
+            with Live(console=printer.console, refresh_per_second=2) as live:
+                while not _all_previous_done(args.jobs, lines, i):
+                    if _any_previous_unclaimed(args.jobs, lines, i):
+                        found_unclaimed = True
+                        break
+                    barrier_waited[i] = barrier_waited.get(i, 0) + 3
+                    dots = "." * ((dot_i % 5) + 1)
+                    msg = _highlighter(Text(f"Barrier '---' at line {i}: Checking every 3s, waited {barrier_waited[i]}s already"))
+                    msg.append(dots)
+                    live.update(msg)
+                    dot_i += 1
+                    time.sleep(3)
+            if found_unclaimed:
+                restart = True
+                break  # restart outer while loop to pick up the unclaimed job
+            barrier_waited.pop(i, None)
+            continue  # barrier cleared, keep scanning for next job
+
+        try:
+            fd = os.open(started_job_flag, os.O_CREAT | os.O_EXCL | os.O_WRONLY)
+        except FileExistsError:
+            continue  # Another worker already claimed this job
+        if os.path.exists(finished_job_flag) or os.path.exists(error_job_flag):
+            os.close(fd)
+            continue
+
+        if barrier_waited:
+            barrier_waited.clear()
+        start_datetime = datetime.datetime.now().strftime("%Y-%m-%d %H:%M:%S")
+        with os.fdopen(fd, 'w') as f:
+            f.write(f"Command: {line.strip()}\n")
+            f.write(f"Started at: {start_datetime}")
+        found_one = True
+        try:
+            printer.information(f"Executing job {i} from '{args.jobs}': {line.strip()}")
+            os.system(f"title Job {i} from '{args.jobs}': {line.strip()}")
+
+            start_time = time.time()
+            exit_code = os.system(line.strip())
+            end_time = time.time()
+
+            if exit_code != 0:
+                raise RuntimeError(f"Command exited with code {exit_code}")
+
+            with open(finished_job_flag, 'w') as f:
+                f.write(f"Command: {line.strip()}\n")
+                f.write(f"Started at:  {start_datetime}\n")
+                f.write(f"Finished at: {datetime.datetime.now().strftime("%Y-%m-%d %H:%M:%S")}\n")
+                f.write(f"Execution time: {end_time - start_time:.2f} seconds (= {(end_time - start_time) / 60 / 60:.2f} hours)\n")
+
+            printer.information(f"Finished job {i} from '{args.jobs}' after {end_time - start_time:.2f} seconds (= {(end_time - start_time) / 60 / 60:.2f} hours).")
+        except Exception as e:
+            printer.error(f"Error while executing job {i}: {e}")
+            with open(error_job_flag, 'w') as f:
+                f.write(f"Command: {line.strip()}\n")
+                f.write(f"Error while executing job {i} from '{args.jobs}': {e}\n")
+                f.write(f"Started at:  {start_datetime}\n")
+                f.write(f"Occurred at: {datetime.datetime.now().strftime("%Y-%m-%d %H:%M:%S")}")
+        break
+
+    if not found_one and not restart:
+        printer.information(f"No more jobs to execute in '{args.jobs}', exiting.")
+        break

From 2a19d6bdedb45e0f7815d8c9446e49007ca8ddef Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 7 Apr 2026 13:54:50 +0200
Subject: [PATCH 03/29] Add log file to Caller.py

---
 Caller.py | 41 ++++++++++++++++++++++++++++++++++++-----
 1 file changed, 36 insertions(+), 5 deletions(-)

diff --git a/Caller.py b/Caller.py
index 21d6b5d..5b28e15 100644
--- a/Caller.py
+++ b/Caller.py
@@ -14,6 +14,17 @@
 
 printer = Printer.getInstance()
 
+
+def _tail(text, n=20):
+    """Return the last n lines of text, or all if fewer."""
+    if not text:
+        return "(no output)"
+    lines = text.rstrip('\n').split('\n')
+    if len(lines) <= n:
+        return text
+    return f"... ({len(lines) - n} lines omitted)\n" + '\n'.join(lines[-n:])
+
+
 parser = argparse.ArgumentParser(description='Calls the exact lines from the given file, can be called multiple times.')
 
 parser.add_argument('jobs', type=str, help='Path to the text-file containing the commands to be called.')
@@ -96,22 +107,41 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
             f.write(f"Command: {line.strip()}\n")
             f.write(f"Started at: {start_datetime}")
         found_one = True
+        log_file = f"{args.jobs}.log{i}"
         try:
             printer.information(f"Executing job {i} from '{args.jobs}': {line.strip()}")
             os.system(f"title Job {i} from '{args.jobs}': {line.strip()}")
 
             start_time = time.time()
-            exit_code = os.system(line.strip())
+            result = subprocess.run(
+                line.strip(),
+                shell=True,
+                stdout=subprocess.PIPE,
+                stderr=subprocess.STDOUT,
+                text=True,
+            )
             end_time = time.time()
 
-            if exit_code != 0:
-                raise RuntimeError(f"Command exited with code {exit_code}")
+            with open(log_file, 'w') as f:
+                f.write(f"Command: {line.strip()}\n")
+                f.write(f"Started at:  {start_datetime}\n")
+                f.write(f"Exit code:   {result.returncode}\n")
+                f.write(f"{'=' * 60}\n")
+                f.write(result.stdout or "")
+
+            if result.returncode != 0:
+                raise RuntimeError(
+                    f"Command exited with code {result.returncode}. "
+                    f"See log: {log_file}\n"
+                    f"Last output:\n{_tail(result.stdout, 20)}"
+                )
 
             with open(finished_job_flag, 'w') as f:
                 f.write(f"Command: {line.strip()}\n")
                 f.write(f"Started at:  {start_datetime}\n")
-                f.write(f"Finished at: {datetime.datetime.now().strftime("%Y-%m-%d %H:%M:%S")}\n")
+                f.write(f"Finished at: {datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S')}\n")
                 f.write(f"Execution time: {end_time - start_time:.2f} seconds (= {(end_time - start_time) / 60 / 60:.2f} hours)\n")
+                f.write(f"Log file: {log_file}\n")
 
             printer.information(f"Finished job {i} from '{args.jobs}' after {end_time - start_time:.2f} seconds (= {(end_time - start_time) / 60 / 60:.2f} hours).")
         except Exception as e:
@@ -120,7 +150,8 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
                 f.write(f"Command: {line.strip()}\n")
                 f.write(f"Error while executing job {i} from '{args.jobs}': {e}\n")
                 f.write(f"Started at:  {start_datetime}\n")
-                f.write(f"Occurred at: {datetime.datetime.now().strftime("%Y-%m-%d %H:%M:%S")}")
+                f.write(f"Occurred at: {datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S')}\n")
+                f.write(f"Log file: {log_file}\n")
         break
 
     if not found_one and not restart:

From 048d03cc7cdf5a7ec7b01b679096813ef6eda9c8 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 7 Apr 2026 14:15:37 +0200
Subject: [PATCH 04/29] Adjust log-handling for Caller.py

---
 Caller.py | 54 +++++++++++++++++++++++++++++++-----------------------
 1 file changed, 31 insertions(+), 23 deletions(-)

diff --git a/Caller.py b/Caller.py
index 5b28e15..7aead2b 100644
--- a/Caller.py
+++ b/Caller.py
@@ -2,6 +2,7 @@
 import datetime
 import os
 import subprocess
+import sys
 import time
 
 from rich.highlighter import ReprHighlighter
@@ -15,14 +16,18 @@
 printer = Printer.getInstance()
 
 
-def _tail(text, n=20):
-    """Return the last n lines of text, or all if fewer."""
-    if not text:
+def _tail_file(filepath, n=20):
+    """Return the last n lines of a file, or all if fewer."""
+    try:
+        with open(filepath, 'r', errors='replace') as f:
+            lines = f.readlines()
+    except OSError:
+        return "(could not read log file)"
+    if not lines:
         return "(no output)"
-    lines = text.rstrip('\n').split('\n')
     if len(lines) <= n:
-        return text
-    return f"... ({len(lines) - n} lines omitted)\n" + '\n'.join(lines[-n:])
+        return ''.join(lines)
+    return f"... ({len(lines) - n} lines omitted)\n" + ''.join(lines[-n:])
 
 
 parser = argparse.ArgumentParser(description='Calls the exact lines from the given file, can be called multiple times.')
@@ -113,27 +118,30 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
             os.system(f"title Job {i} from '{args.jobs}': {line.strip()}")
 
             start_time = time.time()
-            result = subprocess.run(
-                line.strip(),
-                shell=True,
-                stdout=subprocess.PIPE,
-                stderr=subprocess.STDOUT,
-                text=True,
-            )
+            with open(log_file, 'w') as log_f:
+                log_f.write(f"Command: {line.strip()}\n")
+                log_f.write(f"Started at:  {start_datetime}\n")
+                log_f.write(f"{'=' * 60}\n")
+                log_f.flush()
+
+                proc = subprocess.Popen(
+                    line.strip(),
+                    shell=True,
+                    stdout=subprocess.PIPE,
+                    stderr=subprocess.STDOUT,
+                )
+                for raw_line in proc.stdout:
+                    sys.stdout.buffer.write(raw_line)
+                    sys.stdout.buffer.flush()
+                    log_f.write(raw_line.decode(errors='replace'))
+                proc.wait()
             end_time = time.time()
 
-            with open(log_file, 'w') as f:
-                f.write(f"Command: {line.strip()}\n")
-                f.write(f"Started at:  {start_datetime}\n")
-                f.write(f"Exit code:   {result.returncode}\n")
-                f.write(f"{'=' * 60}\n")
-                f.write(result.stdout or "")
-
-            if result.returncode != 0:
+            if proc.returncode != 0:
                 raise RuntimeError(
-                    f"Command exited with code {result.returncode}. "
+                    f"Command exited with code {proc.returncode}. "
                     f"See log: {log_file}\n"
-                    f"Last output:\n{_tail(result.stdout, 20)}"
+                    f"Last output:\n{_tail_file(log_file, 20)}"
                 )
 
             with open(finished_job_flag, 'w') as f:

From 9611e2a41a922e05758795bebd8e62d4c1b36754 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Wed, 8 Apr 2026 13:37:34 +0200
Subject: [PATCH 05/29] Improve KeyboardInterrupt-handling for Caller.py

---
 Caller.py | 13 ++++++++++---
 1 file changed, 10 insertions(+), 3 deletions(-)

diff --git a/Caller.py b/Caller.py
index 7aead2b..8f6b165 100644
--- a/Caller.py
+++ b/Caller.py
@@ -130,11 +130,18 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
                     stdout=subprocess.PIPE,
                     stderr=subprocess.STDOUT,
                 )
+                interrupted = False
                 for raw_line in proc.stdout:
-                    sys.stdout.buffer.write(raw_line)
-                    sys.stdout.buffer.flush()
-                    log_f.write(raw_line.decode(errors='replace'))
+                    decoded = raw_line.decode(errors='replace')
+                    try:
+                        sys.stdout.write(decoded)
+                        sys.stdout.flush()
+                    except KeyboardInterrupt:
+                        interrupted = True
+                    log_f.write(decoded)
                 proc.wait()
+                if interrupted:
+                    raise KeyboardInterrupt
             end_time = time.time()
 
             if proc.returncode != 0:

From c52ac140dc7324c40c11f3176d171685ac7b4df9 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 14 Apr 2026 14:08:40 +0200
Subject: [PATCH 06/29] Save MIP-gap in .sqlite, take LEGO object as parameter

---
 SQLiteWriter.py | 36 +++++++++++++++++++++++++++++-------
 1 file changed, 29 insertions(+), 7 deletions(-)

diff --git a/SQLiteWriter.py b/SQLiteWriter.py
index 7547aeb..f75826e 100644
--- a/SQLiteWriter.py
+++ b/SQLiteWriter.py
@@ -6,13 +6,14 @@
 import pyomo.environ as pyo
 
 from InOutModule.printer import Printer
+from LEGO.LEGO import LEGO
 
 printer = Printer.getInstance()
 
 
-def model_to_sqlite(model: pyo.base.Model, filename: str) -> None:
+def model_to_sqlite(model: pyo.base.ConcreteModel, filename: str) -> None:
     """
-    Save the model to a SQLite database.
+    Save the model to an SQLite database.
     Automatically includes objective decomposition and dual values.
 
     :param model: Pyomo model to save
@@ -61,14 +62,17 @@ def model_to_sqlite(model: pyo.base.Model, filename: str) -> None:
     pass
 
 
-def add_solver_statistics_to_sqlite(filename: str, results, work_units=None) -> None:
+def add_solver_statistics_to_sqlite(filename: str, lego: LEGO) -> None:
     """
-    Add solver statistics (like Gurobi work-units) to an existing SQLite database.
+    Add solver statistics to an existing SQLite database, reading from a LEGO instance.
     :param filename: Path to the SQLite database file
-    :param results: Pyomo solver results object
-    :param work_units: Optional work units value (from Gurobi solver)
+    :param lego: Solved LEGO instance (reads lego.results, lego.work_units, lego.mip_gap)
     :return: None
     """
+    results = lego.results
+    work_units = lego.work_units
+    mip_gap = lego.mip_gap
+
     cnx = sqlite3.connect(filename)
 
     # Extract solver statistics
@@ -106,13 +110,31 @@ def add_solver_statistics_to_sqlite(filename: str, results, work_units=None) ->
                     if value is not None:
                         stats[attr] = float(value) if isinstance(value, (int, float)) else str(value)
 
+        # MIP gap
+        if mip_gap is not None:
+            stats['mip_gap'] = float(mip_gap)
+        else:
+            # Try to derive from bounds as fallback
+            lb = stats.get('lower_bound')
+            ub = stats.get('upper_bound')
+            if lb is not None and ub is not None and isinstance(lb, float) and isinstance(ub, float) and ub != 0:
+                derived_gap = abs(ub - lb) / abs(ub)
+                if derived_gap > 1e-10:
+                    stats['mip_gap'] = derived_gap
+                    printer.information(f"MIP gap derived from bounds: {derived_gap:.6f}")
+                else:
+                    printer.information("MIP gap not stored: bounds are equal (LP or MIP solved to optimality)")
+            else:
+                printer.information("MIP gap not available: not provided and bounds not found in results")
+
         # Create a DataFrame with solver statistics
         if stats:
             df = pd.DataFrame([stats])
             df.to_sql('solver_statistics', cnx, if_exists='replace', index=False)
             cnx.commit()
             work_units_str = f"{stats['work_units']:.2f}" if 'work_units' in stats else 'N/A'
-            printer.information(f"Added solver statistics to SQLite database (work_units: {work_units_str})")
+            mip_gap_str = f"{stats['mip_gap']:.6f}" if 'mip_gap' in stats else 'N/A'
+            printer.information(f"Added solver statistics to SQLite database (work_units: {work_units_str}, mip_gap: {mip_gap_str})")
         else:
             printer.warning("No solver statistics found in results object")
 

From 0064c9729e14f19d42c5ae405b2a406f16c9c780 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 14 Apr 2026 14:49:22 +0200
Subject: [PATCH 07/29] Add README.md and CLAUDE.md including rules

---
 .claude/rules/claude-md-maintenance.md | 23 ++++++
 CLAUDE.md                              | 47 ++++++++++++
 README.md                              | 99 +++++++++++++++++++++++++-
 3 files changed, 168 insertions(+), 1 deletion(-)
 create mode 100644 .claude/rules/claude-md-maintenance.md
 create mode 100644 CLAUDE.md

diff --git a/.claude/rules/claude-md-maintenance.md b/.claude/rules/claude-md-maintenance.md
new file mode 100644
index 0000000..aa5c6e5
--- /dev/null
+++ b/.claude/rules/claude-md-maintenance.md
@@ -0,0 +1,23 @@
+---
+description: Rules for keeping CLAUDE.md and README.md files up to date after code changes in InOutModule/
+---
+
+After making any code changes in `InOutModule/`, update the relevant documentation files.
+
+**`CLAUDE.md`** — context for Claude that is hard to derive from reading the code: non-obvious patterns, gotchas, architectural decisions, naming conventions, and implementation constraints.
+
+**`README.md`** — human-facing documentation: usage instructions, API overview, key concepts, and examples.
+
+**Rules:**
+- Update the file(s) closest to the code that changed.
+- If a change also affects `LEGO/` or root-level code, update the corresponding docs there too (see root `.claude/rules/claude-md-maintenance.md` for the full index).
+- Do not duplicate content between `README.md` and `CLAUDE.md` — if something is in `README.md`, `CLAUDE.md` should reference it, not repeat it.
+- Use prose references (`See README.md for ...`) to point to large human-facing docs, not `@`-imports.
+- Delete stale entries rather than commenting them out.
+
+**Index of documentation files in InOutModule/:**
+
+| File                    | Audience | Contents                                                                |
+|-------------------------|----------|-------------------------------------------------------------------------|
+| `InOutModule/CLAUDE.md` | Claude   | CaseStudy read order, ExcelReader/Writer patterns, SQLiteWriter, Caller |
+| `InOutModule/README.md` | Humans   | CaseStudy API, SQLiteWriter usage, Caller, Excel file format            |
diff --git a/CLAUDE.md b/CLAUDE.md
new file mode 100644
index 0000000..b7e8c04
--- /dev/null
+++ b/CLAUDE.md
@@ -0,0 +1,47 @@
+# CLAUDE.md — InOutModule
+
+This file provides guidance to Claude Code when working with code in this directory.
+See `README.md` for usage, key concepts, and data structure.
+
+## Architecture Notes
+
+### CaseStudy
+
+- Sequential reads happen first (`Global_Parameters`, `Global_Scenarios`, `Power_Parameters`), because subsequent file selection depends on `pEnable*` flags.
+- All remaining files are read in parallel via `ThreadPoolExecutor` — order of assignment is non-deterministic, so no file read may depend on another parallel read.
+- `dPower_WeightsRP` is **computed** from `dPower_Hindex` (counting occurrences per `rp`); if a `Power_WeightsRP.xlsx` file also exists, it is read and compared — a mismatch triggers a warning but uses the **file** value, not the computed one.
+- `merge_single_node_buses()` preserves the `z` (zone) column as a sorted unique union string of all merged zones (e.g. `"R1_R2"`). This is documented in root `CLAUDE.md` as well.
+- `CaseStudy.copy()` is a full `deepcopy` — safe to modify independently.
+- Transition matrices (`rpTransitionMatrixAbsolute`, `rpTransitionMatrixRelativeTo`, `rpTransitionMatrixRelativeFrom`) are computed in the constructor and attached as attributes.
+
+### ExcelReader
+
+- Excel sheets whose name starts with `~` are silently skipped (used to disable scenarios in a multi-sheet file without deleting them).
+- All Excel files have a version specifier in cell `C2` of each sheet. `check_LEGOExcel_version()` warns (or raises, if `fail_on_wrong_version=True`) on mismatch — wrong version can cause silent column misreads.
+- The reader uses `calamine` engine (fast), not `openpyxl`.
+
+### ExcelWriter
+
+- All cell styles, column definitions, and table layouts are declared in `TableDefinitions.xml`, not in Python code. When adding a new output table, define its columns there first.
+- `ExcelWriter.__init__()` parses the XML once and stores resolved objects (`self.columns`, `self.cell_styles`, etc.). Avoid re-instantiating per row.
+
+### SQLiteWriter
+
+- `model_to_sqlite()` automatically calls `add_objective_decomposition_to_sqlite()` and `add_dual_values_to_sqlite()` — these do not need to be called separately.
+- `add_run_parameters_to_sqlite()` stores all run configuration in the `run_parameters` table. Evaluation scripts should read from this table (more reliable than filename parsing).
+- Pyomo component types not handled by the writer emit a `printer.warning()` and are skipped silently — add new `case` branches if new Pyomo types need to be stored.
+
+### Utilities
+
+- `inflowsToCapacityFactors()` joins inflows onto `vresProfiles_df` by dividing by `MaxProd`; generators with missing or zero `MaxProd` are dropped with a warning.
+- `capacityFactorsToInflows()` is the inverse; the `remove_Inflows_from_VRESProfiles_inplace` flag modifies the input DataFrame in place when set.
+
+### Printer
+
+- `Printer` is a singleton — obtain the instance with `Printer.getInstance()`, never call the constructor directly.
+- `set_logfile(path)` redirects all subsequent output to a file (appending). Set to `None` to stop logging.
+
+### Caller
+
+- `Caller.py` is a parallel job runner reading from a text file. It uses sentinel files (`.finished{n}`, `.error{n}`) for barrier synchronization across parallel workers.
+- Lines containing only `---` act as barriers — workers wait until all prior jobs are complete before continuing past the barrier.
diff --git a/README.md b/README.md
index 6c553be..43703ea 100644
--- a/README.md
+++ b/README.md
@@ -1 +1,98 @@
-# InOutModule
\ No newline at end of file
+# InOutModule
+
+Data I/O package for LEGO-Pyomo. Handles reading Excel case study files, writing results to Excel and SQLite, and utility data transformations.
+
+## Key Files
+
+| File                   | Purpose                                                                 |
+|------------------------|-------------------------------------------------------------------------|
+| `CaseStudy.py`         | Loads all Excel input files into a single object; data manipulation     |
+| `ExcelReader.py`       | Low-level Excel parsing (per-file reader functions, version checking)   |
+| `ExcelWriter.py`       | Writes formatted Excel output files, driven by `TableDefinitions.xml`  |
+| `SQLiteWriter.py`      | Exports Pyomo model results to SQLite; stores solver stats & run params |
+| `Utilities.py`         | Data transformations (inflows ↔ capacity factors, Printer helpers)      |
+| `printer.py`           | Singleton console/logfile printer with severity levels                  |
+| `Caller.py`            | Parallel job runner for batch experiments                               |
+| `TableDefinitions.xml` | Declarative column/style definitions used by `ExcelWriter`              |
+| `PypsaReader.py`       | Imports case studies from PyPSA network objects                         |
+| `nrel118-reader.py`    | Converts NREL 118-bus data into LEGO Excel format                       |
+
+## CaseStudy
+
+### Constructor
+
+```python
+CaseStudy(
+    data_folder: str | Path,
+    do_not_scale_units: bool = False,
+    do_not_merge_single_node_buses: bool = False,
+    parallel_read: bool = True,
+    n_jobs: int = 4,
+    # Per-file overrides: pass a DataFrame to skip reading from disk
+    dPower_ThermalGen: pd.DataFrame = None,
+    ...
+)
+```
+
+All Excel files in `data_folder` are read automatically. Any `d*` parameter can be passed directly as a DataFrame to bypass file reading (useful for programmatic construction or testing).
+
+### DataFrame Attributes
+
+| Attribute              | Time dependency | Source file                  |
+|------------------------|-----------------|------------------------------|
+| `dGlobal_Parameters`   | none            | `Global_Parameters.xlsx`     |
+| `dGlobal_Scenarios`    | none            | `Global_Scenarios.xlsx`      |
+| `dPower_Parameters`    | none            | `Power_Parameters.xlsx`      |
+| `dPower_BusInfo`       | none            | `Power_BusInfo.xlsx`         |
+| `dPower_Network`       | none            | `Power_Network.xlsx`         |
+| `dPower_ThermalGen`    | none            | `Power_ThermalGen.xlsx`      |
+| `dPower_VRES`          | none            | `Power_VRES.xlsx`            |
+| `dPower_Storage`       | none            | `Power_Storage.xlsx`         |
+| `dPower_Demand`        | rp + k          | `Power_Demand.xlsx`          |
+| `dPower_VRESProfiles`  | rp + k          | `Power_VRESProfiles.xlsx`    |
+| `dPower_Inflows`       | rp + k          | `Power_Inflows.xlsx`         |
+| `dPower_ImportExport`  | rp + k          | `Power_ImportExport.xlsx`    |
+| `dPower_Hindex`        | rp + k          | `Power_Hindex.xlsx`          |
+| `dPower_WeightsRP`     | rp only         | `Power_WeightsRP.xlsx`       |
+| `dPower_WeightsK`      | k only          | `Power_WeightsK.xlsx`        |
+
+### Key Methods
+
+| Method                         | Description                                              |
+|--------------------------------|----------------------------------------------------------|
+| `copy()`                       | Deep copy — safe to modify independently                 |
+| `equal_to(cs)`                 | Compare all DataFrames with another CaseStudy            |
+| `merge_single_node_buses()`    | Collapse single-bus zones; preserves `z` as union string |
+| `scale_CaseStudy()`            | Applies power and cost scaling factors from parameters   |
+| `get_rpTransitionMatrices()`   | Returns absolute and relative transition matrices        |
+
+## SQLiteWriter
+
+```python
+from InOutModule.SQLiteWriter import model_to_sqlite, add_run_parameters_to_sqlite, add_solver_statistics_to_sqlite
+
+model_to_sqlite(model, "output/results.sqlite")
+add_solver_statistics_to_sqlite("output/results.sqlite", lego)
+add_run_parameters_to_sqlite("output/results.sqlite", zoi="R1", dc_buffer=2)
+```
+
+- `model_to_sqlite()` exports all Pyomo variables, parameters, and sets; automatically appends objective decomposition and dual values.
+- `add_run_parameters_to_sqlite()` creates a `run_parameters` table — evaluation scripts should read from this table rather than parsing filenames.
+
+## Caller (Parallel Job Runner)
+
+```bash
+python InOutModule/Caller.py jobs.txt
+python InOutModule/Caller.py jobs.txt --spawn 4   # open 4 parallel terminal windows
+```
+
+`jobs.txt` is a plain-text file with one shell command per line. Lines containing only `---` act as barriers: all workers wait until every job above the barrier is finished before continuing.
+
+## Excel File Format
+
+All input Excel files follow a versioned multi-sheet format:
+- Each sheet corresponds to one scenario (or `ScenarioA` for deterministic runs).
+- Sheets whose name starts with `~` are skipped.
+- Cell `C2` on each sheet contains a version specifier (e.g., `v0.1.0`); mismatches produce a warning.
+
+See `changelog-LEGOExcels.md` for version history of the Excel format.

From aecee05722c4ecbdc3af5923029cd71fbc052a6e Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 21 Apr 2026 14:29:28 +0200
Subject: [PATCH 08/29] Fix CTRL+C for Caller.py

---
 Caller.py | 22 +++++++++++++---------
 1 file changed, 13 insertions(+), 9 deletions(-)

diff --git a/Caller.py b/Caller.py
index 8f6b165..83ccb75 100644
--- a/Caller.py
+++ b/Caller.py
@@ -130,18 +130,22 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
                     stdout=subprocess.PIPE,
                     stderr=subprocess.STDOUT,
                 )
-                interrupted = False
-                for raw_line in proc.stdout:
-                    decoded = raw_line.decode(errors='replace')
-                    try:
+                try:
+                    for raw_line in proc.stdout:
+                        decoded = raw_line.decode(errors='replace').replace('\r\n', '\n').replace('\r', '\n')
                         sys.stdout.write(decoded)
                         sys.stdout.flush()
-                    except KeyboardInterrupt:
-                        interrupted = True
-                    log_f.write(decoded)
+                        log_f.write(decoded)
+                except KeyboardInterrupt:
+                    # CTRL+C also signals Gurobi, which does a graceful shutdown and saves
+                    # results. Drain remaining output so the log is complete, then let the
+                    # subprocess finish — do not re-raise here.
+                    for raw_line in proc.stdout:
+                        decoded = raw_line.decode(errors='replace').replace('\r\n', '\n').replace('\r', '\n')
+                        sys.stdout.write(decoded)
+                        sys.stdout.flush()
+                        log_f.write(decoded)
                 proc.wait()
-                if interrupted:
-                    raise KeyboardInterrupt
             end_time = time.time()
 
             if proc.returncode != 0:

From f456242510acab939fb4b13baca9fd4700ea58e3 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Wed, 22 Apr 2026 13:54:43 +0200
Subject: [PATCH 09/29] Fix double-newlines in Caller.py logs

---
 Caller.py | 34 ++++++++++++++++------------------
 1 file changed, 16 insertions(+), 18 deletions(-)

diff --git a/Caller.py b/Caller.py
index 83ccb75..56a91ec 100644
--- a/Caller.py
+++ b/Caller.py
@@ -5,12 +5,6 @@
 import sys
 import time
 
-from rich.highlighter import ReprHighlighter
-from rich.live import Live
-from rich.text import Text
-
-_highlighter = ReprHighlighter()
-
 from printer import Printer
 
 printer = Printer.getInstance()
@@ -79,18 +73,22 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
         if line.strip() == "---":
             dot_i = 0
             found_unclaimed = False
-            with Live(console=printer.console, refresh_per_second=2) as live:
-                while not _all_previous_done(args.jobs, lines, i):
-                    if _any_previous_unclaimed(args.jobs, lines, i):
-                        found_unclaimed = True
-                        break
-                    barrier_waited[i] = barrier_waited.get(i, 0) + 3
-                    dots = "." * ((dot_i % 5) + 1)
-                    msg = _highlighter(Text(f"Barrier '---' at line {i}: Checking every 3s, waited {barrier_waited[i]}s already"))
-                    msg.append(dots)
-                    live.update(msg)
-                    dot_i += 1
-                    time.sleep(3)
+            last_len = 0
+            while not _all_previous_done(args.jobs, lines, i):
+                if _any_previous_unclaimed(args.jobs, lines, i):
+                    found_unclaimed = True
+                    break
+                barrier_waited[i] = barrier_waited.get(i, 0) + 3
+                dots = "." * ((dot_i % 5) + 1)
+                msg = f"Barrier '---' at line {i}: Checking every 3s, waited {barrier_waited[i]}s already{dots}"
+                sys.stdout.write(f"\r{msg}{' ' * max(0, last_len - len(msg))}")
+                sys.stdout.flush()
+                last_len = len(msg)
+                dot_i += 1
+                time.sleep(3)
+            if last_len:
+                sys.stdout.write('\n')
+                sys.stdout.flush()
             if found_unclaimed:
                 restart = True
                 break  # restart outer while loop to pick up the unclaimed job

From 823605cf36e3a4e5efcb828f906a8fc9bda45842 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Thu, 23 Apr 2026 14:50:08 +0200
Subject: [PATCH 10/29] Add second check for KeyboardInterrupt

---
 Caller.py | 17 ++++++++++++-----
 1 file changed, 12 insertions(+), 5 deletions(-)

diff --git a/Caller.py b/Caller.py
index 56a91ec..20ef229 100644
--- a/Caller.py
+++ b/Caller.py
@@ -138,11 +138,18 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
                     # CTRL+C also signals Gurobi, which does a graceful shutdown and saves
                     # results. Drain remaining output so the log is complete, then let the
                     # subprocess finish — do not re-raise here.
-                    for raw_line in proc.stdout:
-                        decoded = raw_line.decode(errors='replace').replace('\r\n', '\n').replace('\r', '\n')
-                        sys.stdout.write(decoded)
-                        sys.stdout.flush()
-                        log_f.write(decoded)
+                    try:
+                        for raw_line in proc.stdout:
+                            decoded = raw_line.decode(errors='replace').replace('\r\n', '\n').replace('\r', '\n')
+                            sys.stdout.write(decoded)
+                            sys.stdout.flush()
+                            log_f.write(decoded)
+                    except KeyboardInterrupt:
+                        remaining = proc.stdout.read()
+                        if remaining:
+                            decoded = remaining.decode(errors='replace').replace('\r\n', '\n').replace('\r', '\n')
+                            sys.stdout.write(decoded)
+                            log_f.write(decoded)
                 proc.wait()
             end_time = time.time()
 

From d662c80f04b7c6490aa29063fbfebdd511df240f Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 28 Apr 2026 21:40:29 +0200
Subject: [PATCH 11/29] Add merge_generators to CaseStudy

---
 CLAUDE.md    |   1 +
 CaseStudy.py | 163 +++++++++++++++++++++++++++++++++++++++++++++++++++
 2 files changed, 164 insertions(+)

diff --git a/CLAUDE.md b/CLAUDE.md
index b7e8c04..6767a99 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -11,6 +11,7 @@ See `README.md` for usage, key concepts, and data structure.
 - All remaining files are read in parallel via `ThreadPoolExecutor` — order of assignment is non-deterministic, so no file read may depend on another parallel read.
 - `dPower_WeightsRP` is **computed** from `dPower_Hindex` (counting occurrences per `rp`); if a `Power_WeightsRP.xlsx` file also exists, it is read and compared — a mismatch triggers a warning but uses the **file** value, not the computed one.
 - `merge_single_node_buses()` preserves the `z` (zone) column as a sorted unique union string of all merged zones (e.g. `"R1_R2"`). This is documented in root `CLAUDE.md` as well.
+- `merge_generators()` collapses all generators sharing the same `(tec, i)` into one representative generator with ID `"{i}_{tec}"`. It must be called after construction (scaling is not required). VRESProfiles and Inflows are merged **before** VRES so that the original per-generator `MaxProd` weights are available for the capacity-factor weighted average. Generators in VRESProfiles/Inflows that have no matching entry in `dPower_VRES` (left-join miss) are grouped under `(tec=NaN, i=NaN)` — filter to a single scenario first to avoid mixing scenarios in the groupby.
 - `CaseStudy.copy()` is a full `deepcopy` — safe to modify independently.
 - Transition matrices (`rpTransitionMatrixAbsolute`, `rpTransitionMatrixRelativeTo`, `rpTransitionMatrixRelativeFrom`) are computed in the constructor and attached as attributes.
 
diff --git a/CaseStudy.py b/CaseStudy.py
index 8f1bb5f..9fccb0c 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -627,6 +627,169 @@ def merge_single_node_buses(self, inplace: bool = True) -> typing.Optional[typin
 
         return cs if not inplace else None
 
+    def merge_generators(self, inplace: bool = False) -> Optional['CaseStudy']:
+        """
+        Merge generators of the same technology at the same bus into one representative generator.
+        Affects dPower_ThermalGen, dPower_VRES, dPower_VRESProfiles, and dPower_Inflows.
+        The new generator ID is '{i}_{tec}'.
+
+        :param inplace: If True, modifies the current instance. If False, returns a new instance.
+        :return: None if inplace is True, otherwise a new CaseStudy instance.
+        """
+        cs = self if inplace else self.copy()
+
+        # Save original VRES mapping before any merges (needed for VRESProfiles and Inflows weighting)
+        original_vres_info = None
+        if hasattr(cs, 'dPower_VRES') and cs.dPower_VRES is not None and 'MaxProd' in cs.dPower_VRES.columns:
+            original_vres_info = cs.dPower_VRES[['tec', 'i', 'MaxProd']].copy()
+
+        ### Merge dPower_ThermalGen
+        if hasattr(cs, 'dPower_ThermalGen') and cs.dPower_ThermalGen is not None:
+            df = cs.dPower_ThermalGen.reset_index()
+            groups = ['tec', 'i']
+
+            thermal_simple_agg = {
+                'ExisUnits': 'max',
+                'MaxProd': 'sum',
+                'MinProd': 'min',
+                'RampUp': 'sum',
+                'RampDw': 'sum',
+                'MinUpTime': 'min',
+                'MinDownTime': 'min',
+                'Qmax': 'sum',
+                'Qmin': 'sum',
+                'EnableInvest': 'max',
+                'YearCom': 'min',
+                'YearDecom': 'max',
+                'lat': 'mean',
+                'lon': 'mean',
+            }
+            thermal_weighted_cols = ['InertiaConst', 'FuelCost', 'Efficiency', 'CommitConsumption',
+                                     'OMVarCost', 'StartupConsumption', 'EFOR', 'InvestCost',
+                                     'FirmCapCoef', 'CO2Emis']
+
+            agg_dict = {}
+            skip_cols = set(groups + ['g'] + thermal_weighted_cols)
+            for col in df.columns:
+                if col in skip_cols:
+                    continue
+                agg_dict[col] = thermal_simple_agg.get(col, 'first')
+
+            merged = df.groupby(groups).agg(agg_dict).reset_index()
+
+            for col in thermal_weighted_cols:
+                if col not in df.columns:
+                    continue
+                numer = (df[col] * df['MaxProd']).groupby([df['tec'], df['i']]).sum()
+                denom = df['MaxProd'].groupby([df['tec'], df['i']]).sum()
+                wavg = (numer / denom.replace(0, np.nan)).fillna(df.groupby(groups)[col].mean())
+                wavg.name = col
+                merged = merged.merge(wavg.reset_index(), on=groups, how='left')
+
+            merged['g'] = merged['i'] + '_' + merged['tec']
+            cs.dPower_ThermalGen = merged.set_index('g')
+
+        ### Merge dPower_VRESProfiles (before dPower_VRES so original MaxProd weights are available)
+        if (hasattr(cs, 'dPower_VRESProfiles') and cs.dPower_VRESProfiles is not None
+                and original_vres_info is not None):
+            df = cs.dPower_VRESProfiles.reset_index()
+            vres_cols = original_vres_info.reset_index()[['g', 'tec', 'i', 'MaxProd']]
+            df = df.merge(vres_cols, on='g', how='left')
+
+            groups = ['rp', 'k', 'scenario', 'tec', 'i']
+            key = [df['rp'], df['k'], df['scenario'], df['tec'], df['i']]
+
+            numer = (df['value'] * df['MaxProd']).groupby(key).sum()
+            denom = df['MaxProd'].groupby(key).sum()
+            merged_value = (numer / denom.replace(0, np.nan)).fillna(df.groupby(groups)['value'].mean())
+            merged_value.name = 'value'
+
+            meta_cols = [c for c in ['dataPackage', 'dataSource', 'id'] if c in df.columns]
+            meta = df.groupby(groups)[meta_cols].first().reset_index()
+            merged = meta.merge(merged_value.reset_index(), on=groups, how='left')
+            merged['g'] = merged['i'] + '_' + merged['tec']
+            merged = merged.drop(columns=['tec', 'i'])
+            cs.dPower_VRESProfiles = merged.set_index(['rp', 'k', 'g'])
+
+        ### Merge dPower_Inflows
+        if (hasattr(cs, 'dPower_Inflows') and cs.dPower_Inflows is not None
+                and original_vres_info is not None):
+            df = cs.dPower_Inflows.reset_index()
+            vres_cols = original_vres_info.reset_index()[['g', 'tec', 'i']]
+            df = df.merge(vres_cols, on='g', how='left')
+
+            groups = ['rp', 'k', 'scenario', 'tec', 'i']
+            key = [df['rp'], df['k'], df['scenario'], df['tec'], df['i']]
+
+            merged_value = df['value'].groupby(key).sum()
+            merged_value.name = 'value'
+
+            meta_cols = [c for c in ['dataPackage', 'dataSource', 'id'] if c in df.columns]
+            meta = df.groupby(groups)[meta_cols].first().reset_index()
+            merged = meta.merge(merged_value.reset_index(), on=groups, how='left')
+            merged['g'] = merged['i'] + '_' + merged['tec']
+            merged = merged.drop(columns=['tec', 'i'])
+            cs.dPower_Inflows = merged.set_index(['rp', 'k', 'g'])
+
+        ### Merge dPower_VRES (last, after VRESProfiles and Inflows)
+        if hasattr(cs, 'dPower_VRES') and cs.dPower_VRES is not None:
+            df = cs.dPower_VRES.reset_index()
+            groups = ['tec', 'i']
+
+            vres_simple_agg = {
+                'ExisUnits': 'sum',
+                'EnableInvest': 'max',
+                'Qmax': 'sum',
+                'Qmin': 'sum',
+                'YearCom': 'min',
+                'YearDecom': 'max',
+                'lat': 'mean',
+                'lon': 'mean',
+            }
+            vres_weighted_cols = ['InvestCost', 'OMVarCost', 'FirmCapCoef', 'InertiaConst']
+            special_cols = {'MaxProd', 'MaxInvest'}
+
+            agg_dict = {}
+            skip_cols = set(groups + ['g'] + vres_weighted_cols + list(special_cols))
+            for col in df.columns:
+                if col in skip_cols:
+                    continue
+                agg_dict[col] = vres_simple_agg.get(col, 'first')
+
+            merged = df.groupby(groups).agg(agg_dict).reset_index()
+
+            # Special: newMaxProd = sum(ExisUnits * MaxProd) / sum(ExisUnits); fallback to sum when all units are greenfield
+            if 'MaxProd' in df.columns:
+                total_mw = (df['ExisUnits'] * df['MaxProd']).groupby([df['tec'], df['i']]).sum()
+                total_units = df['ExisUnits'].groupby([df['tec'], df['i']]).sum()
+                new_maxprod = (total_mw / total_units.replace(0, np.nan)).fillna(
+                    df['MaxProd'].groupby([df['tec'], df['i']]).sum()
+                )
+                new_maxprod.name = 'MaxProd'
+                merged = merged.merge(new_maxprod.reset_index(), on=groups, how='left')
+
+            # Special: newMaxInvest = sum(MaxInvest * MaxProd) / newMaxProd
+            if 'MaxInvest' in df.columns and 'MaxProd' in df.columns:
+                invest_mw = (df['MaxInvest'] * df['MaxProd']).groupby([df['tec'], df['i']]).sum()
+                new_maxprod_s = merged.set_index(groups)['MaxProd']
+                new_maxinvest = (invest_mw / new_maxprod_s.replace(0, np.nan)).fillna(0)
+                new_maxinvest.name = 'MaxInvest'
+                merged = merged.merge(new_maxinvest.reset_index(), on=groups, how='left')
+
+            for col in vres_weighted_cols:
+                if col not in df.columns:
+                    continue
+                numer = (df[col] * df['MaxProd']).groupby([df['tec'], df['i']]).sum()
+                denom = df['MaxProd'].groupby([df['tec'], df['i']]).sum()
+                wavg = (numer / denom.replace(0, np.nan)).fillna(df.groupby(groups)[col].mean())
+                wavg.name = col
+                merged = merged.merge(wavg.reset_index(), on=groups, how='left')
+
+            merged['g'] = merged['i'] + '_' + merged['tec']
+            cs.dPower_VRES = merged.set_index('g')
+
+        return None if inplace else cs
+
     # Create transition matrix from Hindex
     def get_rpTransitionMatrices(self, clip_method: str = "none", clip_value: float = 0) -> tuple[pd.DataFrame, pd.DataFrame, pd.DataFrame]:
         rps = sorted(self.dPower_Hindex.index.get_level_values('rp').unique().tolist())

From 8d972dac90cd8b82007bf179415645328ac8d0aa Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Thu, 30 Apr 2026 14:48:25 +0200
Subject: [PATCH 12/29] Add filter_zone to CaseStudy

---
 CLAUDE.md    |  1 +
 CaseStudy.py | 75 ++++++++++++++++++++++++++++++++++++++++++++++++++++
 README.md    | 22 +++++++++------
 3 files changed, 90 insertions(+), 8 deletions(-)

diff --git a/CLAUDE.md b/CLAUDE.md
index 6767a99..fd88fb3 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -12,6 +12,7 @@ See `README.md` for usage, key concepts, and data structure.
 - `dPower_WeightsRP` is **computed** from `dPower_Hindex` (counting occurrences per `rp`); if a `Power_WeightsRP.xlsx` file also exists, it is read and compared — a mismatch triggers a warning but uses the **file** value, not the computed one.
 - `merge_single_node_buses()` preserves the `z` (zone) column as a sorted unique union string of all merged zones (e.g. `"R1_R2"`). This is documented in root `CLAUDE.md` as well.
 - `merge_generators()` collapses all generators sharing the same `(tec, i)` into one representative generator with ID `"{i}_{tec}"`. It must be called after construction (scaling is not required). VRESProfiles and Inflows are merged **before** VRES so that the original per-generator `MaxProd` weights are available for the capacity-factor weighted average. Generators in VRESProfiles/Inflows that have no matching entry in `dPower_VRES` (left-join miss) are grouped under `(tec=NaN, i=NaN)` — filter to a single scenario first to avoid mixing scenarios in the groupby.
+- `filter_zone(zone)` keeps only buses whose `z` column exactly matches the given zone name(s). Cascade order: BusInfo → Network (both endpoints must survive) → ThermalGen / VRES / Storage (by `i`) → Demand (by `i` index level) → VRESProfiles (by surviving VRES gen IDs) → Inflows (by union of surviving VRES + Storage gen IDs) → ImportExport (by `i` index level). `z` is an exact-match filter; merged-bus zone strings like `"R1_R2"` are not matched by `"R1"` alone.
 - `CaseStudy.copy()` is a full `deepcopy` — safe to modify independently.
 - Transition matrices (`rpTransitionMatrixAbsolute`, `rpTransitionMatrixRelativeTo`, `rpTransitionMatrixRelativeFrom`) are computed in the constructor and attached as attributes.
 
diff --git a/CaseStudy.py b/CaseStudy.py
index 9fccb0c..6770e68 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -925,6 +925,81 @@ def filter_scenario(self, scenario_name, inplace: bool = False) -> Optional[Self
 
         return None if inplace else caseStudy
 
+    def filter_zone(self, zone: str | list[str], inplace: bool = False) -> Optional[Self]:
+        """
+        Filters the case study to only include buses in the given zone(s). All generators,
+        network lines, demand entries, and time series profiles connected to buses outside the
+        zone are removed.
+        :param zone: Zone name (value of the 'z' column in Power_BusInfo) or list of zone names to keep.
+        :param inplace: If True, modifies the current instance. If False, returns a new instance.
+        :return: None if inplace is True, otherwise a new CaseStudy instance.
+        """
+        case_study = self if inplace else self.copy()
+
+        zones = [zone] if isinstance(zone, str) else list(zone)
+
+        # Filter BusInfo and derive remaining bus set
+        case_study.dPower_BusInfo = case_study.dPower_BusInfo[case_study.dPower_BusInfo['z'].isin(zones)]
+        remaining_buses = set(case_study.dPower_BusInfo.index)
+
+        # Filter Network: drop lines where either endpoint is outside the zone
+        network_reset = case_study.dPower_Network.reset_index()
+        case_study.dPower_Network = network_reset[
+            network_reset['i'].isin(remaining_buses) & network_reset['j'].isin(remaining_buses)
+            ].set_index(['i', 'j', 'c'])
+
+        # Filter ThermalGen
+        if hasattr(case_study, 'dPower_ThermalGen') and case_study.dPower_ThermalGen is not None:
+            case_study.dPower_ThermalGen = case_study.dPower_ThermalGen[
+                case_study.dPower_ThermalGen['i'].isin(remaining_buses)
+            ]
+
+        # Filter VRES; collect remaining VRES generator IDs for VRESProfiles / Inflows
+        remaining_vres_gens: set = set()
+        if hasattr(case_study, 'dPower_VRES') and case_study.dPower_VRES is not None:
+            case_study.dPower_VRES = case_study.dPower_VRES[
+                case_study.dPower_VRES['i'].isin(remaining_buses)
+            ]
+            remaining_vres_gens = set(case_study.dPower_VRES.index)
+
+        # Filter Storage; collect remaining storage generator IDs for Inflows
+        remaining_storage_gens: set = set()
+        if hasattr(case_study, 'dPower_Storage') and case_study.dPower_Storage is not None:
+            case_study.dPower_Storage = case_study.dPower_Storage[
+                case_study.dPower_Storage['i'].isin(remaining_buses)
+            ]
+            remaining_storage_gens = set(case_study.dPower_Storage.index)
+
+        # Filter Demand
+        demand_reset = case_study.dPower_Demand.reset_index()
+        case_study.dPower_Demand = demand_reset[
+            demand_reset['i'].isin(remaining_buses)
+        ].set_index(['rp', 'k', 'i'])
+
+        # Filter VRESProfiles by remaining VRES generator IDs
+        if hasattr(case_study, 'dPower_VRESProfiles') and case_study.dPower_VRESProfiles is not None:
+            profiles_reset = case_study.dPower_VRESProfiles.reset_index()
+            case_study.dPower_VRESProfiles = profiles_reset[
+                profiles_reset['g'].isin(remaining_vres_gens)
+            ].set_index(['rp', 'k', 'g'])
+
+        # Filter Inflows by remaining VRES + Storage generator IDs
+        if hasattr(case_study, 'dPower_Inflows') and case_study.dPower_Inflows is not None:
+            remaining_gens = remaining_vres_gens | remaining_storage_gens
+            inflows_reset = case_study.dPower_Inflows.reset_index()
+            case_study.dPower_Inflows = inflows_reset[
+                inflows_reset['g'].isin(remaining_gens)
+            ].set_index(['rp', 'k', 'g'])
+
+        # Filter ImportExport by remaining buses
+        if hasattr(case_study, 'dPower_ImportExport') and case_study.dPower_ImportExport is not None:
+            ie_reset = case_study.dPower_ImportExport.reset_index()
+            case_study.dPower_ImportExport = ie_reset[
+                ie_reset['i'].isin(remaining_buses)
+            ].set_index(['hub', 'i', 'rp', 'k'])
+
+        return None if inplace else case_study
+
     def filter_timesteps(self, start: str, end: str, inplace: bool = False, no_weight_k_adjustment: bool = False) -> Optional[Self]:
         """
         Filters each (relevant) dataframe in the case study to only include the timesteps between start and end (both inclusive).
diff --git a/README.md b/README.md
index 43703ea..9670220 100644
--- a/README.md
+++ b/README.md
@@ -8,7 +8,7 @@ Data I/O package for LEGO-Pyomo. Handles reading Excel case study files, writing
 |------------------------|-------------------------------------------------------------------------|
 | `CaseStudy.py`         | Loads all Excel input files into a single object; data manipulation     |
 | `ExcelReader.py`       | Low-level Excel parsing (per-file reader functions, version checking)   |
-| `ExcelWriter.py`       | Writes formatted Excel output files, driven by `TableDefinitions.xml`  |
+| `ExcelWriter.py`       | Writes formatted Excel output files, driven by `TableDefinitions.xml`   |
 | `SQLiteWriter.py`      | Exports Pyomo model results to SQLite; stores solver stats & run params |
 | `Utilities.py`         | Data transformations (inflows ↔ capacity factors, Printer helpers)      |
 | `printer.py`           | Singleton console/logfile printer with severity levels                  |
@@ -58,13 +58,19 @@ All Excel files in `data_folder` are read automatically. Any `d*` parameter can
 
 ### Key Methods
 
-| Method                         | Description                                              |
-|--------------------------------|----------------------------------------------------------|
-| `copy()`                       | Deep copy — safe to modify independently                 |
-| `equal_to(cs)`                 | Compare all DataFrames with another CaseStudy            |
-| `merge_single_node_buses()`    | Collapse single-bus zones; preserves `z` as union string |
-| `scale_CaseStudy()`            | Applies power and cost scaling factors from parameters   |
-| `get_rpTransitionMatrices()`   | Returns absolute and relative transition matrices        |
+| Method                              | Description                                              |
+|-------------------------------------|----------------------------------------------------------|
+| `copy()`                            | Deep copy — safe to modify independently                 |
+| `equal_to(cs)`                      | Compare all DataFrames with another CaseStudy            |
+| `merge_single_node_buses()`         | Collapse single-bus zones; preserves `z` as union string |
+| `merge_generators()`                | Merge generators of the same `(tec, i)` into one         |
+| `scale_CaseStudy()`                 | Applies power and cost scaling factors from parameters   |
+| `get_rpTransitionMatrices()`        | Returns absolute and relative transition matrices        |
+| `filter_scenario(name)`             | Keep only rows matching the given scenario               |
+| `filter_zone(zone)`                 | Keep only buses (and all connected data) in a zone       |
+| `filter_timesteps(start, end)`      | Keep only timesteps in the given `k` range               |
+| `filter_representative_periods(rp)` | Keep only one representative period                      |
+| `apply_kmedoids_aggregation()`      | Temporally aggregate using k-medoids clustering          |
 
 ## SQLiteWriter
 

From 93c19f23a639c9770155da7e4240f81d87c227e3 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Mon, 4 May 2026 13:08:44 +0200
Subject: [PATCH 13/29] Fix clustering if k doesn't start at k0001

---
 Utilities.py | 4 +++-
 1 file changed, 3 insertions(+), 1 deletion(-)

diff --git a/Utilities.py b/Utilities.py
index 2585430..48d3aa7 100644
--- a/Utilities.py
+++ b/Utilities.py
@@ -321,7 +321,9 @@ def _extract_numeric_and_calc_p(df, rp_length):
         """Extract numeric values from rp/k strings and calculate absolute hour."""
         df['rp_num'] = df['rp'].str[2:].astype(int)
         df['k_num'] = df['k'].str[1:].astype(int)
-        df['p'] = (df['rp_num'] - 1) * rp_length + df['k_num']
+        # Normalize k so the first k value maps to 1, regardless of offset (e.g. --limitK k2161-k4320)
+        min_k_num = df['k_num'].min()
+        df['p'] = (df['rp_num'] - 1) * rp_length + (df['k_num'] - min_k_num + 1)
         return df
 
     time_series_tables = [("Power_Demand", case_study.dPower_Demand)]

From 3411ff0adb638d8c3b0f5637c41c5741050c3ce6 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Fri, 8 May 2026 15:04:04 +0200
Subject: [PATCH 14/29] Fix bug in caller when replacement character is printed

---
 Caller.py | 12 ++++++------
 1 file changed, 6 insertions(+), 6 deletions(-)

diff --git a/Caller.py b/Caller.py
index 20ef229..30b826d 100644
--- a/Caller.py
+++ b/Caller.py
@@ -13,7 +13,7 @@
 def _tail_file(filepath, n=20):
     """Return the last n lines of a file, or all if fewer."""
     try:
-        with open(filepath, 'r', errors='replace') as f:
+        with open(filepath, 'r', encoding='utf-8', errors='replace') as f:
             lines = f.readlines()
     except OSError:
         return "(could not read log file)"
@@ -61,7 +61,7 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
 
 
 while True:
-    with open(args.jobs, 'r') as f:
+    with open(args.jobs, 'r', encoding='utf-8') as f:
         lines = f.readlines()
 
     found_one = False
@@ -106,7 +106,7 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
         if barrier_waited:
             barrier_waited.clear()
         start_datetime = datetime.datetime.now().strftime("%Y-%m-%d %H:%M:%S")
-        with os.fdopen(fd, 'w') as f:
+        with os.fdopen(fd, 'w', encoding='utf-8') as f:
             f.write(f"Command: {line.strip()}\n")
             f.write(f"Started at: {start_datetime}")
         found_one = True
@@ -116,7 +116,7 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
             os.system(f"title Job {i} from '{args.jobs}': {line.strip()}")
 
             start_time = time.time()
-            with open(log_file, 'w') as log_f:
+            with open(log_file, 'w', encoding='utf-8') as log_f:
                 log_f.write(f"Command: {line.strip()}\n")
                 log_f.write(f"Started at:  {start_datetime}\n")
                 log_f.write(f"{'=' * 60}\n")
@@ -160,7 +160,7 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
                     f"Last output:\n{_tail_file(log_file, 20)}"
                 )
 
-            with open(finished_job_flag, 'w') as f:
+            with open(finished_job_flag, 'w', encoding='utf-8') as f:
                 f.write(f"Command: {line.strip()}\n")
                 f.write(f"Started at:  {start_datetime}\n")
                 f.write(f"Finished at: {datetime.datetime.now().strftime('%Y-%m-%d %H:%M:%S')}\n")
@@ -170,7 +170,7 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
             printer.information(f"Finished job {i} from '{args.jobs}' after {end_time - start_time:.2f} seconds (= {(end_time - start_time) / 60 / 60:.2f} hours).")
         except Exception as e:
             printer.error(f"Error while executing job {i}: {e}")
-            with open(error_job_flag, 'w') as f:
+            with open(error_job_flag, 'w', encoding='utf-8') as f:
                 f.write(f"Command: {line.strip()}\n")
                 f.write(f"Error while executing job {i} from '{args.jobs}': {e}\n")
                 f.write(f"Started at:  {start_datetime}\n")

From 12399c980462a7d1532d43fba2c1c014381a56cd Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Mon, 11 May 2026 21:04:54 +0200
Subject: [PATCH 15/29] Adjust CaseStudy.to_full_hourly_model to work with
 multiple scenarios

Fix type-annotation for parameter files
---
 CLAUDE.md    |   1 +
 CaseStudy.py | 153 +++++++++++++++++++++++++++++++--------------------
 2 files changed, 95 insertions(+), 59 deletions(-)

diff --git a/CLAUDE.md b/CLAUDE.md
index fd88fb3..880fe18 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -14,6 +14,7 @@ See `README.md` for usage, key concepts, and data structure.
 - `merge_generators()` collapses all generators sharing the same `(tec, i)` into one representative generator with ID `"{i}_{tec}"`. It must be called after construction (scaling is not required). VRESProfiles and Inflows are merged **before** VRES so that the original per-generator `MaxProd` weights are available for the capacity-factor weighted average. Generators in VRESProfiles/Inflows that have no matching entry in `dPower_VRES` (left-join miss) are grouped under `(tec=NaN, i=NaN)` — filter to a single scenario first to avoid mixing scenarios in the groupby.
 - `filter_zone(zone)` keeps only buses whose `z` column exactly matches the given zone name(s). Cascade order: BusInfo → Network (both endpoints must survive) → ThermalGen / VRES / Storage (by `i`) → Demand (by `i` index level) → VRESProfiles (by surviving VRES gen IDs) → Inflows (by union of surviving VRES + Storage gen IDs) → ImportExport (by `i` index level). `z` is an exact-match filter; merged-bus zone strings like `"R1_R2"` are not matched by `"R1"` alone.
 - `CaseStudy.copy()` is a full `deepcopy` — safe to modify independently.
+- `to_full_hourly_model()` processes each enabled scenario from `dGlobal_Scenarios` independently. It is a no-op (returns unchanged) if all scenarios are already hourly (no two p-values share the same 'k' within any scenario). The produced Hindex uses consecutive h0001…hN p-labels and k0001…kN k-labels per scenario, all mapped to rp01 with weight 1.
 - Transition matrices (`rpTransitionMatrixAbsolute`, `rpTransitionMatrixRelativeTo`, `rpTransitionMatrixRelativeFrom`) are computed in the constructor and attached as attributes.
 
 ### ExcelReader
diff --git a/CaseStudy.py b/CaseStudy.py
index 6770e68..1d88a70 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -46,9 +46,9 @@ def __init__(self,
                  do_not_merge_single_node_buses: bool = False,
                  parallel_read: bool = True,
                  n_jobs: int = 4,
-                 global_parameters_file: str = "Global_Parameters.xlsx", dGlobal_Parameters: pd.DataFrame = None,
+                 global_parameters_file: str = "Global_Parameters.xlsx", dGlobal_Parameters: dict = None,
                  global_scenarios_file: str = "Global_Scenarios.xlsx", dGlobal_Scenarios: pd.DataFrame = None,
-                 power_parameters_file: str = "Power_Parameters.xlsx", dPower_Parameters: pd.DataFrame = None,
+                 power_parameters_file: str = "Power_Parameters.xlsx", dPower_Parameters: dict = None,
                  power_businfo_file: str = "Power_BusInfo.xlsx", dPower_BusInfo: pd.DataFrame = None,
                  power_network_file: str = "Power_Network.xlsx", dPower_Network: pd.DataFrame = None,
                  power_thermalgen_file: str = "Power_ThermalGen.xlsx", dPower_ThermalGen: pd.DataFrame = None,
@@ -838,68 +838,103 @@ def to_full_hourly_model(self, inplace: bool) -> Optional['CaseStudy']:
         or return a new `CaseStudy` instance if `inplace` is `False`. The adjustments align the data
         to represent hourly indices and corresponding weights.
 
+        Each enabled scenario in `dGlobal_Scenarios` is processed independently.
+
+        If the case study is already a full hourly model (no two p-values share the same 'k' within
+        any scenario), no adjustment is made and the original instance is returned unchanged.
+
         :param inplace: If `True`, modifies the given instance. If `False`, returns a new `CaseStudy` instance.
         :return: Adjusted `CaseStudy` instance if `inplace` is `False`, otherwise `None`.
         """
         caseStudy = self.copy() if not inplace else self
 
-        # First Adjustment of Hindex (important if the case study was filtered before, to get a coherent p-index)
-        caseStudy.dPower_Hindex = caseStudy.dPower_Hindex.reset_index()
-        for i in caseStudy.dPower_Hindex.index:
-            caseStudy.dPower_Hindex.loc[i, "p"] = f"h{i + 1:0>4}"
-        caseStudy.dPower_Hindex = caseStudy.dPower_Hindex.set_index(["p", "rp", "k"])
-
-        # Adjust Demand
-        adjusted_demand = []
-        for i in caseStudy.dPower_BusInfo.index:
-            for h in caseStudy.dPower_Hindex.index:
-                adjusted_demand.append(["rp01", h[0].replace("h", "k"), i, caseStudy.dPower_Demand.loc[(h[1], h[2], i), "value"], "ScenarioA", None, None, None])
-
-        caseStudy.dPower_Demand = pd.DataFrame(adjusted_demand, columns=["rp", "k", "i", "value", "scenario", "id", "dataPackage", "dataSource"])
-        caseStudy.dPower_Demand = caseStudy.dPower_Demand.set_index(["rp", "k", "i"])
-
-        # Adjust VRESProfiles
-        if hasattr(caseStudy, "dPower_VRESProfiles"):
-            adjusted_vresprofiles = []
-            caseStudy.dPower_VRESProfiles.sort_index(inplace=True)
-            for g in caseStudy.dPower_VRESProfiles.index.get_level_values('g').unique().tolist():
-                for h in caseStudy.dPower_Hindex.index:
-                    adjusted_vresprofiles.append(["rp01", h[0].replace("h", "k"), g, caseStudy.dPower_VRESProfiles.loc[(h[1], h[2], g), "value"], "ScenarioA", None, None, None])
-
-            caseStudy.dPower_VRESProfiles = pd.DataFrame(adjusted_vresprofiles, columns=["rp", "k", "g", "value", "scenario", "id", "dataPackage", "dataSource"])
-            caseStudy.dPower_VRESProfiles = caseStudy.dPower_VRESProfiles.set_index(["rp", "k", "g"])
-
-        # Adjust Inflows
-        if hasattr(caseStudy, "dPower_Inflows"):
-            adjusted_inflows = []
-            caseStudy.dPower_Inflows.sort_index(inplace=True)
-            for g in caseStudy.dPower_Inflows.index.get_level_values('g').unique().tolist():
-                for h in caseStudy.dPower_Hindex.index:
-                    adjusted_inflows.append(["rp01", h[0].replace("h", "k"), g, caseStudy.dPower_Inflows.loc[(h[1], h[2], g), "value"], "ScenarioA", None, None, None])
-            caseStudy.dPower_Inflows = pd.DataFrame(adjusted_inflows, columns=["rp", "k", "g", "value", "scenario", "id", "dataPackage", "dataSource"])
-            caseStudy.dPower_Inflows = caseStudy.dPower_Inflows.set_index(["rp", "k", "g"])
-
-        # Adjust Hindex
-        caseStudy.dPower_Hindex = caseStudy.dPower_Hindex.reset_index()
-        for i in caseStudy.dPower_Hindex.index:
-            caseStudy.dPower_Hindex.loc[i] = f"h{i + 1:0>4}", f"rp01", f"k{i + 1:0>4}", None, None, None, "ScenarioA"
-        caseStudy.dPower_Hindex = caseStudy.dPower_Hindex.set_index(["p", "rp", "k"])
-
-        # Adjust WeightsK
-        caseStudy.dPower_WeightsK = caseStudy.dPower_WeightsK.reset_index()
-        caseStudy.dPower_WeightsK = caseStudy.dPower_WeightsK.drop(caseStudy.dPower_WeightsK.index)
-        for i in range(len(caseStudy.dPower_Hindex)):
-            caseStudy.dPower_WeightsK.loc[i] = f"{caseStudy.dPower_Hindex.index[i][2]}", None, 1, None, None, "ScenarioA"
-        caseStudy.dPower_WeightsK = caseStudy.dPower_WeightsK.set_index("k")
-
-        # Adjust WeightsRP
-        caseStudy.dPower_WeightsRP = caseStudy.dPower_WeightsRP.drop(caseStudy.dPower_WeightsRP.index)
-        caseStudy.dPower_WeightsRP.loc["rp01"] = None, 1, None, None, "ScenarioA"
-
-        if not inplace:
-            return caseStudy
-        else:
-            return None
+        # Check if already hourly: within each scenario no 'k' appears more than once
+        hindex_flat = caseStudy.dPower_Hindex.reset_index()
+        already_hourly = all(
+            not (grp.groupby(['k']).size() > 1).any()
+            for _, grp in hindex_flat.groupby('scenario')
+        )
+        if already_hourly:
+            return None if inplace else caseStudy
+
+        scenario_names = caseStudy.dGlobal_Scenarios.index.tolist()
+
+        all_demand = []
+        all_vresprofiles = []
+        all_inflows = []
+        all_hindex = []
+        all_weightsk = []
+        all_weightsrp = []
+
+        demand_flat = caseStudy.dPower_Demand.reset_index()
+        vres_flat = caseStudy.dPower_VRESProfiles.reset_index() if hasattr(caseStudy, 'dPower_VRESProfiles') and caseStudy.dPower_VRESProfiles is not None else None
+        inflows_flat = caseStudy.dPower_Inflows.reset_index() if hasattr(caseStudy, 'dPower_Inflows') and caseStudy.dPower_Inflows is not None else None
+
+        for scenario in scenario_names:
+            scen_hindex = hindex_flat[hindex_flat['scenario'] == scenario].copy().reset_index(drop=True)
+            num_hours = len(scen_hindex)
+            scen_hindex['new_k'] = [f"k{i + 1:0>4}" for i in range(num_hours)]
+            scen_hindex['new_p'] = [f"h{i + 1:0>4}" for i in range(num_hours)]
+
+            # Demand: merge each original (rp, k) with the buses that have that (rp, k)
+            demand_scen = demand_flat[demand_flat['scenario'] == scenario][['rp', 'k', 'i', 'value']]
+            m = scen_hindex[['rp', 'k', 'new_k']].merge(demand_scen, on=['rp', 'k'])
+            all_demand.append(pd.DataFrame({'rp': 'rp01', 'k': m['new_k'], 'i': m['i'], 'value': m['value'],
+                                            'scenario': scenario, 'id': None, 'dataPackage': None, 'dataSource': None}))
+
+            # VRESProfiles
+            if vres_flat is not None:
+                vres_scen = vres_flat[vres_flat['scenario'] == scenario][['rp', 'k', 'g', 'value']]
+                m = scen_hindex[['rp', 'k', 'new_k']].merge(vres_scen, on=['rp', 'k'])
+                all_vresprofiles.append(pd.DataFrame({'rp': 'rp01', 'k': m['new_k'], 'g': m['g'], 'value': m['value'],
+                                                      'scenario': scenario, 'id': None, 'dataPackage': None, 'dataSource': None}))
+
+            # Inflows
+            if inflows_flat is not None:
+                inflows_scen = inflows_flat[inflows_flat['scenario'] == scenario][['rp', 'k', 'g', 'value']]
+                m = scen_hindex[['rp', 'k', 'new_k']].merge(inflows_scen, on=['rp', 'k'])
+                all_inflows.append(pd.DataFrame({'rp': 'rp01', 'k': m['new_k'], 'g': m['g'], 'value': m['value'],
+                                                 'scenario': scenario, 'id': None, 'dataPackage': None, 'dataSource': None}))
+
+            # Hindex
+            hindex_scen = pd.DataFrame({
+                'p': scen_hindex['new_p'],
+                'rp': 'rp01',
+                'k': scen_hindex['new_k'],
+                'id': None,
+                'dataPackage': None,
+                'dataSource': None,
+                'scenario': scenario,
+            })
+            all_hindex.append(hindex_scen)
+
+            # WeightsK
+            weightsk_scen = pd.DataFrame({
+                'k': scen_hindex['new_k'],
+                'id': None,
+                'pWeight_k': 1,
+                'dataPackage': None,
+                'dataSource': None,
+                'scenario': scenario,
+            })
+            all_weightsk.append(weightsk_scen)
+
+            # WeightsRP
+            all_weightsrp.append({'rp': 'rp01', 'id': None, 'pWeight_rp': 1, 'dataPackage': None, 'dataSource': None, 'scenario': scenario})
+
+        caseStudy.dPower_Demand = pd.concat(all_demand).set_index(['rp', 'k', 'i'])
+
+        if all_vresprofiles:
+            caseStudy.dPower_VRESProfiles = pd.concat(all_vresprofiles).set_index(['rp', 'k', 'g'])
+
+        if all_inflows:
+            caseStudy.dPower_Inflows = pd.concat(all_inflows).set_index(['rp', 'k', 'g'])
+
+        caseStudy.dPower_Hindex = pd.concat(all_hindex).set_index(['p', 'rp', 'k'])
+        caseStudy.dPower_WeightsK = pd.concat(all_weightsk).set_index('k')
+        caseStudy.dPower_WeightsRP = pd.DataFrame(all_weightsrp).set_index('rp')
+
+        return None if inplace else caseStudy
 
     def filter_scenario(self, scenario_name, inplace: bool = False) -> Optional[Self]:
         """

From 47533035bce672ea9f32494e7d17e6297526f08d Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 19 May 2026 15:54:11 +0200
Subject: [PATCH 16/29] Fix issues with potential non-str indices in dataframes

---
 ExcelReader.py | 3 +++
 1 file changed, 3 insertions(+)

diff --git a/ExcelReader.py b/ExcelReader.py
index 6ae4ab9..990a397 100644
--- a/ExcelReader.py
+++ b/ExcelReader.py
@@ -66,6 +66,8 @@ def __read_non_pivoted_file(excel_file_path: str, version_specifier: str, indice
                 df = df[df["excl"].isnull()]  # Only keep rows that are not excluded (i.e., have no value in the "Excl." column)
         else:
             df = df.drop(df.columns[0], axis=1)  # Drop the first column (which is empty)
+
+        df[indices] = df[indices].astype(str)  # Convert index columns to string to prevent issues with integer-indices
         df = df.set_index(indices) if len(indices) > 0 else df
         df["scenario"] = scenario
 
@@ -91,6 +93,7 @@ def __read_pivoted_file(excel_file_path: str, version_specifier: str, indices: l
     df = __read_non_pivoted_file(excel_file_path, version_specifier, [], has_excl_column, keep_excluded_columns, fail_on_wrong_version)
 
     df = df.melt(id_vars=melt_indices + ["scenario"], var_name=pivoted_variable_name, value_name="value")
+    df[indices] = df[indices].astype(str)  # Convert index columns to string to prevent issues with integer-indices (this has to be done here, since we don't pass indices to __read_non_pivoted_file)
     df = df.set_index(indices)
     return df
 

From 9dfc31217b69f8f75e5356ec9bc3a21238620590 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 19 May 2026 16:36:55 +0200
Subject: [PATCH 17/29] Fix issue with non-str indices

---
 Utilities.py | 1 +
 1 file changed, 1 insertion(+)

diff --git a/Utilities.py b/Utilities.py
index 48d3aa7..444e114 100644
--- a/Utilities.py
+++ b/Utilities.py
@@ -224,6 +224,7 @@ def _pivot_technologies(df, value_column, index_cols=None):
     if combined_tech_data is not None:
         # Use right join to keep ALL demand data (even nodes without technology data)
         # Replicates demand for nodes with technology, and preserves demand-only nodes
+        combined_tech_data['i'] = combined_tech_data['i'].astype(str)
         scenario_df = pd.merge(
             combined_tech_data,
             scenario_df,

From dc8cb12a50a34c0c605888ee0662afd61c8e9c11 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 19 May 2026 17:43:44 +0200
Subject: [PATCH 18/29] Filter unused scenarios per default

---
 CaseStudy.py | 9 +++++++++
 1 file changed, 9 insertions(+)

diff --git a/CaseStudy.py b/CaseStudy.py
index 1d88a70..43e73f5 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -44,6 +44,7 @@ def __init__(self,
                  data_folder: str | Path,
                  do_not_scale_units: bool = False,
                  do_not_merge_single_node_buses: bool = False,
+                 do_not_filter_unused_scenarios: bool = False,
                  parallel_read: bool = True,
                  n_jobs: int = 4,
                  global_parameters_file: str = "Global_Parameters.xlsx", dGlobal_Parameters: dict = None,
@@ -65,6 +66,7 @@ def __init__(self,
         self.data_folder = str(data_folder) if str(data_folder).endswith("/") else str(data_folder) + "/"
         self.do_not_scale_units = do_not_scale_units
         self.do_not_merge_single_node_buses = do_not_merge_single_node_buses
+        self.do_not_filter_unused_scenarios = do_not_filter_unused_scenarios
 
         # === SEQUENTIAL READS ===
         if dGlobal_Parameters is not None:
@@ -239,6 +241,13 @@ def __init__(self,
                 printer.warning(f"Executing without 'Power_WeightsRP' (since no file was found at '{self.data_folder + self.power_weightsrp_file}').")
                 self.dPower_WeightsRP = dPower_WeightsRP
 
+        if not self.do_not_filter_unused_scenarios:
+            self.dGlobal_Scenarios = self.dGlobal_Scenarios[self.dGlobal_Scenarios['relativeWeight'] != 0]  # Drop rows in dGlobal_Scenarios where relativeWeight is 0
+            if len(self.dGlobal_Scenarios) == 0:
+                raise ValueError("No scenarios are present in the 'Global_Scenarios' table. Please check if the file exists and contains valid data.")
+            elif len(self.dGlobal_Scenarios) == 1:
+                self.filter_scenario(self.dGlobal_Scenarios.index[0], inplace=True)  # Filter case study to only include actually present scenarios
+
         self.rpTransitionMatrixAbsolute, self.rpTransitionMatrixRelativeTo, self.rpTransitionMatrixRelativeFrom = self.get_rpTransitionMatrices(clip_method=clip_method, clip_value=clip_value)
 
         if not do_not_merge_single_node_buses:

From 15ba89985dd94dd61dcf13ea2499b11fc62cbeae Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 19 May 2026 17:44:01 +0200
Subject: [PATCH 19/29] Fix buses with integer names in merge_generators

---
 CaseStudy.py | 8 ++++----
 1 file changed, 4 insertions(+), 4 deletions(-)

diff --git a/CaseStudy.py b/CaseStudy.py
index 43e73f5..7553da9 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -695,7 +695,7 @@ def merge_generators(self, inplace: bool = False) -> Optional['CaseStudy']:
                 wavg.name = col
                 merged = merged.merge(wavg.reset_index(), on=groups, how='left')
 
-            merged['g'] = merged['i'] + '_' + merged['tec']
+            merged['g'] = merged['i'].astype(str) + '_' + merged['tec']
             cs.dPower_ThermalGen = merged.set_index('g')
 
         ### Merge dPower_VRESProfiles (before dPower_VRES so original MaxProd weights are available)
@@ -716,7 +716,7 @@ def merge_generators(self, inplace: bool = False) -> Optional['CaseStudy']:
             meta_cols = [c for c in ['dataPackage', 'dataSource', 'id'] if c in df.columns]
             meta = df.groupby(groups)[meta_cols].first().reset_index()
             merged = meta.merge(merged_value.reset_index(), on=groups, how='left')
-            merged['g'] = merged['i'] + '_' + merged['tec']
+            merged['g'] = merged['i'].astype(str) + '_' + merged['tec']
             merged = merged.drop(columns=['tec', 'i'])
             cs.dPower_VRESProfiles = merged.set_index(['rp', 'k', 'g'])
 
@@ -736,7 +736,7 @@ def merge_generators(self, inplace: bool = False) -> Optional['CaseStudy']:
             meta_cols = [c for c in ['dataPackage', 'dataSource', 'id'] if c in df.columns]
             meta = df.groupby(groups)[meta_cols].first().reset_index()
             merged = meta.merge(merged_value.reset_index(), on=groups, how='left')
-            merged['g'] = merged['i'] + '_' + merged['tec']
+            merged['g'] = merged['i'].astype(str) + '_' + merged['tec']
             merged = merged.drop(columns=['tec', 'i'])
             cs.dPower_Inflows = merged.set_index(['rp', 'k', 'g'])
 
@@ -794,7 +794,7 @@ def merge_generators(self, inplace: bool = False) -> Optional['CaseStudy']:
                 wavg.name = col
                 merged = merged.merge(wavg.reset_index(), on=groups, how='left')
 
-            merged['g'] = merged['i'] + '_' + merged['tec']
+            merged['g'] = merged['i'].astype(str) + '_' + merged['tec']
             cs.dPower_VRES = merged.set_index('g')
 
         return None if inplace else cs

From 820df55c53a2c5552b2836984f1512a50718dcd8 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Wed, 20 May 2026 17:24:23 +0200
Subject: [PATCH 20/29] Fix integer bus-names

---
 CaseStudy.py | 28 ++++++++++++++++++----------
 1 file changed, 18 insertions(+), 10 deletions(-)

diff --git a/CaseStudy.py b/CaseStudy.py
index 7553da9..fbdcc5a 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -983,26 +983,30 @@ def filter_zone(self, zone: str | list[str], inplace: bool = False) -> Optional[
         zones = [zone] if isinstance(zone, str) else list(zone)
 
         # Filter BusInfo and derive remaining bus set
-        case_study.dPower_BusInfo = case_study.dPower_BusInfo[case_study.dPower_BusInfo['z'].isin(zones)]
+        case_study.dPower_BusInfo = case_study.dPower_BusInfo[case_study.dPower_BusInfo['z'].astype(str).isin(zones)]
         remaining_buses = set(case_study.dPower_BusInfo.index)
 
         # Filter Network: drop lines where either endpoint is outside the zone
         network_reset = case_study.dPower_Network.reset_index()
         case_study.dPower_Network = network_reset[
-            network_reset['i'].isin(remaining_buses) & network_reset['j'].isin(remaining_buses)
-            ].set_index(['i', 'j', 'c'])
+            network_reset['i'].astype(str).isin(remaining_buses) & network_reset['j'].astype(str).isin(remaining_buses)
+            ]
+        case_study.dPower_Network['i'] = case_study.dPower_Network['i'].astype(str)
+        case_study.dPower_Network['j'] = case_study.dPower_Network['j'].astype(str)
+        case_study.dPower_Network['c'] = case_study.dPower_Network['c'].astype(str)
+        case_study.dPower_Network.set_index(['i', 'j', 'c'], inplace=True)
 
         # Filter ThermalGen
         if hasattr(case_study, 'dPower_ThermalGen') and case_study.dPower_ThermalGen is not None:
             case_study.dPower_ThermalGen = case_study.dPower_ThermalGen[
-                case_study.dPower_ThermalGen['i'].isin(remaining_buses)
+                case_study.dPower_ThermalGen['i'].astype(str).isin(remaining_buses)
             ]
 
         # Filter VRES; collect remaining VRES generator IDs for VRESProfiles / Inflows
         remaining_vres_gens: set = set()
         if hasattr(case_study, 'dPower_VRES') and case_study.dPower_VRES is not None:
             case_study.dPower_VRES = case_study.dPower_VRES[
-                case_study.dPower_VRES['i'].isin(remaining_buses)
+                case_study.dPower_VRES['i'].astype(str).isin(remaining_buses)
             ]
             remaining_vres_gens = set(case_study.dPower_VRES.index)
 
@@ -1010,15 +1014,17 @@ def filter_zone(self, zone: str | list[str], inplace: bool = False) -> Optional[
         remaining_storage_gens: set = set()
         if hasattr(case_study, 'dPower_Storage') and case_study.dPower_Storage is not None:
             case_study.dPower_Storage = case_study.dPower_Storage[
-                case_study.dPower_Storage['i'].isin(remaining_buses)
+                case_study.dPower_Storage['i'].astype(str).isin(remaining_buses)
             ]
             remaining_storage_gens = set(case_study.dPower_Storage.index)
 
         # Filter Demand
         demand_reset = case_study.dPower_Demand.reset_index()
         case_study.dPower_Demand = demand_reset[
-            demand_reset['i'].isin(remaining_buses)
-        ].set_index(['rp', 'k', 'i'])
+            demand_reset['i'].astype(str).isin(remaining_buses)
+        ]
+        case_study.dPower_Demand['i'] = case_study.dPower_Demand['i'].astype(str)
+        case_study.dPower_Demand.set_index(['rp', 'k', 'i'], inplace=True)
 
         # Filter VRESProfiles by remaining VRES generator IDs
         if hasattr(case_study, 'dPower_VRESProfiles') and case_study.dPower_VRESProfiles is not None:
@@ -1039,8 +1045,10 @@ def filter_zone(self, zone: str | list[str], inplace: bool = False) -> Optional[
         if hasattr(case_study, 'dPower_ImportExport') and case_study.dPower_ImportExport is not None:
             ie_reset = case_study.dPower_ImportExport.reset_index()
             case_study.dPower_ImportExport = ie_reset[
-                ie_reset['i'].isin(remaining_buses)
-            ].set_index(['hub', 'i', 'rp', 'k'])
+                ie_reset['i'].astype(str).isin(remaining_buses)
+            ]
+            case_study.dPower_ImportExport['i'] = case_study.dPower_ImportExport['i'].astype(str)
+            case_study.dPower_ImportExport.set_index(['hub', 'i', 'rp', 'k'], inplace=True)
 
         return None if inplace else case_study
 

From 0c33cba8e27cee00e48d6c84d29e20311d6d2854 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Thu, 21 May 2026 14:56:17 +0200
Subject: [PATCH 21/29] Extract function to calculate Weights_RP from Hindex

---
 CaseStudy.py | 41 +++++++++++++++++++----------------------
 1 file changed, 19 insertions(+), 22 deletions(-)

diff --git a/CaseStudy.py b/CaseStudy.py
index fbdcc5a..96ab41c 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -203,22 +203,7 @@ def __init__(self,
             self.dPower_WeightsRP = dPower_WeightsRP
         else:
             self.power_weightsrp_file = power_weightsrp_file
-            # Calculate dPower_WeightsRP from Hindex
-            dPower_WeightsRPs = []
-            for scenario in self.dPower_Hindex['scenario'].unique().tolist():
-                # Count occurences of each value in column 'rp' of dPower_Hindex
-                dPower_WeightsRP_scenario = pd.DataFrame(self.dPower_Hindex[self.dPower_Hindex['scenario'] == scenario].reset_index()['rp'].value_counts().sort_index())
-                dPower_WeightsRP_scenario = dPower_WeightsRP_scenario.rename(columns={'count': 'pWeight_rp'})
-                dPower_WeightsRP_scenario['scenario'] = scenario  # Add scenario ID
-
-                # Add other columns with default values
-                dPower_WeightsRP_scenario['id'] = np.nan
-                dPower_WeightsRP_scenario['dataPackage'] = np.nan
-                dPower_WeightsRP_scenario['dataSource'] = np.nan
-
-                dPower_WeightsRPs.append(dPower_WeightsRP_scenario)
-
-            dPower_WeightsRP = pd.concat(dPower_WeightsRPs, ignore_index=False)
+            dPower_WeightsRP = self.calculatePowerWeightsRP(db_id=np.nan, dataPackage=np.nan, dataSource=np.nan)
 
             if os.path.exists(self.data_folder + self.power_weightsrp_file):  # Compare with given file if it exists
                 self.dPower_WeightsRP = ExcelReader.get_Power_WeightsRP(self.data_folder + self.power_weightsrp_file)
@@ -226,19 +211,16 @@ def __init__(self,
                 calculated = dPower_WeightsRP.reset_index().set_index(["rp", "scenario"])
                 fromFile = self.dPower_WeightsRP.reset_index().set_index(["rp", "scenario"])
 
-                # Normalize both to sum to 1 for comparison
-                calc_norm = calculated['pWeight_rp'] / calculated['pWeight_rp'].sum()
-                file_norm = fromFile['pWeight_rp'] / fromFile['pWeight_rp'].sum()
                 # Align indices and fill missing with 0 for comparison
-                combined = pd.concat([calc_norm, file_norm], axis=1, keys=['calculated', 'fromFile']).fillna(0)
+                combined = pd.concat([calculated, fromFile], axis=1, keys=['calculated', 'fromFile']).fillna(0)
                 diff_mask = ~np.isclose(combined['calculated'], combined['fromFile'])
                 if diff_mask.any():
-                    printer.warning(f"Values for 'pWeight_rp' in `{self.data_folder + self.power_weightsrp_file}` do not match the calculated values based on `{self.power_hindex_file}`. Please check if this is intended, using the file `{self.data_folder + self.power_weightsrp_file}` instead of the calculated values.")
+                    printer.warning(f"Values for 'pWeight_rp' in `{self.data_folder + self.power_weightsrp_file}` do not match the calculated values based on `{self.power_hindex_file}`. Please check if this is intended, now using the file `{self.data_folder + self.power_weightsrp_file}` instead of the calculated values.")
                     # Print all differing lines
                     diffs = combined[diff_mask]
                     printer.warning("Differing entries (index -> calculated | fromFile):\n" + diffs.to_string())
             else:  # Use calculated dPower_WeightsRP otherwise
-                printer.warning(f"Executing without 'Power_WeightsRP' (since no file was found at '{self.data_folder + self.power_weightsrp_file}').")
+                printer.warning(f"Executing without 'Power_WeightsRP' (calculating from Power_Hindex, since no file was found at '{self.data_folder + self.power_weightsrp_file}').")
                 self.dPower_WeightsRP = dPower_WeightsRP
 
         if not self.do_not_filter_unused_scenarios:
@@ -799,6 +781,21 @@ def merge_generators(self, inplace: bool = False) -> Optional['CaseStudy']:
 
         return None if inplace else cs
 
+    def calculatePowerWeightsRP(self, db_id, dataPackage, dataSource):
+        dPower_WeightsRPs = []
+        for scenario in self.dPower_Hindex['scenario'].unique().tolist():
+            # Count occurrences of each value in column 'rp' of dPower_Hindex
+            dPower_WeightsRP_scenario = pd.DataFrame(self.dPower_Hindex[self.dPower_Hindex['scenario'] == scenario].reset_index()['rp'].value_counts().sort_index())
+            dPower_WeightsRP_scenario = dPower_WeightsRP_scenario.rename(columns={'count': 'pWeight_rp'})
+            dPower_WeightsRP_scenario['scenario'] = scenario  # Add scenario ID
+            dPower_WeightsRPs.append(dPower_WeightsRP_scenario)
+
+        dPower_WeightsRP = pd.concat(dPower_WeightsRPs, ignore_index=False)
+        dPower_WeightsRP['id'] = db_id
+        dPower_WeightsRP['dataPackage'] = dataPackage
+        dPower_WeightsRP['dataSource'] = dataSource
+        return dPower_WeightsRP
+
     # Create transition matrix from Hindex
     def get_rpTransitionMatrices(self, clip_method: str = "none", clip_value: float = 0) -> tuple[pd.DataFrame, pd.DataFrame, pd.DataFrame]:
         rps = sorted(self.dPower_Hindex.index.get_level_values('rp').unique().tolist())

From f064f6006d9a72b48c9d58a1b1322676d07e8100 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Thu, 21 May 2026 17:47:23 +0200
Subject: [PATCH 22/29] Implement shift_transition_matrix

Fix bug in calculating PowerWeights_RP
---
 CLAUDE.md    |  1 +
 CaseStudy.py | 68 +++++++++++++++++++++++++++++++++++++++++++++++++++-
 2 files changed, 68 insertions(+), 1 deletion(-)

diff --git a/CLAUDE.md b/CLAUDE.md
index 880fe18..58d5f1a 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -13,6 +13,7 @@ See `README.md` for usage, key concepts, and data structure.
 - `merge_single_node_buses()` preserves the `z` (zone) column as a sorted unique union string of all merged zones (e.g. `"R1_R2"`). This is documented in root `CLAUDE.md` as well.
 - `merge_generators()` collapses all generators sharing the same `(tec, i)` into one representative generator with ID `"{i}_{tec}"`. It must be called after construction (scaling is not required). VRESProfiles and Inflows are merged **before** VRES so that the original per-generator `MaxProd` weights are available for the capacity-factor weighted average. Generators in VRESProfiles/Inflows that have no matching entry in `dPower_VRES` (left-join miss) are grouped under `(tec=NaN, i=NaN)` — filter to a single scenario first to avoid mixing scenarios in the groupby.
 - `filter_zone(zone)` keeps only buses whose `z` column exactly matches the given zone name(s). Cascade order: BusInfo → Network (both endpoints must survive) → ThermalGen / VRES / Storage (by `i`) → Demand (by `i` index level) → VRESProfiles (by surviving VRES gen IDs) → Inflows (by union of surviving VRES + Storage gen IDs) → ImportExport (by `i` index level). `z` is an exact-match filter; merged-bus zone strings like `"R1_R2"` are not matched by `"R1"` alone.
+- `shift_transition_matrix(positions, inplace=True, seed=42)` cyclically shifts every row of `rpTransitionMatrixAbsolute` right by `positions` (via `np.roll`, negative shifts left), normalizes, then rebuilds Hindex by Markov chain sampling from the shifted row distributions — anchoring the first period of each scenario to its original RP assignment. Recomputes `dPower_WeightsRP` and all three TM attributes from the new Hindex. The RNG seed defaults to 42.
 - `CaseStudy.copy()` is a full `deepcopy` — safe to modify independently.
 - `to_full_hourly_model()` processes each enabled scenario from `dGlobal_Scenarios` independently. It is a no-op (returns unchanged) if all scenarios are already hourly (no two p-values share the same 'k' within any scenario). The produced Hindex uses consecutive h0001…hN p-labels and k0001…kN k-labels per scenario, all mapped to rp01 with weight 1.
 - Transition matrices (`rpTransitionMatrixAbsolute`, `rpTransitionMatrixRelativeTo`, `rpTransitionMatrixRelativeFrom`) are computed in the constructor and attached as attributes.
diff --git a/CaseStudy.py b/CaseStudy.py
index 96ab41c..18a1beb 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -212,7 +212,7 @@ def __init__(self,
                 fromFile = self.dPower_WeightsRP.reset_index().set_index(["rp", "scenario"])
 
                 # Align indices and fill missing with 0 for comparison
-                combined = pd.concat([calculated, fromFile], axis=1, keys=['calculated', 'fromFile']).fillna(0)
+                combined = pd.concat([calculated["pWeight_rp"], fromFile["pWeight_rp"]], axis=1, keys=['calculated', 'fromFile']).fillna(0)
                 diff_mask = ~np.isclose(combined['calculated'], combined['fromFile'])
                 if diff_mask.any():
                     printer.warning(f"Values for 'pWeight_rp' in `{self.data_folder + self.power_weightsrp_file}` do not match the calculated values based on `{self.power_hindex_file}`. Please check if this is intended, now using the file `{self.data_folder + self.power_weightsrp_file}` instead of the calculated values.")
@@ -794,6 +794,10 @@ def calculatePowerWeightsRP(self, db_id, dataPackage, dataSource):
         dPower_WeightsRP['id'] = db_id
         dPower_WeightsRP['dataPackage'] = dataPackage
         dPower_WeightsRP['dataSource'] = dataSource
+
+        scenario_sums = self.dPower_WeightsK.groupby('scenario')['pWeight_k'].sum()
+        dPower_WeightsRP['pWeight_rp'] = dPower_WeightsRP['pWeight_rp'] / dPower_WeightsRP['scenario'].map(scenario_sums)
+
         return dPower_WeightsRP
 
     # Create transition matrix from Hindex
@@ -1158,6 +1162,68 @@ def shift_ks(self, shift: int, inplace: bool = False) -> Optional[Self]:
 
         return None if inplace else case_study
 
+    def shift_transition_matrix(self, positions: int, inplace: bool = True, seed: int = 42) -> Optional['CaseStudy']:
+        """
+        Adjust the transition matrix by shifting the probabilities by <positions> positions, resampling Hindex from the adjusted
+        target distribution. dPower_WeightsRP and all three transition-matrices are recomputed from the new Hindex.
+
+        :param positions: Number of positions to shift to the right (negative shifts to the left).
+        :param inplace: If True, modifies the current instance. If False, returns a new instance.
+        :param seed: Random seed to guarantee reproduceability.
+        :return: None if inplace is True, otherwise a new CaseStudy instance.
+        """
+        cs = self if inplace else self.copy()
+
+        rps = cs.rpTransitionMatrixAbsolute.index.tolist()
+        n = len(rps)
+
+        target_probs: dict[str, np.ndarray] = {}
+
+        printer.information(f"Adjusting transition matrix, shifting it by {positions} positions")
+        for rp in rps:
+            c = cs.rpTransitionMatrixAbsolute.loc[rp].values.astype(float)
+            total = c.sum()
+            if total == 0:
+                target_probs[rp] = np.ones(n) / n
+            else:
+                c_shifted = np.roll(c, positions)
+                target_probs[rp] = c_shifted / total
+
+        # Rebuild Hindex: sample a new RP sequence per scenario
+        hindex_flat = cs.dPower_Hindex.reset_index()
+        new_parts = []
+        rng = np.random.default_rng(seed)
+
+        for scenario in hindex_flat['scenario'].unique().tolist():
+            sc = hindex_flat[hindex_flat['scenario'] == scenario].copy()
+            sc = sc.sort_values(['p'])
+            n_ks_per_rp = len(cs.dPower_WeightsK['scenario'] == scenario)
+
+            period_labels = sc['p'].tolist()
+            first_rp = sc['rp'].iloc[0]
+            n_periods = len(period_labels)
+
+            # Sample new RP sequence, anchoring the first period to the original
+            new_rp_seq = [first_rp for _ in range(n_ks_per_rp)]
+            for _ in range(n_periods - 1):
+                probs = target_probs[new_rp_seq[-1]]
+                next_rp = str(rng.choice(rps, p=probs))
+                new_rp_seq.extend([next_rp for _ in range(n_ks_per_rp)])
+
+            sc['rp'] = sc['p'].map(dict(zip(period_labels, new_rp_seq)))
+            new_parts.append(sc)
+
+        new_hindex = pd.concat(new_parts, ignore_index=False)
+        cs.dPower_Hindex = new_hindex.set_index(['p', 'rp', 'k'])
+
+        # Recompute WeightsRP from new Hindex
+        cs.dPower_WeightsRP = cs.calculatePowerWeightsRP(np.nan, np.nan, np.nan)
+
+        # Recompute actual TM from new Hindex (approximates the target distributions)
+        cs.rpTransitionMatrixAbsolute, cs.rpTransitionMatrixRelativeTo, cs.rpTransitionMatrixRelativeFrom = cs.get_rpTransitionMatrices()
+
+        return None if inplace else cs
+
     def apply_kmedoids_aggregation(self, number_rps: int, rp_length: int = 24,
                                    cluster_strategy: Literal["aggregated", "disaggregated"] = "aggregated",
                                    capacity_normalization: Literal["installed", "maxInvestment"] = "maxInvestment",

From 27e3207c606d731208391490bd0a65f14d7dacf4 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Thu, 21 May 2026 18:42:13 +0200
Subject: [PATCH 23/29] Implement plot_transition_matrix

---
 Utilities.py | 97 ++++++++++++++++++++++++++++++++++++++++++++++++++++
 1 file changed, 97 insertions(+)

diff --git a/Utilities.py b/Utilities.py
index 444e114..fab01a1 100644
--- a/Utilities.py
+++ b/Utilities.py
@@ -3,6 +3,7 @@
 import typing
 from typing import TYPE_CHECKING, Literal, Dict
 
+import matplotlib.pyplot as plt
 import numpy as np
 import pandas as pd
 import tsam.timeseriesaggregation as tsam
@@ -619,3 +620,99 @@ def apply_kmedoids_aggregation(
         inplace=inplace,
         verbose=verbose
     )
+
+
+def plot_transition_matrix(tm: pd.DataFrame, title: str | None = None, output: str | None = None):
+    """Plot a transition matrix as a Blue-tinted coloured table.
+
+    Cell text shows absolute counts and row-normalised percentages; cell colour
+    intensity encodes the row-normalised probability (darker = more likely).
+    Row and column sums are appended in bold.
+
+    :param tm: Transition matrix DataFrame (index and columns are RP labels).
+               Typically ``cs.rpTransitionMatrixAbsolute``.
+    :param title: Optional subtitle shown below the main heading.
+    :param output: If given, save the figure to this path instead of displaying it.
+    """
+    labels = list(tm.index)
+    n = len(labels)
+    data = tm.values.astype(float)
+
+    row_totals = data.sum(axis=1)  # shape (n,)
+    col_totals = data.sum(axis=0)  # shape (n,)
+    grand_total = data.sum()
+
+    rel = data / np.where(row_totals == 0, 1, row_totals)[:, np.newaxis]
+
+    fig_w = max(4.0, 0.9 * (n + 1) + 1.5)
+    fig_h = max(2.5, 0.7 * (n + 1) + 1.5)
+    fig, ax = plt.subplots(figsize=(fig_w, fig_h))
+    ax.set_axis_off()
+
+    cmap = plt.cm.Blues
+
+    cell_text = [
+        [f"{int(data[r, c])}\n{rel[r, c]:.0%}" for c in range(n)] + [f"{int(row_totals[r])}"]
+        for r in range(n)
+    ]
+    cell_text.append([f"{int(col_totals[c])}" for c in range(n)] + [f"{int(grand_total)}"])
+
+    cell_colors = [
+        [cmap(rel[r, c]) for c in range(n)] + ["white"]
+        for r in range(n)
+    ]
+    cell_colors.append(["white"] * (n + 1))
+
+    sum_labels = labels + ["Sum"]
+
+    tbl = ax.table(
+        cellText=cell_text,
+        rowLabels=sum_labels,
+        colLabels=sum_labels,
+        cellColours=cell_colors,
+        loc="center",
+        cellLoc="center",
+    )
+    tbl.auto_set_font_size(False)
+    tbl.set_fontsize(8)
+    tbl.scale(1, 1.6)
+
+    for r in range(n + 2):
+        tbl[r, n].get_text().set_fontweight("bold")
+    for c in range(-1, n + 1):
+        tbl[n + 1, c].get_text().set_fontweight("bold")
+
+    heading = "Transition matrix (count / row %)"
+    if title:
+        heading += f"\n{title}"
+    ax.set_title(heading, fontsize=10, pad=8)
+    fig.tight_layout()
+
+    # Force layout so cell bounding boxes are finalised
+    fig.canvas.draw()
+    renderer = fig.canvas.get_renderer()
+
+    def _to_axes(xd, yd):
+        return ax.transAxes.inverted().transform((xd, yd))
+
+    bb0 = tbl[n + 1, 0].get_window_extent(renderer)
+    bb_rc = tbl[n + 1, n].get_window_extent(renderer)
+    bb_rl = tbl[n + 1, -1].get_window_extent(renderer)
+    x0_ax, y_ax = _to_axes(bb_rl.x0, bb0.y1)
+    x1_ax, _ = _to_axes(bb_rc.x1, bb0.y1)
+    ax.plot([x0_ax, x1_ax], [y_ax, y_ax], transform=ax.transAxes,
+            color="black", linewidth=2, clip_on=False, zorder=10)
+
+    bb_hdr = tbl[0, n].get_window_extent(renderer)
+    bb_bot = tbl[n + 1, n].get_window_extent(renderer)
+    x_ax, y0_ax = _to_axes(bb_hdr.x0, bb_bot.y0)
+    _, y1_ax = _to_axes(bb_hdr.x0, bb_hdr.y1)
+    ax.plot([x_ax, x_ax], [y0_ax, y1_ax], transform=ax.transAxes,
+            color="black", linewidth=2, clip_on=False, zorder=10)
+
+    if output:
+        fig.savefig(output, dpi=150, bbox_inches="tight")
+        print(f"Saved to {output}")
+    else:
+        plt.show()
+    plt.close(fig)

From d908cc9f483b4e9a9a63113bfb49506a9607e2b9 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Thu, 21 May 2026 19:50:29 +0200
Subject: [PATCH 24/29] Fix possibility of RPs not occuring at all after
 shifting TM

---
 CaseStudy.py | 16 ++++++++++++++++
 1 file changed, 16 insertions(+)

diff --git a/CaseStudy.py b/CaseStudy.py
index 18a1beb..0838f15 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -1211,6 +1211,22 @@ def shift_transition_matrix(self, positions: int, inplace: bool = True, seed: in
                 new_rp_seq.extend([next_rp for _ in range(n_ks_per_rp)])
 
             sc['rp'] = sc['p'].map(dict(zip(period_labels, new_rp_seq)))
+
+            # Ensure every RP appears at least once; if not, replace the most frequent period with the missing one
+            rp_count = sc.groupby('rp')['p'].count() / n_ks_per_rp
+            missing_rps = [rp for rp in rps if rp_count.get(rp, 0) == 0]
+            while missing_rps:
+                missing_rp = missing_rps.pop(0)
+                max_rp = rp_count.idxmax()
+                if rp_count[max_rp] <= 1:
+                    raise ValueError(f"It seems like there are more RPs than periods in the case study. Check your data and settings.")
+
+                first_occurence_of_max_rp = sc.rp.eq(max_rp).idxmax()
+                sc.loc[sc.index[first_occurence_of_max_rp:first_occurence_of_max_rp + n_ks_per_rp], 'rp'] = missing_rp
+
+                rp_count.loc[max_rp] -= 1
+                rp_count.loc[missing_rp] = 1
+
             new_parts.append(sc)
 
         new_hindex = pd.concat(new_parts, ignore_index=False)

From 2b82a25b228e0b65b3f91652102224f4697e215b Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Fri, 22 May 2026 13:04:41 +0200
Subject: [PATCH 25/29] Implmeent perturb_transition_matrix in CaseStudy

---
 CLAUDE.md    |   1 +
 CaseStudy.py | 111 ++++++++++++++++++++++++++++++++++++---------------
 2 files changed, 80 insertions(+), 32 deletions(-)

diff --git a/CLAUDE.md b/CLAUDE.md
index 58d5f1a..c6a2005 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -14,6 +14,7 @@ See `README.md` for usage, key concepts, and data structure.
 - `merge_generators()` collapses all generators sharing the same `(tec, i)` into one representative generator with ID `"{i}_{tec}"`. It must be called after construction (scaling is not required). VRESProfiles and Inflows are merged **before** VRES so that the original per-generator `MaxProd` weights are available for the capacity-factor weighted average. Generators in VRESProfiles/Inflows that have no matching entry in `dPower_VRES` (left-join miss) are grouped under `(tec=NaN, i=NaN)` — filter to a single scenario first to avoid mixing scenarios in the groupby.
 - `filter_zone(zone)` keeps only buses whose `z` column exactly matches the given zone name(s). Cascade order: BusInfo → Network (both endpoints must survive) → ThermalGen / VRES / Storage (by `i`) → Demand (by `i` index level) → VRESProfiles (by surviving VRES gen IDs) → Inflows (by union of surviving VRES + Storage gen IDs) → ImportExport (by `i` index level). `z` is an exact-match filter; merged-bus zone strings like `"R1_R2"` are not matched by `"R1"` alone.
 - `shift_transition_matrix(positions, inplace=True, seed=42)` cyclically shifts every row of `rpTransitionMatrixAbsolute` right by `positions` (via `np.roll`, negative shifts left), normalizes, then rebuilds Hindex by Markov chain sampling from the shifted row distributions — anchoring the first period of each scenario to its original RP assignment. Recomputes `dPower_WeightsRP` and all three TM attributes from the new Hindex. The RNG seed defaults to 42.
+- `perturb_transition_matrix(randomness, inplace=True, seed=42)` interpolates each row between its original distribution and a random draw: `new_prob = (1-randomness)*orig_prob + randomness*random_draw`, where `random_draw` is `n` uniform values normalized to sum to 1. `randomness=0.0` leaves the matrix unchanged; `randomness=1.0` replaces it fully. Rebuilds Hindex identically to `shift_transition_matrix` (Markov chain sampling, first period anchored, missing-RP correction loop). Uses a single RNG instance for both the random draws and the chain sampling.
 - `CaseStudy.copy()` is a full `deepcopy` — safe to modify independently.
 - `to_full_hourly_model()` processes each enabled scenario from `dGlobal_Scenarios` independently. It is a no-op (returns unchanged) if all scenarios are already hourly (no two p-values share the same 'k' within any scenario). The produced Hindex uses consecutive h0001…hN p-labels and k0001…kN k-labels per scenario, all mapped to rp01 with weight 1.
 - Transition matrices (`rpTransitionMatrixAbsolute`, `rpTransitionMatrixRelativeTo`, `rpTransitionMatrixRelativeFrom`) are computed in the constructor and attached as attributes.
diff --git a/CaseStudy.py b/CaseStudy.py
index 2cc017b..1ab153f 100644
--- a/CaseStudy.py
+++ b/CaseStudy.py
@@ -1203,42 +1203,14 @@ def shift_ks(self, shift: int, inplace: bool = False) -> Optional[Self]:
 
         return None if inplace else case_study
 
-    def shift_transition_matrix(self, positions: int, inplace: bool = True, seed: int = 42) -> Optional['CaseStudy']:
-        """
-        Adjust the transition matrix by shifting the probabilities by <positions> positions, resampling Hindex from the adjusted
-        target distribution. dPower_WeightsRP and all three transition-matrices are recomputed from the new Hindex.
-
-        :param positions: Number of positions to shift to the right (negative shifts to the left).
-        :param inplace: If True, modifies the current instance. If False, returns a new instance.
-        :param seed: Random seed to guarantee reproduceability.
-        :return: None if inplace is True, otherwise a new CaseStudy instance.
-        """
-        cs = self if inplace else self.copy()
-
-        rps = cs.rpTransitionMatrixAbsolute.index.tolist()
-        n = len(rps)
-
-        target_probs: dict[str, np.ndarray] = {}
-
-        printer.information(f"Adjusting transition matrix, shifting it by {positions} positions")
-        for rp in rps:
-            c = cs.rpTransitionMatrixAbsolute.loc[rp].values.astype(float)
-            total = c.sum()
-            if total == 0:
-                target_probs[rp] = np.ones(n) / n
-            else:
-                c_shifted = np.roll(c, positions)
-                target_probs[rp] = c_shifted / total
-
-        # Rebuild Hindex: sample a new RP sequence per scenario
-        hindex_flat = cs.dPower_Hindex.reset_index()
+    def _resample_hindex_from_target_probs(self, target_probs: dict[str, np.ndarray], rps: list[str], rng: np.random.Generator) -> None:
+        hindex_flat = self.dPower_Hindex.reset_index()
         new_parts = []
-        rng = np.random.default_rng(seed)
 
         for scenario in hindex_flat['scenario'].unique().tolist():
             sc = hindex_flat[hindex_flat['scenario'] == scenario].copy()
             sc = sc.sort_values(['p'])
-            n_ks_per_rp = len(cs.dPower_WeightsK['scenario'] == scenario)
+            n_ks_per_rp = len(self.dPower_WeightsK['scenario'] == scenario)
 
             period_labels = sc['p'].tolist()
             first_rp = sc['rp'].iloc[0]
@@ -1271,7 +1243,82 @@ def shift_transition_matrix(self, positions: int, inplace: bool = True, seed: in
             new_parts.append(sc)
 
         new_hindex = pd.concat(new_parts, ignore_index=False)
-        cs.dPower_Hindex = new_hindex.set_index(['p', 'rp', 'k'])
+        self.dPower_Hindex = new_hindex.set_index(['p', 'rp', 'k'])
+
+    def shift_transition_matrix(self, positions: int, inplace: bool = True, seed: int = 42) -> Optional['CaseStudy']:
+        """
+        Adjust the transition matrix by shifting the probabilities by <positions> positions, resampling Hindex from the adjusted
+        target distribution. dPower_WeightsRP and all three transition-matrices are recomputed from the new Hindex.
+
+        :param positions: Number of positions to shift to the right (negative shifts to the left).
+        :param inplace: If True, modifies the current instance. If False, returns a new instance.
+        :param seed: Random seed to guarantee reproduceability.
+        :return: None if inplace is True, otherwise a new CaseStudy instance.
+        """
+        cs = self if inplace else self.copy()
+
+        rps = cs.rpTransitionMatrixAbsolute.index.tolist()
+        n = len(rps)
+
+        target_probs: dict[str, np.ndarray] = {}
+
+        printer.information(f"Adjusting transition matrix, shifting it by {positions} positions")
+        for rp in rps:
+            c = cs.rpTransitionMatrixAbsolute.loc[rp].values.astype(float)
+            total = c.sum()
+            if total == 0:
+                target_probs[rp] = np.ones(n) / n
+            else:
+                c_shifted = np.roll(c, positions)
+                target_probs[rp] = c_shifted / total
+
+        rng = np.random.default_rng(seed)
+        cs._resample_hindex_from_target_probs(target_probs, rps, rng)
+
+        # Recompute WeightsRP from new Hindex
+        cs.dPower_WeightsRP = cs.calculatePowerWeightsRP(np.nan, np.nan, np.nan)
+
+        # Recompute actual TM from new Hindex (approximates the target distributions)
+        cs.rpTransitionMatrixAbsolute, cs.rpTransitionMatrixRelativeTo, cs.rpTransitionMatrixRelativeFrom = cs.get_rpTransitionMatrices()
+
+        return None if inplace else cs
+
+    def perturb_transition_matrix(self, randomness: float, inplace: bool = True, seed: int = 42) -> Optional['CaseStudy']:
+        """
+        Adjust the transition matrix by interpolating each row between its original distribution and a random draw:
+        new_prob = (1 - randomness) * orig_prob + randomness * random_draw
+        Resamples Hindex from the adjusted target distribution. dPower_WeightsRP and all three transition-matrices
+        are recomputed from the new Hindex.
+
+        :param randomness: Interpolation factor [0.0, 1.0]. 0.0 leaves the matrix unchanged; 1.0 replaces it fully with random draws.
+        :param inplace: If True, modifies the current instance. If False, returns a new instance.
+        :param seed: Random seed to guarantee reproducibility.
+        :return: None if inplace is True, otherwise a new CaseStudy instance.
+        """
+        if not 0.0 <= randomness <= 1.0:
+            raise ValueError(f"randomness must be in [0.0, 1.0], got {randomness}")
+
+        cs = self if inplace else self.copy()
+
+        rps = cs.rpTransitionMatrixAbsolute.index.tolist()
+        n = len(rps)
+
+        target_probs: dict[str, np.ndarray] = {}
+        rng = np.random.default_rng(seed)
+
+        printer.information(f"Adjusting transition matrix, perturbing it with randomness={randomness}")
+        for rp in rps:
+            c = cs.rpTransitionMatrixAbsolute.loc[rp].values.astype(float)
+            total = c.sum()
+            if total == 0:
+                target_probs[rp] = np.ones(n) / n
+            else:
+                orig_prob = c / total
+                raw = rng.random(n)
+                random_draw = raw / raw.sum()
+                target_probs[rp] = (1 - randomness) * orig_prob + randomness * random_draw
+
+        cs._resample_hindex_from_target_probs(target_probs, rps, rng)
 
         # Recompute WeightsRP from new Hindex
         cs.dPower_WeightsRP = cs.calculatePowerWeightsRP(np.nan, np.nan, np.nan)

From 467ab074bed73f177b18bdd17ea49798524155a5 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 26 May 2026 17:38:14 +0200
Subject: [PATCH 26/29] Fix import for Caller.py so it can be called from the
 root

---
 Caller.py | 5 ++++-
 1 file changed, 4 insertions(+), 1 deletion(-)

diff --git a/Caller.py b/Caller.py
index 30b826d..248b389 100644
--- a/Caller.py
+++ b/Caller.py
@@ -5,7 +5,10 @@
 import sys
 import time
 
-from printer import Printer
+if __name__ == "__main__":
+    sys.path.insert(0, os.path.dirname(os.path.dirname(os.path.abspath(__file__))))
+
+from InOutModule.printer import Printer
 
 printer = Printer.getInstance()
 

From c0ecd5d99a8d6184ac5f813b90994ab2a0b6b7a2 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Mon, 1 Jun 2026 16:09:21 +0200
Subject: [PATCH 27/29] Fix setting Caller-title if a lot of jobs are processed

---
 Caller.py | 3 ++-
 1 file changed, 2 insertions(+), 1 deletion(-)

diff --git a/Caller.py b/Caller.py
index 248b389..0124040 100644
--- a/Caller.py
+++ b/Caller.py
@@ -1,4 +1,5 @@
 import argparse
+import ctypes
 import datetime
 import os
 import subprocess
@@ -116,7 +117,7 @@ def _any_previous_unclaimed(jobs_file, lines, barrier_index):
         log_file = f"{args.jobs}.log{i}"
         try:
             printer.information(f"Executing job {i} from '{args.jobs}': {line.strip()}")
-            os.system(f"title Job {i} from '{args.jobs}': {line.strip()}")
+            ctypes.windll.kernel32.SetConsoleTitleW(f"Job {i} from '{args.jobs}': {line.strip()}")  # type: ignore[attr-defined] - it is Windows-only anyway, so we can ignore the fact that this attribute doesn't exist on other platforms
 
             start_time = time.time()
             with open(log_file, 'w', encoding='utf-8') as log_f:

From ea282c58ec027f16b4beb5a377ddbddc77cec2d5 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Wed, 17 Jun 2026 23:45:12 +0200
Subject: [PATCH 28/29] Fix spacing in Transition Matrix plots

---
 CLAUDE.md    | 1 +
 Utilities.py | 9 +++++----
 2 files changed, 6 insertions(+), 4 deletions(-)

diff --git a/CLAUDE.md b/CLAUDE.md
index c6a2005..f4fca6a 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -40,6 +40,7 @@ See `README.md` for usage, key concepts, and data structure.
 
 - `inflowsToCapacityFactors()` joins inflows onto `vresProfiles_df` by dividing by `MaxProd`; generators with missing or zero `MaxProd` are dropped with a warning.
 - `capacityFactorsToInflows()` is the inverse; the `remove_Inflows_from_VRESProfiles_inplace` flag modifies the input DataFrame in place when set.
+- `plot_transition_matrix()` passes `bbox=[0, 0, 1, 1]` to `ax.table()` so the table is forced to fill the entire Axes — without it, matplotlib sizes table rows from font metrics rather than the Axes height, leaving large dead space above/below when the figure is taller than the table's natural size. `fig_w`/`fig_h` are tuned per-cell-inch (not just a generic min-size heuristic) since they now directly determine the rendered cell aspect ratio.
 
 ### Printer
 
diff --git a/Utilities.py b/Utilities.py
index fab01a1..859011e 100644
--- a/Utilities.py
+++ b/Utilities.py
@@ -644,8 +644,9 @@ def plot_transition_matrix(tm: pd.DataFrame, title: str | None = None, output: s
 
     rel = data / np.where(row_totals == 0, 1, row_totals)[:, np.newaxis]
 
-    fig_w = max(4.0, 0.9 * (n + 1) + 1.5)
-    fig_h = max(2.5, 0.7 * (n + 1) + 1.5)
+    title_h = 0.55 if title else 0.35
+    fig_w = max(4.0, 0.95 * (n + 1) + 0.3)
+    fig_h = 0.4 * (n + 2) + title_h
     fig, ax = plt.subplots(figsize=(fig_w, fig_h))
     ax.set_axis_off()
 
@@ -672,10 +673,10 @@ def plot_transition_matrix(tm: pd.DataFrame, title: str | None = None, output: s
         cellColours=cell_colors,
         loc="center",
         cellLoc="center",
+        bbox=[0, 0, 1, 1],
     )
     tbl.auto_set_font_size(False)
     tbl.set_fontsize(8)
-    tbl.scale(1, 1.6)
 
     for r in range(n + 2):
         tbl[r, n].get_text().set_fontweight("bold")
@@ -686,7 +687,7 @@ def plot_transition_matrix(tm: pd.DataFrame, title: str | None = None, output: s
     if title:
         heading += f"\n{title}"
     ax.set_title(heading, fontsize=10, pad=8)
-    fig.tight_layout()
+    fig.tight_layout(pad=0.3)
 
     # Force layout so cell bounding boxes are finalised
     fig.canvas.draw()

From 8b1f53a75d152645e5ff8c7d5f417befaf316da5 Mon Sep 17 00:00:00 2001
From: "Felix C. A. Auer" <10127354+FelixCAAuer@users.noreply.github.com>
Date: Tue, 7 Jul 2026 22:37:32 +0200
Subject: [PATCH 29/29] Create plots for both the column- and row-normalized
 Transition Matrices

Add white text for sections which are otherwise too dark
---
 CLAUDE.md    |   2 +-
 Utilities.py | 186 +++++++++++++++++++++++++++++----------------------
 2 files changed, 107 insertions(+), 81 deletions(-)

diff --git a/CLAUDE.md b/CLAUDE.md
index f4fca6a..e8974d3 100644
--- a/CLAUDE.md
+++ b/CLAUDE.md
@@ -40,7 +40,7 @@ See `README.md` for usage, key concepts, and data structure.
 
 - `inflowsToCapacityFactors()` joins inflows onto `vresProfiles_df` by dividing by `MaxProd`; generators with missing or zero `MaxProd` are dropped with a warning.
 - `capacityFactorsToInflows()` is the inverse; the `remove_Inflows_from_VRESProfiles_inplace` flag modifies the input DataFrame in place when set.
-- `plot_transition_matrix()` passes `bbox=[0, 0, 1, 1]` to `ax.table()` so the table is forced to fill the entire Axes — without it, matplotlib sizes table rows from font metrics rather than the Axes height, leaving large dead space above/below when the figure is taller than the table's natural size. `fig_w`/`fig_h` are tuned per-cell-inch (not just a generic min-size heuristic) since they now directly determine the rendered cell aspect ratio.
+- `plot_transition_matrix()` emits **two** figures per call — a row-normalised ("from" perspective) and a column-normalised ("to" perspective) version — via the inner `_render()` helper. When `output` is given, `os.path.splitext` inserts a `-rowNorm` / `-colNorm` suffix before the extension (e.g. `MK-foo.png` → `MK-foo-rowNorm.png`, `MK-foo-colNorm.png`). Colour intensity and the per-cell percentage encode the respective normalisation; absolute counts and the row/col/grand totals are identical between the two. Cell text flips from black to white where the Blue fill is too dark to read on (perceived luminance `0.299·R + 0.587·G + 0.114·B < 0.5`). It passes `bbox=[0, 0, 1, 1]` to `ax.table()` so the table is forced to fill the entire Axes — without it, matplotlib sizes table rows from font metrics rather than the Axes height, leaving large dead space above/below when the figure is taller than the table's natural size. `fig_w`/`fig_h` are tuned per-cell-inch (not just a generic min-size heuristic) since they now directly determine the rendered cell aspect ratio.
 
 ### Printer
 
diff --git a/Utilities.py b/Utilities.py
index 859011e..712835a 100644
--- a/Utilities.py
+++ b/Utilities.py
@@ -1,5 +1,6 @@
 from __future__ import annotations
 
+import os
 import typing
 from typing import TYPE_CHECKING, Literal, Dict
 
@@ -623,16 +624,20 @@ def apply_kmedoids_aggregation(
 
 
 def plot_transition_matrix(tm: pd.DataFrame, title: str | None = None, output: str | None = None):
-    """Plot a transition matrix as a Blue-tinted coloured table.
+    """Plot a transition matrix as row- and column-normalised Blue-tinted tables.
 
-    Cell text shows absolute counts and row-normalised percentages; cell colour
-    intensity encodes the row-normalised probability (darker = more likely).
-    Row and column sums are appended in bold.
+    Two figures are produced: one normalised by row sum ("from" perspective) and
+    one normalised by column sum ("to" perspective). In each, cell text shows
+    absolute counts and the respective normalised percentage; cell colour
+    intensity encodes that probability (darker = more likely). Row and column
+    sums are appended in bold.
 
     :param tm: Transition matrix DataFrame (index and columns are RP labels).
                Typically ``cs.rpTransitionMatrixAbsolute``.
     :param title: Optional subtitle shown below the main heading.
-    :param output: If given, save the figure to this path instead of displaying it.
+    :param output: If given, save the figures to this path with ``-rowNorm`` /
+                   ``-colNorm`` inserted before the extension, instead of
+                   displaying them.
     """
     labels = list(tm.index)
     n = len(labels)
@@ -642,78 +647,99 @@ def plot_transition_matrix(tm: pd.DataFrame, title: str | None = None, output: s
     col_totals = data.sum(axis=0)  # shape (n,)
     grand_total = data.sum()
 
-    rel = data / np.where(row_totals == 0, 1, row_totals)[:, np.newaxis]
-
-    title_h = 0.55 if title else 0.35
-    fig_w = max(4.0, 0.95 * (n + 1) + 0.3)
-    fig_h = 0.4 * (n + 2) + title_h
-    fig, ax = plt.subplots(figsize=(fig_w, fig_h))
-    ax.set_axis_off()
-
-    cmap = plt.cm.Blues
-
-    cell_text = [
-        [f"{int(data[r, c])}\n{rel[r, c]:.0%}" for c in range(n)] + [f"{int(row_totals[r])}"]
-        for r in range(n)
-    ]
-    cell_text.append([f"{int(col_totals[c])}" for c in range(n)] + [f"{int(grand_total)}"])
-
-    cell_colors = [
-        [cmap(rel[r, c]) for c in range(n)] + ["white"]
-        for r in range(n)
-    ]
-    cell_colors.append(["white"] * (n + 1))
-
-    sum_labels = labels + ["Sum"]
-
-    tbl = ax.table(
-        cellText=cell_text,
-        rowLabels=sum_labels,
-        colLabels=sum_labels,
-        cellColours=cell_colors,
-        loc="center",
-        cellLoc="center",
-        bbox=[0, 0, 1, 1],
-    )
-    tbl.auto_set_font_size(False)
-    tbl.set_fontsize(8)
-
-    for r in range(n + 2):
-        tbl[r, n].get_text().set_fontweight("bold")
-    for c in range(-1, n + 1):
-        tbl[n + 1, c].get_text().set_fontweight("bold")
-
-    heading = "Transition matrix (count / row %)"
-    if title:
-        heading += f"\n{title}"
-    ax.set_title(heading, fontsize=10, pad=8)
-    fig.tight_layout(pad=0.3)
-
-    # Force layout so cell bounding boxes are finalised
-    fig.canvas.draw()
-    renderer = fig.canvas.get_renderer()
-
-    def _to_axes(xd, yd):
-        return ax.transAxes.inverted().transform((xd, yd))
-
-    bb0 = tbl[n + 1, 0].get_window_extent(renderer)
-    bb_rc = tbl[n + 1, n].get_window_extent(renderer)
-    bb_rl = tbl[n + 1, -1].get_window_extent(renderer)
-    x0_ax, y_ax = _to_axes(bb_rl.x0, bb0.y1)
-    x1_ax, _ = _to_axes(bb_rc.x1, bb0.y1)
-    ax.plot([x0_ax, x1_ax], [y_ax, y_ax], transform=ax.transAxes,
-            color="black", linewidth=2, clip_on=False, zorder=10)
-
-    bb_hdr = tbl[0, n].get_window_extent(renderer)
-    bb_bot = tbl[n + 1, n].get_window_extent(renderer)
-    x_ax, y0_ax = _to_axes(bb_hdr.x0, bb_bot.y0)
-    _, y1_ax = _to_axes(bb_hdr.x0, bb_hdr.y1)
-    ax.plot([x_ax, x_ax], [y0_ax, y1_ax], transform=ax.transAxes,
-            color="black", linewidth=2, clip_on=False, zorder=10)
-
-    if output:
-        fig.savefig(output, dpi=150, bbox_inches="tight")
-        print(f"Saved to {output}")
-    else:
-        plt.show()
-    plt.close(fig)
+    row_rel = data / np.where(row_totals == 0, 1, row_totals)[:, np.newaxis]
+    col_rel = data / np.where(col_totals == 0, 1, col_totals)[np.newaxis, :]
+
+    def _render(rel: np.ndarray, pct_label: str, out_path: str | None):
+        title_h = 0.55 if title else 0.35
+        fig_w = max(4.0, 0.95 * (n + 1) + 0.3)
+        fig_h = 0.4 * (n + 2) + title_h
+        fig, ax = plt.subplots(figsize=(fig_w, fig_h))
+        ax.set_axis_off()
+
+        cmap = plt.cm.Blues
+
+        cell_text = [
+            [f"{int(data[r, c])}\n{rel[r, c]:.0%}" for c in range(n)] + [f"{int(row_totals[r])}"]
+            for r in range(n)
+        ]
+        cell_text.append([f"{int(col_totals[c])}" for c in range(n)] + [f"{int(grand_total)}"])
+
+        cell_colors = [
+            [cmap(rel[r, c]) for c in range(n)] + ["white"]
+            for r in range(n)
+        ]
+        cell_colors.append(["white"] * (n + 1))
+
+        sum_labels = labels + ["Sum"]
+
+        tbl = ax.table(
+            cellText=cell_text,
+            rowLabels=sum_labels,
+            colLabels=sum_labels,
+            cellColours=cell_colors,
+            loc="center",
+            cellLoc="center",
+            bbox=[0, 0, 1, 1],
+        )
+        tbl.auto_set_font_size(False)
+        tbl.set_fontsize(8)
+
+        for r in range(n + 2):
+            tbl[r, n].get_text().set_fontweight("bold")
+        for c in range(-1, n + 1):
+            tbl[n + 1, c].get_text().set_fontweight("bold")
+
+        # Flip to white text where the Blue shade is too dark to read black on.
+        # (data cells are table rows 1..n, cols 0..n-1; row 0 is the header)
+        for r in range(n):
+            for c in range(n):
+                red, green, blue, _ = cmap(rel[r, c])
+                luminance = 0.299 * red + 0.587 * green + 0.114 * blue
+                if luminance < 0.5:
+                    tbl[r + 1, c].get_text().set_color("white")
+
+        heading = f"Transition matrix (count / {pct_label})"
+        if title:
+            heading += f"\n{title}"
+        ax.set_title(heading, fontsize=10, pad=8)
+        fig.tight_layout(pad=0.3)
+
+        # Force layout so cell bounding boxes are finalised
+        fig.canvas.draw()
+        renderer = fig.canvas.get_renderer()
+
+        def _to_axes(xd, yd):
+            return ax.transAxes.inverted().transform((xd, yd))
+
+        bb0 = tbl[n + 1, 0].get_window_extent(renderer)
+        bb_rc = tbl[n + 1, n].get_window_extent(renderer)
+        bb_rl = tbl[n + 1, -1].get_window_extent(renderer)
+        x0_ax, y_ax = _to_axes(bb_rl.x0, bb0.y1)
+        x1_ax, _ = _to_axes(bb_rc.x1, bb0.y1)
+        ax.plot([x0_ax, x1_ax], [y_ax, y_ax], transform=ax.transAxes,
+                color="black", linewidth=2, clip_on=False, zorder=10)
+
+        bb_hdr = tbl[0, n].get_window_extent(renderer)
+        bb_bot = tbl[n + 1, n].get_window_extent(renderer)
+        x_ax, y0_ax = _to_axes(bb_hdr.x0, bb_bot.y0)
+        _, y1_ax = _to_axes(bb_hdr.x0, bb_hdr.y1)
+        ax.plot([x_ax, x_ax], [y0_ax, y1_ax], transform=ax.transAxes,
+                color="black", linewidth=2, clip_on=False, zorder=10)
+
+        if out_path:
+            fig.savefig(out_path, dpi=150, bbox_inches="tight")
+            print(f"Saved to {out_path}")
+        else:
+            plt.show()
+        plt.close(fig)
+
+    for rel_matrix, norm_suffix, norm_label in (
+        (row_rel, "rowNorm", "row %"),
+        (col_rel, "colNorm", "column %"),
+    ):
+        target = None
+        if output:
+            base, ext = os.path.splitext(output)
+            target = f"{base}-{norm_suffix}{ext}"
+        _render(rel_matrix, norm_label, target)